回归建模的基础与要领(Ⅳ)——样品状态与相互间关系

胡良平

(1.军事科学院研究生院，北京 100850；2.世界中医药学会联合会临床科研统计学专业委员会，北京 100029*通信作者：胡良平，E-mail：lphu812@sina.com)

1 概 述

回归分析是研究因变量如何依赖自变量变化而变化的规律的重要统计分析方法之一，然而，回归分析的基本要素涉及两个方面，其一，变量状态及相互关系；其二，样品(测定变量取值的对象)状态及相互关系。因篇幅所限，本文仅讨论前述的“第二个要素”。

2 样品状态

2.1 单个体型样品

通常，适合于采用回归分析的数据结构中的每个“样品”对应着“一个个体”。若受试对象或调查对象是“人”，则每个“人”在统计学上被称为一个“样品”。例如，从30例某病患者的血液样品中测得“载脂蛋白A1、载脂蛋白B、载脂蛋白E、载脂蛋白C、低密度脂蛋白中胆固醇”的含量[1]，其数据结构见表1。

表1 30例患者载脂蛋白和低密度脂蛋白中胆固醇含量的测量结果

又例如，调查某地区某时间段内685例年龄≥70岁老年人的“一般情况、与健康有关的各项指标的取值和生活质量”所得到的资料[2]。显然，若受试对象或调查对象是某种“动物”，则每只“动物”在统计学上也被称为一个“样品”。由此可知，以每个“样品”为观察单位时，就可以称其为“单个体型样品”。

2.2 多个体型样品

有时，观察单位不是一个个体，而是由具有相近条件的多个个体组成。例如，某试验研究药物剂量与有效率之间的关系，数据结构如表2所示[1]。

表2 药物不同剂量与有效数

在表2中，每个剂量组的“全部动物”被视为一个“观察单位”，共有10个剂量组。显然，每个“观察单位”有多只动物。

又例如，某棉纺厂为减轻试验工作量，拟用较易测定的每毫克重纤维的根数x估计测定工作量较大的原棉单纤维强力y。研究者收集到的试验资料见表3。其中，mi为在第i个试验点xi上进行的独立重复试验次数；yi实际上是第i个试验点上mi个yij(j=1,2,…,mi)的算术平均值(注：若未求平均值，就可求出方差)[3]。

表3 每毫克重纤维的根数x与原棉单纤维强力y之间关系的测定结果

在表3中，各“编号”代表一批试验或被称为一个“观察单位”，若对其进行回归分析，各“编号”对应的数据被称为“样品”。则每个“编号”由重复试验次数不等的多个样品构成，被称为“多个体型样品”。

3 样品间相互关系

3.1 概述

与变量间相互关系相比，样品间相互关系比较难理解，因为样品间关系需要借助“几何图形”呈现出来才便于直观判断。通常，可通过在二维直角坐标系中的全部(x，y)散点分布情况，用目测法得出全部“样品”或“试验点”间实际存在的相互关系(简称为“几何方法”)；然而，当自变量数目≥2时，要在高维空间中直接呈现全部样品间相互关系非常困难。

解决前述困难的办法是：在二维空间中，找到合适的统计处理方法(简称“代数方法”)，从而建立起“几何方法”与“代数方法”之间的联系。由此，可将“代数方法”推广到高维空间中去研究“样品”间的相互关系。

根据数学理论和实践结果，上面提及的“代数方法”可归结为给“每个样品”一个“权重系数”，它的作用是反映每个样品在计算中的“分量”或“作用大小”。也就是说，“权重系数”大的“样品”要比“权重系数”小的“样品”发挥更大的作用。对于同一个回归分析资料，选取不同的依据来构造“权重系数”并据此来构建回归模型，其精确度是不同的。因此，可将全部可能的“依据”都用来构造“权重系数”，从而可构建出多种不同的回归模型。于是，可从中选出“最精准的回归模型”。

3.2 样品的同质性与异常点

3.2.1 样品的同质性

