第三十课Spearman等级相关分析一、秩相关的Spearman等级相关分析前面介绍了使用非参数方法比较总体的位置或刻度参数,我们同样也可以用非参数方法 比较两总体之间的相关问题秩相关(rank correlation)又称等级相关,它是一种分析和等级 间是否相关的方法适用于某些不能准确地测量指标值而只能以严重程度、名次先后、反应 大小等定出的等级资料,也适用于某些不呈正态分布或难于判断分布的资料设和分别为和各自在变量X和变量Y中的秩,如果变量X与变量Y之间存在着正相关, 那么X与Y应当是同时增加或减少,这种现象当然会反映在(,)相应的秩(,)上反之, 若(,)具有同步性,那么(,)的变化也具有同步性因此:d =工d 2 =X(R -Q )2 (3o.i)i i ii=1 i=1具有较小的数值如果变量X与变量Y之间存在着负相关,那么X与Y中一个增加时,另一 个在减小,具有较大的数值既然由(,)构成的样本相关系数反映了 X与Y之间相关与否的 信息,那么在参数相关系数的公式r(X,Y)中以和分别代替和,不是同样地反映了这种信息吗?基于这种想法,Charles Spearman秩相关系数r (R, Q)应运而生:sr (R, Q)=工(R - 1 工 R )(Q - 1 工 Q )』 n 』 』 n i(R -1 工R )2工(Q -1 工Q )2i n i i n i(30.2)r (R, Q)与r(X, Y)形式上完全一致,但在r (R, Q)中的秩,不管X与Y取值如何,总s s是只取1到之间的数值,因此它不涉及X与Y总体其他的内在性质,例如,秩相关不需要总 体具有有限两阶矩的要求。
由于:ii=1i=1 i 2》R =Hq = 1 + 2 + …+ n = n(n +1)Hr 2 =工 Q 2 = 12 + 22 + …+ n 2 = n(n + 1)(2n +1)i ii=1 i=1因此,公式(30.2)可以化简为:6工(R - Q )2 r = 1 — i 匚s n(n 2 -1)6工d 2i n(n 2 — 1)(30.3)显然在=时,秩相关系数达到最大值+1又因为:工(R - Q )2 =i i工R2 + 工Q2 -2工RQ =i i i in(n + 1)(2n +31) - 2 工 RQi i而工RQ在每对+= n +1时达到最小值,最小值求法为:i i工(n +1)2 =工 R 2 +Z Q 2 + 2 工 RQi i i i所以,最小的i in(n +1)2 n(n + 1)(2n +1)2 - 6最大的工(R 一 Q )2为:i i2n(n + 1)(2n +1) , .、 n(n2 -1)-n(n +1)2 =3 3故秩相关系数的最小值为1—2=—1在原假设和不相关的情况为真时,即秩相关系数为0时,的期望值为0,样本的方差为1 — r 2 s 2 = s (30.4)rs n — 2自由度为n - 2且分布关于零点对称。
当10时,的样本分布可以标准化为近似的t分布:t == XzL = r 口 〜t (n - 2)s 1 一 r 2 s 1 — r2 (30.5)r I o ss 、;n-2例30・1某公司想要知道是否职工期望成为好的销售员而实际上就能有好的销售记录为 了调查这个问题,公司的副总裁仔细地查看和评价了公司10个职工的初始面试摘要、学科成 绩、推荐信等材料,最后副总裁根据他们成功的潜能给出了单独的等级评分二年后获得了 实际的销售记录,得到了第二份等级评分,见表30.1中的第1到4列所示统计问题为是否 职工的销售潜能与开始二年的实际销售成绩一致表30.1职工的销售潜能与销售成绩的秩相关分析职工编号潜能等级销售成绩成绩等级d = R — Qi i i124001112436031137300524412956—525562807—11633504—11710200100089260811982209—11105385239工d 2 =i44Spearman秩相关系数r (R,Q)的计算过程见表30.1中的第5到6列所示,最后计算结s果为二 1 纟 二 1 - 6(44) 二 0.7333n(n2 -1) 10(100 -1)表明潜能与成绩之间是较强的正相关,高的潜能趋向于好的成绩。
秩相关系数r (R,Q)原假 s设为0的t检验统计量为:t = 0.7333 —二 3.051 - (0.7333)2查表自由度为8, t=3.05的双侧p=0.0158在0.05显著水平上,t分布的上临界点为2.30,由 于3.05>2.30,因此,拒绝秩相关系数为0的原假设,接受潜能与成绩之间存在秩相关二、Corr相关过程Corr相关过程用于计算变量之间的相关系数,包括Pearson(皮尔逊)的乘积矩相关和加 权乘积矩相关还能产生三个非参数的关联测量:Spearman的秩相关,Kendall的tau-b和 Hoeffding的相关性度量D该过程也可以计算偏相关等一些单变量的描述性统计量1. Corr过程说明proc co proccorrdata=数据集〈选项〉;var 变量列表;with 变量列表; partial变量列表; weight 变量; freq 变量; By 变量列表; run ; proc corr语句调用corr过程,且是唯一必需的语句如果只使用proc corr这一条语句,过程计算输入数据集中所有数值变量之间的相关系数其余语句是供选择的2. proc corr语句的选项• 0"切=数据集名 产生含有Pearson相关系数的一个新数据集。
• 0%霑=数据集名 产生含有Spearman等级相关系数的一个新数据集• 0%饮=数据集名——产生含有KendallT b相关系数的一个新数据集• outh-数据集名 产生含有Hoeffding D统计量的一个新数据集• pearson 求计算通常的pearson乘积矩相关系数,是缺省值• hoeffding 求计算并输出Hoeffding的D统计量• kendall——求计算并输出KendallT b相关系数b• spearman 求计算并输出Spearman等级相关系数• vardef=df^weigh1\wgtwdf- 指定计算方差时的除数:df (自由度n-1),weight或wgt (权重之和),n (观察数),wdf (权重之和一1)缺省值为dfC”——计算协方差一方差矩阵 sscp 求输出平方和与交叉积和• csscp 求输出偏差平方和与交叉积和• best=数^一对每个变量输出指定个数的绝对值最大的相关系数• noprint 止所有打印输出• noprob 止输出同这些相关有联系的显著性概率• nosimple 对原始数据执行标准方差分析• rank—— 求按绝对值从高到低的次序对每个变量输出相关系数。
• nocorr 制Pearson相关的计算及输出• nomiss—— 带有某一变量缺失值的观测值从所有计算中除去• nosimple 不输出每个变量的简单描述性统计量3. var语句该语句列出要计算相关系数的变量例如,var a b c;则计算a和b, a和c, b和c三对变 量之间的相关系数4. with语句为了得到变量间的特殊组合的相关系数,该语句和var语句联合使用用var语句列出的 变量在输出相关阵的上方,而用with语句列出的变量竖在相关阵左边例如,var a b;with x y z;则生成x和a, y和a, z和a, x和b, y和b, z和b5. partial 语句为了计算Pearson偏相关,Spearman偏秩相关,Kendal 1偏tau-b,用该语句给出偏出去(即 固定)的变量名6. weight 语句为了计算加权的乘积矩相关系数,用该语句给出权数变量名该语句仅用于Pearson相关7. freq语句当规定freq语句时,输入数据集中的每个观察假定代表n个观察,其中n是该观察中freq 变量中的值观察的总数规定为freq变量值的和8. by语句使用by语句能够获得用by变量定义的分组观察的独立分析结果。
三、实例分析例30・1的SAS程序如下:data study.persons ;input x y @@;y=400-y;2 400 4 360 7 300 1 295 6 2803 350 10 200 9 260 8 220 5 385proc corr data二study.persons spearman;var x;with y;run;程序说明:建立输入数据集persons,要注意实际数据所表示的等级次序大小与SAS系 统中自动给出的等级次序大小的不同输入变量x,获得从1到10的数据,表示潜能等级从 最高到最低,而输入变量y,获得从最大销售额400到最小销售额220,转换销售成绩等级应 该是从高到低,即从1到10但在SAS系统中把销售成绩数值从小到大按等级值从1到10 给予因此,需要颠倒变量x或变量y中数值大小的次序,本程序用最大销售额400减去原 来的销售额实现次序颠倒,即语句y=400—y等级相关与一般参数相关一样仍然调用corr过 程,只需要在选择项中指定为何种等级相关,我们选择计算spearman秩相关系数var语句 列出要计算相关系数的第一个变量x,with语句必须要与var语句联合使用,列出的要计算相 关系数的第二变量y。
主要结果如表30.2所示表30・2用corr过程进行多样本输出结果Correlation Analysis1 'WITH' Variables: Y1 'VAR' Variables: XSimple StatisticsVariableNMeanStd DevMedianMinimumMaximumY1095.00000067.905163102.5000000200.000000X105.5000003.0276505.5000001.00000010.000000Spearman Correlation Coefficients / Prob > |R| under Ho: Rho=0 / N = 10 XY 0.733330.0158结果说明:Spearman等级相关系数为0.73333,是一个比较大的正相关系数这个相关 系数为0的原假设检验结果是p=0.0158v0.05,因此,我们拒绝相关系数为0的原假设,接受 了这个0.73333等级相关系数结论为销售潜能的高低与销售成绩好坏之间存在明显的正相 关性第三十一课一元线性回归分析回归分析是一种统计分析方法,它利用两个或两个以上变量之间的关系,由一个或几个 变量来预测另一个变量。
在SAS/STAT中有多个进行回归的过程,如REG、GLM 等, REG过 程常用于进行一般线性回归模型分析四、回归模型1. 基本概念回归模型是一种正规工具,它表示统计关系中两个基本的内容:①用系统的形式表示因 变量Y随一个或几个自变量X变化的趋势;②表现观察值围绕统计关系曲线的散布情况这 两个特点是由下列假设决定的:•在与抽样过程相联系的观察值总体中,对应于每一个X值,存在Y的一个概率分布; 这些概率分布的均值以一些系统的方式随X变化•图31.1是用透视的方法来显示回归曲线Y对给定X具有概率分布这一概念总是 与统计关系中的经验分布形式上相对应;同样,描述概率分布的均值与X之间关系的回归 曲线,与统计关系中Y系统地随X变化的一般趋势相对应统计关系线20019018017018015014013012011010090807080504030201010 15 20 25 30 35 40 45 50 55 60 65 70 75 80 85 90批量X图31.1线性回归模型的图示在回归模型中,X称为“自变量”,Y称为“因变量”;这只是传统的称法,并不表明在 给定的情况下Y因果地依赖于X,无论统计关系多么密切,回归模型不一定是因果关系,在 某些应用中,比如我们由温度表水银柱高度(自变量)来估计温度(因变量)时,自变量实 际上依赖于因变量。
此外,回归模型的自变量可以多于一个2. 回归模型的构造1) 自变量的选择构造回归模型时必须考虑到易处理性,所以在有关的任何问题中,回归模型只能(或只 应该)包括有限个自变量或预测变量2) 回归方程的函数形式选择回归方程函数形式与选择自变量紧密相关有时有关理论可能指出适当的函数形式 然而,通常我们预先并不能知道回归方程的函数形式,要在收集和分析数据后,才能确定函 数形式我们经常使用线性和二次回归函数来作为未知性质回归方程的最初近似值图31. 2(a) 表示复杂回归函数可以由线性回归函数近似的情况,图31.2(b)表示复杂回归函数可以由两个 线性回归函数分段近似的情况图31.2用线性回归函数近似复杂回归函数3) 模型的范围在建立回归模型时,通常需要限制模型的自变量或因变量取值的区间范围,这个范围由 调查设计和已掌握数据的情况决定4) 回归分析的运用回归分析主要有3个目的:描述,控制和预测五、未指定误差项分布的回归模型1. 模型的正规表述现在我们只限于一个自变量的基本回归模型,且回归函数是线性的,可表述如下:Y = a + BX* + s * (31.1)其中,是第次观测或试验中因变量的取值,和是参数,为第次观测或试验中自变量的取 值,是随机误差项,其基本假设应该满足三个条件:• 均值E() = 0•方差 差Var (s ) =o2t• 协方差Cov(& , s )二0,当详j时。
即对所有的详j,与互不相关i j模型(31.1)称为简单模型,参数是线性的,自变量也是线性的所谓“简单”,是因 为它只有一个自变量,“参数线性”是指没有参数具有指数形式,或者被另一个参数相乘或 相除,“自变量线性”是指这个自变量是一次的参数和自变量都是线性的模型称为一阶模 型2. 模型的重要特点第次观察中的观察值包括2部分:常数项空墮]和随机项的和所以,是随机变量 因为E()=0,这样:E(Y)二 a + pX + E(& )二 a + pX (31.2)t t t t其中,吓阪]是常数因此,当第次试验中取为时,相应的来自一个概率分布,其均 值是:E(Y) =a + pX I (31.3)t t所以,模型(31.1 )的回归函数是:(31.4)E (Y) =a + px这样对任何给定的,回归函数把水平与的概率分布均值联系起来在第次试验中,的观察值超过或低于回归函数值的部分为误差项部分假设误差项具有 相同的方差,则相应的的方差为:Var(Y)二 b 2t这是因为:Var(Y)二 Var(a + PX +£ )二 Var(£ )二b2t t t t无论自变量取值如何,模型(31.1 )总是假设的概率分布具有相同的方差,且假设误差 项互不相关。
因此,任何一次试验的结果对其他各次试验的误差项都没有影响,相应的与也 互不相关总之,模型(31.1)的含义为:对所有水平的来说,因变量观察值都来自均值 E(Y) =a +阪]、方差的概率分布此外,任何两个观察值与是互不相关的六、■ =J最【小二乘估计法1.观测数据图设有一组T期间内关于二变量和的样本观测值(,)(t=1, 2,…,N),在和之间存 在着函数关系,如果将这些观测数据,在2维平面上用图来表示,只要数据至少有3个以上, 那么所有的点大概不可能都在一条直线上以被认为在X和Y之间成立的未知回归直线:Y=a+ pX为中心,观测点总是适当地散布在其周围未知回归直线和各观测点的垂直方向的间隔就是 上节引进的概率误差项由于a和卩的数值未知,因此,不能准确地知道与各观测点对应的 概率误差项的值大致来说,可以认为回归直线是从散布在平面上的各观测点的中央穿过的 直线根据所给的观测数据来估计这条直线的位置(a和p的值),是我们需要解决的主要问 题2.误差二乘和的1:小化估计回归直线的方式(规则)有各种各样的考虑但是,对于确定^和卩的值时,要使所有的观测点和直线的“距离”从整体来说为最小这个一般的规则,大概无论谁也没有异议。
意见的分歧在于究竟要用什么尺度来衡量各观测点和回归直线的“距离”也就是说,即使 都承认上述的一般规则,但由于按什么标准来测定“点和线的距离”的看法不同,推导出的 估计方式也是多种多样的假定估计出的直线为:Y = a * + p * X (31.5)则同X=X对应的估计直线上的点是J* +卩*x观测点(x, yt)同估计直线垂直方向的间 隔:(31.1.6)e 二 y - (a * + p * x )叫做残差(residual)这里将各观测点看作是已经观测完毕的一对已知数组,用小写字母 来表示)应当注意的是误差项和残差的区别:误差项是未知回归直线同观测点的间隔,而 残差是已知的估计直线同观测点的间隔为了便于讨论,我们暂且将测量点和直线之间距离的“评价函数”限定为残差的函数 对照我们的常识,要求评价函数满足以下各条件:1)残差可能为正也可能为负,但不管是正的残差还是负的残差,只要其绝对值相等, 用与直线的离差这一标准来衡量,就应当完全平等地评价2)评价函数必须是各残差绝对值的非减函数把评价函数记为V5霜,…,J将以上两条件用数学方式表现,可得:(31.7)av(31.8)同时,为了方便起见,除以上2个条件外,暂且再追加以下2个条件。
3) N个观测点都具有同等资格即和(tzs)作为评价函数的变量应得到同样的对待 这一条件同各期误差项的方差为一定值的假定有着密切的关系将条件(3)用数学方式表…,小的任意重新排列叽…r有:现,可得,对于(1, 2,4)我们已经假定时期不同的概率误差项相互之间不相关因此,评价函数中各的作用 最好是相互无关的将这一叙述用数学方式表示,可得:二 0, t 丰 s(31.9)a 2Vde det根据以上的讨论,备择的评价函数被限定在相当狭的范围内,作为满足资格的函数,例 如可以考虑:(31.10)当k为偶数时,绝对值的符号就失去意义残差是回归系数的估计值(a*, P* )的函数因此,如果给定了观测数据(x,y),则可以把V看作是以和为变量的二变量函数从而可t t以考虑确定能使V为最小的和的值当然,使V的值为最小的和的值要依存于N个观测数据 当k= 1时,评价函数式(31.10)是残差绝对值的总和就某种意义来说,这一评价函数 在直观上也许是最容易理解的通过使它为最小来确定和的方式,叫做最小绝对离差估计法 (least absolute deviation estimation method)当k=2时,评价函数是残差的平方和。
确定能使这一评价函数为最小的和的方式,便是 最小二乘法(least squares method)令 k= 2,将式(31.6)代入式(31.10),可得:V =乙0 —a*—卩*x 丄 (31.11)t tt^1把样本观测值看作已知数,从而可以把V当作和的函数来考虑,利用解决最大最小问题 的方法,令V对和的偏导数为零,可以推导出关于和的二元联立一次方程组为:6Vda—卩(31.12)dVxt^y —a *tt(31.13)这一联立方程叫做正规方程式,其解如下:迓(x — x)G — y )t t- 迓(x — x )2tt 1(31.14)(31.15)x _ 捍迓 xt ‘y _ N 迓 ytt_1 t_1(31.16)在求解时,利用了下列恒等式:迓(x — x)2 _Z^x 2 ——迓x t t N { ■t _1 t _1 ' t _1因为,V的驻点(使偏导数同时为0的和的值)只有唯一的一个,而且通过增大和的值, 可以使V无限增大,所以正规方程的解的确给出了V的最小值于是,可知最小二乘估计量 是:为(x — x )(y — y )t tP _迓(x — x )2t(31.17)(31.18)3. 最=Jt=1迓ytt=1(31.19)【小二乘估计量的平均值和方差我们已经相当详细地论述了关于“估计量的优劣”问题的一般理论。
从18世纪由高斯 (Gauss)发明的所谓最小二乘法直到今天仍得到如此广泛的实际运用这一事实来看,最小二 乘估计法理论应具有某些特别的优点如前所述,最小二乘法并不是“确定使T个观测点与 回归直线之间的距离就整体来说为最小的直线位置”的独一无二的方法,它只不过是多种方 法中的一个罢了尽管如此,最小二乘法还能够绝对地凌驾于其他任何方法之上,一直被应 用于现实数据的分析,这并不仅仅是由于计算简单,而且还有其他合适的理由一一理论上的 根据事实上,在计算技术有了非常大的进步的今天,计算简便已经不再具有那样大的价值 了以下,我们首先来研究一下最小二乘估计量的性质将|Y =a + BXt田代入估计量(31.17)和估计量(31.18),并作以下变形:(31.20)为(x 一 X、=P+ — l -迓(X - X)2tt=1t=1于是,和的期望值分别为:E Ce )=eX (x - x )t 为(x — x )2tt=1(31.21)(31.22)(31.23)从而很简单地证明了和分别是e和卩的无偏估计量这样,最小二乘估计量顺利地通过了 第一道关卡既然已表明最小二乘估计量具有无偏性,那么下一个问题就是估计量的方差的 大小。
我们暂且先根据方差的公式进行形式上的推导根据前面的假定:Var(8 ) =Q2和Cov(8 ,8 ) = 0,由定义得t t sVar(3 )= E (3 —卩)(31.24)按照同样的方法也可以推导出:Var心)=E(d -a》=c21(31.25)+N兰(x - X》tt=1这里顺便再计算一下和的协方差:Cov(31.26)C,0)= Ea-a)k 一卩)=-XC2 兰(x -x》tt=1从式(31.25)和式(31.26)可知,估计量的方差与样本的大小大致成反比同时,解释变 量在较广的范围内分布得越散,估计量的方差就越小估计量的方差越小即意味着估计值的精度越高当lim迓(x - X》=8时,和都是一致估计量 N ” t=1 t七、检验与预测从最小二乘估计表达式(31.17)和(31.18)知,只要给出了 N组数据(x , y ), i = 1,2,…,N,i i总可将它们代入这两个表达式获得a和卩的估计,从而写出回归方程但这个回归方程是否有 意义呢?需要有个检验准则为作检验,首先要建立假设我们求回归方程的目的是要去反 映随变化的一种统计规律,那么如果0=0,从式(31.4)可知,不管如何变化,不会随之而改变, 在这种情况下求出的回归方程是无意义的。
所以,检验回归方程是否有意义的问题转化为检 验下列假设是否为真:H °: 0 = 0 (31.27)常用的方法有F检验和t检验方法1. F检验这一方法类似于第三章所介绍的方差分析的想法,也是从观察值的偏差平方和分解入手我们观察到的y ,y ,y的差异可以用总偏差平方和表示:1 2 NTSS =X (y - y)2, df = N -1 (31.28)i Ti=1造成这一差异的原因有如下两个方面:一是由于假设0 =0不真,从而对不同的值,随而变化我们可以用下列偏差平方和来表示由此引起的差异:RSS =迟(y - y)2,df = 1 (31.29)i R称为回归平方和其中,y = & + [3x = y - 0x + [3xi i i=y + 0(x — x)所以,公式(31.29)i又可以写成:RSS =迟(y — y )2i(31.30)(31.31)i=1=艺[0 (x — x )]2 = 0 2 艺(x — x)2i i ii=1 i=1根据公式(31.24)可知,其期望值:E(RSS) = E02 •迓(x — x)2ii=1=[(E0 )2 + Var (0 )]迓(x — x )2ii=1=02(x - x)2 + ◎ 2ii=1这便表明,RSS中除了误差波动外,还反映了由于0所引起的数据间的差异。
二是由其他一切随机因素引起的差异,它可以用残差平方和:ESS 仝(y — y )2,df = N — 2 (31.31)i i Ei=1表示由于可以证明:ESS/g2 〜咒2(N — 2) (31.32)于是有:E(ESS) = (N — 2)g 2 (31.33)所以,其自由度为N-2利用公式工(y — y ) = 0,工(y — y )x = 0,从而有下列平方和分解式:i i i i iTSS =工(y - y)2 =工(y - y + y - y)21 1 1 1(31.34)=Z(y - y )2+乙(y - y)2i i i=ESS + RSS由于在0= 0为真时,RSS与ESS /(N — 2)都是的无偏估计,因而采用F统计量:RSS / g 2/1ESS / g 2/( N — 2)RSSESS /(N — 2)〜F (1, N — 2)(31.35)来检验原假设0 =0是否为真2. t检验由公式(31.24)和(31.25)知(31.36)1 x 2a 〜N(a,a2[ + — ])N 乙(———)2i(31.37)在原假设卩=0为真时,- 〜N(0,1),但其中未知,常用& 2 = ESS/(N — 2)CT / J》(X — X)2i去代替,根据公式(31.32)和又与独立,从而在卩二0时有:0 a / 込(x — x£t = = i 〜t (N — 2)a / L(Xi— x)2 :空/(N — 2)a 2实质上,对于一元回归方程t检验与F检验是等价的,因为只要将公式(31.30)中的RSS(31.38)代入到公式(31.35沖去,就不难发现t2二F。
我们同样可以得到原假设Q二0为真时的t统计量:a a L1/N + X 2 /L (x — X)Tt = = i 〜t (N — 2)& J1/N + X2 / 乙(x — X)2i:空 /(N — 2)a 2(31.39)3.利用回归方程作预测当求得回归方程y二a+併后,并经检验,方程是显著的,则可将该回归方程用于预测所谓预测是指当取某一个具体值时,对相应的取值所作的推断由模型知y 0 =a + PX0 +s是一个随机变量,要预测随机变量的取值是不可能的,只能预测其期望值E(yo)根据前面公式(31.24)、(31.25)和(31.26)可知,在x = X处的回归值是亍a +札,且:y o~ N(E (y o),Var (y 0))(31.40)其中:E (y ) = a + px0 0(31.41)Var(Y0)=X 2C 2—X )2t2T(31.42)1 (x X) 2其中,未知,用&2 = ESS /(N _ 2)去代替,设杠杆率h =诂+ —0 N 工(x — X)2t,所以预测均值的预测区间为(—t 时 + 5)0 01/2 k 0 0 a/2 0(31.43)其中,,o/2的自由度为N - 2。
注意在SAS系统model语句中的elm选项是按公式(31.43) 来计算的然而在x = x0时,随机变量的取值与预测均值总会有一定的偏离,我们根据公式(31.43) 不难求出j _ y的均值E(j j )和方差Var(y _ Y ),且它符合正态分布,故有:0 0 0 0 0 0y 0 _ y 0~ N(°,(x — x )20—工(x — x )2t其中,未知,用c2 = ESS/(N — 2)去代替,所以y — y的预测区间为:( 11)b 0 - y 0) - 八(1+ 讣 2,( y 0- y 0) + 八(1+ 讣丿(31.44)(31.45)其中,to/2的自由度为N - 2注意在SAS系统model语句中的eli选项是按公式(31.44) 来计算的从方差Var(y — y )表达式中我们可以看到,当取值离均值越近,预测精度就越好,当 0 0取值离均值越远,预测精度就越差,其预测区间两头呈喇叭状因此,我们要特别注意取值 应该在样本数据最小的和最大的之间,否则预测很不可靠八、回归诊断回归诊断主要用于检验关于回归假设是否成立,以及检验模型形式是否错误,否则我们 通过最小二乘法求得的回归方程就缺乏理论依据。
这些检验主要探究的问题为:•残差是否为随机性、是否为正态性、是否不为异方差•高度相关的自变量是否引起了共线性•模型的函数形式是否错误或在模型中是否缺少重要的自变量•样本数据中是否存在异常值1.残差图分析所谓残差图就是以残差e =y -y为纵坐标,某一个合适的自变量为横坐标的散点图 t t t残差中包含了许多有关数据和模型的信息,它是研究回归诊断最基本及最重要的统计量残 差图分析的基本思想是,在回归模型的假设中,我们总是假定误差项是独立的正态分布随机 变量,且均值为零和方差相等为如果模型适合于观察到的数据,那么残差作为误差的无偏 估计,应基本反映误差的假设习性即残差图应该在零点附近对称地密布,越远离零点的地 方就疏散,则在形象上似有正态趋势,常认为模型与数据拟合得很好如图31. 3所示,是残 差的各种可能出现情况C.i LXx xXX x * X X x x XXx x Xx XX XX xXG)正常=正态分布(b)异常点’可疑偏离很大点(c)异方差’残差随区而増大Xx Xxxx xv X x X* X XXXXXX x yX X X(d)异方差’残差随X而増减@)非随机性二残差非线性趋势(f)非随机性’残差线性趋势图31.3残差的主要几种类型若残差图呈现如图31.3 (a)所示的形式,残差是随机的且不表示出一定的趋势与形式, 我们认为建立的回归模型应诊断为无甚大问题。
更进一步的诊断应该采用学生化残差鉴别是 否正态性一个简单的思想就是,如果模型假设正确的话,残差就应该是误差的良好估计, 那么残差全体构成的直方图应当与正态曲线很相似我们可以求出估计残差的方差Var® ),t且符合正态分布:(31.46)(31.47)Q ~ N (0,(1 - h )孚t t N — 2那么学生化残差:八 /\& y — yt 二 〜N(0,1)...Var(Q ) ...;(1 — h )ESS/(N — 2)t t则遵循标准正态分布在实际中,学生化残差常与配合作图,会有更好的直观判断效果 若残差图呈现如图31.3 (b)所示的形式,有一个对既定模型偏离很大的观察数据点, 称为异常点如果怀疑异常点是由于记录数据中发生的错误或者在测量过程中采用了拙劣的 技巧,我们理应从数据集中删除,重新回归模型但对异常点的处理须持谨慎态度,因为异 常点的出现可能代表了相当重要的某些数据,它恰好成为我们探究某些事先不清楚或许是更 为重要的因素的线索在SAS系统的reg回归过程中用来度量异常点影响大小的统计量是COOKD统计量,计算方法请参阅SAS/STAT软件使用手册若残差图呈现如图31.3 (c)所示的形式,残差随的增大而增大。
如图31.3 (d)所示的 形式,残差随的增大而先增后减,则蕴含着残差乃至误差对于不同的观察值具有不同的方差 变化,称为异方差在这种场合应该考虑在回归之前对数据或进行变换,实现方差稳定后再 拟合回归模型原则上,当误差方差变化不太快时取变换,当误差方差变化较快时取变 换log y或ln y,当误差方差变化很快时取变换1/y当然,还存在着不少其他变换,如著. y 九一1名的Box-Cox幕变换 一九若残差图呈现如图31.3 (e)所示的形式,显示了模型本身具有非线性趋势,或者提示人 们在模型中是否忽略了若干重要的变量如图31.3 (f)所示的形式,显示了模型本身具有线 性趋势同样表示了模型的错误选定2. 共线性回归研究中很容易发生模型中两个或两个以上的自变量高度相关,从而引起最小二乘估 计可能很不精确高度相关的自变量以及由它们所引起的估计问题合在一起称之为共线性(collinearity) 问题为什么共线性会引起参数估计可能很不精确呢?主要原因是最小二乘法所利用的数据信 息,如果存在共线性,就可能已经被其他的自变量说明了大部分,因此用剩余的少量数据估 计参数,将产生估计参数的方差很大,置信区间也会很大,假设检验也使人缺乏信任感。
在 实际中,最常见的问题是一些重要的自变量很可能由于在假设检验中t值不显著而被不恰当 地剔除了共线性诊断问题就是要找出哪些变量间存在共线性关系SAS系统的reg过程中 提供了特征值法、条件指数collin和方差膨胀因子vf 请参阅SAS/STAT软件使用手册3. 误差的独立性在回归诊断中,有一个非常重要的回归模型假设需要诊断和检验,那就是回归模型中的 误差项的独立性如果误差项不独立,那么我们对回归模型的许多处理,包括误差项估计、 假设检验等都将没有推导依据由于残差是误差的合理估计,因此检验统计量通常是建立在 残差的基础上检验误差独立性的最常用方法,是对残差的一阶自相关性进行Durbin-Watson 检验原假设H°:误差项是相互独立的,备选假设H1:误差项是相关的检验统计量为:DW =迓(e - e )2 /ESSt t—1t=2我们可以通过简单不等式证明:0 < X(e — e )2 < 2(迓e2 + 迓e2 ) < 4迓e2 = 4ESStt=2t —1 t t—1t=2 t=2tt=1因此,DW统计量应满足:0
在给定显著水平下,我们可以查Durbin-Watson表得到不能拒绝独立性原假设的区间DW
缺省 SLSTAY =0.1• SELECTION=STEPWISE SLENTRY =入选水平 SLSTAY=i 除水平逐步法(STEPWISE):按前进法进入变量,再对模型内所有变量检验,看是否有新因 变量引入而对模型的贡献变得不显著的变量,若有就剔除,若无则保留,直至方程内所有的 变量均显著,显然逐步法有两个水平,即选入水平和剔除水平,而且剔除水平应低于选入水 平缺省 SLENTRY =0.15 SLSTAY =0.1在上述三种方法的使用中,若要求打印出每一次选入或剔除变量进行模型拟合时的所有 统计量,可以加选DETAILS o(2)・NOINT——表示拟合无常数项(截距)的回归模型与屏幕输出有关的选项有:• CORRB――输出参数估计的相关阵• STB――输出标准化偏回归系数矩阵• P——输出个体观测值、预测值及残差若已选了 CLI、CLM、R,则无需该选项• R――输出每个个体观测值、残差及标准误差• CLM――输出每个观测值因变量期望值的95%的上、下限• CLI――输出每个个体观测值的95%的上、下限与残差分析有关的选项有:• VIF——输出变量间相关性的方差膨胀系数(Variance Inflation Factor) ,VIF越大, 说明由于共线性存在,使方差变大。
• COLLIN——输出条件数(Condition index),它表示最大的本征性与每个自变量 本征值之比的平方根一般情况下,条件数越大越可能存在共线性• TOL——表示共线性水平的容许值,TOL(Tolerance Value)越小说明其可用别 的自变量解释的部分多,自然可能与别的自变量存在共线性关系• DW 输出 Durbin-Watson 统计量3. 其他选择语句注意,这部分的语句可以在REG过程被被激活后,以交互式方式运行• OUTPUT语句——建立弘S的输出结果数据集语句格式为:OUTPUT OUT=SAS数据集名关键字名=输出数据集中的变量名其中关键字名为需要的统计量名,它们有P (预测值)、R (残差)、L95M (期望值的 95%的下限)、U95M (期望值的95%的上限)、L95 (个体预测值的95%的下限)、U95 (个 体预测值的95%的上限)、STDP (期望值的标准误差)、STDR (残差的标准误差)、STDI (预测值的标准误差)、STUDENT (学生化残差)、COOKD (COOK氏D值)• PLOT语句——绘制两变量的散点图语句格式为:PLOT X*Y / 选项• ADD变量名列表——向模型中增加变量• DELETE变量名列表——删除原拟合模型中的有关变量• REFIT——重新拟合模型• PRINT~ 出有关模型的相关信息七、应用举例例31・1广告花费X与销售额Y的回归模型。
大多数公司最终会询问关于花费在广告上 的费用对公司产品销售额的影响程度由于广告需要一定的时间才能达到它的效应,同时它 的效应也不是永久持续的,它的影响也许仅仅延续开头的一段时期假设公司相信销售额与 当月以及前两个月内所花的广告费有较密切的关系,即意味着:Yt与Xt,Xti, Xt2有密切的 t t t-1 t - 2关系假设它们之间存性关系,建立模型为:Yt=P o+P !Xt+P 2 Xt-i+P 3 X"t我们现在有某公司15个月内有关广告花费X与销售额Y的数据,如表31.1所示表31・1广告额与销售额月t月销售额Yt月广告花费Xt12945280242954003564545046995590583456506969575071104589081239510009137451050101509512001116445125012177951350131914514601420495150015218451650程序如下:LIBNAME STUDY "D:\SASDATA\MYDIR";DATA STUDY.AAA06;INPUT SALES ADV @@;ADVLAG1=LAG1(ADV);ADVLAG2=LAG2(ADV);OUTPUT;CARDS;2945 280 4295 400 5645 450 6995 590 8345 6509695 750 11045 890 12395 1000 13745 1050 15095 120016445 1250 17795 1350 19145 1460 20495 1500 21845 1650PROC REG DATA=STUDY.AAA06 OUTEST=STUDY.AAA06OUT;MODEL SALES=ADV ADVLAG1 ADVLAG2;RUN;输出的结果见表31.2。
表31.2 回归分析的第一次结果The SAS SystemModel: M0DEL1Dependent Variable: SALESAnalysis of VarianceSum ofMeanSourceDFSquaresSquareF ValueProb〉FModel3331387422.75110462474.253232.2360.0001Error9307577.2513734175.25015C Total12331695000Root MSE184.86549R-square0.9991Dep Mean13745.00000Adj R-sq0.9988C.V.1.34497Parameter EstimatesParameterStandardT for H0:VariableDFEstimateErrorParameter=0Prob > |T|INTERCEP1522.130666 :372.394170171.4020.1944ADV13.6814841.778890922.0700.0684ADVLAG114.9658061.466149043.3870.0080ADVLAG2 1 5.199508 1.90840150运行后,得到的最小二乘回归形式为:Y*=552・1+3・681X*+4・966X .+5.200X _t t t-1 t-2进一步统计分析,按显著性水平=0.05需要剔除ADV变量,并要求绘制残差图,要再 次提交下列程序: DELETE ADV;Pri nt ;PLOT R. * P. / SYMB0L='*'; RUN; 再一次的输出的结果见表31.3。
表31.3回归分析的第二次结果The SAS SystemModel: M0DEL1Dependent Variable: SALESAnalysis of VarianceSum ofMeanSourceDFSquaresSquareF ValueProb〉FModel2331241050.48165620525.243648.4350.0001Error10453949.5245545394.95245C Total12331695000Root MSE213.06091R-square0.9986Dep Mean13745.00000 ,Adj R-sq0.9984C.V.1.55010Parameter EstimatesParameterStandardT for H0:VariableDFEstimateErrorParameter=0Prob > |T|INTERCEP11161.641018239.517605774.8500.0007ADVLAG115.8732331.612430143.6420.0045ADVLAG217.9446401.581319785.0240.0005将Xt从模型中剔除而重新建立模型,得到的估计方程为:Y=1161・6+5・873X* +7.945X ,t t-1 t-2注意到新的估计方程式并不是从原来的方程式中单纯剔除Xt项而获得,新的估计方程式 也必须重新进行完全的F检验与部分的T检验,结果显示余下的变量应当保留在模型中。