阻滞增长模型--Logistic模型

写在前面

做了一个作业，觉得挺有意思，有空就发上来了。

阻滞增长模型–Logistic模型的部分推导、思路，并给出参考代码。
参考书籍：《数学模型》第四版，姜启源 谢金星 叶俊 编
MATLAB版本2021b
点我下载源代码及原始数据
提取码：p3uu

一、思路

阻滞作用体现在对人口增长率r的影响上，使得r随着人口数量x的增加而下降。于是有：

对r(x)做假设，设r(x)是x的线性函数，即，
这里r称固有增长率；引入人口容量xm，当x=xm时人口不在增长，于是有，

结合上述三式，有，

可解得，
式5即为Logistic阻滞增长模型。
取数据得第一个时间点得值为x0，即x0=3.9。
由式4可得，

等号左边项可由实际数据计算得出，左端项与右端x呈线性关系，可通过线性拟合得r，s,进而求得x_m。
至此，式5中固有人口增长率r、人口容量x_m、人口初始值x_0均已求得。

二、原始数据

点我下载mat格式原始数据
提取码：75d5
数据在参考书P164

三、源代码

// function_xt.m
function output = function_xt(xm, x0, r, t )
%UNTITLED2 此处提供此函数的摘要
% Logistic模型
% 输入：
% xm:人口容量
% r :固有增长率
% x0:初始人口数
% t :时间点
% 输出
% x（t），时间点对应的人口数
output = xm / ( 1 + ( xm/x0 - 1 ) * exp(-r*t) );
end
1
2
3
4
5
6
7
8
9
10
11
12
13

// main
%%
% 运行环境：matlab2021b
%
% @9490，2022-3-29
%
clear
close all

%% 加载数据并进行 固有增长率、人口容量 的计算
% 加载原始数据
load('data.mat')
% 求微分
dxdt = diff(Data.value);
% 线性拟合中 y 的计算
y = dxdt./Data.value(1:end-1);
% 进行线性拟合
% 返回值p(1)为系数，p(2)为截距
p = polyfit(Data.value(1:end-1), y, 1);
% 对应赋值
s = -p(1);
r = p(2);
xm = r/s;
x0 = Data.value(1); % 取第一个数据点为x0

%% 获取模型关键值后进行计算
for t = 0:1:20
    yhat(t+1) = function_xt(xm,x0,r,t);
end
%% 画图，返回关键值
plot(Data.t(1:end-1), Data.value(1:end-1),'blue')
hold on
plot(Data.t(1:end-1), yhat, '*')
set(gca, 'XLim',[Data.t(1) Data.t(end)]);
xlabel('时间（年）','FontSize',24)
ylabel('人数（百万）','FontSize',24)
legend('原始数据','模型计算','FontSize',24)
title('阻滞增长模型拟合图形','FontSize',24)

fprintf('人口容量xm是：%f\n',xm)
fprintf('固有增长率r是：%f\n',r)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41

四、结果

人口容量xm是：294.385993
固有增长率r是：0.324519