半参数纵向数据的Logistic模型

cij
=
lnççèæ
nij
yij -
+ yij
0.5 + 0.5
÷÷øö,
i =1,2,L , m, j =1,2,L , ni ,
１，３，１８，３３，３７，４１，４２，４７，５９，６２，６７，７５，７６，８３，８９。协变 量 ｘｉｊ１～Ｎ （μ１，τ２１），ｘｉｊ２～Ｎ （μ２，τ２２），其中 μ１＝２，τ２１＝０．０１， μ２＝３．５，τ２２＝０．０１。
１ 参数估计
同一般的线性回归一样，首先我们要利用观测数据 对模型中的参数 βｐ 及函数ｇ（ｔｉｊ）作出估计。函数ｇ（ｔｉｊ）是 关于时间ｔｉｊ 的函数，函数ｇ（ｔｉｊ）可由ＴＴｉｊβｑ 逼近，即设
（２）
其中Ｔｉｊ＝（Ｂ１（ｔｉｊ），Ｂ２（ｔｉｊ），…，Ｂｑ（ｔｉｊ））Ｔ 是基函数向量。 设对变量ｘｉｊ＝（Ｘｉｊ１，Ｘｉｊ２，…，Ｘｉｊｐ－１）Ｔ 给定了 ｍ×ｎｉ 组值，
和 Ｖ（１）代人式（１６），求 β（２）。 ④重复步骤 ３。直到 ｔ＝ｔ０ 时，对给定的误差限 ε，有
设 β（ｔ）为第 ｔ 次迭代所得的 β值，代入上式可得 Ｕ（ｔ） 和 Ｖ（ｔ）。在 β（ｔ）处将 ｈ（β）按 ｐ＋ｑ 元函数的 Ｔａｙｌｏｒ 公式展 开，并取至二次项（记为 Ｑ（ｔ）（β）），有
这时 β的最大似然估计 ⑤由（１５）的第二个等式得，在迭代过程中，我们可 同时得到信息矩阵的估计为
第 ２１ 卷第 ２ 期 ２００８ 年 ４ 月
四川理工学院学报（ 自然科学版） ＪＯＵＲＮＡＬ ＯＦ ＳＩＣＨＵＡＮ ＵＮＩＶＥＲＳＩＴＹ ＯＦ ＳＣＩＥＮＣＥ ＆ ＥＮＧＩＮＥＥＲＩＮＧ（ ＮＡＴＵＲＡＬ ＳＣＩＥＮＣＥ ＥＤＩＴＩＯＮ）
文章编号： １６７３－ １５４９（ ２００８） ０２－ ０００８－ ０４
从而可得β! 的协方差矩阵的估计为 （１７）
它是 Ｎｅｗｔｏｎ－Ｒａｐｈｓｏｎ 迭代计算的一个很有用的副 产品。 ３． ３ 初值 β（０） 的确定
理论上，虽然初值 β（０）可以任意给定（如取 β（０）＝（０，０ …，０）Ｔ），但是，初值对于达到指定精度 ε所需要的迭代 次数 ｔ０ 有较大的影响。因此，初值 β（０）的确定有时是很重 要的。下面介绍一种常用的迭代初值 β（０）的选取方法，即 利用所谓修正的经验 Ｌｏｇｉｓｔｉｃ 变换确定 β（０），具体方法 如下： 令
从而对数似然函数为 由于（
（７）
（１０）
β的信息矩阵是上述 （ｐ＋ｑ）２ 个二阶偏导数所组成
（ｐ＋ｑ）×（ｐ＋ｑ）方阵的负矩阵的期望，记为 Ｉ（β）。Ｉ（β）中
各元素与随机变量
的观测值
无关，故
１０
四川理工学院学报（ 自然科学版）
２００８ 年 ４ 月
由式（８）及式（１０）可得当 β＝β（ｔ）时，
（１１） 根据最大似然估计理论，β的最大似然估计β! 的协 方差矩阵为
（Ｙ１１＝ｙ１１，Ｙ１２ ＝ｙ１２， … ，
Ｙ２１＝ｙ２１，
），
的似然函数为
则得似然方程
（９） 关于 β求解方程组（９），便得 β的最大似然估计β! 。 但式（９）关于 β是非线性方程组，所以我们要用迭代方 法求解。
２ 信息矩阵及其估计
下面我们在（８）式的基础上，求对数似然函数关于 参数的二阶偏导数，并给出 β的信息矩阵及其估计式。
参 考 文 献： ［１］ Ｄａｌｔｏｎ Ｆ Ａｎｄｒａｄｅ，Ｈｅｌｉｔｏｎ Ｒ Ｔａｖａｒｅｓ． Ｉｔｅｍ ｒｅｓｐｏｎｓｅ ｔｈｅｏｒｙ
ｆｏｒ ｌｏｎｇｉｔｕｄｉｎａｌ ｄａｔａ：ｐｏｐｕｌａｔｉｏｎ ｐａｒａｍｅｔｅｒ ｅｓｔｉｍａｔｉｏｎ ［Ｊ］． Ｊｏｕｒｎａｌ ｏｆ Ｍｕｌｔｉｖａｒｉａｔｅ Ａｎａｌｙｓｉｓ，２００５，９５：１－ ２２． ［２］ Ｑｉｎ Ｇｕｏｙｏｕ，Ｚｈｕ Ｚｈｏｎｇｙｉ．Ｒｏｂｕｓｔ ｅｓｔｉｍａｔｉｏｎ ｉｎ ｇｅｎｅｒａｌｉｚｅｄ ｓｅｍｉｐａｒａｍｅｔｒｉｃ ｍｉｘｅｄ ｍｏｄｅｌｓ ｆｏｒ ｌｏｎｇｉｔｕｄｉｎａｌ ｄａｔａ ［Ｊ］． Ｊｏｕｒｎａｌ ｏｆ Ｍｕｌｔｉｖａｒｉａｔｅ Ａｎａｌｙｓｉｓ．（２００７）， ｄｏｉ： １０．１０１６／ｊ．ｊｍｖａ． ２００７． ［３］ 张尧庭．定性资料的统计分析［Ｍ］．桂林：广西师范大学出版 社， １９９１． ［４］ 范金城，梅长林．数据分析［Ｍ］．北京：科学出版社，２００２． ［５］ 田 萍，薛留根．纵向数据半参数回 归 模 型 估 计 的 强 相 合 性 ［Ｊ］．工程数学学报，２００６，２３（２）：３６９－ ３７２． ［６］ 孙 孝 前，尤 进 红．纵 向 数 据 半 参 数 建 模 中 的 迭 代 加 权 偏 样 条最小二乘估计［Ｊ］．中国科学，２００３，３３（５）：４７０－ ４８０． ［７］ 余 火 军，朱 仲 义．纵 向 数 据 模 型 均 值 参 数 和 方 差 参 数 的 影 响分析［Ｊ］．应用概率统计，２００５，２１（１）：４０－ ５２． ［８］ 林 金 官，韦 博 成．非 线 性 纵 向 数 据 模 型 中 方 差 和 自 相 关 系 数的齐性检验［Ｊ］．应用数学学报，２００４，２７（３）：４６６－ ４８０． ［９］ 钱 伟 民，柴 根 象．纵 向 数 据 混 合 效 应 模 型 的 统 计 分 析［Ｊ］．数 学年刊， ２００４，（２５）：４４５－ ４５６．
（１２） 若利用方程（９）求得 β的最大似然估计β! ，代入 π（Ｄｉｊ）的 表达式便得 π（Ｄｉｊ）的估计
其中
，这时，
从而可得
的估计为
。
３ 似然方程的 Ｎｅｗｔ ｏｎ－ Ｒａｐｈｓ ｏｎ 迭代求解
３． １ Ｎｅｗｔ ｏｎ－ Ｒａｐｈｓ ｏｎ 迭代法的一般描述
设 ｈ（β）是关于
，的 ｐ＋ｑ 元函
数，我们要求β! ，使
我们考虑纵向数据半参数回归模型
且
（３）
在我们研究的问题中，我们可能更关心因变量是一 个二值变量，即ｙｉｊ 只取 ０ 与 １ 两个值时，Ｙｉｊ＝１ 的概率 Ｐ＝ Ｐｒ （ｙｉｊ＝１）。因此我们建立如下的纵向数据的半 参 数 Ｌｏｇｉｓｔｉｃ 模型。
（１）
其中 （ｘ!ｉｊ，ｔｉｊ）∈Ｒｐ×Ｒ 是已知的设计点列，βｐ 是 ｐ 维未 知参数，ｇ（·）是定义在区间 Ｉ 上的未知回归函数。Ｅ（εｉｊ） ＝０，０＜Ｖａｒ（εｉｊ）＝σ２＜∞，当ｉ≠ｊ 时，εｉｋ 与 εｉｒ 相互独立。其 中 ｙｉｊ 是第ｉ个体在时间ｔｉｊ 的响应值。记ｎｉ 为有界正整数， Ｎ＝ｍｎｉ。
则初值 β（０）可取为 （１８）
４ 随机模拟
我们从下面的半参数 Ｌｏｇｉｓｔｉｃ 回归模型产生纵向 数据，
其中
，随机变量 εｉｊ 假设为独立
同wk.baidu.com布，且服从正态分布
。β和 σ２ 的真实值为
β＝ （－１，１，－０．５，１．５），σ２＝０．０１。我们观测 的 个 体 数 为
ｍ＝２０，每个个体进行了 ｎｉ＝ １５ 次观测，ｎｉ 的取值时间为
我们利用上面模型产生的随机数据，根据纵向数据 Ｌｏｇｉｓｔｉｃ 模型参数最大似然估计的 Ｎｅｗｔｏｎ－Ｒａｐｈｓｏｎ迭 代 法 编 程 ，我 们 选 两 组 初 值 ，一 组 是 β（０）＝ （－０．５，０．５， ０．５，０．５），另一组是 β（０）＝（－１．５，０．２，０．２，１），经过 １５ 次 迭代都收敛到 β＝（－０．９９４３，１．０２１１，－０．５３６５，１．６０７４）， 则得到 β的估计值 β＝（－０．９９４３，１．０２１１，－０．５３６５，１．６０７４）。 可以看出，在迭代次数不多的情况下，收敛到 β估计值； 而且偏差比较小。
上面的模型我们考虑的ｙｉｊ 是多个输出或连续输出 的情形，但在有些情况下，我们考虑ｙｉｊ 只取两个值。如在 医药和其它行业的试验中，我们可以用ｙｉｊ＝０ 表示第 ｉ个 体在ｔｉｊ 时刻正常，用ｙｉｊ＝１ 表示不正常。因此我们可以建立 如下纵向数据的半参数模型
（４） 其中（!ｘ ｉｊ，ｔｉｊ）∈Ｒｐ×Ｒ 是已知的设计点列，βｐ 是ｐ 维未知参 数，ｇ（·）是定义在区间 Ｉ 上的未知回归函数。Ｅ（εｉｊ）＝０， ０＜Ｖａｒ（εｉｊ）＝σ２＜∞，当ｉ≠ｊ 时，εｉｋ 与 εｉｒ 相互独立。π（ｘｉｊ， ｔｉｊ）是随机变量 Ｙｉｊ 取 １ 时的概率大小。
半参数纵向数据的 Ｌｏｇｉ ｓ ｔ ｉ ｃ 模型
２Ｖ００ｏ８ｌ．２年１ Ｎ４ｏ月．２ Ａｐｒ．２００８
魏强强， 魏立力
（ 宁夏大学数学计算机学院， 银川 ７５００２１）
摘 要： 以纵向数据的半参数模型为基础， 建立了半参数纵向数据的 Ｌｏｇｉ ｓ ｔ ｉ ｃ 回 归 模 型 ， 对 此 模 型 中 的 参 数
，
（１３）
令
（１５）
将式（１５）代入（１４）可得 Ｌｏｇｉｓｔｉｃ 模型参数 β的最大似 然估计的 Ｎｅｗｔｏｎ－Ｒａｐｈｓｏｎ 迭代式为：
（１６）
具体计算步骤如下： ①给定初值 β（０），计算 π（０）（Ｄｉｊ），并计算 Ｕ（０）和 Ｖ（０）。 ②将 β（０），Ｕ（０）和 Ｖ（０）代人式（１６），求 β（１）。 ③由 β（１）计算 π（１）（Ｄｉｊ），并计算 Ｕ（１）和 Ｖ（１）。将 β（１），Ｕ（１）
（５） 故对于 bl (l = 0,1,2,L , p + q -1)
令
，
其中
再令
，
则相应的 Ｌｏｇｉｓｔｉｃ 模型为
其中
。令
或
９ （８）
（６）
其中
为未知参数。
设在 ｘ＝ｘｉｊ 的 ｎｉｊ 次独立观测中，事件 Ａ 发生了 ｙｉｊ 次
（即 Ｙｉｊ＝ｙｉｊ），ｉ＝１，…，ｍ， ｊ＝１，…，ｎｉ。则由式（５）和（６）知
令
可解得 β（看作 β的第 ｔ＋１ 次迭代所得的值，记为 β（ｔ＋１）） 为
（１４） 这便是求解问题（１３）的一般 Ｎｅｗｔｏｎ－Ｒａｐｈｓｏｎ 迭代公式。 ３． ２ Ｌｏｇｉ ｓ ｔ ｉ ｃ 模型参数的最大似然估计的 Ｎｅｗｔ ｏｎ－ Ｒａｐｈｓ ｏｎ 迭代法
这时，对数似然函数（７）为其目标函数，即
进行了参数估计， 讨论了它的信息矩阵， 给出了似然方程的 Ｎｅｗｔ ｏｎ－ Ｒａｐｈｓ ｏｎ 迭代求解过程。
关键词： Ｌｏｇｉ ｓ ｔ ｉ ｃ 回归； 纵向数据； Ｎｅｗｔ ｏｎ－ Ｒａｐｈｓ ｏｎ 迭代
中图分类号： Ｏ２１２．４
文献标识码： Ａ
引言与模型
纵向数据是指对每一组个体在不同时间进行观测 而得到的由截面和时间序列融合在一起的数据，其特点 是将截面数据和时间序列数据结合在一起，既能分析出 个体随时间变化趋势，又能分析总体的变化趋势。近年 来，对纵向数据各种模型的研究引起国内外统计学者的 广泛关注。本文用 Ｌｏｇｉｓｔｉｃ 回归对纵向数据建立模型， 可用此模型分析纵向数的内部结构，及个体间的差异。
B1(t11) B1(t12)
L
B1(t1n1 ) B1(t21)
L B1(tmnm )
B2(t11) L B2(t12) L
LL
B2(t1n1 ) L B2(t21) L
LL B2(tmnm ) L
Bq(t11) ù
Bq(t12)
ú ú
Lú
Bq
(t1n1
)
ú ú
Bq(t21)
ú ú
Lú
Bq(tmnm )úû
Ｌｏｇｉｓ ｔｉｃ Ｍｏｄｅ ｌ ｏｆ Ｌｏｎｇｉｔｕｄｉｎａ ｌ Ｄａ ｔａ
ＷＥＩ Ｑｉａｎｇ－ｑｉａｎｇ， ＷＥＩ Ｌｉ－ｌｉ （Ｃｏｌｌｅｇｅ ｏｆ Ｍａｔｈｅｍａｔｉｃｓ ａｎｄ Ｃｏｍｐｕｔｅｒ Ｓｃｉｅｎｃｅ， Ｎｉｎｇｘｉａ Ｕｎｉｖｅｒｓｉｔｙ， Ｙｉｎｃｈｕａｎ ７５００２１， Ｃｈｉｎａ）
第 ２１ 卷第 ２ 期
魏强强等： 半参数纵向数据的 Ｌｏｇｉ ｓ ｔ ｉ ｃ 模型
１１
é x110
ê ê
x120
êL
X
=
ê ê
x1n10
ê ê
x210
êL
êëxmnm0
x111 L x121 L LL
x1n11 L x211 L LL xmnm1 L
x11p-1 x12p-1 L
x1n1p-1 x21p-1 L xmnmp-1
收稿日期： ２００７－１０－１６ 基金项目： 国家自然科学基金项目（ ６０６６３００３） 作者简介： 魏强强（ １９８１－） ， 男， 甘肃静宁人， 硕士生， 主要从事统计学方面的研究。
第 ２１ 卷第 ２ 期
魏强强等： 半参数纵向数据的 Ｌｏｇｉ ｓ ｔ ｉ ｃ 模型
ｘｉｊ＝（Ｘｉｊ１，Ｘｉｊ２，…，Ｘｉｊｐ－１）Ｔ，ｉ＝１，…，ｍ，ｊ＝１，…，ｎｉ，对于其中的 ｉｊ 则 组直，共独立观测了ｎｉｊ 次，ｉ＝１，…，ｍ，ｊ＝１，…，ｎｉ。令 Ｙｉｊ 为 在对 ｘｉｊ 的ｎｉｊ 次观测中事件 Ａ 发生的次数，以 π（ｘｉｊ，ｔｉｊ） 记在 ｘ＝ｘｉｊ 且ｔ＝ｔｉｊ 时事件 Ａ 发生的概率，则 Ｙｉｊ 服从参数 为ｎｉｊ 和 π（ｘｉｊ，ｔｉｊ）的二项分布。即