新型冠状病毒传播规律离散微观模型（结果与实际情况一致）-附源码

新型冠状病毒传播规律离散微观模型

0 前言
1、模型
- 1.1 个体病毒感染概率计算模型
- 1.2 个体移动模型
- 1.3 求解过程
2 实例分析
- 2.1 正常活动传播模拟
- 2.2 隔离模拟
- 2.3 频繁活动模拟
3 结论和展望
- 3.1 结论
- 3.2 展望
4 特别说明
5 计算程序Python代码（参数需要自行调）

0 前言

大年三十，全家自我隔离，各路关于新型冠状病毒的消息刷屏，有喜有忧。老裴完全没有心思过年，也无法工作，但自己又无能为力。老裴一直以科学计算工作者自诩，希望能做点什么。跟有相关工作经验的同学交流得知，目前的模型都是基于连续的常微分方程，该种模型理论性强、模型简单，但不够直观。于是老裴设计了一个离散的微观模型，通过实例表明，所有人应果断拥护和执行ZF的相关决定，
1、毫不犹豫实行自我隔离
2、发现不适症状及时就医
由于老裴水平有限，不当之处请多指教。

1、模型

所提出模型包括两个部分，分别为个体病毒感染概率计算模型和个体移动模型。

1.1 个体病毒感染概率计算模型

用平面坐标系上的点P(xi,yi)描述个体，假设任意两个体P(x_i,y_i)和Q(x_j,y_j)之间的距离为d_ij，如果Q是病毒携带者，Q向P传播病毒并使P感染的几率为p_ij，二者的关系用式（1）描述：

其中:
d_0i为个体P(x_i,y_i)的临界距离；
w_ij反映个体P和Q的亲密度，其值取0~1；
p_0i为w_ij=0时P在临界距离d_0i内被感染的几率，也称基础感染几率；
α_i为大于1的系数。
个体的临界距离d_0i大小以及其在临界距离内被感染的几率r_0i反映了个体自身的抵抗能力，w_ij反映个体之间的亲密度指数，α_i反映影响范围。
需要注意的是，此处所指的距离是广义距离，并不是指现实中两个个体的实际距离，可以认为是两个体多方面的综合距离，包括个体之间的时空关系、亲密度关系、抵抗能力等，同时p_ij (d_ij)≠p_ji (d_ij)。
以上公式反应了个体距离与被感染病毒几率之间的关系，如图1所示为具体实例（p_0i=0.05，d_0i=5，α_i=2，w_ij=0和0.5）。

该模型假设个体在临界范围d_0i内的感染几率是一定的，超过该范围后感染几率将逐渐减小，达到一定值后，感染几率迅速减小为0。个体的亲密度越高，感染几率越大。

1.2 个体移动模型

假设任意个体P(x_i,y_i)在平面内随机游走（Random Walk），用个体的移动步长d_step来描述个体的活跃度，值越大表明越活跃。随机游走相关知识这里不再赘述，请读者自行查阅相关资料。
初始时刻，假设每个个体之间的距离为固定值d_f，且满足p(d_f )≈0。如图2所示。

1.3 求解过程

整个求解过程如下：

Step 1：确定模型参数，主要包括：

个体数量NI；
个体初始化距离d_f；
个体临界距离矩阵
个体在临界距离内被携带者感染的0亲密度概率矩阵，也称为基础感染几率矩阵
每个个体与其他个体的亲密度矩阵
个体的概率计算模型底数矩阵
随机游走时间步长t_step；
和移动距离步长矩阵d_step；
个体游走持续时间t₁.

Step 2：初始化个体位置；
Step 3：放置传染源个体；

Step 4：在每个时间步长t_step内更新个体位置并根据感染概率模型更新被感染的个体，直到达到时间t₁。

2 实例分析

接下来将逐一模拟多种不同的情况。需要注意的是，由于模型中很多参数的随机性，如果要达到较好的结果，可能需要不断的调整参数，做一个快乐的调参侠。

2.1 正常活动传播模拟

该种情形模拟的是病毒发展初期，已有人感染了病毒，但未引起重视，包括感染者在内，所有人正常活动。通过设置不同的参数观察传播过程。

基本参数
N = 20
NI = N*N  #个体数量为400
D0 = 2.   #临界距离
P0 = 0.3  #个体基础感染几率
W = 0    #亲密度
df = 6.   #初始距离
t_step = 0.01  #时间步长
d_step = 0.05 #移动步长
m = Diffusion(N)
m.init_pos(df)   #初始化位置
m.init_para(D0,P0,W,t_step = t_step,d_step = d_step)   #模型参数设置
m.place_source(153)  #放置感染源个体
for t in[1,2,3,4,5,6,7]: m.walk(t)m.plot()

接下来的分析是在以上参数的基础上进行修改。
以上模型个体数量为400，每个个体的临界距离d₀相同，基础感染几率p₀=0.3也相同，所有个体之间亲密度w=0，每个个体的移动步长也相等，默认感染概率模型底数为e。
分别计算移动步长d_step = 0.05和d_step = 0.1，基础感染几率p₀=0.3，投放1个病源时间增量t=[1,2,3,4,5,6,7]各节点上感染过程图如下。
其中：
红点代表未感染个体；
蓝色星星为被感染个体。
绿点为隔离个体。

接下来，将基础感染几率降低为p₀=0.05，相同时间节点上的感染过程图如下

将不同参数下各时间点被感染总数回执在下图中

可以看到，时间t=5内感染总数基本一致，随着时间的增加，高基础感染率p₀=0.3和高移动步长d_step = 0.1增长趋势更快。
这也意味着，早期的扩散增长与参数大小无太大的关系，及时发现及时隔离

