在R中模拟OLS回归并存储残差的方法是否正确？

残差计算逻辑的正确性 #

你写的残差计算公式residuals <- y - (beta[1,i] + beta[2,i] * x)本身是正确的，OLS残差的定义就是观测值减去回归拟合值，这里的拟合值计算完全符合定义。

存在的遗漏问题 #

当前代码最大的问题是残差没有被完整存储：

• 每次循环中你都用residuals变量覆盖前一次迭代的残差，最终只会保留第1000次迭代的残差，完全丢失了前999次的结果。
• 如果你需要存储1000次迭代的所有残差（对应n×1000的结构），需要提前初始化一个残差存储容器，比如矩阵或者列表。

修正后的代码示例 #

Data <- read.csv("cons.csv")

ols <- function(y,x){
  ones <- rep(1,length(x))
  X = matrix(c(ones, x), ncol = 2)
  beta = solve(t(X) %*% X) %*% t(X) %*% as.matrix(y)
  return(beta)
}

b0 = 0.42
b1 = 0.95
n = 198
x = log(Data$di[-199])
xT = log(Data$di[199])

beta <- matrix(0, 2, 1000)
# 初始化残差矩阵：n行（观测数），1000列（迭代次数）
residuals_mat <- matrix(0, nrow = n, ncol = 1000)

for (i in 1:1000) {
  set.seed(i) 
  u <- rnorm(n, mean = 0, sd = 1) 
  y <- b0 + b1 * x + u
  beta[, i] <- ols(y, x) 
  # 将当前迭代的残差存入对应列
  residuals_mat[, i] <- y - (beta[1,i] + beta[2,i] * x)
}

额外说明 #

• 如果更偏好列表存储（比如每次迭代的残差作为列表的一个元素），也可以初始化residuals_list <- vector("list", 1000)，然后在循环中赋值residuals_list[[i]] <- y - (beta[1,i] + beta[2,i] * x)。
• 你的ols函数可以简化，比如直接用lm函数的结果提取系数，但如果是手动实现OLS的话当前写法没问题。

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