将曲线拟合到 curve_fit 的幂律分布不起作用
Fitting a curve to a power-law distribution with curve_fit does not work
我试图找到一条曲线拟合我的数据,该曲线在视觉上似乎具有幂律分布。
我希望利用 scipy.optimize.curve_fit,但无论我尝试什么函数或数据规范化,我都会收到 RuntimeError(未找到参数或溢出)或一条与我的数据不相符的曲线,甚至远程.请帮我弄清楚我在这里做错了什么。
%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
df = pd.DataFrame({
'x': [ 1000, 3250, 5500, 10000, 32500, 55000, 77500, 100000, 200000 ],
'y': [ 1100, 500, 288, 200, 113, 67, 52, 44, 5 ]
})
df.plot(x='x', y='y', kind='line', style='--ro', figsize=(10, 5))
def func_powerlaw(x, m, c, c0):
return c0 + x**m * c
target_func = func_powerlaw
X = df['x']
y = df['y']
popt, pcov = curve_fit(target_func, X, y)
plt.figure(figsize=(10, 5))
plt.plot(X, target_func(X, *popt), '--')
plt.plot(X, y, 'ro')
plt.legend()
plt.show()
输出
---------------------------------------------------------------------------
RuntimeError Traceback (most recent call last)
<ipython-input-243-17421b6b0c14> in <module>()
18 y = df['y']
19
---> 20 popt, pcov = curve_fit(target_func, X, y)
21
22 plt.figure(figsize=(10, 5))
/Users/evgenyp/.virtualenvs/kindle-dev/lib/python2.7/site-packages/scipy/optimize/minpack.pyc in curve_fit(f, xdata, ydata, p0, sigma, absolute_sigma, check_finite, bounds, method, **kwargs)
653 cost = np.sum(infodict['fvec'] ** 2)
654 if ier not in [1, 2, 3, 4]:
--> 655 raise RuntimeError("Optimal parameters not found: " + errmsg)
656 else:
657 res = least_squares(func, p0, args=args, bounds=bounds, method=method,
RuntimeError: Optimal parameters not found: Number of calls to function has reached maxfev = 800.
作为回溯状态,在没有找到固定点(终止算法)的情况下达到了函数计算的最大次数。您可以使用选项 maxfev
增加最大数量。对于此示例,设置 maxfev=2000
足以成功终止算法。
但是,解决方案并不令人满意。这是由于算法为变量选择了一个(默认的)初始估计,对于这个例子来说,这是不好的(所需的大量迭代就是一个指标)。提供另一个初始化点(通过简单的反复试验找到)可以很好地拟合,而无需增加 maxfev
.
两个拟合以及与数据的视觉比较如下所示。
x = np.asarray([ 1000, 3250, 5500, 10000, 32500, 55000, 77500, 100000, 200000 ])
y = np.asarray([ 1100, 500, 288, 200, 113, 67, 52, 44, 5 ])
sol1 = curve_fit(func_powerlaw, x, y, maxfev=2000 )
sol2 = curve_fit(func_powerlaw, x, y, p0 = np.asarray([-1,10**5,0]))
你的 func_powerlaw
严格来说不是幂律,因为它有一个加法常数。
一般来说,如果您想快速直观地评估幂律关系,您会
plot(log(x),log(y))
或
loglog(x,y)
两者都应该给出一条直线,尽管它们之间存在细微差别(特别是关于曲线拟合)。
所有这些都没有附加常数,这打乱了幂律关系。
如果您想要根据对数对数尺度(通常需要)来拟合加权数据的幂律,您可以使用下面的代码。
import numpy as np
from scipy.optimize import curve_fit
def powlaw(x, a, b) :
return a * np.power(x, b)
def linlaw(x, a, b) :
return a + x * b
def curve_fit_log(xdata, ydata) :
"""Fit data to a power law with weights according to a log scale"""
# Weights according to a log scale
# Apply fscalex
xdata_log = np.log10(xdata)
# Apply fscaley
ydata_log = np.log10(ydata)
# Fit linear
popt_log, pcov_log = curve_fit(linlaw, xdata_log, ydata_log)
#print(popt_log, pcov_log)
# Apply fscaley^-1 to fitted data
ydatafit_log = np.power(10, linlaw(xdata_log, *popt_log))
# There is no need to apply fscalex^-1 as original data is already available
return (popt_log, pcov_log, ydatafit_log)
我试图找到一条曲线拟合我的数据,该曲线在视觉上似乎具有幂律分布。
我希望利用 scipy.optimize.curve_fit,但无论我尝试什么函数或数据规范化,我都会收到 RuntimeError(未找到参数或溢出)或一条与我的数据不相符的曲线,甚至远程.请帮我弄清楚我在这里做错了什么。
%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
df = pd.DataFrame({
'x': [ 1000, 3250, 5500, 10000, 32500, 55000, 77500, 100000, 200000 ],
'y': [ 1100, 500, 288, 200, 113, 67, 52, 44, 5 ]
})
df.plot(x='x', y='y', kind='line', style='--ro', figsize=(10, 5))
def func_powerlaw(x, m, c, c0):
return c0 + x**m * c
target_func = func_powerlaw
X = df['x']
y = df['y']
popt, pcov = curve_fit(target_func, X, y)
plt.figure(figsize=(10, 5))
plt.plot(X, target_func(X, *popt), '--')
plt.plot(X, y, 'ro')
plt.legend()
plt.show()
输出
---------------------------------------------------------------------------
RuntimeError Traceback (most recent call last)
<ipython-input-243-17421b6b0c14> in <module>()
18 y = df['y']
19
---> 20 popt, pcov = curve_fit(target_func, X, y)
21
22 plt.figure(figsize=(10, 5))
/Users/evgenyp/.virtualenvs/kindle-dev/lib/python2.7/site-packages/scipy/optimize/minpack.pyc in curve_fit(f, xdata, ydata, p0, sigma, absolute_sigma, check_finite, bounds, method, **kwargs)
653 cost = np.sum(infodict['fvec'] ** 2)
654 if ier not in [1, 2, 3, 4]:
--> 655 raise RuntimeError("Optimal parameters not found: " + errmsg)
656 else:
657 res = least_squares(func, p0, args=args, bounds=bounds, method=method,
RuntimeError: Optimal parameters not found: Number of calls to function has reached maxfev = 800.
作为回溯状态,在没有找到固定点(终止算法)的情况下达到了函数计算的最大次数。您可以使用选项 maxfev
增加最大数量。对于此示例,设置 maxfev=2000
足以成功终止算法。
但是,解决方案并不令人满意。这是由于算法为变量选择了一个(默认的)初始估计,对于这个例子来说,这是不好的(所需的大量迭代就是一个指标)。提供另一个初始化点(通过简单的反复试验找到)可以很好地拟合,而无需增加 maxfev
.
两个拟合以及与数据的视觉比较如下所示。
x = np.asarray([ 1000, 3250, 5500, 10000, 32500, 55000, 77500, 100000, 200000 ])
y = np.asarray([ 1100, 500, 288, 200, 113, 67, 52, 44, 5 ])
sol1 = curve_fit(func_powerlaw, x, y, maxfev=2000 )
sol2 = curve_fit(func_powerlaw, x, y, p0 = np.asarray([-1,10**5,0]))
你的 func_powerlaw
严格来说不是幂律,因为它有一个加法常数。
一般来说,如果您想快速直观地评估幂律关系,您会
plot(log(x),log(y))
或
loglog(x,y)
两者都应该给出一条直线,尽管它们之间存在细微差别(特别是关于曲线拟合)。
所有这些都没有附加常数,这打乱了幂律关系。
如果您想要根据对数对数尺度(通常需要)来拟合加权数据的幂律,您可以使用下面的代码。
import numpy as np
from scipy.optimize import curve_fit
def powlaw(x, a, b) :
return a * np.power(x, b)
def linlaw(x, a, b) :
return a + x * b
def curve_fit_log(xdata, ydata) :
"""Fit data to a power law with weights according to a log scale"""
# Weights according to a log scale
# Apply fscalex
xdata_log = np.log10(xdata)
# Apply fscaley
ydata_log = np.log10(ydata)
# Fit linear
popt_log, pcov_log = curve_fit(linlaw, xdata_log, ydata_log)
#print(popt_log, pcov_log)
# Apply fscaley^-1 to fitted data
ydatafit_log = np.power(10, linlaw(xdata_log, *popt_log))
# There is no need to apply fscalex^-1 as original data is already available
return (popt_log, pcov_log, ydatafit_log)