在对计量资料进行t检验或方差分析时，统计学教科书上都会明确交代：资料必须满足“独立性”“正态性”和“方差齐性”三个前提条件；而在对资料进行相关与回归分析，尤其是进行多重回归分析时，统计学教科书上则很少提及极其重要的“前提条件”，即所有样品对于全部变量应满足“同质性”。其含义是：所有样品或个体在全部变量上的“取值规律”是基本相同的。例如：研究某地区某时间段内正常成年人的体重是如何随身高变化而变化的依赖关系时，当所有被观测个体(或称为“样品”)在(身高，体重)两个变量上的取值对应的“数据点”沿一条直线(或曲线)变化趋势随机地散布，没有偏离“绝大多数样品”所在“区域”较远或很远的“数据点”，就称该资料中的所有“个体或样品”具有较好的“同质性”。

3.2.2 异常点

在前面的(身高，体重)例子中，若其中包含了少数特体型个体(例如，身高约为2.3 m，但体重约为50 kg；体重约为250 kg，但身高仅1.6 m；身高约为1.0 m，但体重为80 kg)，那么，这少数特体型人与绝大多数正常成年人就不是“同质的”。于是，那些“特体型个体”在统计学上就被称为“异常点(即异常的个体)”。之所以说它们是“异常点”，是因为当采用“几何方法”呈现时，它们所处的“空间位置”会偏离其他“数据点”所在的“变化区域”。在二维直角坐标系中，绘制出资料的散布图，数据点的分布情况将一览无余，见图1。

图1 四组数据的散布图

在图1中，从左上角至右下角的对角线上有数据点较粗的4幅小图，其中，最后两幅小图中均各有一个“异常点”。

在进行直线回归分析时，若存在“异常点”，但分析者对其视而不见，就很容易得出错误的结果和结论；在进行多重回归分析时，若存在“异常点”，分析者也同样容易“误入歧途”。在SAS的“REG过程”中，可以通过“学生化残差”和“Cook’s D距离统计量”来进行“异常点诊断”，淘汰掉资料中“严重的异常点”，将有助于提高回归模型的拟合质量。当一个资料中存在较大比例的“异常点”且又不适合将它们全部删除时，需要找到对“异常点”有一定“耐受性”的回归建模方法，常称为“稳健回归分析法”[4]。经过比较发现：“分位数回归分析法”[5]比“参数法中的诸多稳健回归分析法”[4]更加“稳健”。

3.3 可用于构造“权重系数”的常见“依据”[3，6]

3.3.1 某变量在各样品上取值的平方的倒数

选取某个变量(因变量或自变量)，其在每个“样品”上会有不同的取值。若在大多数样品上的取值比较接近，而在少数样品上的取值“非常大”，那么，就这个“变量”而言，可能意味着：取值“非常大”的那几个“样品”很可能是“异常点”。若取该“变量平方的倒数”为“权重系数”，则取值“非常大”的“样品”的“权重系数”就比较“小”，从而，它们在回归系数估计中所起的作用就相应地变小了。也就是说，这实际上是间接削弱了“可疑异常点”的影响，使回归系数的估计趋于“稳健”。

3.3.2 因变量在各样品上残差平方的倒数

也可以这样做：先不盲目地寻找任何“依据”来构建“权重系数”，而采取通常的方法构建回归模型，再利用此回归模型计算出各“样品”上因变量的“预测值”；进而可计算出各“样品”上的“残差”。于是，可以求出各“样品”上“残差平方的倒数”。估计回归系数时，以各“样品”上“残差平方的倒数”为“权重系数”。道理如前所述，此处从略。

3.3.3 某变量在每个样品上全部取值的方差的倒数

若在某个变量取每个特定值的条件下都进行了多次重复试验，就会获得因变量的多个观测值，于是，就可计算出多个因变量观测值的方差，进而可计算出“因变量方差的倒数”。若选取各样品上“因变量方差的倒数”为“权重系数”，则“方差大”的样品上的“权重系数”就很小，故它们在回归系数估计中所起的作用就相对变小了。

