数值诊断
from datascience import *
path_data = '../../../assets/data/'
import numpy as np
from scipy import stats
import matplotlib
matplotlib.use('Agg')
%matplotlib inline
import matplotlib.pyplot as plots
plots.style.use('fivethirtyeight')除了可视化之外，我们还可以使用残差的数值属性来评估回归的质量。我们不会从数学上证明这些性质，而是通过计算来观察它们，并看看它们能告诉我们关于回归的哪些信息。
下面列出的所有事实对所有形状的散点图都成立，无论它们是否为线性。
family_heights = Table.read_table(path_data + 'family_heights.csv')
heights = family_heights.select('midparentHeight', 'childHeight')
heights = heights.relabel(0, 'MidParent').relabel(1, 'Child')
dugong = Table.read_table(path_data + 'dugongs.csv')
dugong = dugong.move_to_start('Length')
hybrid = Table.read_table(path_data + 'hybrid.csv')def standard_units(x):
return (x - np.mean(x))/np.std(x)
def correlation(table, x, y):
x_in_standard_units = standard_units(table.column(x))
y_in_standard_units = standard_units(table.column(y))
return np.mean(x_in_standard_units * y_in_standard_units)
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):
a = slope(table, x, y)
b = intercept(table, x, y)
return a * table.column(x) + b
def residual(table, x, y):
return table.column(y) - fit(table, x, y)
def scatter_fit(table, x, y):
table.scatter(x, y, s=15)
plots.plot(table.column(x), fit(table, x, y), lw=4, color='gold')
plots.xlabel(x)
plots.ylabel(y)
def residual_plot(table, x, y):
x_array = table.column(x)
t = Table().with_columns(
x, x_array,
'residuals', residual(table, x, y)
)
t.scatter(x, 'residuals', color='r')
xlims = make_array(min(x_array), max(x_array))
plots.plot(xlims, make_array(0, 0), color='darkblue', lw=4)
plots.title('Residual Plot')
def regression_diagnostic_plots(table, x, y):
scatter_fit(table, x, y)
residual_plot(table, x, y) heights = heights.with_columns(
'Fitted Value', fit(heights, 'MidParent', 'Child'),
'Residual', residual(heights, 'MidParent', 'Child')
)残差图显示无趋势
对于每一个线性回归，无论好坏，残差图都不显示趋势。总体上是平坦的。换句话说，残差与预测变量不相关。
你可以在上述所有残差图中看到这一点。我们还可以计算每种情况下预测变量与残差之间的相关性。
correlation(heights, 'MidParent', 'Residual')-2.719689807647064e-16这看起来不像零，但这是一个非常小的数字，除了由于计算的舍入误差外，它等于 0。这里再次显示，精确到 10 位小数。负号是由于上述的舍入造成的。
round(correlation(heights, 'MidParent', 'Residual'), 10)-0.0dugong = dugong.with_columns(
'Fitted Value', fit(dugong, 'Length', 'Age'),
'Residual', residual(dugong, 'Length', 'Age')
)
round(correlation(dugong, 'Length', 'Residual'), 10)0.0残差的平均值
无论散点图的形状如何，残差的平均值都是 0。
这类似于以下事实：如果你取任意一组数字，计算它们与平均值的偏差，这些偏差的平均值为 0。
在上述所有残差图中，你都看到了穿过图中心的水平线 0。这就是这个事实的可视化体现。
作为数值示例，以下是在基于父母中位身高对子女身高进行回归时的残差平均值。
round(np.mean(heights.column('Residual')), 10)0.0儒艮年龄对体长的回归中残差的平均值也是如此。除舍入误差外，残差的均值为 0。
round(np.mean(dugong.column('Residual')), 10)0.0残差的标准差
无论散点图的形状如何，残差的标准差都是响应变量标准差的一个分数。这个分数是 $\sqrt{1-r^2}$。
$$ \mbox{SD of residuals} ~=~ \sqrt{1 - r^2} \cdot \mbox{SD of }y $$
我们很快就会看到这如何度量回归估计的准确性。但首先，让我们通过例子来确认它。
在子女身高和父母中位身高的例子中，残差的标准差约为 3.39 英寸。
np.std(heights.column('Residual'))3.3880799163953426这与 $\sqrt{1-r^2}$ 乘以响应变量的标准差相同：
r = correlation(heights, 'MidParent', 'Child')
np.sqrt(1 - r**2) * np.std(heights.column('Child'))3.388079916395342对于混合动力汽车里程对加速度的回归也是如此。相关 $r$ 为负（约 -0.5），但 $r^2$ 为正，因此 $\sqrt{1-r^2}$ 是一个分数。
r = correlation(hybrid, 'acceleration', 'mpg')
r-0.5060703843771186hybrid = hybrid.with_columns(
'fitted mpg', fit(hybrid, 'acceleration', 'mpg'),
'residual', residual(hybrid, 'acceleration', 'mpg')
)
np.std(hybrid.column('residual')), np.sqrt(1 - r**2)*np.std(hybrid.column('mpg'))(9.43273683343029, 9.43273683343029)现在让我们看看残差的标准差如何衡量回归的好坏。记住残差的平均值为 0。因此，残差的标准差越小，残差越接近 0。换句话说，如果残差的标准差很小，那么回归中误差的总体大小就很小。
极端情况是当 $r=1$ 或 $r=-1$ 时。在这两种情况下，$\sqrt{1-r^2} = 0$。因此残差的平均值和标准差都为 0，因此所有残差都等于 0。回归线完美地完成了估计。正如我们在本章前面看到的，如果 $r = \pm 1$，散点图就是一条完美的直线，与回归线相同，因此回归估计中确实没有误差。
但通常 $r$ 不在极端值。如果 $r$ 既不是 $\pm 1$ 也不是 0，那么 $\sqrt{1-r^2}$ 是一个真分数，回归估计误差的大致总体大小介于 0 和 $y$ 的标准差之间。
最坏的情况是当 $r = 0$ 时。此时 $\sqrt{1-r^2} = 1$，残差的标准差等于 $y$ 的标准差。这与以下观察一致：如果 $r=0$，则回归线是 $y$ 平均值处的一条水平线。在这种情况下，回归的均方根误差就是与 $y$ 平均值的均方根偏差，即 $y$ 的标准差。在实际应用中，如果 $r = 0$，则两个变量之间没有线性关联，因此使用线性回归没有好处。
解释 $r$ 的另一种方式
我们可以将上述结果改写为：无论散点图的形状如何，
$$ \frac{\mbox{SD of residuals}}{\mbox{SD of }y} ~=~ \sqrt{1-r^2} $$
一个互补的结果是：无论散点图的形状如何，拟合值的标准差是 $y$ 观测值标准差的一个分数。这个分数是 $\vert r \vert$。
$$ \frac{\mbox{SD of fitted values}}{\mbox{SD of }y} ~=~ \vert r \vert $$
要理解这个分数的来源，注意拟合值都位于回归线上，而 $y$ 的观测值是散点图中所有点的纵坐标，变化更大。
scatter_fit(heights, 'MidParent', 'Child')Scatterplot with 'MidParent' on the x-axis and 'Child' on the y-axis. The data points are in dark blue and have a weak to moderate positive association. A gold line is drawn through the center of the data with a positive slope.拟合值的范围大约从 64 到 71，而所有子女的身高变化要大得多，范围大约从 55 到 80。
为了从数值上验证这个结果，我们只需计算等式两边。
correlation(heights, 'MidParent', 'Child')0.32094989606395924以下是拟合值的标准差与出生体重观测值的标准差的比率：
np.std(heights.column('Fitted Value'))/np.std(heights.column('Child'))0.32094989606395957该比率等于 $r$，证实了我们的结果。
绝对值是从哪里来的？首先注意，标准差不能为负，标准差的比率也不能为负。那么当 $r$ 为负时会发生什么？燃油效率和加速度的例子将告诉我们。
correlation(hybrid, 'acceleration', 'mpg')-0.5060703843771186np.std(hybrid.column('fitted mpg'))/np.std(hybrid.column('mpg'))0.5060703843771186两个标准差的比率为 $\vert r \vert$。
表达这个结果的一个更标准的方式是回顾：
$$ \mbox{variance} ~=~ \mbox{mean squared deviation from average} ~=~ \mbox{SD}^2 $$
因此，将我们结果的两边平方，
$$ \frac{\mbox{variance of fitted values}}{\mbox{variance of }y} ~=~ r^2 $$