流行病学

警告

pyro.contrib.epidemiology 中的代码正在开发中。 此代码不保证保持向后兼容性。

pyro.contrib.epidemiology 提供了一种用于一类随机离散时间离散计数分区模型的建模语言。该模块实现了黑盒推断（包括随机变分推断和哈密顿蒙特卡洛）、潜在变量的预测以及未来轨迹的预测。

例如，使用情况请参见以下教程：

基础隔室模型

class CompartmentalModel(compartments, duration, population, *, approximate=())[source]

基础类：abc.ABC

离散时间离散值随机分区模型的抽象基类。

派生类必须实现方法 initialize() 和 transition()。派生类可以选择性地实现 global_model()、compute_flows() 和 heuristic()。

示例用法：

# First implement a concrete derived class.
class MyModel(CompartmentalModel):
    def __init__(self, ...): ...
    def global_model(self): ...
    def initialize(self, params): ...
    def transition(self, params, state, t): ...

# Run inference to fit the model to data.
model = MyModel(...)
model.fit_svi(num_samples=100)  # or .fit_mcmc(...)
R0 = model.samples["R0"]  # An example parameter.
print("R0 = {:0.3g} ± {:0.3g}".format(R0.mean(), R0.std()))

# Predict latent variables.
samples = model.predict()

# Forecast forward.
samples = model.predict(forecast=30)

# You can assess future interventions (applied after ``duration``) by
# storing them as attributes that are read by your derived methods.
model.my_intervention = False
samples1 = model.predict(forecast=30)
model.my_intervention = True
samples2 = model.predict(forecast=30)
effect = samples2["my_result"].mean() - samples1["my_result"].mean()
print("average effect = {:0.3g}".format(effect))

一个示例工作流程是在寻找良好的模型结构和先验时使用更便宜的近似推理，然后在模型合理时转向更准确但更昂贵的推理。

  1. 从.fit_svi(guide_rank=0, num_steps=2000)开始，以便在寻找好模型时进行低成本推理。

  2. 此外，通过使用.fit_svi(guide_rank=None, num_steps=5000)移动到低秩多元正态指南来推断长程相关性。

  3. 可选地，通过转向更昂贵（但仍然通过矩匹配近似）的.fit_mcmc(num_quant_bins=1, num_samples=10000, num_chains=2)来推断非高斯后验。

  4. 可选地通过使用更昂贵的基于枚举的算法 .fit_mcmc(num_quant_bins=4, num_samples=10000, num_chains=2) 来改善小计数的拟合（推荐使用GPU）。

Variables

样本 (dict) – 后验样本的字典。

Parameters

  • compartments (list) – 一个包含隔间名称的字符串列表。

  • duration (int) – 此模型中的离散时间步数。

  • population (int 或 torch.Tensor) – 单区域模型中的总人口或区域模型中每个区域的人口张量。

  • 近似 (元组) – 在transition()中应提供点近似值的隔间名称，例如，如果您指定approximate=("I")，则state["I_approx"]将是state["I"]的连续值非枚举点估计。近似值有助于降低计算成本。近似值是连续值的，支持范围为(-0.5, population + 0.5)。

property time_plate

用于时间维度的pyro.plate。

property region_plate

根据此模型是否为区域性的.is_regional，可能是pyro.plate或一个简单的ExitStack。

property full_mass

一个包含全局随机变量名称的单个元组的列表。

property series

每次时间步长采样的样本站点名称的冻结集合。

global_model()[source]

采样并返回任何全局参数。

Returns

一个任意的参数对象（例如 None 或一个元组）。

abstract initialize(params)[source]

返回每个隔间的初始计数。

Parameters

params – 由global_model()返回的全局参数。

Returns

一个将隔间名称映射到初始值的字典。

Return type

dict

abstract transition(params, state, t)[来源]

动态的前向生成过程。

这输入一个当前的state并随机地就地更新该状态。

请注意，此方法在多种不同的解释下被调用，包括批处理和向量化解释。在generate()期间，此方法被调用以生成单个样本。在heuristic()期间，此方法被调用以生成一批样本用于SMC。在fit_mcmc()期间，此方法以向量化形式（在时间上向量化）和顺序形式（针对单个时间步）被调用；两种形式都枚举离散的潜在变量。在predict()期间，此方法被调用以预测一批样本，条件是时间间隔[0:duration]的后验样本。

