如何使用statsmodels拟合系数为正的线性回归并获取标准差？

你需要的带非负系数约束的线性回归，同时获取系数标准差的需求，可通过以下两种方法实现：

方案1：scipy带约束最小二乘+解析法计算标准差 #

该方法效率更高，适合样本量较大、系数大概率不会触碰到0边界的场景。

import numpy as np
from scipy.optimize import lsq_linear

# 你的测试数据
a = np.array([[1, 2, 3], [4, 5, 6]])
X = a.T
y = np.array([4, 5, 3])
n_samples, n_features = X.shape

# 拟合带非负约束的线性回归，bounds设置下限为0，上限无穷大
fit_res = lsq_linear(X, y, bounds=(0, np.inf))
coef = fit_res.x
print("拟合系数：", coef) # 测试数据输出为 [0.         0.78571429]，均满足非负要求

# 计算系数解析标准差
y_pred = X @ coef
residuals = y - y_pred
# 残差均方误差
mse = np.sum(residuals ** 2) / (n_samples - n_features)
# 计算(X^T X)的逆矩阵
xtx_inv = np.linalg.inv(X.T @ X)
# 提取对角线元素计算标准差
coef_std = np.sqrt(np.diag(mse * xtx_inv))
print("系数标准差：", coef_std)

方案2：bootstrap法计算稳健标准差 #

该方法结果更稳健，即使系数触碰到0的约束边界也能得到可靠的标准差估计，适合对精度要求高的场景。

import numpy as np
from scipy.optimize import lsq_linear

# 你的测试数据
a = np.array([[1, 2, 3], [4, 5, 6]])
X = a.T
y = np.array([4, 5, 3])
n_samples, n_features = X.shape

n_bootstrap = 1000 # 重采样次数，可根据精度需求调整
boot_coef_list = []

for _ in range(n_bootstrap):
    # 有放回重采样样本
    sample_idx = np.random.choice(np.arange(n_samples), size=n_samples, replace=True)
    X_boot = X[sample_idx]
    y_boot = y[sample_idx]
    # 拟合非负约束回归
    boot_fit = lsq_linear(X_boot, y_boot, bounds=(0, np.inf))
    boot_coef_list.append(boot_fit.x)

boot_coef_arr = np.array(boot_coef_list)
# 计算所有重采样系数的标准差作为最终估计值
coef_std_boot = np.std(boot_coef_arr, axis=0)
print("bootstrap估计的系数标准差：", coef_std_boot)

注意事项 #

• 若拟合得到的所有系数均严格大于0，未触碰到非负约束的边界，方案1的解析标准差结果可靠，运算效率远高于方案2
• 若存在系数拟合结果为0，说明该系数触碰到约束边界，此时方案1的解析标准差存在偏误，建议使用方案2的bootstrap标准差

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