Bonjour, je voudrais vous partager un morceau de code, pour ceux qui s’intéressent à la simulation d'épidémies.

C'est pas le plus beau code du monde puisque j'ai appris à coder sur une TI 83+, mais ça fonctionne correctement, et je pense que vous pourriez facilement l'optimiser et l'améliorer.

N'hésitez donc pas à modifier, améliorer ou juste vous amusez avec ce programme.

Je peux répondre aux questions sur le code si il y a des zones d'ombre.

N'oubliez pas d'installer matplotlib (pip install matplotlib) si ce n'est pas déjà fait.

Bien à vous, Pierre.



Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
 
#from random import *
import matplotlib.pyplot as plt
from tkinter import *
from tkinter import messagebox
from matplotlib.lines import Line2D
 
 
 
def start_simulation(_): #Function that will run the simulation after the start button or enter is pressed.
    font = {'size': 8}
    plt.rc('font', **font)
    plt.figure(1, figsize=(10, 6))
    plt.xlabel('Time (in days)')
    plt.ylabel('Population')
 
#Setting variables :
    simulation_length = int(simulation_length_field.get())
    increment_1=1
    infected_total = 0
    immunized = 0
    desimmunized = 0
    death_total = 0
    icu_in_the_same_time = 0
    icu_total = 0
    contagiousness = float(contagiousness_field.get()) / 100
    number_of_contact_per_day = float(contact_per_day_field.get())
    healthy = int(population_size_field.get())
    infected = int(initially_infected_field.get())
    time_before_immunization = int(disease_duration_field.get())
    death_rate = float(death_rate_field.get()) / 100
    icu_rate = float(icu_rate_field.get()) / 100
    mean_icu_stay_length = int(icu_duration_field.get())
    table_immunization = simulation_length * [0]
    table_icu = time_before_immunization * [0]
    death_multiplier = int(death_multiplier_v.get())
    icu_multiplier = int(icu_multiplier_v.get())
    display_healthy = int(display_healthy_v.get())
    immunization_length = int(immunisation_duration_field.get())
    display_death = int(display_death_v.get())
    display_icu = int(diplay_icu_v.get())
    dynamic_graph = int(dynamic_display_v.get())
    isolate_infected = int(infected_detection_v.get())
    time_before_detection = int(sick_days_before_detection_v.get())
    table_isolated = time_before_detection * [0]
    isolation_rate = int(proportion_isolated_v.get()) / 100
    isolated_total = 0
    healthy = healthy - infected
    infected_initially=infected
    healthy_initially = healthy
    table_infected = [infected]
    growth_rate = 0
    death_before_reducing_contact=int(death_before_quarantine_v.get())
    reduce_contact = int(quarantine_v.get())
 
    # find and fix errors before starting graph :
    if time_before_detection > time_before_immunization and isolate_infected == 1:
        time_before_detection = time_before_immunization
        messagebox.showinfo('Error', 'Delay before infected detection must be lower than the disease duration. Detection time has been set to ' + str(time_before_immunization) + '.')
 
    # enable static view if ticked :
    if dynamic_graph==0 :
        plt.axis([0, simulation_length, 0, healthy])
 
    # enable/disable death and/or ICU multiplier, to have a better overview when plotted :
    if death_multiplier == 1:
        death_multiply = 10
    else:
        death_multiply = 1
 
    if icu_multiplier ==1:
        icu_multiply = 10
    else:
        icu_multiply = 1
 
    # allow line 132 to run properly if isolation of infected is disabled :
    if isolate_infected == 0:
        table_isolated = simulation_length * [0]
 
    # Start of the main loop, which calculate data for each day :
    for actual_day in range(0,simulation_length):
        # infected_per_day=0  #reset the count of infected for each new day (only used if slow calculation is used)
        increment_1 = increment_1+1  # allows line 182 to 187 to run properly, avoid unwanted early epidemic termination.
 
        # reduces the number of contact per day to 1, when the total number of death is superior to the chosen limit:
        if death_total > death_before_reducing_contact and reduce_contact==1:
            number_of_contact_per_day = 1
 
        # This secondary loop is deactivated by default, it could be interesting if you plan on doing simulation of small populations (<50 000), to have a unique simulation at each run
        # It calculate each day and for each subject if they are infected or not by running a random function
        # But it is way too slow for big populations
        # If you decide to activate it, delete the triple quote at the beginning and at the end of the paragraph and add some to the next paragraph, also remove the # before infected_per_day=0 (line 81) and before random import at the beginning.
        '''for subject in range(0,healthy):
            probability_of_being_contaminated = number_of_contact_per_day*contagiousness*infected/(healthy+infected+immunized)
            if random() <= probability_of_being_contaminated:
               infected_per_day = infected_per_day + 1
               healthy = healthy -1'''
 
        # Add triple quote here if needed (cf line 91)
        # Calculation of the daily new infected :
        probability_of_being_contaminated = number_of_contact_per_day * contagiousness * infected / (
                healthy + infected + immunized)
        infected_per_day =  int((healthy * probability_of_being_contaminated))
        healthy = healthy - int((healthy * probability_of_being_contaminated)) # Add triple quote here if needed (cf line 91)
 
 
        table_infected.append(infected_per_day) # Saving the infected of the day in a table in order to remove them after healing time has past
 
        # Increment the infected counters
        infected=infected + infected_per_day
        infected_total = infected_total + infected_per_day
 
        # Growth rate calculation & displaying title + growth rate on the graph :
        if table_infected[actual_day-1] != 0:
            growth_rate=table_infected[actual_day]/table_infected[actual_day-1]
        graph_title = 'Epidemic Simulation,  Growth factor = ' + str(round(growth_rate,2))
        plt.suptitle(graph_title)
 
        # Detection and isolation of the infected :
        if isolate_infected == 1 and actual_day >= time_before_detection:
            infected_isolated = int(isolation_rate * table_infected[actual_day - time_before_detection])
            table_infected[actual_day - time_before_detection] = table_infected[actual_day - time_before_detection] - infected_isolated
            infected = infected - infected_isolated
            table_isolated.append(infected_isolated)
            isolated_total = isolated_total + infected_isolated
 
        plt.pause(0.05) # For live display
 
        # Infected become immunized after the time_before_immunization, or go to ICU or die :
        if actual_day >= time_before_immunization:
            infected = infected - table_infected[(actual_day - time_before_immunization)]
            daily_death = death_rate * (table_infected[actual_day - time_before_immunization])
            daily_icu = icu_rate*(table_infected[actual_day - time_before_immunization])
            immunized = immunized + table_infected[(actual_day - time_before_immunization)] - daily_death - daily_icu
 
            #same job for the infected in quarantine (become immunized, go to ICU or die) :
            if isolate_infected == 1 and actual_day >= time_before_detection:
                isolated_icu = icu_rate * (table_isolated[actual_day - time_before_immunization])
                daily_icu = daily_icu +isolated_icu
                isolated_death = death_rate * (table_isolated[actual_day - time_before_immunization])
                daily_death = daily_death + isolated_death
                immunized = immunized + table_isolated[actual_day - time_before_immunization] - isolated_death - isolated_icu
                isolated_total = isolated_total - table_isolated[actual_day - time_before_immunization]
 
            # If the immunization is not permanent, immunized subject get back to the healthy state and can be infected again :
            if immunization_length !=0:
                table_immunization[actual_day] = (table_infected[(actual_day - time_before_immunization)] - daily_death - daily_icu) + (table_isolated[actual_day - time_before_immunization])
                desimmunized = round(table_immunization[actual_day - immunization_length])
            if actual_day >= immunization_length != 0:
                healthy = healthy + desimmunized
                immunized = immunized - desimmunized
 
            # ICU count, adding daily_ICU to ICU_table in order to be able to retrieve and switch ICU subject from ICU to Immunized, Death count :
            icu_in_the_same_time = icu_in_the_same_time + daily_icu
            icu_total = icu_total + daily_icu
            table_icu.append(daily_icu)
            death_total = death_total + daily_death
 
            # When the mean ICU stay length is gone, subject from ICU switch to Immunized subjects :
            if actual_day >= mean_icu_stay_length:
                icu_in_the_same_time = icu_in_the_same_time - table_icu[actual_day - mean_icu_stay_length]
                immunized = immunized + table_icu[actual_day - mean_icu_stay_length]
 
            #Display of the legend on the graph :
            infected_patch = Line2D([0], [0], marker='o', color='w', label='Actual infected : ' + str(round(infected)), markerfacecolor='red', markersize=7)
            infected_per_day_patch = Line2D([0], [0], marker='^', color='w', label='Infected per day : ' + str(round(infected_per_day)), markerfacecolor='red', markersize=7)
            infected_total_patch = Line2D([0], [0], marker='*', color='w',label='Inf total : ' + str(round(infected_total)) + ' (' + str(round(infected_total / (infected_initially + healthy_initially) * 100, 1)) + ' %)', markerfacecolor='r', markersize=11)
            immunized_patch = Line2D([0], [0], marker='o', color='w', label='Immunized : ' + str(round(immunized)) + ' (' + str(round(immunized / (infected_initially + healthy_initially) * 100, 1)) + ' %)',markerfacecolor='green', markersize=7)
            death_patch = Line2D([0], [0], marker='o', color='w', label='Actual deaths : ' + str(round(daily_death)), markerfacecolor='black', markersize=7)
            death_total_patch = Line2D([0], [0], marker='*', color='w', label='Deaths total : ' + str(round(death_total)) + ' (' + str(round(death_total / (infected_initially + healthy_initially) * 100, 1)) + ' %)', markerfacecolor='black', markersize=11)
            icu_patch = Line2D([0], [0], marker='o', color='w', label='Actual ICU : ' + str(round(daily_icu)), markerfacecolor='blue', markersize=7)
            icu_total_patch = Line2D([0], [0], marker='*', color='w', label='ICU total : ' + str(round(icu_total)) + ' (' + str(round(icu_total / (infected_initially + healthy_initially) * 100, 1)) + ' %)', markerfacecolor='blue', markersize=11)
            healthy_patch = Line2D([0], [0], marker='o', color='w', label='Healthy not immunized', markerfacecolor='purple', markersize=7)
            isolated_patch = Line2D([0], [0], marker='o', color='w', label='Actual isolated : ' + str(round(isolated_total)), markerfacecolor='orange', markersize=7)
            info_prct_patch = Line2D([0], [0], marker='o', color='w', label='(% are on total population)', markerfacecolor='white',markersize=7)
 
            #Display the healthy subject on the legend if ticked, else not :
            if display_healthy == 1:
                plt.legend(handles=[infected_patch, infected_per_day_patch, infected_total_patch, immunized_patch, icu_patch, icu_total_patch, death_patch, death_total_patch, healthy_patch, isolated_patch,info_prct_patch ])
            else:
                plt.legend(handles=[infected_patch, infected_per_day_patch, infected_total_patch, immunized_patch, icu_patch, icu_total_patch,death_patch, death_total_patch, isolated_patch,info_prct_patch])
 
 
        else: # If the time of immunization/death/ICU (time_before_immunization) is not yet passed, this legend is shown. It's because the calculations of immunized, death and ICU can start only after this delay.
            infected_patch = Line2D([0], [0], marker='o', color='w', label='Inf : ' + str(round(infected-infected_initially)), markerfacecolor='r', markersize=7)
            infected_per_day_patch = Line2D([0], [0], marker='^', color='w', label='Inf per day : ' + str(round(infected_per_day)), markerfacecolor='r', markersize=7)
            infected_total_patch = Line2D([0], [0], marker='*', color='w', label='Inf total : ' + str(round(infected_total)) +' ---> ' + str(round(infected_total / (infected_initially + healthy_initially)*100, 1)) + ' %', markerfacecolor='r', markersize=11)
            plt.legend(handles=[infected_patch, infected_per_day_patch, infected_total_patch])
 
        # If the epidemic is gone, there's no more cases, no more virus to spread, a pop-up window appears :
        if infected == 0 and increment_1 > 10:
            result_end_question = messagebox.askquestion('Epidemic is gone', 'Do you want to quit all?')
            if result_end_question == 'yes':
               quit()
            else:
               break
 
        # Plot the different curves on the graph, sometimes with conditions (if the checkbox is ticked or not..)
        plt.scatter(actual_day, infected, color='red', s=1)
        plt.scatter(actual_day, immunized, color='green', s=1)
        if display_death == 1:
            plt.scatter(actual_day, death_total * death_multiply, color='black', s=1)
        if display_icu == 1:
            plt.scatter(actual_day, icu_in_the_same_time * icu_multiply, color='blue', s=1)
        if isolate_infected == 1:
            plt.scatter(actual_day, isolated_total, color='orange', s=1)
        if display_healthy == 1:
            plt.scatter(actual_day, healthy, color='purple', s=1)
 
    plt.show()
 
 
