logistic模型在中国人口预测应用

经济管理学院 08金融 孔繁哲 0865147113

Logistic模型在中国人口增长中的应用

1. 统计数据

本文数据来源于国家统计局的《中国统计年鉴2009》，取每年的全国总人口数，时间段为从1990年到2008年。

年份 （x） 人口数（y） 114333 1990（0） 115823 1991（1） 117171 1992（2） 118517 1993（3） 119850 1994（4） 121121 1995（5） 122389 1996（6） 123626 1997（7） 124761 1998（8） 125786 1999（9） 2000（10） 126743 2001（11） 127627 2002（12） 128453 2003（13） 129227 2004（14） 129988 2005（15） 130756 2006（16） 131448 2007（17） 132129 2008（18） 132802 注：人口数单位为：万人。取1990年为基准年第0年，以后年份以此类推。

2. 模型建立

2.1 Logistic模型

Logistic模型是1938年Verhulst．Pearl在修正非密度方程时提出的，又称自我抑制性方程。用群体中的发病的普遍率或严重度表示病害数量（x），将环境最大容量k定为1（100），则其微分方程为：

dxrx(1x) dt 式中r为速率参数，来源于实际调查是观察到的症状明显的病害，范.德.普朗克（1963）将r称作表观侵染速率，该方程与指数模型的主要不同之处，是方程的右边增加了（1-x）修正因子，使模型包含着我抑制作用。 2.2 模型应用

本文将模型用于预测我国人口增长。此处用x表示时间t时刻我国人口数量，r表示人口增长率，r随着x的增长而下降。假设我国总人口容量为xm，r为x的减函数，r(x)rsx，r，s>0。其中，r为固有增长率，即当人口为Xo时的增长率。设t时刻人口为x(t)，其中x(0)x0。则有：

x(t+t)-x(t)=r(x)x(t)t

令t0得方程：

dxr(x)x(t),x(0)x0 （1） dt又r(x)rsx （2）

当xxm时，人口不再增长，即r(xm)0代人（2）得sr/xm，于是（2）

x变为 r()r(1) （3）

xm将（3）代入（1）得

dxxrx(1),x(0)x0 （4） dtxmx)作为修正因子表示资源xm方程（4）右端因子rx表示人口增长趋势，因子(1环境等对人口增长的阻滞作用。对方程（4）求解得：

x(t)xm （5）

xm1(-1)e-rtx03. 参数估计

确定模型后需要对方程中的参数进行估计，这里需要估计的参数为人口最大容量xm和人口增长速率r。根据第一部分给出的数据，结合模型通过Matlab进行关于参数的最小二乘估计，从而确定各参数量。最小二乘估计是利用各观测点到模型的垂直距离最小来确定的，即求下式最小值：

xmE(xm,r)=(yi)2

xi1(m1)ertix0n

4. 数据分析

通过Matlab解决最小二乘法曲线拟合问题，求得：xm=14.2190（亿人），r=0.0685，同时求得的最小二乘的残差值为e=0.0011。

把xm和r值代入（5）式预测2010年人口数量为：y=13.3911（亿人） 通过matlab得到的拟合曲线和观测散点图如下图所示：

代码附录：

（1）population

function y2010=population

x=[0 1 2 3 4 5 6 7 8 9 0 11 12 13 14 15 16 17 18];

y=[11.4333 11.5823 11.7171 11.8517 11.9850 12.1121 12.2389 12.3626 12.4761 12.5786 12.6743 12.7627 12.8453 12.9227 12.9988 13.0756 13.1448 13.2129 13.2802];

x0=[20 0.003];

[a,resnorm]=lsqcurvefit(@logistic1,x0,x,y) x1=linspace(0,25,100); b=11.4333;

y1=a(1)./(1+(a(1)/b-1).*exp(-a(2)*x1)); plot(x1+1990,y1); hold

plot(x+1990,y,'r*') xlabel('年份/年'); ylabel('人口/亿人'); title('中国人口预测');

legend('拟合曲线','实际值');

y2010=a(1)./(1+(a(1)/b-1).*exp(-a(2)*20));

（2）logistic1

function y=logistic1(a,x) b=11.4333;

y=a(1)./(1+(a(1)/b-1).*exp(-a(2)*x))