Parameters

  • params – 由global_model()返回的全局参数。

  • state (dict) – 一个将隔间名称映射到当前张量值的字典。这应该就地更新。

  • t (int 或 slice) – 一个时间索引。在推理过程中，t 可能是一个切片（用于向量化推理）或一个整数时间索引。在预测过程中，t 将是整数时间索引。

finalize(params, prev, curr)[source]

对于依赖于整个时间序列的似然函数的可选方法。

这应该仅用于跨时间耦合状态的不可分解似然。可分解的似然应该添加到transition()方法中，从而使其能够在heuristic()初始化中使用。由于此方法仅在最后一个时间步之后调用，因此它不用于heuristic()初始化。

警告

目前不支持潜在变量。

Parameters

  • params – 由global_model()返回的全局参数。

  • prev (dict) –

  • curr (dict) – 字典将隔间名称映射到整个时间序列的张量。这两个参数偏移了1步，从而使得计算通量的时间序列变得容易。对于量化推断，这使用了近似的点估计，因此用户必须在__init__()中请求任何需要的时间序列，例如，如果似然依赖于I和E时间序列，可以通过调用super().__init__(..., approximate=("I", "E"))来实现。

compute_flows(prev, curr, t)[source]

计算在时间步长 t 前后，给定隔间人口之间的流动。

默认实现假设顺序流终止于名为“R”的隐式隔间。例如，如果：

compartment_names = ("S", "E", "I")

默认实现在时间步长 t = 9 时计算：

flows["S2E_9"] = prev["S"] - curr["S"]
flows["E2I_9"] = prev["E"] - curr["E"] + flows["S2E_9"]
flows["I2R_9"] = prev["I"] - curr["I"] + flows["E2I_9"]

对于更复杂的流程（非顺序、分支、循环、重复等），用户可以重写此方法。

Parameters

  • state (dict) – 一个将隔间名称映射到当前张量值的字典。这应该就地更新。

  • t (int 或 slice) – 一个时间索引。在推理过程中，t 可能是一个切片（用于向量化推理）或一个整数时间索引。在预测过程中，t 将是整数时间索引。

Returns

一个将流程名称映射到张量值的字典。

Return type

dict

generate(fixed={})[source]

从先验生成数据。

Pram dict fixed

一个参数字典，用于条件化。 这些必须是顶层的无父节点，即没有 上游的随机依赖。

Returns

一个字典，将样本站点名称映射到采样值。

Return type

dict

fit_svi(*, num_samples=100, num_steps=2000, num_particles=32, learning_rate=0.1, learning_rate_decay=0.01, betas=(0.8, 0.99), haar=True, init_scale=0.01, guide_rank=0, jit=False, log_every=200, **options)[source]

运行随机变分推断以生成后验样本。

这将运行SVI，并在完成后设置.samples属性。

这种近似推理方法对于快速迭代概率模型非常有用。

Parameters

  • num_samples (int) – 从训练好的指南中抽取的后验样本数量。默认为100。

  • num_steps (int) – SVI 的步数。

  • num_particles (int) – 每步的SVI粒子数量。

  • learning_rate (int) – ClippedAdam 优化器的学习率。

  • learning_rate_decay (int) – ClippedAdam 优化器的学习率。请注意这是整个计划中的衰减，而不是每步的衰减。

  • betas (tuple) – ClippedAdam 优化器的动量参数。

  • haar (bool) – 是否使用Haar小波重参数化器。

  • guide_rank (int) – 自动正态引导的秩。如果为零（默认值），则使用AutoNormal引导。如果为正整数或None，则使用AutoLowRankMultivariateNormal引导。如果为字符串“full”，则使用AutoMultivariateNormal引导。后两者需要更多的num_steps来拟合。

  • init_scale (float) – AutoLowRankMultivariateNormal 引导的初始比例。

  • jit (bool) – 是否使用jit编译的ELBO。

  • log_every (int) – 记录svi损失的频率。

  • heuristic_num_particles (int) – 传递给 heuristic() 作为 num_particles。默认值为1024。

Returns

SVI损失的时序（有助于诊断收敛性）。

