LTS回归与M估计稳健性的比较
作者：王彤 何大卫
单位：山西医科大学卫生统计教研室(030001)
关键词：稳健回归；失效点；LTS
中国卫生统计990205 【提 要】 目的 探讨可以同时抵抗回归数据中出现的Y与X空间的多个异常点的稳健回归估计方法。方法 介绍了M估计与LTS估计，引入失效点的概念以讨论其稳健性尺度，并对两个模拟数据进行了分析比较。结果 LTS估计优于M估计，失效点较高。结论 LTS估计可同时抵抗X与Y方面的异常，具有极高的稳健性。
Comparison of M-estimator and LTS Estimator in Regression.
Wang Tong,He Dawei.
Shan Xi Medical University(030001),Taiyuan
, 百拇医药
【Abstract】 Objective We will explore some robust regression estimates to resist multiple outliers from both sides.Methods We introduce M estimator and LTS estimator,with breakdown point to discuss their robustness.Also,we compare two simulated data as demonstration.Results LTS estimator is superior to M estimator,and has higher breakdown point.Conclusion LTS regression can resist outliers from both sides,and it is more robust.
【Key words】 Robust regression,Breakdown point,LTS
, http://www.100md.com
在线性回归分析中，使用经典的最小二乘法(LS估计)求解具有许多优良性质，但此时一般都要求模型满足GM假设。在大量的医学回归实例中，许多数据并不一定满足既定的假设，于是就可以通过“残差分析”或“影响分析”对数据进行回归诊断，用以识别其中的离群点或影响点。诊断之后对识别出的异常点应审慎处理，一般不宜简单地剔除(除非经核查证实为过失误差造成)，如果此时出于专业上的要求，需将此数据拟合既定模型，则应采用一些可以抵抗异常点影响作用的稳健估计方法〔1〕。本文将对M估计及LTS估计两种稳健回归技术作一比较，并通过两个模拟例子说明其各自的稳健性程度，以期促进稳健回归技术在医学领域中的应用。
稳健性尺度及M估计
为评价某种估计方法抵抗异常值扰动作用的能力，在Hodges的研究工作基础上，Hampel提出了失效点BP(Break-down Point)的概念。由于起初它是一个渐近结果，不便计算，故而Donoho和Huber定义了对于有限样本的BP为〔2〕：
, http://www.100md.com
其中T为参数估计值，n为样本X的含量，β(m;T,X)表示将X中任意替换m个点后两次参数估计值之差的上确界。这样，有限样本的BP就表示用该法作参数估计时，允许数据中存在的不会使估计值失效的异常点数目的最小比例，当数据中的异常点比例超过BP时，估计值会变得很不稳定。某法的BP值越高，其稳健性越好，可容忍的异常点数目就越多，而且这里的稳健性对于回归来说，应同时不受X与Y两方面异常的影响。
由以上定义可知，LS估计的BP值为
，大样本情况下趋近于零，其估计原则是使各点残差平方和最小，故而它很不稳健。如果用增长速度较慢的函数ρ代替平方和函数，就可减轻异常点的影响。如Huber的M估计，其目标函数为：
Q(β)=Σρ(Yi-Xiβ)
使其极小化，求导得到估计方程：Σψ(ri)Xi=0。式中ri=Yi-Xiβ，ψ为ρ的导函数。为使M估计满足尺度同变性(Scale equivariant)，应将ri除以一个稳健的尺度估计值s，然后将标准化残差ri/s代入估计方程迭代求解。s通常取值为“绝对离差中位数”MAD(Median absolute deviation)〔3〕,其值为1.4826×Sm,Sm=med(|ri-med(ri)|)，med表示中位数，1.4826=1/Φ-1(0.75)是假定误差正态分布情况下的一致性校正因子。
, http://www.100md.com
由于目标函数增长较慢，M估计比LS估计要稳健一些，但它只能抵抗Y方面异常，而仍然会受X方面异常点(即杠杆点)的影响，故而严格地说，其失效点仍为1/n。
LTS估计
为了从X与Y两方面同时抵抗异常点的扰动，下面介绍一种更为稳健的LTS(Least Trimmed Squares)回归估计。它的目标函数不是使所有各点的残差平方和最小，而是定义为使升序排列在前一半的残差平方和取最小时的估计值，这样就可使LTS估计对排序靠后的残差值较大的异常点稳健。其解可通过以下无放回的重复抽样算法来实现：
1.确定重复抽样次数SN。
2.从n个观察值中随机抽取P+1个不同的样本点(假定回归方程不通过原点，P为X变量个数)，代入模型，通过解P+1元方程组，得到一个回归参数向量估计值b。
, http://www.100md.com
3.将b代回方程，计算所有n个点的残差平方和
(i=1,2,…，n)，由小到大排序为(R2)1<(R2)2<…<(R2)n，并对排在前一半位置的各点残差平方求和，得到
。
4.令j=j+1，重复以上步骤至SN=j，在所得到的SN个ST中取最小值，将此时对应的b作为β的LTS估计。
这样得到的LST估计可以同时抵抗X与Y两方面异常的影响，且其BP=50%，达到了很高的稳健性。下面通过模拟数据来比较这两个方法。
模拟数据结果
, http://www.100md.com
以下模拟数据由SAS程序产生，M估计与LTS估计用SAS中的“交互式矩阵语言”(IML)实现。
1.表1数据中前20个观测由模型Y=8-2.1X+e产生，其中X服从正态分布N(2.8，1.852)，e服从正态分布N(0，0.22)，后边加入的5个观测为Y方向上的异常值，它们服从均数为(2.82，19.50)，标准差为0.2的双变量正态分布。
表1 Y方向异常点的模拟数据 序号
Y
X
序号
Y
X
, 百拇医药 序号
Y
X
1
5.7141
0.96428
10
0.0114
3.82292
19
5.7791
1.20331
2
, 百拇医药
6.1544
0.87827
11
1.3081
3.17510
20
2.8934
2.52784
3
-0.3644
3.96026
12
6.2338
, 百拇医药
0.90553
21
19.4144
2.75822
4
5.1700
1.30844
13
-1.7811
4.70899
22
19.6156
2.56923
, 百拇医药
5
9.6299
-0.94079
14
0.5898
3.60698
23
19.5580
2.80399
6
6.3929
0.69015
15
, 百拇医药
-2.9446
5.27164
24
19.3509
3.22539
7
7.0465
0.46168
16
0.4342
3.58672
25
19.4290
, 百拇医药
2.72924
8
2.3733
2.68818
17
3.1189
2.47573
9
3.9577
1.94393
18
-4.1691
5.84596
, http://www.100md.com
表2为对此数据LS，M及LTS的估计值，其中“基本点模拟模型”是指没有异常点干扰的前20对基本数据组成的模型。表2 Y方向异常数据的几种估计方法比较 回归参数
LS估计
M估计
LTS估计
基本点
模拟模型
α
10.3512
7.9882
8.0128
8.0000
, http://www.100md.com
β
-1.6442
-2.0514
-2.0930
-2.1000
此结果说明，当数据中仅存在Y方面异常点时，LS估计会受到很大影响，M估计可在一定程度上抵抗其影响，而LTS估计非常稳健。
2.表3数据中前20个观测由模型Y=8-2.1X+e产生，其中X服从正态分布N(2.8，2.852)，e服从正态分布N(0，0.22)，后边加入的5个观测在X与Y方向上均为异常，它们服从均数为(10.20，19.50)，标准差为0.2的双变量正态分布。表3 X与Y方向异常点的模拟数据 序号
, http://www.100md.com
Y
X
序号
Y
X
序号
Y
X
1
13.4357
-2.6180
10
2.2621
, 百拇医药 2.9511
19
-0.0145
3.7567
2
7.3808
0.2061
11
8.1841
-0.0436
20
-4.8288
6.0448
, 百拇医药
3
-9.0345
8.2139
12
-0.7296
4.2653
21
19.8492
10.3570
4
9.1248
-0.4925
13
, 百拇医药
-2.2529
5.1118
22
19.4543
10.5364
5
-3.9024
5.5795
14
12.4265
-2.2128
23
19.4951
, http://www.100md.com
10.1587
6
-1.9019
4.8249
15
-11.4241
9.3575
24
19.3384
9.8128
7
13.3114
-2.5200
, http://www.100md.com
16
-7.1410
7.2220
25
19.2779
10.0949
8
-4.3471
5.9405
17
-6.0822
6.7026
9
, http://www.100md.com
4.6383
1.8593
18
4.3972
1.7845
表4为对此数据LS，M及LTS的估计值，“基本点模拟模型”含义同上。表4 X与Y方向异常数据的几种估计方法比较 回归参数
LS估计
M估计
LTS估计
基本点
模拟模型
, 百拇医药 α
5.0780
5.0780
8.0048
8.0000
β
-0.0516
-0.0516
-2.0744
-2.1000
此结果说明，当数据在X与Y方面均出现异常点时，LS估计会受到很大影响，M估计此时也很不稳健，而LTS估计却不受其影响，结果非常稳健。讨 论
, http://www.100md.com
1.在LTS估计的重复抽样算法中，对抽样数SN的确定可以有两种考虑。其一是当n和P都不太大时，可穷尽一切抽样组合，取SN=CPn，此时得到的是LTS估计的唯一解；其二是当n和P比较大时，为减小计算时间的开销，可从一种概率角度考虑：从SN个样本中得到至少一个只包括稳定观察值(其中无异常点)的样本的概率为：1-[1-(1-e)P+1]SN，e为观察值中可能存在的异常点的比例，当取最大值e=50%，令此概率接近1(如95%)时，即可求得所需的抽样数SN，此时得到的是LTS估计的近似解。本文中的估计结果是在抽样数仅为100时得到的，近似已经很好。
2.在描述复杂生命现象的医学回归数据中，可能出现的异常点数目是未知的，我们往往不能先验地判断其可能的数目，也不好预知是在Y空间或X空间有异常。当X方向存在异常值时，经典的M估计就会失效，故而作稳健回归分析时应同时考虑抵抗两方面的异常。从失效点的角度考虑，该法的BP值要尽可能的高。由以上模拟数据分析结果可以看到，LTS估计可同时抵抗X与Y方向多个异常点的影响，是一种非常稳健的回归估计技术。
, http://www.100md.com
3.LTS回归的误差估计与M估计中的s类似，只是其中的SM由排序在前一半的残差平方和开根来代替，在此基础上就可以进一步进行有关的统计推断。另外，与其他一些用在简单线性回归分析中的稳健方法(如Theil回归)比较，LTS估计可用于多元回归过程〔4〕。
4.与LTS很类似的还有LMS(Least Median of Squares)估计，二者都满足“仿射等积性”(Affine Equivariat)，失效点都高达50%，同样可用于解决多个异常点诊断中复杂的“掩盖”(Masking)问题〔5〕，但对于稳健估计而言，LTS估计的“渐近有效性”(Assymptotic Efficency)及收敛速度都较LMS估计更好。有关结果笔者将另文给出。
(本研究得到SAS公司赞助，特此致谢)
*本题目为全国统计科学研究项目，山西省归国人员基金项目
, http://www.100md.com
参考文献
1.Myers,R.H.Classical and modern regression analysis with applications(2nd ed.).Boston:PWS-Kent.1990,128.
2.Lopuhaa H.P.Highly efficiency estimators of multivariate location with high breakdown point.The Annals of Statistics,1992,20(1):398
3.Joseph W.M et al.The use and interpretation of residuals based on robust estimation.The Journal of American Statistical Association,1993,88(424):1254
4.Rousseeuw P.J & Leroy A.M.Robust regression and outlier detection.Toronto:John Wiley & Sons 1987,43
5.王彤，何大卫.线性回归中多个异常点的稳健诊断及其医学应用.中国卫生统计1998，15(5)：1, http://www.100md.com
单位：山西医科大学卫生统计教研室(030001)
关键词：稳健回归；失效点；LTS
中国卫生统计990205 【提 要】 目的 探讨可以同时抵抗回归数据中出现的Y与X空间的多个异常点的稳健回归估计方法。方法 介绍了M估计与LTS估计，引入失效点的概念以讨论其稳健性尺度，并对两个模拟数据进行了分析比较。结果 LTS估计优于M估计，失效点较高。结论 LTS估计可同时抵抗X与Y方面的异常，具有极高的稳健性。
Comparison of M-estimator and LTS Estimator in Regression.
Wang Tong,He Dawei.
Shan Xi Medical University(030001),Taiyuan
, 百拇医药
【Abstract】 Objective We will explore some robust regression estimates to resist multiple outliers from both sides.Methods We introduce M estimator and LTS estimator,with breakdown point to discuss their robustness.Also,we compare two simulated data as demonstration.Results LTS estimator is superior to M estimator,and has higher breakdown point.Conclusion LTS regression can resist outliers from both sides,and it is more robust.
【Key words】 Robust regression,Breakdown point,LTS
, http://www.100md.com
在线性回归分析中，使用经典的最小二乘法(LS估计)求解具有许多优良性质，但此时一般都要求模型满足GM假设。在大量的医学回归实例中，许多数据并不一定满足既定的假设，于是就可以通过“残差分析”或“影响分析”对数据进行回归诊断，用以识别其中的离群点或影响点。诊断之后对识别出的异常点应审慎处理，一般不宜简单地剔除(除非经核查证实为过失误差造成)，如果此时出于专业上的要求，需将此数据拟合既定模型，则应采用一些可以抵抗异常点影响作用的稳健估计方法〔1〕。本文将对M估计及LTS估计两种稳健回归技术作一比较，并通过两个模拟例子说明其各自的稳健性程度，以期促进稳健回归技术在医学领域中的应用。
稳健性尺度及M估计
为评价某种估计方法抵抗异常值扰动作用的能力，在Hodges的研究工作基础上，Hampel提出了失效点BP(Break-down Point)的概念。由于起初它是一个渐近结果，不便计算，故而Donoho和Huber定义了对于有限样本的BP为〔2〕：
, http://www.100md.com
其中T为参数估计值，n为样本X的含量，β(m;T,X)表示将X中任意替换m个点后两次参数估计值之差的上确界。这样，有限样本的BP就表示用该法作参数估计时，允许数据中存在的不会使估计值失效的异常点数目的最小比例，当数据中的异常点比例超过BP时，估计值会变得很不稳定。某法的BP值越高，其稳健性越好，可容忍的异常点数目就越多，而且这里的稳健性对于回归来说，应同时不受X与Y两方面异常的影响。
由以上定义可知，LS估计的BP值为
，大样本情况下趋近于零，其估计原则是使各点残差平方和最小，故而它很不稳健。如果用增长速度较慢的函数ρ代替平方和函数，就可减轻异常点的影响。如Huber的M估计，其目标函数为：Q(β)=Σρ(Yi-Xiβ)
使其极小化，求导得到估计方程：Σψ(ri)Xi=0。式中ri=Yi-Xiβ，ψ为ρ的导函数。为使M估计满足尺度同变性(Scale equivariant)，应将ri除以一个稳健的尺度估计值s，然后将标准化残差ri/s代入估计方程迭代求解。s通常取值为“绝对离差中位数”MAD(Median absolute deviation)〔3〕,其值为1.4826×Sm,Sm=med(|ri-med(ri)|)，med表示中位数，1.4826=1/Φ-1(0.75)是假定误差正态分布情况下的一致性校正因子。
, http://www.100md.com
由于目标函数增长较慢，M估计比LS估计要稳健一些，但它只能抵抗Y方面异常，而仍然会受X方面异常点(即杠杆点)的影响，故而严格地说，其失效点仍为1/n。
LTS估计
为了从X与Y两方面同时抵抗异常点的扰动，下面介绍一种更为稳健的LTS(Least Trimmed Squares)回归估计。它的目标函数不是使所有各点的残差平方和最小，而是定义为使升序排列在前一半的残差平方和取最小时的估计值，这样就可使LTS估计对排序靠后的残差值较大的异常点稳健。其解可通过以下无放回的重复抽样算法来实现：
1.确定重复抽样次数SN。
2.从n个观察值中随机抽取P+1个不同的样本点(假定回归方程不通过原点，P为X变量个数)，代入模型，通过解P+1元方程组，得到一个回归参数向量估计值b。
, http://www.100md.com
3.将b代回方程，计算所有n个点的残差平方和
(i=1,2,…，n)，由小到大排序为(R2)1<(R2)2<…<(R2)n，并对排在前一半位置的各点残差平方求和，得到。4.令j=j+1，重复以上步骤至SN=j，在所得到的SN个ST中取最小值，将此时对应的b作为β的LTS估计。
这样得到的LST估计可以同时抵抗X与Y两方面异常的影响，且其BP=50%，达到了很高的稳健性。下面通过模拟数据来比较这两个方法。
模拟数据结果
, http://www.100md.com
以下模拟数据由SAS程序产生，M估计与LTS估计用SAS中的“交互式矩阵语言”(IML)实现。
1.表1数据中前20个观测由模型Y=8-2.1X+e产生，其中X服从正态分布N(2.8，1.852)，e服从正态分布N(0，0.22)，后边加入的5个观测为Y方向上的异常值，它们服从均数为(2.82，19.50)，标准差为0.2的双变量正态分布。
表1 Y方向异常点的模拟数据 序号
Y
X
序号
Y
X
, 百拇医药 序号
Y
X
1
5.7141
0.96428
10
0.0114
3.82292
19
5.7791
1.20331
2
, 百拇医药
6.1544
0.87827
11
1.3081
3.17510
20
2.8934
2.52784
3
-0.3644
3.96026
12
6.2338
, 百拇医药
0.90553
21
19.4144
2.75822
4
5.1700
1.30844
13
-1.7811
4.70899
22
19.6156
2.56923
, 百拇医药
5
9.6299
-0.94079
14
0.5898
3.60698
23
19.5580
2.80399
6
6.3929
0.69015
15
, 百拇医药
-2.9446
5.27164
24
19.3509
3.22539
7
7.0465
0.46168
16
0.4342
3.58672
25
19.4290
, 百拇医药
2.72924
8
2.3733
2.68818
17
3.1189
2.47573
9
3.9577
1.94393
18
-4.1691
5.84596
, http://www.100md.com
表2为对此数据LS，M及LTS的估计值，其中“基本点模拟模型”是指没有异常点干扰的前20对基本数据组成的模型。表2 Y方向异常数据的几种估计方法比较 回归参数
LS估计
M估计
LTS估计
基本点
模拟模型
α
10.3512
7.9882
8.0128
8.0000
, http://www.100md.com
β
-1.6442
-2.0514
-2.0930
-2.1000
此结果说明，当数据中仅存在Y方面异常点时，LS估计会受到很大影响，M估计可在一定程度上抵抗其影响，而LTS估计非常稳健。
2.表3数据中前20个观测由模型Y=8-2.1X+e产生，其中X服从正态分布N(2.8，2.852)，e服从正态分布N(0，0.22)，后边加入的5个观测在X与Y方向上均为异常，它们服从均数为(10.20，19.50)，标准差为0.2的双变量正态分布。表3 X与Y方向异常点的模拟数据 序号
, http://www.100md.com
Y
X
序号
Y
X
序号
Y
X
1
13.4357
-2.6180
10
2.2621
, 百拇医药 2.9511
19
-0.0145
3.7567
2
7.3808
0.2061
11
8.1841
-0.0436
20
-4.8288
6.0448
, 百拇医药
3
-9.0345
8.2139
12
-0.7296
4.2653
21
19.8492
10.3570
4
9.1248
-0.4925
13
, 百拇医药
-2.2529
5.1118
22
19.4543
10.5364
5
-3.9024
5.5795
14
12.4265
-2.2128
23
19.4951
, http://www.100md.com
10.1587
6
-1.9019
4.8249
15
-11.4241
9.3575
24
19.3384
9.8128
7
13.3114
-2.5200
, http://www.100md.com
16
-7.1410
7.2220
25
19.2779
10.0949
8
-4.3471
5.9405
17
-6.0822
6.7026
9
, http://www.100md.com
4.6383
1.8593
18
4.3972
1.7845
表4为对此数据LS，M及LTS的估计值，“基本点模拟模型”含义同上。表4 X与Y方向异常数据的几种估计方法比较 回归参数
LS估计
M估计
LTS估计
基本点
模拟模型
, 百拇医药 α
5.0780
5.0780
8.0048
8.0000
β
-0.0516
-0.0516
-2.0744
-2.1000
此结果说明，当数据在X与Y方面均出现异常点时，LS估计会受到很大影响，M估计此时也很不稳健，而LTS估计却不受其影响，结果非常稳健。讨 论
, http://www.100md.com
1.在LTS估计的重复抽样算法中，对抽样数SN的确定可以有两种考虑。其一是当n和P都不太大时，可穷尽一切抽样组合，取SN=CPn，此时得到的是LTS估计的唯一解；其二是当n和P比较大时，为减小计算时间的开销，可从一种概率角度考虑：从SN个样本中得到至少一个只包括稳定观察值(其中无异常点)的样本的概率为：1-[1-(1-e)P+1]SN，e为观察值中可能存在的异常点的比例，当取最大值e=50%，令此概率接近1(如95%)时，即可求得所需的抽样数SN，此时得到的是LTS估计的近似解。本文中的估计结果是在抽样数仅为100时得到的，近似已经很好。
2.在描述复杂生命现象的医学回归数据中，可能出现的异常点数目是未知的，我们往往不能先验地判断其可能的数目，也不好预知是在Y空间或X空间有异常。当X方向存在异常值时，经典的M估计就会失效，故而作稳健回归分析时应同时考虑抵抗两方面的异常。从失效点的角度考虑，该法的BP值要尽可能的高。由以上模拟数据分析结果可以看到，LTS估计可同时抵抗X与Y方向多个异常点的影响，是一种非常稳健的回归估计技术。
, http://www.100md.com
3.LTS回归的误差估计与M估计中的s类似，只是其中的SM由排序在前一半的残差平方和开根来代替，在此基础上就可以进一步进行有关的统计推断。另外，与其他一些用在简单线性回归分析中的稳健方法(如Theil回归)比较，LTS估计可用于多元回归过程〔4〕。
4.与LTS很类似的还有LMS(Least Median of Squares)估计，二者都满足“仿射等积性”(Affine Equivariat)，失效点都高达50%，同样可用于解决多个异常点诊断中复杂的“掩盖”(Masking)问题〔5〕，但对于稳健估计而言，LTS估计的“渐近有效性”(Assymptotic Efficency)及收敛速度都较LMS估计更好。有关结果笔者将另文给出。
(本研究得到SAS公司赞助，特此致谢)
*本题目为全国统计科学研究项目，山西省归国人员基金项目
, http://www.100md.com
参考文献
1.Myers,R.H.Classical and modern regression analysis with applications(2nd ed.).Boston:PWS-Kent.1990,128.
2.Lopuhaa H.P.Highly efficiency estimators of multivariate location with high breakdown point.The Annals of Statistics,1992,20(1):398
3.Joseph W.M et al.The use and interpretation of residuals based on robust estimation.The Journal of American Statistical Association,1993,88(424):1254
4.Rousseeuw P.J & Leroy A.M.Robust regression and outlier detection.Toronto:John Wiley & Sons 1987,43
5.王彤，何大卫.线性回归中多个异常点的稳健诊断及其医学应用.中国卫生统计1998，15(5)：1, http://www.100md.com