
Logit模型
-
2023年2月26日发(作者:环境友好型)Logistic回归模型
1Logistic回归模型的基本知识
1.1Logistic模型简介
主要应用在研究某些现象发生的概率p,比如股票涨还是跌,公司成功或失败的概率,以及讨论概率
p与那些因素有关。显然作为概率值,一定有10p,因此很难用线性模型描述概率p与自变量的关
系,另外如果p接近两个极端值,此时一般方法难以较好地反映p的微小变化。为此在构建p与自变量关
系的模型时,变换一下思路,不直接研究p,而是研究p的一个严格单调函数)(pG,并要求)(pG在p接
近两端值时对其微小变化很敏感。于是Logit变换被提出来:
p
p
pLogit
1
ln)((1)
其中当p从10时,)(pLogit从,这个变化范围在模型数据处理上带来很大的方便,
解决了上述面临的难题。另外从函数的变形可得如下等价的公式:
X
T
X
T
T
e
e
pX
p
p
pLogit
1
1
ln)((2)
模型(2)的基本要求是,因变量(y)是个二元变量,仅取0或1两个值,而因变量取1的概率
)|1(XyP就是模型要研究的对象。而T
k
xxxX),,,,1(
21
,其中
i
x表示影响y的第i个因素,它可以
是定性变量也可以是定量变量,T
k
),,,(
10
。为此模型(2)可以表述成:
k
x
k
x
k
x
k
x
kke
e
pxx
p
p
110
110
1101
1
ln(3)
显然pyE)(,故上述模型表明
)(1
)(
ln
yE
yE
是
k
xxx,,,
21
的线性函数。此时我们称满足上面条件
的回归方程为Logistic线性回归。
Logistic线性回归的主要问题是不能用普通的回归方式来分析模型,一方面离散变量的误差形式服从伯
努利分布而非正态分布,即没有正态性假设前提;二是二值变量方差不是常数,有异方差性。不同于多元线
性回归的最小二乘估计法则(残差平方和最小),Logistic变换的非线性特征采用极大似然估计的方法寻求最佳
的回归系数。因此评价模型的拟合度的标准变为似然值而非离差平方和。
定义1称事件发生与不发生的概率比为优势比(比数比oddsratio简称OR),形式上表示为
OR=k
x
k
xe
p
p
110
1
(4)
定义2Logistic回归模型是通过极大似然估计法得到的,故模型好坏的评价准则有似然值来表征,称
-2
ˆ
ln()L为估计值
ˆ
的拟合似然度,该值越小越好,如果模型完全拟合,则似然值
ˆ
()L为1,而拟合似然
度达到最小,值为0。其中
ˆ
()lnL表示
ˆ
的对数似然函数值。
定义3记
)
ˆ
(Var
为估计值
ˆ
的方差-协方差矩阵,2
1
)]
ˆ
([)
ˆ
(VarS
为
ˆ
的标准差矩阵,则称
ki
S
w
ii
i
i
,,2,1,]
ˆ
[2
(5)
为
i
ˆ
的Wald统计量,在大样本时,
i
w近似服从)1(2分布,通过它实现对系数的显著性检验。
定义4假定方程中只有常数项
0
,即各变量的系数均为0,此时称
2
0
ˆˆ
2[ln()ln()]LL(6)
为方程的显著性似然统计量,在大样本时,2近似服从)(2k分布。
1.2Logistic模型的分类及主要问题
根据研究设计的不同,Logistic回归通常分为成组资料的非条件Logistic回归和配对资料的条件
Logistic回归两种大类。还兼具两分类和多分类之分,分组与未分组之分,有序与无序变量之分。具体如
下:
两分类非条件Logistic回归:分组数据的Logistic回归,未分组数据的Logistic回归;
多分类非条件Logistic回归:无序变量Logistic回归,无序变量Logistic回归;
条件Logistic回归:1:1型、1:M型和M:N型Logistic回归。
关于Logistic回归,主要研究的内容包括:
1.模型参数的估计及检验
2.变量模型化及自变量的选择
3.模型评价和预测问题
4.模型应用
2Logistic模型的参数估计及算法实现
2.1两分类分组数据非条件Logistic回归
因变量(反应变量)分为两类,取值有两种,设事件发生记为y=1,不发生记为y=0,设自变量
T
k
xxxX),,,(
21
是分组数据,取有限的几个值;研究事件发生的概率)|1(XyP与自变量X的关
系,其Logistic回归方程为:
kk
xx
XyP
XyP
110)|0(
)|1(
ln或
k
x
k
x
k
x
k
x
e
e
XyP
110
110
1
)|1(
例2.1.1分组数据[1]在一次住房展销会上,与房地产商签订初步购房意向书的有n=325人,在随后的3个月
时间内,只有一部分顾客购买了房屋。购买房屋的顾客记为1,否则记为0。以顾客的年家庭收入(万元)作为
自变量X,对数据统计后如表2.1.1所示,建立Logistic回归模型。
表2.1.1购房分组数据
序号
年家庭收入
X(万元)
签订意
向人数
实际购
买人数
11.5258
22.53213
33.55826
44.55222
55.54320
66.53922
77.52816
88.52112
99.51510
例2.1.2药物疗效数据[2]为考察某药物疗效,随机抽取220例病人并分配到治疗组和对照组,治疗组采用治
疗药物,对照组采用安慰剂。治疗一段时间后观察病人的疗效,得到表2.1.2数据。设y为疗效指标(y=1有
效,y=0无效),
1
x为治疗组指标(1为治疗组,0为对照组),
2
x为年龄组指标(1为>45岁,0为其他)。
表2.1.2药物疗效数据
序号
治疗分组
1
x年龄分组
2
x有疗效无效合计
111321850
210402060
301213152
400184058
上述两个例子数据都是经过统计加工后的分组数据,对此类数据进行Logistic回归,首先要明确应变
量对应事件的发生概率如何确定和进行Logit变换,其次才能建立Logistic回归。为便于数据处理,我们将
此类数据的格式作个约定,排列格式为(组序号,自变量X,该组事件发生数,该组总例数)。
表2.1.3分组数据的标准格式
表2.1.1改造表
序
号
年家庭收入
X(万元)
实际购买
人数
i
m
签订意向
总人数
i
n
11.5825
22.51332
33.52658
44.52252
55.52043
66.52239
77.51628
88.51221
99.51015
表2.1.2改造表
序
号
治疗分
组
1
x
年龄分
组
2
x
有效例
数
i
m
观察例
数
i
n
1113250
2104060
3012152
4001858
经过改造后,可得我们关心的事件的发生的频率为
ni
n
m
p
i
,,2,1,
i
i
该组总例数
该组发生事件数
。其中
n
为
分组数,然后作Logit变换,即
i
i
iip
p
pLogitp
1
ln)(
~
。变换后的数据,形式上已经可以采用一般的线
性回归的处理方式来估计回归参数了。此时方程变为:
k
j
ijji
nixp
1
0
,,2,1,
~
当然这样处理并没有解决异方差性,当
i
n较大时,
i
p
~
的近似方差为:
)(,
)1(
1
)
~
(
ii
iii
i
yE
n
pD
(7)
所以选择权重
nippn
iiii
,,2,1),1(,最后采用加权最小二乘法估计参数。
注意,分组数据的Logistic回归只适用于大样本分组数据,对小样本的为分组数据不适用,并且以组数
n
为回归拟合的样本量,明显降低了拟合精度,在实际应用中必须谨慎。
求解算法及步骤:
1.依据分组数据的标准格式,计算频率
i
p、Logit变换
i
p
~
和权重
i
2.构建加权最小二乘估计:
n
i
k
j
ijjiiii
n
i
k
j
ijjii
xyxy
11
2
0
11
2
0
)(min)(min
(8)
令
iii
yy*,T
ikiiiii
xxX),,,(
1
*,T
k
),,,(
10
则方程又变成一般的线性回归模型:
n
i
i
T
i
Xy
1
2**)(min(9)
3.构造增广矩阵
21
****][
kk
TTYXXX利用消去法得]
ˆ
)
ˆ
([VarI矩阵,得到估计
ˆ
其中
2,1KK
I为残差平方和SE,回归方差
1
ˆ2
kn
SE
各系数检验采用)1(~
ˆ
ˆ
knt
I
t
ii
i
i
总平方和
n
i
n
i
i
n
i
ii
ii
y
yST
1
1
2
1
2
2
)(
)(
,回归平方和SESTSR
总平方和求解相当于拟合
ii
y*
0
*方程的残差平方和,故得上式ST
所以方程的检验为)1,(~
)1/(
/
knkF
knSE
kSR
F
例2.1.1的求解过程如下(由LLLStat统计软件计算):
表2.1.4数据Logit变换及权重
家庭年收入x实际购买mi签订意向ni比例pi逻辑变换Logit权重ni*pi(1-pi)
1.5000008250.320000-0.7537725.440000
2.5.406250-0.3794907.718750
3.5.448276-0.20763914.344828
4.5.423077-0.31015512.692308
5.5.465116-0.13976210.697674
6.5.5641030.2578299.589744
7.5.5714290.2876826.857143
8.5.5714290.2876825.142857
9.5.6666670.6931473.333333
表2.1.5回归模型基本信息
总样本9
求解方法加权最小二乘
仅常数项beta0-0.095029
方程F统计量51.982160
F分布自由度1,7
方程检验p值0.000176
总平方和8.798294
回归平方和7.754112
残差平方和1.044181
表2.1.6分组Logistic回归系数检验
序号均值回归系数系数标准误t统计量自由度df检验P值
常数项2.837815-0.8488820.113578-7.47399470.000056
家庭年收入x14.9011400.1493230.0207117.20986570.000056
表2.1.7
1][XXT
0.086479-0.014517
-0.0145170.002876
本例Logistic模型的回归方程:x
e
x
e
p
i
149323.0848882.0
149323.0848882.0
1
ˆ
对于多分类无序自变量的Logistic回归,即某个自变量为m个水平的名义变量(如治疗方法
A,B,C),只需要引入m-1(2个)个哑变量,然后采用上述方法进行分析。
例2.1.3研究三种治疗方法对不同性别病人的治疗效果[2],数据如表2.1.4
表2.1.4性别和治疗法对某病治愈情况的影响
性别治疗方法
有效
i
m
无效
总例数
i
n
男
A78
28
106
B10111112
C6846114
女
A40545
B54559
C34640
由于治疗方法有三种,没有等级关系,所以属于无序的名义变量,故引入两个哑变量
32
,xx分别代表A
和B疗法,其中
0,1
32
xx表示方法A,1,0
32
xx表示方法B,0,0
32
xx表示方法C,将上述数
据转化成标准格式,得表2.1.5。
表2.1.5性别和治疗法对某病治愈情况的影响
性别
1
x
2
x
3
x有效
i
m总例数
i
n
1
1078106
101101112
10068114
0104045
0015459
0003440
对于分类数据,也可以采用极大似然法进行参数估计,具体见2.2节最后部分内容。
2.2两分类未分组(连续)非条件Logistic回归
应变量y取值为0和1,设事件发生记为y=1,否则为0,设自变量T
k
xxxx),,,(
21
,n组观测数据
记为
),,,,(
21iikii
yxxx,ni,,2,1。记T
ikiii
xxxX),,,,1(
21
,1
0
i
x,则
i
y与
ikii
xxx,,,
21
的
Logistic回归模型是:
i
X
T
i
X
T
ik
x
ki
x
ik
x
ki
x
ikkiii
e
e
e
e
xxfyE
1
1
)()(
110
110
110
(10)
易知,
i
y是均值为
i
的0-1型分布,其分布律为
i
y
i
i
y
ii
yf1)1()(,niy
i
,,2,1;1,0
则
n
yyy,,,
21
的似然函数和对数似然函数分别为:
n
i
i
y
i
i
y
i
L
1
1)1(
n
i
i
i
i
i
n
i
iiii
yyyL
11
)]1ln(
1
ln[)]1ln()1(ln[ln
代入
ik
x
ki
x
ik
x
ki
x
ie
e
110
110
1
,得
n
i
i
X
T
i
T
i
n
i
ik
x
ki
x
ikkii
eXy
exxyL
1
1
110
110
)]1ln([
)]1ln()([ln
(11)
记)(ln)(LLL,选取T
k
),,,(
10
的估计T
k
)
ˆ
,,
ˆ
,
ˆ
(
ˆ
10
使得)(LL达到极
大,这就是Logistic回归模型的极大似然估计,该过程的求解需要采用牛顿迭代法。
构造得分函数
kg
LL
F
g
g
,,2,1,0,
)(
)(
,共k+1个非线性方程组,令其=0求解,其中
kg
e
ex
xyF
n
i
i
X
T
i
X
T
ig
igig
,,2,1,0,]
1
[)(
0
(12)
构造信息矩阵
khg
hg
LL
I
gh
,,2,1,0,,
)(
)(
2
,即)(LL二阶导矩阵的负矩阵,其中
khg
e
exx
I
n
i
i
X
T
i
X
T
ihig
gh
,,2,1,0,,
)1(
)(
0
2
(13)
很明显
)()(
hggh
II,故)(I是一个对称矩阵。
求解算法及步骤:
1.根据公式(12)计算得分函数)(
g
F,公式(13)计算信息矩阵)(
gh
I
给定初值)0,,0,0(00,k=1和精度,可取0.000001
2.采用牛顿迭代式1kk,)()]([111kkFI,通过以下方式求解。
构造增广矩阵
)(1kIF=))()((11kkFI,通过对IF矩阵作k+1次ij消去变换求解
若
k
g
g
0
2||||或者
k
g
g
0
||||||或者
|}{|max
0
g
kg
,则转3
否则k=k+1,继续执行第2步
3.此时k就是回归系数的数值估计
ˆ
,k就是迭代次数,消去变换后的IF矩阵的前
11kk子阵就是方差-协方差矩阵的估计阵
11
)()
ˆ
(
kkgh
VVar=V,下面给出检验有关计算:
计算Wald统计量
gg
g
gV
W
2ˆ
,近似服从)1(2分布,检验p值))1((2
gg
WPp
标准误
ggg
VES).(.,g
g
eORˆ
)(,kg,,1,0
例2.2.1公共交通调查数据[1]在一次关于公共交通的社会调查中,调查项目为“是乘坐公共汽车上下班,还
是骑自行车上下班”。因变量y=1表示乘坐公共汽车,y=0表示骑自行车。自变量
1
x是年龄,作为连续变
量;
2
x是月收入(元);
3
x是性别,
3
x=1表示男性,
3
x=0表示女性。调查对象为工薪族群体,数据如表
2.2.1所示。
表2.2.1公共交通社会调查
序号年龄
1
x月收入
2
x性别
3
x交通y
11885000
221120000
32385001
42395001
528120001
63185000
736150001
842100001
94695001
1048120000
1155180001
1256210001
1358180001
141885010
1520100010
1625120010
1727130010
1828150010
193095011
2032100010
2133180010
2233100010
2338120010
2441150010
2545180011
2648100010
2752150011
2856180011
以下计算结果采用LLLStat1.0软件得到:
表2.2.2主要计算结果
序号均值回归系数系数标准误wald统计量自由度df检验p值OR=Exp(B)
常数项0.535714-3.6550162.0912233.05476610.0805010.025861
年龄1273.2142860.0821680.0521192.48551610.1148991.085639
月收入0.4642860.0015170.0018650.66146610.4160431.001518
性别36.107143-2.5018441.1578184.66917510.0307090.081934
表2.2.3Logistic模型基本信息
总样本28
求解方法极大似然法&Newton迭代
迭代次数(仅beta0)7(4)
-2LogLikelihood(Beta)25.970652
仅常数项beta0-0.143101
-2LogLikelihood(beta0)38.673263
方程Wald值(相减)12.702611
方程自由度4
方程检验p值0.012824
对于例2.1.3分组数据的极大似然估计法,主要过程如下:
n
i
i
m
i
n
i
i
m
i
m
n
i
i
CL
1
)1(
n
i
ii
i
i
i
m
n
n
i
iiiii
m
n
nmCmnmCLi
i
i
i
11
)]1ln(
1
ln[ln)]1ln()(ln[lnln
代入
ik
x
ki
x
ik
x
ki
x
ie
e
110
110
1
,得
n
i
i
X
T
ii
T
i
m
n
enXmCLi
i
1
)]1ln([lnln
则有
g
g
LL
F
)(
)(
kg
e
ex
nxm
n
i
i
X
T
i
X
T
ig
iigi
,,2,1,0,]
1
[
1
hg
LL
I
gh
)(
)(
2
khg
e
exxnn
i
i
X
T
i
X
T
ihigi,,2,1,0,,
)1(1
2
;
其中
i
m,
i
n分别表示分组i中事件发生次数和总观察数,如表2.1.4和2.1.5所示。然后可采用Newton-
Raphson迭代法进行求解。由LLLStat计算得到如下结果。
表2.2.4性别和疗法对某病治愈的影响(未分组Logistic似然估计法)
序号均值回归系数系数标准误wald统计量自由度df检验P值
常数项1.0000001.4183990.29869022.55051310.000002
性别0.500000-0.9616180.29979710.28847210.001339
治疗A0.3333330.5847450.2641084.90196610.026826
治疗B0.3333331.5607630.31596124.40099310.000001
表2.2.5回归系数方差矩阵V(beta)(信息矩阵I(Beta)的逆矩阵)
0.089215-0.072957-0.029931-0.030097
-0.0729570.089878-0.0000780.000128
-0.029931-0.0000780.0697530.029993
-0.0300970.0001280.0299930.099831
2.3条件Logistic回归[2,3]
条件Logistic回归是配对设计(病例-对照)中常用的一种统计分析方法,通过配对方法收集资料:每一配
对组可包括一个病例和一个或多个对照,有1:1型、1:m型配对。假设收集了如下数据:
表2.3.1n个1:m配对组,k个协变量的比例资料
配对组号
病例组0X第1对照组1X…
第m个对照组mX
10
1
0
12
0
11
,,,
k
xxx1
1
1
12
1
11
,,,
k
xxx
…m
k
mmxxx
11211
,,,
20
2
0
22
0
21
,,,
k
xxx1
2
1
22
1
21
,,,
k
xxx
…m
k
mmxxx
22221
,,,
……………
n00
2
0
1
,,,
nknn
xxx11
2
1
1
,,,
nknn
xnxx
…m
nk
m
n
m
n
xxx,,,
21
配对资料用配对的方法来控制影响因素的干扰,并且每个配对组都可以建立一个Logistic回归方程:
nixxpLogit
kk
i,,2,1,)(
110
为此需要估计的参数有n个常数项n
0
1
0
,,
和k个回归系数
k
,,
1
,配对数越多估计的参数就越
多,但是一般的数据量难以支撑这样的估计,故一般的Logistic回归不适合配对资料。不过在参数估计时,
常数项会被消去,所以方程组减少了n个常数项n
0
1
0
,,
的估计,复杂度大大降低。对于回归参数的估计
采用条件似然函数替代一般的似然函数进行。
对于第i个配对组而言,共有m+1个观察对象,记为
m
BBBA,,,,
21
,其中仅有一例发病,且正好是
病例组A发病,而对照组均没有发病的条件概率
i
p(类似Bayes概率)可以表示成:
m
j
mjm
m
i
BBBAPBBBAP
BBBAP
p
1
121
21
)()(
)(
(14)
其中
)(
21m
BBBAP
=
)|0()|0()|1(1100m
i
m
iiiii
XyPXyPXyP
,而
j
i
X
T
j
i
X
T
j
i
j
i
e
e
XyP
1
)|1(,
j
i
X
T
j
i
j
i
e
XyP
1
1
)|0(,mj,,2,1(15)
故n个配对组的条件似然函数表示为:
n
i
m
j
i
X
j
i
X
T
n
i
m
j
i
X
j
i
X
T
n
i
m
k
m
kj
j
j
i
X
T
k
i
X
T
k
i
X
T
i
X
T
m
j
j
i
X
T
i
X
T
i
X
T
m
j
j
i
X
T
i
X
T
i
X
T
e
e
e
e
e
e
e
e
e
e
e
e
L
1
1
1
)
0
(
1
1
)
0
(
1
1
1
0
1
0
0
1
0
0
]1[
1
1
1
1
11
1
1
1
1
1
1
1
)(
(16)
则对数似然函数)()(LnLLL为:
n
i
m
j
i
X
j
i
X
T
eLnLLL
1
1
)
0
()1ln()()(
(17)
令
)(0
i
j
i
j
i
XXD,它是一个与第i个样本点有关的k维向量,j
ig
D表示向量中的第g个元素,
则有如下得分函数和信息矩阵:
g
g
LL
F
)(
)(=
n
i
m
j
j
i
D
T
m
j
j
i
D
T
j
ig
e
eD
1
1
1
1
hg
LL
I
gh
)(
)(
2
=khg
e
eDeD
e
eDD
m
j
j
i
D
T
m
j
j
i
D
T
j
ih
m
j
j
i
D
T
j
ig
n
i
m
j
j
i
D
T
m
j
j
i
D
T
j
ih
j
ig
,,2,1,]
)1(1
[
2
1
11
1
1
1
注意此时的T
k
),,,(
21
,没有常数项
0
。至此(17)式中的参数可采用Newton-Raphson迭
代法求解了,初值依然取为0向量。不过该方程的求解已经相对复杂多了。
方程似然度检验和回归系数的wald检验同非条件Logistic回归。
例2.3.1研究肥胖、口服避孕药雌激素与子宫内膜癌的关系,随机抽取20名患者,对于每名患者,在
随机抽取年龄相近的正常人作为对照。检测患者与正常人的肥胖程度和雌激素服用情况[3]。
表2.3.1肥胖和雌激素与子宫内膜癌关系病例-对照研究数据
配对组
病例
肥胖
病例
雌激素
对照1
肥胖
对照1
雌激素
对照2
肥胖
对照2
雌激素
1110000
2111101
3110111
4010001
5001001
6110010
7110111
8111011
9101111
10010100
11010110
12010101
13110011
14110010
15110001
16010101
17010010
18111001
19100101
20110100
例2.3.1求解的主要结果,由LLLStat软件计算得到:
表2.3.2条件Logistic回归系数检验
序号均值(病例)回归系数系数标准误wald统计量自由度df检验P值
肥胖0.6500001.8239140.54719211.11039010.000859
雌激素0.8500001.5896210.45054412.44836710.000419
表2.3.3条件Logistic回归模型基本信息
样本量20
求解方法极大似然+牛顿迭代
迭代次数45
-2LogLikelihood(Beta)33.306763
-2LogLikelihood(0)43.944492
方程Wald值(相减)10.637728
方程自由度2
方程检验p值0.004898
2.4多分类有序反应变量Logistic回归
在实际应用中,经常遇到反应变量为多分类有序变量的情况,例如评价指标分为差、中、良、优
等,各等级之间是有序的。这种资料的Logistic回归分析通常称为比例比数模型(累积概率模型)[4],它
需要拟合m-1(m为水平或等级个数)个Logistic回归模型。
有序累积概率Logistic模型:
1,,2,1;,,2,1,
1
)|(
mjni
e
e
XjyP
i
X
T
j
a
i
X
T
j
a
ii
或(18)
1,,2,1,
)|(1
)|(
ln
1
1
mjX
XkyP
XkyP
i
T
j
j
k
ii
j
k
ii
(19)
有序累积概率模型参数的极大似然估计就是寻找参数使得联合概率实现最大化,由于观测之间相互
独立,联合概率被分解成边缘概率之积。而观测到
jy
i
的概率就是累积概率之差:
)|1()|()|(
iiiiii
XjyPXjyPXjyP
第i个观测值对应似然值的贡献取决于观测到哪一个j值,因此对于次序响应的每个j值,取所有
jy
i
的观测之的乘积,有似然函数:
n
i
m
j
ij
d
ii
XjyPL
11
)|(,其中若jy
i
,则
1
ij
d
,否则0
ij
d
并且对于任一个观测
i
X而言,只有一个等级事件发生,即
m
j
ii
XjyP
1
1)|(,故有(19)式。其对数
似然函数如下(对于分组数据,似然函数变为:
n
i
m
j
ij
d
i
n
ii
XjyPL
11
)|(,
i
n分组中各分类例数)。
])
11
ln()
1
1ln(
1
ln[
)|(lnln
1
1
2
1
1
1
1
1
1
1
11
n
i
m
j
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
ij
i
X
T
m
a
i
X
T
m
a
im
i
X
T
a
i
X
T
a
i
n
i
m
j
iiij
e
e
e
e
d
e
e
d
e
e
d
XjyPdL
(20)
其中:
mj
mj
j
e
e
e
e
e
e
e
e
XjyP
i
X
T
m
a
i
X
T
m
a
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
i
X
T
a
i
X
T
a
11
1
1
1
11
1
)|(
1
1
1
1
1
1
(21)
然后就可以通过极大似然法,就上Newton-Raphson方法加以求解参数,,,
11m
aa了,注意的是
121
m
aaa。下面给出具体推导,,,
11m
aa求解的详细过程。对(20)式进行化简,可得
)}1ln()]1ln([
])1ln()1ln()ln([{ln
11
11
1
1
2
11
i
X
T
m
a
im
i
X
T
a
i
T
i
n
i
m
j
i
X
T
j
a
i
X
T
j
a
j
a
j
a
i
T
ij
edeXad
eeeeXdL
(22)
n
i
i
X
T
a
i
X
T
a
aa
a
i
i
X
T
a
i
e
e
ee
e
d
e
d
a
L
1
1
1
12
1
2
1
1
1
)]
1
(
1
1
[
ln
(23)
n
i
i
X
T
m
a
i
X
T
m
a
im
i
X
T
m
a
i
X
T
m
a
m
a
m
a
m
a
im
me
e
d
e
e
ee
e
d
a
L
1
1
1
1
1
21
1
1
1
]
1
)
1
([
ln
(24)
2,,2,)]
1
()
1
([
ln
1
1
1
1
mg
e
e
ee
e
d
e
e
ee
e
d
a
Ln
i
i
X
T
g
a
i
X
T
g
a
g
a
g
a
g
a
ig
i
X
T
g
a
i
X
T
g
a
g
a
g
a
g
a
ig
g
(25)
kg
e
e
e
e
d
e
e
d
e
dx
Ln
i
m
j
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
ij
i
X
T
m
a
i
X
T
m
a
im
i
X
T
a
iig
g
,,2,1
,])
11
1(
11
1
[
ln
1
1
2
1
1
1
1
1
1
(26)
n
i
i
X
T
a
i
X
T
a
aa
aa
i
i
X
T
a
i
X
T
a
i
e
e
ee
e
d
e
e
d
aa
L
1
2
1
1
2
12
21
2
2
1
1
1
11
2
)]
)1(
)(
(
)1(
[
ln
(27)
n
i
i
X
T
m
a
i
X
T
m
a
im
i
X
T
m
a
i
X
T
m
a
m
a
m
a
m
a
m
a
im
mme
e
d
e
e
ee
e
d
aa
L
1
2
1
1
2
1
1
2
21
21
1
11
2
]
)1(
)
)1(
)(
([
ln
(28)
2,,2
)]
)1(
)(
()
)1(
)(
([
ln
1
2
2
1
1
1
2
2
1
1
2
mg
e
e
ee
e
d
e
e
ee
e
d
aa
Ln
i
i
X
T
g
a
i
X
T
g
a
g
a
g
a
g
a
g
a
ig
i
X
T
g
a
i
X
T
g
a
g
a
g
a
g
a
g
a
ig
gg
(29)
2,,2,1,
)(
ln
1
2
1
1
1
1
2
mg
ee
e
d
aa
Ln
i
g
a
g
a
g
a
g
a
ig
gg
(30)
khmg
e
e
ddx
a
Ln
i
i
X
T
g
a
i
X
T
g
a
igigih
hg
,,2,1;1,,1,
)1(
)(
ln
1
2
1
2
(31)
khg
e
e
e
e
d
e
e
d
e
e
dxx
Ln
i
m
j
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
i
X
T
j
a
ij
i
X
T
m
a
i
X
T
m
a
im
i
X
T
a
i
X
T
a
iihig
hg
,,2,1,
)]
)1()1(
(
)1()1(
[
ln
1
1
2
2
1
1
2
2
1
1
2
1
1
1
2
(32)
由此构建信息矩阵),(aI和),(aF,并可迭代求解了。注:若为分组数据,上述每项乘以
i
n。
例2.4.1研究性别和两种治疗方法对某种疾病疗效的影响[3],将疗效分成效果显、有效和无效三个等级,
根据试验调查,得到如下资料。
表2.4.1性别和两种治疗方法对某种疾病疗效的影响
性别治疗方法显著有效无效合计
女
新药165627
传统671932
男
新药52714
传统101011
表2.4.2多分类有序反应变量数据格式
行号性别治疗方法频数疗效等级
111161
21152
31163
41061
51072
610193
70151
80122
90173
100011
110002
1200103
计算结果,由LLLStat统计软件给出:
表2.4.3回归系数方差矩阵V(beta)(信息矩阵I(Beta)的逆矩阵)
0.3747330.324880-0.257757-0.192823
0.3248800.323782-0.244457-0.169612
-0.257757-0.2444570.2894880.069404
-0.192823-0.1696120.0694040.236257
表2.4.4有序分类因变量Logistic回归系数检验
序号回归系数系数标准误wald统计量自由度df检验P值
常数项a1-2.6935760.61215519.36137710.000011
常数项a2-1.8120400.56901810.14105910.001450
性别1.0523520.5380413.82552810.050477
治疗方法2.1872720.48606320.24980010.000007
表2.4.5有序分类因变量Logistic回归模型基本信息
样本分组数12
求解方法极大似然+牛顿迭代
迭代次数17
注意:该结果与SAS,DPS不一致。
Poisson回归模型
1简介
一般情况下,单位容积水中的细菌数,单位时间内某些事件发生的次数,单位面积上降落的灰尘的
颗粒数等,都可以用Poisson分布来描述。一般Poisson分布描述成随机变量)(~PY,概率分布律为:
,2,1,0,
!
)(y
y
eyYP
y
易知EY,通常可能受到众多因素的影响,不妨假设这些因素为
k
xxx,,,
21
(自变量,协变
量),令),,,,1(
21k
xxxX,对于分组数据,Poisson分布的期望发生数假设为[7]
:
i
X
T
i
ik
x
ki
x
iiii
enenXyE
110)|((1)
其中T
k
),,,(
10
为回归参数,
i
n为第i组的总观测数。回归模型的似然函数为Poisson分
布条件下各个格子概率函数的乘积,因此Poisson分布的极大似然函数和对数似然函数具体形式分别为:
n
i
i
i
y
i
n
i
i
n
i
i
i
y
i
i
n
i
iy
e
y
epL
1
1
11
!!
n
i
n
i
iii
n
i
i
yyL
111
)!ln(lnln
代入i
X
T
ii
en,得
n
i
n
i
i
y
j
i
X
T
ii
T
iii
n
i
n
i
i
y
j
i
X
T
i
i
X
T
ii
jenXyny
jenenyL
111
111
ln])ln([
ln]ln([ln
(2)
令
n
i
i
X
T
igiigi
g
g
eXnXy
L
F
1
][
)(ln
)(
(3)
n
i
i
X
T
ihigi
hg
gh
eXXn
L
I
1
2)(ln
)(
,khg,,1,0,(4)
则可采用Newton-Raphson迭代法求解参数T
k
),,,(
10
的极大似然估计了。
对于仅有常数项的Poisson模型,其估计值为
n
i
i
n
i
i
n
y
1
1
0
ln
ˆ
,用于计算对数似然比。
2.案例分析
例1[3]Doll和Hill(1966)研究英国男性医生患冠心病与抽烟、年龄关系。由于死亡与追踪人数和追踪时间
有关,故用追踪人数和追踪时间的乘积(人年)作为观察单位数。假定其目标变量(死亡人数)近似服从Poisson
分布,其调查取样共74588调查单位,死亡598例。主要研究因素有抽烟(1为抽烟,0为不抽烟);调查对
象年龄分成4组(35-44岁,45-54岁,55-64岁,65-74岁),此为多分类变量,需要设置三个变量加以区分,可将
其中一个年龄组作为参照组,不妨取35-44岁,计算时不考虑对照组信息。
表1英国男性医生患冠心病与抽烟、年龄关系
分组抽烟34-44岁45-54岁55-64岁65-74岁死亡数总例数
111
218
312
413
5
63
7
8
由LLLStat软件计算的如下结果:
表2回归系数方差矩阵V(beta)(信息矩阵I(Beta)的逆矩阵)
0.040354-0.013325-0.028763-0.028467-0.028496
-0.0133250.016227-0.000790-0.001151-0.001115
-0.028763-0.0007900.0380710.0294680.029466
-0.028467-0.0011510.0294680.0337670.029491
-0.028496-0.0011150.0294660.0294910.034161
表3分组Poisson回归系数检验