3.3.4 各样品上重复试验次数的倒数

当各样品上有次数不等的重复试验时，可以取“重复试验次数的倒数”为“权重系数”。因为重复试验次数较多的样品上“因变量”的不同观测值的数目可能就会多一些，也就间接反映了该样品上“因变量的方差可能比较大”。于是，取“重复试验次数的倒数”为“权重系数”，相当于取“因变量方差的倒数”为“权重系数”。道理同上，此处不再赘述。

3.3.5基于“观测权重”与“抽样权重”构造“综合权重”[7]

对抽样调查资料进行回归建模时，选取合适的“权重系数”至关重要。这里可能涉及到多个有关的概念，如观测权重、抽样权重。

观测权重是基于综合评价中权重系数的思想，在回归分析中引入反映各个体或观测对总体的重要性的度量，表示在其他观测不变的情况下，该观测的变化对结果的影响程度。常用的有经验权重法、试验次数权重法、贡献权重法等。

抽样权重是在抽样研究中，为反映所抽取样本中各个观测在总体中的重要程度，或样本中各个观测代表总体中个体的数目。抽样权重的大小与抽样方法有关，分为基础抽样权重、调整抽样权重与总抽样权重。例如，将某省或州划分为3个地区(即“层”)，各层总样本量、抽取的样本量和抽样权重的计算方法和结果见表4。

表4 美国两个州各地区(层)中农场的数目、抽样数目和抽样权重

表4中，抽样权重=农场数目/抽样数目，例如，第1行上被抽取的3个农场中的每一个代表全部农场(100个)中的1/3，即33.333个。

综合权重是在对随机抽样所得数据进行统计分析时，不仅考虑抽样权重，还考虑观测权重，计算各个观测对结果的总的重要程度。其计算方法是：综合权重=观测权重×抽样权重。

在SAS/STAT的“SURVEREG过程”[8]中，若采用分层随机抽样，则选择“抽样权重”作为“权重系数”。

3.4 实例演示

3.4.1 无重复试验的回归分析问题

3.4.1.1 问题与数据结构

【例1】某公司对12次投标情况进行研究。设投标规模为x(单位：百万美元)，企业准备投标的费用为y(单位：千美元)。具体数据见表5。试建立y关于x的回归方程[4]。

表5 某地某年投标规模x与企业准备投标费用y的数据

注：表4数据摘自文献《应用线性回归模型》(约翰·内特, 威廉·沃塞曼, 迈克尔·H·库特纳, 著, 张勇, 王国明, 赵秀珍, 译. 北京: 中国统计出版社, 1990: 178)

若绘制出(x，y)的散布图，各散点随自变量x的增加，y的离散度也变大，为节省篇幅，绘制散布图的SAS程序和散布图均省略。

下面，先不考虑“权重系数”拟合直线回归模型，然后再选取不同的“依据”构建“权重系数”，并据此构建直线回归模型。基于模型的拟合优度(R2、误差等)确定最合适的“权重系数”。

3.4.1.2 所需要的SAS程序

data a1;

inputxy@@;

/*以下产生三个新变量，分别代表不同的权重系数*/

wx=1/(x**2);

wy=1/(y**2);

wxy=1/(x*y);

cards;

2.1315.5 1.2111.111.0062.6 6.0035.4 5.6024.9 6.9128.1 2.9715.0 3.3523.210.3942.0 1.1010.0 4.3620.0 8.0047.5

;

run;

/*以下程序采用普通最小平方法1创建含截距项的直线回归模型*/

proc reg data=a1;

modely=x/r;

quit;

/*以下程序采用普通最小平方法2创建不含截距项的直线回归模型*/

proc reg data=a1;

modely=x/ nointr;

quit;

/*以下程序采用加权最小平方法1创建直线回归模型*/

proc reg data=a1;

modely=x/r;

weightwx;

quit;