# This is the first window of the program :
window = Tk()
window.title('EPIDEMIC SIMULATION')
 
label_1 = Label()
label_1.grid(column=0, rowspan=4) #this is just a space...
 
#Every left sliders are declared here :
contagiousness_v = StringVar()
contagiousness_v.set(5)
contagiousness_field = Scale(window, variable = contagiousness_v, orient='horizontal', from_=0, to=20,
                             resolution=0.1, tickinterval=0, length=350,
                             label='% Contagiousness')
contagiousness_field.grid(column=0, rowspan=4, padx=20)
 
contact_per_day_v = StringVar()
contact_per_day_v.set(4)
contact_per_day_field = Scale(window, variable = contact_per_day_v, orient='horizontal', from_=0, to=30,
                              resolution=1, tickinterval=0, length=350,
                              label='N contact per day')
contact_per_day_field.grid(column=0, rowspan=4)
 
population_size_v = StringVar()
population_size_v.set(1000000)
population_size_field = Scale(window, variable = population_size_v, orient='horizontal', from_=100000, to=10000000,
                              resolution=100000, tickinterval=0, length=800,
                              label='Population size')
population_size_field.grid(column=0, rowspan=4, columnspan=4, sticky='W', padx=20)
 
initially_infected_v = StringVar()
initially_infected_v.set(500)
initially_infected_field = Scale(window, variable = initially_infected_v, orient='horizontal', from_=1, to=10000,
                                 resolution=1, tickinterval=0, length=800, sliderlength=15,
                                 label='N initially infected')
