R语言生成指定相关系数三变量条件分布的实现方法

现有代码偏差原因 #

你当前代码计算得到的cor(y,x)和目标值不符，核心原因是噪声项的方差未做校准：
当x ~ N(0,1)、噪声项e ~ N(0,1)时，若构造y = r_yx * x + e，则Var(y) = r_yx² + 1，x与y的协方差为r_yx，实际相关系数为r_yx / sqrt(r_yx² + 1)。以目标值0.9为例，实际相关系数仅约为0.67，自然存在明显偏差。

实现方案 #

要生成x、y、z三个变量，精准匹配你设定的三个两两相关系数，最稳定通用的方法是基于Cholesky分解多元正态相关矩阵，步骤如下：

• 按变量顺序构造3*3的目标相关矩阵，对角线固定为1，非对角线位置填入你设定的相关系数
• 先校验矩阵合法性：相关矩阵必须满足半正定要求，否则不存在符合条件的多元分布，需要调整相关系数取值
• 对合法的相关矩阵做Cholesky分解得到变换矩阵，将独立生成的标准正态随机数乘变换矩阵，即可得到符合目标相关系数的变量，样本量越大，实际计算的相关系数和目标值偏差越小

可直接运行的R代码 #

# 自定义目标相关系数
r_yx <- 0.9
r_zx <- -0.5
r_yz <- 0.9
n <- 10000

# 构造顺序为x/y/z的相关矩阵
cor_mat <- matrix(c(
  1,    r_yx, r_zx,
  r_yx, 1,    r_yz,
  r_zx, r_yz, 1
), nrow = 3, byrow = TRUE)

# 合法性校验
eigen_vals <- eigen(cor_mat)$values
if (any(eigen_vals < -1e-8)) {
  stop("当前设置的相关系数不合法，不存在对应分布，请调整取值")
}

# Cholesky分解变换
chol_mat <- chol(cor_mat)
set.seed(123) # 固定随机种子用于复现，不需要可删除
random_mat <- matrix(rnorm(n * 3), ncol = 3)
res <- random_mat %*% chol_mat

# 提取三个变量
x <- res[, 1]
y <- res[, 2]
z <- res[, 3]

# 校验结果，输出值会非常接近设定的目标值
cor(x, y)
cor(x, z)
cor(y, z)

逐变量生成的修正逻辑 #

如果你希望保留逐变量生成的写法，只需要对每一步的残差方差做校准：

• 生成标准正态x后，构造y时将残差标准差设为sqrt(1 - r_yx^2)，即y = r_yx * x + sqrt(1 - r_yx^2) * rnorm(n)，此时x和y的相关系数会精准等于r_yx
• 后续生成z时，需要同时根据r_zx、r_yz两个目标值推导z对x、y的回归系数，再校准残差方差，推导过程相对繁琐，更推荐直接使用上述矩阵分解方法。

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