ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

数模-微分方程(人口预测之马尔萨斯模型和阻滞增长模型)

2022-05-08 10:35:40  阅读:517  来源: 互联网

标签:人口预测 xm 模型 拟合 数模 exp year x0 population


模型:

image
image

代码:

%% Malthus模型(马尔萨斯模型)
clear;clc
x = dsolve('Dx=r*x','x(0)=x0','t')    % x = dsolve('Dx=r*x','x(t0)=x0','t')
% x = x0*exp(r*t)
% 怎么把上面这个式子中的x0和r替换成确定的值?
x0 = 100; 
r = 0.1;
subs(x)

% 初始人口为1000,画出50年且增长率分别为0.5%,1%,1.5% 一直到5%的人口变化曲线
x0 = 1000;
for i = 1:10
    r = 0.005*i;
    xx = subs(x);
    fplot(xx,[0,50])   % fplot函数可以绘制表达式的图形
    hold on
end


%% 阻滞增长模型(logistic模型)
clear;clc
x = dsolve('Dx=r*(1-x/xm)*x','x(t0)=x0','t')   % 化简后和书上的结果一样
% mupad  % 把计算出来的结果粘贴过去可以得到直观的表达式
% 高版本可以使用实时脚本
t0 = 0;
x0 = 1000;
xm = 10000;
r = 0.05;
xx = subs(x);    %  10000/(exp(log(9) - t/20) + 1)
fplot(xx,[0,200])  

% 如果不会用上面的fplot函数怎么办?
t = 0:200;
x = 10000 ./ (exp(log(9) - t/20) + 1);
plot(t,x,'-')

image

image

例题

image

求解过程

image
image
image
image

代码

createFit.m

function [fitresult, gof] = createFit(year, population)
[xData, yData] = prepareCurveData( year, population );

% Set up fittype and options.
ft = fittype( 'xm/(1+(xm/3.9-1)*exp(-r*(t-1790)))', 'independent', 't', 'dependent', 'x' );
opts = fitoptions( 'Method', 'NonlinearLeastSquares' );
opts.Display = 'Off';
opts.StartPoint = [0.2 500];

% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft, opts );

% Plot fit with data.
figure( 'Name', 'untitled fit 1' );
h = plot( fitresult, xData, yData );
legend( h, 'population vs. year', 'untitled fit 1', 'Location', 'NorthEast' );
% Label axes
xlabel year
ylabel population
grid on

code.m

clear;clc
year = 1790:10:2000;
population = [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,106.5,123.2,131.7,150.7,179.3,204.0,226.5,251.4,281.4];
cftool  % 拟合工具箱
% (1) X data 选择 year
% (2) Y data 选择 population
% (3) 拟合方式选择:Custom Equation (自定义方程)
% (4) 修改下方的方框为:x = f(t) = xm/(1+(xm/3.9-1)*exp(-r*(t-1790)))
% (5) 左边的result一栏最上面显示:Fit computation did not converge:即没有找到收敛解,右边的拟合图形也表明拟合结果不理想
% (6) 点击Fit Options,修改非线性最小二乘估计法拟合的初始值(StartPoint), r修改为0.02,xm修改为500
% (7) 此时左边的result一览得到了拟合结果:r = 0.02735, xm = 342.4
% (8) 依次点击拟合工具箱的菜单栏最左边的文件—Generate Code(导出代码到时候可以放在你的论文附录),可以得到一个未命名的脚本文件
% (9) 在这个打开的脚本中按快捷键Ctrl+S,将这个文件保存到当前文件夹。
% (10) 在现在这个文件中调用这个函数得到参数的拟合值和预测的效果
[fitresult, gof] = createFit(year, population) 
t = 2001:2030;
xm = 342.4;   
r =  0.02735;
predictions = xm./(1+(xm./3.9-1).*exp(-r.*(t-1790)));  % 计算预测值(注意这里要写成点乘和点除)
figure(2)
plot(year,population,'o',t,predictions,'.')  % 绘制预测结果图
disp(predictions)  % 预测的数值

image

标签:人口预测,xm,模型,拟合,数模,exp,year,x0,population
来源: https://www.cnblogs.com/jgg54335/p/15181619.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有