initially_infected_field.grid(column=0, rowspan=4, columnspan=4, sticky='W', padx=20)
 
death_rate_v = StringVar()
death_rate_v.set(1)
death_rate_field = Scale(window, variable = death_rate_v, orient='horizontal', from_=0, to=100,
                         resolution=0.5, tickinterval=0, length=350,
                         label='Death rate in %, when infected')
death_rate_field.grid(column=0, rowspan=4)
 
icu_rate_v = StringVar()
icu_rate_v.set(4)
icu_rate_field = Scale(window, variable = icu_rate_v, orient='horizontal', from_=0, to=100,
                       resolution=0.5, tickinterval=0, length=350,
                       label='ICU rate in %, when infected')
icu_rate_field.grid(column=0, rowspan=4)
 
icu_duration_v = StringVar()
icu_duration_v.set(12)
icu_duration_field = Scale(window, variable = icu_duration_v, orient='horizontal', from_=1, to=20,
                           resolution=1, tickinterval=0, length=350,
                           label='ICU mean duration of stay, in days')
icu_duration_field.grid(column=0, rowspan=4)
 
simulation_length_v = StringVar()
simulation_length_v.set(100)
simulation_length_field = Scale(window, variable = simulation_length_v, orient='horizontal', from_=50, to=1000,
                                resolution=50, tickinterval=0, length=350,
                                label='Simulation length, in days')
