数学建模：Malthus人口模型

clear all 
clc
dsolve('Dx=r*x','x(0)=x0');%指数增长模型，也称为Malthus人口模型
%得到 x（t）=x0*exp(r*t)

%美国人口数据
t=0:1:21;
x=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76.0,92.0,105.7,122.8,131.7,150.7,179.3,203.2,226.5,248.7,281.4];

%方法一,直接用人口数据和最小二乘法
y=log(x);
x0=6.0496;
a=log(x0);
r1=(y-a)/t;
disp("r="+ r1 +"/10年")
y1=x0*exp(r1*t);
subplot(2,2,1);
plot(t,x,'*r',t,y1,'-b')
%方法二
xx(1)=(-3*x(1)+4*x(2)-x(3))/2;
xx(22)=(x(20)-4*x(21)+3*x(22))/2;
for k=2:21
    xx(k)=(x(k+1)-x(k-1))/2;
end
r=xx./x;
ar=sum(r)/22
y2=x(1)*exp(ar*t);
subplot(2,2,2)
plot(t,x,'*r',t,y2,'-g')
ar=ar*ones(size(t));
%改进的指数增长模型
subplot(2,2,3)
plot(t,r,'*r',t,ar,'b-')
dsolve('Dx=(r0+r1*t)*x','x(0)=x0','t')

p1=polyfit(t,r,1)%对增长率r产生拟合参数
y3=polyval(p1,t);
hold on
plot(t,y3,'-b')
r0=p1(2);
r1=p1(1);
y4=x(1)*exp((r1*t.^2)/2 + r0*t);
subplot(2,2,4)
plot(t,x,'*r',t,y4,'-b')