Return type

list

fit_mcmc(**options)[source]

运行NUTS推断以生成后验样本。

这使用NUTS内核来运行 MCMC，并在完成后设置.samples 属性。

当num_quant_bins > 1时，使用基于渐近精确枚举的模型，而当num_quant_bins == 1时，使用更经济的矩匹配近似模型。

Parameters

  • **options – 传递给 MCMC 的选项。其余的选项被提取出来并具有特殊含义。

  • num_samples (int) – 通过mcmc绘制的后验样本数量。默认为100。

  • max_tree_depth (int) – (默认值 5). NUTS 核的最大树深度。

  • full_mass – NUTS核的质量矩阵的规范。默认为全局随机变量的完整质量。

  • arrowhead_mass (bool) – 是否将 full_mass 视为箭头矩阵的头部，而不是简单地视为一个块。默认为 False。

  • num_quant_bins (int) – 如果大于1，则通过本地枚举使用渐近精确的推理方法，枚举的量化箱数量为此值。 如果等于1，则使用连续值的松弛近似推理。 请注意，计算成本在num_quant_bins中是指数级的。 默认值为1，用于松弛推理。

  • haar (bool) – 是否使用Haar小波重参数化器。 默认为True。

  • haar_full_mass (int) – 包含在全质量矩阵中的低频Haar分量的数量。如果 haar=False 则忽略此参数。默认值为10。

  • heuristic_num_particles (int) – 传递给 heuristic() 作为 num_particles。默认值为1024。

Returns

用于诊断的MCMC对象，例如 MCMC.summary()。

Return type

MCMC

predict(forecast=0)[来源]

预测潜在变量并可选地进行向前预测。

这只能在fit_mcmc()之后运行，并且绘制与传递给fit_mcmc()相同的num_samples。

Parameters

forecast (int) – 向前预测的时间步数。

Returns

一个字典，将样本站点名称（或隔间名称）映射到一个张量，该张量的第一维度对应于样本批处理。

Return type

dict

heuristic(num_particles=1024, ess_threshold=0.5, retries=10)[source]

找到所有潜在变量的初始可行猜测，与观察到的数据一致。这是必要的，因为并非所有假设都是可行的，HMC需要从一个可行的解决方案开始才能继续。

默认实现尝试使用SMCFilter从前提出发寻找一个可行的状态。然而，在SMC表现不佳的情况下，例如在高维模型中，可以覆盖此方法。

Parameters

  • num_particles (int) – 用于SMC的粒子数量。

  • ess_threshold (float) – SMC的有效样本量阈值。

Returns

一个将样本站点名称映射到张量值的字典。

Return type

dict

示例模型

简单的SIR模型

class SimpleSIRModel(population, recovery_time, data)[source]

易感-感染-恢复模型。

要自定义此模型，我们建议分叉并编辑此类。

这是一个随机的离散时间离散状态模型，包含三个部分：“S”代表易感者，“I”代表感染者，“R”代表康复者（康复者是隐含的：R = population - S - I），其转换过程为S -> I -> R。

Parameters

  • 人口 (int) – 总人口 population = S + I + R.

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> I转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

简单的SEIR模型

class SimpleSEIRModel(population, incubation_time, recovery_time, data)[source]

易感-暴露-感染-恢复模型。

要自定义此模型，我们建议分叉并编辑此类。

这是一个具有四个隔室的随机离散时间离散状态模型：“S”表示易感者，“E”表示暴露者，“I”表示感染者，“R”表示康复者（康复者是隐含的：R = population - S - E - I），其转换过程为S -> E -> I -> R。

Parameters

  • 人口 (int) – 总 population = S + E + I + R.

  • incubation_time (float) – 平均潜伏时间（在状态 E中的持续时间）。必须大于1。

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> E转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

简单的SEIRD

class SimpleSEIRDModel(population, incubation_time, recovery_time, mortality_rate, data)[source]

易感-暴露-感染-康复-死亡模型。

要自定义此模型，我们建议分叉并编辑此类。

这是一个具有四个隔室的随机离散时间离散状态模型：“S”表示易感者，“E”表示暴露者，“I”表示感染者，“D”表示死亡个体，“R”表示康复个体（康复个体是隐含的：R = population - S - E - I - D），具有转换S -> E -> I -> R和I -> D。

