条件Logistic回归分析(Conditional Logistic Regression Analysis)——R软件实现

关键词：R语言; R软件; 条件logistic回归; 配对logistic回归; 条件逻辑回归; 配对逻辑回归

一、案例介绍

某肾内科医师拟探究急性肾损伤的危险因素，回顾性收集了109例在院内发生急性肾损伤患者的性别(gender)、年龄(age)、体质指数(BMI)、血肌酐(Cr)和血清乳酸(Serum Lactate)，并根据性别和年龄进行1:1配对，收集了109例未发生肾损伤患者的相关信息，进行配对病例对照研究。“sAKI”为二分类变量，代表急性肾损伤发生情况(“发生”记录为“1”，“未发生”记录为“0”)；“BMI”为有序多分类变量，代表体质指数(“偏瘦”为“1”，“中等”为“2”，“超重”为“3”)；“Cr”和“Lactate”都为连续变量，分别表示血肌酐和血清乳酸含量；“subclass”表示对子数。部分数据见图1。本文案例可从“附件下载”处下载。

二、问题分析

本案例的分析目的是探讨某肾内科医师拟探究急性肾损伤的危险因素，采用配对设计，研究多个因素对二分类因变量的影响，可以采用条件logistic回归分析。但需要满足7个条件：

条件1：因变量为二分类变量。本研究中因变量是是否发生急性肾损伤“是”和“否”，为二分类变量，该条件满足。

条件2：至少有1个自变量。自变量可以是分类变量也可以是连续变量。本研究中有多个自变量，类型各异，该条件满足。

条件3：观察变量为配对设计或具有相关性，即不满足独立性。本研究中，两组患者是根据性别和年龄进行1:1配对，因此属于相关样本，不满足独立性。该条件满足。

条件4：因变量对子数为自变量个数的10~15倍(EPV原则)，最好＞30对，自变量的参照水平组不应少于30或50例。该条件需要通过软件分析后判断。

条件5：自变量之间无多重共线性。该条件需要通过软件分析后判断。

条件6：自变量不存在显著的异常值。该条件需要通过软件分析后判断。

条件7：数据未出现完全分离或拟完全分离现象。该条件需要通过软件分析后判断。

三、软件操作及结果解读

(一) 导入数据

mydata <- read.csv("条件logistic回归.csv")  #导入CSV数据
View(mydata)  #查看数据

在数据栏目中可以查看全部数据情况，数据集中共有7个变量和218个观察数据，7个变量分别为急性肾损伤(sAKI)、性别(Gender)、年龄(Age)、体质指数(BMI)、血肌酐(Cr)、血清乳酸(Lactate)及对子数(subclass)。

如果数据集较大也可使用如下命令查看数据框结构：

str(mydata)  #查看数据框结构

(二) 适用条件判断

1. 条件4判断(因变量样本例数)

(1) 软件操作

##预后变量计数##

table(mydata$sAKI) #计算样本数

(2) 结果解读

由图3可见，“1”代表“发生急性肾损伤”，“0”代表“未发生急性肾损伤”，两组均为109例，即109对。根据“因变量对子数为自变量个数的10~15倍(EPV原则) ”，本案例可纳入7~11个自变量进行多因素条件logistic回归分析。

2. 条件4判断(自变量样本例数)

为了增加结果的可读性，给部分分类变量增加标签，之后再逐一计算分类变量各类别的因变量例数。

(1) 软件操作

##给变量添加标签label ##
mydata$Gender <- factor(mydata$Gender,levels = c(1,0),labels = c("男", "女"))
mydata$BMI <- factor(mydata$BMI,levels = c(3,2,1),labels = c("超重", "中等","偏瘦"))

##变量计数##
table(mydata$sAKI,mydata$Gender)
table(mydata$sAKI,mydata$BMI)

(2) 结果解读

