import numpy as np from scipy.optimize import linprog # 生成数据 np.random.seed(2021) p = 3 n = 50 X = np.random.normal(size=(n, p)) e = np.random.chisquare(5, size=n) y = X.dot(np.array([1, 2, 3])) + e # 定义线性规划问题 c = np.ones(n) A_ub = np.vstack([np.hstack([X, -np.ones((n, 1))]), np.hstack([-X, -np.ones((n, 1))])]) b_ub = np.hstack([y, -y]) bounds = [(None, None)]*p + [(0, None)] res = linprog(c, A_ub=A_ub, b_ub=b_ub, bounds=bounds, method='simplex') # 计算结果 b_median = res.x[:p] print('中位数回归系数为：', b_median)

这段代码实现了中位数回归，其思路是将线性回归问题转化为线性规划问题，并通过线性规划求解得到中位数回归系数。具体来说，该问题的目标函数为 $\sum_{i=1}^n w_i$, 其中 $w_i$ 是一些权重系数，代表第 $i$ 个样本对目标函数的贡献。同时，为了约束模型，需要加入一些约束条件，比如回归系数的范围等。在这里，我们将回归系数的范围限制在 $(\text{None}, \text{None})$ 和 $(0, \text{None})$ 之间，这个约束在 bounds 中给出。

通过求解该线性规划问题，我们得到的中位数回归系数为 b_median。