由于过渡不是简单的线性连续，此模型实现了一个自定义的 compute_flows() 方法。

Parameters

  • 人口 (int) – 总人口 population = S + E + I + R + D.

  • incubation_time (float) – 平均潜伏时间（在状态 E中的持续时间）。必须大于1。

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • mortality_rate (float) – 导致死亡的感染比例。 必须在开区间 (0, 1) 内。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> E转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

过度分散的SIR

class OverdispersedSIRModel(population, recovery_time, data)[来源]

将SimpleSIRModel推广到具有过度分散分布的情况。

要自定义此模型，我们建议分叉并编辑此类。

这增加了一个单一的全局过度分散参数，用于控制转移和观察分布的过度分散。有关分布详情，请参见 binomial_dist() 和 beta_binomial_dist()。有关先前结合过度分散分布的工作，请参见 [1,2,3,4]。

参考文献：

[1] D. Champredon, M. Li, B. Bolker. J. Dushoff (2018)

“预测埃博拉合成流行病的两种方法” https://www.sciencedirect.com/science/article/pii/S1755436517300233

[2] Carrie Reed et al. (2015)

“基于美国人口监测数据估计流感疾病负担” https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4349859/

[3] A. Leonard, D. Weissman, B. Greenbaum, E. Ghedin, K. Koelle (2017)

“从病原体深度测序数据估计传播瓶颈大小，应用于人类甲型流感病毒” https://jvi.asm.org/content/jvi/91/14/e00171-17.full.pdf

[4] A. Miller, N. Foti, J. Lewnard, N. Jewell, C. Guestrin, E. Fox (2020)

“移动趋势提供了SARS-CoV-2传播变化的领先指标” https://www.medrxiv.org/content/medrxiv/early/2020/05/11/2020.05.07.20094441.full.pdf

Parameters

  • 人口 (int) – 总人口 population = S + I + R.

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> I转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

过度分散的SEIR

class OverdispersedSEIRModel(population, incubation_time, recovery_time, data)[source]

将SimpleSEIRModel推广到具有过度分散分布的情况。

要自定义此模型，我们建议分叉并编辑此类。

这增加了一个单一的全局过度分散参数，用于控制转移和观察分布的过度分散。有关分布详情，请参见 binomial_dist() 和 beta_binomial_dist()。有关先前结合过度分散分布的工作，请参见 [1,2,3,4]。

参考文献：

[1] D. Champredon, M. Li, B. Bolker. J. Dushoff (2018)

“预测埃博拉合成流行病的两种方法” https://www.sciencedirect.com/science/article/pii/S1755436517300233

[2] Carrie Reed et al. (2015)

“基于美国人口监测数据估计流感疾病负担” https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4349859/

[3] A. Leonard, D. Weissman, B. Greenbaum, E. Ghedin, K. Koelle (2017)

“从病原体深度测序数据估计传播瓶颈大小，应用于人类甲型流感病毒” https://jvi.asm.org/content/jvi/91/14/e00171-17.full.pdf

[4] A. Miller, N. Foti, J. Lewnard, N. Jewell, C. Guestrin, E. Fox (2020)

“移动趋势提供了SARS-CoV-2传播变化的领先指标” https://www.medrxiv.org/content/medrxiv/early/2020/05/11/2020.05.07.20094441.full.pdf

Parameters

  • 人口 (int) – 总 population = S + E + I + R.

  • incubation_time (float) – 平均潜伏时间（在状态 E中的持续时间）。必须大于1。

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> E转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

超级传播SIR

class SuperspreadingSIRModel(population, recovery_time, data)[源代码]

通过添加超级传播效应来推广SimpleSIRModel。

要自定义此模型，我们建议分叉并编辑此类。

该模型通过假设每个感染者感染BetaBinomial多个易感个体来解释超级传播（过度分散的个体繁殖数），其中BetaBinomial分布作为过度分散的二项分布，适应了更标准的负二项分布，后者作为过度分散的泊松分布[1,2]，适用于有限人口的情况。为了保持马尔可夫结构，我们遵循[2]并假设单个个体的所有感染发生在该个体进行I -> R转换的单个时间步长上。也就是说，尽管SimpleSIRModel假设感染者在每个感染时间步长（平均超过tau个步长）期间感染Binomial(S,R/tau)多个易感个体，但该模型假设他们仅在恢复前的最后一个时间步长感染BetaBinomial(k,…,S)多个易感个体。

