[PYTHON] [Introduction à la simulation] J'ai essayé de jouer en simulant une infection corona ♬ Partie 2

Cette fois, il s'agit d'une suite de Simulation résumée le 18 mars. En fait, j'avais l'impression d'avoir terminé, mais quand j'y ai pensé, j'ai dû approfondir les deux points suivants. ① Transmission d'infection qui ne peut pas être expliquée par le modèle SIR ② Propagation lorsqu'elle se produit dans un cluster En d'autres termes, dans (1), en regardant la transition du nombre d'infections dans la région de transition dépendante de la densité de particules précédente, le pic du nombre d'infections est divisé en plusieurs pics, ce qui ne peut pas être expliqué par le modèle SIR. En outre, les infections réelles se produisent souvent en grappes et ne semblent pas être transmises de manière monotone. Donc, cette fois, nous avons implémenté ce qui suit. ・ Transmission de l'infection dans la région de densité critique ・ Transmission lorsque la distribution des porteurs de susceptibilité est répartie en grappes Le code est presque le même que le précédent, mais il a été étendu pour donner de la vitesse à la distribution et aux valeurs initiales des détenteurs de susceptibilité. Je vais mettre le code en bonus.

• Cette fois, je vais omettre la simulation pour donner de la vitesse.

・ Transmission de l'infection dans la région de densité critique

Lorsque les détenteurs de susceptibilité sont répartis aléatoirement et uniformément comme dans le temps précédent ① Gamme d'infection; Rc ② Probabilité d'infection; p ③ Densité des porteurs de susceptibilité; d ④ Vitesse de déplacement du support de susceptibilité; mr Peut être modifié pour contrôler la transmission et le taux d'infection.

Transmission typique

Le chiffre de droite est bleu; nombre de porteurs sensibles, rouge; infecté, vert; guérisseur, et ce graphique est une image familière du modèle SIR que vous voyez souvent. Paramètres: ①Rc = 20x20, ②p = 0,05, ③d = 1000 / 200x200, ④mr = 1.0

Transmission critique

Si la probabilité d'infection est abaissée comme suit, la transmission de l'infection ne se produira pas de sitôt. La transmission de l'infection juste avant peut être réalisée avec les paramètres suivants, et ce sera comme l'animation Gif ci-dessous.

Paramètres: ①Rc = 20x20, ②p = 0,0075, ③d = 1000 / 200x200, ④mr = 1.0

Comme on peut le voir sur cette figure, comme le montre la figure de droite, cette propagation est dans une région qui ne peut plus être représentée par le modèle SIR, et on peut voir que la propagation se produit de manière probabiliste. Sur cette figure, la hausse est exponentielle, mais après cela, l'infection se propage et guérit en tirant lentement la queue.

Paramètres: ①Rc = 13x13, ②p = 0,03, ③d = 1000 / 200x200, ④mr = 1.0

Une transmission critique de l'infection pourrait être obtenue même si la gamme d'infection était réduite, comme suit. Dans cette animation Gif, même au début, l'augmentation exponentielle disparaît, elle augmente lentement, diminue minable et augmente à nouveau. Ainsi, la transmission de l'infection dans cette région dépend de la distribution des porteurs sensibles et peut être considérée comme un processus presque probabiliste.

・ Transmission lors de la distribution en clusters

Dans la distribution réelle, on suppose que la densité de population fluctue généralement dans chaque village ou ville, ou que la distribution fluctue à chaque rassemblement ou rassemblement de personnes. Donc, j'ai essayé de propager l'infection lorsque la distribution des porteurs de sensibilité est distribuée dans un cluster.

Lorsque les grappes sont densément présentes

Cette fois, par souci de simplicité, nous nous sommes d'abord penchés sur la transmission en grappes réparties sur une grille 10x10.

Paramètres: ①Rc = 20x20, ②p = 0,5, ③d = 100x10 / 200x200, ④mr = 1,5

La signification de d = 100x10 / 100x100 signifie que 100 grappes sont placées dans la zone de 200x200 et 10 supports de sensibilité sont placés à chaque point de la grille comme indiqué dans la figure ci-dessous. Avec ce paramètre, la transmission de l'infection est une région qui peut être exprimée par le modèle SIR, qui est identique à la distribution uniforme.

Paramètres: ①Rc = 20x20, ②p = 0,5, ③d = 100x10 / 200x200, ④mr = 1,25

Lorsque la vitesse de la marche aléatoire du porteur de sensibilité diminue, la probabilité de contact avec la personne infectée diminue et la transmission entre les grappes se produit en fonction de l'infection par le porteur de la susceptibilité comme suit, et ce qui suit n'est pas possible. Une propagation régulière se produit.

