如何使用Python和NumPy计算R方？

Question

我正在使用Python和Numpy来计算任意阶数的最佳拟合多项式。我传入一组x值和y值，以及我想要拟合的多项式的阶数（比如线性、二次等）。

到这里都没问题，但我还想计算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

Answer 1

来自 yanl（又一个库）的 sklearn.metrics 里有一个叫 r2_score 的函数；

from sklearn.metrics import r2_score

coefficient_of_dermination = r2_score(y, p(x))

Answer 2

虽然回复得很晚，但如果有人需要一个现成的函数来解决这个问题，可以参考以下内容：

scipy.stats.linregress

也就是说：

slope, intercept, r_value, p_value, std_err = scipy.stats.linregress(x, y)

就像@Adam Marples的回答中提到的那样。

Answer 3

根据numpy.polyfit的说明，它用于进行线性回归。具体来说，使用numpy.polyfit时，如果你设置一个数字'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