参考文献

[1] J. O. Lloyd-Smith, S. J. Schreiber, P. E. Kopp, W. M. Getz (2005)

“超级传播及个体差异对疾病出现的影响” https://www.nature.com/articles/nature04153.pdf

[2] Lucy M. Li, Nicholas C. Grassly, Christophe Fraser (2017)

“使用病原体系统发育和发病率时间序列量化传播异质性” https://academic.oup.com/mbe/article/34/11/2982/3952784

Parameters

  • 人口 (int) – 总人口 population = S + I + R.

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> I转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

超级传播SEIR

class SuperspreadingSEIRModel(population, incubation_time, recovery_time, data, *, leaf_times=None, coal_times=None)[来源]

通过添加超级传播效应来推广SimpleSEIRModel。

要自定义此模型，我们建议分叉并编辑此类。

该模型通过假设每个感染者感染BetaBinomial数量的易感个体来解释超级传播（过度分散的个体繁殖数），其中BetaBinomial分布作为过度分散的二项分布，适应了更标准的负二项分布，后者作为过度分散的泊松分布[1,2]，适用于有限人口的情况。为了保持马尔可夫结构，我们遵循[2]并假设单个个体的所有感染发生在该个体进行I -> R转换的单个时间步长上。也就是说，尽管SimpleSEIRModel假设感染者在每个感染时间步长（平均超过tau个步长）期间感染Binomial(S,R/tau)数量的易感个体，但该模型假设他们仅在恢复前的最后一个时间步长感染BetaBinomial(k,…,S)数量的易感个体。

该模型还为观察到的系统发育数据添加了一个可选的似然性，以合并时间的形式提供。这些数据以一对(leaf_times, coal_times)的形式提供，分别表示基因组测序和谱系合并的时间。我们使用CoalescentRateLikelihood将这些数据纳入模型，其中基础合并率是从S和I种群计算得出的。这种似然性在时间上是独立的，并保留了推理所需的马尔可夫性质。

参考文献

[1] J. O. Lloyd-Smith, S. J. Schreiber, P. E. Kopp, W. M. Getz (2005)

“超级传播及个体差异对疾病出现的影响” https://www.nature.com/articles/nature04153.pdf

[2] Lucy M. Li, Nicholas C. Grassly, Christophe Fraser (2017)

“使用病原体系统发育和发病率时间序列量化传播异质性” https://academic.oup.com/mbe/article/34/11/2982/3952784

Parameters

  • 人口 (int) – 总 population = S + E + I + R.

  • incubation_time (float) – 平均潜伏时间（在状态 E中的持续时间）。必须大于1。

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> E转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

异质SIR

class HeterogeneousSIRModel(population, recovery_time, data)[源代码]

通过允许Rt和rho随时间变化，对SimpleSIRModel进行了泛化。

要自定义此模型，我们建议分叉并编辑此类。

在这个模型中，响应率 rho 在三段中是分段常数，其值未知。繁殖数 Rt 是一个常数 R0 与一个因子 beta 的乘积，该因子在对数空间中通过布朗运动漂移。 rho 和 Rt 都可以作为时间序列使用。

Parameters

  • 人口 (int) – 总人口 population = S + I + R.

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> I转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

稀疏SIR

class SparseSIRModel(population, recovery_time, data, mask)[source]

将SimpleSIRModel推广到允许稀疏观察到的感染。

要自定义此模型，我们建议分叉并编辑此类。

该模型允许在不均匀的时间间隔内观察累积感染情况。为了保持马尔可夫结构（从而便于推理），该模型添加了一个辅助隔间O，表示每个时间点的完全观察到的累积观察次数。在观察到的时间点（当mask[t] == True时），O必须与提供的数据完全匹配；在观察到的时间点之间，O随机插值提供的数据。

该模型展示了如何实现自定义的compute_flows()方法。在此模型中需要自定义方法，因为S隔间的居民可以过渡到I和O隔间，从而允许重复。

