使用statsmodels做单变量逻辑回归遇LinAlgError(Singular matrix)求解

问题场景 #

使用statsmodels库（样本量N=153），针对每组连续预测变量x与结局变量y，在校正BMI（连续）、性别、年龄协变量后构建单变量逻辑回归模型，以分析变量关联。但部分x-y组合触发错误：LinAlgError("Singular matrix")。

已排除的潜在原因：

• 预测变量高度相关：同一x搭配其他y可正常运行，排除该可能
• 变量列值单一：所有预测变量、协变量及结局变量均无此情况
• 因变量方差过低：所有y的0/1占比均不低于15%
• 完全分离：statsmodels针对完全分离会触发PerfectSeparationError，且数据集无此问题
• 更换求解器：statsmodels的Logit无法更换求解器，sklearn实现不适用当前场景

示例代码：

for y in Y:
    for x in X_df.columns:
        X = Covariates
        X[x] = X_df[x]
        X = sm.add_constant(X)
        model = sm.Logit(y, X)
        result = model.fit()
        
# Y.iloc[:, 0] = 运行正常
# Y.iloc[:, 1] = 运行正常
# Y.iloc[:, 2] = 触发LinAlgError("Singular matrix")

解决方案 #

1. 排查特定y子集下的共线性 #

同一x搭配其他y正常，但在触发错误的y对应样本子集中，x可能与某个协变量出现完全共线性。可通过方差膨胀因子（VIF）排查：

from statsmodels.stats.outliers_influence import variance_inflation_factor
import pandas as pd

# 针对触发错误的y子集提取数据
y_problem = Y.iloc[:, 2]
valid_idx = ~y_problem.isna()  # 处理缺失值（如有）
X_check = Covariates.loc[valid_idx].copy()
X_check[x] = X_df.loc[valid_idx, x]
X_check = sm.add_constant(X_check)

# 计算VIF值
vif_data = pd.DataFrame()
vif_data["特征"] = X_check.columns
vif_data["VIF"] = [variance_inflation_factor(X_check.values, i) for i in range(X_check.shape[1])]
print(vif_data)

若VIF值远大于10，说明存在严重共线性，可结合业务逻辑移除相关协变量，或对变量做中心化处理。

2. 中心化连续变量+正则化拟合 #

连续变量的量纲差异或均值偏移可能导致矩阵奇异，中心化处理可缓解数值问题；同时使用正则化拟合压缩系数：

# 中心化协变量和预测变量
Covariates_centered = Covariates.copy()
Covariates_centered['BMI'] = Covariates_centered['BMI'] - Covariates_centered['BMI'].mean()
Covariates_centered['年龄'] = Covariates_centered['年龄'] - Covariates_centered['年龄'].mean()

X_df_centered = X_df.copy()
X_df_centered = X_df_centered.apply(lambda col: col - col.mean())

# 带正则化的逻辑回归拟合
for y in Y:
    for x in X_df_centered.columns:
        X = Covariates_centered.copy()
        X[x] = X_df_centered[x]
        X = sm.add_constant(X)
        model = sm.Logit(y, X)
        # 使用L1正则化，alpha可根据情况调整
        result = model.fit_regularized(method='l1', alpha=0.1)

3. 检查样本子集的异常值/缺失值 #

针对触发错误的y，检查对应样本中是否存在极端异常值或隐性缺失：

y_problem = Y.iloc[:, 2]
combined_data = pd.concat([Covariates, X_df, y_problem], axis=1)

# 查看描述统计
print(combined_data.describe())
# 检查缺失值
print(combined_data.isna().sum())

若存在异常值，可考虑剔除或截断处理；缺失值需根据情况做填充或删除。

4. 改用GLM模型（带logit链接） #

statsmodels的GLM支持更多求解器，可替换Logit为GLM并指定求解器：

for y in Y:
    for x in X_df.columns:
        X = Covariates.copy()
        X[x] = X_df[x]
        X = sm.add_constant(X)
        # 使用二项分布+logit链接的GLM
        model = sm.GLM(y, X, family=sm.families.Binomial())
        # 指定bfgs求解器，可尝试其他如lbfgs
        result = model.fit(method='bfgs')

内容的提问来源于stack exchange，提问作者JED HK