Python中非线性最小二乘曲线拟合的R²值计算

Question

我已经写了一个代码，我想计算一个非线性拟合的R²值，通过给定的力-深度关系。在

我对给定的x和y数据使用的代码是：

   ydata=npy.array(N)
   xdata=npy.array(depth)

   #Exponential Law

   def func1 (x,Bo,k):
       return Bo*npy.exp(x*k)
   popt, pcov, infodict, mesg, ier = curve_fit(func1, xdata, ydata,p0=(1.0,0.2),full_output=True)        
   Bo=popt[0]
   k=popt[1]

   SSR1=sum((func1(ydata,Bo,k)-xdata.mean())**2)
   SST1=sum((xdata-func1(ydata,Bo,k))**2)
   rsquared1=SSR1/SST1

指数定律和：

   #Power Law

   def func2(a,Bp,z):
       return Bp*a**z
   popt2, pcov2=curve_fit(func2,xdata,ydata,p0=(1,0.2),bounds=([-npy.inf,0],npy.inf))

   Bp=popt2[0]
   z=popt2[1]

   residuals2 = func2(ydata,Bp,z)-xdata.mean()
   fres2=sum(residuals2**2)
   ss_tot2=sum((xdata-func2(ydata,Bp,z))**2)
   rsquared2=(fres2/ss_tot2)

为了幂律。 根据rsquared=SSR/SST，这个值应该在0和1之间。不幸的是，对于一些值，我得到一个比1稍大的rsquared。在

r平方大于1的值的一个例子是：

扩展数据（深度）： [0。246810121416182022242628 303638404244464850525456第58条。]

ydata（力）： [0。000004.4 8条。2036 30.8 12.4 5.8 3.2 4。3.8 54.6 15.6 37.2 39.6 76.8 81.2 111。142.4 76.8 107.2 151.8 131.4]

我很感激你的每一次帮助

Answer 1

来自Wikipedia

Values of R2 outside the range 0 to 1 can occur where it is used to measure the agreement between observed and modeled values and where the "modeled" values are not obtained by linear regression and depending on which formulation of R2 is used.

我认为在你的例子中，你可能只是把x和{}放错了方向。。。根据this回答你的数据，我知道

import numpy as npy
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

#Exponential Law
def func1(x,Bo,k):
   return Bo*npy.exp(x*k)

#Power Law
def func2(a, Bp, z):
   return Bp*npy.power(a, z)


def get_rsq(f, y, popt):

    ss_res = npy.dot((y - func1(x, *popt)),(y - func1(x, *popt)))
    ymean = npy.mean(y)
    ss_tot = npy.dot((y-ymean),(y-ymean))
    return 1-ss_res/ss_tot



x = npy.array([0., 2., 4., 6., 8., 10., 12., 14., 16., 18., 20., 22., 24., 26., 28., 30., 36., 38., 40., 42., 44., 46., 48., 50., 52., 54., 56., 58.])
y = npy.array([0., 0., 0., 0., 0., 0., 4.4, 8., 20., 36., 30.8, 12.4, 5.8, 3.2, 4., 3.8, 54.6, 15.6, 37.2, 39.6, 76.8, 81.2, 111., 142.4, 76.8, 107.2, 151.8, 131.4])


popt, pcov = curve_fit(func1, x, y, p0=(1.0,0.2))
popt2, pcov2 = curve_fit(func2, x, y,p0=[0.0008,3.0],bounds=([-npy.inf,0],npy.inf))


plt.figure()
plt.plot(x, y, 'ko', label="Data")
plt.plot(x, func1(x, *popt), 'r-', label="Exponential Law")
plt.plot(x, func2(x, *popt), 'b-', label="Power Law")
plt.legend()
plt.show()

print "Mean R Exponential:",  get_rsq(func1, y, popt)
print "Mean R Power:",  get_rsq(func2, y, popt2)

实验结果显示

^{pr2}$

幂律拟合失败得很厉害（-4.64462440385e+140），我想这是基于类似的questions）的，我猜你为什么要加边界。我的scipy是0.17之前的版本，所以不能测试，但也许你在这里也会有更多的运气。在