Paramètres: ①Rc = 20x20, ②p = 0.5, ③d = 100x10 / 200x200, ④mr = 1.09

Et la vitesse autour de 1,09 est la limite de savoir si elle est critique et se propage. Et lorsque la vitesse est inférieure à cela, la propagation disparaît. En d'autres termes, il ne sera pas infecté. Avant l'infection, la transmission s'arrête à la première personne infectée ou à au moins l'un des groupes voisins.

Lorsque les grappes sont rares

Que faire si les grappes sont plus grandes et séparées? Autrement dit, l'image est lorsque de grands villages sont dispersés.

Paramètres: ①Rc = 20x20, ②p = 0,5, ③d = 4x250 / 200x200, ④mr = 10

La signification de d = 4x250 / 100x100 signifie que 200x200 est parsemé de grappes contenant chacune 250 porteuses de susceptibilité. Dans ce cas, il s'effondre un peu à cause de la distribution des grappes, mais il devient un seul pic pour le moment.

Paramètres: ①Rc = 20x20, ②p = 0,5, ③d = 4x250 / 200x200, ④mr = 7,5

Dans le cas de cette distribution de clusters, cette zone de paramètres donne une propagation critique au sens de propagation entre clusters.

Lorsque la grappe est clairsemée et que la personne infectée initiale est inégalement répartie sur le côté inférieur et que la distribution de la vitesse est inégale

Cela ressemble à un cadre compliqué lorsqu'il est écrit en lettres, mais cela semble possible en réalité.

Paramètres: ①Rc = 20x20, ②p = 0,5, ③d = 4x250 / 200x200, ④mr = 7,5, ⑤mx = 1, my = 0,5

Avec le dernier paramètre (5), la vitesse est augmentée dans la direction x et la vitesse dans la direction y est divisée par deux.

Paramètres: ①Rc = 20x20, ②p = 0,0075, ③d = 4x250 / 200x200, ④mr = 7,5, ⑤mx = 1, my = 0,75

Si vous placez la première personne infectée dans le village en bas à gauche depuis le début, ce sera comme suit. Et la probabilité d'infection est également fixée à une petite valeur de 0,0075 (la valeur de transmission critique de distribution uniforme). Même si la vitesse dans la direction y est de 0,75, on constate que les porteurs de susceptibilité appartenant au haut village sont à peine infectés. Dans ce domaine, on constate que la transition du nombre d'infections est très structurée et compliquée.

Discussion

Ce qui est important, c'est qu'il existe des zones où le nombre de belles infections n'augmente pas ou ne diminue pas comme on le voit dans le modèle SIR. Et en cas de propagation critique ① Si la vitesse de déplacement des personnes infectées et des porteurs de sensibilité est faible, l'infection ne sera pas transmise ② Si le taux d'infection est abaissé, la transmission de l'infection ne se produira pas. ③ Si la densité des porteurs de sensibilité est faible, il ne sera pas infecté. Par conséquent, il est dangereux de mélanger (voyager) des personnes infectées avec des porteurs sensibles. Le taux d'infection peut être réduit à 0 par un masque (personne infectée) ou un lavage des mains (porteur de sensibilité). Eh bien, certaines personnes sont asymptomatiques, donc tous les masques doivent être portés. De plus, pratiquez surtout des mesures pour éviter trois miels (surpeuplés, serrés, fermés), comme ne pas cueillir parce que l'infection peut être évitée simplement en abaissant la densité en dessous d'un certain niveau. Il semble que l'infection puisse être évitée en pratiquant les mesures qui se disent dans la rue. Enfin, bien que cela ne soit pas montré dans les calculs, je comprends que le concept de soi-disant immunité de masse n'a presque aucun sens car il change lorsque la densité et la probabilité d'infection sont simplement modifiées dans ces calculs.

Résumé

・ Le modèle SIR ne semble pas tenir dans la région de propagation critique ・ La distribution structurelle des porteurs de sensibilité entraîne des infections structurelles ・ Lorsque les grappes sont réparties dans les villages, il est peu probable que la transmission de l'infection entre les villages se produise. ・ En cas de propagation critique, elle ne se propage pas si la vitesse est faible. ・ Le concept d'immunité de masse n'a pas de sens (changements en fonction de la probabilité d'infection, de la plage d'infection, de la distribution des porteurs de sensibilité, etc.)

・ Si la personne infectée a une vitesse, le taux d'augmentation de l'infection augmentera naturellement, mais la simulation n'est pas terminée.

prime

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import random
import time

