Is it possible to get the trace value of this matrix? I want to find a general formula for the trace of the following $N\times N$ matrix raised to the power of $d$, where $d \in \mathbb{N}$.
$$
\begin{bmatrix}
 N-2 & 1 & 0 & \cdots & 0 & 0 \\
 1 & N-3 & 1 & \cdots & 0 & 0 \\
 0 & 1 & N-3 & \cdots & 0 & 0 \\
 \vdots  & \vdots & \vdots & \ddots & \vdots & \vdots \\
 0 & 0 & 0 & \cdots & N-3 & 1 \\
 0 & 0 & 0 & \cdots & 1 & N-2     
 \end{bmatrix}^d
$$
For example when $N=4$ and $d=2$, it will look like this
$$
\begin{bmatrix}
 2 & 1 & 0 & 0 \\
 1 & 1 & 1 & 0 \\
 0 & 1 & 1 & 1 \\
 0 & 0 & 1 & 2 \\   
 \end{bmatrix}^2
=
\begin{bmatrix}
 5 & 3 & 1 & 0 \\
 3 & 3 & 2 & 1 \\
 1 & 2 & 3 & 3 \\
 0 & 1 & 3 & 5 \\   
 \end{bmatrix},
$$
so the value of trace for $N=4$ and $d=2$ is $5+3+3+5=16$.
So far I found the solutions for $N=2$ and $N=3$ and $N=4$, but I can't find any patterns from it.
$$
N=2 : a_d=(-1)^d+1
$$
$$
N=3 : a_d=2^d+(-1)^d+1
$$
$$
N=4 : a_d=3^d+(1+\sqrt2)^d+(1-\sqrt2)^d+1
$$
 A: I can't put the full source code in the comments so I'll just put it here.
I calculated the value in two different ways: one with your formula, and one by just calculating everything.
If you run this code, you'll see that there's difference in $d>=2N$ by varying the $N$ value.
P = 998244353

def matrix_mult(A, B, N):
    temp = [[0] * N for _ in range(N)]
    for i in range(N):
        for j in range(N):
            for k in range(N):
                temp[i][j] += (A[i][k] * B[k][j]) % P
    return temp

def Mtrace(N, d):
    A =[[0 for _ in range(N)] for _ in range(N)]

    for i in range(1, N):
        A[i][i] = N-3
        A[i-1][i] = 1
        A[i][i-1] = 1

    A[0][0] = N-2
    A[N-1][N-1] = N-2

    C = [[0 for _ in range(N)] for _ in range(N)]
    for i in range(N):
        C[i][i] = 1

    l = [N]
    for i in range(d):
        sum = 0
        C = matrix_mult(A, C, N)
        for i in range(N): 
            sum += C[i][i]
        l.append(sum%P)

    return l


def comb(n, k):
    k = min(k, n - k)
    if k == 0:
        return 1

    top = 1
    bot = 1

    for i in range(k):
        top = top * (n - i) % P
        bot = bot * (i + 1) % P

    inv = pow(bot, P-2, P)

    return top * inv % P

def trace(N,d):
    if N==2:
        return 0 if d&1 else 2

    s1 = 0
    for i in range(0, d+1, 2):
        s1 += comb(d,i)*comb(i,i//2)*pow(N-3, d-i)
        s1 %= P

    s1 = (s1 * N) % P

    s2 = pow(N-1, d, P) - pow(N-5, d, P)
    if s2 % 2 == 1:
        s2 = (s2 + P) // 2
    else:
        s2 //= 2
    
    return (s1+s2)%P


N=20
list = Mtrace(N, 101)
for d in range(1, 100):
    sum = trace(N, d)
    print(f'd={d} | diff : {sum-list[d]} | formula : {sum} | matrix : {list[d]}')
sum = Mtrace(N,d)

Here's the new code based on your new edits, and it's perfect!!
P = 998244353

def matrix_mult(A, B, N):
    temp = [[0] * N for _ in range(N)]
    for i in range(N):
        for j in range(N):
            for k in range(N):
                temp[i][j] += (A[i][k] * B[k][j]) % P
    return temp

def Mtrace(N, d):
    A =[[0 for _ in range(N)] for _ in range(N)]

    for i in range(1, N):
        A[i][i] = N-3
        A[i-1][i] = 1
        A[i][i-1] = 1

    A[0][0] = N-2
    A[N-1][N-1] = N-2

    C = [[0 for _ in range(N)] for _ in range(N)]
    for i in range(N):
        C[i][i] = 1

    l = [N]
    for i in range(d):
        sum = 0
        C = matrix_mult(A, C, N)
        for i in range(N): 
            sum += C[i][i]
        l.append(sum%P)

    return l


def comb(n, k):
    k = min(k, n - k)
    if k == 0:
        return 1

    top = 1
    bot = 1

    for i in range(k):
        top = top * (n - i) % P
        bot = bot * (i + 1) % P

    inv = pow(bot, P-2, P)

    return top * inv % P

def trace(N,d):
    if N==2:
        return 0 if d&1 else 2

    s1 = 0
    for i in range(d//2+1):
        tmp = i//N
        ss1=0
        for j in range(-tmp, tmp+1):
            ss1 += comb(2*i, i-N*j)
        s1 += comb(d,2*i)*ss1*pow(N-3, d-2*i)
        s1 %= P

    s1 = (s1 * N) % P

    s2 = pow(N-1, d, P) - pow(N-5, d, P)
    if s2 % 2 == 1:
        s2 = (s2 + P) // 2
    else:
        s2 //= 2
    
    return (s1+s2)%P


N=20
list = Mtrace(N, 101)
for d in range(1, 100):
    sum = trace(N, d)
    print(f'd={d} | diff : {sum-list[d]} | formula : {sum} | matrix : {list[d]}')
sum = Mtrace(N,d)
```

