纵向数据的稳健估计及案例分析_胡倩.pdf

第 43 卷 第 2 期 高 师 理 科 学 刊 Vol.43 No.2 2023 年 2 月 Journal of Science of Teachers College and University Feb.2023 文章编号：1007-9831（2023）02-0030-05 纵向数据的稳健估计及案例分析 胡倩1，罗文塔2（1.贵阳信息科技学院 信息工程系，贵州 贵阳 550025；2.贵州宏信创达工程检测咨询有限公司，贵州 贵阳 550016）摘要：纵向数据是对一组观测个体按时间或空间顺序重复跟踪监测而得，即对每一个体在不同时间或不同实验条件下进行多次测量，而多次重复观测之间一般具有相关性，会影响对应估计的准确性.研究一般线性模型下纵向数据的理论推导和实例分析，发现模型的最小二乘估计、极大似然估计以及约束极大似然估计的参数估计结果差别不大甚至基本相等，验证了纵向数据估计的稳健性.关键词：纵向数据；稳健估计；相关性；参数估计 中图分类号：O212 文献标识码：A doi：10.3969/j.issn.1007-9831.2023.02.007 Robust estimation of longitudinal data and case analysis HU Qian1，LUO Wenta2 （1.Department of Information Engineering，Guiyang Institute of Information Science and Technology，Guiyang 550025，China；2.Guizhou Hongxin Chuangda Engineering Detection&Consultation Co.Ltd，Guiyang 550016，China）AbstractAbstract：Longitudinal data is obtained by repeated tracking and monitoring of a group of observation individuals in time or space sequence，that is，multiple measurements are made for each individual at different times or under different experimental conditions However，there is a correlation between multiple repeated observations，which will affect the accuracy of corresponding estimationThe theoretical derivation and case analysis of longitudinal data under the general linear model is studied，and it is found that the parameter estimation results of least squares estimation，maximum likelihood estimation and constrained maximum likelihood estimation of the model have little difference or even almost equal，which verifies the robustness of longitudinal data estimation Key wordsKey words：longitudinal data；robust estimation；relevance；parameter estimation 纵向数据是指对一组个体按时间顺序或空间顺序追踪重复测得的数据，对每一个体在不同时间或不同实验条件下多次测量，所得的数据兼有时间序列和截面数据的特点纵向数据在医学、生物学、社会学、经济学以及心理学等领域极为常见 纵向数据的参数回归分析方法是早期研究的主要方法，随着研究的不断发展，已有许多学者将纵向数据研究应用于不同的模型中1-4一般线性模型往往假定误差项为多元正态分布，均值为零向量，协方差阵为分块对角阵，进一步按协方差阵可细分为均匀相关、指数相关和一步相关等，模型可用极大似然法或加权最小二乘法以及广义估计方程的方法进行估计或统计推断；广义线性模型520可将连续型反应变量的研究推广至离散型，如 Logistic 边缘模型和泊松回归模型等，对于反应变量是分类数据的情形同样适用.由于在纵向数据中对同一个体的多次重复观察之间往往具有相关性，如何处理这种个体内的相关性便成为纵向数据分析中不可回避的问题.长期以来，实际工作者和统计学家都在关注纵向数据的统计推断问 收稿日期：2022-07-01 作者简介：胡倩（1995-），女，贵州独山人，助教，硕士，从事概率论与数理统计研究E-mail： 第 2 期 胡倩，等：纵向数据的稳健估计及案例分析 31 题，Larid6等在一族广泛的参数模型下研究了纵向数据，提出了经验贝叶斯方法和极大似然方法估计模型参数Liang7等于 1986 年提出了经典的广义估计方程（Generalized estimating equations，GEE）方法来处理纵向数据中广义线性模型的参数统计推断，还进一步将协方差阵建模，假定各观测对象的组内相关性不变，也就是说相关矩阵对所有观测对象一致，根据实际情况将相关矩阵建模为讨厌参数的函数，并证明了即使对相关矩阵建模错误 GEE 方法仍能得到相合估计但是，这种建模错误会大大影响参数估计的效率，针对这个问题，文献8对 GEE 方法做了扩展，基于广义矩估计（Generalized method of moments，GMM）的思想，提出了二次推断函数（Quadratic inference function，QIF）方法近些年纵向数据研究问题越发成为学者们的研究热点关晓妮9等在广义线性模型下，研究了纵向数据下的稳健二次推断函数估计，通过改进得到了模型参数有效且稳健的二次推断函数估计（Effective and robust quadratic inferential function，ERQIF），并验证了其相合性芦飞10基于指数平方损失函数，针对纵向数据提出了均值和广义自回归参数的模型参数的稳健估计方法，并且其估计量的渐近正态性是可验证的.综上，可以看出学者们对于纵向数据的研究一直持续深入，且取得了较好的研究成果，但是对于纵向数据中参数估计的稳健性研究还有进一步的探究空间本文结合已有研究以及现有理论基础，研究了一般线性模型下的纵向数据，分别使用最小二乘估计、极大似然估计以及约束极大似然估计来计算模型的参数估计，对结果加以分析整合，验证纵向数据的参数估计稳健性.1 模型结构 1.1 纵向数据的线性模型 假定每个个体分别有n次观测，即假定每个个体的观测次数相同（平衡数据的情形）设第i个个体的第j次观测为,1,;1,ijYimjn=，而与ijY对应的p个协变量的观测值分别为1,ijijpXX，并且满足线性回归模型 T11ijijijppijijijYXX=+=+X （1）式中：()T1,ijijijpXX=X；ij为随机误差；()T1,p=为未知的回归系数.令 11121T11121222T12,iii piiiiii piiiininininininpXXXYXXXXYXXXX|=|YX 则模型（1）可写为 iii=+YX （2）更进一步地，如果令 11111,NNpNmmmYXYX|=|YX 式中：N为总观测次数，且满足Nmn=，则模型可进一步写为=+YX （3）在纵向数据中，通常假定不同个体之间的观测是相互独立的，而同一个体内部的不同观测是相关的，并且满足()()2,Var(1,)iiiiEim=YXYV 则Y的均值和协方差矩阵分别为()()2,VarE=YXYV 式中：V为分块对角矩阵，共包含有m个分块对角元 32 高 师 理 科 学 刊 第 43 卷 1.2 加权最小二乘估计 设W是一个对称正定矩阵，最小化随机误差的加权平方和()()()TS=-WYXW YX 将()SW 关于 求导，可以得到估计方程TT=X WXX WY，求解估计方程可得到加权最小二乘估计()1TT-=WX WXX WY 并且容易得到加权最小二乘估计W 的均值和协方差矩阵分别为()()1TTE-=WX WXX WX()()()()11TTTVar,Var=-=WX WXX WYWX X WX 容易看出，对于任意的加权矩阵W，W 都是 的无偏估计.当()1Var-=WVY时，可以进一步得到()()()11T1T1T1,Var-=VVX VXX V YX VX 加权最小二乘估计始终都是回归系数向量 的无偏估计，特别地，如果加权矩阵为单位矩阵，即有=WI，则此时加权最小二乘估计就变成了普通最小二乘估计()1TT-=IX XX Y，其协方差矩阵为()()()11TTVar-=TIX XX VX X X 虽然加权最小二乘估计具有计算简单的优点，但在某些实际情况下，一旦给定响应变量的分布，则加权最小二乘法将不再是最优选择因此，很多实际案例中会考虑用极大似然估计取而代之.由于在确定相关结构矩阵和参数后，得到加权最小二乘估计与极大似然估计两者结果等价因此，此处不再单独对极大似然估计进行介绍.2 标准误差稳健估计 为了推导 的稳健方法，定义加权最小二乘估计()1TT-=?WX WXX WY （4）进一步可得到估计的方差矩阵 ()()1TTTT1-=|WRX WXX W V WX X WX X WX （5）式中：V是估计的方差矩阵无论真正的协方差结构如何，V与V是一致的.因此，可以得到(),WMVNWR?（6）将1-W称为方差工作矩阵，以便将其与真实方差矩阵V区分开来.这与参数模型方法之间的关键区别在于，W的错误选择仅影响对 的推断效率，而不影响其有效性.特别是，无论V的真实形式如何，置信区间和式（6）得出的假设检验都是渐近正确的.假设在g个实验处理组的第h组中，对实验个体中的时点数据进行测量将完整的测量集写为,1,;1,;1,hijYhgim jn=，则得到均值响应的饱和模型为()(1,)hijhjE Yhg=，并且协方差结构的饱和模型是分块对角矩阵()Var=VY，其中所有的非零分块阵等于正定的nn阶矩阵0V.当每个个体的观测时间不同时，以下的修改形式的估计也同样可行，即所需的方差矩阵V仍然是分块对角矩阵，但是与个体内的测量组相对应的非零分块阵在个体之间不再是恒定的.对第i个个体的测量集合的nn方差矩阵为0iV，这些测量的均值向量为i.用i 来估计i，得到0iV的估计()()T0iiiii=-VYY （7）因此，对V的稳健估计则是V 3 案例分析 数据来自美国疾控与预防中心（Centers for Disease Control and Prevention，CDC）的一个儿童铅汞中毒治疗研究（Treatment of Lead and Mercury Poisoning in Children，TLC）523，属于连续型数据.根据有关调查第 2 期 胡倩，等：纵向数据的稳健估计及案例分析 33 结果，铅会对儿童的认知功能造成损害，当儿童血液中的铅含量高于 10 g/dL 时就会存在一定的风险.从医学角度来看，蟹合剂可以用于对铅中毒进行治疗，但是需要注射并且住院治疗，给患者的治疗带来不便.随着医学的发展，1990 年出现了一种新型的口服药，为琥巯酸（Succimer）为了研究琥巯酸能否有效降低血液中的铅含量，TLC 实验组进行了一次随机试验.试验对象是在中心注册的血液中铅含量在 20 40 g/dL 的 1233 个月的儿童，共 100 名将这些儿童随机分成组，一组服用安慰剂，另一组服用琥巯酸他们接受了 3 次治疗，其中每次治疗的一个疗程为 26 d，并且被跟踪研究年.随机挑选 10 名儿童在基线水平、服药后的第 1 周、第 4 周和第 6 周血液中铅含量的数据（见表）.表 1 TLC10 名儿童次测量中血液中铅含量 g/dL 编号 组别 基线水平 第周 第周 第周 编号 组别 基线水平 第周 第周 第周 1 P 30.80 26.90 25.80 23.80 6 A 20.40 05.40 04.50 11.90 2 A 26.50 14.80 19.50 21.00 7 P 28.60 20.80 19.20 18.40 3 A 25.80 23.00 19.10 23.20 8 P 33.70 31.60 28.50 25.10 4 P 24.70 24.50 22.00 22.50 9 P 19.70 14.90 15.30 14.70 5 A 20.40 02.80 03.20 09.40 10 P 31.10 31.20 29.20 30.10 注：P 为安慰剂组；A 为琥巯酸组 由表可以看出，服用安慰剂的儿童血液中的铅含量最高可达 33.70 g/dL，而服用了琥巯酸的儿童血液中的铅含量普遍都在 26.5 g/dL 以下，这表面琥巯酸对于儿童血液中铅含量有抑制减轻的作用.为了进一步验证该结论，计算每组儿童在 4 次测量中血液中铅含量的均值和标准差，结果见表 2.表 2 组儿童次测量中血液中铅含量的均值和标准差 琥巯酸组 安慰剂组 测量时间 均值/gdL-1 标准差 均值/gdL-1 标准差 基线水平 26.54 5.02 26.27 5.02 第周 13.52 7.67 24.66 5.46 第周 15.51 7.85 24.07 5.75 第周 20.76 9.25 23.65 5.64 由表 2 可以看出，随着测量时间的推移，琥巯酸组（A）中血液铅含量的均值和方差呈现的是逐渐递增的趋势，说明服用了琥巯酸组的儿童血液中的药效开始逐渐发挥作用，且铅含量均值逐渐增大接近基线水平，标准差都比基线水平时要大.而安慰剂组（P）中血液铅含量的均值和标准差的趋势波动都表现得较为平缓，说明效果不显著，甚至可以说没有影响.为了比较不同方法下的参数估计，使用几个常用的估计方法分别对问题进行研究.3.1 最小二乘估计 调用 R 软件中的 lm 函数，即可以得到参数的最小二乘估计，运行结果见图，并绘制组儿童血液中铅含量的相关图（见图）.图 最小二乘估计运行结果 图 组儿童血液铅含量相关性图 由图 12 可以看出，10 名儿童血液中铅含量的平均值为 26.17 g/dL，与表 1 中组的均值（琥巯酸组 26.54，安慰剂组 26.27）接近，表明最小二乘估计的效果较为良好在相关图 2 中也显示了组的次测量结果之间存在着较强的相关性.34 高 师 理 科 学 刊 第 43 卷 3.2 极大似然估计与约束极大似然估计 调用 R 软件中的 nlme 包里面的 gls 函数，并在 method 处指定 method=“ML”得到参数的极大似然估计，指定 method=“REML”即可以得到参数的约束极大似然估计.根据运行结果和计算，参数的极大似然估计与约束极大似然估计以及最小二乘估计的结果是一样的，只是极大似然估计的标准误差和比最小二乘估计的标准误差要小一些，因此要更精确.为了更好地体现数据的分布情况，绘制组儿童血液中铅含量的箱线图（见图）图 组儿童血液铅含量箱线图 由图可以看出，第周、第周和第周的均值相对于基线水平的均值而言呈现的是下降趋势，并且从第周开始，直到第周的数据波动较大.4 结语 通过案例分析可以看出，对于一般线性模型下的纵向数据研究而言，模型参数的最小二乘估计、极大似然估计以及约束极大似然估计得到的参数估计结果差别不大甚至基本相等，换句话说，一般线性模型下的纵向数据估计具有稳健性，唯一仅有的差别只是不同的参数估计下的标准误差不同一般而言，标准误差越小，则说明该参数估计方法越有效，因此在处理数据时根据不同参数估计标准差大小的对比，选择合适参数估计方法可提高估计的有效性.参考文献：1 Wang S J，Qian L F，Carroll R J Generalized empirical likelihood methods for analyzing longitudinal dataJ Biometrika，2010，97：79-93.2 You J，Chen G，Zhou YBlock empirical likelihood for longitudinal partially linear regression modelsJCanadian Journal of Statistics，2006，34：79-96.3 Lin D F，Ying Z Semiparametric and nonparametric regression analysis for longitudinal dataJ Journal of the American Statistical Association，2001，96：103-126 4 Xue L G，Zhu L XEmpirical likelihood semiparametric regression analysis for longitudinal dataJBiometrika，2007，94：921-937 5 王友乾纵向数据分析M北京：高等教育出版社，2015 6 Laird N，Ware JRandom-e ects models for longitudinal dataJBiometrics，1982，38（4）：963-974 7 Liang K Y，Zeger S LLongitudinal Data Analysis Using Generalized Linear ModelsJBiometrika，1986，73：13-22 8 Qu A，Lindsy B G，Li BImproving generalized estimating equations using quadratic inference functionsJBiometrika，2000，87：823-836 9 关晓妮，黄彬 纵向数据下广义线性模型的稳健二次推断函数估计J 北京化工大学学报（自然科学版），2018，45（2）：100-104 10 芦飞纵向数据下均值-协方差模型的估计理论D北京：北京工业大学，2020