2.2 隔离模拟

现在将感染源附近的进行自我隔离，及令其活动步长d_step=0（绿点），基础感染率仍设为p₀=0.3，其他个体的基础感染率仍为p₀=0.3。移动步长d_step = 0.1。感染过程图如下

可以看到，自我隔离后的个体被感染的几率明显降低。

如果将隔离者的基础感染几率设置为p₀=0.0，将不会感染，而且将会将低传播速度，因为隔离减少了传播路径。
这也意味着，自我隔离是降低传播的重要途径

2.3 频繁活动模拟

现在提高部分个体的活跃度，将1/3个体的活动步长设置为d_step=0.5，感染过程图如下

可以看到，随着个体活跃度的增加，感染速度增加。但隔离的个体依然不会受到感染。

3 结论和展望

3.1 结论

结论和现在的实际情况完全一致：

1、自我隔离
2、如若有感染征兆，及时就医

3.2 展望

目前的分析比较简单，调整参数很花时间，但时间有限。有可能希望后期能跟进一些工作：
1、更多感染几率模型的考虑，包括临界距离和模型底数的调整
2、超级感染者的传播路径
3、社区传播模拟（个体间设置亲密度系数）

4 特别说明

1、因为是概率模型，每次计算结果都不尽相同，但差距在可控范围内。时间仓促，老裴没有时间调出最优参数，发现最符合实际的规律。还可以做更多的分析工作，但老裴现在实在没有心思。
2、允许任何人在此模型基础上开展完善和改进，调出最有参数，但必须有引用说明
3、该列将收录至老裴新书《Python科学计算基础》一书

5 计算程序Python代码（参数需要自行调）

#导入所需的库
import numpy as np
import matplotlib.pyplot as plt#创建扩散模型类
class Diffusion:#初始化def __init__(self,N):self.N = Nself.NI = N**2self.indivduals = Noneself.df = Noneself.NoP = 0#初始化个体位置def init_pos(self,df):self.df = dfN = self.Nx,y = np.linspace(0,(N-1)*df,N),\np.linspace(0,(N-1)*df,N)vx,vy = np.meshgrid(x,y)self.ix,self.iy = vx.flatten(),vy.flatten()self.iix,self.iiy = vx.flatten(),vy.flatten()#设置模型参数def init_para(self,D0,P0,W,d_step = 0.05,t_step = 0.01,A=np.e):self.D0 = np.asarray(D0)   self.P0 = np.asarray(P0)   self.W = np.asarray(W)   self.A = np.asarray(A)     d0 = np.max(self.D0)p0 = np.min(self.P0)w = np.min(self.W)a = np.min(self.A)self.t_step = t_stepself.d_step = d_stepdf = self.dfself.p_min = (1+w)*p0/(a**(df-d0))self.patients = np.zeros(self.NI,dtype = int)#放置感染源def place_source(self,source):self.patients[source] = 1#个体游走def walk(self,t):t_step = self.t_stepd_step = self.d_stepNI = self.NIstep_d = d_step*self.dfn_step = int(t//t_step)moves = np.random.random_integers(1,4,n_step*NI).reshape(n_step,NI)ix,iy = self.ix,self.iyiix,iiy = self.iix,self.iiyUP,DOWN,LEFT,RIGHT = 1,2,3,4for step in range(n_step):this = moves[step]ix += np.where(this == RIGHT,d_step,0.)ix -= np.where(this == LEFT,d_step,0.)iy += np.where(this == UP,d_step,0.)iy -= np.where(this == DOWN,d_step,0.)if np.any(self.patients):self.update_patients()self.NoP = np.nonzero(self.patients)[0].shape[0]#更新感染者def update_patients(self):df = self.dfcurrent_patients = np.nonzero(self.patients)[0]for patient in current_patients:Dx = self.ix[patient] - self.ixDy = self.iy[patient] - self.iyDxy = np.hypot(Dx,Dy)Dxy = np.where(Dxy>df,df,Dxy)Dij = Dxy-self.D0WP0 = (1+self.W)*self.P0WP1 = (1+self.W)*self.P0/(self.A**(Dij))Pij0 = np.where(Dij<=0,WP0,0.)Pij1 = np.where(Dij>0,WP1,0.)Pij = Pij0 + Pij1Pr = np.random.uniform(self.p_min,1,self.NI)DP = Pr - Pijinfected = np.argwhere(DP<0)self.patients[infected] = 1def plot(self):ix,iy = self.ix,self.iypatients = self.patientsx_max,x_min = np.max(ix),np.min(ix)y_max,y_min = np.max(iy),np.min(iy)x_max,x_min = x_max+10,x_min-10y_max,y_min = y_max+10,y_min-10infected = np.where(patients == 1.)noninfected = np.where(patients != 1.)iix = ix[infected]iiy = iy[infected]nix = ix[noninfected]niy = iy[noninfected]plt.axis([x_min,x_max,y_min,y_max])plt.scatter(iix,iiy,marker = "*",color = "b")plt.scatter(nix,niy,marker = ".",color = "r")plt.show()if __name__ == "__main__":N = 20NI = N*ND0 = 2.P0 = 0.3W = 0df = 6.t_step = 0.01d_step = 0.05m = Diffusion(N)m.init_pos(df)m.init_para(D0,P0,W,t_step = t_step,d_step = d_step)m.place_source(153)numbers = []m.plot()for t in [1,2,3,4,5,6,7]:m.walk(t)m.plot()numbers.append(m.NoP)