Parameters

  • 人口 (int) – 总人口 population = S + I + R.

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • data (可迭代对象) – 累计观察感染的时间序列。 每当 mask[t] == True 时，data[t] 对应一个观察值；否则 data[t] 可以是任意的，例如 NAN。

  • mask (iterable) – 布尔时间序列，表示每个时间步是否进行了观察。应满足 len(mask) == len(data)。

未知起始SIR

class UnknownStartSIRModel(population, recovery_time, pre_obs_window, data)[source]

通过允许未知的首次感染日期，对SimpleSIRModel进行了推广。

要自定义此模型，我们建议分叉并编辑此类。

该模型展示了：

  1. 如何整合来自外部来源的自发感染；

  2. 如何通过在transition()中支持预测来整合时变分段rho。

  3. 如何重写predict()方法来计算额外的统计信息。

Parameters

  • 人口 (int) – 总人口 population = S + I + R.

  • recovery_time (float) – 平均恢复时间（在状态 I中的持续时间）。必须大于1。

  • pre_obs_window (int) – 在开始 data 之前可能发生初始感染的时间步数。必须为正数。

  • data (可迭代对象) – 新观察到的感染时间序列。每个时间步长在0和S -> I转换次数之间是二项分布的。这允许假阴性但不允许假阳性。

区域SIR

class RegionalSIRModel(population, coupling, recovery_time, data)[source]

将SimpleSIRModel推广到同时建模多个区域，并在区域之间具有弱耦合。

要自定义此模型，我们建议分叉并编辑此类。

区域通过一个coupling矩阵耦合，矩阵的条目在[0,1]范围内。 全1矩阵等同于一个单一区域。单位矩阵等同于一组独立区域。这不需要是对称的， 但对称矩阵可能更符合物理实际。每个时间步长S2I的新感染预期数量是二项分布的， 其均值为：

E[S2I] = S (1 - (1 - R0 / (population @ coupling)) ** (I @ coupling))
       ≈ R0 S (I @ coupling) / (population @ coupling)  # for small I

因此，在一个几乎完全易感的人群中，一个单一的感染者平均会感染大约R0个新个体，这与coupling无关。

该模型展示了：

  1. 如何使用population向量创建区域模型。

  2. 如何使用self.region_plate来建模同质参数（这里是R0）和具有层次结构的异质参数（这里是rho）。

  3. 如何在transition()中使用state["I_approx"]来近似耦合区域。

Parameters

  • 人口 (torch.Tensor) – 每个区域的人口的张量，定义 population = S + I + R.

  • 耦合 (torch.Tensor) – 成对耦合矩阵。条目应在 [0,1] 范围内。

  • recovery_time (float) – 平均恢复时间（在状态 I 中的持续时间）。 必须大于1。

  • data (可迭代) – 时间 x 区域大小的新观察感染张量。每个时间步是介于0和S -> I转换次数之间的二项分布向量。这允许假阴性但不允许假阳性。

异质区域SIR

class HeterogeneousRegionalSIRModel(population, coupling, recovery_time, data)[来源]

通过允许Rt和rho随时间变化，推广了RegionalSIRModel。

要自定义此模型，我们建议分叉并编辑此类。

在这个模型中，响应率 rho 随时间和地区变化， 而繁殖数 Rt 随时间变化但在各地区之间共享。 这两个参数根据学习的漂移率按照转换后的布朗运动漂移。

该模型展示了如何对层次潜在时间序列进行建模，而不是对分区变量进行建模。

Parameters

  • 人口 (torch.Tensor) – 每个区域的人口张量，定义 population = S + I + R.

  • 耦合 (torch.Tensor) – 成对耦合矩阵。条目应在 [0,1] 范围内。

  • recovery_time (float) – 平均恢复时间（在状态 I 中的持续时间）。 必须大于1。

  • data (可迭代) – 时间 x 区域大小的新观察到的感染张量。每个时间步是介于0和S -> I转换次数之间的二项分布向量。这允许假阴性但不允许假阳性。

分布

set_approx_sample_thresh(thresh)[source]

实验性的上下文管理器/装饰器，用于临时设置全局默认值 Binomial.approx_sample_thresh，从而降低从 Binomial, BetaBinomial, ExtendedBinomial, ExtendedBetaBinomial, 以及由 infection_dist() 返回的分布中进行采样的计算复杂度。

