파이썬과 Numpy를 사용하여 r-제곱을 계산하려면 어떻게 해야 합니까?
저는 파이썬과 Numpy를 사용하여 임의의 정도의 가장 적합한 다항식을 계산하고 있습니다.x 값 목록, y 값 및 적합할 다항식의 정도(선형, 2차 등)를 전달합니다.
이 정도는 효과가 있지만, 저는 또한 r(상관계수)과 r-제곱(결정계수)도 계산하려고 합니다.Excel의 가장 적합한 추세선 기능과 계산한 r 제곱 값을 비교하고 있습니다.이를 사용하여 선형 최적 적합도(도는 1)에 대한 r-제곱을 올바르게 계산하고 있습니다.그러나 차수가 1보다 큰 다항식에서는 함수가 작동하지 않습니다.
Excel은 이것을 할 수 있습니다.Numpy를 사용하여 고차 다항식에 대한 r 제곱을 어떻게 계산합니까?
제 기능은 다음과 같습니다.
import numpy
# Polynomial Regression
def polyfit(x, y, degree):
results = {}
coeffs = numpy.polyfit(x, y, degree)
# Polynomial Coefficients
results['polynomial'] = coeffs.tolist()
correlation = numpy.corrcoef(x, y)[0,1]
# r
results['correlation'] = correlation
# r-squared
results['determination'] = correlation**2
return results
매우 늦은 답변입니다. 하지만 누군가 이를 위해 준비된 기능이 필요할 경우:
예.
slope, intercept, r_value, p_value, std_err = scipy.stats.linregress(x, y)
@아담 마플스의 대답처럼.
yanl(아직 다른 라이브러리)로부터.sklearn.metrics기능이 있습니다.
from sklearn.metrics import r2_score
coefficient_of_dermination = r2_score(y, p(x))
바보 같은 놈한테서.다중 적합 문서, 적합 선형 회귀 분석입니다.특히, 바보야.정도가 'd'인 폴리핏은 평균 함수를 사용한 선형 회귀 분석을 적합시킵니다.
E(y|x) = p_d * x**d + p_{d-1} * x **(d-1) + ...p_1 * x + p_0
따라서 해당 적합치에 대한 R-제곱을 계산하면 됩니다.선형 회귀에 대한 위키백과 페이지는 자세한 내용을 제공합니다.당신은 몇 가지 방법으로 계산할 수 있는 R^2에 관심이 있습니다, 아마도 가장 쉬운 것은
SST = Sum(i=1..n) (y_i - y_bar)^2
SSReg = Sum(i=1..n) (y_ihat - y_bar)^2
Rsquared = SSReg/SST
여기서 'y_bar'를 y의 평균으로 사용하고 'y_ihat'을 각 점의 적합치로 사용합니다.
저는 numpy에 대해 잘 모르기 때문에(저는 보통 R에서 일합니다), 아마도 R 제곱을 계산하는 더 쉬운 방법이 있을 것입니다. 하지만 다음은 정확해야 합니다.
import numpy
# Polynomial Regression
def polyfit(x, y, degree):
results = {}
coeffs = numpy.polyfit(x, y, degree)
# Polynomial Coefficients
results['polynomial'] = coeffs.tolist()
# r-squared
p = numpy.poly1d(coeffs)
# fit values, and mean
yhat = p(x) # or [p(z) for z in x]
ybar = numpy.sum(y)/len(y) # or sum(y)/len(y)
ssreg = numpy.sum((yhat-ybar)**2) # or sum([ (yihat - ybar)**2 for yihat in yhat])
sstot = numpy.sum((y - ybar)**2) # or sum([ (yi - ybar)**2 for yi in y])
results['determination'] = ssreg / sstot
return results
저는 x와 y가 배열과 같은 이것을 성공적으로 사용하고 있습니다.
참고: 선형 회귀 분석 전용
def rsquared(x, y):
""" Return R^2 where x and y are array-like."""
slope, intercept, r_value, p_value, std_err = scipy.stats.linregress(x, y)
return r_value**2
를 추천할 목적으로 했습니다.numpy.corrcoef어리석게도 원래 질문이 이미 사용되고 있다는 것을 깨닫지 못합니다.corrcoef그리고 사실은 고차 다항식 적합치에 대해 묻고 있었습니다.저는 통계 모형을 사용하여 다항식 r 제곱 질문에 실제 솔루션을 추가했고, 주제에서 벗어나 있지만 누군가에게 잠재적으로 유용한 원래 벤치마크를 남겼습니다.
statsmodels에는 계기 있이니다를 이 있습니다.r^2직접 다항식 적합의 두 가지 방법이 있습니다.
import statsmodels.api as sm
import statsmodels.formula.api as smf
# Construct the columns for the different powers of x
def get_r2_statsmodels(x, y, k=1):
xpoly = np.column_stack([x**i for i in range(k+1)])
return sm.OLS(y, xpoly).fit().rsquared
# Use the formula API and construct a formula describing the polynomial
def get_r2_statsmodels_formula(x, y, k=1):
formula = 'y ~ 1 + ' + ' + '.join('I(x**{})'.format(i) for i in range(1, k+1))
data = {'x': x, 'y': y}
return smf.ols(formula, data).fit().rsquared # or rsquared_adj
다음과 같은 이점을 더욱 활용하기 위해statsmodels 합니다. 이 요약은 Jupyter/JupyterRich HTML 할 수 .나는 Python 노트북.결과 개체는 유용한 통계 메트릭에 대한 액세스를 제공할 뿐만 아니라rsquared.
model = sm.OLS(y, xpoly)
results = model.fit()
results.summary()
아래는 다양한 선형 회귀 r^2 방법을 벤치마킹한 저의 원래 답변입니다.
질문에 사용된 코르코프 함수는 상관 계수를 계산합니다.r단일 선형 회귀에 대해서만, 그래서 그것은 질문을 다루지 않습니다.r^2고차 다항식 적합치의 경우.회귀의 인 계산 회귀라는 것을 r.
def get_r2_numpy_corrcoef(x, y):
return np.corrcoef(x, y)[0, 1]**2
1000개의 랜덤(x, y) 포인트에 대해 여러 방법을 비교한 결과는 다음과 같습니다.
- Python( Pure Python()
r계산)- 1,000개의 루프(루프당 최고 3: 1.59ms)
- Numpy 폴리핏(n차 다항식 적합치에 적용 가능)
- 1000개 루프, 3개 루프당 최적: 326µs
- (직접 Numpy Manual (직접))
r계산)- 10000개 루프, 3개 루프당 62.1µs 최적
- Numpy Corrcoef (으)로 표시)
r계산)- 10000개 루프(루프당 56.6µs), 최적의 3개 루프
- Scipy를 사용한 분석)
r결과로서)- 1000개 루프, 루프당 3:676µs 중 최고치
- 통계분석 모형(일차 다항식 및 기타 여러 적합치)
- 1,000개의 루프(루프당 3:422µs)
코르코프 방법은 numpy 방법을 사용하여 r^2를 "수동으로" 계산하는 것을 아슬아슬하게 능가합니다.이는 폴리핏 방법보다 >5배 빠르고 scipy.lin 회귀보다 ~12배 빠릅니다.Numpy가 당신에게 하는 일을 강화하기 위해, 그것은 순수한 파이썬보다 28배 더 빠릅니다.같은 에 대해 잘 에 다른 그 할 입니다. 이 numba와 pypy라는 것을 할 수 합니다. 하지만 저는 이것이 저에게 충분히 설득력이 있다고 생각합니다.corrcoef 가적도다니구입한합장에ating▁for다를 계산하는 데 가장 좋은 입니다.r단순 선형 회귀 분석의 경우.
여기 제 벤치마킹 코드가 있습니다.제가 주피터 노트북(아이피톤 노트북이라고 하기는 힘들지만...)에서 복사 붙여넣기를 했기 때문에 도중에 뭔가 고장이 났다면 사과드립니다.%timeit magic 명령에는 IPython이 필요합니다.
import numpy as np
from scipy import stats
import statsmodels.api as sm
import math
n=1000
x = np.random.rand(1000)*10
x.sort()
y = 10 * x + (5+np.random.randn(1000)*10-5)
x_list = list(x)
y_list = list(y)
def get_r2_numpy(x, y):
slope, intercept = np.polyfit(x, y, 1)
r_squared = 1 - (sum((y - (slope * x + intercept))**2) / ((len(y) - 1) * np.var(y, ddof=1)))
return r_squared
def get_r2_scipy(x, y):
_, _, r_value, _, _ = stats.linregress(x, y)
return r_value**2
def get_r2_statsmodels(x, y):
return sm.OLS(y, sm.add_constant(x)).fit().rsquared
def get_r2_python(x_list, y_list):
n = len(x_list)
x_bar = sum(x_list)/n
y_bar = sum(y_list)/n
x_std = math.sqrt(sum([(xi-x_bar)**2 for xi in x_list])/(n-1))
y_std = math.sqrt(sum([(yi-y_bar)**2 for yi in y_list])/(n-1))
zx = [(xi-x_bar)/x_std for xi in x_list]
zy = [(yi-y_bar)/y_std for yi in y_list]
r = sum(zxi*zyi for zxi, zyi in zip(zx, zy))/(n-1)
return r**2
def get_r2_numpy_manual(x, y):
zx = (x-np.mean(x))/np.std(x, ddof=1)
zy = (y-np.mean(y))/np.std(y, ddof=1)
r = np.sum(zx*zy)/(len(x)-1)
return r**2
def get_r2_numpy_corrcoef(x, y):
return np.corrcoef(x, y)[0, 1]**2
print('Python')
%timeit get_r2_python(x_list, y_list)
print('Numpy polyfit')
%timeit get_r2_numpy(x, y)
print('Numpy Manual')
%timeit get_r2_numpy_manual(x, y)
print('Numpy corrcoef')
%timeit get_r2_numpy_corrcoef(x, y)
print('Scipy')
%timeit get_r2_scipy(x, y)
print('Statsmodels')
%timeit get_r2_statsmodels(x, y)
7월 28일/21일 벤치마크 결과. (파이썬 3.7, 누피 1.19, 스키피 1.6, 통계 모델 0.12)
Python
2.41 ms ± 180 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
Numpy polyfit
318 µs ± 44.3 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
Numpy Manual
79.3 µs ± 4.05 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
Numpy corrcoef
83.8 µs ± 1.37 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
Scipy
221 µs ± 7.12 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
Statsmodels
375 µs ± 3.63 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
다음은 Python과 Numpy의 가중치 r 제곱을 계산하는 함수입니다(대부분의 코드는 sklearn에서 옵니다).
from __future__ import division
import numpy as np
def compute_r2_weighted(y_true, y_pred, weight):
sse = (weight * (y_true - y_pred) ** 2).sum(axis=0, dtype=np.float64)
tse = (weight * (y_true - np.average(
y_true, axis=0, weights=weight)) ** 2).sum(axis=0, dtype=np.float64)
r2_score = 1 - (sse / tse)
return r2_score, sse, tse
예:
from __future__ import print_function, division
import sklearn.metrics
def compute_r2_weighted(y_true, y_pred, weight):
sse = (weight * (y_true - y_pred) ** 2).sum(axis=0, dtype=np.float64)
tse = (weight * (y_true - np.average(
y_true, axis=0, weights=weight)) ** 2).sum(axis=0, dtype=np.float64)
r2_score = 1 - (sse / tse)
return r2_score, sse, tse
def compute_r2(y_true, y_predicted):
sse = sum((y_true - y_predicted)**2)
tse = (len(y_true) - 1) * np.var(y_true, ddof=1)
r2_score = 1 - (sse / tse)
return r2_score, sse, tse
def main():
'''
Demonstrate the use of compute_r2_weighted() and checks the results against sklearn
'''
y_true = [3, -0.5, 2, 7]
y_pred = [2.5, 0.0, 2, 8]
weight = [1, 5, 1, 2]
r2_score = sklearn.metrics.r2_score(y_true, y_pred)
print('r2_score: {0}'.format(r2_score))
r2_score,_,_ = compute_r2(np.array(y_true), np.array(y_pred))
print('r2_score: {0}'.format(r2_score))
r2_score = sklearn.metrics.r2_score(y_true, y_pred,weight)
print('r2_score weighted: {0}'.format(r2_score))
r2_score,_,_ = compute_r2_weighted(np.array(y_true), np.array(y_pred), np.array(weight))
print('r2_score weighted: {0}'.format(r2_score))
if __name__ == "__main__":
main()
#cProfile.run('main()') # if you want to do some profiling
출력:
r2_score: 0.9486081370449679
r2_score: 0.9486081370449679
r2_score weighted: 0.9573170731707317
r2_score weighted: 0.9573170731707317
f_i는 적합치의 예측 값이고, y_{av}는 관측 데이터의 평균 y_i는 관측 데이터 값입니다.w_i는 각 데이터 점(일반적으로 w_i=1)에 적용되는 가중치입니다.SSE는 오차로 인한 제곱합이고 SST는 전체 제곱합입니다.
관심이 있는 경우 R: https://gist.github.com/dhimmel/588d64a73fa4fef02c8f (미러)의 코드입니다.
다음은 y와 y_hat가 판다 계열이라고 가정할 때 실제 및 예측 값으로부터 R^2를 계산하는 매우 간단한 파이썬 함수입니다.
def r_squared(y, y_hat):
y_bar = y.mean()
ss_tot = ((y-y_bar)**2).sum()
ss_res = ((y-y_hat)**2).sum()
return 1 - (ss_res/ss_tot)
R-제곱은 선형 회귀 분석에만 적용되는 통계량입니다.
기본적으로 선형 회귀 분석을 통해 데이터의 변동을 설명할 수 있는 정도를 측정합니다.
따라서 각 결과 변수의 평균으로부터의 총 제곱 편차인 "총 제곱합"을 계산합니다.
여기서 y_bar는 그들의 평균입니다.
그런 다음 적합치가 평균과 얼마나 다른지 나타내는 "회귀 제곱합"을 계산합니다.
그리고 그 둘의 비율을 찾아보세요.
이제 다항식 적합치를 위해 해야 할 일은 해당 모델의 y_hats를 연결하는 것뿐이지만, r-제곱이라고 부르는 것은 정확하지 않습니다.
여기 제가 찾은 링크가 그것을 조금 말해줍니다.
r 제곱에 대한 위키백과 기사는 단순한 선형 회귀가 아닌 일반 모델 적합에 사용될 수 있음을 시사합니다.
numpy 모듈 사용(python3에서 테스트됨):
import numpy as np
def linear_regression(x, y):
coefs = np.polynomial.polynomial.polyfit(x, y, 1)
ffit = np.poly1d(coefs)
m = ffit[0]
b = ffit[1]
eq = 'y = {}x + {}'.format(round(m, 3), round(b, 3))
rsquared = np.corrcoef(x, y)[0, 1]**2
return rsquared, eq, m, b
rsquared, eq, m, b = linear_regression(x,y)
print(rsquared, m, b)
print(eq)
출력:
0.013378252355751777 0.1316331351105754 0.7928782850418713 y = 0.132x + 0.793
참고: r² ≠ R²
r²를 "결정 계수"라고 합니다.
R²는 Pearson 계수의 제곱입니다.
공식적으로 r²로 통합되는 R²는 최소 제곱 적합이므로 아마도 원하는 것일 것입니다. R²는 단순한 합의 부분보다 더 좋습니다.Numpy는 Pearson이 사실상의 상관 계수라고 가정하는 "corcoef"라고 부르는 것을 두려워하지 않습니다.
이 코드를 직접 실행하면 다항식을 찾을 수 있으며, 더 많은 설명이 필요할 경우 아래에 주석을 달 수 있는 R-값을 찾을 수 있습니다.
from scipy.stats import linregress
import numpy as np
x = np.array([1,2,3,4,5,6])
y = np.array([2,3,5,6,7,8])
p3 = np.polyfit(x,y,3) # 3rd degree polynomial, you can change it to any degree you want
xp = np.linspace(1,6,6) # 6 means the length of the line
poly_arr = np.polyval(p3,xp)
poly_list = [round(num, 3) for num in list(poly_arr)]
slope, intercept, r_value, p_value, std_err = linregress(x, poly_list)
print(r_value**2)
scipy.stats.lin회귀 소스에서.평균 제곱합 방법을 사용합니다.
import numpy as np
x = np.array(x)
y = np.array(y)
# average sum of squares:
ssxm, ssxym, ssyxm, ssym = np.cov(x, y, bias=1).flat
r_num = ssxym
r_den = np.sqrt(ssxm * ssym)
r = r_num / r_den
if r_den == 0.0:
r = 0.0
else:
r = r_num / r_den
if r > 1.0:
r = 1.0
elif r < -1.0:
r = -1.0
언급URL : https://stackoverflow.com/questions/893657/how-do-i-calculate-r-squared-using-python-and-numpy
'programing' 카테고리의 다른 글
| 성능 - 스프링 부팅 - 서버 응답 시간 (0) | 2023.07.21 |
|---|---|
| 오라클에 내장된 해시 기능이 있습니까? (0) | 2023.07.21 |
| 부울 ID == 참 대 참 (0) | 2023.07.21 |
| REST API용 spring-boot-starter-web과 spring-boot-starter-data-rest의 차이점 (0) | 2023.07.21 |
| OracleDB에서 "예상 NUMBER but got BINARY"를 가져오지 않고 Long 값의 null을 쿼리하려면 어떻게 해야 합니까? (0) | 2023.07.21 |