[PYTHON] [Introduction to simulation] I tried playing by simulating corona infection ♬

In the field of numerical calculation, there have been simulations of free-moving particles and collective motion of Brownian-moving particles for a long time. This time, I would like to apply this to the simulation of coronavirus infection. There is PSO as a similar simulation, but this time I referred to these python codes. 【reference】 -About particle swarm optimization (PSO) parameters ・ Saved Particle Swarm Optimization (PSO)

What i did

・ Code explanation ・ Infection rate dependence ・ Particle density dependence ・ Particle motion dependence

・ Code explanation

The entire code is below. ・ Collive_particles / snow.py Explain the main parts of the code. The libraries used are as follows

import numpy as np
import matplotlib.pyplot as plt
import random
import time

The definition of the initial value is as follows

PARTICLE_NO = 1000 #Number of particles
ITERATION = 200 #Maximum number of loops: Stops when the number of infected people reaches 0
MIN_X, MIN_Y = -100.0, -100.0 #Minimum range at the start of search
MAX_X, MAX_Y = 100.0, 100.0 #Maximum range at the start of search
recovery=30 #Healed after a certain period of time
p=0.03 #probability of infecion

The graph display, as you can see from the above results, ax1: Particle position information ax2: Blue; Normal number, Red; Infection number, Green; Frequency graph of cure number Is displayed as follows. Arguments are graph number, position information, elapsed time, r, g, b values

def plot_particle(sk,positions,elt,r,g,b):
    el_time = time.time()-start
    fig, (ax1, ax2) = plt.subplots(2, 1, sharey=False,figsize=(8, 16))
        
    for j in range(0,PARTICLE_NO):
        x=positions[j]["x"]
        y=positions[j]["y"]
        c=positions[j]["c"]
        s = 5**2  #Particle size
        ax1.scatter(x, y, s, c, marker="o")
    ax1.set_xlim([MIN_X, MAX_X])
    ax1.set_ylim([MIN_Y, MAX_Y])
    ax1.set_xlabel("x")
    ax1.set_ylabel("y")
    ax1.set_title("{:.2f}:InfectionRate;{:.2f} %".format(el_time,(PARTICLE_NO-b[-1])/PARTICLE_NO*100)) #Cumulative infection rate
    ind = np.arange(len(elt))  # the x locations for the groups
    width = 0.3       # the width of the bars
    ax2.set_ylim([0, PARTICLE_NO])
    ax2.set_title("{:.2f}:red_{} green_{} blue_{}".format(el_time,r[-1],g[-1],b[-1]))
    rect1 = ax2.bar(ind, b,width, color="b")
    rect2 = ax2.bar(ind+width, g, width, color="g")
    rect3 = ax2.bar(ind+2*width, r,width, color="r")
    plt.pause(0.1)
    plt.savefig('./fig/fig{}_.png'.format(sk)) 
    plt.close()

The main logic this time was pushed into the particle position information update function. That is, the attributes of the particle object are defined as follows.

position.append({"x": new_x, "y": new_y, "c": p_color, "t": t_time,"flag":s})

Use the following function to find the time change of the above variables. Infection is judged by the following formula. It is evaluated whether or not it is infected with a certain probability p, assuming that it has come into contact with the inside of a circle with a radius of 20. In other words, ** this radius and the probability of infection p determine the strength of infection **.

if (x-x0[k])**2+(y-y0[k])**2 < 400 and random.uniform(0,1)<p:

Now that we have defined the particle object as above, we are extracting and assigning each variable. Another thing, for the sake of simplicity, the healing is ** automatically cured after a certain period of time **. You can cure it with probability here as well, but since it only complicates the event here, it does not need to be refined. In addition, ** particle movement (movement of people) was set to random walk ** on average. It is an approximation that people do not move so much. On the other hand, it is possible to make free movement like a simulation of molecular motion, but I have not done it. And ** the assumption that infection does not infect a cured person **.

#Particle position update function
def update_position(positions):
    x0 = []
    y0 = []
    for i in range(PARTICLE_NO):
        c=positions[i]["c"]
        t_time = positions[i]["t"]  #Initial value 0, infection time at infection
        k_time = time.time()-start  #elapsed time
        s = positions[i]["flag"]    #No infection 0, infection: 1
        if s == 1 and c == "red":   #If infected
            if k_time-t_time>recovery:  #Healed after a certain period of time
                c = "blue"
                positions[i]["c"] = "green"
                positions[i]["flag"] = 1   #However, the infection history remains
        if c == "red":  #Get location information if infected red
            x0.append(positions[i]["x"])
            y0.append(positions[i]["y"])
    position = []
    for j in range(PARTICLE_NO):
        x=positions[j]["x"]
        y=positions[j]["y"]
        c=positions[j]["c"]
        s = positions[j]["flag"]
        t_time = positions[j]["t"]
        for k in range(len(x0)):
            if (x-x0[k])**2+(y-y0[k])**2 < 400 and random.uniform(0,1)<p:
                if s ==0:
                    c = "red"
                    t_time = time.time()-start
                    s = 1
                    positions[j]["flag"]=s
                else:
                    continue
        vx = 1*random.uniform(-1, 1)
        vy = 1*random.uniform(-1, 1)
        new_x = x + vx
        new_y = y + vy
        p_color = c
        s=s
        
        position.append({"x": new_x, "y": new_y, "c": p_color, "t": t_time,"flag":s})
    return position, x0

The main function is as follows The point is that one particle is infected by default; red, flag; 1, etc. Other particles are randomly arranged without infection. There is an idea to place the position of the first particle in the center, but I wanted to see the transmission of infection from various positions, so I chose an arbitrary position. The count_brg () function simply counts.

