14 分析
在上一章节，我们使用仿真来预测易感人群中的传染病的影响，并设计干预措施来让传染病的影响最小化。
在这一章节，我们将使用分析来研究参数 和
与仿真结果的关系。
14.1 无量纲化
13.3节中的图表明SIR模型中的参数 和
之间有一个关系确定了仿真的结果，即感染学生的比例。让我们想想这种关系可能是什么。
- 当
超过
时，这意味着日接触率（或其他时间单位）多于日治愈率。
- 作为另一种选择，我们可以考虑比值
，即每个康复病例所对应的接触数量。因为分子和分母的单位相同，所以这个比值是无量纲的，即没有单位。
在建模和仿真项目中，使用无量纲的单位描述现实系统通常是一个实用的技巧。事实上，它非常有用以至于它有一个名字：无量纲化 (http://modsimpy.com/nondim)。
所以我们先尝试第二个选择，在本章的 notebook 中你可以将第一个选项作为练习来探索。
14.2 探索结果
假设我们有一个扫描框（SweepFrame），其中每个 值对应一行，每个
值对应一列。 扫描框中的每一个元素都是在仿真中使用给定的参数对得到的学生感染的比例。
我们可以像这样输出扫描框中的值：
for gamma in frams.columns:
column = frame[gamma]
for beta in column.index:
frac_infected = column[beta]
print(beta, gamma, frac_infected)
这是我们看到的第一个，一个 for 循环嵌套在另一个for循环中的例子：
- 每次外部循环运行的时候，它都从数据框（DataFrame）的列中选择一个
- 每次内部的循环运行的时候，它都从外部循环选出的列中选出一个
在上一个章节的例子中，框（frame）有4列，每列为 的一个取值，有11行，每行对应
的一个取值。因此，这些循环输出44行，每个输出对应一个参数对。
下面的函数封装了上述的循环，并绘制了受感染的分数与比率 的关系图：
图 14.1: 总感染比例关于接触数的函数
def plot_sweep_frame(frame):
for gamma in frame.columns:
series = frame[gamma]
for beta in series.index:
frac_infected = series[beta]
plot(beta/gamma, frac_infected, 'ro')
图14.1显示结果至少大致落在一条曲线上。这意味着我们可以只用一个参数，即 来预测将会被感染的学生的比例。我们不再需要分别知道
和
的值。
14.3 接触数
从11.2节回忆起，在给定的一天中，新增的感染数量为 ，康复的数量为
。如果我们将这些数相除，结果为
，即每个治愈新增的感染数（占总人口的比例）。
当一种新的疾病被易感人群感染后， 的值接近于1，因此每个患者感染的人数为
。该比率被称为“接触数”或“基本传染数”（请查看http://modsimpy.com/contact）。按照惯例，它通常表示为
，但是在SIR模型的语境中，这种表示方式令人混乱，因此我们将改用
。
上一节中的结果表明， 与感染总数之间存在关联。 我们可以通过分析11.3节中的微分方程来推导这种关系：
我们使用与之前将感染率除以治愈率以得到无量纲的 的相同的方式，将
除以
来得到一个比率的比例：
将一个微分方程除以另一个微分方程并不是一个易理解的做法，但是在这个例子中这个做法十分有用，因为这个它使我们得到了一个无关于时间的 、
和
之间的关系。理论上，这个等式让通过观测真实疫情推导出
成为可能。
这里是推导步骤。我们将等式两边同时乘以 ：
然后两边同时积分：
其中 是积分里的常数。将等式重新整理：
现在让我们看看是否可以确定 是什么。在流行开始时，如果感染的比例很小，并且每个人都易感，那么我们可以使用近似值
和
来计算
。
由于 ，我们得到
。
现在，在流行结束时，我们假设 ，
是一个未知量，简写为
。现在我们知道：
求解 ，我们得到：
通过将 和
关联起来，这个方程式可以通过数据来估计
，并可能预测未来流行病的趋势。
14.4 分析并仿真
让我们将分析结果与仿真结果进行比较。我将为 的创建一个值的数组
s_inf_array = linspace(0.0001, 0.9999, 31)
并计算 的对应值：
c_array = log(s_inf_array) / (s_inf_array - 1)
为了得到总体感染数，我们计算 和
的差距，并将它存储到序列（Series）中:
frac_infected = 1 - s_inf_array
frac_infected_series = Series(frac_infected, index=c_array)
图 14.2： 总感染分数关于接触数的函数展现了模拟和分析的结果
从10.2节回想，一个序列（Series）对象包含一个索引和一个相应的值序列。 在这种情况下，索引是c_array，值来自frac_infected。
现在我们可以绘制结果：
plot(frac_infected_series)
图14.2将本章的分析结果与14.1节的仿真结果进行了比较。当接触数大于1时，分析和仿真结果一致。当接触数小于1时，它们不一致：分析认为没有感染，但仿真中会出现一小部分感染。
出现差异的原因是，仿真将时间看作离散的时间序列，而分析则将时间视为连续的量。 换句话说，这两种方法实际上基于不同的模型。 那么哪个模型更好呢？
可能两者都不好。当接触数较小的时候，传染病的早期传播依赖于环境的详细信息。如果我们幸运，最初的患者，即“零号病人”，没有传染任何人，传染病也没有流行起来。如果我们不幸，零号病人有大量密友，或者可能在一家餐厅工作（并且未遵守安全食品处理程序）。
对于接近或小于1的接触数，我们可能需要更精细的模型。但是对于更高的接触数，SIR模型可能就足够了。
14.5 估计接触数
图14.2显示，如果我们知道接触数，我们就可以计算感染的比例。但我们也可以用另一种方式来解读图像，在一场流行病结束时，如果我们能估计出曾经感染过的人口比例，我们就可以用它来估计接触数。
好吧，在理论上我们是可以的。但在实际上，由于曲线的形状，估计结果可能并不好。当接触数接近2时，曲线相当陡峭，这意味着 的微小变化会导致感染数量的巨大变化。如果我们观察到的总体感染的分数在20%到80%之间的任何位置，我们都只能得出
接近2的结论。
另一方面，对于较大的接触数，几乎整个人口都受到了感染，因此曲线几乎是平坦的。 在这种情况下，我们将无法精确估计 ，因为事实上任何大于3的值都会产生相同的结果。幸运的是，这在现实世界中不太可能发生。 几乎没有什么流行病影响到接近90％的人口。
因此SIR模型具有局限性。 但是它可以洞察传染病的传播进程，尤其是群体免疫现象。 正如我们在第12章中看到的那样，如果我们知道模型的参数，则可以使用它来评估可能的干预措施。 正如我们在本章中看到的那样，我们也许可以使用早期爆发的数据来估计参数。
在继续之前，您可能需要阅读本章的 notebook ，chap14.ipynb 并进行练习。有关下载和运行代码的说明，请参见0.4节。
本书的中文翻译由南开大学医学院智能医学工程专业2018级、2019级的师生完成，方便后续学生学习《Python仿真建模》课程。翻译人员（排名不分前后）：薛淏源、金钰、张雯、张莹睿、赵子雨、李翀、慕振墺、许靖云、李文硕、尹瀛寰、沈纪辰、迪力木拉、樊旭波、商嘉文、赵旭、连煦、杨永新、樊一诺、刘志鑫、彭子豪、马碧婷、吴晓玲、常智星、陈俊帆、高胜寒、韩志恒、刘天翔、张艺潇、刘畅。
整理校订由刘畅完成，如果您发现有翻译不当或者错误，请邮件联系changliu@nankai.edu.cn。
本书中文版不用于商业用途，供大家自由使用。
未经允许，请勿转载。
