数模(6)：Leslie矩阵人口模型

上期中介绍了两种利用非线性函数拟合人口与物种增长趋势的方法。这两种方法都可以用于对人口与物种增长的总体趋势进行预测，但预测不够精细。我们知道在正常社会条件或自然条件下，生育率与死亡率是与群体的年龄构成息息相关的。我们需要对整个群体按年龄进行层次划分，构建与年龄相联系的人口模型。典型的例子就是Leslie矩阵模型。

Leslie矩阵介绍

我们把整个社会中的人群按年龄等距分成n组，每组中该年的人口总数为ai,i=1,2,...,na_i,i=1,2,...,nai,i=1,2,...,n，每组人口的每年的普遍存活率为ci,i=1,2,...,n−1c_i,i=1,2,...,n-1ci,i=1,2,...,n−1(设最后一组下一年全部死亡)，每组人口的每年普遍生育率为bi,i=1,2,...,nb_i,i=1,2,...,nbi,i=1,2,...,n，则下一年每组中的人口总数ai′,i=1,2,...,na'_i,i=1,2,...,nai′,i=1,2,...,n就满足递推关系式{ai′=ai−1ci−1,i=2,3,...,na1′=∑i=1naibi\begin{cases}a'_i=a_{i-1}c_{i-1},i=2,3,...,n\\a'_1=\sum_{i=1}^{n}a_ib_i\end{cases}{ai′=ai−1ci−1,i=2,3,...,na1′=∑i=1naibi