/*以下程序采用加权最小平方法2创建直线回归模型*/

proc reg data=a1;

modely=x/r;

weightwy;

quit;

/*以下程序采用加权最小平方法3创建直线回归模型*/

proc reg data=a1;

modely=x/r;

weightwxy;

quit;

/*以下程序采用普通最小平方法创建直线回归模型，提取各样品点上的残差*/

proc reg data=a1 noprint;

modely=x/ nointr;

output out=aaa residual=resid;

quit;

/*以下程序为求取各样品点上残差平方的倒数*/

data a2;

set aaa;

wr=1/resid**2;

run;

/*以下程序采用加权最小平方法4创建含截距项的直线回归模型*/

proc reg data=a2;

modely=x/r;

weightwr;

quit;

/*以下程序采用加权最小平方法5创建不含截距项的直线回归模型*/

proc reg data=a2;

modely=x/ nointr;

weightwr;

quit;

【SAS程序说明】

以上SAS程序很长，各段SAS程序之前都有“注释语句”，这些注释语句解释了其后面程序的作用，此处不再赘述。

3.4.1.3 主要计算结果汇总

“普通最小平方法”和“加权最小平方法”拟合直线回归模型的参数估计结果见表6。

表6 普通与加权最小平方法拟合直线回归模型的参数估计值等内容比较

在表6中，“普通最小平方法1”对应的结果中，截距项无统计学意义；“普通最小平方法2”中就没有包含截距项。

在表6中，“加权最小平方法”有5种，其中，前4种对应的“权重系数”分别为“自变量x的平方的倒数”“因变量y的平方的倒数”“自变量x与因变量y的乘积的倒数”“因变量y的残差平方的倒数”，而“加权最小平方法5”与“加权最小平方法4”的“权重系数”相同，都是“因变量y的残差平方的倒数”，它们的区别在于“是否保留截距项”。

淘汰掉“普通最小平方法1”和“加权最小平方法4”的结果之后，还有5种方法对应的结果，分别为“普通最小平方法2”和“加权最小平方法1”“加权最小平方法2”“加权最小平方法3”“加权最小平方法5”。那么，这5种方法对应的结果哪一个相对更好？

若从假设检验的P值来看，很难分辨出孰优孰劣。但可以依据参数“标准误”大小进行比较，标准误小者为好。由此可知，“加权最小平方法5”给出的斜率的“标准误0.05057”最小，故该法相对其他4种更好。

下面，再列出上述各种方法对应的“R2”“调整R2”“均方根误差”和“预测残差平方和，简称PRESS”，见表7。

表7 普通与加权最小平方法拟合直线回归模型的拟合优度内容比较

续表1：

加权最小平方法10.91510.90660.8896710.32205加权最小平方法20.89970.88970.164170.34635加权最小平方法30.90680.89750.383971.90498加权最小平方法40.90870.89961.04916109.89517加权最小平方法50.99890.99881.0028913.00709

由统计学基础知识可知：对于直线回归模型而言，“R2”和“调整R2”的数值越大越好；而“均方根误差”和“PRESS(预测残差平方和)”越小越好。由此可知，在表7中，最后一行的结果是最好的。

3.4.1.4 本例小结

结合表6和表7的结果以及比较得出的结论可知：本例以“因变量y的残差平方的倒数”构建“权重系数”，并采取“加权最小平方法”拟合直线回归模型且不含截距项(因为截距项无统计学意义，需删除)，其拟合效果最佳。

3.4.2 多重线性回归分析问题

3.4.2.1 问题与数据结构

【例2】沿用前面的“表1资料”，设y代表“低密度脂蛋白”，x1～x4分别代表表1中第2列至第5列上的4种“载脂蛋白”，试建立y依赖4个自变量的多重线性回归模型。

3.4.2.2 所需要的SAS程序

data a1;

input idx1-x4y@@;

wy=1/y**2;

cards;