def main():
    #Initial position of each particle,speed, 
    position = []
    velocity = []  #Speed is not used this time
    #Initial position,Initial speed
    position.append({"x": random.uniform(MIN_X, MAX_X), "y": random.uniform(MIN_Y, MAX_Y), "c": "red", "t":0, "flag":1})
    for s in range(1,PARTICLE_NO):
        position.append({"x": random.uniform(MIN_X, MAX_X), "y": random.uniform(MIN_Y, MAX_Y), "c": "blue", "t": 0, "flag":0})
    sk = 0    
    red=[]
    green=[]
    blue=[]
    elapsed_time = []
    while sk < ITERATION:
        position, x0 = update_position(position)
        r,g,b = count_brg(position)
        red.append(r)
        green.append(g)
        blue.append(b)
        el_time=time.time()-start
        elapsed_time.append(el_time)
        plot_particle(sk,position,elapsed_time,red,green,blue)
        if x0==[]:
            break
        sk += 1

・ Infection rate dependence

First, how does the state of infection transmission and the cumulative infection rate as a whole change when the infection rate p is changed? In the above example, the infection rate is p = 30%, the cumulative infection rate is 100%, and it can be seen that the infection has definitely occurred and is transmitted like a tsunami. The peak infections were 404 and 67.5 seconds. In this case, it seems that a group velocity such as infection transmission rate can be defined, but I somehow stopped it. On the other hand, when the infection rate p = 5%, the cumulative infection rate dropped to 92.20%, and the state of infection transmission became conspicuous in blue as shown below, the infection peak was suppressed to 235, and the peak position extended to 178 seconds. Furthermore, when the critical infection rate of whether or not infection transmission occurs = 3%, the results are as follows. In other words, the transmission of infections such as the tsunami disappears, and it is an image of terrifying transmission. The peak infection was as low as 117 and as long as 260 sec. Cumulative infection rate decreased to 52%. In the case of this infection rate, this calculation result was also completed in the middle, but it sometimes disappeared at the initial stage with almost no infection. In other words, ** If you can reduce the infection rate below a certain level (3% or less in this simulation) by using hand washing and masks well, it is possible to eliminate the transmission of infection. ** **

These results are considered to indicate that ** to some extent, control of communication methods such as masks, hugs, kisses and cough etiquette can prolong but reduce the cumulative infection rate **. it can. And if you look closely at this simulation, ** the blue particles (uninfected person) in the area where there is a person who has healed once are protected from infection **, and the transmission of infection does not spread inward. That is, ** populations in this area are robust with respect to infection and rarely avoid the risk of infection. ** That is, it can be said that ** herd immunity was acquired **. This herd immunity can be simply rephrased as ** preventing transmission of infection because healers substantially reduce the density of uninfected individuals in the population **.

・ Particle density dependence

Another interest is whether the common sense of ** not gathering **, which usually causes event cancellations or rally cancellations, is the correct answer. Here, under the conditions of an infection probability of 30% and a recovery time of 30 seconds, the particle density was changed to see the difference in behavior. The results are as follows when the particle density dependence of the infection rate is shown in the table.

That is, at a density of 60 or less, infection transmission does not occur, and infection transmission begins to occur from around 80 / 200X200, so that the cumulative infection rate begins to increase and reaches about 100% at around 200. As an example, the simulation when the transition region is 140 is shown below. The characteristic was a single peak when the above density was sufficiently high, but at this density, the peak of the infected person is gentler, and some peaks are reflected, reflecting the density fluctuation of the particle distribution. To show. And, the state of infection transmission is struggling to overcome the area with low density, and it can be seen that ** sometimes the infection cannot be transmitted if there is too much space **. In other words, the policy of not gathering not only has the effect of lowering the cumulative infection rate by lowering the overall density **, but also ** by acting with the policy of not gathering locally, the probability of infection is increased. It means that it can be lowered **.

・ Particle motion dependence

Next, let's verify whether ** not going out ** is correct. For the time being, we looked at how this effect affects infection transmission and cumulative infection rate when particle motion is large. When the particle density is 120, the results are shown in the table below.

From this assessment, there is almost no infection at this density when stationary. On the other hand, it can be seen that the cumulative infection rate increases when there is particle motion, and the cumulative infection rate tends to increase as the intensity increases. It can be imagined that this is simply overcoming the gap that hinders the transmission of infection by exercise as described above. A typical example of an actual simulation is as follows. As you can imagine from this simulation, you can see that the infection has spread over a slightly large space. In other words, it can be seen that ** not going out ** is also an act of lowering the risk of infection. This moving particle model can also be simulated by considering various situations that imitate the real world, but this time it is up to this point.

Summary

・ I tried playing with a simulation of corona infection ・ ** The conclusion is that if people do not gather and the infection rate can be reduced to below a certain level (3% or less in this simulation) by using hand washing and masks well, it is possible to prevent the transmission of infection by itself **. Obtained ・ It was found that if the infection rate is reduced by not gathering, masking, washing hands, etc., it takes time to end the infection, but the cumulative infection rate can be reduced. ・ As a verification of the effect of not going out, it was found that the cumulative infection rate increases when exercising. ・ At high density, transmission is similar to the spread of ripples that drop water droplets. ・ In the case of low density but transmission of infection, the area where healers and uninfected persons coexist expands, but this area is robust against infection and can be said to be a herd immunity state.

・ If this simulation is performed on a larger scale, it can be expanded to imitate an actual city, so I will try to calculate with a little higher firepower. By the way, I used Jetson-nano this time, but I could calculate up to 1000 without any problem.