我想用一些平均向量和协方差矩阵来模拟来自多元正态分布的两组数据。我想要 N(=3000) 模拟,每个 k(=100) 变量的 n(=30) 个值。所以会有两组 100 个变量,有 30 个观察值。然后我想在这两组之间进行测试。我想做 N 次同样的过程。所以我应该得到每个变量的 N p 个值,即 p 将是一个 2000X100 的矩阵。我怎么能在 r 中做到这一点?我有 ar 代码,但它不起作用。
library('MASS')
N=3000
n=30
k=100
cov <- matrix(.5,ncol=k,nrow=k);diag(cov) <- 1
mu1=rep(0,k)
mu2=rep(0,k)
for (j in 1:N){
x=mvrnorm(n, mu1, cov, tol = 1e-6, empirical = FALSE)
y=mvrnorm(n, mu2, cov, tol = 1e-6, empirical = FALSE)
xm<- apply(x,2,mean)
ym<- apply(y,2,mean)
xv<- apply(x,2,sd)
yv<- apply(y,2,sd)
s=sqrt((xv^2*(n-1)+yv^2*(n-1))/(2*n-2))
t=(xm-ym)/(s*sqrt(2/n))
p=2*(1-pt(abs(t),df=2*n-2))
}