PARTICLE_NO = 1000 #Nombre de particules
ITERATION = 200 #Nombre maximum de boucles S'arrête lorsque le nombre de personnes infectées atteint 0
MIN_X, MIN_Y = -100.0, -100.0 #Portée minimale au début de la recherche
MAX_X, MAX_Y = 100.0, 100.0 #Portée maximale au début de la recherche
recovery=30 #Guéri après une certaine période de temps
p=0.5 #0.03 #probability of infecion
rc=100 #121 #169 #225 Gamme d'infection Rayon du cercle^2

start = time.time()

def plot_particle(sk,positions,elt,r,g,b):
    #fig, ax = plt.subplots()
    el_time = time.time()-start
    fig, (ax1, ax2) = plt.subplots(1, 2, sharey=False,figsize=(8*2, 8))
        
    for j in range(0,PARTICLE_NO):
        x=positions[j]["x"]
        y=positions[j]["y"]
        c=positions[j]["c"]
        s = 5**2
        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))
    
    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") #, bottom=b)
    rect3 = ax2.bar(ind+2*width, r,width, color="r") #, bottom=b)
    plt.pause(0.1)
    plt.savefig('./fig/fig{}_.png'.format(sk)) 
    plt.close()
    
#Fonction de mise à jour de la position des particules
def update_position(positions,velocity):
    x0 = []
    y0 = []
    for i in range(PARTICLE_NO):
        c=positions[i]["c"]
        t_time = positions[i]["t"]  #Valeur initiale 0, heure de l'infection au moment de l'infection
        k_time = time.time()-start  #temps écoulé
        s = positions[i]["flag"]    #Aucune infection 0, infection: 1
        if s == 1 and c == "red":   #Si infecté
            if k_time-t_time>recovery:  #Guéri après une certaine période de temps
                #print("inside",i,s,c,k_time-t_time)
                c = "blue"
                positions[i]["c"] = "green"
                positions[i]["flag"] = 1   #Cependant, les antécédents d'infection restent
        if c == "red":  #Obtenir des informations de localisation en cas d'infection rouge
            x0.append(positions[i]["x"])
            y0.append(positions[i]["y"])
            #print("4",i,s,c,t_time)
    #print(x0,y0)   
    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 < rc 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 = velocity[j]["x"]+1.085*random.uniform(-1, 1) #Le coefficient est la grandeur de la motilité des particules
        vy = velocity[j]["y"]+1.085*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})
        velocity.append({"x": vx, "y": vy})

    return position, velocity, x0

def count_brg(position):
    r=0
    g=0
    b=0
    for j in range(len(position)):
        if position[j]["c"] == "red":
            r += 1
        elif position[j]["c"] == "green":
            g += 1
        else:
            b += 1
    return r,g,b        

def main():
    #Mesure de l'heure de début
    #start = time.time()
    xy_min, xy_max = -32, 32
    #Position initiale de chaque particule,la vitesse, personal best,meilleurs paramètres globaux et espace de recherche
    position = []
    velocity = []  #Augmentation de la vitesse d'utilisation
    
    #Position initiale,Vitesse initiale
    #position.append({"x": random.uniform(MIN_X, MAX_X), "y": random.uniform(MIN_Y, MAX_Y), "c": "red", "t":0, "flag":1})
    position.append({"x": 0, "y": 0, "c": "red", "t":0, "flag":1}) #Placer une première personne infectée au milieu (0,0)
    velocity.append({"x": 0, "y": 0}) #La vitesse initiale de la personne infectée est de 0
    for i in range(0,10): #Disposer 10 villages dans la direction x
        for j in range(0,10): #Disposer les villages avec le maillage 10 dans la direction xy
            for k in range(0,10): #Répartition de 10 porteurs de susceptibilité par village
                s=k+j*10+i*100;
                position.append({"x": 10+(-100+i*20)+random.uniform(MIN_X/100, MAX_X/100), "y":10+(-100+j*20)+ random.uniform(MIN_Y/100, MAX_Y/100), "c": "blue", "t": 0, "flag":0})
                velocity.append({"x": 0, "y": 0})
    print(len(position))
    sk = 0
    red=[]
    green=[]
    blue=[]
    elapsed_time = []
    while sk < ITERATION:
        position,velocity, x0 = update_position(position,velocity) ######
        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)
        #print("{:.2f}:red_{} green_{} blue_{}".format(el_time,r,g,b))
        plot_particle(sk,position,elapsed_time,red,green,blue)
        if x0==[]:
            break
        sk += 1

    #Mesure du temps terminée
    process_time = time.time() - start
    print("time:", process_time)
    
if __name__ == '__main__':
    main()