第十六讲 残差分析

变量变换
• 线性模型假设 E(y|x)的线性性, 和误差方差齐性： E(y|x)=a+b’x var(y|x)=常数 • 我们知道（y，x）联合正态分布时，该假设是正确的。 • 实际问题中，若x，y都是连续变量，通常对x或y或两者做 变换，使得变换之后(x,y)近似服从正态分布。 • 若某些自变量是因子（x1)，其它自变量(x2)是连续型，那 么, 变换，使得(y,x2)|x1～正态。
变换的一般原则
• 总的原则是变换后每个变量都比较对称、均衡，换言之联合分布接近 正态。变换包括
– – – – Log变换 或 Box-Cox变换， 连续变量离散化， 有次序的因子变量的连续化， 无次序因子变量的合并
•
log 原则 如果一个非负变量的取值不在一个尺度或量级(magnitude)上,则取对数 后分析可能是有益的. 如果一个变量的取值在一个尺度或量级内,任何变换都可能无益. 不容易确定何种变换时，采用Box-Cox变换。
(b) 若x, z不独立 ⇒ var( y | x)一般依赖于x，除非 var( z | x)不依赖于x ⇒ var( y | x) = 常数 若( z , x) ~ 正态, 则
(3) δ = y-{α + β ' x} = a + b' x + c' z − {α + β ' x} = a − α + (b − β )' x + c' z
δ与x是否相关可通过 (1), (2)部分地检查，但一般无 法完全验证。
注：通常我们只有“工作模型”而不知道完全模型，只能通过 工作模型的残差探讨其拟合好坏。
残差分析
• 残差分析：拟合线性回归模型之后，通过 分析残差特征，检查拟合的好坏，即检查 数据是否满足模型假设
(i) 响应变量均值函数是否是线性 (ii) 误差方差是否为常数 (iii) 误差项是否与自变量不相关
即使都是基于条件分布即线性模型的统计推断行的但有时未明确写在给定自变量条件下进的期望和方差时都是在计算向量模型中是随机向量机变量在矩阵在总体模型中是一个随比如区分模型有些符号未特意模型多数时候是样本有时线性模型是总体定其具体含义交视具体情况而表示独立不相关正符号我们不加区分地使用
第十六讲. 残差分析
残差
解决方案： 变量变换。
(iii)误差项与自变量相关：通常表现为非线性或/和误差方差 非常数 例：真实模型：y=a+bx+cz+error, z=ux2, u~poisson(1) E(z|x)=x2 nonlinear in x, var(z|x)=x4 拟合：y~x
解决方案： 变量变换或增加 高阶项。
残差图： 纵坐标为残差，横坐标为拟合值（或每一个自变量）
模型正确时，残差与拟合值（或每一个自变量）不相关。
(i) 响应变量均值函数非线性（误差项与自变量相关） 例：真实模型：y=1+x^2+error 拟合：y=a+bx+error
残差
解决方案： 增加x的高阶项 或做变量变换。
(ii) 误差方差不是常数（误差项与自变量相关） 例：真实模型：y=a+bx+cz+error, z=x*u, u~poisson(1) E(z|x)=ux linear in x, var(z|x)=x2 拟合：y~x
• Income是人均收入，为了可比，将驾照人 数(Drivers)和汽油消耗总数(FuelC)人均 化： Drivers= Drivers/Pop FuelC= FuelC/Pop MPC=MPC/Pop • 对Drivers，Income， Miles，MPC作log变 换
> summary(m1) Call: lm(formula = FuelC ~ log(Drivers) + log(Income) + log(Miles) + log(MPC) + Tax, data = fuel01)
= E ( a + b' x + c ' z | x ) = a + b' x + c ' E ( z | x )
( a) 若z与x独立 ⇒ E ( z | x) = E(z) = 常数 E ( y|x) = (a + c' E ( z )) + b' x = ˆ α + β ' x，β = b
(b) 若z与不x独立 ⇒ E ( z | x)与 x有关,
Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 0.470488 0.771341 0.610 0.544956 log(Drivers) 0.247722 0.112181 2.208 0.032365 * log(Income) -0.014599 0.070636 -0.207 0.837190 log(Miles) 0.062485 0.012793 4.884 1.35e-05 *** log(MPC) 0.046964 0.012819 3.664 0.000652 *** Tax -0.005165 0.001804 -2.863 0.006353 ** --Residual standard error: 0.05717 on 45 degrees of freedom Multiple R-squared: 0.6282, Adjusted R-squared: 0.5869 F-statistic: 15.21 on 5 and 45 DF, p-value: 9.641e-09
对两个变量取对数, LS拟合得回归直线 : log( BrainWt ) = 2.14 + 0.75 log( BodyWt ) ⇔ BrainWt = 8.5 × BodyWt 3/4
回忆第七讲 - 幂次律 很多有效的有机系统服从3/4幂次律： - 加油站数目与城市规模 - 脑重量与体重
注意： 回归方程非常显著(p值 < 10-15 ). 如果数据来自于各种动物个体, 而不是平均, 那么 结果不会如此显著 （ecological correlation)
Boxcox变换
变换后的残差图较好 （注意y轴尺度不同）
例：汽车汽油消耗量与税率的关系
• alr3中的数据fuel2001收集了美国2001年51 个州的汽油消耗量数据，变量包括
Drivers: Number of Licensed drivers in the state FuelC: Gasoline sold for road use (1000s of gal.) Income: Per capita personal income (year 2000) Miles: Miles of Federal-aid highway miles in the state MPC: Estimated miles driven per capita Pop: Population age 16 and over Tax: Gasoline state tax rate, cents per gallon State:State name
评价工作模型是否正确，需要检查 1）E( y|x) = x的线性函数, 且β = b？ （ （2）var( y|x) = 常数？ (3) δ = y-(α + β ' x)与x不相关？
（1）E( y|x) = x的线性函数, 且β = b？
E ( y|x) = E (E ( y | x, z ) | x )
BodyWt (kg) 3.385 0.48 1.35 464.983 36.328 27.66 14.831 1.04 4.19 0.425 0.101 0.92 1 0.005 0.06 3.5 2 1.7 2547.07 0.023
a=lm(BrainWt~BodyWt, data=brains) plot(brains[,1], resid(a), xlab="BrainWt", ylab="residual")
若( z, x) ~ 正态, 或 E ( z | x) = ξ + ω ' x为x的线性函数，则 E ( y|x) = (a + c' ξ ) + (b + ωc)' x = ˆ α + β ' x 仍为线性, 但 α = a + c ' ξ , β = b + ωc ≠ b
若E ( z | x)为x的非线性函数, 则 E ( y|x)不是x的线性函数;
残差
残差图无法检测出问题的一种情况是: 未包括进模型y~x 的变量z与x联合正态（线性相关），此时 x的相关系数估 计的有偏性无法通过残差图发现.
完全模型： y = a + b' x + c ' z + ε , ε ~ （0，σ 2），ε ⊥ ( x, z ),
工作模型： y = α + β ' x + δ ,δ ~ （0，τ 2）
错误指定模型【从总体模型考虑】
完全模型： y = a + b' x + c ' z + ε , ε ~ （0，σ 2），ε ⊥ ( x, z ), ⇔ E( y|x,z ) = a + b' x + c' z, var( y | x, z ) = σ 2 目标是估计、推断 b。
工作模型： y = α + β ' x + δ ,δ ~ （0，τ 2），δ ⊥ x, ⇔ E(y | x) = α + β ' x, var( y | x) = τ 2
所以E ( y | x) = α + β ' x, var( y | x) = var(c' ε 1 + ε ) = 常数，且δ = c' ε 1 + ε ⊥ x 即工作模型满足（ 1）、（2） , 可以通过残差图的检查， 但其斜率β = b + ϕ ' c ≠ b这一事实无法通过残差图检测出来。
这一事实说明, 即使通过各种手段(包括变量变换, 引入连续变量的高阶项,因子变量的交互作用等 等), 得到良好的残差图即良好的拟合效果,也不 能保证回归系数的无偏性. 唯一的策略是数据的采集过程中保证与y,x都有 关的因素都采集到并包括进回归模型.
假设 （z , x) ~ 正态：E ( z | x) = c + dx, var( z | x) = 常数
则完全模型 y = a + b' x + c' z + ε = a + b' x + c' (φ + ϕx + ε 1 ) + ε = (a + c'φ )+ (b + ϕ ' c)' x + (c' ε 1 + ε ) = ˆ α + β'x +δ
（2）var( y|x) = 常数？
var( y | x) = E (var( y | x, z ) | x) + var(E ( y | x, z ) | x)
= E (σ 2 | x) + var(a + b' x + c' z | x) = σ 2 + c' var( z | x)c
(a) 若x, z独立 ⇒ var( z | x) = 常数 ⇒ var( y | x) = 常数，不依赖于x
一般地,只要λ不是远离0, 取log变换.
boxcox(BrainWt~BodyWt, data=brains)
λ ≈ 0, 对BrainWt做 log变换
boxcox(BodyWt~BrainWt, data=brains)
λ ≈ 0, 对BodyWt做 log变换
例：各国生育率与人均GDP关系
Box-Cox变换
u > 0, 连续变量 Box - Cox变换为 : u → u ( λ ) = u λ -1
λ
当λ = 0时, u ( λ ) = log(u)
λ由极大似然准则确定,
使用 library(MASS)中的函数boxcox可确定最优的λ :
确定响应变量y的最优变换 : boxcox( y ~ x + z + ...) 确定自变量x的最优变换 : boxcox( x ~ y + z + ...)
例1. 动物的脑重量与体重. 62 种动物的平均数据，图示如下:
BrainWt (g) Arctic_fox Owl_monkey Beaver Cow Gray_wolf Goat Roe_deer Guinea_pig Vervet Chinchilla Ground_squirrel Arctic_ground_squirrel African_giant_pouched_rat Lesser_short-tailed_shrew Star-nosed_mole Nine-banded_armadillo Tree_hyrax N._American_opossum Asian_elephant Big_brown_bat 44.5 15.499 8.1 423.012 119.498 114.996 98.199 5.5 57.998 6.4 4 5.7 6.6 0.14 1 10.8 12.3 6.3 4603.17 0.3
• 变量： PPgdp：人均GDP Fertility：平均生育个数 Purban：城市人口比率
PPgdp Fertility Purban Algeria 1784 2.80 58 Armenia 687 1.15 67 Australia 18788 1.70 91 Austria 23260 1.28 67 Azerbaijan 695 2.10 52 …