Proof that a(n) = 9*9^n + 16*8^n - 24*6^n + 10*5^n - 12*4^n + 6*3^n.

--- Setup ---

For fixed width w, encode each row as a w-tuple a = (a_0,...,a_{w-1}) with entries in {0,...,4}. Define the 5^w x 5^w transfer matrix T_w by T_w[a,b] = 1 if consecutive rows are compatible (every 2x2 subblock sums to 8), else 0. Compatibility forces b_{j+1} = 8-a_j-a_{j+1}-b_j for j=0,...,w-2; given b_0, all other entries are determined. Then a(n, w) = 1^T * T_w^n * 1 counts valid (n+1)×w arrays, and a(n) = a(n, n+1).

Define row statistics pi(a) = (a_0, min_E, max_E, min_O, max_O), where E = {a_{2j} : j >= 1} and O = {a_{2j+1} : j >= 0} (with sentinels for empty sets). There are 1155 possible statistics states, independent of w. Define the 1155-state statistics transfer matrix T_stats: from state s, b_0 is constrained to [lo(s), hi(s)], and each b_0 yields a unique successor state. Define the lift J_w: R^{Stats} -> R^{rows} by (J_w f)(a) = f(pi(a)). Define the column-extension operators A_odd, A_even on R^{Stats}: (A_p f)(s) = sum_{c=0}^{4} f(extend(s, c, p)), where extend updates the relevant min/max statistic with c. Let f_w denote the statistics distribution of all width-w rows.

--- Lemma 1 (statistics quotient) ---

T_w * J_w = J_w * T_stats for all w >= 1.

Proof: For any row b and any valid predecessor a with first entry a_0, the recurrence a_{j+1} = 8-a_j-b_j-b_{j+1} gives by induction:

  a_{2k}   = a_0 + b_0 - b_{2k}                for k >= 1
  a_{2k+1} = 8 - a_0 - b_0 - b_{2k+1}   for k >= 0

(Base case a_1 = 8-a_0-b_0-b_1; inductive steps follow by substitution.) Therefore:

  min_E(a) = a_0 + b_0 - max_E(b),  max_E(a) = a_0 + b_0 - min_E(b)
  min_O(a) = 8 - a_0 - b_0 - max_O(b),  max_O(a) = 8 - a_0 - b_0 - min_O(b)

These depend only on a_0 and pi(b). Hence the multiset {pi(a) : a->b} depends only on pi(b), which gives (T_w * J_w * f)(b) = (J_w * T_stats * f)(b) for all f. QED.

--- Lemma 2 (annihilator) ---

Let p(x) = (x-1)*(x-2)*(x-3)*(x-4)*(x-5). Then p(T_w) * 1 = 0 for every w >= 1.

Proof: An exact integer computation on the fixed 1155-state matrix T_stats confirms p(T_stats) * 1 = 0. By Lemma 1, p(T_w) * 1 = p(T_w) * (J_w * 1) = J_w * (p(T_stats) * 1) = 0. QED.

Since p has 5 distinct roots, the annihilator method gives:

  a(n, w) = sum_{k=1}^{5} alpha_k(w) * k^n

where alpha_k(w) = (1^T * Q_k(w)) / prod_{k'!=k}(k-k'), Q_k(w) = prod_{k'!=k}(T_w - k'*I)*1.

The computed values for w=1,...,5:
  w=1: alpha_1=0, alpha_2=0, alpha_3=0, alpha_4=0, alpha_5=5
  w=2: alpha_1=2, alpha_2=4, alpha_3=6, alpha_4=8, alpha_5=5
  w=3: alpha_1=24, alpha_2=36, alpha_3=36, alpha_4=24, alpha_5=5
  w=4: alpha_1=194, alpha_2=220, alpha_3=150, alpha_4=56, alpha_5=5
  w=5: alpha_1=1320, alpha_2=1140, alpha_3=540, alpha_4=120, alpha_5=5

--- Lemma 3 (width recurrence) ---

Define Q_k' = Q_k / denom_k (rational vectors in R^{1155}), so that alpha_k(w) = f_w . Q_k'. Let u_w = (alpha_1(w), ..., alpha_{5}(w)). There exists a fixed 5x5 rational upper-triangular matrix M1 with diagonal (5, 4, 3, 2, 1) such that u_{w+1} = M1 * u_w for all w >= 1. Explicitly:

  row k=1: [5  1  2/3  1/2  2/5]
  row k=2: [0  4  4/3  1  4/5]
  row k=3: [0  0  3  3/2  6/5]
  row k=4: [0  0  0  2  8/5]
  row k=5: [0  0  0  0  1]

Proof: For each k and column parity p, define the error vector:
  e_{k,p} = A_p^T * Q_k' - sum_j M1[k,j] * Q_j'  (a rational vector in R^{1155}).

The identity u_{w+1} = M1 * u_w at width w is equivalent to f_w . e_{k,p} = 0, where p is the parity of the column added at step w.

