R与医学统计学·系列 | 之logistic回归

除了变量之间的线性相关关系，临床研究设计中也会常用logistic回归来进行两个或者多个因素直接的相关性分析，常用于进行疾病的病因分析。

Logistic回归按照反应变量的类型可分为：

二分类反应变量的logistic回归、多分类有序反应变量的logistic回归、多分类无序反应变量的logistic回归；

按照研究设计的类型可分为：

研究对象未经过匹配的非条件logistic回归、研究对象经匹配的条件logistic回归

今天主要介绍最为常见的二分类反应变量的logistic回归。

【基本概念】

1.基本原理

Logistic回归起源于生物统计学，是为了解决线性回归在处理二分类问题时的局限性而诞生的，主要通过Logit函数将事件发生的概率映射到线性组合上，使得我们可以使用线性模型来处理分类问题。logistic函数可以用下图中的式子表示（左、右分别为一元和多元logistic回归）。X为自变量，β为自变量的系数。可以看出，这个式子的取值范围一定在0～1之间。

本质上，logistic函数计算的是概率，我们要做的是找到一个阈值（如下图），如果计算概率大于这个阈值，则预测为阳性（Y=1），否则为阴性（Y=0）。

2.β系数（Beta Coefficient）

β系数是Logistic回归模型中的关键参数，它量化了自变量每变化一个单位时对结果变量对数几率的影响。在模型中，β系数的符号（正或负）指示了自变量与结果变量之间的关系方向：正相关或负相关。

正β系数：表示自变量的增加与结果发生概率的增加相关。

负β系数：表示自变量的增加与结果发生概率的减少相关。

3.标准误（Standard Error）

标准误是β系数估计的不确定性度量，反映了估计值的可靠性。较小的标准错意味着β系数的估计更为精确，而较大的标准错则意味着估计值可能存在较大的变异。

4.Wald检验与P值

Wald检验是一种统计检验方法，用于评估模型中β系数是否显著不同于0。P值是该检验的结果，提供了β系数显著性的概率水平。在大多数情况下，P值小于0.05被认为是统计学上显著的，意味着可以拒绝零假设，认为自变量对结果有显著影响。

5.Odds Ratio（OR值）

Odds Ratio（OR值）是Logistic回归分析中的一个重要指标，它是β系数的指数化（e^β），表示自变量每变化一个单位时，结果发生几率相对于不发生几率的比值。OR值大于1表明自变量增加结果发生的概率，而OR值小于1则表明自变量减少结果发生的概率。

6.95%置信区间（Confidence Interval, CI）

95%置信区间是OR值估计的不确定性范围，提供了对OR值估计精度的一种度量。如果95%置信区间相反，如果置信区不包含1，这通常意味着自变量对结果的影响是统计学上显著的。间包含1，则没有足够的证据表明自变量对结果有显著影响。

【Logistic回归模型的应用条件】

1. 因变量是二分类变量：Logistic回归专门用于处理二分类问题，如判断病人是否患病。

2. 自变量可以是连续或分类变量：自变量可以是多种类型的数据，不局限于某一种。

3. 假设线性关系：Logistic回归假设自变量与因变量的Logit变换之间存在线性关系。

4.独立性：观测值之间应该相互独立。

【模型检验】

常用的假设检验包括wald检验和似然比检验。评估指标: 拟合优度指标、信息测量指标和广义R2指标。假设检验: 检查多元线性回归假设是否满足（线性关系、独立性、正态性、同方差性）。

【R实操】

示例数据使用MASS包中的biopsy数据集。

1. 数据读取与模型创建

为了确保数据质量，我们需要处理缺失值并转换目标变量。

Null deviance 表示仅包含截距项的模型的偏差，Residual deviance 表示包含所有自变量的模型的偏差。偏差越低，模型越好。AIC（赤池信息准则）值越低，模型越好。

2. 模型预测与评估

使用测试集进行预测，并评估模型的性能（通过混淆矩阵，我们可以评估模型的性能。）。

01.Accuracy（准确率）：0.9755，准确率表示正确预测的比例。这里模型的准确率为97.55%，说明模型整体预测效果很好。

02.95% 置信区间为 (0.9437, 0.992)，表示我们有95%的信心准确率落在这个区间内。

03.No Information Rate为0.652，表示如果我们总是预测多数类（良性），得到的准确率。P-Value [Acc > NIR] < 2e-16 表明模型的准确率显著高于No Information Rate。

04.Kappa（卡帕系数）：0.9458，卡帕系数考虑了偶然一致的影响，范围从 -1 到 1。0.9458 表示模型具有非常高的一致性。

05.Mcnemar's Test P-Value：1，该值用于比较分类器是否在不同类别上有不同的错误率。P值为1表示分类器在不同类别上没有显著不同的错误率。

06.Sensitivity（敏感性/召回率）：0.9850，敏感性表示真实为正类（良性）且被正确预测为正类的比例。0.9850 表示98.50%的良性样本被正确预测为良性。

07.Specificity（特异性）：0.9577，特异性表示真实为负类（恶性）且被正确预测为负类的比例。0.9577 表示95.77%的恶性样本被正确预测为恶性。

08.Pos Pred Value（阳性预测值）：0.9776，阳性预测值表示被预测为正类（良性）中真实为正类的比例。0.9776 表示97.76%的预测为良性样本实际是良性。

