r - 生成具有给定相关性的随机变量 :

标签 r statistics simulation genetics

我想生成 2 个连续的随机变量 Q1 , Q2 (数量性状,每个都是正常的)和2个二元随机变量Z1 , Z2 (二元特征)在所有可能的对之间具有给定的成对相关性。

(Q1,Q2):0.23 
(Q1,Z1):0.55 
(Q1,Z2):0.45 
(Q2,Z1):0.4 
(Q2,Z2):0.5 
(Z1,Z2):0.47 

请帮我在 R 中生成这样的数据。

最佳答案

这很粗糙,但可能会让您朝着正确的方向开始。

library(copula)

options(digits=3)
probs <- c(0.5,0.5)
corrs <- c(0.23,0.55,0.45,0.4,0.5,0.47)  ## lower triangle

模拟相关值(前两个定量,后两个转换为二进制)
sim <- function(n,probs,corrs) {
    tmp <- normalCopula( corrs, dim=4 , "un")
    getSigma(tmp) ## test
    x <- rCopula(1000, tmp)
    x2 <- x
    x2[,3:4] <- qbinom(x[,3:4],size=1,prob=rep(probs,each=nrow(x)))
    x2
}

测试观察到的和目标相关性之间的 SSQ 距离:
objfun <- function(corrs,targetcorrs,probs,n=1000) {
    cc <- try(cor(sim(n,probs,corrs)),silent=TRUE)
    if (is(cc,"try-error")) return(NA)
    sum((cc[lower.tri(cc)]-targetcorrs)^2)
}

当输入 corrs=target 时,看看事情有多糟糕:
cc0 <- cor(sim(1000,probs=probs,corrs=corrs))
cc0[lower.tri(cc0)]
corrs
objfun(corrs,corrs,probs=probs) ## 0.112

现在尝试优化。
opt1 <- optim(fn=objfun,
              par=corrs,
              targetcorrs=corrs,probs=c(0.5,0.5))
opt1$value     ## 0.0208

在 501 次迭代后停止并显示“超出最大迭代次数”。这永远不会很好地工作,因为我们试图在随机目标函数上使用确定性爬山算法......
cc1 <- cor(sim(1000,probs=c(0.5,0.5),corrs=opt1$par))
cc1[lower.tri(cc1)]
corrs

也许尝试模拟退火?
opt2 <- optim(fn=objfun,
              par=corrs,
              targetcorrs=corrs,probs=c(0.5,0.5),
              method="SANN")

它似乎并没有比以前的值好多少。两个可能的问题(留给读者作为练习)(1)我们已经指定了一组与我们选择的边缘分布不可行的相关性,或者(2)目标函数表面的误差正在进入方式 - 为了做得更好,我们必须对更多的重复进行平均(即增加 n )。

关于r - 生成具有给定相关性的随机变量 :,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/23846435/

相关文章:

r - 带有 knitr 字幕的循环

statistics - 我如何计算这些统计数据?

R:计算积分累计和变化?

python - 在遍历列表时如何有效地删除元素?

Java - 如何让类 1 和类 2 互相了解?

r - 获取R中字符串第一个大写字母的索引?

R - 图中的 "max subgraphs"

algorithm - Netlogo,创建避障算法

R:使用 ggplot2/ggmap 的世界地图 - 如何加载 png 图像作为 map

statistics - 访问前瞻性测试统计数据