Step 1 (candidate subspaces). Let f_w denote the statistics distribution of all width-w rows. By definition of A_p, f_{w+1} = A_p(f_w) where p is the parity of column w+1 (odd if w is odd, even if w is even) — this is exactly what A_p means as a column-extension operator on statistics distributions. Define the candidate subspaces:
  V_odd  = span{f_1, f_3, ..., f_{29}}  (15 specific vectors computed by exact Gaussian elimination over Q)
  V_even = span{f_2, f_4, ..., f_{30}}  (15 specific vectors computed similarly)
These are linear subspaces defined from computed data; their invariance under the dynamics is established in Step 2. These are candidate bases; Step 3 shows that all subsequent f_w remain in their alternating spans.

Step 2 (invariance). An exact rational computation (checking the image of each basis vector under the 1155-dimensional matrix representation of A_odd and A_even) verifies:
  (i)  A_odd(v) in V_even  for each v in V_odd   (15 checks)
  (ii) A_even(v) in V_odd  for each v in V_even  (15 checks)

By linearity: (i) gives A_odd(V_odd) ⊆ V_even, and (ii) gives A_even(V_even) ⊆ V_odd.

Step 3 (orbit confinement). We now prove that all f_w lie in V_odd (odd w) or V_even (even w); no assumption of this was made in Step 1. Recall that by definition of the lifted statistics dynamics, f_{w+1} = A_{parity(w+1)}(f_w) as vectors in R^{1155}. By construction, f_{w+1} = A_odd(f_w) for odd w, and f_{w+1} = A_even(f_w) for even w. By induction on w:
  Base cases (exact): f_1 in V_odd; f_2 = A_odd(f_1) in A_odd(V_odd) ⊆ V_even [by (i)].
  If f_{2m+1} in V_odd: f_{2m+2} = A_odd(f_{2m+1}) in A_odd(V_odd) ⊆ V_even [by (i)].
  If f_{2m}   in V_even: f_{2m+1} = A_even(f_{2m}) in A_even(V_even) ⊆ V_odd [by (ii)].
Hence f_w in V_odd for odd w, and f_w in V_even for even w, for all w >= 1. Thus the parity of w determines the operator independently of the subspace membership, so Step 2 applies at every step without circularity.


Step 4 (orthogonality). Exact rational arithmetic verifies e_{k,odd} . v = 0 for all 15 basis vectors v of V_odd, and e_{k,even} . v = 0 for all 15 basis vectors v of V_even (all k = 1,...,5). By linearity, e_{k,odd} perp V_odd and e_{k,even} perp V_even.

Since f_w in V_odd (odd w) or V_even (even w) by Step 3, and e_{k,odd} perp V_odd and e_{k,even} perp V_even by Step 4, f_w . e_{k,p} = 0 for all w >= 1. Therefore u_{w+1} = M1 * u_w for all w >= 1. QED.

Since M1 is upper triangular with diagonal (5, 4, 3, 2, 1), the k-th coordinate alpha_k(w) satisfies a linear recurrence in w. The active roots for each k are determined by the initial conditions: starting from alpha_5(w) = 5 (constant, root {1}), and propagating upward through the upper-triangular structure. Concretely:

  alpha_1(w) = 1*3^w + -2*4^w + 1*5^w
  alpha_2(w) = 2*2^w + -4*3^w + 2*4^w
  alpha_3(w) = 3*1^w + -6*2^w + 3*3^w
  alpha_4(w) = -8*1^w + 4*2^w
  alpha_5(w) = 5*1^w

--- Conclusion ---

Since a(n) = a(n, n+1), substituting w = n+1 (valid for all n >= 0; the OEIS lists terms starting at n=1):

  a(n) = sum_k alpha_k(n+1) * k^n
       = sum_k sum_j c_{k,j} * rho_{k,j}^(n+1) * k^n
       = sum_k sum_j (c_{k,j} * rho_{k,j}) * (rho_{k,j} * k)^n

The contributing cross-terms:
  k=1, rho=3: 1*3^(n+1)*1^n = 3*3^n
  k=1, rho=4: -2*4^(n+1)*1^n = -8*4^n
  k=1, rho=5: 1*5^(n+1)*1^n = 5*5^n
  k=2, rho=2: 2*2^(n+1)*2^n = 4*4^n
  k=2, rho=3: -4*3^(n+1)*2^n = -12*6^n
  k=2, rho=4: 2*4^(n+1)*2^n = 8*8^n
  k=3, rho=1: 3*1^(n+1)*3^n = 3*3^n
  k=3, rho=2: -6*2^(n+1)*3^n = -12*6^n
  k=3, rho=3: 3*3^(n+1)*3^n = 9*9^n
  k=4, rho=1: -8*1^(n+1)*4^n = -8*4^n
  k=4, rho=2: 4*2^(n+1)*4^n = 8*8^n
  k=5, rho=1: 5*1^(n+1)*5^n = 5*5^n

Collecting by eigenvalue (rho*k):
  a(n) = 9*9^n + 16*8^n - 24*6^n + 10*5^n - 12*4^n + 6*3^n

This formula is verified against 30 computed terms and all 6 OEIS terms.