由图4可知，“BMI”水平为“1”偏瘦时，因变量的例数<30，如果该变量在多因素分析过程中进入模型，应注意避免例数较少的水平被选为参照。

3. 条件5判断(多重共线性诊断)

详见本章后文。

4. 条件6判断(异常值检测)

(1) 软件操作

##异常值检验-cook距离##
fit<-glm(sAKI~BMI+Cr+Lactate,data=mydata,family = "binomial") #建立二项logistic回归方程
cook<-cooks.distance(fit) #计算cook距离
cook[cooks.distance(fit)>0.5] #显示cook距离>0.5的个案编号和cook值
max(cook) #显示最大cook距离

(2) 结果解读

使用库克距离来判断强影响点是否为因变量的异常值点。一般认为当D<0.5时不是异常值点，当D>0.5时认为是异常值点。图5结果显示，本案例最大库克距离为0.354<0.5，提示不存在显著异常值，本研究数据满足条件6。

5. 条件7判断(完全分离检测)

完全分离，指某一个自变量本身或者某几个自变量的线性组合，对因变量的预测结果与实际情况完全一致，常表现为OR值无穷大。通过图4可见并不存在这种情况。该条件满足。

(三) 条件logistic回归分析

1. 软件操作

##建立包含所有自变量的条件logistic回归模型##
library(survival) #加载包“survival”
条件logistic分析：clogit(因变量~自变量+strata(配对变量),data,方法)
clogit_full<-clogit(sAKI~BMI+Cr+Lactate+strata(subclass),data=mydata,method="exact")
summary(clogit_full) #查看模型参数

2. 结果解读

图6显示了条件logistic回归分析的结果，可知，模型整体Likelihood=66.30，P<0.001，提示模型有统计学意义，说明模型纳入的自变量中至少有一个有统计学意义。

BMI偏瘦的OR=0.201 (95%CI：0.060-0.680, P=0.010)，说明偏瘦发生急性肾损伤的风险是超重的0.201倍(将0.201取倒数，即为超重组发生急性肾损伤的风险是偏瘦组的4.98倍)；BMI中等的OR=0.690 (95%CI：0.290-1.641, P=0.401>0.05)，说明BMI中等相比BMI超重的患者，急性肾损伤风险未发生明显变化。Cr组和Lacate组的OR值分别为1.011 (95%CI：1.006~1.017, P<0.001)和1.169 (95%CI：1.029~1.329, P=0.017)，均大于1，说明急性肾损伤风险随着Cr和Lacate的增高而增加，Cr每升高一个单位，发生急性肾损伤的风险增加0.011倍；Lacate每升高一个单位，发生急性肾损伤的风险增加0.169倍。

(四) 适用条件判断(补充条件5判断:多重共线性诊断)

1. 软件操作

##共线性诊断##

library(car)
vif(clogit_full) #计算vif

2. 结果解读 

图7中计算了自变量的方差膨胀因子(variance inflation factor，VIF)。可见，所有自变量的VIF均<10，提示自变量之间不存在严重共线性问题。满足条件5。

四、结论

本研究采用条件logistic回归模型探讨发生急性肾损伤的危险因素。通过分析可知因变量对子数和自变量个数满足需求，变量之间不存在严重共线性和异常值，数据不存在完全分离现象。满足条件logistic回归分析条件。

条件logistic回归分析结果提示，所建立的模型有统计学意义(χ²=66.30，P<0.001)。三个自变量经过分析后发现均有统计学意义，其中BMI偏瘦发生急性肾损伤的风险是超重的0.201倍 (95%CI：0.060~0.680, P=0.010)；急性肾损伤风险随着Cr的升高而增加，Cr每升高一个单位发生风险增加0.011倍 (95%CI：0.006~0.017；P<0.001)；急性肾损伤风险随着Lacate的升高而增加，Lacate每升高一个单位发生风险增加0.169倍 (95%CI：0.029~0.329, P=0.017)。