simulation_length_field.grid(column=0, rowspan=4)
 
disease_duration_v = StringVar()
disease_duration_v.set(10)
disease_duration_field = Scale(window, variable = disease_duration_v, orient='horizontal', from_=1, to=30,
                               resolution=1, tickinterval=0, length=350,
                               label='Disease duration, in days (before dying or being immunized)')
disease_duration_field.grid(column=0, rowspan=4)
 
immunization_duration_v = StringVar()
immunization_duration_v.set(0)
immunisation_duration_field = Scale(window, variable = immunization_duration_v, orient='horizontal', from_=0, to=200,
                                    resolution=5, tickinterval=0, length=350,
                                    label='Immunization duration, in days (0 = forever)')
immunisation_duration_field.grid (column=0, rowspan=4)
 
label_2 = Label()
label_2.grid(column=0, rowspan=4) # Just a space..
 
# Every checkbox are declared here + Start button :
infected_detection_v = IntVar()
infected_detection_v.set(0)
infected_detection_checkbox = Checkbutton(window, text='Detect INF after ', variable=infected_detection_v)
infected_detection_checkbox.grid(column=2, row=7, sticky='W')
 
sick_days_before_detection_v = StringVar()
sick_days_before_detection_v.set(5)
sick_days_before_detection_field = Scale(window, variable = sick_days_before_detection_v, orient='horizontal',
                                         from_=1, to=30, resolution=1, tickinterval=0, length=35,
                                         sliderlength=5, width=3)