这对于从非常大的total_count中进行采样非常有用。

这是由 CompartmentalModel 内部使用的。

Parameters

thresh (int 或 float.) – 新的临时阈值。

set_approx_log_prob_tol(tol)[源代码]

实验性的上下文管理器/装饰器，用于临时设置全局默认值 Binomial.approx_log_prob_tol 和 BetaBinomial.approx_log_prob_tol，从而降低评分 Binomial 和 BetaBinomial 分布的计算复杂度。

这是由 CompartmentalModel 内部使用的。

Parameters

tol (int 或 float.) – 新的临时 tolold。

binomial_dist(total_count, probs, *, overdispersion=0.0)[源代码]

返回一个Beta-Binomial分布，这是一个根据参数overdispersion进行过度分散的Binomial分布版本，通常设置在0.1到0.5的范围内。

这对于（1）拟合相对于二项分布过度分散的实际数据，以及（2）放宽大群体模型以改进推断非常有用。特别是overdispersion参数在随机模型中设定了相对不确定性的下限，使得随着群体增加，系统趋向于一个具有有限随机性的无标度动力系统，这与基于二项分布的SDEs在大群体极限下收敛到确定性ODEs形成对比。

此参数化满足以下属性：

  1. 方差在overdispersion中单调增加。

  2. overdispersion = 0 导致二项分布。

  3. overdispersion 限制了相对不确定性 std_dev / (total_count * p * q)，其中 probs = p = 1 - q，并且当 total_count → ∞ 时，相对不确定性趋近于一个渐近线。这与二项分布形成对比，二项分布的相对不确定性趋近于零。

  4. 如果 X ~ binomial_dist(n, p, overdispersion=σ) 那么在大 总体极限 n → ∞ 下，缩放后的随机变量 X / n 在分布上收敛到 LogitNormal(log(p/(1-p)), σ)。

为了实现这些属性，我们设置 p = probs, q = 1 - p，并且：

concentration = 1 / (p * q * overdispersion**2) - 1
Parameters

  • total_count (int 或 torch.Tensor) – 伯努利试验的次数。

  • probs (float 或 torch.Tensor) – 事件概率。

  • 过度分散 (float 或 torch.tensor) – 过度分散的量，在半开区间 [0,2) 内。默认为零。

beta_binomial_dist(concentration1, concentration0, total_count, *, overdispersion=0.0)[source]

返回一个Beta-Binomial分布，这是一个根据额外参数overdispersion（通常在0.1到0.5范围内设置）的通常Beta-Binomial分布的过度分散版本。

Parameters

  • concentration1 (float 或 torch.Tensor) – Beta 分布的第一个浓度参数 (alpha)。

  • concentration0 (float 或 torch.Tensor) – Beta 分布的第二个浓度参数 (beta)。

  • total_count (float 或 torch.Tensor) – 伯努利试验的次数。

  • 过度分散 (float 或 torch.tensor) – 过度分散的量，在半开区间 [0,2) 内。默认为零。

infection_dist(*, individual_rate, num_infectious, num_susceptible=inf, population=inf, concentration=inf, overdispersion=0.0)[源代码]

创建一个Distribution来表示在离散时间步长内新感染的数量。

这将返回泊松分布、负二项分布、二项分布或贝塔二项分布，具体取决于population和concentration是否为有限值。在Pyro模型中，population通常是有限的。在极限情况下population → ∞和num_susceptible/population → 1，二项分布收敛到泊松分布，贝塔二项分布收敛到负二项分布。在极限情况下concentration → ∞，负二项分布收敛到泊松分布，贝塔二项分布收敛到二项分布。

过度分散的分布（当concentration < ∞时返回的负二项分布和贝塔二项分布）对于建模超级传播者个体非常有用[1,2]。有限支持的分布（二项分布和负二项分布）在小群体中以及在截断或审查成本较高的概率编程系统中非常有用[3]。

参考文献

[1] J. O. Lloyd-Smith, S. J. Schreiber, P. E. Kopp, W. M. Getz (2005)

“超级传播及个体差异对疾病出现的影响” https://www.nature.com/articles/nature04153.pdf

[2] Lucy M. Li, Nicholas C. Grassly, Christophe Fraser (2017)

