20250720

import numpy as np

import matplotlib.pyplot as plt

from scipy.integrate import odeint

import pandas as pd

import matplotlib.dates as mdates

plt.rcParams['font.family'] = 'Malgun Gothic'

plt.rcParams['axes.unicode_minus'] = False

def sirv_model_unique(y, t, N, beta_func, gamma, reinfection_rate, vaccination_schedule):

U, S, I, R, C = y

beta = beta_func(t)

vac_effect = 0

for start_day, rate in vaccination_schedule:

if t >= start_day:

vac_effect += rate

new_infections = beta * S * I / N

first_infections = new_infections * (U / S if S > 0 else 0)

dUdt = -first_infections

dSdt = -new_infections + reinfection_rate * R - vac_effect * S

dIdt = new_infections - gamma * I

dRdt = gamma * I - reinfection_rate * R + vac_effect * S

dCdt = first_infections

return dUdt, dSdt, dIdt, dRdt, dCdt

def beta_time_dependent(t):

gamma = 1/14

# t는 0부터 시작, 0 = 2020-01-01

if t <= 365:

R0 = 2.5 # 초기 감염 확산

elif t <= 600:

R0 = 1.2 # 사회적 거리두기 등으로 감소

elif t <= 800:

R0 = 4.0 # 델타 변이 등장, 급증

elif t <= 1000:

R0 = 1.5 # 델타 이후 완화

elif t <= 1200:

R0 = 5.0 # 오미크론 등장, 급격한 확산

else:

R0 = 1.0 # 오미크론 이후 완화

return R0 * gamma

def run_covid_simulation_and_analysis():

N = 50_000_000 # 인구 수

I0, R0, C0 = 10, 0, 10 # 초기 감염자, 회복자, 누적 감염자

S0 = N - I0 - R0

U0 = N - I0

gamma = 1/14

reinfection_rate = 1/180

vaccination_schedule = [(500, 0.002), (900, 0.003)] # 백신 접종 시작일과 효과

start_date = pd.to_datetime('2020-01-01')

end_date = pd.to_datetime('2024-06-30')

total_days = (end_date - start_date).days

t = np.linspace(0, total_days, total_days + 1)

y0 = U0, S0, I0, R0, C0

ret = odeint(sirv_model_unique, y0, t,

args=(N, beta_time_dependent, gamma, reinfection_rate, vaccination_schedule))

U, S, I, R, C = ret.T

U_pct, S_pct, I_pct, R_pct, C_pct = (U / N) * 100, (S / N) * 100, (I / N) * 100, (R / N) * 100, (C / N) * 100

dates = start_date + pd.to_timedelta(t, unit='D')

plt.figure(figsize=(14, 7))

plt.plot(dates, I_pct, 'r', linewidth=2, label='현재 감염자 (I%)')

plt.plot(dates, R_pct, 'g', linewidth=2, label='회복자 (R%)')

plt.plot(dates, C_pct, 'm--', linewidth=2, label='감염 경험자 (C%)')

plt.plot(dates, U_pct, 'c--', linewidth=2, label='아직 안 걸린 사람 (U%)')

plt.xlabel('날짜')

plt.ylabel('인구 대비 비율 (%)')

plt.title('현실 반영 코로나 SIR+ 모델 시뮬레이션')

plt.legend()

plt.grid()

plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[1, 6]))

plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%Y-%m'))

plt.xticks(rotation=45)

plt.tight_layout(rect=[0, 0, 1, 0.95])

plt.show()

if __name__ == "__main__":

run_covid_simulation_and_analysis()