python拟合sir模型_如何在python中实现自动模型确定和两种状态模型拟合？

Currently, I've been doing model fitting in Prism manually for all my data. It's quite tedious and time consuming. I wonder if there is any way to improve the efficiency in data analysis. I'm familiar with Python so I'd like to start with python to think of a better workflow. Your help is much appreciated.

Two questions:

How to do two state model fitting in python (figure 1)? In my case, it will has an initial linear increase followed by horizontal plateau state. I wish to the way of implementation in python and the way to automatically detect the turning point where the changes occurs (Ideally I can get the time where it occurs and the slope as well)

Another case is when the state become exponentially or polynomial. Is there any way in python to automatically tell which model is the best.

解决方案

Scipy provides a least square curve fit method that supports custom defined functions. Here is an example for the first model:

import numpy as np

from scipy.optimize import curve_fit

import matplotlib.pyplot as plt

#custom fit function - first slope steeper than second slope

def two_lin(x, m1, n1, m2, n2):

return np.min([m1 * x + n1, m2 * x + n2], axis = 0)

#x/y data points

x = np.asarray([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10])

y = np.asarray([2, 4, 8, 12, 14, 18, 20, 21, 22, 23, 24])

#initial guess for a steep rising and plateau phase

start_values = [3, 0, 0, 3]

#curve fitting

fit_param, pcov = curve_fit(two_lin, x, y, p0 = start_values)

#output of slope/intercept for both parts

m1, n1, m2, n2 = fit_param

print(m1, n1, m2, n2)

#calculating sum of squared residuals as parameter for fit quality

r = y - two_lin(x, *fit_param)

print(np.sum(np.square(r)))

#point, where the parts intersect

if m1 != m2:

x_intersect = (n2 - n1) / (m1 - m2)

print(x_intersect)

else:

print("did not find two linear components")

#plot data and fit function

x_fit = np.linspace(-1, 11, 100)

plt.plot(x, y, 'o', label='data')

plt.plot(x_fit, two_lin(x_fit, *fit_param), '--', label='fit')

plt.axis([-2, 12, 0, 30])

plt.legend()

plt.show()

More information about scipy.optimize.curve_fit can be found in the reference guide. For polynomials, numpy provides standard functions with numpy.polyfit and numpy.poly1d, but you still have to provide the expected degree.

The sum of squared residuals can be used to compare the accuracy of different fit functions.