我这里的代码执行得非常糟糕。当我在滑块上改变东西时,我的速度几乎不超过10 fps。虽然我对matplotlib不是很精通,但有人能指出我做错了什么以及如何修复它吗
注意:我正在处理大量数据,在最坏的情况下大约3*100000点。。。 也不确定这是否是必要的,但我是在'TkAgg'后端运行
这是我的代码(这是一个绘制和运行SIR流行病学数学模型的代码):
import matplotlib.pyplot as plt
from matplotlib.widgets import Slider, Button
import matplotlib.patches as patches
p = 1 #population
i = 0.01*p #infected
s = p-i #susceptible
r = 0 #recovered/removed
a = 3.2 #transmission parameter
b = 0.23 #recovery parameter
initialTime = 0
deltaTime = 0.001 #smaller the delta, better the approximation to a real derivative
maxTime = 10000 #more number of points, better is the curve generated
def sPrime(oldS, oldI, transmissionRate): #differential equations being expressed as functions to
return -1*((transmissionRate*oldS*oldI)/p) #calculate rate of change between time intervals of the
#different quantities i.e susceptible, infected and recovered/removed
def iPrime(oldS, oldI, transmissionRate, recoveryRate):
return (((transmissionRate*oldS)/p)-recoveryRate)*oldI
def rPrime(oldI, recoveryRate):
return recoveryRate*oldI
maxTimeInitial = maxTime
def genData(transRate, recovRate, maxT):
global a, b, maxTimeInitial
a = transRate
b = recovRate
maxTimeInitial = maxT
sInitial = s
iInitial = i
rInitial = r
time = []
sVals = []
iVals = []
rVals = []
for t in range(initialTime, maxTimeInitial+1): #generating the data through a loop
time.append(t)
sVals.append(sInitial)
iVals.append(iInitial)
rVals.append(rInitial)
newDeltas = (sPrime(sInitial, iInitial, transmissionRate=a), iPrime(sInitial, iInitial, transmissionRate=a, recoveryRate=b), rPrime(iInitial, recoveryRate=b))
sInitial += newDeltas[0]*deltaTime
iInitial += newDeltas[1]*deltaTime
rInitial += newDeltas[2]*deltaTime
if sInitial < 0 or iInitial < 0 or rInitial < 0: #as soon as any of these value become negative, the data generated becomes invalid
break #according to the SIR model, we assume all values of S, I and R are always positive.
return (time, sVals, iVals, rVals)
fig, ax = plt.subplots()
plt.subplots_adjust(bottom=0.4, top=0.94)
plt.title('SIR epidemiology curves for a disease')
plt.xlim(0, maxTime+1)
plt.ylim(0, p*1.4)
plt.xlabel('Time (t)')
plt.ylabel('Population (p)')
initialData = genData(a, b, maxTimeInitial)
susceptible, = ax.plot(initialData[0], initialData[1], label='Susceptible', color='b')
infected, = ax.plot(initialData[0], initialData[2], label='Infected', color='r')
recovered, = ax.plot(initialData[0], initialData[3], label='Recovered/Removed', color='g')
plt.legend()
transmissionAxes = plt.axes([0.125, 0.25, 0.775, 0.03], facecolor='white')
recoveryAxes = plt.axes([0.125, 0.2, 0.775, 0.03], facecolor='white')
timeAxes = plt.axes([0.125, 0.15, 0.775, 0.03], facecolor='white')
transmissionSlider = Slider(transmissionAxes, 'Transmission parameter', 0, 10, valinit=a, valstep=0.01)
recoverySlider = Slider(recoveryAxes, 'Recovery parameter', 0, 10, valinit=b, valstep=0.01)
timeSlider = Slider(timeAxes, 'Max time', 0, 100000, valinit=maxTime, valstep=1, valfmt="%i")
def updateTransmission(newVal):
newData = genData(newVal, b, maxTimeInitial)
susceptible.set_ydata(newData[1])
infected.set_ydata(newData[2])
recovered.set_ydata(newData[3])
r_o.set_text(r'$R_O$={:.2f}'.format(a/b))
fig.canvas.draw_idle()
def updateRecovery(newVal):
newData = genData(a, newVal, maxTimeInitial)
susceptible.set_ydata(newData[1])
infected.set_ydata(newData[2])
recovered.set_ydata(newData[3])
r_o.set_text(r'$R_O$={:.2f}'.format(a/b))
fig.canvas.draw_idle()
def updateMaxTime(newVal):
global susceptible, infected, recovered
newData = genData(a, b, int(newVal.item()))
del ax.lines[:3]
susceptible, = ax.plot(newData[0], newData[1], label='Susceptible', color='b')
infected, = ax.plot(newData[0], newData[2], label='Infected', color='r')
recovered, = ax.plot(newData[0], newData[3], label='Recovered/Removed', color='g')
transmissionSlider.on_changed(updateTransmission)
recoverySlider.on_changed(updateRecovery)
timeSlider.on_changed(updateMaxTime)
resetAxes = plt.axes([0.8, 0.025, 0.1, 0.05])
resetButton = Button(resetAxes, 'Reset', color='white')
r_o = plt.text(0.1, 1.5, r'$R_O$={:.2f}'.format(a/b), fontsize=12)
def reset(event):
transmissionSlider.reset()
recoverySlider.reset()
timeSlider.reset()
resetButton.on_clicked(reset)
plt.show()
使用适当的ODE解算器,如
scipy.integrate.odeint
来提高速度。然后可以对输出使用较大的时间步长。使用隐式解算器,如odeint
或solve_ivp
和method="Radau"
,精确解中作为边界的坐标平面也将作为数值解中的边界,因此值永远不会变为负值减少打印数据集以匹配打印图像的实际分辨率。 从300点到1000点的差异可能仍然是可见的,从1000点到5000点没有可见的差异,甚至可能没有实际的差异
matplotlib使用缓慢的python迭代,通过场景树将其图像绘制为对象。如果要绘制的对象超过10000个,则绘制速度会非常慢,因此最好将详细信息的数量限制为该数量
ODE解算器的代码
为了解决ODE,我使用了solve_ivp,但如果使用odeint,则没有什么区别
简化了绘图更新过程的代码
可以删除许多重复的代码。我还添加了一条线,以便时间轴随maxTime变量变化,这样就可以真正放大
中间和周围的代码保持不变
相关问题 更多 >
编程相关推荐