解决计算Mahalanobis距离时Lapack routine dgesv奇异矩阵错误

先处理原问题：修复奇异矩阵，让Mahalanobis距离可计算 #

1. 移除重复/完全线性相关的变量 #

你的数据里存在重复列（示例中的两个Q3），这是协方差矩阵奇异的典型诱因。先清理这类变量：

# 检查并移除重复列
duplicated_cols <- duplicated(t(all_data_EFA))
all_data_clean <- all_data_EFA[, !duplicated_cols]

# 或者检测完全相关的变量（相关系数绝对值=1）
cor_matrix <- cor(all_data_EFA)
high_cor_pairs <- which(abs(cor_matrix) == 1 & row(cor_matrix) < col(cor_matrix), arr.ind = TRUE)
# 保留每对相关变量中的第一个，移除其余
cols_to_remove <- unique(high_cor_pairs[,2])
all_data_clean <- all_data_EFA[, -cols_to_remove]

清理后重新计算Mahalanobis距离，通常能解决奇异问题。

2. 使用正则化协方差矩阵 #

如果数据存在共线性但不能删除变量，可以给协方差矩阵添加小对角扰动使其可逆，推荐用corpcor包的收缩协方差估计：

library(corpcor)
# 计算收缩后的协方差矩阵
cov_shrunk <- cov.shrink(all_data_EFA)
# 基于正则化矩阵计算马氏距离
distance <- mahalanobis(all_data_EFA, colMeans(all_data_EFA), cov = cov_shrunk)

3. 降维后计算马氏距离 #

通过PCA将原变量转换为不相关的主成分，主成分的协方差矩阵为对角矩阵，不会出现奇异问题：

# 执行PCA，保留解释95%方差的主成分
pca_result <- prcomp(all_data_EFA, scale. = TRUE)
cum_var <- cumsum(pca_result$sdev^2 / sum(pca_result$sdev^2))
n_components <- which(cum_var >= 0.95)[1]
pca_scores <- pca_result$x[, 1:n_components]

# 在主成分空间计算马氏距离
distance_pca <- mahalanobis(pca_scores, colMeans(pca_scores), cov = cov(pca_scores))

替代的多元异常值检测方法 #

如果不想纠结马氏距离，这些方法更适配高维问卷数据：

1. 稳健马氏距离 #

用稳健估计的均值和协方差矩阵（不受异常值干扰，且不易奇异），依赖robustbase包：

library(robustbase)
# 计算稳健均值与协方差
robust_cov <- covMcd(all_data_EFA)
# 计算稳健马氏距离
distance_robust <- mahalanobis(all_data_EFA, robust_cov$center, robust_cov$cov)
# 筛选p值<=0.001的异常值
mah_p_robust <- pchisq(distance_robust, df = ncol(all_data_EFA), lower.tail = FALSE)
outliers_robust <- which(mah_p_robust <= 0.001)

2. 孤立森林 #

针对高维数据设计的异常检测算法，无需计算协方差矩阵，用isotree包：

library(isotree)
# 训练孤立森林模型
iso_model <- isolation.forest(all_data_EFA, ntrees = 100)
# 获取异常得分（得分越高越可能是异常值）
anomaly_scores <- predict(iso_model, all_data_EFA, type = "score")
# 筛选前0.1%的高得分样本（对应原阈值0.001）
outliers_iso <- which(anomaly_scores >= quantile(anomaly_scores, 0.999))

3. 局部离群因子（LOF） #

基于局部密度差异识别异常值，适合有聚类结构的数据，依赖dbscan包：

library(dbscan)
# 计算LOF值（k取5个邻居，可按需调整）
lof_scores <- lof(all_data_EFA, k = 5)
# LOF值大于1表示样本密度低于邻居，值越大异常程度越高
outliers_lof <- which(lof_scores >= 2)

内容的提问来源于stack exchange，提问作者Cristina Procentese