最小二乘回归
在前面的章节中，我们推导了穿过“橄榄球形”散点图的回归线的斜率和截距公式。事实证明，最小二乘线的斜率和截距与之前推导的公式相同，无论散点图的形状如何。
我们在关于《小妇人》的例子中看到了这一点，但让我们在一个散点图明显不是橄榄球形的例子中确认一下。对于这些数据，我们再次感谢佛罗里达大学 Larry Winner 教授丰富的数据档案。《国际运动科学杂志》上的一项2013 年研究考察了大学铅球运动员，研究了力量与铅球距离之间的关系。研究对象包括 28 名女性大学运动员。力量是通过运动员在赛季前“1RM 高翻”中举起的最大重量（以公斤为单位）来衡量的。距离（以米为单位）是运动员的个人最好成绩。
from datascience import *
%matplotlib inline
path_data = '../../../assets/data/'
import matplotlib.pyplot as plots
plots.style.use('fivethirtyeight')
import numpy as npdef standard_units(any_numbers):
"Convert any array of numbers to standard units."
return (any_numbers - np.mean(any_numbers))/np.std(any_numbers)
def correlation(t, x, y):
return np.mean(standard_units(t.column(x))*standard_units(t.column(y)))
def slope(table, x, y):
r = correlation(table, x, y)
return r * np.std(table.column(y))/np.std(table.column(x))
def intercept(table, x, y):
a = slope(table, x, y)
return np.mean(table.column(y)) - a * np.mean(table.column(x))
def fit(table, x, y):
"""Return the height of the regression line at each x value."""
a = slope(table, x, y)
b = intercept(table, x, y)
return a * table.column(x) + bshotput = Table.read_table(path_data + 'shotput.csv')shotputWeight Lifted | Shot Put Distance
37.5 | 6.4
51.5 | 10.2
61.3 | 12.4
61.3 | 13
63.6 | 13.2
66.1 | 13
70 | 12.7
92.7 | 13.9
90.5 | 15.5
90.5 | 15.8
... (18 rows omitted)shotput.scatter('Weight Lifted')Scatterplot with 'Weight Lifted' on the x-axis and 'Shot Put Distance' on the y-axis. There is a clear positive correlation, with more data points that have higher x values.这不是一个橄榄球形的散点图。事实上，它似乎有轻微的非线性成分。但如果我们坚持使用直线进行预测，仍然存在一条在所有直线中最佳的选择。
我们为橄榄球形散点图推导的回归线斜率和截距公式给出了以下值。
slope(shotput, 'Weight Lifted', 'Shot Put Distance')0.09834382159781997intercept(shotput, 'Weight Lifted', 'Shot Put Distance')5.959629098373952即使散点图不是橄榄球形的，使用这些公式是否仍然合理？我们可以通过找到最小化 mse 的直线的斜率和截距来回答这个问题。
我们将定义函数 shotput_linear_mse，它接受任意斜率和截距作为参数，并返回相应的 mse。然后对 shotput_linear_mse 应用 minimize 将返回最佳的斜率和截距。
def shotput_linear_mse(any_slope, any_intercept):
x = shotput.column('Weight Lifted')
y = shotput.column('Shot Put Distance')
fitted = any_slope*x + any_intercept
return np.mean((y - fitted) ** 2)minimize(shotput_linear_mse)array([0.09834382, 5.95962911])这些值与使用我们的公式得到的值相同。总结如下：
无论散点图的形状如何，都存在唯一一条最小化估计均方误差的直线。它被称为回归线，其斜率和截距由下式给出：
$$ \mathbf{\mbox{slope of the regression line}} ~=~ r \cdot \frac{\mbox{SD of }y}{\mbox{SD of }x} $$
$$ \mathbf{\mbox{intercept of the regression line}} ~=~ \mbox{average of }y ~-~ \mbox{slope} \cdot \mbox{average of }x $$
fitted = fit(shotput, 'Weight Lifted', 'Shot Put Distance')
shotput.with_column('Best Straight Line', fitted).scatter('Weight Lifted')The previous graph is reproduced with the addition of gold data points labeled 'Best Straight Line.' There is one data point per x value of the dark blue data points. We see that the gold points if connected would make a line with a positive slope. One notable deviation of the gold data points from the dark blue data points is for the smallest x value; it has an actual y value of just over 6 but a predicted y value of just under 10.非线性回归
上图加强了我们之前的观察，即散点图有点弯曲。因此，拟合一条曲线比拟合一条直线更好。该研究假设举起重量与铅球距离之间存在二次关系。那么，让我们使用二次函数作为预测变量，看看能否找到最佳的二次函数。
我们必须在所有二次函数中找到最佳二次函数，而不是在所有直线中找到最佳直线。最小二乘法使我们能够做到这一点。
这个最小化的数学原理很复杂，仅通过检查散点图不容易看出。但数值最小化与线性预测变量时一样简单！我们可以再次使用 minimize 来获得最佳的二次预测变量。让我们看看它是如何工作的。
回顾一下，二次函数的形式为
$$ f(x) ~=~ ax^2 + bx + c $$
其中 $a$、$b$、$c$ 为常数。
为了根据最小二乘准则找到基于举起重量预测距离的最佳二次函数，我们首先编写一个函数，该函数接受三个常数作为参数，使用上述二次函数计算拟合值，然后返回均方误差。
该函数称为 shotput_quadratic_mse。注意其定义类似于 lw_mse，只是拟合值基于二次函数而非线性函数。
def shotput_quadratic_mse(a, b, c):
x = shotput.column('Weight Lifted')
y = shotput.column('Shot Put Distance')
fitted = a*(x**2) + b*x + c
return np.mean((y - fitted) ** 2)现在我们可以像之前一样使用 minimize 来找到最小化均方误差的常数。
best = minimize(shotput_quadratic_mse)
bestarray([-1.04004838e-03, 2.82708045e-01, -1.53182115e+00])我们对举起 $x$ 公斤的运动员的铅球距离的预测约为
$$ -0.00104x^2 ~+~ 0.2827x - 1.5318 $$
米。例如，如果运动员可以举起 100 公斤，预测距离为 16.33 米。在散点图上，这位于 100 公斤附近垂直条带的中心附近。
(-0.00104)*(100**2) + 0.2827*100 - 1.531816.3382以下是所有 Weight Lifted 值的预测值。你可以看到它们大致穿过了散点图的中心。
x = shotput.column(0)
shotput_fit = best.item(0)*(x**2) + best.item(1)*x + best.item(2)shotput.with_column('Best Quadratic Curve', shotput_fit).scatter(0)The previous graph is reproduced with different gold data points labeled 'Best Quadratic Line.' There is one data point per x value of the dark blue data points. We see that the gold points if connected would make a curve that increases rapidly until it begins to even out on the right hand side of the graph. The gold data points better follow the dark blue data points, including the data point with the smallest x value.注意： 此处拟合二次曲线是因为原始研究中提出了这一建议。但值得注意的是，在图的最右端，二次曲线似乎接近峰值，之后曲线将开始下降。因此，对于能够举起远高于我们数据集中重量值的新运动员，我们可能不希望使用此模型。