R语言蒙特卡洛仿真存储结果到矩阵的实现及报错解决

• 循环语法硬错：原代码写for (i in iter)，R的for循环需要遍历一个向量，iter是值为2000的单元素向量，这个写法只会执行1次迭代（i=2000），完全达不到2000轮仿真的效果。
• 存储索引方向颠倒：预分配的results.matrix为2000行×12列结构，行对应单轮仿真的输出，列对应该轮要存储的各个指标。原代码results.matrix[,1]是固定取矩阵第一列，将当前轮次的计算结果往整列写入，也就是要把长度为11个的p值向量，填充到长度为2000的列空间里——2000无法被11整除，直接触发"number of items to replace is not a multiple of replacement length"报错。
• 冗余计算：单独执行summary(glm0)只会把模型摘要打印到控制台，没有存入对象，后续取系数p值时会重新计算一遍模型摘要，平白浪费计算资源。
• 隐藏适配错误：原代码生成的Y是二分类结果，但glm没有指定逻辑斯蒂回归的族参数；如果是按默认线性回归拟合，p值列为Pr(>|t|)是正确的，如果后续改逻辑斯蒂回归，需要对应取Pr(>|z|)列，取错列会得到空值。
• 无列名的矩阵存储后无法溯源：没有给预分配的矩阵设置列名，后续拿到结果很容易混淆每列对应的指标含义。

核心逻辑遵循蒙特卡洛仿真的标准操作：提前预分配固定大小的存储矩阵，给每列绑定明确的指标名称，每轮循环计算完成后，将结果按索引写入对应行，避免动态扩矩阵的性能损耗。
修复后的可直接运行代码：

iter <- 2000
# 提前定义所有待存储的指标，共12列，和原设定的列数匹配
col_indicators <- c(
  "p_intercept",
  "p_x1", "p_x2", "p_x3", "p_x4", "p_x5", "p_x6", "p_x7",
  "p_age", "p_gender", "p_inkomst", "model_converge"
)
# 预分配带列名的结果矩阵
results.matrix <- matrix(
  ncol = length(col_indicators),
  nrow = iter,
  dimnames = list(NULL, col_indicators) # 列绑定名称，避免后续混淆
)

# 修正循环遍历逻辑：遍历1到iter的完整序列
for (i in 1:iter) {
  # 生成仿真因变量
  Ylink <- (1 * (age - mean(age)) + (1.5 * (gender == "1")) + 1.2 * (inkomst - mean(inkomst)))
  Yprob <- 1/(1 + exp(-Ylink))
  Y <- rbinom(2000, size = 1, prob = Yprob)
  df <- data.frame(x1, x2, x3, x4, x5, x6, x7, age, gender, inkomst, Y)

  # 拟合模型，直接把摘要结果存入对象
  # * 如果是拟合二分类逻辑斯蒂模型，补充参数 family = binomial(link = "logit")
  glm0 <- glm(Y ~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + age + gender + inkomst, data = df)
  glm0_result <- summary(glm0)

  # 将当前轮次结果写入矩阵的第i行，和列名一一对应
  # * 如果是逻辑斯蒂模型，把列名改为 "Pr(>|z|)" 即可
  results.matrix[i, 1:11] <- glm0_result$coef[, "Pr(>|t|)"]
  results.matrix[i, 12] <- glm0$converged # 存储模型是否收敛的状态
}

• 绝对不要在循环里用rbind动态拼接结果，2000轮规模下动态扩矩阵的速度会比预分配慢数十倍，预分配是R做仿真的通用高效方案。
• 牢记R矩阵的索引规则：矩阵对象[行索引, 列索引]，行位置对应仿真轮次，列位置对应单个指标，索引顺序搞反是新手最容易触发长度不匹配报错的原因。
• 正式跑全量仿真前，先把iter改成2-3轮小规模测试，确认每轮写入的结果长度和列数完全匹配，避免跑很久才发现索引错误浪费时间。

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