Use Python to simulate the spread of plague Gif images on a map.

Source: Internet
Author: User
Tags vmin

Use Python to simulate the spread of plague Gif images on a map.

Inspired by Jason's "Almost Looks Like Work", I will show some virus transmission models. Note that this model does not reflect the actual situation, so do not mistakenly think it is a terrible infectious disease in West Africa. On the contrary, it should be seen as a fictitious zombie outbreak. Let's enter the topic.

This is the SIR model. The letters S, I, and R indicate the differences between individuals in the botnet epidemic.

  • S represents the susceptible group, that is, the number of potential changes in healthy individuals.
  • I represents the number of infected groups, that is, the number of botnets.
  • R indicates the number of removed botnets that quit the game due to death, or the number of botnets that are transferred back to humans after infection. But there is no cure for botnets, so we should not be fooled by ourselves (if we want to apply the SIR model to the flu infection, there is still a cure ).
  • For beta and gamma ):
  • Beta indicates the degree to which the disease is contagious. Once bitten, It is infected.
  • Gamma (gamma) indicates the speed from botnets to death, depending on the average working rate of botnets. Of course, this is not a perfect model. Please be patient with me.
  • S' =? Beta IS tells us the rate at which healthy people turn into botnets, And S' IS the derivative of time.
  • I '= βis? Gamma I tells us how infected people increase, and how quickly the zombie enters the removal state (puntion ).
  • R' = gamma I is only added with (gamma I), which is negative in the preceding equation.

The Space Distribution of S/I/R is not taken into account in the above model. Let's fix it!

One way is to split the Swedish and Nordic countries into a grid. Each unit can infect Neighboring units, as described below:

For the unit, and for the four units around it. (Do not get tired because of diagonal units. We need our brains not to be eaten ).

Initialize some things.
 

import numpy as npimport mathimport matplotlib.pyplot as plt  %matplotlib inlinefrom matplotlib import rcParamsimport matplotlib.image as mpimgrcParams['font.family'] = 'serif'rcParams['font.size'] = 16rcParams['figure.figsize'] = 12, 8from PIL import Image

Appropriate beta and gamma values can destroy more than half of the country
 

beta = 0.010gamma = 1

Do you still remember the definition of derivative? When the derivative is known and assume that the △t is very small, it can be used to approximate the next value of the prediction function. We have already declared U' (t ).

Initialize some things.
 

import numpy as npimport mathimport matplotlib.pyplot as plt  %matplotlib inlinefrom matplotlib import rcParamsimport matplotlib.image as mpimgrcParams['font.family'] = 'serif'rcParams['font.size'] = 16rcParams['figure.figsize'] = 12, 8from PIL import Image

Appropriate beta and gamma values can destroy more than half of the country
 

beta = 0.010gamma = 1

Do you still remember the definition of derivative? When the derivative is known and assume that the △t is very small, it can be used to approximate the next value of the prediction function. We have already declared U' (t ).

The code for this method is as follows:
 

def euler_step(u, f, dt):  return u + dt * f(u)

We need the function f (u ). The friendly numpy provides simple array operations. I may review it in another article, because they are too powerful and require more explanation, but now they can achieve the effect:

 def f(u):  S = u[0]  I = u[1]  R = u[2]   new = np.array([-beta*(S[1:-1, 1:-1]*I[1:-1, 1:-1] +              S[0:-2, 1:-1]*I[0:-2, 1:-1] +              S[2:, 1:-1]*I[2:, 1:-1] +              S[1:-1, 0:-2]*I[1:-1, 0:-2] +              S[1:-1, 2:]*I[1:-1, 2:]),           beta*(S[1:-1, 1:-1]*I[1:-1, 1:-1] +              S[0:-2, 1:-1]*I[0:-2, 1:-1] +              S[2:, 1:-1]*I[2:, 1:-1] +              S[1:-1, 0:-2]*I[1:-1, 0:-2] +              S[1:-1, 2:]*I[1:-1, 2:]) - gamma*I[1:-1, 1:-1],           gamma*I[1:-1, 1:-1]          ])   padding = np.zeros_like(u)  padding[:,1:-1,1:-1] = new  padding[0][padding[0] < 0] = 0  padding[0][padding[0] > 255] = 255  padding[1][padding[1] < 0] = 0  padding[1][padding[1] > 255] = 255  padding[2][padding[2] < 0] = 0  padding[2][padding[2] > 255] = 255   return padding

Import the population density map of the Nordic countries and perform downsampling to get the results quickly.
 

from PIL import Imageimg = Image.open('popdens2.png')img = img.resize((img.size[0]/2,img.size[1]/2))img = 255 - np.asarray(img)imgplot = plt.imshow(img)imgplot.set_interpolation('nearest')

Population density map of Nordic countries (excluding Denmark)

The s matrix, that is, the susceptible individual, should be similar to the population density. The initial value of infected people is 0. We use Stockholm as the first source of infection.
 

S_0 = img[:,:,1]I_0 = np.zeros_like(S_0)I_0[309,170] = 1 # patient zero

Because no one has died, the matrix is also set to 0.
 

R_0 = np.zeros_like(S_0)

Then initialize the simulation duration.
 

T = 900             # final timedt = 1             # time incrementN = int(T/dt) + 1        # number of time-stepst = np.linspace(0.0, T, N)   # time discretization # initialize the array containing the solution for each time-stepu = np.empty((N, 3, S_0.shape[0], S_0.shape[1]))u[0][0] = S_0u[0][1] = I_0u[0][2] = R_0

We need to customize a color table to display the infection matrix on the map.
 

import matplotlib.cm as cmtheCM = cm.get_cmap("Reds")theCM._init()alphas = np.abs(np.linspace(0, 1, theCM.N))theCM._lut[:-3,-1] = alphas

Sit down and enjoy it...

 for n in range(N-1):  u[n+1] = euler_step(u[n], f, dt)

Let's make another rendering and make it a gif. Everyone loves gifs!
 

from images2gif import writeGif keyFrames = []frames = 60.0 for i in range(0, N-1, int(N/frames)):  imgplot = plt.imshow(img, vmin=0, vmax=255)  imgplot.set_interpolation("nearest")  imgplot = plt.imshow(u[i][1], vmin=0, cmap=theCM)  imgplot.set_interpolation("nearest")  filename = "outbreak" + str(i) + ".png"  plt.savefig(filename)  keyFrames.append(filename) images = [Image.open(fn) for fn in keyFrames]gifFilename = "outbreak.gif"writeGif(gifFilename, images, duration=0.3)plt.clf()

Contact Us

The content source of this page is from Internet, which doesn't represent Alibaba Cloud's opinion; products and services mentioned on that page don't have any relationship with Alibaba Cloud. If the content of the page makes you feel confusing, please write us an email, we will handle the problem within 5 days after receiving your email.

If you find any instances of plagiarism from the community, please send an email to: info-contact@alibabacloud.com and provide relevant evidence. A staff member will contact you within 5 working days.

A Free Trial That Lets You Build Big!

Start building with 50+ products and up to 12 months usage for Elastic Compute Service

  • Sales Support

    1 on 1 presale consultation

  • After-Sales Support

    24/7 Technical Support 6 Free Tickets per Quarter Faster Response

  • Alibaba Cloud offers highly flexible support services tailored to meet your exact needs.