Given a regular language $L$, consider some DFA accepting $L$, let $A$ be its transfer matrix ($A_{ij}$ is the number of edges leading from state $i$ to state $j$), let $x$ be the characteristic vector of the initial state, and let $y$ be the characteristic vector of the accepting states. Then
$$ s_L(n) = x^T A^n y. $$
Jordan's theorem states that over the complex numbers, $A$ is similar to a matrix with blocks of one of the forms
$$ \begin{pmatrix} \lambda \end{pmatrix}, \begin{pmatrix} \lambda & 1 \\ 0 & \lambda \end{pmatrix}, \begin{pmatrix} \lambda & 1 & 0 \\ 0 & \lambda & 1 \\ 0 & 0 & \lambda \end{pmatrix}, \begin{pmatrix} \lambda & 1 & 0 & 0 \\ 0 & \lambda & 1 & 0 \\ 0 & 0 & \lambda & 1 \\ 0 & 0 & 0 & \lambda \end{pmatrix}, \ldots $$
If $\lambda \neq 0$, then the $n$th powers of these blocks are
$$ \begin{pmatrix} \lambda^n \end{pmatrix}, \begin{pmatrix} \lambda^n & n\lambda^{n-1} \\ 0 & \lambda^n \end{pmatrix}, \begin{pmatrix} \lambda^n & n\lambda^{n-1} & \binom{n}{2} \lambda^{n-2} \\ 0 & \lambda^n & n\lambda^{n-1} \\ 0 & 0 & \lambda^n \end{pmatrix}, \begin{pmatrix} \lambda^n & n\lambda^{n-1} & \binom{n}{2}\lambda^{n-2} & \binom{n}{3}\lambda^{n-3} \\ 0 & \lambda^n & n\lambda^{n-1} & \binom{n}{2}\lambda^{n-2} \\ 0 & 0 & \lambda^n & n\lambda^{n-1} \\ 0 & 0 & 0 & \lambda^n \end{pmatrix}, \ldots $$
Here's how we got to these formulas: write the block as $B = \lambda + N$. Successive powers of $N$ are successive secondary diagonals of the matrix. Using the binomial theorem (using the fact that $\lambda$ commutes with $N$),
$$ B^n = (\lambda + n)^N = \lambda^n + n \lambda^{n-1} N + \binom{n}{2} \lambda^{n-2} N^2 + \cdots. $$
When $\lambda = 0$, the block is nilpotent, and we get the following matrices (the notation $[n = k]$ is $1$ if $n=k$ and $0$ otherwise):
$$ \begin{pmatrix} [n=0] \end{pmatrix}, \begin{pmatrix} [n=0] & [n=1] \\ 0 & [n=0] \end{pmatrix}, \begin{pmatrix} [n=0] & [n=1] & [n=2] \\ 0 & [n=0] & [n=1] \\ 0 & 0 & [n=0] \end{pmatrix}, \begin{pmatrix} [n=0] & [n=1] & [n=2] & [n=3] \\ 0 & [n=0] & [n=1] & [n=2] \\ 0 & 0 & [n=0] & [n=1] \\ 0 & 0 & 0 & [n=0] \end{pmatrix} $$
Summarizing, every entry in $A^n$ is either of the form $\binom{n}{k} \lambda^{n-k}$ or of the form $[n=k]$, and we deduce that
$$ s_L(n) = \sum_i p_i(n) \lambda_i^n + \sum_j c_j [n=j], $$
for some complex $\lambda_i,c_j$ and complex polynomials $p_i$. In particular, for large enough $n$,
$$ s_L(n) = \sum_i p_i(n) \lambda_i^n. $$
This is the precise statement of the result.
We can go on and obtain asymptotic information about $s_L(n)$, but this is surprisingly non-trivial. If there is a unique $\lambda_i$ of largest magnitude, say $\lambda_1$, then
$$ s_L(n) = p_1(n) \lambda_1^n (1 + o(1)). $$
Things get more complicated when there are several $\lambda$s of largest magnitude. It so happens that their angle must be rational (i.e. up to magnitude, they are roots of unity). If the LCM of the denominators is $d$, then the asymptotics of $s_L$ will very according to the remainder of $n$ modulo $d$. For some of these remainders, all $\lambda$s of largest magnitude cancel, and then the asymptotics "drops", and we have to iterate this procedure. The interested reader can check the details in Flajolet and Sedgewick's Analytic Combinatorics, Theorem V.3. They prove that for some $d$, integers $p_0,\ldots,p_{d-1}$ and reals $\lambda_0,\ldots,\lambda_{d-1}$,
$$ s_L(n) = n^{p_{n\pmod{d}}} \lambda_{n\pmod{d}}^n (1 + o(1)). $$
