给定y的测量误差，如何计算线性回归系数的误差及官方实现方法

我手头有一组声速测量的实验室数据，x和y的对应值如下：

x       y
0       0
1     212
2     426
3     640
4     858
5    1074
6    1290
7    1506
8    1722
9    1939

已知每个y的测量值存在±2的误差，比如x=1时，y的实际取值范围是210到214。我想搞清楚这个误差会对线性回归的系数（斜率和截距）产生多大影响。

之前用sklearn的LinearRegression时，设置fit_intercept=False的情况下，只要分别计算y-2和y+2序列的回归系数，再求差值就能得到误差范围，但现在需要保留截距项（也就是允许x=0时y不为0），想知道有没有官方实现的方法（不限于sklearn）。

可行的官方实现方案 #

• 用statsmodels计算参数置信区间
  statsmodels是专门做统计建模的库，它的线性回归模型会直接给出参数的置信区间，其中就包含了测量误差带来的影响。示例代码：

  import statsmodels.api as sm
  import numpy as np
  
  # 准备数据
  x = np.array([0,1,2,3,4,5,6,7,8,9]).reshape(-1,1)
  y = np.array([0,212,426,640,858,1074,1290,1506,1722,1939])
  # 添加常数项（对应截距）
  x_with_intercept = sm.add_constant(x)
  
  # 拟合模型
  model = sm.OLS(y, x_with_intercept).fit()
  # 查看参数的95%置信区间，这里的区间已经考虑了测量误差的影响
  print(model.conf_int())
  

  输出的置信区间就是截距和斜率在统计意义下的误差范围，完全符合需求。

• 用scipy的linregress计算统计指标
  scipy的linregress函数会返回斜率、截距以及它们的标准误差，你可以用标准误差乘以临界值（比如95%置信度用1.96）得到误差范围：

  from scipy.stats import linregress
  import numpy as np
  
  x = np.array([0,1,2,3,4,5,6,7,8,9])
  y = np.array([0,212,426,640,858,1074,1290,1506,1722,1939])
  
  result = linregress(x, y)
  # 计算95%置信度下的斜率和截距误差范围
  slope_error = result.stderr * 1.96
  intercept_error = result.intercept_stderr * 1.96
  print(f"斜率误差范围: ±{slope_error:.2f}")
  print(f"截距误差范围: ±{intercept_error:.2f}")
  

• Bootstrap抽样法（sklearn也支持）
  如果你想更直观地模拟误差的影响，可以用bootstrap方法：多次随机扰动y值（在±2范围内），每次拟合回归模型，最后统计系数的分布。sklearn的resample函数可以实现：

  from sklearn.linear_model import LinearRegression
  from sklearn.utils import resample
  import numpy as np
  
  x = np.array([0,1,2,3,4,5,6,7,8,9]).reshape(-1,1)
  y = np.array([0,212,426,640,858,1074,1290,1506,1722,1939])
  
  slopes = []
  intercepts = []
  n_iterations = 1000
  
  for _ in range(n_iterations):
      # 给每个y添加±2范围内的随机误差
      y_perturbed = y + np.random.uniform(-2, 2, size=y.shape)
      model = LinearRegression(fit_intercept=True)
      model.fit(x, y_perturbed)
      slopes.append(model.coef_[0])
      intercepts.append(model.intercept_)
  
  # 计算95%置信区间（取2.5%和97.5%分位数）
  slope_ci = np.percentile(slopes, [2.5, 97.5])
  intercept_ci = np.percentile(intercepts, [2.5, 97.5])
  print(f"斜率95%置信区间: [{slope_ci[0]:.2f}, {slope_ci[1]:.2f}]")
  print(f"截距95%置信区间: [{intercept_ci[0]:.2f}, {intercept_ci[1]:.2f}]")
  

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