09.Neg Pred Value（阴性预测值）：0.9714，阴性预测值表示被预测为负类（恶性）中真实为负类的比例。0.9714 表示97.14%的预测为恶性样本实际是恶性。

10.Prevalence（患病率）：0.6520，患病率表示数据集中正类（良性）的比例。0.6520 表示65.20%的样本为良性。

11.Detection Rate（检测率）：0.6422，检测率表示正确检测到的正类样本的比例。0.6422 表示64.22%的样本被正确预测为良性。

12.Detection Prevalence（检测流行率）：0.6569，检测流行率表示模型预测为正类的比例。0.6569 表示65.69%的样本被预测为良性。

13.Balanced Accuracy（平衡准确率）：0.9714，平衡准确率是敏感性和特异性的平均值，0.9714 表示模型在正负类上都表现良好。

通过这些指标，我们可以看出模型在分类良性和恶性方面表现非常好。高准确率、高敏感性和高特异性表明模型在检测和区分良性和恶性样本方面非常有效。卡帕系数接近1，也表明模型具有高度一致性。

3. Logistic模型可视化结果

3.1. 绘制预测概率分布

首先，我们可以绘制模型预测的概率分布，以便观察预测结果的整体分布情况。

这个图显示了Logistic回归模型预测的概率分布，其中每个柱状条代表预测的概率值，颜色区分了实际类别（良性和恶性）。

从这个图中，我们可以得出以下结论：

1.模型在预测极端概率值（接近0或1）时，准确性较高。

2.实际类别为良性的样本大多集中在预测概率低的区域（接近0）。

3.实际类别为恶性的样本大多集中在预测概率高的区域（接近1）。

4.在中间概率区域，预测的准确性较低，样本数量较少。

因此，该模型在区分良性和恶性样本方面表现良好，但在概率预测接近中间值时存在一定的不确定性。

3.2. 绘制ROC曲线

ROC曲线是评估分类模型性能的重要工具，AUC值可以量化模型的区分能力。

3.3 绘制混淆矩阵热图

混淆矩阵热图可以直观地显示分类模型的性能。

这张图展示了模型的混淆矩阵热图，直观地显示了模型在预测良性（Benign）和恶性（Malignant）样本时的表现。通过颜色和数值，我们可以清晰地看到模型的预测结果与实际情况的对比。

从这张混淆矩阵热图中可以看出：

1.模型在预测良性和恶性样本方面表现非常好，大部分样本被正确分类。

2.少数样本被错误分类，其中良性样本被错误预测为恶性（3个），恶性样本被错误预测为良性（2个）。

3.4 绘制列线图（Nomogram）

列线图（Nomogram）是一种直观的工具，用于展示Logistic回归模型的预测结果，特别适用于医学领域的风险预测和决策支持。我们可以使用rms包来绘制列线图。

列线图的每个变量都有一个独立的刻度，表示其对结果的影响。通过将变量的实际值在刻度上找到对应的点，然后将所有变量的点相加，我们可以得到一个总分，再通过总分找到对应的预测概率。

我们使用上面这个的列线图来解释如何使用列线图进行预测。列线图是一种直观工具，可以帮助我们通过简单的加分法计算个体的总分，然后将总分转换为预测风险。

3.4.1选择变量值

首先，确定每个变量的取值。在本例中，我们假设某个个体的变量值如下：

· V1 = 5

· V2 = 6

· V3 = 4

· V4 = 3

· V5 = 8

· V6 = 7

· V7 = 2

· V8 = 9

· V9 = 6

3.4.2. 查找每个变量的分值

根据列线图中的各变量刻度，查找每个变量取值对应的分值：

· V1 = 5：在V1刻度上找到5对应的分值（大约为42分）。

· V2 = 6：在V2刻度上找到6对应的分值（大约为9分）。

· V3 = 4：在V3刻度上找到4对应的分值（大约为21分）。

· V4 = 3：在V4刻度上找到3对应的分值（大约为11分）。

· V5 = 8：在V5刻度上找到8对应的分值（大约为2分）。

· V6 = 7：在V6刻度上找到7对应的分值（大约为44分）。

· V7 = 2：在V7刻度上找到2对应的分值（大约为11分）。

· V8 = 9：在V8刻度上找到9对应的分值（大约为64分）。

· V9 = 6：在V9刻度上找到6对应的分值（大约为40分）。

3.4.2. 计算总分

将所有变量的分值相加，得到总分：

42+9+21+11+2+44+11+64+40 = 244 分

3.4.3. 查找线性预测值

在Total Points刻度线找到总分244对应的位置，然后在下方的Linear Predictor刻度线上找到对应的线性预测值（大约为4）。

3.4.4. 查找风险概率

根据线性预测值在Risk刻度线上找到对应的风险概率。这里，线性预测值为4，对应的风险概率大约为0.98（98%）。

通过上述步骤，我们可以得出以下结论：

1.该个体的总分为244。

2.根据总分，对应的线性预测值为4。

3.对应的风险概率为0.98（98%）。

更多干货分享↓↓↓

临床常用的十大计数公式，用这个工具免费查

【手把手教学】如何用STRING和Cytoscape绘制蛋白质互作网络图

常规转录组入门分析——数据清洁：counts、FPKM、RPKM、TPM、CPM傻傻分不清？