用for循环实现传染病传播的动画情节(Matlab)

Question

我是Matlab的初学者，我试图用Matlab来模拟传染病的传播。然而，我遇到了一些问题。

首先，我定义了需要填充的矩阵及其初始状态：

diseasematrix=zeros(20,20);
inirow=10;
inicol=10;
diseasematrix(inirow,inicol)=1;         % The first place where a sick person is
infectionmatrix=zeros(20,20);           % Infected people, initially all 0
healthymatrix=round(rand(20,20)*100);   % Initial healthy population (randomly)
Rate=0.0001;                            % Rate of spread

现在，我想用for循环绘制一个图表，显示疾病的传播情况。但我被困在这里..。

for t=1:365
    Zneighboursum=zeros(size(diseasematrix));
    out_ZT = calc_ZT(Zneighboursum, diseasematrix);
    infectionmatrix(t) = round((Rate).*(out_ZT));
    diseasematrix(t) = diseasematrix(t-1) + infectionmatrix(t-1);
    healthymatrix(t) = healthymatrix(t-1) - infectionmatrix(t-1);
    imagesc(diseasematrix(t));
    title(sprintf('Day %i',t));
    drawnow;
end

这基本上是说，感染矩阵是根据循环中的公式计算的，疾病矩阵是通过将上一阶段的病人与前一次感染者相加来计算的。残存的健康人是通过减去前一时间的健康人与感染者的比例来计算的。变量out_ZT是我创建的一个函数：

function [ZT] = calc_ZT(Zneighboursum, diseasematrix)
Zneighboursum = Zneighboursum + circshift(diseasematrix,[1 0]);
Zneighboursum = Zneighboursum + circshift(diseasematrix,[0 1]);
ZT=Zneighboursum;
end

这是为了量化中央细胞周围的病人人数。

然而，结果不是我想要的。情节不会动态发展，而且这些值似乎也不正确。有谁可以帮我？

提前感谢！

Answer 1

代码有几个问题：

1. (Rate).*(out_ZT)错了。因为第一个是标量，第二个是矩阵，而.*要求两者都是相同大小的矩阵。所以一个*就行了。
2. infectionmatrix，diseasematrix，healthymatrix都是二维矩阵，为了将它们保存在内存中，需要有一个三维矩阵。但是既然你以后不使用你储存的东西，你就可以重写旧的了。
3. 在infectionmatrix中存储整数，因为使用round()计算整数。这将结果始终设置为零。
4. Rate的值太低，无法看到任何结果。所以我把它增加到0.01了
5. (只是警告点)在代码中没有使用healthymatrix。

函数的代码很好，所以根据我所看到的进行调试之后，下面是代码：

diseasematrix=zeros(20,20);
inirow=10;
inicol=10;
diseasematrix(inirow,inicol)=1;         % The first place where a sick person is
infectionmatrix=zeros(20,20);           % Infected people, initially all 0
healthymatrix=round(rand(20,20)*100);   % Initial healthy population (randomly)
Rate=0.01;          

for t=1:365
    Zneighboursum=zeros(size(diseasematrix));
    out_ZT = calc_ZT(Zneighboursum, diseasematrix);
    infectionmatrix = (Rate*out_ZT);
    diseasematrix = diseasematrix + infectionmatrix;
    healthymatrix = healthymatrix - infectionmatrix;
    imagesc(diseasematrix);
    title(sprintf('Day %i',t));
    drawnow;
end

Answer 2

有几个问题：

1)如果您想保存一个3D矩阵，您将需要一个3D向量：

所以你必须用myvariable(:,:,t);代替myvariable(:,:,t);

2)你为什么使用round？如果舍入值< 0.5，结果将为0。所以在你的循环中什么都不会改变。

3)需要定义边界条件(t=1)，然后以t= 2启动循环。

diseasematrix=zeros(20,20);
inirow=10;
inicol=10;
diseasematrix(inirow,inicol)=1;         % The first place where a sick person is
infectionmatrix =zeros(20,20);           % Infected people, initially all 0
healthymatrix=round(rand(20,20)*100);   % Initial healthy population (randomly)
Rate=0.01;  % Rate of spread

for t=2:365
    Zneighboursum=zeros(size(diseasematrix,1),size(diseasematrix,2));
    out_ZT = calc_ZT(Zneighboursum, diseasematrix(:,:,t-1));
    infectionmatrix(:,:,t) = (Rate).*(out_ZT);
    diseasematrix(:,:,t) = diseasematrix(:,:,t-1) + infectionmatrix(:,:,t-1);
    healthymatrix(:,:,t) = healthymatrix(:,:,t-1) - infectionmatrix(:,:,t-1);
    imagesc(diseasematrix(:,:,t));
    title(sprintf('Day %i',t));
    drawnow;
end

重要：circshift克隆你的矩阵，以处理边界效应。