用Python制作在地圖上模擬瘟疫擴散的Gif圖
更新時間：2015年03月31日 15:30:03 作者：Max Berggren
這篇文章主要介紹了如何用Python制作在地圖上模擬瘟疫擴散的Gif圖,其中用到了歐拉公式等數(shù)學知識、需要一定的算法基礎，需要的朋友可以參考下
受杰森的《Almost Looks Like Work》啟發(fā)，我來展示一些病毒傳播模型。需要注意的是這個模型并不反映現(xiàn)實情況，因此不要誤以為是西非可怕的傳染病。相反，它更應該被看做是某種虛構的僵尸爆發(fā)現(xiàn)象。那么，讓我們進入主題。
這就是SIR模型，其中字母S、I和R反映的是在僵尸疫情中，個體可能處于的不同狀態(tài)。
- S 代表易感群體，即健康個體中潛在的可能轉變的數(shù)量。
- I 代表染病群體，即僵尸數(shù)量。
- R 代表移除量，即因死亡而退出游戲的僵尸數(shù)量，或者感染后又轉回人類的數(shù)量。但對與僵尸不存在治愈者，所以我們就不要自我愚弄了（如果要把SIR模型應用到流感傳染中，還是有治愈者的）。
- 至于β(beta)和γ(gamma):
- β(beta)表示疾病的傳染性程度，只要被咬就會感染。
- γ(gamma)表示從僵尸走向死亡的速率，取決于僵尸獵人的平均工作速率，當然，這不是一個完美的模型，請對我保持耐心。
- S′=?βIS告訴我們健康者變成僵尸的速率，S′是對時間的導數(shù)。
- I′=βIS?γI告訴我們感染者是如何增加的，以及行尸進入移除態(tài)速率（雙關語）。
- R′=γI只是加上(gamma I)，這一項在前面的等式中是負的。
上面的模型沒有考慮S/I/R的空間分布，下面來修正一下！
一種方法是把瑞典和北歐國家分割成網格，每個單元可以感染鄰近單元，描述如下：
其中對于單元，和是它周圍的四個單元。（不要因為對角單元而腦疲勞，我們需要我們的大腦不被吃掉）。
初始化一些東東。
import numpy as np import math import matplotlib.pyplot as plt %matplotlib inline from matplotlib import rcParams import matplotlib.image as mpimg rcParams['font.family'] = 'serif' rcParams['font.size'] = 16 rcParams['figure.figsize'] = 12, 8 from PIL import Image
適當?shù)腷eta和gamma值就能夠摧毀大半江山
beta = 0.010 gamma = 1
還記得導數(shù)的定義么？當導數(shù)已知，假設Δt很小的情況下，經過重新整理，它可以用來近似預測函數(shù)的下一個取值，我們已經聲明過u′(t)。
初始化一些東東。
import numpy as np import math import matplotlib.pyplot as plt %matplotlib inline from matplotlib import rcParams import matplotlib.image as mpimg rcParams['font.family'] = 'serif' rcParams['font.size'] = 16 rcParams['figure.figsize'] = 12, 8 from PIL import Image
適當?shù)腷eta和gamma值就能夠摧毀大半江山
beta = 0.010 gamma = 1
還記得導數(shù)的定義么？當導數(shù)已知，假設Δt很小的情況下，經過重新整理，它可以用來近似預測函數(shù)的下一個取值，我們已經聲明過u′(t)。
這種方法叫做歐拉法，代碼如下：
def euler_step(u, f, dt): return u + dt * f(u)
我們需要函數(shù)f(u)。友好的numpy提供了簡潔的數(shù)組操作。我可能會在另一篇文章中回顧它，因為它們太強大了，需要更多的解釋，但現(xiàn)在這樣就能達到效果：
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
導入北歐國家的人口密度圖并進行下采樣，以便較快地得到結果
from PIL import Image
img = 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')
北歐國家的人口密度圖（未包含丹麥）
S矩陣，也就是易感個體，應該近似于人口密度。感染者初始值是0，我們把斯德哥爾摩作為第一感染源。
S_0 = img[:,:,1] I_0 = np.zeros_like(S_0) I_0[309,170] = 1 # patient zero
因為還沒人死亡，所以把矩陣也置為0.
R_0 = np.zeros_like(S_0)
接著初始化模擬時長等。
T = 900 # final time dt = 1 # time increment N = int(T/dt) + 1 # number of time-steps t = np.linspace(0.0, T, N) # time discretization # initialize the array containing the solution for each time-step u = np.empty((N, 3, S_0.shape[0], S_0.shape[1])) u[0][0] = S_0 u[0][1] = I_0 u[0][2] = R_0
我們需要自定義一個顏色表，這樣才能將感染矩陣顯示在地圖上。
import matplotlib.cm as cm
theCM = cm.get_cmap("Reds")
theCM._init()
alphas = np.abs(np.linspace(0, 1, theCM.N))
theCM._lut[:-3,-1] = alphas
下面坐下來欣賞吧…
for n in range(N-1): u[n+1] = euler_step(u[n], f, dt)
讓我們再做一下圖像渲染，把它做成gif，每個人都喜歡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()
相關文章
這篇文章主要介紹了python如何實現(xiàn)騰訊滑塊驗證碼識別，幫助大家更好的理解和學習使用python，感興趣的朋友可以了解下2021-04-04
RSA加密算法是一種強大的公鑰加密算法,安全性很高,這里我們來看一下Python使用Pycrypto庫進行RSA加密的方法詳解,需要的朋友可以參考下2016-06-06
OpenCV是計算機視覺領域一個圖像和視頻處理庫,用于各種圖像和視頻分析,如面部識別和檢測,車牌閱讀,照片編輯,高級機器人視覺,光學字符識別等等,下面這篇文章主要給大家介紹了關于openCV入門學習基礎教程第一篇的相關資料,需要的朋友可以參考下2022-11-11
python django框架中使用FastDFS分布式文件系統(tǒng)的安裝方法
這篇文章主要介紹了python-django框架中使用FastDFS分布式文件系統(tǒng)的安裝方法，本文圖文并茂給大家介紹的非常詳細，具有一定的參考借鑒價值 ,需要的朋友可以參考下2019-06-06
pycharm無法安裝第三方庫的問題及解決方法以scrapy為例(圖解)
這篇文章主要介紹了pycharm無法安裝第三方庫的解決辦法以scrapy為例,本文通過圖文并茂的形式給大家介紹的非常詳細，對大家的學習或工作具有一定的參考借鑒價值，需要的朋友可以參考下2020-05-05
這篇文章主要介紹了python 讀取二進制 顯示圖片案例，具有很好的參考價值，希望對大家有所幫助。一起跟隨小編過來看看吧2020-04-04
這篇文章主要介紹了python?time模塊計算時間之間的差距,本文通過實例代碼給大家介紹的非常詳細，對大家的學習或工作具有一定的參考借鑒價值，需要的朋友可以參考下2023-05-05
今天小編就為大家分享一篇對pytorch中的梯度更新方法詳解，具有很好的參考價值，希望對大家有所幫助。一起跟隨小編過來看看吧2019-08-08
這篇文章主要為大家介紹了Python Flask框架模塊安裝級使用介紹，有需要的朋友可以借鑒參考下，希望能夠有所幫助，祝大家多多進步，早日升職加薪2023-03-03
python生成器generator:深度學習讀取batch圖片的操作
這篇文章主要介紹了python生成器generator:深度學習讀取batch圖片的操作，具有很好的參考價值，希望對大家有所幫助。如有錯誤或未考慮完全的地方，望不吝賜教2021-05-05