sick_days_before_detection_field.grid(column=2, row=7)
 
text_2 = Label(window, text ='days sick and isolate')
text_2.grid(column=2, row=7, sticky='E')
 
proportion_isolated_v=StringVar()
proportion_isolated_v.set(50)
proportion_isolated_field = Scale(window, variable = proportion_isolated_v, orient='horizontal', from_=1, to=100,
                                  resolution=5, tickinterval=0, length=60, sliderlength=10, width=3)
proportion_isolated_field.grid(column=3, row=7, sticky='W')
 
text_3 = Label(window, text ='% of them')
text_3.grid(column=3, row=7)
 
quarantine_v = IntVar()
quarantine_v.set(0)
quarantine_checkbox = Checkbutton(window, text='Reduce to 1 contact per day when reaching ', variable=quarantine_v)
quarantine_checkbox.grid(column=2, row=11, sticky='W')
 
death_before_quarantine_v=StringVar()
death_before_quarantine_v.set(100)
death_before_quarantine_field = Scale(window, variable = death_before_quarantine_v, orient='horizontal',
                                      from_=1, to=5000, resolution=50, tickinterval=0, length=120,
                                      sliderlength=10, width=3)
death_before_quarantine_field.grid(column=3, row=11, sticky='W')
 
text_1 = Label(window, text ='deaths')
text_1.grid(column=3, row=11, sticky='E', padx=20)
 
display_death_v = IntVar()
display_death_v.set(1)
display_death_checkbox = Checkbutton(window, text='Display Deaths', variable=display_death_v)
display_death_checkbox.grid(column=2, row=23, sticky='W')
 
death_multiplier_v = IntVar()
death_multiplier_v.set(1)
death_multiplier_checkbox = Checkbutton(window, text='Display Deaths x10', variable=death_multiplier_v)
death_multiplier_checkbox.grid(column=2, row=23, sticky='E')
 
diplay_icu_v = IntVar()
diplay_icu_v.set(1)
display_icu_checkbox = Checkbutton(window, text='Display ICU', variable=diplay_icu_v)
display_icu_checkbox.grid(column=2, row=27, sticky='W')
 
icu_multiplier_v = IntVar()
icu_multiplier_v.set(1)
icu_multiplier_checkbox = Checkbutton(window, text='Display ICU x10', variable=icu_multiplier_v)
icu_multiplier_checkbox.grid(column=2, row=27, sticky='E')
 
display_healthy_v = IntVar()
display_healthy_v.set(1)
display_healthy_checkbox = Checkbutton(window, text='Display Healthy persons', variable=display_healthy_v)
display_healthy_checkbox.grid(column=3, row = 39, sticky='W', padx=20)
 
dynamic_display_v = IntVar()
dynamic_display_v.set(0)
dynamic_display_checkbox = Checkbutton(window, text='Dynamic display', variable=dynamic_display_v)
dynamic_display_checkbox.grid(column=3, row=40, sticky='W', padx=20)
 
start_button = Button(window, text='Start simulation')
start_button.grid(column=3, row=43)
start_button.bind('<Button-1>', start_simulation)
window.bind('<Return>', start_simulation)
 
text_4 = Label(window, text ='@pirog, 2020')
text_4.grid(column=2, row=44)
 
window.mainloop()