2022年2022年广义线性模型 .pdf
1 广义线性模型广义线性模型 *(Nelder 和 Wedderburn ,1972)除了正态分布,也允许反应分布,以及模型结构中的一定程度的非线性。GLM 具有基本结构g( i)=Xi,其中 i E(Yi) ,g 是光滑单调 链接函数 ,Xi 是模型矩阵的第i 行,X 和 是未知参数的向量。此外,GLM 通常会做出 Yi 是独立的和 Yi 服从一些指数族分布的假设。指数族分布包括许多对实际建模有用的分布,如泊松分布,二项分布,伽马分布和正态分布。 GLM 的综合参考文献是McCullagh 和 Nelder (1989) , 而 Dobson(2001)提供了一个全面的介绍。因为广义线性模型是以 “ 线性预测器 ”X 的形式详细说明的,所以线性模型的许多一般想法和概念通过一些修改而继续存在到广义线性模型中。除了必须选择的链接函数和分布之外,基本模型公式与线性模型公式基本相同。当然,如果恒等函数被选择作为链接以及正态分布,那么普通线性模型将作为特例被恢复。然而,泛化是以某种成本为代价的:现在的模型拟合必须要迭代完成,而且用于推理的分布结果是近似的,并且由大样本限制结果证明是正确的而不是精确的。但在深入探讨这些问题之前,请考虑几个简单的例子。i=cexp(bti), 例 1:在疾病流行的早期阶段, 新病例的发生率通常会随着时间以指数方式增加。因此,如果 i是第 ti 天的新病例的预期数量,则该形式的模型为请注意, “ 广义” 和“ 一般” 线性模型之间存在区别-后一个术语有时用于指除简单直线以外的所有线性模型。可能是合适的,其中c 和 b 是未知参数。通过使用对数链路,这样的模型可以变成 GLM 形式log( i)=log(c)+bti=0+ti1(根据 0=logc 和 1=b 的定义) 。 请注意, 模型的右侧现在在参数中是线性的。反应变量是每天新病例的数量,因为这是一个计数,所以泊松分布可能是一个合理的可以尝试的分布。因此,针对这种情况的GLM 使用泊松反应分布,对数链路和线性预测器 0+ti 1。, 例 2:狩猎动物捕获猎物的速度yi 往往随着猎物密度xi 的增加而增加,但最终会趋于平衡,当捕食者捕获尽可能多的猎物时。对于这种情况一个合适的模型可能是名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 1 页,共 35 页 - - - - - - - - - 2 其中 a是未知参数,表示最大捕获率,h 是未知参数,表示捕获率为最大速率一半时的猎物密度。很显然,这个模型在其参数中是非线性的,但是通过使用倒数链路,右边的参数可以是线性的:(这里 01/a 和 1h/a) 。在这种情况下, 猎物捕获率的标准差可能与平均速率大致成比例,建议使用Gamma分布作为反应,并完成模型设定。们不限于示例的简单直线形式,但可以有对于线性模型可能具有的线性预测器的任何结构。2.1GLMs 的理论GLM 的估计和推理基于最大似然估计理论,尽管可能性的最大化需要迭代最小二乘法,与 1.8.6 节的方法相关。本部分首先介绍指数分布族,它允许发展一种通用方法来最大化GLM 的可能性。然后基于似然理论的一般结果 (在本章末尾推导)来讨论 GLM 的推论。在本节中,区分反应数据y 和 Y 的观察值的随机变量有时很有用,所以它们在符号上具有区别性:对于估计和估计值,它尚未完成。名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 2 页,共 35 页 - - - - - - - - - 3 2.1.1 指数族分布f(y)=expy-b( )/a()+c(y,),GLM 中的反应变量可以来自指数族的任何分布。如果一个分布的概率密度函数或概率质量函数可以写成,则该分布属于指数族分布名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 3 页,共 35 页 - - - - - - - - - 4 其中 b,a和 c 是任意函数, 是任意的 “ 尺度” 参数,并且 被称为分布的 “ 典范参数 ” (在 GLM 上下文中, 将完全依赖于模型参数 ,但是目前没有必要做这个明确) 。例如,很容易看出,正态分布是指数族的一员,因为, . 表 2.1 给出了在 R 中为 GLM 实施的指数族成员的类似分解。用 a,b 和 可以得到指数族分布的均值和方差的一般表达式。给定一个特定的 y,的对数似然性仅仅只是将logf (y)视为 的一个函数。那是. 将 l 作为一个随机变量来处理, 通过用随机变量Y 替换特定的观测值y,可以评估?l/? 的期望值:使用 E(?l/? )=0 这个一般结果,(在 取真值时,参见2.4 节中的( 2.14) )和重新排列意味着E(Y)=b0( ) 。 (2.1)即任何指数族随机变量的均值由bw.r.t.的一阶导数给出。 ,其中 b 的形式取决于特定的分布。该等式是将GLM 的模型参数 与指数族的典范参数联系起来的关键。在 GLM 中,参数 决定了反应变量的均值,并且通过( 2.1) ,它们决定了每个反应观测值的典范参数。, 再次对似然性微分处理, 并将其插入到一般结果中,E(?2l/? 2)=-E(?l/? )2(衍生物在真实值下计算,参见结果( 2.16) ,第 2.4 节) ,重新安排第二个有用的一般结果:var(Y)=b00( )a( ).a原则上可以是 的任何函数,并且当与 GLM 一起工作时,如果 是已知的,处理任何形式的a 都是没有困难的。然而,当未知,事情就会变得很尴尬,除名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 4 页,共 35 页 - - - - - - - - - 5 非我们可以写出( )=/ ,其中 是一个已知常数。事实上,这种限制形式涵盖了所有有实际意义的案例(见表2.1) 。a( )=/ 允许基于正态分布的模型中的不等方差的可能性,但是在大多数情况下,仅为 1.因此,我们现在有var(Y)=b00( ) / .(2.2) 在随后的章节中,我们通常会将var(Y)视为 E (Y)的函数,这会很方便,并且由于 和 通过(2.1)链接,我们总能定义一个函数V( ) =b00 ( )/ ,使得 var(Y)=V( ) 。表 2.1 列出了几个这样的功能。2.1.2 拟合广义线性模型回想一下, GLM 模拟独立反应变量的n 阶向量 Y,其中 E (Y) ,通过g( i)=Xi和Yi f i(yi),其中 f i(yi)表示指数族分布,典范参数i由 i(通过方程 2.1)决定,因此最终由 决定。给定 Y 的一个观测值向量y,的最大似然估计是可能的。由于Yi 是相互独立的, 的似然函数是n L( )=Yfi(yi), i=1 因此 的对数似然函数是, 右边对 的依赖是通过 i对 的依赖。请注意,函数a,b 和 c 可能随着 i 而变化-例如,允许不同的二项分母ni,对于每个二项式反应的观测值,或对于正常反应的不同方差 (但在常数内已知)。另一方面,对于所有的 i,假设 是相同的。正如前一节所讨论的那样,对于实际工作来说,只考虑可以写出ai( )=/ i的情况就足够了,其中i 是已知常数(通常为1) ,在这种情况下. 通过偏导 lw.r.t 来最大化收益。 的每个元素,令结果表达式为零并求出 。, 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 5 页,共 35 页 - - - - - - - - - 6 通过这个链式法则, 所以微分( 2.1) ,我们可以得到, 然后推出. 将(2.1)和( 2.2)代入最后一个方程,意味着求解 的方程是(2.3) , (2.4) 然而,如果权重 V(i)事先已知且与 独立,那么这些方程就正好是为了通过非线性加权最小二乘寻找而必须解出的方程。在这种情况下,最小二乘的目标是其中 i非线性地依赖于 ,但权重 V(i)被视为固定的。要找到最小平方估计值,包括求解 ?S/?j=0 ?j,但当 V(i)项被视为固定时,这个方程组很容易被看作是( 2.3) 。这种对应立即提出了一种求解(2.3)的迭代方法。令k 表示在第 k 次迭代处的估计参数向量,并且令 k 和 k 为具有元素 ik=Xik和 ik=g -1 (ik)的向量,其中g-1( )是链路的反函数。从参数估计开始,0 ,迭代以下步骤,直到 k 的序列收敛:1.计算当前 k 隐含的 V(ik )项。2.为了得到 k+1 (V(ik )被视为固定的而不是作为的函数),给出这些估计,使用第 1.8.6节的方法来最小化( 2.4) ) 。3.设置 k 到 k+1 实际上,这种方法比需要的要慢。第2步本身涉及迭代,但在已经收敛之前实际上将非线性最小二乘方法迭代到收敛没有多大意义。因此,第2 步通常被替换为:2.使用 k 作为初始值,仅执行第1.8.6节中给出的求解( 2.4)的迭代方法的一次迭代来获得 k+1 。名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 6 页,共 35 页 - - - - - - - - - 7 应用这种方法会产生一个相当紧凑和简洁的方案。为了看到这个,让我们以矩阵形式写出非线性最小二乘问题。定义对角矩阵Vk ,其中 Vkii=V (ik)(2.4) 按照第 1.8.6节的方法, 被它在 k 附近的一阶泰勒展开替代,所以因此,没有进一步的近似Jij=?i/? j| ?k.Now . 由“伪数据”的定义对角线权重矩阵. 因此,下面的步骤可迭代到收敛1.使用当前的 k 和 k 计算伪数据 zk和迭代权重 Wk 。2.将对于 的平方和最小化,以便获得 k+1 , 因此 k+1=X k+1和 k+1 。将 k 增加 1。收敛的 解决了( 2.3) ,因此是 ? 的最大似然估计。该算法在大多数实际情况下趋于一致,但也有例外(例如,二项数据的不良或过于灵活的模型)。请注意,要开始迭代,我们只需要0 和 0 的值,不需要 0 。因此,迭代通常通过设置 0i=yi和 i0=g (0i )开始,并根据需要轻微调整0i ,以避免 无 限 的i0( 例 如 , 如 果yi=0且 有 对 数 链 路 )。 该 方 法 被 称 为IterativelyRe-weightedLeastSquares(IRLS) , 原因很明显,在此背景下,归因于 Nelder和 Wedderburn (1972) 。2.1.3IRLS 的目标是对数似然的二次近似IRLS 迭代中的工作线性模型不仅仅是寻找参数的最大似然估计的手段。在一个加性常数内名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 7 页,共 35 页 - - - - - - - - - 8 (在收敛时)也是模型在附近的对数似然的二次近似。显然,第一个衍生工具是关于对数似然和S 之间的 j匹配:实际上它们都是零。S 的二阶导数矩阵为-XWX/ ,并且这被证明与对数似然的期望的二阶导数矩阵匹配,并且因此在大样本极限中由大数定律来匹配二阶导数矩阵本身。为了证明这一点,首先将u 定义为关于模型参数的对数似然的导数向量,因此 ui=?l/? i,然后将( 2.3)中的导数以矩阵向量形式重写为u=XTG- 1V- 1(y- )/ .?注意,如果 V(i)被视为 的函数,算法不会最小化(2.4) ,因为在这种情况下,令导数为零不会产生(2.3) 。换句话说,最大似然与具有平均方差关系的最小二乘基本上不同。E(uuT) =XTG- 1V- 1E(Y- )(Y- )TV-1G- 1X/ 2=XTG-1V-1VV-1G-1X/ =XTWX / 因为 E(Y-T) (Y- )T=V。根据一般似然结果( 2.19) ,在 2.4.2 节中,-E(uu)也是对数似然的期望二阶导数矩阵。导数的这种对应关系足以证明S 是 附近的对数似然的二次近似,并且由于MLE 的一致性,它们在真实参数值附近。2.1.4AICforGLMs 通过对可能性进行直接比较的模型选择存在如下问题:如果将冗余的参数添加到正确的模型中,可能性几乎总是增加(并且从不减小),因为额外的参数让模型更接近数据,即使这只意味着对数据的组成部分进行“噪声建模”。正如在线性模型的情况下,如果我们能够根据它们拟合数据均值而不是数据y 的能力来选择模型,这个问题会得到缓解。在GLM 背景下,合理的方法是根据模型最大化l(; )而不是l( ;y)的能力来选择模型,但是为了做到这一点,我们必须能够估计出 l(; ) 。实际上这个估计很明确。从2.1.3节我们已经有, 并且因为 y=时,这也必须成立. 然后,这个论据导致在第1.8.5节的( 1.15) (只修改权重)产生了估计量tr(A) 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 8 页,共 35 页 - - - - - - - - - 9 l( ?;y)-tr( A)+n/2 A=X(XTWX )-1XTWandhencetr( A)=p,(可识别的)模型参数的数量。因此,在模型之间进行选择时, 我们会选择哪个模型具有最高的l( )-p 值,这相当于选择了Akaike 信息准则的最低值( Akaike,1973) ,AIC=2-l ( )+p。前面的论点假设已知。如果不是,那么就需要一个的估计来计算AIC,结果 AIC 中的惩罚项 p 将变为 p+1。这种归纳在 2.4.7 节中是合理的。2.1.5 的大样本分布GLM 的分布结果并不精确,反而是基于大样本近似,利用包括一致性在内的最大似然估计的一般性质(见2.4 节) 。根据最大似然估计量的一般性质,我们认为,在大样本限制下,? N(,I- 1), 其中 I=E(uuT)是模型参数的信息矩阵,u 是关于模型参数的对数似然函数的导数向量,因此ui=?l/? i(见 2.4 节(2.20)和(2.19) ) 。在第 2.1.3节中显示 E(uT)=XTWX/ ,因此在大样本限制中? N(,( XTWX )- 1).对于具有已知尺度参数的分布,可以直接使用此结果来查找参数的置信区间,但是如果尺度参数未知(例如对于正态分布),则必须估计它,并且区间必须基于一个合适的 t 分布。尺度参数估计在第2.1.7节中介绍。2.1.6 通过假设检验比较模型考虑检验H0:g( )=X00相反H1:g( )=X11, 其中 是反应向量 Y 的期望, Y 的元素是来自指数族分布的相同成员的独立随机变量,其中X0? X1。如果我们有反应向量的观测值y,则可以执行广义似然比检验。令 l(0 )和 l(1 )为两个模型的最大似然率。如果H0 为真,则在大样本限制中,2l(?1)-l(?0) 2p1- p0, (2.5) sothattwicethe?whichshouldstrictlybeamaximumlikelihoodestimate ,其中 pi 是模型 i 中(可识别的)参数( i)的数量(参见 2.4.5 节和 2.4.6 节推导该结果)。如果原假设是假的,则模型1 倾向于具有比模型0 高得多的可能性,从而两倍于严格应该是最大似然估计的 ?,或者在大样本限制下倾向于MLE 的估计值。名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 9 页,共 35 页 - - - - - - - - - 10 对数似然对于相关的2 分布的一致性差异太大。如果可以计算相关模型的对数似然性,那么近似结果(2.5)仅仅是有用的。在由 IRLS 估算的 GLM 情况下,只有当尺度参数已知时才是这种情况。因此,结果可以与泊松和二项模型一起直接使用,但不能与正态 , 伽马或反高斯分布 (其中尺度参数未知)一起使用。稍后将简短地讨论在后面这些情况下做什么。偏差在实践中使用 GLM 时,在普通的线性建模中, 用与残差平方和类似的方式来解释数量是有用的。Xi=1h i 这 个 数 量 是模型的偏差,被定义为thetermsinsidethesummation其中 l(max )表示饱和模型的最大可能性:每个数据点具有一个参数的模型。在给定数据的情况下,l(max )是可能性具有的最大值,并且可以通过简单地设置=y并计算出可能性来计算。 ? 和 分别表示关于饱和模型和兴趣模型的典范参数的最大似然估计。注意如何将偏差定义为与独立。表 2.1 列出了单个数据对偏差的贡献,对于若干分布-这些是偏差定义中求和中的项(术语)。scaleddeviance, 与偏差相关的是调整偏差,D?=D/,这取决于尺度参数。 对于二项分布和泊松分布, 其中 =1 ,偏差和比例偏差是相同的,但通常情况并非如此。通过广义似然比检验结果(2.5) ,我们可以预计,如果模型是正确的,那么近似D? 2n- p, (2.8) tojustify(2.8)asalargesampleapproximationundermanycircumstances 在大样本限制中。实际上,这样的论点是捏造的,因为证明(2.5)的限制论证依赖于模型中参数的数量保持不变,而样本容量趋于无穷大,但饱和模型具有D = = 2l(?max)-l(?) n 2iyi( ?i- ?i)-b( ?i)+b( ?i), (2.6) (2.7) 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 10 页,共 35 页 - - - - - - - - - 11 与数据一样多的参数。渐近结果可用于表2.1中的一些分布,以在许多情况下将(2.8)证明为大样本近似值(参见当然,对于正态分布和一致性链接,我们使用第1 章的结果。这对于正态分布情形下是精确的。但是,请注意,它完全打破了二进制数据的二项式。鉴于偏差的定义,很容易看出,本节开始的似然比检验可以通过重新表达两次对数似然比统计量作为。然后在H0 下(2.9) (在大样本极限内), 其中 Di* 是具有 pi 个可识别参数的模型i 的偏差。 但是,这只有在尺度参数已知时才有用,以便可以计算D*。与未知 的模型比较在 H0 下我们有近似的结果and, 并且, f如果被看作是渐近独立的,这就意味着, 在大样本限制下(当然,在普通线性模型特例中结果确实如此)。F 的有用性质是它可以在不知道的情况下进行计算,这可以从比率收益的顶部和底部取消,在 H0 下,近似结果. (2.10) 这个结果的优点是,当未知时,它可用于基于模型比较的假设检验。缺点是对于 Di 的可疑的分布假设以及它所基于的独立性近似。当然, 一个明显的替代方法是使用估计值来获得每个模型的估计值Di*=Di ,然后使用( 2.9)进行假设检验。然而, 如果我们使用估计(2.11) 来达到这个目的,则很容易看出它只是(n-p1) F,所以我们的测试将完全等同于使用F 比率结果( 2.10) ,但是使用 Fp1-p0,作为参考分布。显然直接使用(2.10)是一种更保守的方法,因此通常是首选:它至少在估计尺度参数时考虑到了不确定性。2.1.7 ?andPearson sstatistic正如我们已经看到的那样, 参数 的 MLE 可以在不知道尺度参数的情况下获得,但是在这些参数未知的情况下,通常必须进行估计。近似结果(2.8)提供了一个明显的估计量。一个 2n -p 随机变量的期望是n-p,所以等同于观测值名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 11 页,共 35 页 - - - - - - - - - 12 达到我们的近似预期值?D=D/?(n-p).(2.11) 第二个估计量基于Pearson统计量,定义为显然, X2/是一组零均值,单位方差,随机变量的平方和,具有n-p 个自由度,这表明如果模型适当,则近似X2/ ? x2n-p:这个近似值是有根据的。将观测到的 Pearson统计量设置为我们得到的预期值?=X?2/(n-p).请注意,它很明确的表明这一点, 其中 W 和 z是 IRLS 权重和伪数据,在收敛时进行评估。2.1.8 典范链接函数典范链接 gc 用于分布,是链接函数,使得gc(i)=i,其中 i 是分布的典范参数。例如,对于泊松分布,典范链接是对数函数(其他示例见表2.1) 。典范链接的使用意味着i=Xi (其中 Xi 是 X 的第 i 行) 。典范链接倾向于具有一些很好的性质, 例如确保 保持在反应变量的范围内,但它们也具有更多微妙精细的优点,其中之一在此处被推导出。回想一下,似然最大化涉及对每个j的对数似然性微分,并令结果为零,以获得方程组但是,如果正在使用典范链接,则?i/? j=Xij,如果像通常情况, wi=1?i,这个方程组可简化为XTy- XT?= 0, i.e.toXTy=XT?. 即,简化为现在考虑 X 包含 1 列的情况:任意其他加权求和,其中权重由模型矩阵P 给出,这意味着该系统中的一个方程简单地为iyi=ii。?RecallthatifZi:i=1.nareasetofi.i.d.N(0,1)r.v.sthenPZi2 2n. 回想一下,如果Zi :i=1.n 是一组独立同分布, N (0,1) r.v.然后是 PZi2-2n 。列(或这些的线性组合)在原始数据和拟合值之间保存。其中一个实际结果是,对于任何具有截距项和典型联系的GLM ,残差总和为零:这种“ 观察到的无偏性 ” 是一个令人放心的特性。 结果的另一个实际用途是在使名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 12 页,共 35 页 - - - - - - - - - 13 用对数线性模型的分类数据分析中,其中它通过规定模型提供了一种方法,可以保存任何模型中保存的内置于研究设计中的总数。2.1.9 残差模型检查可能是应用统计建模中最重要的部分。在普通线性模型情况下,这是基于对模型残差的检验,模型残差包含数据中的所有信息,而不是由模型的系统部分来解释。残差检验也是GLM 案例模型检验的主要手段,但在这种情况下,残差标准化是必要的,并且更加困难。对于 GLM 而言,不仅仅只检查原始残差的主要原因是难以检查原始残差检验假设的均值方差关系的有效性。例如,如果采用泊松模型,则残差的方差应与拟合值(i)的大小成正比增加。然而,如果将原始残差与拟合值作图,则需要一种杰出的能力来判断残差的可变性是否与均值成比例增加,而不是均值的平方根或平方。出于这个原因,通常将 GLM 残差标准化,以这种方式,如果模型假设正确,则标准化残差应该具有大致相等的方差,并且尽可能与来自普通线性模型的残差相似(尽管见第 6.5 节中图 6.9 的替代绘图方法)。皮尔森残差根据拟合模型,将残差标准化的最明显的方法是将它们除以与其标准差成比例的量。这导致了皮尔森残差, 如果模型正确, 它应该具有近似为零的均值和方差 。当与拟合值或任何协变量(无论是否包含在模型中)进行绘图时,这些残差不应显示均值或方差的任何趋势。 “ 皮尔森残差 ” 这个名字与相关的皮尔森残差的平方和给出了2.1.7 节中讨论的皮尔森统计量。请注意,皮尔森残差是来自收敛IRLS 方法的工作线性模型的残差除以收敛IRLS 权重的平方根。异常残差在实践中,皮尔森残差的分布在零附近可以是非常不对称的,因此它们的行为不像可能希望的那样接近普通线性模型残差。异常残差在这方面通常是优选的。异常残差是通过注意到偏差对GLM 起到与残差平方和对普通线性模型起相同作用来得到:确实,对于普通线性模型来说,偏差是残差平方和。在普通线性模型情况下,偏差由平方残差和组成。这就是具有适当标示的偏差组成部分的平方根的残差。因此,把 di 写成第 i 个数据所贡献的偏差的一部分(即(2.7)中求和的第i项) ,我们有名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 13 页,共 35 页 - - - - - - - - - 14 并且通过与普通线性模型类比,我们可以定义. 根据需要,这些“异常残差” 的平方和给出了偏差本身。现在如果计算了所有参数已知的模型的偏差,则(2.8)将变成 D*? x2n,这可能表明对于单个数据di? x12,可推出。当然, (2.8)不能合理地应用于单个数据,但对于良好的拟合模型来说,这表明我们可能预期异常残差具有类似于 N(0,1)随机变量的行为,特别是在对(2.8)预计是一个合理的近似值的情况下。2.1.10拟极大似然迄今为止, GLM 的处理假定反应变量的分布是指数族的已知成员。如果有充分的理由假设该反应变量遵循特定的分布,则它对该分布上的基本模型很有吸引力,但是在许多情况下,反应分布的性质并不十分清楚,并且可能它只能详细说明反应变量的方差与其均值之间的关系。也就是说,函数V( )可以被指定,但是好不了多少。然后出现的问题是, 是否有可能开发用于拟合和推理GLM 的理论,从仅指定平均方差关系的位置开始。事实证明,基于拟极大似然的概念,开发出令人满意的方法是可能的。考虑一个具有均值 i 和已知方差函数V(i)的随机变量的观测值yi。然后对于给定yi 的 i的对数拟极大似然定义为(2.12) 正如我们将会看到的,这个函数的关键特征是它具有li 的许多有用特性,即对应于单个观测的对数似然,但只需要 V 的知识而不是 Yi 的全部分布。 假设数据是独立随机变量的观测值,我们可以定义所有反应数据的平均向量的对数拟极大似然或定义 的任意参数向量为. q 的关键特征是,为了推导GLM ,它表现得与对数似然函数非常相似,但只需要方差函数就可定义它。例如,考虑获得GLM 参数 的最大拟似然参数估计。对q 关于 j产量微分, 所以参数估计是方程的解名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 14 页,共 35 页 - - - - - - - - - 15 但这正是方程组( 2.3) ,必须解决这个问题才能找到GLM 的 m.l.e.s。因此,最大拟似然参数估计可以通过通常的GLMIRLS 方法找到,在任何情况下只需要知道 V( ) 。此外,对数拟似然性与参数估计值的大样本分布结果(参见2.1.5节)的对数似然性具有足够的性质,并且对于参数的最大拟似然估计量值也适用。类似地,当对数似然度 l 被对数拟似然函数 q代替时, 2.1.6节的大样本分布结果支持用GLM进行假设检验。第2.4.8节提供了这些断言的理论基础。请注意,饱和模型的对数拟似然始终为零,因此GLM 的拟偏差很简单Dq=-2q( ?) .很显然,关于残差和尺度参数估计的讨论也从似然性中延续到拟似然的情况,但不超过用 q 代替 l。名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 15 页,共 35 页 - - - - - - - - - 16 拟似然性方法的实际应用要求对(2.12)中的积分进行评估,但对于大多数实际有用的 V 的形式,这是可能的: McCullagh 和 Nelder(1989)给出了例子,或者在R 中你可以输入例如以获得所执行的任何特定均值方差关系的qi 形式。对于表 2.1 中对应于指数族分布的均值方差关系,拟偏差的形式精确对应于该族的偏差形式。拟似然性的一个主要实际用途是提供一种来模拟比泊松或二项分布(具有其固定尺度参数)更可变的计数数据的方法:拟似然性方法假定未知。这种 “ 过度分散” 的数据在实践中很常见。另一个实际用途是提供一种对具有平均方差关系的数据建模的方法,其中没有明显的指数族分布:例如,预期方差与均值成比例的连续数据。2.2GLM 的几何图形GLM 和 GLM 拟合的几何图形比普通线性模型的几何图形更难,因为用于判断模型拟合的可能性通常不意味着拟合可以通过模型和数据之间的欧几里得距离来判断。图 2.1 用一个含有 2 个参数的 GLM 的 3 个数据拟合 Gamma分布和对数链接的例子说明了GLMs 的几何情况。 1.4 节的平坦模板子空间现在被一个弯曲的“ 模型流形 ” 所取代,该模型流形由模型可预测的所有可能的拟合值向量组成。由于模型流形与数据之间的欧几里德距离不再被用来进行拟合度的测量,因此必须采用不同的方法来说明估计的几何结构。图 2.1 右端板的黑线显示了所有反应变量的组合,这些变量产生了相同的估计模型。注意这些线条通常并不平行,并且通常不与模型流形正交。要充分理解图2.1,它可能有助于考虑一些不同的含有2 个参数模型的图形。1.对于一个普通的线性模型, 模型流形将是一个平面, 所有的相等线拟合将是正交的(因此相互平行) 。2.对于假定正态分布(但是非一致性链接)的GLM,相等线拟合将与它们遇到的模型流形的(切线空间)正交。3.对于拟合 4个数据的含有 2个参数的模型,相等线拟合将变成相等平面拟合。一般来说,图2.1 所示的几何图形适用于任何GLM 。有了更多的数据,相等线拟合就变成了 n-p 维相等面的拟合, 其中 n 和 p 分别是数据和参数的数量: 对于任何固定的 ,方程( 2.3)给出了定义这样一个平面的y 的限制。请注意这些平面名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 16 页,共 35 页 - - - - - - - - - 17 x 图 2.1GLM 的几何图形。左端板展示了所示的三个x,y 数据的广义线性模型E(y) =exp(0+1x)的最佳拟合,假设每个yi 是由模型给出的带有均值的Gamma 分布随机变量的一个观测值。右端板展示了使用此模型作为示例的GLM拟合的几何图形。所示的单位立方体表示向量(y1,y2,y3)T 定义单点空间。灰色表面根据模型显示所有可能的预测值(在单位立方体内),即它表示所有可能的(1 ,2 ,3 )T 值。由于参数0 和 1允许在所有可能的值上变化,这就是相应模型 “ 拟合值 ” 追踪的表面: “ 模型流形 ” 。连续线从立方体的一个面开始,在另一个面上离开,是等效拟合线:位于这样一条线上的反应数据(y1,y2,y3)T 的值各自导致相同的0 ,1的最大似然估计,即相同的(1 ,2 ,3 )T。注意等效拟合线既不相互平行,也不与模型流形正交。可以相交 -稍后将返回的点。对于离散反应数据,图片没有什么不同,尽管只有在可能性连续概括下(通常可以通过在概率函数中用合适的伽玛函数代替阶乘)才能严格得到相等线拟合。只有正态分布才能与模型流形有正交的相等线/平面拟合。对于其他分布,相等线/平面拟合有时可能彼此平行,但决不会全部与模型流形正交。2.2.1IRLE 的几何图形通过考虑一个参数模型对2 个反应数据的拟合, IRLS 估计算法的几何图形是最容易理解的。图 2.2 说明了这种模型的几何结构:在这种情况下,一个具有对数链接和 Gamma错误的 GLM ,1.20.00.20.40.60.81.00.20.30.40.50.60.70.8y名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 17 页,共 35 页 - - - - - - - - - 18 图 2.2GLM 的几何图形。 E (yi)i=20exp(-xi )其中 yi? Gamma和 i=1,2。左端板展示了拟合2 个 x,y 数据的模型(连续线)的最大似然估计值,显示为?。右端板显示了拟合几何图形。15 15 的正方形是空间 2 中的一部分,其中( y1,y2)定义为单点。 粗体曲线是 “ 模型流形 ” :它根据模型包含所有可能的点 (1 ,2 )(即,当 变化时, (1 ,2 )可以追踪该曲线)。细线是相等线拟合的例子。位于这些线之一上的所有点(y1,y2)共享相同的的 MLE ,即( 1 ,2 ) :这个MLE 是相等拟合线切割模型流形的地方。对于 =.1 ,.2,.3,.4,.5,.6,.7,.8,.9,1,1.2,1.5,2,3,4 绘制相等线拟合。(=.1 ,.7 和 2 用实线表示,当=2 时候的线接近图形的底部。当=.1 时候的线在该图的绘图区外,但出现在随后的图中。)但是可以通过链接和分布假设的任意组合来为GLM 构建类似的图片。现在拟合 GLM 的关键问题是模型流形不平坦, 并且相等线拟合与它们遇到的模型流形不正交。 IRLS 方法对拟合问题进行线性转换和重新调整,因此在当前的的估计值时, 模型流形和相等线拟合是正交的,并且在重新调整的空间中, 当前的估计值的位置由X 乘以当前的 估计值给出。这种重新调整会产生一个拟合问题,可以将其视为局部线性,从而可以用最小二乘法更新 估计值。图 2.3 说明了 IRLS 步骤如何涉及形成伪数据并对其进行加权,从而将拟合问题有效地转化为可以通过线性最小二乘法近似解决的问题。该图说明了 IRLS 步骤中涉及的转换,这些转换重复进行,直到IRLS 方法被迭代到收敛。(a) (b) 2.01.50.00.51.08101214161820 xy15051005101512名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 18 页,共 35 页 - - - - - - - - - 19 1 1 (c) (d) 1 1 图 2.3 基于图 2.2 所示的 GLM 的 IRLS 估计的几何图形。(a)显示拟合问题的几何图形 -模型流形是粗黑色曲线,等拟合线是细线(如图2.2) ,数据是?和拟合值的当前估计值 k 。 (b)问题以当前拟合值为新中心(yi 被 yi-i 代替) 。 (c)问题被重新线性调整,使得X 的列现在跨越模型流形的切线空间。切线空间用灰线表示(这个步骤用g0(ik) (yi-ik) )代替 yi-ik。(d)该问题是线性转换的,因此现在由Xk给出位置。对于大多数GLM,现在必须通过将相对于每个轴的分量乘以Wi来再次调整问题, 其中 Wi 是迭代权重:这将确保通过的相等估计线与切线空间正交。在当前的例子中,这些权重都是 1,因此所需的正交性已经成立。现在对于转换后的问题,在模型流形的附近,可以用切线空间来近似,其中相等线拟合近似正交:因此,可以通过找到转换数据的最小二乘投影来获得和 的更新估计,到切线空间(灰线) 。1.50.00.51.02.02.53.00.00.51.01.52.02.53.023.00.50.01.01.52.02.50.00.