2016-11-23 10 views
0

私のデータにべき乗則曲線を当てはめる際に問題が発生しました。私は2つのデータセットしている:bins1とカーブフィッティングnumpy.linalg.lstsqを使用して、中に細かい演技bins2

bins1(私はその後、べき乗則方程式を得るためにnp.exp(coefs[0])*x**coefs[1]を使用)

一方、bins2は奇妙な演技されます

どちらのデータも、excelが私に見せているものとは違う方程式を持っています(悪化したR-二乗)。累乗法則フィッティングscipy、numpy not working

import numpy as np 
import matplotlib.pyplot as plt 
bins1 = np.array([[6.769318871738219667e-03, 
      1.306418618130891773e-02, 
      1.912138120913448383e-02, 
      2.545189874466026111e-02, 
      3.214689891729670401e-02, 
      4.101898933375244805e-02, 
      5.129862592803200588e-02, 
      6.636505322669797313e-02, 
      8.409809827572585494e-02, 
      1.058164348650862258e-01, 
      1.375849753230810046e-01, 
      1.830664031837437311e-01, 
      2.682454535427478137e-01, 
      3.912508246490400410e-01, 
      5.893271848997768680e-01, 
      8.480213305038615257e-01, 
      2.408136266017391058e+00, 
      3.629192766488219313e+00, 
      4.639246557509275171e+00, 
      9.901792214343277720e+00], 
      [8.501658465758301112e-04, 
       1.562697718429977012e-03, 
       1.902062808421856087e-04, 
       4.411817741488644959e-03, 
       3.409236963162485048e-03, 
       1.686099657013027898e-03, 
       3.643231240239608402e-03, 
       2.544120616413291154e-04, 
       2.549036204611017029e-02, 
       3.527340723977697573e-02, 
       5.038482027310990652e-02, 
       5.617932487522721979e-02, 
       1.620407270423956103e-01, 
       1.906538999080910068e-01, 
       3.180688368126549093e-01, 
       2.364903188268162038e-01, 
       3.267322385964683273e-01, 
       9.384571074801122403e-01, 
       4.419747716107813029e-01, 
       9.254710022316929852e+00]]).T 
bins2 = np.array([[6.522512685133712192e-03, 
       1.300415548684437199e-02, 
       1.888928895701269539e-02, 
       2.509905819337970856e-02, 
       3.239654633369139919e-02, 
       4.130706234846069635e-02, 
       5.123820846515786398e-02, 
       6.444380072984744190e-02, 
       8.235238352205621892e-02, 
       1.070907072127811749e-01, 
       1.403438221033725120e-01, 
       1.863115065963684147e-01, 
       2.670209758710758163e-01, 
       4.003337413814173074e-01, 
       6.549054078382223754e-01, 
       1.116611087124244062e+00, 
       2.438604844718367914e+00, 
       3.480674117919704269e+00, 
       4.410201659398489404e+00, 
       6.401903059926267403e+00], 
      [1.793454543936148608e-03, 
       2.441092334386309615e-03, 
       2.754373929745804715e-03, 
       1.182752729942167062e-03, 
       1.357797177773524414e-03, 
       6.711673916715021199e-03, 
       1.392761674092503343e-02, 
       1.127957613093066511e-02, 
       7.928803089359596004e-03, 
       2.524609593305639915e-02, 
       5.698702885370290905e-02, 
       8.607729156137132465e-02, 
       2.453761830112021203e-01, 
       9.734443815196883176e-02, 
       1.487480479168299119e-01, 
       9.918002699934079791e-01, 
       1.121298151253063535e+00, 
       1.389239135742518227e+00, 
       4.254082922056571237e-01, 
       2.643453492951096440e+00]]).T 

bins = bins1 #change to bins2 to see results for bins2 

def fit(x,a,m): # power-law fit (based on previous studies) 
    return a*(x**m) 

coefs= np.linalg.lstsq(np.vstack([np.ones(len(bins[:,0])), np.log(bins[:,0]), bins[:,0]]).T, np.log(bins[:,1]))[0] # calculating fitting coefficients (a,m) 
y_predict = fit(bins[:,0],np.exp(coefs[0]),coefs[1]) # prediction based of fitted model 
model_plot = plt.loglog(bins[:,0],bins[:,1],'o',label="error") 
fit_line = plt.plot(bins[:,0],y_predict,'r', label="fit") 
plt.ylabel('Y (bins[:,1])') 
plt.xlabel('X (bins[:,0])') 
plt.title('model') 
plt.legend(loc='best') 
plt.show(model_plot,fit_line) 

def R_sqr (y,y_predict): # calculating R squared value to measure fitting accuracy 
    rsdl = y - y_predict 
    ss_res = np.sum(rsdl**2) 
    ss_tot = np.sum((y-np.mean(y))**2) 
    R2 = 1-(ss_res/ss_tot) 
    R2 = np.around(R2,decimals=4) 
    return R2 

R2= R_sqr(bins[:,1],y_predict) 
print ('(R^2 = %s)' % (R2)) 

bins1に適合式[X]、[Y]:python: y = 0.337*(x)^1.223 (R^2 = 0.7773), excel: y = 0.289*(x)^1.174 (R^2 = 0.8548)

bins2に適合式[X]、[ここで

コード(およびデータ)でありますy]]:python: y = 0.509*(x)^1.332 (R^2 = -1.753), excel: y = 0.311*(x)^1.174 (R^2 = 0.9116)

これらは30個のうち2つのサンプルデータセットです。このデータは私のデータにランダムに表示され、一部は「-150」の周りにR-二乗されています。
私は実際に悪い "scissy" curve_fitしかし、私は良い結果を得ていない!

誰でも、どのようにPythonでExcelのようなフィット感を得るのを知っていますか?

答えて

1

ログスペースに変換されていないYを使用してR-平方和を計算しようとしています。以下の変更により、妥当なR二乗値が得られます。

R2 = R_sqr(np.log(bins[:,1]), np.log(y_predict)) 
+0

驚くばかり!ありがとう – Sina