Enter your search text in the box above

Select a result to preview

VIF方差膨胀因子

#garden/🌱Seedling

VIF（方差膨胀因子）

核心问题：多重共线性

当多个自变量之间高度相关时，线性模型会出现多重共线性问题。

直觉理解：假设用 BIO1（年均温）和 BIO10（最暖季均温）同时预测等位基因频率，两者相关性可能高达 0.98。模型无法分辨"到底是哪个变量在起作用"，结果是两个变量的回归系数都变得不稳定、标准误极大。

数学定义

对第 j 个变量：

VIFj=11−Rj2

其中 Rj2 是将第 j 个变量作为因变量，其余所有变量作为自变量回归得到的决定系数。

Rj2 含义 VIF值
0 与其他变量完全不相关 1
0.9 90% 方差能被其他变量解释 10
0.99 高度冗余 100

阈值选择

阈值 适用场景
VIF < 3 严格，生态学/景观基因组常用
VIF < 5 一般推荐
VIF < 10 宽松，社会科学常用

vifstep 逐步剔除逻辑（R: usdm包）

library(usdm)
vif_result <- vifstep(as.data.frame(env_scaled), th = 10)

逐步迭代流程：

  1. 计算所有变量的 VIF
  2. 找到 VIF 最大的变量
  3. 若 > 阈值，剔除它
  4. 用剩余变量重新计算 VIF
  5. 重复，直到所有变量 VIF ≤ 阈值

为什么逐步而非一次性剔除？
剔除一个变量后，其他变量的 VIF 会改变。A 和 B 高度相关，剔除 A 后 B 的 VIF 可能从 12 降到 3，就无需再剔除。

在景观基因组（GEA）中的应用

因变量 vs 自变量

VIF 筛的是气候变量之间的共线性，与基因组数据无关。

哪些方法需要 VIF 筛选？

方法 是否需要VIF 原因
RDA ✅ 必须 多变量，所有环境变量同时进入模型
LFMM2 建议 单变量，但减少冗余检验
BayPass 建议 单变量，但减少冗余检验
GEMMA/单变量GWAS ❌ 不需要 每次一个变量，不存在共线性

Sun et al. (iMeta 2026) 用 GEMMA+LFMM+BayPass 单变量方法，直接保留了全部 28 个变量，并非不做 VIF，而是单变量方法本身不需要。

WorldClim BIO 变量的典型共线性结构

温度类（BIO1,2,4,5,6,7,8,9,10,11）→ 内部高度相关
降水类（BIO12,13,14,15,16,17,18,19）→ 内部相关
温度 vs 降水 → 相对独立

VIF 筛选后通常覆盖：年均温、温度季节性、极端温度、年均降水、降水季节性 ← 对应牛的热耐受、季节适应、极端气候、水分代谢。

重要注意

VIF 筛选的是自变量之间的共线性，与因变量无关。筛掉的变量不代表对等位基因频率没影响——只是其信息已被保留变量涵盖，保留反而干扰模型估计。

统计上无法证明因果方向，气候→基因组的方向是基于生物学先验知识，需后续功能验证（eQTL、ATAC-seq等）支持因果推断。

实操代码（本项目：多源环境变量 GEA）

完整流程

library(usdm)

# 读入环境变量总表
env <- read.csv("breeds_master.csv", row.names=1)

# 排除非环境列（品种信息和坐标不参与VIF计算）
non_env <- c("Num","Species","Group","PaperGroup","NewData",
             "Location","Longitude_E","Latitude_N")

# 排除Outgroup和n<5的品种（VIF在和GEA相同的品种集合上计算）
exclude_breeds <- rownames(env)[env$PaperGroup == "Outgroup" | env$Num < 5]
env_vars <- env[!rownames(env) %in% exclude_breeds,
                !colnames(env) %in% non_env]
cat("品种数：", nrow(env_vars), "，变量数：", ncol(env_vars), "\n")

# VIF逐步过滤（阈值=10）
# ⚠️ 运行时间约1-2分钟（130个变量）
# ⚠️ 警告 "essentially perfect fit" 属正常，不影响结果
vif_result <- vifstep(env_vars, th=10)
cat("保留变量数：", nrow(vif_result@results), "\n")
print(vif_result@results$Variables)

# 提取过滤后的变量
env_filtered <- env_vars[, vif_result@results$Variables]

# 对因NA的沿海品种进行均值填补
# （SOILW、PEVPR坐标落在NCEP海洋格点，无法提取）
for(col in c("PEVPR_dist","SOILW_max","SOILW_dist")){
  if(col %in% colnames(env_filtered)){
    env_filtered[is.na(env_filtered[,col]), col] <- mean(env_filtered[,col], na.rm=TRUE)
  }
}

# 加回坐标（供pRDA或地图可视化使用）
env_for_analysis <- cbind(
  env[!rownames(env) %in% exclude_breeds, c("Longitude_E","Latitude_N")],
  env_filtered
)

本项目实际结果

参数 数值
输入品种数 180（排除后127）
输入变量数 130（WorldClim BIO×19 + NCEP×64 + CRU×24 + CHELSA×3 + GPCC×4 + PDSI×4 + elevation）
VIF过滤后 32个
运行时间 约1分钟

保留的32个变量：
wc2.1_2.5m_bio_2, wc2.1_2.5m_bio_8, wc2.1_2.5m_bio_9, wc2.1_2.5m_bio_14, wc2.1_2.5m_bio_18, wc2.1_2.5m_bio_19, OLR_dist, PEVPR_dist, PRATE_mean, RHUM_dist, SNOW_mean, SOILW_max, SOILW_dist, SRFP_dist, TCDC_mean, TCDC_dist, UWND_mean, UWND_dist, VWND_mean, VWND_dist, scd, vpd_dist, PRE_dist, PDSI_mean, PDSI_dist, POTTMP_max, WSPD_dist, WET_dist, CLD_max, CLD_dist, FRS_dist, DTR_dist

注意事项

相关笔记