日期
Jan 8, 2025
SVEAICFR模型
1 概述
1.1 模型假设
1.人群划分：在本模型中，人群首先被分为两大类：非接种疫苗组(1组S1)和接种疫苗组(2组V2)，和3个年龄组（a：＜18岁，b:18-64岁，c:≥65岁）。每个大类又细分为七个类别：易感者(S), 暴露者(E)，有症状感染者(I)，无症状感染者(A)，重症感染者(C)，死亡的感染者(F)，以及康复者(R)。
2.人口变化：模型假设人群分布是同质均匀的，未考虑人群出生、死亡对疾病传播的影响。
3.传播情况：该传染病通过易感人群接触有症状感染者、无症状感染者以β、κβ (0 < κ < 1)的速率实现传播。
4.死亡：假设该传染病死亡只发生在重症患者中。
5.疫苗覆盖率：假设人群中该传染病的疫苗覆盖率为n，流感疫苗接种后，人群易感性、传染性、重症率和死亡率会减弱。
6.人口流动性：使用20XX年XX迁徙指数模拟人口流动状态。
7.康复者永久免疫：康复者永久免疫该传染病，无再感染风险。
8.分为三个年龄组人群：三个年龄组人群之间会相互影响和传播。
1.2 模型描述
将人群分为三个年龄组（<18，18-64和≥65岁人群）和两大类：非接种疫苗组(1组S1)和接种疫苗组(2组V2)，以20XX年XX月至20XX年XX月在中国31个省市自治区（香港、澳门和台湾除外）的某传染病病例负担（假设的验证指标：VI）和接种情况作为基准，构建可用的SVEAICFR模型（易感者-接种者-暴露者-无症状感染者-感染者-重症者-致命者-康复者/移除者）。
同时上调各省疫苗接种率，预测不同疫苗接种情境、不同人群策略下的疾病负担变化。
1.3 模型仓室及参数含义
| 参数 | 定义 | 单位 |
|---|---|---|
| β | 人传人传播系数 | |
| k | 无症状感染者相对症状感染者的传播能力系数 | |
| p | 无症状感染者的比例 | % |
| q | 症状感染到重症转化率 | % |
| f | 重症感染者病死率 | % |
| ωS | 潜伏期速率 | day-1 |
| ωA | 潜隐期速率 | day-1 |
| γA | 无症状感染者恢复期速率 | day-1 |
| γI | 症状感染者恢复期速率 | day-1 |
| γI2 | 重症状感染者恢复期速率 | day-1 |
| α | 症状感染到重症的速率 | day-1 |
| Mt | 各省每日人口流入和流出之差 | 人 |
| Pt | 当前时刻各省各自然史状态人数占总人口的比例 | % |
| ni | 各组疫苗覆盖率 | % |
| VES | 疫苗接种的降低人群易感性效果 | % |
| VEI | 疫苗接种的降低人群传染性效果 | % |
| VEC | 疫苗接种的降低人群重症率效果 | % |
| VEF | 疫苗接种的降低人群死亡率效果 | % |
2 模型流程图及方程
$$
\begin{aligned}
\frac{dS1i}{dt} &= \frac{-\beta ij S1i[(kA1i + I1i) + (kA2i + I2i)(1 - VEIi)]}{Ni} + M_t P_t \\
\frac{dE1i}{dt} &= \frac{\beta ij S1i[(kA1i + I1i) + (kA2i + I2i)(1 - VESi)] }{Ni} - p\omega_A E1i - (1-p)\omega_S E1i + M_t P_t \\
\frac{dI1i}{dt} &= (1-p)\omega_S E1i - (1-q)\gamma_I I1i - q\alpha I1i + M_t P_t \\
\frac{dA1i}{dt} &= p\omega_A E1i - \gamma_A A1i + M_t P_t \\
\frac{dC1i}{dt} &= q\alpha I1i - fC1i - (1-f)\gamma_{I2} C1i \\
\frac{dF1i}{dt} &= fC1i \\
\frac{dV2i}{dt} &= \frac{-\beta ij V2i(1 - VESi)[(I1i + kA1i) + (I2i + kA2i)(1 - VEIi)]}{Ni} + M_t P_t \\
\frac{dE2i}{dt} &= \frac{\beta ij V2i(1 - VESi)[(I1i + kA1i) + (I2i + kA2i)(1 - VEIi)]}{Ni} - p\omega_A E2i(1 - VEIi) - (1-p)\omega_S E2i(1 - VEIi) + M_t P_t \\
\frac{dI2i}{dt} &= (1-p)\omega_S E2i (1 - VEIi) - q\alpha I2i (1 - VECi) - (1-q)\gamma_I I2i + M_t P_t \\
\frac{dA2i}{dt} &= p\omega_A E2i (1 - VEIi) - \gamma_A A2i + M_t P_t \\
\frac{dC2i}{dt} &= q\alpha I2i (1 - VECi) - fC2i (1 - VEFi) - (1-f) \gamma_{I2} C2i \\
\frac{dF2i}{dt} &= fC2i (1 - VEFi) \\
\frac{dRi}{dt} &= (1-f) \gamma_{I2} C1i + (1-q) \gamma_I I1i + \gamma_A A1i + \gamma_A A2i + (1-q) \gamma_I I2i + (1-f) \gamma_{I2} C2i \\
Ni &= S1i + E1i + I1i + C1i + F1i + A1i + V2i + E2i + I2i + A2i + C2i + F2i + Ri
\end{aligned}
$$
3 适用疾病流行条件
适用于存在疾病病程中存在易感者-暴露者-无症状感染者-感染者-重症者-致命者-康复者/移除者七种状态，无/有相应疫苗实施接种皆可，疫苗接种后某传染病的传染性、人群易感性、人群重症率和人群死亡率会减弱的疾病，例如流感。
4 代码（Python）
4.1 条件设定
本模型模拟了全国31个省市自治区（香港、澳门和台湾除外）的某传染病流行趋势情况，并通过情景模拟：
1.使用全人群策略，将各省市自治区全人群疫苗接种率在基线基础上提高0%至60%，预测传染病发病、重症、死亡数的疾病负担变化。
2.使用高危人群策略，将各省市自治区的<18岁和≥65岁人群疫苗接种率在基线基础上提高0%至60%，预测传染病发病、重症、死亡数的疾病负担变化。
适用于存在疾病病程中存在易感者-暴露者-无症状感染者-感染者-重症者-致命者-康复者/移除者七种状态，无/有相应疫苗实施接种，疫苗接种后某传染病的传染性、人群易感性、人群重症率和人群死亡率会减弱的疾病，例如流感。
4.2 代码
1. Basic
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import pandas as pd
import json
import warnings
warnings.filterwarnings("ignore", category=FutureWarning)
index_to_province = {
0: 'Hebei',
1: 'Xizang',
2: 'Hainan',
3: 'Liaoning',
4: 'Heilongjiang',
5: 'Henan',
6: 'Sichuan',
7: 'Qinghai',
8: 'Shandong',
9: 'Guizhou',
10: 'Xinjiang',
11: 'Gansu',
12: 'Guangdong',
13: 'Shanghai',
14: 'Hunan',
15: 'Beijing',
16: 'Yunnan',
17: 'Shanxi',
18: 'Hubei',
19: 'Chongqing',
20: 'Fujian',
21: 'Ningxia',
22: 'Tianjin',
23: 'Shaanxi',
24: 'Anhui',
25: 'Jiangsu',
26: 'Inner Mongolia',
27: 'Zhejiang',
28: 'Jiangxi',
29: 'Guangxi',
30: 'Jilin'
}
2. Data and params
#Baidu migration data
flow_matrix_number_3d = np.load(r"C:\Users\xieyulun\Desktop\案例汇报\seir\\flow_matrix_number_3d_2023.npy")
sliced_data1 = flow_matrix_number_3d[211:365]*50000
sliced_data2 = flow_matrix_number_3d[0:63]*50000
sliced_data = np.concatenate((sliced_data1, sliced_data2), axis=0)
daily_outflow_np = sliced_data.sum(axis=2)
daily_inflow_np = sliced_data.sum(axis=1)
diff_flow = daily_inflow_np-daily_outflow_np
diff_flow_sum = diff_flow.sum(axis=0)
daily_inflow_np.shape
daily_inflow = daily_inflow_np
daily_outflow = daily_outflow_np
#Initial parameter andn values
para = pd.read_excel(r"C:\Users\xieyulun\Desktop\案例汇报\SVEAICFR\\初值.xlsx")
Nas = np.array(para['Polulation_<18'])
Nbs = np.array(para['Polulation_18-64'] )
Ncs = np.array(para['Polulation_65+'] )
I0as = np.array(para['<18 cases'])
I0bs = np.array(para['18-64 cases'] )
I0cs = np.array(para['65+ cases'])
n1s = np.array(para['VC1']*0.01)
n2s = np.array(para['VC2']*0.01)
n3s = np.array(para['VC3']*0.01)
VES1 = 0.65
VES2 = 0.50
VES3 = 0.42
VEI1 = 0.8
VEI2 = 0.8
VEI3 = 0.8
VEC1 = 0.69
VEC2 = 0.49
VEC3 = 0.50
VEF1 = 0.29
VEF2 = 0.29
VEF3 = 0.29
p = 0.500
kappa= 0.333
sigma = 0.53
delta = 0.83
gamma_A = 0.25
gamma_I = 0.31
q = 0.003
f = 0.17
gamma_I2 = 0.15
alpha = 0.23
days = 216
3. Model
def model_with_migration(y, t, Na,Nb,Nc, beta1, beta2,beta3,beta12, beta21,beta31,beta13,beta32,beta23,\
sigma, delta, gamma_A, gamma_I, kappa, p, q, f, n1,n2,n3, alpha, gamma_I2, VEI1, VES1, VEC1, VEF1, \
VEI2, VES2, VEC2, VEF2,VEI3, VES3, VEC3, VEF3,inflow, outflow):
S1a, E1a, I1a, A1a, R1a, C1a, F1a, V2a, E2a, I2a, A2a, R2a, C2a, F2a, Ia, Ca, Fa, \
S1b, E1b, I1b, A1b, R1b, C1b, F1b, V2b, E2b, I2b, A2b, R2b, C2b, F2b, Ib, Cb, Fb, \
S1c, E1c, I1c, A1c, R1c, C1c, F1c, V2c, E2c, I2c, A2c, R2c, C2c, F2c, Ic, Cc, Fc = y
#group a
dS1adt = - beta1 * S1a * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Na \
- beta21 * S1a * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Na \
- beta31 * S1a * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Na \
+ (inflow - outflow) * (S1a / Na)
dE1adt = beta1 * S1a * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Na \
+ beta21 * S1a * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Na \
+ beta31 * S1a * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Na \
- sigma * E1a * (1-p) - delta * E1a * p \
+ (inflow - outflow) * (E1a / Na)
dI1adt = sigma * E1a * (1-p) - (1-q)* gamma_I * I1a + (inflow - outflow) * (I1a / Na) - q * alpha *I1a
dA1adt = delta * E1a * p - gamma_A * A1a + (inflow - outflow) * (A1a / Na)
dC1adt = q * alpha * I1a - (1-f) * gamma_I2 * C1a - f * C1a
dF1adt = f * C1a
dR1adt = (1-q)* gamma_I * I1a + gamma_A * A1a + (1-f) * gamma_I2 * C1a + (inflow - outflow) * (R1a / Na)
dV2adt = - beta1 *(1-VES1)* V2a * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Na \
- beta21 *(1-VES1)* V2a * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Na \
- beta31 *(1-VES1)* V2a * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Na \
+ (inflow - outflow) * (V2a / Na)
dE2adt = beta1 *(1-VES1)* V2a * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Na \
+ beta21 * (1-VES1)* V2a * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Na \
+ beta31 *(1-VES1)* V2a * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Na \
- sigma * E2a * (1-p) * (1-VEI1) - delta * E2a * p * (1-VEI1) \
+ (inflow - outflow) * (E2a / Na)
dI2adt = sigma * E2a * (1-p) * (1-VEI1) - (1-q)* gamma_I * I2a - q * alpha *I2a + (inflow - outflow) * (I2a / Na)
dA2adt = delta * E2a * p * (1-VEI1) - gamma_A * A2a + (inflow - outflow) * (A2a / Na)
dC2adt = q * alpha * I2a * (1-VEC1) - (1-f) * gamma_I2 * C2a - f * C2a * (1-VEF1)
dF2adt = f * C2a * (1-VEF1)
dR2adt = (1-q)* gamma_I * I2a + gamma_A * A2a + (inflow - outflow) * (R2a / Na) + (1-f) * gamma_I2 * C2a
#新发病例，重症，死亡数
dIadt = sigma * E1a * (1-p) + delta * E1a * p + sigma * E2a * (1-p) * (1-VEI1) + delta * E2a * p * (1-VEI1)
dCadt = q * alpha * I1a + q * alpha * I2a * (1-VEC1)
dFadt = f * C1a + f * C2a * (1-VEF1)
# Group b
dS1bdt = - beta2 * S1b * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nb \
- beta12 * S1b * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nb \
- beta32 * S1b * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nb \
+ (inflow - outflow) * (S1b / Nb)
dE1bdt = beta2 * S1b * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nb \
+ beta12 * S1b * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nb \
+ beta32 * S1b * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nb \
- sigma * E1b * (1-p) - delta * E1b * p \
+ (inflow - outflow) * (E1b / Nb)
dI1bdt = sigma * E1b * (1-p) - (1-q)* gamma_I * I1b + (inflow - outflow) * (I1b / Nb) - q * alpha *I1b
dA1bdt = delta * E1b * p - gamma_A * A1b + (inflow - outflow) * (A1b / Nb)
dC1bdt = q * alpha * I1b - (1-f) * gamma_I2 * C1b - f * C1b
dF1bdt = f * C1b
dR1bdt = (1-q)* gamma_I * I1b + gamma_A * A1b + (1-f) * gamma_I2 * C1b + (inflow - outflow) * (R1b / Nb)
dV2bdt = - beta2 * (1-VES2)*V2b * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nb \
- beta12*(1-VES2)* V2b * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nb \
- beta32 *(1-VES2)* V2b * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nb \
+ (inflow - outflow) * (V2b / Nb)
dE2bdt = beta2 *(1-VES2)* V2b * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nb \
+ beta12 *(1-VES2)* V2b * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nb \
+ beta32 *(1-VES2)* V2b * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nb \
- sigma * E2b * (1-p) - delta * E2b * p \
+ (inflow - outflow) * (E2b / Nb)
dI2bdt = sigma * E2b * (1-p) * (1-VEI2) - (1-q)* gamma_I * I2b - q * alpha *I2b + (inflow - outflow) * (I2b / Nb)
dA2bdt = delta * E2b * p * (1-VEI2) - gamma_A * A2b + (inflow - outflow) * (A2b / Nb)
dC2bdt = q * alpha * I2b * (1-VEC2) - (1-f) * gamma_I2 * C2b - f * C2b * (1-VEF2)
dF2bdt = f * C2b * (1-VEF2)
dR2bdt = (1-q)* gamma_I * I2b + gamma_A * A2b + (inflow - outflow) * (R2b / Nb) + (1-f) * gamma_I2 * C2b
#新发病例，重症，死亡数
dIbdt = sigma * E1b * (1-p) + delta * E1b * p + sigma * E2b * (1-p) * (1-VEI2) + delta * E2b * p * (1-VEI2)
dCbdt = q * alpha * I1b + q * alpha * I2b * (1-VEC2)
dFbdt = f * C1b + f * C2b * (1-VEF2)
# Group c
dS1cdt = - beta3 * S1c * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nc \
- beta13 * S1c * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nc \
- beta23 * S1c * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nc \
+ (inflow - outflow) * (S1c / Nc)
dE1cdt = beta3 * S1c * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nc \
+ beta13 * S1c * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nc \
+ beta23 * S1c * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nc \
- sigma * E1c * (1-p) - delta * E1c * p \
+ (inflow - outflow) * (E1c / Nc)
dI1cdt = sigma * E1c * (1-p) - (1-q)* gamma_I * I1c + (inflow - outflow) * (I1c / Nc) - q * alpha *I1c
dA1cdt = delta * E1c * p - gamma_A * A1c + (inflow - outflow) * (A1c / Nc)
dC1cdt = q * alpha * I1c - (1-f) * gamma_I2 * C1c - f * C1c
dF1cdt = f * C1c
dR1cdt = (1-q)* gamma_I * I1c + gamma_A * A1c + (1-f) * gamma_I2 * C1c + (inflow - outflow) * (R1c / Nc)
dV2cdt = - beta3* (1-VES3)* V2c * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nc \
- beta13 *(1-VES3)* V2c * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nc \
- beta23 *(1-VES3)* V2c * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nc \
+ (inflow - outflow) * (V2c / Nc)
dE2cdt = beta3 * (1-VES3)* V2c * ((I1c + kappa*A1c) + (I2c + kappa*A2c) * (1-VEI3)) / Nc \
+ beta13 *(1-VES3)* V2c * ((I1a + kappa*A1a) + (I2a + kappa*A2a) * (1-VEI1)) / Nc \
+ beta23 *(1-VES3)* V2c * ((I1b + kappa*A1b) + (I2b + kappa*A2b) * (1-VEI2)) / Nc \
- sigma * E2c * (1-p) - delta * E2c * p \
+ (inflow - outflow) * (E2c / Nc)
dI2cdt = sigma * E2c * (1-p) * (1-VEI3) - (1-q)* gamma_I * I2c - q * alpha *I2c + (inflow - outflow) * (I2c / Nc)
dA2cdt = delta * E2c * p * (1-VEI3) - gamma_A * A2c + (inflow - outflow) * (A2c / Nc)
dC2cdt = q * alpha * I2c * (1-VEC3) - (1-f) * gamma_I2 * C2c - f * C2c * (1-VEF3)
dF2cdt = f * C2c * (1-VEF3)
dR2cdt = (1-q)* gamma_I * I2c + gamma_A * A2c + (inflow - outflow) * (R2c / Nc) + (1-f) * gamma_I2 * C2c
#新发病例，重症，死亡数
dIcdt = sigma * E1c * (1-p) + delta * E1c * p + sigma * E2c * (1-p) * (1-VEI3) + delta * E2c * p * (1-VEI3)
dCcdt = q * alpha * I1c + q * alpha * I2c * (1-VEC3)
dFcdt = f * C1c + f * C2c * (1-VEF3)
Na= S1a+E1a+I1a+A1a+R1a+C1a+F1a+V2a+E2a+I2a+A2a+R2a+C2a+F2a
Nb= S1b+E1b+I1b+A1b+R1b+C1b+F1b+V2b+E2b+I2b+A2b+R2b+C2b+F2b
Nc= S1c+E1c+I1c+A1c+R1c+C1c+F1c+V2c+E2c+I2c+A2c+R2c+C2c+F2c
return dS1adt, dE1adt, dI1adt, dA1adt, dR1adt, dC1adt, dF1adt, \
dV2adt, dE2adt, dI2adt, dA2adt, dR2adt, dC2adt, dF2adt, \
dIadt, dCadt, dFadt, \
dS1bdt, dE1bdt, dI1bdt, dA1bdt, dR1bdt, dC1bdt, dF1bdt, \
dV2bdt, dE2bdt, dI2bdt, dA2bdt, dR2bdt, dC2bdt, dF2bdt, \
dIbdt, dCbdt, dFbdt, \
dS1cdt, dE1cdt, dI1cdt, dA1cdt, dR1cdt, dC1cdt, dF1cdt, \
dV2cdt, dE2cdt, dI2cdt, dA2cdt, dR2cdt, dC2cdt, dF2cdt, \
dIcdt, dCcdt, dFcdt
#VI data
VI_plus = pd.read_excel(r"C:\Users\xieyulun\Desktop\案例汇报\SVEAICFR\\VI.xlsx",index_col = 0)
4. Verification
beta1s = np.array(para['Betas'])
seiar_results_with_migration = []
for idx, (beta1, Na, Nb, Nc, I0xa, I0xb, I0xc, n1, n2, n3) in enumerate(zip(beta1s, Nas, Nbs, Ncs, I0as, I0bs, I0cs, n1s, n2s, n3s)):
beta2 = beta1 * 0.94
beta3 = beta1 * 0.55
beta12 = beta1 * 0.9
beta32 = beta1 * 0.5
beta21 = beta1 * 0.9
beta31 = beta1 * 0.5
beta13 = beta1 * 0.9
beta23 = beta1 * 0.9
E1a0, I1a0, A1a0, R1a0, C1a0, F1a0, E2a0, I2a0, A2a0, R2a0, C2a0, F2a0 = 0, I0xa * (1 - n1), 0, 0, 0, 0, 0, I0xa * n1, 0, 0, 0, 0
I0a = I0xa
C0a = 0
F0a = 0
S1a0 = Na * (1 - n1) - E1a0 - I1a0 - A1a0 - R1a0 - C1a0 - F1a0
V2a0 = Na * n1 - E2a0 - I2a0 - A2a0 - R2a0 - C2a0 - F2a0
E1b0, I1b0, A1b0, R1b0, C1b0, F1b0, E2b0, I2b0, A2b0, R2b0, C2b0, F2b0 = 0, I0xb * (1 - n2), 0, 0, 0, 0, 0, I0xb * n2, 0, 0, 0, 0
I0b = I0xb
C0b = 0
F0b = 0
S1b0 = Nb * (1 - n2) - E1b0 - I1b0 - A1b0 - R1b0 - C1b0 - F1b0
V2b0 = Nb * n2 - E2b0 - I2b0 - A2b0 - R2b0 - C2b0 - F2b0
E1c0, I1c0, A1c0, R1c0, C1c0, F1c0, E2c0, I2c0, A2c0, R2c0, C2c0, F2c0 = 0, I0xc * (1 - n3), 0, 0, 0, 0, 0, I0xc * n3, 0, 0, 0, 0
I0c = I0xc
C0c = 0
F0c = 0
S1c0 = Nc * (1 - n3) - E1c0 - I1c0 - A1c0 - R1c0 - C1c0 - F1c0
V2c0 = Nc * n3 - E2c0 - I2c0 - A2c0 - R2c0 - C2c0 - F2c0
y0 = S1a0, E1a0, I1a0, A1a0, R1a0, C1a0, F1a0, V2a0, E2a0, I2a0, A2a0, R2a0, C2a0, F2a0, I0a, C0a, F0a, \
S1b0, E1b0, I1b0, A1b0, R1b0, C1b0, F1b0, V2b0, E2b0, I2b0, A2b0, R2b0, C2b0, F2b0, I0b, C0b, F0b, \
S1c0, E1c0, I1c0, A1c0, R1c0, C1c0, F1c0, V2c0, E2c0, I2c0, A2c0, R2c0, C2c0, F2c0, I0c, C0c, F0c
t = np.linspace(0, days, days + 1)
def get_daily_migration(day, idx):
if day < days:
return daily_inflow[int(day), idx], daily_outflow[int(day), idx]
else:
return 0, 0
ret = odeint(lambda y, t: model_with_migration(y, t, Na,Nb,Nc, beta1, beta2,beta3,beta12, beta21,beta31,beta13,beta32,beta23,\
sigma, delta, gamma_A, gamma_I, kappa, p, q, f, n1,n2,n3, alpha, gamma_I2, VEI1, VES1, VEC1, VEF1, \
VEI2, VES2, VEC2, VEF2,VEI3, VES3, VEC3, VEF3,
*get_daily_migration(t, idx)), y0, t)
S1a, E1a, I1a, A1a, R1a, C1a, F1a, V2a, E2a, I2a, A2a, R2a, C2a, F2a, Ia, Ca, Fa, \
S1b, E1b, I1b, A1b, R1b, C1b, F1b, V2b, E2b, I2b, A2b, R2b, C2b, F2b, Ib, Cb, Fb, \
S1c, E1c, I1c, A1c, R1c, C1c, F1c, V2c, E2c, I2c, A2c, R2c, C2c, F2c, Ic, Cc, Fc = ret.T
seiar_results_with_migration.append({"Province": index_to_province[idx], 'Ia': Ia.tolist(), 'Ca': Ca.tolist(), 'Fa': Fa.tolist(),
'Ib': Ib.tolist(), 'Cb': Cb.tolist(), 'Fb': Fb.tolist(), 'Ic': Ic.tolist(), 'Cc': Cc.tolist(), 'Fc': Fc.tolist()})
moni_df = pd.DataFrame(seiar_results_with_migration)
moni_df.to_excel(r"C:\Users\xieyulun\Desktop\案例汇报\SVEAICFR\\模拟数据.xlsx", index=False)
5. Verification Plot
import os
plt.rc('font', family='Times New Roman')
plt.rcParams.update({'font.size': 15})
plt.rcParams.update({'font.weight': 'bold'})
fig, axs = plt.subplots(8, 4, figsize=(25, 32), sharex=True)
# 设置共享 x 轴格式和间隔
date_format = "%Y-%m" # 年-月格式
xticks_interval = 30 # 每 30 天显示一个刻度
total_cases = []
for i in range(31):
choose_number = i
example_province_seiar_migration = seiar_results_with_migration[choose_number]
province_name = index_to_province[choose_number].split('省')[0]
row = i // 4 # 计算子图行
col = i % 4 # 计算子图列
ax = axs[row, col] # 获取当前子图
ax2 = ax.twinx() # 创建共享 x 轴的第二个坐标轴
# 计算每个省份的每日abc三组的新发总病例
y = [sum(x) for x in zip(example_province_seiar_migration["Ia"], example_province_seiar_migration["Ib"], example_province_seiar_migration["Ic"])]
diff_y = [y[i] - y[i - 1] for i in range(1, len(y))]
diff_y = [y[0]] + diff_y
date_range = pd.date_range(start='2023-07-31', end='2024-03-03', freq='D')
#重采样为每周的新发病例数的变化
cases_df = pd.DataFrame(diff_y, index=date_range, columns=['Cases'])
weekly_cases_change = cases_df.resample('W').sum()
ax.set_xlim(pd.to_datetime('2023-07-31'), pd.to_datetime('2024-03-03')) # 设置 x 轴范围
ax.bar(weekly_cases_change.index, weekly_cases_change['Cases'].values, color='#ffd700', label='Model fit cases', width=6) # 绘制柱状图
fontdict = {'fontsize': 15, 'fontweight': 'bold', 'family': 'Times New Roman'}
ax.set_xlabel('Time', fontdict=fontdict)
ax.set_ylabel('Number of cases', fontdict=fontdict)
ax.set_title(f'{province_name}', fontdict=fontdict)
ILI_weekly = VI_plus[province_name]
ax2.plot(ILI_weekly, color='#FFA07A', label='VI', linewidth=2) # 绘制第二个坐标轴上的折线图
ax2.set_ylabel('Validation index', fontdict=fontdict)
# 调整 X 轴刻度
if row == 7: # 仅最后一行设置 X 轴刻度
xticks = pd.date_range(start='2023-07-31', end='2024-03-03', freq=f'{xticks_interval}D')
ax.set_xticks(xticks)
ax.set_xticklabels([tick.strftime(date_format) for tick in xticks], rotation=45, ha='right')
# 隐藏第32个子图的位置（即最后一个子图）
fig.delaxes(axs[7, 3])
# 合并图例
handles1, labels1 = ax.get_legend_handles_labels() # 获取 ax1 的图例句柄和标签
handles2, labels2 = ax2.get_legend_handles_labels() # 获取 ax2 的图例句柄和标签
# 合并图例句柄和标签
handles = handles1 + handles2
labels = labels1 + labels2
# 在全局图形中显示合并后的图例，并将其放置到指定位置
fig.legend(handles, labels, loc='lower right', ncol=2, fontsize=13, bbox_to_anchor=(0.9, 0.05), frameon=False)
plt.tight_layout() # 自动调整子图布局
plt.savefig(r"C:\\Users\xieyulun\Desktop\\案例汇报\SVEAICFR\\Plot_VI.pdf", dpi=450, bbox_inches='tight')
plt.show() # 显示图形
5 Vaccine Change（干预措施，全人群提高疫苗接种率）
seiar_results_with_migration = []
for rx in range(0,60):
for idx, (beta1,Na, Nb,Nc, I0xa, I0xb, I0xc,n1,n2,n3) in enumerate(zip(beta1s,Nas,Nbs,Ncs, I0as,I0bs,I0cs,n1s,n2s,n3s)):
n1 = n1 + rx*0.01
n2 = n2 + rx*0.01
n3 = n3 + rx*0.01
beta2 = beta1*0.94
beta3 = beta1*0.55
beta12 = beta1*0.9
beta32 = beta1*0.5
beta21 = beta1*0.9
beta31 = beta1*0.5
beta13 = beta1*0.9
beta23 = beta1*0.9
#Group a
E1a0, I1a0, A1a0, R1a0, C1a0, F1a0, E2a0, I2a0, A2a0, R2a0, C2a0, F2a0 = 0, I0xa*(1-n1), 0, 0 , 0, 0 ,0, I0xa*n1, 0, 0 , 0, 0
I0a = I0xa
C0a = 0
F0a = 0
S1a0 = Na*(1-n1) - E1a0 - I1a0 - A1a0 - R1a0 - C1a0 - F1a0
V2a0 = Na*n1 - E2a0 - I2a0 - A2a0 - R2a0 - C2a0 - F2a0
#Group b
E1b0, I1b0, A1b0, R1b0, C1b0, F1b0, E2b0, I2b0, A2b0, R2b0, C2b0, F2b0 = 0, I0xb*(1-n2), 0, 0 , 0, 0 ,0, I0xb*n2, 0, 0 , 0, 0
I0b = I0xb
C0b = 0
F0b = 0
S1b0 = Nb*(1-n2) - E1b0 - I1b0 - A1b0 - R1b0 - C1b0 - F1b0
V2b0 = Nb*n2 - E2b0 - I2b0 - A2b0 - R2b0 - C2b0 - F2b0
#Group c
E1c0, I1c0, A1c0, R1c0, C1c0, F1c0, E2c0, I2c0, A2c0, R2c0, C2c0, F2c0 = 0, I0xc*(1-n3), 0, 0 , 0, 0 ,0, I0xc*n3, 0, 0 , 0, 0
I0c = I0xc
C0c = 0
F0c = 0
S1c0 = Nc*(1-n3) - E1c0 - I1c0 - A1c0 - R1c0 - C1c0 - F1c0
V2c0 = Nc*n3 - E2c0 - I2c0 - A2c0 - R2c0 - C2c0 - F2c0
y0 = S1a0, E1a0, I1a0, A1a0, R1a0, C1a0, F1a0, V2a0, E2a0, I2a0, A2a0, R2a0, C2a0, F2a0, I0a, C0a, F0a, \
S1b0, E1b0, I1b0, A1b0, R1b0, C1b0, F1b0, V2b0, E2b0, I2b0, A2b0, R2b0, C2b0, F2b0, I0b, C0b, F0b, \
S1c0, E1c0, I1c0, A1c0, R1c0, C1c0, F1c0, V2c0, E2c0, I2c0, A2c0, R2c0, C2c0, F2c0, I0c, C0c, F0c
t = np.linspace(0, days, days+1)
def get_daily_migration(day, idx):
if day < days:
return daily_inflow[int(day), idx], daily_outflow[int(day), idx]
else:
return 0, 0
ret = odeint(lambda y, t: model_with_migration(y, t, Na,Nb,Nc, beta1, beta2,beta3,beta12, beta21,beta31,beta13,beta32,beta23,\
sigma, delta, gamma_A, gamma_I, kappa, p, q, f, n1,n2,n3, alpha, gamma_I2, VEI1, VES1, VEC1, VEF1, \
VEI2, VES2, VEC2, VEF2,VEI3, VES3, VEC3, VEF3,
*get_daily_migration(t, idx)), y0, t)
S1a, E1a, I1a, A1a, R1a, C1a, F1a, V2a, E2a, I2a, A2a, R2a, C2a, F2a, Ia, Ca, Fa, \
S1b, E1b, I1b, A1b, R1b, C1b, F1b, V2b, E2b, I2b, A2b, R2b, C2b, F2b, Ib, Cb, Fb, \
S1c, E1c, I1c, A1c, R1c, C1c, F1c, V2c, E2c, I2c, A2c, R2c, C2c, F2c, Ic, Cc, Fc = ret.T
S1a, E1a, I1a, A1a, R1a, C1a, F1a, V2a, E2a, I2a, A2a, R2a, C2a, F2a, Ia, Ca, Fa = S1a.tolist(), E1a.tolist(), I1a.tolist(), A1a.tolist(), R1a.tolist(), C1a.tolist(), F1a.tolist(), V2a.tolist(), E2a.tolist(), I2a.tolist(), A2a.tolist(), R2a.tolist(), C2a.tolist(), F2a.tolist(), Ia.tolist(), Ca.tolist(), Fa.tolist()
S1b, E1b, I1b, A1b, R1b, C1b, F1b, V2b, E2b, I2b, A2b, R2b, C2b, F2b, Ib, Cb, Fb = S1b.tolist(), E1b.tolist(), I1b.tolist(), A1b.tolist(), R1b.tolist(), C1b.tolist(), F1b.tolist(), V2b.tolist(), E2b.tolist(), I2b.tolist(), A2b.tolist(), R2b.tolist(), C2b.tolist(), F2b.tolist(), Ib.tolist(), Cb.tolist(), Fb.tolist()
S1c, E1c, I1c, A1c, R1c, C1c, F1c, V2c, E2c, I2c, A2c, R2c, C2c, F2c, Ic, Cc, Fc = S1c.tolist(), E1c.tolist(), I1c.tolist(), A1c.tolist(), R1c.tolist(), C1c.tolist(), F1c.tolist(), V2c.tolist(), E2c.tolist(), I2c.tolist(), A2c.tolist(), R2c.tolist(), C2c.tolist(), F2c.tolist(), Ic.tolist(), Cc.tolist(), Fc.tolist()
Ia_sum = Ia[-1]
Ib_sum = Ib[-1]
Ic_sum = Ic[-1]
Ca_sum = Ca[-1]
Cb_sum = Cb[-1]
Cc_sum = Cc[-1]
Fa_sum = Fa[-1]
Fb_sum = Fb[-1]
Fc_sum = Fc[-1]
seiar_results_with_migration.append({'rx':rx,'n1':n1,'n2':n2,'n3':n3,'Province': index_to_province[idx],'未成年人累积病例数':Ia_sum,'中年人累积病例数' :Ib_sum, '老年人累积病例数':Ic_sum,
'未成年人累积重症数':Ca_sum,'中年人累积重症数' :Cb_sum, '老年人累积重症数':Cc_sum,
'未成年人累积死亡数':Fa_sum,'中年人累积死亡数' :Fb_sum, '老年人累积死亡数':Fc_sum,})
res_df = pd.DataFrame(seiar_results_with_migration)
res_df.to_excel(r"C:\\Users\xieyulun\Desktop\\案例汇报\SVEAICFR\\全人群策略疫苗接种率提高.xlsx")
5.1 Vaccine Change Plot
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
# 读取Excel文件
file_path = r"C:\\Users\xieyulun\\Desktop\\案例汇报\\SVEAICFR\\全人群策略疫苗接种率提高.xlsx"
df = pd.read_excel(file_path)
# 定义需要相加的列
columns_to_sum = {
'infected cases': ['未成年人累积病例数', '中年人累积病例数', '老年人累积病例数'],
'critical cases': ['未成年人累积重症数', '中年人累积重症数', '老年人累积重症数'],
'fatality cases': ['未成年人累积死亡数', '中年人累积死亡数', '老年人累积死亡数']
}
# 初始化结果字典
results = {key: [] for key in columns_to_sum.keys()}
# 对每个RX情景进行处理
for rx in df['rx'].unique():
rx_df = df[df['rx'] == rx]
for result_key, column_keys in columns_to_sum.items():
total = rx_df[column_keys].sum().sum() # 先对列求和，再对行求和
results[result_key].append(total)
# 将结果转换为DataFrame
results_df = pd.DataFrame(results, index=df['rx'].unique())
# 定义颜色、字体
colors = ['#87CEEB', '#98FB98', '#FF1493']
plt.rc('font', family='Times New Roman')
plt.rcParams.update({'font.size': 15})
plt.rcParams.update({'font.weight': 'bold'})
# 创建一个包含三个子图的图形，子图排成一行
fig, axes = plt.subplots(1, 3, figsize=(18, 6))
# 绘制条形图
bar_width = 0.2
index = np.arange(len(results_df.index))
# 为每个指标绘制一张图
for i, (metric, values) in enumerate(results_df.items()):
ax = axes[i]
ax.bar(index, values, bar_width, color=colors[i], label=metric)
# 添加图例
ax.legend()
# 添加标题和标签
ax.set_title(f' {metric} under Variou Vaccine Rates')
ax.set_xlabel('Increase Vaccine Rates (Universal Strategy)')
ax.set_ylabel('Total Cases')
# 调整 X 轴刻度，间隔显示一个刻度
xticks_to_show = index[::3] # 每隔两个显示一个刻度
ax.set_xticks(xticks_to_show)
ax.set_xticklabels(results_df.index[::3]) # 间隔两个显示rx情景名称
# 去掉上边框和右边框
ax.spines['top'].set_visible(False)
ax.spines['right'].set_visible(False)
# 显示图表
plt.tight_layout() # 自动调整子图布局
plt.savefig(r"C:\\Users\xieyulun\Desktop\\案例汇报\SVEAICFR\\Universal Strategy.pdf", dpi=450, bbox_inches='tight')
plt.show()
6 Vaccine Change（干预措施，高危人群：未成年组和老年组提高疫苗接种率）
seiar_results_with_migration = []
for rx in range(0,60):
for idx, (beta1,Na, Nb,Nc, I0xa, I0xb, I0xc,n1,n2,n3) in enumerate(zip(beta1s,Nas,Nbs,Ncs, I0as,I0bs,I0cs,n1s,n2s,n3s)):
n1 = n1+ rx*0.01
n3 = n3 + rx*0.01
beta2 = beta1*0.94
beta3 = beta1*0.55
beta12 = beta1*0.9
beta32 = beta1*0.5
beta21 = beta1*0.9
beta31 = beta1*0.5
beta13 = beta1*0.9
beta23 = beta1*0.9
#Group a
E1a0, I1a0, A1a0, R1a0, C1a0, F1a0, E2a0, I2a0, A2a0, R2a0, C2a0, F2a0 = 0, I0xa*(1-n1), 0, 0 , 0, 0 ,0, I0xa*n1, 0, 0 , 0, 0
I0a = I0xa
C0a = 0
F0a = 0
S1a0 = Na*(1-n1) - E1a0 - I1a0 - A1a0 - R1a0 - C1a0 - F1a0
V2a0 = Na*n1 - E2a0 - I2a0 - A2a0 - R2a0 - C2a0 - F2a0
#Group b
E1b0, I1b0, A1b0, R1b0, C1b0, F1b0, E2b0, I2b0, A2b0, R2b0, C2b0, F2b0 = 0, I0xb*(1-n2), 0, 0 , 0, 0 ,0, I0xb*n2, 0, 0 , 0, 0
I0b = I0xb
C0b = 0
F0b = 0
S1b0 = Nb*(1-n2) - E1b0 - I1b0 - A1b0 - R1b0 - C1b0 - F1b0
V2b0 = Nb*n2 - E2b0 - I2b0 - A2b0 - R2b0 - C2b0 - F2b0
#Group c
E1c0, I1c0, A1c0, R1c0, C1c0, F1c0, E2c0, I2c0, A2c0, R2c0, C2c0, F2c0 = 0, I0xc*(1-n3), 0, 0 , 0, 0 ,0, I0xc*n3, 0, 0 , 0, 0
I0c = I0xc
C0c = 0
F0c = 0
S1c0 = Nc*(1-n3) - E1c0 - I1c0 - A1c0 - R1c0 - C1c0 - F1c0
V2c0 = Nc*n3 - E2c0 - I2c0 - A2c0 - R2c0 - C2c0 - F2c0
y0 = S1a0, E1a0, I1a0, A1a0, R1a0, C1a0, F1a0, V2a0, E2a0, I2a0, A2a0, R2a0, C2a0, F2a0, I0a, C0a, F0a, \
S1b0, E1b0, I1b0, A1b0, R1b0, C1b0, F1b0, V2b0, E2b0, I2b0, A2b0, R2b0, C2b0, F2b0, I0b, C0b, F0b, \
S1c0, E1c0, I1c0, A1c0, R1c0, C1c0, F1c0, V2c0, E2c0, I2c0, A2c0, R2c0, C2c0, F2c0, I0c, C0c, F0c
t = np.linspace(0, days, days+1)
def get_daily_migration(day, idx):
if day < days:
return daily_inflow[int(day), idx], daily_outflow[int(day), idx]
else:
return 0, 0
ret = odeint(lambda y, t: model_with_migration(y, t, Na,Nb,Nc, beta1, beta2,beta3,beta12, beta21,beta31,beta13,beta32,beta23,\
sigma, delta, gamma_A, gamma_I, kappa, p, q, f, n1,n2,n3, alpha, gamma_I2, VEI1, VES1, VEC1, VEF1, \
VEI2, VES2, VEC2, VEF2,VEI3, VES3, VEC3, VEF3,
*get_daily_migration(t, idx)), y0, t)
S1a, E1a, I1a, A1a, R1a, C1a, F1a, V2a, E2a, I2a, A2a, R2a, C2a, F2a, Ia, Ca, Fa, \
S1b, E1b, I1b, A1b, R1b, C1b, F1b, V2b, E2b, I2b, A2b, R2b, C2b, F2b, Ib, Cb, Fb, \
S1c, E1c, I1c, A1c, R1c, C1c, F1c, V2c, E2c, I2c, A2c, R2c, C2c, F2c, Ic, Cc, Fc = ret.T
S1a, E1a, I1a, A1a, R1a, C1a, F1a, V2a, E2a, I2a, A2a, R2a, C2a, F2a, Ia, Ca, Fa = S1a.tolist(), E1a.tolist(), I1a.tolist(), A1a.tolist(), R1a.tolist(), C1a.tolist(), F1a.tolist(), V2a.tolist(), E2a.tolist(), I2a.tolist(), A2a.tolist(), R2a.tolist(), C2a.tolist(), F2a.tolist(), Ia.tolist(), Ca.tolist(), Fa.tolist()
S1b, E1b, I1b, A1b, R1b, C1b, F1b, V2b, E2b, I2b, A2b, R2b, C2b, F2b, Ib, Cb, Fb = S1b.tolist(), E1b.tolist(), I1b.tolist(), A1b.tolist(), R1b.tolist(), C1b.tolist(), F1b.tolist(), V2b.tolist(), E2b.tolist(), I2b.tolist(), A2b.tolist(), R2b.tolist(), C2b.tolist(), F2b.tolist(), Ib.tolist(), Cb.tolist(), Fb.tolist()
S1c, E1c, I1c, A1c, R1c, C1c, F1c, V2c, E2c, I2c, A2c, R2c, C2c, F2c, Ic, Cc, Fc = S1c.tolist(), E1c.tolist(), I1c.tolist(), A1c.tolist(), R1c.tolist(), C1c.tolist(), F1c.tolist(), V2c.tolist(), E2c.tolist(), I2c.tolist(), A2c.tolist(), R2c.tolist(), C2c.tolist(), F2c.tolist(), Ic.tolist(), Cc.tolist(), Fc.tolist()
Ia_sum = Ia[-1]
Ib_sum = Ib[-1]
Ic_sum = Ic[-1]
Ca_sum = Ca[-1]
Cb_sum = Cb[-1]
Cc_sum = Cc[-1]
Fa_sum = Fa[-1]
Fb_sum = Fb[-1]
Fc_sum = Fc[-1]
seiar_results_with_migration.append({'rx':rx,'n1':n1,'n2':n2,'n3':n3,'Province': index_to_province[idx],'未成年人累积病例数':Ia_sum,'中年人累积病例数' :Ib_sum, '老年人累积病例数':Ic_sum,
'未成年人累积重症数':Ca_sum,'中年人累积重症数' :Cb_sum, '老年人累积重症数':Cc_sum,
'未成年人累积死亡数':Fa_sum,'中年人累积死亡数' :Fb_sum, '老年人累积死亡数':Fc_sum,})
res_df = pd.DataFrame(seiar_results_with_migration)
res_df.to_excel(r"C:\\Users\xieyulun\Desktop\\案例汇报\SVEAICFR\\高危人群策略疫苗接种率提高.xlsx")
6.1 Vaccine Change Plot
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
# 读取Excel文件
file_path = r"C:\\Users\xieyulun\\Desktop\\案例汇报\\SVEAICFR\\高危人群策略疫苗接种率提高.xlsx"
df = pd.read_excel(file_path)
# 定义需要相加的列
columns_to_sum = {
'infected cases': ['未成年人累积病例数', '中年人累积病例数', '老年人累积病例数'],
'critical cases': ['未成年人累积重症数', '中年人累积重症数', '老年人累积重症数'],
'fatality cases': ['未成年人累积死亡数', '中年人累积死亡数', '老年人累积死亡数']
}
# 初始化结果字典
results = {key: [] for key in columns_to_sum.keys()}
# 对每个RX情景进行处理
for rx in df['rx'].unique():
rx_df = df[df['rx'] == rx]
for result_key, column_keys in columns_to_sum.items():
total = rx_df[column_keys].sum().sum() # 先对列求和，再对行求和
results[result_key].append(total)
# 将结果转换为DataFrame
results_df = pd.DataFrame(results, index=df['rx'].unique())
# 定义颜色、字体
colors = ['#87CEEB', '#98FB98', '#FF1493']
plt.rc('font', family='Times New Roman')
plt.rcParams.update({'font.size': 15})
plt.rcParams.update({'font.weight': 'bold'})
# 创建一个包含三个子图的图形，子图排成一行
fig, axes = plt.subplots(1, 3, figsize=(18, 6))
# 绘制条形图
bar_width = 0.2
index = np.arange(len(results_df.index))
# 为每个指标绘制一张图
for i, (metric, values) in enumerate(results_df.items()):
ax = axes[i]
ax.bar(index, values, bar_width, color=colors[i], label=metric)
# 添加图例
ax.legend()
# 添加标题和标签
ax.set_title(f' {metric} under Variou Vaccine Rates')
ax.set_xlabel('Increase Vaccine Rates (High-Risk Population Strategy)')
ax.set_ylabel('Total Cases')
# 调整 X 轴刻度，间隔显示一个刻度
xticks_to_show = index[::3] # 每隔两个显示一个刻度
ax.set_xticks(xticks_to_show)
ax.set_xticklabels(results_df.index[::3]) # 间隔两个显示rx情景名称
# 去掉上边框和右边框
ax.spines['top'].set_visible(False)
ax.spines['right'].set_visible(False)
# 显示图表
plt.tight_layout() # 自动调整子图布局
plt.savefig(r"C:\\Users\xieyulun\Desktop\\案例汇报\SVEAICFR\\High-Risk Population Strategy.pdf", dpi=450, bbox_inches='tight')
plt.show()
7 结果
提高疫苗接种率可以显著降低疾病负担，总体而言，全人群策略在基线基础上提高XX%左右的疫苗接种率，疾病负担变化达到平稳阈值；高危人群策略在基线基础上提高XX%左右的疫苗接种率，疾病负担变化达到平稳阈值。 选择哪种疫苗策略，需要考虑不同情景下的成本效益问题。