The new coronavirus infection (COVID-19) is spreading from China to Europe and the Americas and is rampant. In Japan as well, as of March 27, there are 1,254 infected people in Japan, and the governor of Tokyo is requesting restrictions on weekend outings in Tokyo as well. However, compared to the G7 countries where medical care seems to be relatively well-developed, the number of people in Japan is small, with 85594 in the United States, 80589 in Italy, 43938 in Germany, 9155 in France, 11658 in the United Kingdom, and 4043 in Canada. It stands out. In this article, as one hypothesis to explain this reason, we examined whether it could be explained by ** a mathematical model focusing on cluster countermeasures **. The flow of examination is the problems of the conventional model, the examination of the new model, the simulation, and the consideration.
First of all, from the conclusion,
In the article below, I have calculated using a model called the SEIR model, which calculates the transition of infectious diseases.
The SEIR model is simply S: Susceptible for infectious diseases, E: Exposed for infectious diseases, I: Infectious, R: Infectious diseases It is a model that expresses how the distribution of a certain group changes among the four people who recovered from the disease and acquired immunity (Recovered) by a differential equation. It seems that it is often used in the trend analysis of infectious diseases, and it often appears in the explanation materials of the Ministry of Health and Welfare in Japan, the Ministry of Health in the United Kingdom, and the US CDC. However, if you try to calculate the population group S based on a large group such as the whole country or the whole prefecture, it will not be a digit at all **. For example, in the calculation of the SEIR model examined in the above article, assuming a population of 1000 people, it peaks in 60 days, and at the peak, more than 30% of the population becomes infected, but in terms of population ratio Even if you look at Italy, where the number of infected people is extremely high, the current number of infected people is only about 0.1% of the population of 60.48 million. There are various possible reasons, but the following are assumed, for example.
So, let's think about ** Is there a model that can explain the facts more from the actually observed data **?
Therefore, let's check the phenomenon that was actually observed. The figure below is the result calculated by Let's examine the convergence time from the global trend of the effective reproduction number of the new coronavirus.
The figures in the graph are calculated based on a certain standard for the effective basic reproduction number (** below, referred to as Rp. **). The dotted line is an approximate expression (log-linear) that approximates the graph. The trends that can be read from this graph are as follows.
So ** let's obediently model this observed trend **. In particular,
If you write it in a mathematical formula, it will look like this.
Rp(t) = Rp(t_0) \cdot 2^{-\frac{t-t_0}{T_h}}, T_h=7.5[days]
Incorporate this effective reproduction number Rp (t) into the SEIR model to make it a ** cluster-based SEIR model **.
\begin{eqnarray}
\frac{dRp}{dt} &=& - \frac{ln2}{T_h}Rp \\
\frac{dS}{dt} &=& -\frac{1}{ip} Rp \cdot I \\
\frac{dE}{dt} &=& -\frac{dS}{dt} - \frac{1}{lp} E \\
\frac{dI}{dt} &=& \frac{1}{lp} E - \frac{1}{ip} I \\
\frac{dR}{dt} &=& \frac{1}{ip} I
\end{eqnarray}
The point is that, unlike the original SEIR model, it does not depend on the initial value of ** S (t) **. In the following, I would like to examine the effect of cluster countermeasures on this model, that is, the effect of discovering clusters and discovering and isolating infected people originating from clusters **.
Let's calculate the above cluster-based SEIR model in Python. Import the library.
import numpy as np
import matplotlib.pyplot as plt
Define the ODE. Here, in order to calculate the effect of cluster countermeasures, after a certain time $ T_c $, it is calculated with $ Rp = 0 $.
#define differencial equation of seir model
def seir_eq7(v, t, keys):
Th = keys['Th']
lp = keys['lp']
ip = keys['ip']
Tc = keys['Tc']
#
Rp = v[0];
s = v[1];
e = v[2];
i = v[3];
r = v[4];
#
if t >= Tc:
Rp = 0
#
dRp = - np.log(2)/Th * Rp
ds = - Rp/ip * i
de = - ds - (1/lp) * e
di = (1/lp)*e - (1/ip)*i
dr = (1/ip)*i
return [dRp, ds, de, di, dr]
Define a function to solve the ODE.
def my_odeint(deq, ini_state, tseq, keys):
sim = None
v = np.array(ini_state).astype(np.float64)
dt = (tseq[1] - tseq[0])*1.0
for t in tseq:
dv = deq(v,t, keys)
v = v + np.array(dv) * dt
if sim is None:
sim = v
else:
sim = np.vstack((sim, v))
return sim
This is the code to simulate by changing the initial cluster value $ Rp (t_0) $ in various ways. The parameters are $ lp = 5, ip = 8, Th = 7.5 $.
# case 7.1
def calcsim(Rp0, keys):
# solve seir model
# Rp S E I R
ini_state=[Rp0, 0, 0, 1, 0]
t_max=180
dt=0.01
t=np.arange(0,t_max,dt)
#
sim = my_odeint(seir_eq7, ini_state, t, keys)
#
plt.rcParams["font.size"] = 12
fig, ax = plt.subplots(figsize=(10,5))
ax.plot(t,sim[:,[2,3,4]]) # extract Rp, E I R
ax.set_xticks(np.linspace(0,t_max,19))
ax.grid(which='both')
ax.legend([ 'Exposed', 'Infected','Recoverd'])
ax.set_xlabel('date')
ax.set_ylim(0,)
plt.show()
print("R(tmax):{}".format(sim[-1,4]))
# try different Rp0
keys = {'lp':5, 'ip':8, 'Th':7.5, 'Tc':200 }
calcsim(100, keys)
calcsim(30, keys)
calcsim(20, keys)
calcsim(10, keys)
# try different Tc
keys = {'lp':5, 'ip':8, 'Th':7.5, 'Tc':20 }
calcsim(10, keys)
keys = {'lp':5, 'ip':8, 'Th':7.5, 'Tc':10 }
calcsim(10, keys)
This is the code to simulate by changing the cluster countermeasure time $ T_c $ in various ways.
# case 7.2
def calcTctoRp(Rp0):
# Rp S E I R
ini_state=[Rp0, 0, 0, 1, 0]
t_max=180
dt=0.01
t=np.arange(0,t_max,dt)
#
lp = 5
ip = 8
keys = {'lp':lp, 'ip':ip, 'Th':7.5, 'Tc':(5+0) }
#
rslt = []
for i in range(0,60):
keys['Tc'] = lp + i
sim = my_odeint(seir_eq7, ini_state, t, keys)
rslt.append([keys['Tc'], sim[-1:,4]]) # (i, R(tmax))
#
rslt = np.array(rslt)
ymax = max(rslt[:,1])
plt.rcParams["font.size"] = 12
fig, ax = plt.subplots(figsize=(10,5))
#
ax.plot([0 , 0],[0,ymax],'m:')
ax.plot([lp,lp],[0,ymax],'b:')
ax.plot([lp+ip,lp+ip],[0,ymax],'g:')
#
ax.plot( rslt[:,0], rslt[:,1], 'r')
ax.legend([ '0 day', 'lp days','lp + ip days','Total infected'], loc='lower right')
ax.grid(which='both')
ax.set_xlabel('cluster shutdown date since cluster occured')
ax.set_ylabel('R(tmax)')
#
for tat in [10,13,20,30,40]:
idx = [i for i in range(len(rslt[:,0])) if rslt[i,0] >= tat][0]
print("R_fin with Tc:{} is {}".format(rslt[idx,0], rslt[idx,1]))
#
ax.set_xlim(-1,)
ax.set_ylim(0,)
plt.show()
calcTctoRp(100)
calcTctoRp(30)
calcTctoRp(20)
calcTctoRp(10)
Let's take a look at the calculation result. The result of changing some patterns of the initial reproduction number $ Rp (t_0) $ and changing the cluster countermeasure time $ T_c $ by fixing $ Rp (t_0) $ is shown.
Eventually 113 people will be infected.
Eventually 3493 people will be infected.
Eventually 58,437 people will be infected.
When the initial reproduction number Rp (t_0) = 10, the cluster countermeasure time Tc is changed and verified. Of course, the sooner the cluster measures are taken, the more effective it is expected. In the above simulation, 113 people were eventually infected without clustering measures.
If the infected person can be quarantined within 10 days of the cluster outbreak, the final infected person can be reduced from 113 to 21.
If the infected can be quarantined within 20 days of the cluster, the final number of infected can be reduced from 113 to 64.
Let's graph the relationship between the cluster countermeasure time Tc and the final number of infected people. Dotted lines are drawn at three points: when a cluster occurs, when the incubation period lp has elapsed, and when the incubation period + infection period (lp + ip) has elapsed.
Compared to when Rp (t_0) = 10, the effect of cluster measures is relatively increased.
Compared to when Rp (t_0) = 20, the effect of cluster measures is relatively increased.
From the above, the following trends can be derived from the simulation regarding the transition of the number of infected persons and the effect of the cluster countermeasure time based on the cluster-based SEIR model.
Even in the cluster-based SEIR model that does not assume the initial value of + S, the curve has almost the same shape as the normal SEIR model.
I referred to the following page.