2017-01-20 6 views
2

제로 오더를 가진 두 번째 차수 함수에 적합하려고합니다.2 차 함수를 피팅하는 데 제로 인터 셉션을 강제로 적용하는 방법은 무엇입니까? (Python)

y**2 = 14.29566 * np.pi * x 

또는

y = np.sqrt(14.29566 * np.pi * x) 

: 지금 나는 그것을 그릴 때, 나는 기능에 의해 포인트 출력 세트에 맞추려고 오전 Y-INT> 0으로 선을 얻고있다 D = 3.57391553 인 두 개의 데이터 세트 x 및 y로 변환합니다. 내 피팅 루틴은 다음과 같습니다

mod_ols = sm.OLS(y,x) 
res_ols = mod_ols.fit() 

하지만 선형 함수와 달리 2 차 함수의 계수를 생성하는 방법을 이해하지 않으며 방법 :

z = np.polyfit(x,y,2)       # Generate curve coefficients 
p = np.poly1d(z)        # Generate curve function 

xp = np.linspace(0, catalog.tlimit, 10)  # generate input values 

plt.scatter(x,y) 
plt.plot(xp, p(xp), '-') 
plt.show() 

가 나는 또한 statsmodels.ols를 사용하여 시도 y-int를 0으로 설정합니다. 비슷한 모양의 다른 게시물에서 선형 맞춤으로 0의 y-int를 처리하는 것을 보았습니다. 그러나 2 차 함수로이를 수행하는 방법을 파악할 수 없었습니다. 지금

플롯 : I want to anchor the curve at (0,0)

데이터 : 난 당신이 OLS과 다항식 회귀 라인에 맞게 원하는 제대로 이해하면

x = [0., 0.00325492, 0.00650985, 0.00976477, 0.01301969, 0.01627462, 0.01952954, 0.02278447, 
    0.02603939, 0.02929431, 0.03254924, 0.03580416, 0.03905908, 0.04231401] 
y = [0., 0.38233801, 0.5407076, 0.66222886, 0.76467602, 0.85493378, 0.93653303, 1.01157129, 
    1.0814152, 1.14701403, 1.20905895, 1.26807172, 1.32445772, 1.3785393] 
+0

다항식 설계 행렬을 만들고 상수 열 x [:, 없음] ** np.arange (1, k + 1)을 제거하고 OLS를 사용하려면 np.vander를 사용하십시오. – user333700

답변

1

, 당신이 시도 할 수있는 우리가 볼 수있는 (다음과 같은 피팅 된 다항식의 차수가 증가하면, 모델은 데이터에 점점 더 잘 맞는다) :

xp = np.linspace(0, 0.05, 10)  # generate input values 
pred1 = p(xp) # prediction with np.poly as you have done 

import pandas as pd 
data = pd.DataFrame(data={'x':x, 'y':y}) 

# let's first fit 2nd order polynomial with OLS with intercept 
olsres2 = sm.ols(formula = 'y ~ x + I(x**2)', data = data).fit() 
print olsres2.summary() 


          OLS Regression Results 
============================================================================== 
Dep. Variable:      y R-squared:      0.978 
Model:       OLS Adj. R-squared:     0.974 
Method:     Least Squares F-statistic:      243.1 
Date:    Sat, 21 Jan 2017 Prob (F-statistic):   7.89e-10 
Time:      04:16:22 Log-Likelihood:     20.323 
No. Observations:     14 AIC:       -34.65 
Df Residuals:      11 BIC:       -32.73 
Df Model:       2 
Covariance Type:   nonrobust 
============================================================================== 
       coef std err   t  P>|t|  [95.0% Conf. Int.] 
------------------------------------------------------------------------------ 
Intercept  0.1470  0.045  3.287  0.007   0.049  0.245 
x    52.4655  4.907  10.691  0.000  41.664 63.267 
I(x ** 2) -580.4730 111.820  -5.191  0.000  -826.588 -334.358 
============================================================================== 
Omnibus:      4.803 Durbin-Watson:     1.164 
Prob(Omnibus):     0.091 Jarque-Bera (JB):    2.101 
Skew:       -0.854 Prob(JB):      0.350 
Kurtosis:      3.826 Cond. No.      6.55e+03 
============================================================================== 

