2019-nCoV 疫情传播模拟

武汉必胜！中国必胜！

2020 新年之际，新型冠状病毒疫情在全国各地陆续出现，尤其武汉人民深受其害。为了支援疫区，同时打发被“软禁”在家的时光，决定研究一下目前的疫情，看看是否能发现一些有用的结论。

在网上浏览一下，对于疫情数据的可视化以及建模分析，目前已经做了大量的工作，并且已经建立实时监控疫情的系统。于是准备从另一个方面入手，从本质上模拟疫情的传播，并对比理想状态下的传播于实际情况的区别。

下面是模拟思路和结论分析。

理想状态下建模

首先在理想状态下，我们设定一个“游戏规则”。从第一个感染了病毒的人开始，在随机的 N 个时间点传染给随机的 N 个人。每个感染的人经过随机的 M 日的潜伏期后病发并隔离治疗，不再传染给别人（有点绕，其实就是对现实情况做一个简化）。

根据这个规则，可以分析出病人对象所包含的属性：

class patient:
    '''
    ----------Attributions----------
    id:             [int]按照感染顺序确定的一个递增的id；
    period:         [int]所在的迭代周期，以天为单位，第一个感染者出现为周期 1；
    infect_num:     [int]该病人传染的人数，是一个服从正态分布的随机数；
    attack_interval:[int]潜伏期天数，是一个服从正态分布的随机数；
    infect_period:  [list]传染给别人的时间，潜伏期内均匀分布；
    if_attack:      [bool]是否已经发病；
    '''

    def __init__(self, id, i):
        self.id = id
        self.period = i
        self.infect_num = int(np.abs(np.random.normal(3,2,1)[0]))
        self.attack_interval = int(np.abs(np.random.normal(7,2,1)[0]))
        self.infect_period = list(self.get_infect_period())
        self.if_attack = False
    
    def get_infect_period(self):
        if self.attack_interval:
            return np.random.randint(0, self.attack_interval, self.infect_num)
        else:
            return np.array([0])

接下来模拟传播的过程，为了提高可读性，这里直接按照实际逻辑来写：

def iteration(id, i, patients):
     '''
     单次迭代，模拟一天之中发生的事情
     '''

     new_patients = []
     # 迭代当前患者
     for p in patients:
        # 发病隔离
        if p.if_attack:
            continue
        # 传染新患者
        else:
            if i-p.period >= p.attack_interval - 1:
                p.if_attack = True
            n = sum(np.array(p.infect_period) == (i-p.period))
            for _ in range(n):
                id += 1
                new_p = patient(id, i)
                new_patients.append(new_p)

     # 返回参数
     for new_p in new_patients:
         patients.append(new_p)

     return id, patients

最后进行迭代，模拟疫情的发展趋势：

if __name__ == "__main__":
    id = 1
    p = patient(id,1)   # 第一个感染者
    patients = [p]      # 感染者列表
    numbers = [1]       # 累计感染人数
    attack = [0]        # 累计发病人数

    for i in range(1, 30):
        # 迭代
        id, patients = iteration(id, i, patients)
        # 记录参数
        a = sum([p.if_attack for p in patients])
        attack.append(a)
        numbers.append(len(patients))

模拟疫情趋势

如下图所示，由于上述的简化模型中包含随机的参数，所以每次模拟的结果不尽相同。可以看到，在不采取其他措施的情况下（除了发病隔离），疫情将呈指数式爆发。对结果影响较大的是前几个迭代周期的病人对象，如果最先出现的感染者传染给了较多的人，那么疫情将以更快的速度爆发（如下面的中图）。相反如果最初的病人被及时隔离，那么疫情可能被扼制在摇篮中，没有机会传播（如下面的右图）。

另外从图中可以发现，在统计意义上，感染人数和发病人数具有一个固定的倍数关系（在本例的假设下，这个比值约为 4:1）。这是由于感染病毒后不会立即发病，存在一个较长的潜伏期。这就意味着，在备受关注的确诊病例数据后面，还隐藏这几倍甚至更多的已经感染的人群。所以要抑制疫情，必须控制好这群潜伏的感染者，早发现早隔离。目前社会上呼吁避免外出、避免接触也是出于这点考虑。

上面找到了一个关键问题，就是需要控制带病潜伏的人群。下面考虑一个理想的状态，如果所有人都居家隔离，断绝一切接触，疫情会如何发展？

修改程序隔绝一切传播，得到的模拟结果如下：

可以看到，疫情还会经历一个快速发展的过程，直到增速减缓最后停止增加。也就是说，当确诊病例的增速开始明显减小的时候，就可以期望疫情即将消散。

总结

本文得出两个（非严格的）结论

• 存在大量潜伏期感染者，减少外出接触对疫情控制有极其重要的意义。
• 当确诊病例增速减小，疫情将趋于结束。 参考疫情实时动态

本例存在的缺陷

• 理想状态下，忽略大量实际因素。例如复工影响，治愈率变化，等。

• 本例假设的传染率和潜伏期的概率分布分别为正态分布 N(3, 4) 和 N(7, 4)：

• 
  使用更准确的统计数据进行模拟将得出更贴近实际的结果。 参考 Early Transmission Dynamics in Wuhan, China, of Novel Coronavirus–Infected Pneumonia

• ……