11731067.014.713721391326.417.816231981126.916.713441181387.115.71885139948.613.6138617516012.120.3215713115411.221.517181581419.729.614891581377.418.2197101321517.517.2113111621106.015.91451214411310.142.881131621377.220.7185141691298.516.7157151291386.310.11971616614811.533.4156171851186.017.5156181551216.120.4154191751114.127.2144201361109.426.090211531338.516.9215221101499.524.718423160865.310.8118241121238.016.6127251471108.518.4137262041226.121.0126271311026.613.4130281701278.424.7135291731238.719.01883013213113.829.2122

;

run;

/*不加权且保留截距项，三种常规筛选方法所得结果相同，仅留一种*/

proc reg data=a1;

modely=x1-x4/selection=backward sls=0.05r;

title1'此处创建的是模型1';

quit;

/*不加权且不保留截距项，三种常规筛选方法所得结果相同，仅留一种*/

/*数据集aaa中包含各样品点上的残差resid1*/

proc reg data=a1;

modely=x1-x4/noint selection=backward sls=0.05r;

output out=aaa residual=resid1;

title1'此处创建的是模型2';

quit;

/*以下程序基于aaa数据集，用残差平方的倒数作为权重系数*/

/*这实际上就是做了一次加权多重线性回归分析*/

data a2;

set aaa;

wr1=1/resid1**2;

run;

proc reg data=a2;

modely=x1-x4/noint selection=backward sle=0.05r;

weightwr1;

title1 '此处创建的是模型3';

quit;

/*加权且保留截距项，三种常规筛选方法所得结果相同，仅留一种*/

proc reg data=a1;

modely=x1-x4/selection=backward sls=0.05r;

weightwy;

title1'此处创建的是模型4';

quit;

/*加权且不保留截距项，三种常规筛选方法所得结果相同，仅留一种*/

/*数据集bbb中包含各样品点上的残差resid2*/

proc reg data=a1;

modely=x1-x4/noint selection=backward sls=0.05r;

weightwy;

output out=bbb residual=resid2;

title1'此处创建的是模型5';

quit;

/*以下程序基于bbb数据集，用残差平方的倒数作为权重系数*/

/*这实际上就是做了第二次加权多重线性回归分析*/

data b2;

set bbb;

wr2=1/resid2**2;

run;

proc reg data=b2;

modely=x1-x4/noint selection=backward sle=0.05r;

weightwr2;

title1'此处创建的是模型6';

quit;

3.4.2.3 主要计算结果汇总

“普通最小平方法”和“加权最小平方法”拟合多重线性回归模型的参数估计结果见表8。

表8 普通与加权最小平方法拟合多重线性回归模型的参数估计值等内容比较

由表8可知：模型1和模型4均不够理想，因为它们都包含了无统计学意义的截距项；模型2与模型3具有可比性，但模型3中参数的标准误小于模型2中参数的标准误，稍好一些；同理，模型6比模型5更好。

那么模型3与模型6哪一个更好？为回答这个问题，需要列出与“拟合优度”有关统计量的计算结果，见表9。

表9 普通与加权最小平方法拟合多重线性回归模型的拟合优度内容比较

在表9中，NI代表模型中自变量的个数，Cp值越接近模型中自变量的个数，表明模型对资料的拟合度越好。模型3与模型6相比，R2和调整R2都比较接近；模型6的均方根误差小于模型3的均方根误差；特别是Cp值，说明模型6优于模型3。

3.4.2.4 本例小结

模型1和模型2都是基于普通最小平方法建模，前者保留截距项，后者不保留截距项；模型3仅采取残差平方的倒数为“权重系数”，进行了一次“加权最小平方法构建回归模型”；模型4和模型5都是基于因变量y平方的倒数为“权重系数”，进行了第一次“加权最小平方法构建回归模型”；而模型6在模型5的基础上，又基于残差平方的倒数为“权重系数”，进行了第二次“加权最小平方法构建回归模型”。最终的结论是：基于两次加权回归分析得到的模型6优于其他模型。其回归模型为：