pred2 = olsres2.predict(data) # predict with the fitted model 
# fit 3rd order polynomial with OLS with intercept 
olsres3 = sm.ols(formula = 'y ~ x + I(x**2) + I(x**3)', data = data).fit() 
pred3 = olsres3.predict(data) # predict 

plt.scatter(x,y) 
plt.plot(xp, pred1, '-r', label='np.poly') 
plt.plot(x, pred2, '-b', label='ols.o2') 
plt.plot(x, pred3, '-g', label='ols.o3') 
plt.legend(loc='upper left') 
plt.show() 
# now let's fit the polynomial regression lines this time without intercept 
olsres2 = sm.ols(formula = 'y ~ x + I(x**2)-1', data = data).fit() 
print olsres2.summary() 
           OLS Regression Results 
============================================================================== 
Dep. Variable:      y R-squared:      0.993 
Model:       OLS Adj. R-squared:     0.992 
Method:     Least Squares F-statistic:      889.6 
Date:    Sat, 21 Jan 2017 Prob (F-statistic):   9.04e-14 
Time:      04:16:24 Log-Likelihood:     15.532 
No. Observations:     14 AIC:       -27.06 
Df Residuals:      12 BIC:       -25.79 
Df Model:       2 
Covariance Type:   nonrobust 
============================================================================== 
       coef std err   t  P>|t|  [95.0% Conf. Int.] 
------------------------------------------------------------------------------ 
x    65.8170  3.714  17.723  0.000  57.726 73.908 
I(x ** 2) -833.6787 109.279  -7.629  0.000  -1071.777 -595.580 
============================================================================== 
Omnibus:      1.716 Durbin-Watson:     0.537 
Prob(Omnibus):     0.424 Jarque-Bera (JB):    1.341 
Skew:       0.649 Prob(JB):      0.511 
Kurtosis:      2.217 Cond. No.       118. 
============================================================================== 
pred2 = olsres2.predict(data) 
# fit 3rd order polynomial with OLS without intercept 
olsres3 = sm.ols(formula = 'y ~ x + I(x**2) + I(x**3) -1', data = data).fit() 
pred3 = olsres3.predict(data) 

plt.scatter(x,y) 
plt.plot(xp, pred1, '-r', label='np.poly') 
plt.plot(x, pred2, '-b', label='ols.o2') 
plt.plot(x, pred3, '-g', label='ols.o3') 
plt.legend(loc='upper left') 
plt.show() 

enter image description here

1

3,215,895,913,이 코드는 질문에 매우 직접적인 대답을 제공합니다. 나는 patsy 언어로 먼저 실험했음을 인정해야합니다 (https://patsy.readthedocs.io/en/latest/ 참조).

여기의 핵심은 모델의 선언입니다 : np.power (y, 2) ~ x - 1. 절편을 생략하고 y의 사각형을 x에 회귀한다는 것을 의미합니다. 일단 계수가 x으로 나누어지면 pi로 나누어 원하는 숫자를 얻습니다.

x = [0., 0.00325492, 0.00650985, 0.00976477, 0.01301969, 0.01627462, 0.01952954, 0.02278447, 0.02603939, 0.02929431, 0.03254924, 0.03580416, 0.03905908, 0.04231401] 
y = [0., 0.38233801, 0.5407076, 0.66222886, 0.76467602, 0.85493378, 0.93653303, 1.01157129, 1.0814152, 1.14701403, 1.20905895, 1.26807172, 1.32445772, 1.3785393] 

data = { 'x': x, 'y': y } 

import statsmodels.formula.api as smf 
import numpy as np 

model = smf.ols(formula = 'np.power(y,2) ~ x - 1', data = data) 
results = model.fit() 
print (results.params) 

predicteds = results.predict(data) 

import matplotlib.pyplot as plt 

plt.scatter(x,y) 
plt.plot(x, [_**2 for _ in y], 'go') 
plt.plot(x, predicteds, '-b') 
plt.show() 

회귀 계수는 44.911147입니다. 다음은 그래픽 결과입니다.

plot

OK, 나는 의심 해요. ;)