“使用病原体系统发育和发病率时间序列量化传播异质性” https://academic.oup.com/mbe/article/34/11/2982/3952784

[3] Lawrence Murray et al. (2018)

“延迟采样和概率程序的自动Rao-Blackwell化” https://arxiv.org/pdf/1708.07787.pdf

Parameters

  • individual_rate – 在大人口限制下，每个感染个体每个时间步的平均感染数，等于 R0 / tau 其中 R0 是基本再生数，tau 是感染的平均持续时间。

  • num_infectious – 当前时间步的感染个体数量，有时为 I，有时为 E+I。

  • num_susceptible – 当前时间步易感个体的数量 S。默认情况下，这是一个无限大的群体。

  • population – 种群中的个体总数。 默认情况下，这是一个无限种群。

  • concentration – 在超级传播者的过度分散模型中的浓度或分散参数 k [1,2]。默认值为最小方差 concentration = ∞。

  • 过度分散 (float 或 torch.tensor) – 过度分散的量，在半开区间 [0,2) 内。默认为零。

class CoalescentRateLikelihood(leaf_times, coal_times, duration, *, validate_args=None)[来源]

基础类：object

实验性 这不是一个Distribution，而是作为CoalescentTimesWithRate的转置版本，使得rate_grid的元素独立，从而与plate和poutine.markov兼容。对于非批处理输入，以下都是等效的似然：

# Version 1.
pyro.sample("coalescent",
            CoalescentTimesWithRate(leaf_times, rate_grid),
            obs=coal_times)

# Version 2. using pyro.plate
likelihood = CoalescentRateLikelihood(leaf_times, coal_times, len(rate_grid))
with pyro.plate("time", len(rate_grid)):
    pyro.factor("coalescent", likelihood(rate_grid))

# Version 3. using pyro.markov
likelihood = CoalescentRateLikelihood(leaf_times, coal_times, len(rate_grid))
for t in pyro.markov(range(len(rate_grid))):
    pyro.factor("coalescent_{}".format(t), likelihood(rate_grid[t], t))

第三个版本对于例如 SMCFilter 是有用的，其中 rate_grid 可能是 按顺序计算的。

Parameters

  • leaf_times (torch.Tensor) – 采样事件时间的张量，即 系统发育中的叶节点。这些可以是任意实数，具有 任意顺序和重复。

  • coal_times (torch.Tensor) – 一个合并时间的张量。这些表示沿着尾随维度的大小为 leaf_times.size(-1) - 1 的集合，并且应该沿着该维度排序。

  • duration (int) – 速率网格的大小，rate_grid.size(-1)。

__call__(rate_grid, t=slice(None, None, None))[来源]

计算[1]方程7-9在一个或所有时间点的可能性。

参考文献

[1] A. Popinga, T. Vaughan, T. Statler, A.J. Drummond (2014)

“使用贝叶斯合并推断推断流行病学动态：确定性和随机模型的优点” https://arxiv.org/pdf/1407.1792.pdf

Parameters

  • rate_grid (torch.Tensor) – 基础合并率的张量（合并的成对率）。例如，在一个简单的SIR模型中，这可能是 beta S / I。最右边的维度是时间，这个张量表示在时间上是分段恒定的（一批）速率。

  • time (int 或 slice) – 可选的时间索引，用于切片输入，如 rate_grid[..., t] 这可以是顺序模型的整数 或向量化模型的 slice(None)。

Returns

似然 p(coal_times | leaf_times, rate_grid)，或对应于单个时间步长的似然部分。

Return type

torch.Tensor

bio_phylo_to_times(tree, *, get_time=None)[source]

从系统发育中提取合并的摘要统计信息，适合与CoalescentRateLikelihood一起使用。

Parameters

  • tree (Bio.Phylo.BaseTree.Clade) – 一个系统发育树。

  • get_time (callable) – 可选的函数，用于提取每个子Clade的时间点。如果不存在，时间将通过累积的.branch_length计算。

Returns

一对 Tensor s (leaf_times, coal_times) 其中 leaf_times 是采样事件的时间（系统发育树中的叶节点），而 coal_times 是合并事件的时间（系统发育二叉树中的叶节点）。

Return type

tuple