Python 中 scipy 的 curve_fit 曲线拟合
curve fitting by curve_fit from scipy in Python
我正在尝试用指数函数拟合我的数据
import numpy as np
def exponentional(k, alpha, k0, c):
return k0 * np.exp(k *-alpha) + c
我使用了 scipy.optimize、
中的 curve_fit
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
uniq_deg = [2,...,103,..,203,...,307,...,506]
normalized_deg_dist = [0.99,...,0.43,..0.12,..,0.04,..., 0.01]
popt, pcov = curve_fit(exponentional, uniq_deg, normalized_deg_dist,
p0 = [1,0.00001,1,1], maxfev = 6000)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.semilogy(uniq_deg, normalized_deg_dist, 'bo', label = 'Real data')
ax.semilogy(uniq_deg,[exponentional(d,*popt) for d in uniq_deg], 'r-', label = 'Fit')
ax.set_xlabel('Degree' )
ax.set_ylabel('1-CDF degree')
ax.legend(loc='best')
ax.set_title(f'Degree distribution in {city}')
plt.show()
导致:
看起来不太合适。
我哪里错了?
最后没用curve_fit
。我使用了 https://mathworld.wolfram.com/LeastSquaresFittingExponential.html
中指数拟合的定义
而且我还需要通过 power-law 拟合一些其他数据,我也这样做了。
#%%
def fit_powerlaw(xs, ys):
S_lnx_lny = 0.0
S_lnx_S_lny = 0.0
S_lny = 0.0
S_lnx = 0.0
S_lnx2 = 0.0
S_ln_x_2 = 0.0
n = len(xs)
for (x,y) in zip(xs, ys):
S_lnx += np.log(x)
S_lny += np.log(y)
S_lnx_lny += np.log(x) * np.log(y)
S_lnx_S_lny = S_lnx * S_lny
S_lnx2 += np.power(np.log(x),2)
S_ln_x_2 = np.power(S_lnx,2)
#end
b = (n * S_lnx_lny - S_lnx_S_lny ) / (n * S_lnx2 - S_ln_x_2)
a = (S_lny - b * S_lnx) / (n)
return (np.exp(a), b)
#%%
def fit_exp(xs, ys):
S_x2_y = 0.0
S_y_lny = 0.0
S_x_y = 0.0
S_x_y_lny = 0.0
S_y = 0.0
for (x,y) in zip(xs, ys):
S_x2_y += x * x * y
S_y_lny += y * np.log(y)
S_x_y += x * y
S_x_y_lny += x * y * np.log(y)
S_y += y
#end
a = (S_x2_y * S_y_lny - S_x_y * S_x_y_lny) / (S_y * S_x2_y - S_x_y * S_x_y)
b = (S_y * S_x_y_lny - S_x_y * S_y_lny) / (S_y * S_x2_y - S_x_y * S_x_y)
return (np.exp(a), b)
第一张图表用于指数拟合,第二张图表用于幂律。
我觉得这个结果很有说服力。
我正在尝试用指数函数拟合我的数据
import numpy as np
def exponentional(k, alpha, k0, c):
return k0 * np.exp(k *-alpha) + c
我使用了 scipy.optimize、
中的curve_fit
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
uniq_deg = [2,...,103,..,203,...,307,...,506]
normalized_deg_dist = [0.99,...,0.43,..0.12,..,0.04,..., 0.01]
popt, pcov = curve_fit(exponentional, uniq_deg, normalized_deg_dist,
p0 = [1,0.00001,1,1], maxfev = 6000)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.semilogy(uniq_deg, normalized_deg_dist, 'bo', label = 'Real data')
ax.semilogy(uniq_deg,[exponentional(d,*popt) for d in uniq_deg], 'r-', label = 'Fit')
ax.set_xlabel('Degree' )
ax.set_ylabel('1-CDF degree')
ax.legend(loc='best')
ax.set_title(f'Degree distribution in {city}')
plt.show()
导致:
看起来不太合适。
我哪里错了?
最后没用curve_fit
。我使用了 https://mathworld.wolfram.com/LeastSquaresFittingExponential.html
而且我还需要通过 power-law 拟合一些其他数据,我也这样做了。
#%%
def fit_powerlaw(xs, ys):
S_lnx_lny = 0.0
S_lnx_S_lny = 0.0
S_lny = 0.0
S_lnx = 0.0
S_lnx2 = 0.0
S_ln_x_2 = 0.0
n = len(xs)
for (x,y) in zip(xs, ys):
S_lnx += np.log(x)
S_lny += np.log(y)
S_lnx_lny += np.log(x) * np.log(y)
S_lnx_S_lny = S_lnx * S_lny
S_lnx2 += np.power(np.log(x),2)
S_ln_x_2 = np.power(S_lnx,2)
#end
b = (n * S_lnx_lny - S_lnx_S_lny ) / (n * S_lnx2 - S_ln_x_2)
a = (S_lny - b * S_lnx) / (n)
return (np.exp(a), b)
#%%
def fit_exp(xs, ys):
S_x2_y = 0.0
S_y_lny = 0.0
S_x_y = 0.0
S_x_y_lny = 0.0
S_y = 0.0
for (x,y) in zip(xs, ys):
S_x2_y += x * x * y
S_y_lny += y * np.log(y)
S_x_y += x * y
S_x_y_lny += x * y * np.log(y)
S_y += y
#end
a = (S_x2_y * S_y_lny - S_x_y * S_x_y_lny) / (S_y * S_x2_y - S_x_y * S_x_y)
b = (S_y * S_x_y_lny - S_x_y * S_y_lny) / (S_y * S_x2_y - S_x_y * S_x_y)
return (np.exp(a), b)
第一张图表用于指数拟合,第二张图表用于幂律。 我觉得这个结果很有说服力。