多重共线性检测—相关性系数矩阵和方差膨胀系数(VIF)分析学习

多重共线性检测—相关性系数矩阵和方差膨胀系数(VIF)分析学习

多重共线性是在多元线性回归分析中经常遇到的一个问题，它发生在当两个或更多的预测变量（自变量）在统计模型中高度相关。本文将详细介绍多重共线性的概念、检测方法以及处理策略，并通过TCGA数据集演示具体的检测过程。

多重共线性的概念

多重共线性（Multicollinearity）是在多元线性回归分析中经常遇到的一个问题，它发生在当两个或更多的预测变量（自变量）在统计模型中高度相关。在这种情况下，这些变量之间的关系会影响模型对各个变量影响的准确评估，从而导致以下几个问题：

1. 参数估计的不准确：由于变量之间的高度相关性，模型中的系数（参数）估计可能会变得非常不稳定。微小的数据变化可能导致估计值大幅波动。
2. 系数解释困难：当变量彼此相关时，很难区分单个变量对因变量影响的独立效应。这会使得模型的解释变得复杂和模糊。
3. 统计显著性的问题：共线性往往会增加参数的标准误差，这可能导致统计检验无法拒绝原假设，即使在实际中某些变量是有影响的。

检测多重共线性的常用方法

1. 方差膨胀因子（VIF）：这是最常用的检测多重共线性的量化方法。一般认为，VIF值大于5或者10表明存在严重的多重共线性，需要进一步处理。
2. 容忍度（Tolerance）：这是VIF的倒数，较低的容忍度值（通常小于0.1）表明高共线性。
3. 相关性系数矩阵：检查预测变量之间的相关系数。高度相关（例如，相关系数大于0.8或小于-0.8）可能指示共线性。

处理多重共线性的策略

1. 移除变量：如果某些变量之间存在高共线性，可以考虑从模型中移除一些变量。
2. 合并变量：将相关的变量合并为一个新变量，例如，通过计算几个相关变量的平均值。
3. 增加样本量：有时增加样本量可以帮助减少共线性带来的影响，因为更多的数据提供了更多的信息，有助于更准确地估计模型参数。
4. 正则化方法：使用如岭回归（Ridge Regression）或LASSO这类引入惩罚项的技术可以有效地处理共线性问题。

步骤流程

1. 导入数据

使用TCGA数据集作为示例数据。

rm(list = ls())
load(”./est_data.Rdata“)
dat_test <- est_data
head(dat_test)[1:5,1:5]
#                                ID  OS    EVA1A    TNFAIP6
# TCGA-CR-7374-01A TCGA-CR-7374-01A   0 -1.28317273 -1.6738964
# TCGA-CV-A45V-01A TCGA-CV-A45V-01A   1 -0.85167331 -0.1011465
# TCGA-CV-7102-01A TCGA-CV-7102-01A   1  1.41296836  0.4649019
# TCGA-MT-A67D-01A TCGA-MT-A67D-01A   0 -0.06444696  0.8510312
# TCGA-P3-A6T4-01A TCGA-P3-A6T4-01A   1 -0.76278103 -0.3687546
dim(est_data)
# [1] 493  40  

准备好一个含有结局变量和基因/临床参数(自变量)的数据。

2. 多重共线性检验

这里使用两种方法进行多重共线性检测，容忍度检测是VIF的倒数就不再演示了。

相关性系数矩阵

变量之间的相关性值一般认为要小于0.8，否则认为系数之间存在"严重"的共线性。如果变量较多时，不建议使用相关性系数矩阵，不然就会出现下图这样看不清的情况。

#相关性系数矩阵
colnames(dat_test)
library(PerformanceAnalytics) 
chart.Correlation(dat_test[,c(4:41)], histogram=TRUE, method="pearson")  

方差膨胀因子

一般宽松一点会把标准设定为10，严格一点会把标准设定为5。

# 方差膨胀因子
library(car)
colnames(dat_test)
# 找出除了ID,OS以外的所有变量
variables <- setdiff(names(dat_test), c("ID", "OS"))
# 把所有的变量用+连起来,并创建模型
e <- paste(variables, collapse = " + ")
full_formula <- as.formula(paste("OS ~", e)) #构建函数的时候一定要包含响应变量，其中的符号为“~”
M <- lm(full_formula, data = est_data)
vif_values <- vif(M)
#做图
library(ggplot2)
vif_data <- data.frame(Variable = names(vif_values), VIF = vif_values)
# 使用 ggplot2 绘制 VIF 值的条形图
png("VIF.png",width = 2000, height = 2800, res = 300)
ggplot(vif_data, aes(x = reorder(Variable, VIF), y = VIF, fill = VIF)) +
    geom_bar(stat = "identity") +
    theme_minimal() +
    labs(title = "VIF Values", x = "Variables", y = "Variance Inflation Factor (VIF)") +
    geom_hline(yintercept = 5, linetype = "dashed", color = "red",size = 3) +
    #coord_flip()+  # 翻转坐标轴，使得条形图水平显示
    theme(axis.text.x = element_text(angle = 45, hjust = 1),
          axis.text = element_text(size = 10),  # 调整轴标签字体大小
          axis.title = element_text(size = 16),  # 调整轴标题字体大小
          plot.title = element_text(size = 20))  # 调整图标题字体大小
dev.off()  

参考资料

1. 医学统计学，主编：孙振球/徐勇勇
2. 医学和生信笔记： https://mp.weixin.qq.com/s/yk8Y51ilnnMwdaVcK457KQ
3. 宇哥读文献： https://mp.weixin.qq.com/s/-jyP-QicDJvqTlLfx9008w
4. 观科研： https://mp.weixin.qq.com/s/icQ8Amfhvx6rwQFhD44Mnw

注：若对内容有疑惑或者有发现明确错误的朋友，请联系后台(欢迎交流)。更多内容可关注公众号：生信方舟

本文原文来自腾讯云开发者社区