[proofplan]
By the [Hilbert–Serre](/theorems/2952) theorem, the Poincaré series $P(M, T) = \sum_{n \ge 0} \ell(M_n) T^n$ has the form $f(T)/(1 - T)^d$ where $d = d(M)$ and $f \in \mathbb{Z}[T]$ satisfies $f(1) \neq 0$. We expand $(1 - T)^{-d}$ using the binomial series, convolve with the coefficients of $f$, and show that for $n \ge \deg f$, the coefficient of $T^n$ equals a polynomial in $n$ of degree $d - 1$ with leading coefficient $f(1)/(d-1)! \neq 0$. Uniqueness follows from the fact that a polynomial is determined by its values on all sufficiently large integers.
[/proofplan]
[step:Write the Poincaré series in the form $f(T)/(1-T)^d$ with $f(1) \neq 0$]
By the [Hilbert–Serre](/theorems/2952) theorem, since $k_1 = \cdots = k_s = 1$:
\begin{align*}
P(M, T) = \frac{f(T)}{(1 - T)^s}
\end{align*}
for some $f(T) \in \mathbb{Z}[T]$. Write $d = d(M)$, which by definition is the order of the pole of $P(M, T)$ at $T = 1$. Factor out as many copies of $(1 - T)$ from $f(T)$ as possible: write $f(T) = (1 - T)^{s - d} h(T)$ where $h(T) \in \mathbb{Z}[T]$ and $h(1) \neq 0$. Then:
\begin{align*}
P(M, T) = \frac{h(T)}{(1 - T)^d},
\end{align*}
with $h(1) \neq 0$ and $d \ge 0$. Relabel $h$ as $f$ for notational simplicity: $P(M, T) = f(T)/(1 - T)^d$ with $f \in \mathbb{Z}[T]$, $f(1) \neq 0$.
[guided]
The [Hilbert–Serre](/theorems/2952) theorem gives $P(M, T) = f(T)/\prod_{i=1}^s (1 - T^{k_i})$. Since $k_i = 1$ for all $i$, the denominator is $(1 - T)^s$.
The quantity $d(M)$ is defined as the order of the pole of $P(M, T)$ at $T = 1$. The denominator contributes a pole of order $s$ at $T = 1$, and the numerator $f(T)$ may have $(1 - T)$ as a factor, which cancels some poles. If $f(T) = (1 - T)^a g(T)$ with $g(1) \neq 0$, then the pole order is $s - a$, so $d = s - a$. After cancellation, $P(M, T) = g(T)/(1 - T)^d$ with $g(1) \neq 0$.
Why must $d \ge 0$? Because $P(M, T) = \sum_{n \ge 0} \ell(M_n) T^n$ is a power series with non-negative integer coefficients; if $d$ were negative, $P(M, T)$ would be a polynomial vanishing at $T = 1$ to order $|d|$, which would force $M_n = 0$ for all $n$ eventually, making $d = 0$ by convention.
[/guided]
[/step]
[step:Expand using the binomial series and extract the coefficient of $T^n$]
The binomial series for the negative power gives:
\begin{align*}
\frac{1}{(1 - T)^d} = \sum_{j \ge 0} \binom{j + d - 1}{d - 1} T^j,
\end{align*}
where $\binom{j + d - 1}{d - 1} = \frac{(j+d-1)(j+d-2) \cdots (j+1)}{(d-1)!}$ is a polynomial in $j$ of degree $d - 1$ with leading coefficient $1/(d-1)!$.
Write $f(T) = \sum_{i=0}^{N} a_i T^i$ where $N = \deg f$. The coefficient of $T^n$ in $P(M, T) = f(T) \cdot (1 - T)^{-d}$ is the Cauchy product:
\begin{align*}
\ell(M_n) = \sum_{i=0}^{\min(n, N)} a_i \binom{n - i + d - 1}{d - 1}.
\end{align*}
For all $n \ge N$, the upper limit of summation equals $N$, so:
\begin{align*}
\ell(M_n) = \sum_{i=0}^{N} a_i \binom{n - i + d - 1}{d - 1}.
\end{align*}
[guided]
The generating function identity $P(M, T) = f(T)/(1-T)^d$ means we need to expand $(1-T)^{-d}$ as a power series and multiply by $f(T)$. The binomial series gives:
\begin{align*}
(1 - T)^{-d} = \sum_{j=0}^{\infty} \binom{j + d - 1}{d - 1} T^j.
\end{align*}
This is the standard identity for negative binomial coefficients: $\binom{-d}{j}(-1)^j = \binom{j+d-1}{d-1}$.
Writing $f(T) = \sum_{i=0}^{N} a_i T^i$, the coefficient of $T^n$ in the product $f(T) \cdot (1-T)^{-d}$ is obtained by the Cauchy convolution:
\begin{align*}
[T^n] P(M,T) = \sum_{i=0}^{N} a_i \cdot [T^{n-i}] (1-T)^{-d} = \sum_{i=0}^{N} a_i \binom{n - i + d - 1}{d-1},
\end{align*}
provided $n \ge N$, so that $n - i \ge 0$ for all $0 \le i \le N$. For $n < N$, some of the binomial coefficients $\binom{n-i+d-1}{d-1}$ have $n - i < 0$, which complicates the formula. But for $n \ge N$, the formula is a clean sum.
[/guided]
[/step]
[step:Identify the sum as a polynomial of degree $d - 1$ in $n$ with non-vanishing leading coefficient]
Define:
\begin{align*}
\operatorname{HP}_M: \mathbb{Q} &\to \mathbb{Q}, \\
n &\mapsto \sum_{i=0}^{N} a_i \binom{n - i + d - 1}{d - 1}.
\end{align*}
Each binomial coefficient $\binom{n - i + d - 1}{d - 1}$ is a polynomial in $n$ of degree $d - 1$. Indeed:
\begin{align*}
\binom{n - i + d - 1}{d - 1} = \frac{(n - i + d - 1)(n - i + d - 2) \cdots (n - i + 1)}{(d-1)!},
\end{align*}
which is a product of $d - 1$ linear factors in $n$ divided by $(d-1)!$. The leading term (the coefficient of $n^{d-1}$) is $1/(d-1)!$ for every $i$.
Therefore $\operatorname{HP}_M(n)$ is a polynomial in $n$ of degree at most $d - 1$. The leading coefficient is:
\begin{align*}
\frac{1}{(d-1)!} \sum_{i=0}^{N} a_i = \frac{f(1)}{(d-1)!}.
\end{align*}
Since $f(1) \neq 0$, the degree of $\operatorname{HP}_M$ is exactly $d - 1$.
[guided]
Why is $\operatorname{HP}_M$ a polynomial in $n$? Each summand $a_i \binom{n - i + d - 1}{d-1}$ is a polynomial in $n$: it is $a_i$ times a product of $d-1$ linear functions of $n$ (namely $n - i + 1, n - i + 2, \ldots, n - i + d - 1$), divided by $(d-1)!$. A finite sum of polynomials in $n$ is a polynomial in $n$, so $\operatorname{HP}_M \in \mathbb{Q}[n]$.
The degree is at most $d - 1$ since each summand has degree $d - 1$. Could cancellation reduce the degree? The leading coefficient of $a_i \binom{n - i + d - 1}{d-1}$ as a polynomial in $n$ is $a_i/(d-1)!$, regardless of $i$ (the shift by $i$ affects lower-order terms but not the leading term). Summing over $i$:
\begin{align*}
\text{leading coefficient of } \operatorname{HP}_M = \frac{1}{(d-1)!} \sum_{i=0}^{N} a_i = \frac{f(1)}{(d-1)!}.
\end{align*}
Since $f(1) \neq 0$ by construction, the leading coefficient is nonzero, confirming $\deg(\operatorname{HP}_M) = d - 1$.
[/guided]
[/step]
[step:Conclude: $\ell(M_n) = \operatorname{HP}_M(n)$ for large $n$, and uniqueness]
By the computation in the previous steps, for all $n \ge N = \deg f$:
\begin{align*}
\ell(M_n) = \operatorname{HP}_M(n),
\end{align*}
where $\operatorname{HP}_M \in \mathbb{Q}[T]$ has degree $d(M) - 1$.
For uniqueness: suppose $Q \in \mathbb{Q}[T]$ also satisfies $\ell(M_n) = Q(n)$ for all sufficiently large $n$. Then $\operatorname{HP}_M(n) = Q(n)$ for all sufficiently large $n$, so the polynomial $\operatorname{HP}_M - Q$ has infinitely many roots (all sufficiently large integers). A nonzero polynomial has only finitely many roots, so $\operatorname{HP}_M - Q = 0$, i.e., $\operatorname{HP}_M = Q$.
[/step]