R语言手动实现最小二乘法出现参数估计异常的问题

嗨Marco，你的这个问题很好地展现了数值优化和解析解之间的差异，以及参数敏感性对数值结果的影响。我来帮你拆解问题根源和解决办法：

1. 截距a偏差的核心原因 #

你已经通过3D图发现SSE对a的变化不敏感，但更关键的点是：最小二乘的截距a和斜率b是高度相关的。

当你用数值方法遍历a和b的组合时，找到的是当前离散网格下SSE最小的点，但这个点不一定是全局最优的解析解——因为SSE的“谷底”在a方向上非常平缓，哪怕b的取值有一点点偏差（比如你的possible.b.vals是步长0.01的序列，并非精准的解析解b值），对应的最优a值就会偏离真实值。

举个简单例子：假设真实的b是27.84，如果你取的b是27.839，为了让SSE最小，对应的a就会比真实值大一点——斜率稍小的话，截距需要往上调才能更好拟合数据，反之亦然。当你扩大取值范围时，b的候选值离真实值更远的概率变大，对应的a偏差自然也就更明显了。

2. 数值方法的优化方向 #

（1）优先提升斜率b的分辨率 #

因为SSE对b的敏感性远高于a，你可以先精准锁定最优b，再固定b去寻找最优a：

# 先缩小b的步长，精准定位最优b
possible.b.vals <- seq(27.83, 27.85, by=0.0001)
# 固定a的合理范围，计算每个b对应的最小SSE对应的a

（2）用向量化操作替代循环提升效率 #

你的for循环在参数组合量大时效率极低，用apply()或向量化运算可以大幅提速：

# 用apply替代循环计算SSE
possible.ab.SSE$SSE <- apply(possible.ab, 1, function(x) {
  sum((longley$GNP - (x[1] + x[2] * longley$Employed))^2)
})

（3）使用专业数值优化函数 #

其实不需要手动遍历所有组合，R内置的optim()函数可以高效寻找SSE最小值：

# 定义SSE计算函数
sse_fun <- function(params) {
  a <- params[1]
  b <- params[2]
  sum((longley$GNP - (a + b * longley$Employed))^2)
}

# 传入初始猜测值
init_guess <- c(-1430, 27.8)
# 调用optim寻找最小值
result <- optim(init_guess, sse_fun)
# 查看最优参数
result$par
# 结果会非常接近lm()的解析解

optim()采用了高效的优化算法（如Nelder-Mead），比暴力遍历快得多，精度也更高。

3. 为什么解析解能精准估计？ #

解析解是通过对SSE求偏导并令其为0，直接解出a和b的最优值，不存在离散网格的精度限制，也不需要遍历所有可能组合，所以能得到精准结果。而数值方法是在离散参数空间中寻找近似最优解，必然会受到网格分辨率、参数相关性的影响。

总结一下：你的问题本质是离散数值搜索的精度限制，加上参数相关性放大了截距的偏差。通过提升关键参数（b）的分辨率、优化计算效率或者使用专业优化函数，就能让数值方法的结果更接近解析解啦。

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