1. 程式人生 > 其它 >數模-微分方程(SIS模型)

數模-微分方程(SIS模型)

SIS模型

程式碼

fun1.m

function dx=fun1(t,x)   % 大家可以修改裡面的引數,來看結果的變化
    global TOTAL_N   % 定義總人數為全域性變數
    beta = 0.1;  % 易感染者與已感染者接觸且被傳染的強度
    alpha = 0.06; % 由感染狀態I恢復為易感者狀態S的恢復率
    dx = zeros(2,1);  % x(1)表示S  x(2)表示I
    dx(1) = alpha*x(2) - beta*x(1)*x(2)/TOTAL_N;  
    dx(2) = beta*x(1)*x(2)/TOTAL_N - alpha*x(2);
end

code.m

%% 最簡單的SIS模型
clc;clear
global TOTAL_N   % 定義總人數為全域性變數(可以在子函式中使用)
TOTAL_N = 1000;  % 總人數
i0 = 1; % 初始時刻患者(已感染者)的人數
[t,x]=ode45('fun1',[1:500],[TOTAL_N-i0 i0]); 
x = round(x);  % 對x進行四捨五入(人數為整數)
figure(1)
plot(t,x(:,1),'r-',t,x(:,2),'b-','Linewidth',1.5)   % x的第一列是易感染者S的數量,x的第二列是患者I的數量
legend('易感染者S','患者I')

結果

SIS模型擴充套件

程式碼

fun2.m

function dx=fun2(t,x)   % 大家可以修改裡面的引數,來看結果的變化
    global TOTAL_N   % 定義總人數為全域性變數
    beta = 0.1;  % 易感染者與已感染者接觸且被傳染的強度
    alpha = 0.02; % 由感染狀態I恢復為易感者狀態S的恢復率
    if t > 200
        alpha = alpha * 10; % 第200期後引入了新的醫療裝備,使得恢復率alpha增加為原來的10倍
    end
    dx = zeros(2,1);  % x(1)表示S  x(2)表示I
    dx(1) = alpha*x(2) - beta*x(1)*x(2)/TOTAL_N;  
    dx(2) = beta*x(1)*x(2)/TOTAL_N - alpha*x(2);
end

code.m

%% 考慮某種使得恢復率alpha增加的因素(例如建立醫院、升級醫療裝備等)
% 第200期後引入了新的醫療裝備,使得恢復率alpha增加為原來的10倍
clc;clear
global TOTAL_N   % 定義總人數為全域性變數(可以在子函式中使用)
TOTAL_N = 1000;  % 總人數
i0 = 1; % 初始時刻患者(已感染者)的人數
[t,x]=ode45('fun2',[1:500],[TOTAL_N-i0 i0]); 
x = round(x);  % 對x進行四捨五入(人數為整數)
figure(2)
plot(t,x(:,1),'r-',t,x(:,2),'b-','Linewidth',1.5)   % x的第一列是易感染者S的數量,x的第二列是患者I的數量
legend('易感染者S','患者I')
axis([0 500 0 1000])  % 設定座標軸範圍,x軸為0-500 y軸為0-1000

結果