2016-07-06 4 views
0

ダブルループを作成する方法を知りたいと思います。 私のコードでは、私は1000サンプル(各サンプルサイズ:25)の重回帰を行います 次に、null仮説を持つ1000個のサンプルのそれぞれのt検定値を作成します:sample = 'real' 。私は、モンテカルロシミュレーション(ベータ3 =回帰の第3の係数の値)から「実際の」ベータ3の値を知っています。 しかし、コードはこれまでのところ動作します。 サンプルサイズ50,100,250,500,1000(各サンプルサイズは1000回)について同じ手順を実行したいと思います。 ループでこの目標を実現する方法を教えてください。あなたが私を助けることができれば嬉しいです!ここに私のコードを見ることができます:ダブルループの作成方法は?

n <- 25 
B <- 1000 
beta3 <- 1.01901 #'real' beta3 value 

t.test.values <- rep(NA, B) 
for(rep in 1:B){ 

##data generation 
    d1 <- runif(25, 0, 1) 
    d2 <- rnorm(25, 0, 1) 
    d3 <- rchisq(25, 1, ncp=0) 
    x1 <- (1 + d1) 
    x2 <- (3 * d1 + 0.6 * d2) 
    x3 <- (2 * d1 + 0.6 * d3) 
    exi <- rchisq(25, 5, ncp = 0) 
    y <- beta0 + beta1*x1 + beta2*x2 + beta3*x3 + exi 

## estimation 
    lmobj  <- lm(y ~ x1 + x2 + x3)   

## extraction 
    betaestim <- coefficients(lmobj)[2:4] 
    betavar <- vcov(lmobj)[2:4, 2:4] 

## t-test 
    t.test.values[rep] <- (betaestim[3] - beta3)/sqrt((betavar)[9]) 

    } 
+2

_4ヶ月_ 06JUL16に元の質問を投稿し、@bouncyballの回答を受け入れた後、あなたは17NOV16で大量にQを変更しました。変更を元に戻し、新しい質問を提出してください。 – Uwe

答えて

0

私は結果を格納するためにdata.frameを使用することができます。また、beta0beta1、またはbeta2の値が含まれていないので、ちょうどプレースホルダの値を使用しました。

n <- c(50,100,250,500,1000) #how big are our sample sizes? 
B <- 1000 
beta3 <- 1.01901 #'real' beta3 value 
#other beta values (note that these were not included in your question) 
beta1 <- 2 
beta2 <- 4 
beta0 <- 6 

iter <- 1 

#initialize our results data.frame 
result_df <- data.frame(sample_size = numeric(length(n) * B), 
         t.test.values = numeric(length(n) * B) 
         ) 

for(size in n){ 

for(rep in 1:B){ 

    ##data generation 
    d1 <- runif(size, 0, 1) 
    d2 <- rnorm(size, 0, 1)  
    d3 <- rchisq(size, 1, ncp=0)  
    x1 <- (1 + d1)  
    x2 <- (3 * d1 + 0.6 * d2)  
    x3 <- (2 * d1 + 0.6 * d3)  
    exi <- rchisq(size, 5, ncp = 0)  
    y <- beta0 + beta1*x1 + beta2*x2 + beta3*x3 + exi 

    ## estimation 
    lmobj  <- lm(y ~ x1 + x2 + x3)   

    ## extraction 
    betaestim <- coefficients(lmobj)[2:4] 
    betavar <- vcov(lmobj)[2:4, 2:4] 

    ## store our values 
    result_df[iter, 1] <- size 

    result_df[iter, 2] <- (betaestim[3] - beta3)/sqrt((betavar)[9]) 

    iter = iter + 1 #iterate 

    } 
} 

限り、あなたは(私はここiterを使用しました)反復処理を追跡するために何かを使用して、二重のforループは悪くないです。 data.frameを正しいサイズに初期化してください。より多くのシミュレーションを行う予定がある場合は、replicate関数と*apply関数クラスを調べると役に立ちます。

関連する問題