该式可写成矩阵乘向量的形式：
a′⃗=(b1b2...bn−1bnc10...000c2...00⋮⋮⋮00...cn−10)(a1,a2,...,an)T\vec{a'}= \left( \begin{matrix} b_1&b_2&...&b_{n-1}&b_n\\ c_1&0&...&0&0\\ 0&c_2&...&0&0\\ \vdots&\vdots&&\vdots\\ 0&0&...&c_{n-1}&0 \end{matrix} \right) (a_1,a_2,...,a_n)^T a′=⎝⎜⎜⎜⎜⎜⎛b1c10⋮0b20c2⋮0............bn−100⋮cn−1bn000⎠⎟⎟⎟⎟⎟⎞(a1,a2,...,an)T

该式中左边的矩阵就是Leslie矩阵。

Leslie矩阵性质

Leslie矩阵有唯一的单重正特征值λ1\lambda_1λ1，对应的特征向量x⃗1=(1,b1/λ1,c1c2/λ12,...,c1c2...cn−1/λ1n−1)T\vec x_1=(1,b_1/\lambda_1,c_1c_2/\lambda_1^2,...,c_1c_2...c_{n-1}/\lambda_1^{n-1})^Tx1=(1,b1/λ1,c1c2/λ12,...,c1c2...cn−1/λ1n−1)T

证明：设n阶的该矩阵为Ln，n阶的特征多项式为Pn，则有
Pn=∣λI−Ln∣=∣λ−b1−b2...−bn−1−bn−c1λ...000−c2...00⋮⋮⋮00...−cn−1λ∣P_n=|\lambda I-L_n|= \left| \begin{matrix} \lambda-b_1&-b_2&...&-b_{n-1}&-b_n\\ -c_1&\lambda&...&0&0\\ 0&-c_2&...&0&0\\ \vdots&\vdots&&\vdots&&\\ 0&0&...&-c_{n-1}&\lambda \end{matrix} \right| Pn=∣λI−Ln∣=∣∣∣∣∣∣∣∣∣∣∣λ−b1−c10⋮0−b2λ−c2⋮0............−bn−100⋮−cn−1−bn00λ∣∣∣∣∣∣∣∣∣∣∣
=>Pn=λ∣λ−b1−b2...−bn−2−bn−1−c1λ...000−c2...00⋮⋮⋮00...−cn−2λ∣+cn−1∣λ−b1−b2...−bn−3−bn−1−c1λ...000−c2...00⋮⋮⋮00...−cn−20∣=>P_n=\lambda \left| \begin{matrix} \lambda-b_1&-b_2&...&-b_{n-2}&-b_{n-1}\\ -c_1&\lambda&...&0&0\\ 0&-c_2&...&0&0\\ \vdots&\vdots&&\vdots&&\\ 0&0&...&-c_{n-2}&\lambda \end{matrix} \right|+c_{n-1} \left| \begin{matrix} \lambda-b_1&-b_2&...&-b_{n-3}&-b_{n-1}\\ -c_1&\lambda&...&0&0\\ 0&-c_2&...&0&0\\ \vdots&\vdots&&\vdots&&\\ 0&0&...&-c_{n-2}&0 \end{matrix} \right| =>Pn=λ∣∣∣∣∣∣∣∣∣∣∣λ−b1−c10⋮0−b2λ−c2⋮0............−bn−200⋮−cn−2−bn−100λ∣∣∣∣∣∣∣∣∣∣∣+cn−1∣∣∣∣∣∣∣∣∣∣∣λ−b1−c10⋮0−b2λ−c2⋮0............−bn−300⋮−cn−2−bn−1000∣∣∣∣∣∣∣∣∣∣∣
=>Pn=λPn−1+cn−1(−bn−1)(c1c2...cn−2)(−1)n−2(−1)n−2=>Pn=λPn−1−bn−1c1c2...cn−1=>P_n=\lambda P_{n-1}+c_{n-1}(-b_{n-1})(c_1c_2...c_{n-2})(-1)^{n-2}(-1)^{n-2}=> P_n=\lambda P_{n-1}-b_{n-1}c_1c_2...c_{n-1} =>Pn=λPn−1+cn−1(−bn−1)(c1c2...cn−2)(−1)n−2(−1)n−2=>Pn=λPn−1−bn−1c1c2...cn−1
=>Pn=λPn−1−βn−1=>Pn=λn−β1λn−1−β2λn−2−...−βn=>=>P_n=\lambda P_{n-1}-\beta_{n-1}=> P_n=\lambda^n-\beta_1\lambda^{n-1}-\beta_2\lambda^{n-2}-...-\beta_n=> =>Pn=λPn−1−βn−1=>Pn=λn−β1λn−1−β2λn−2−...−βn=>
0=λn−β1λn−1−β2λn−2−...−βn=>Pn=λn−β1λn−1−β2λn−2−...−βn=>0=\lambda^n-\beta_1\lambda^{n-1}-\beta_2\lambda^{n-2}-...-\beta_n=>P_n=\lambda^n-\beta_1\lambda^{n-1}-\beta_2\lambda^{n-2}-...-\beta_n=> 0=λn−β1λn−1−β2λn−2−...−βn=>Pn=λn−β1λn−1−β2λn−2−...−βn=>
1=β1λ−1+β2λ−2+...+βnλ−n1=\beta_1\lambda^{-1}+\beta_2\lambda^{-2}+...+\beta_n\lambda^{-n} 1=β1λ−1+β2λ−2+...+βnλ−n右边的函数是单调连续减函数，且λ\lambdaλ无穷大时趋近0、λ\lambdaλ趋近于0时趋近正无穷，所以有唯一正特征根λ1\lambda_1λ1，对应的特征向量为x⃗1=(1,b1/λ1,c1c2/λ12,...,c1c2...cn−1/λ1n−1)T\vec x_1=(1,b_1/\lambda_1,c_1c_2/\lambda_1^2,...,c_1c_2...c_{n-1}/\lambda_1^{n-1})^Tx1=(1,b1/λ1,c1c2/λ12,...,c1c2...cn−1/λ1n−1)T

所有负的特征值都满足∣λ∣<λ1|\lambda|<\lambda_1∣λ∣<λ1，称λ1\lambda_1λ1为严格优势特征值

证明:设有特征值满足∣λ∣≥λ1=>λ≥−λ1|\lambda|\geq\lambda_1=>\lambda\geq-\lambda_1∣λ∣≥λ1=>λ≥−λ1，则有其依然满足1=β1λ−1+β2λ−2+...+βnλ−n1=\beta_1\lambda^{-1}+\beta_2\lambda^{-2}+...+\beta_n\lambda^{-n}1=β1λ−1+β2λ−2+...+βnλ−n ，而 1=β1λ1−1+β2λ1−2+...+βnλ1−n≥β∣λ−1∣+β∣λ−2∣+...+β∣λ−n∣>βλ−1+βλ−2+...+βλ−n1=\beta_1\lambda_1^{-1}+\beta_2\lambda_1^{-2}+...+\beta_n\lambda_1^{-n} \geq\beta|\lambda^{-1}|+\beta|\lambda^{-2}|+...+\beta|\lambda^{-n}|>\beta\lambda^{-1}+\beta\lambda^{-2}+...+\beta\lambda^{-n}1=β1λ1−1+β2λ1−2+...+βnλ1−n≥β∣λ−1∣+β∣λ−2∣+...+β∣λ−n∣>βλ−1+βλ−2+...+βλ−n，矛盾

对于任意人口分布向量x⃗\vec xx，其迭代k次后的结果有lim⁡k−>+∞x⃗(k)λ1k=cx⃗1\displaystyle \lim_{k->+∞} \frac{\vec x^{(k)}}{\lambda_1^k}=c\vec x_1k−>+∞limλ1kx(k)=cx1(c为常数)，即迭代了无穷多次时，人口的分布比例趋近于特征向量x⃗1\vec x_1x1，而人口增长率趋近于特征值λ1\lambda_1λ1

证明：仅对可化为对角阵的情况进行证明（一般情况需要用到约旦标准型）。lim⁡k−>+∞x⃗(k)λ1k=lim⁡k−>+∞Lkx⃗(0)λ1k=lim⁡k−>+∞(Pdiag(λ1,λ2,...,λn)P−1)kx⃗(0)λ1k=lim⁡k−>+∞Pdiag(λ1k,λ2k,...,λnk)P−1x⃗(0)λ1k=lim⁡k−>+∞Pdiag(1,λ2k/λ1k,...,λnk/λ1k)P−1x⃗(0)\displaystyle \lim_{k->+∞} \frac{\vec x^{(k)}}{\lambda_1^k}=\displaystyle \lim_{k->+∞} \frac{L^k\vec x^{(0)}}{\lambda_1^k}=\displaystyle \lim_{k->+∞} \frac{(Pdiag(\lambda_1,\lambda_2,...,\lambda_n)P^{-1})^k\vec x^{(0)}}{\lambda_1^k}=\displaystyle \lim_{k->+∞} \frac{Pdiag(\lambda_1^k,\lambda_2^k,...,\lambda_n^k)P^{-1}\vec x^{(0)}}{\lambda_1^k}=\displaystyle \lim_{k->+∞} Pdiag(1,\lambda_2^k/\lambda_1^k,...,\lambda_n^k/\lambda_1^k)P^{-1}\vec x^{(0)}k−>+∞limλ1kx(k)=k−>+∞limλ1kLkx(0)=k−>+∞limλ1k(Pdiag(λ1,λ2,...,λn)P−1)kx(0)=k−>+∞limλ1kPdiag(λ1k,λ2k,...,λnk)P−1x(0)=k−>+∞limPdiag(1,λ2k/λ1k,...,λnk/λ1k)P−1x(0)，由于λ1\lambda_1λ1为严格优势特征值，有原式=lim⁡k−>+∞Pdiag(1,0,...,0)P−1x⃗(0)=(x⃗1,x⃗2,...,x⃗n)diag(1,0,...,0)(x⃗1′,x⃗2′,...,x⃗n′)(a1,a2,...,an)T=cx⃗1原式=\displaystyle \lim_{k->+∞} Pdiag(1,0,...,0)P^{-1}\vec x^{(0)}=(\vec x_1,\vec x_2,...,\vec x_n)diag(1,0,...,0)(\vec x'_1,\vec x'_2,...,\vec x'_n)(a_1,a_2,...,a_n)^T=c\vec x_1原式=k−>+∞limPdiag(1,0,...,0)P−1x(0)=(x1,x2,...,xn)diag(1,0,...,0)(x1′,x2′,...,xn′)(a1,a2,...,an)T=cx1

总结

列出Leslie矩阵，我们即可对人口年龄分布进行迭代。且无论一开始的人口分布向量如何，人口比例在迭代无数次之后总趋近于特征向量x⃗1\vec x_1x1。而人口增长率趋近于特征值λ1\lambda_1λ1，说明特征值λ1\lambda_1λ1可以用于预测人口增长速度，对于计生有重要意义。