1

私はRの初心者で、シミュレーション研究を行っており、自分がやりたいことを1つのサンプルで作成することができました。しかし、自分がやったことをどのように再現すればよいのかわかりません。

これまでに書いたプログラムは次のとおりです。

I <- 500       # number of observations
J <- 18        # total number of items
K <- 6         # number of testlets
JK <-3         # number of items within a testlet
response <- matrix(0, I, J)  # null binary (0, 1) response matrix 
unit <- matrix(1, JK, 1)     # unit vector

set.seed(1234)

# Multidimensional 3-pl model
pij <- function(a,b,c,theta,gamma) {c+(1-c)*(1/(1+exp(-1.7*a*(theta-b-gamma))))}

# Assigning a and b parameter values
a <- c(.8,.9,.7,.8,.9,.7,.8,.9,.7,.8,.9,.7,.8,.9,.7,.8,.9,.7)
b <-c(1,0,-1.5,1,0,-1.5,1,0,-1.5,1,0,-1.5,1,0,-1.5,1,0,-1.5)
# Assigning c-parameter, each 3 items (c-parameter & testlet effect)
#(small&small, small&large, large&small, large&large, mixed&small, mixed&large)
c <- c(.2,.2,.2,.2,.2,.2,.5,.5,.5,.5,.5,.5,.2,.33,.5,.2,.33,.5)    

theta <- rnorm(I, 0, 1)   # random sampling theta-values from normal dist. M=0, SD=1

gamma1 <- rnorm(I, 0, .2)  # small testlet effect: random sampling gamma from normal dist. M=0, SD=.2
gamma2 <- rnorm(I, 0, 1)   # large testlet effect: random sampling gamma from normal dist. M=0, SD=1
gamma3 <- rnorm(I, 0, .2)  # small testlet effect: random sampling gamma from normal dist. M=0, SD=.2
gamma4 <- rnorm(I, 0, 1)   # large testlet effect: random sampling gamma from normal dist. M=0, SD=1
gamma5 <- rnorm(I, 0, .2)  # small testlet effect: random sampling gamma from normal dist. M=0, SD=.2
gamma6 <- rnorm(I, 0, 1)   # large testlet effect: random sampling gamma from normal dist. M=0, SD=1

# implementing that the testlet effect is same for the items within a testlet
gamma1T <- gamma1 %*% t(unit)
gamma2T <- gamma2 %*% t(unit)
gamma3T <- gamma3 %*% t(unit)
gamma4T <- gamma4 %*% t(unit)
gamma5T <- gamma5 %*% t(unit)
gamma6T <- gamma6 %*% t(unit)

gammaT <- matrix(c(gamma1T, gamma2T, gamma3T, gamma4T, gamma5T, gamma6T), I, J)  # getting all the gammas together in a large matrix

# Generating data using the information above
for(i in 1:I) {
  for(j in 1:J) {
    response[i, j] <- ifelse(pij(a=a[j], b=b[j], c=c[j], theta=theta[i], gamma=gammaT[i,j]) < runif(1), 0, 1)
  }
}

したがって、データセット「応答」を取得します。私がやりたいのは、これを複製して、たとえば1000個の「応答」データセットを取得することです。これは、「シータ」と「ガンマ」のランダムサンプリングを複製することで実現できると思いますが、実際にこれを行うことは考えていません。

よろしくお願いします。

半条。

4

2 に答える 2

4

Stedyのアドバイスは、1つのことを除けば、合理的です。forループでシードをインクリメントしないでください。

Stedyの提案を理解しているのと同じように、シミュレーションごとset.seed(i)にループ内で呼び出され、反復ごとに増分されます。この戦略は、生成されたシーケンス間の相関により、(潜在的に大きな)バイアスを導入することが知られています。fori

代わりに、シードを最初、つまりforループの前に1回設定します。たとえば、現在のUnix時間をシードとして使用したり、乱数を含むファイル(random.orgなど)から読み取ることができます。また、シードを結果とともに保存するようにしてください。たとえば、ログファイルに出力します。以前の一連の複製の正確な結果を再度再現する場合は、対応するシードを設定するだけです。

他の人が結果を正確に複製できるようにする場合は、使用したRバージョン(および場合によってはオペレーティングシステム)も指定する必要があります(RNGの実装は異なる場合があるため)。

余談ですが、シミュレーションレプリケーションは驚異的並列タスクです。つまり、マルチコアマシン(rparallelなど)を使用している場合は、レプリケーションを簡単に並列で実行できます。ただし、この場合、乱数に関する特別な注意が必要です(たとえば、このペーパーを参照してください)。

于 2012-02-22T11:59:34.850 に答える
1

ローカル変数を取得して関数にします。次に、for()ループを作成し、関数を呼び出し、set.seed()関数が呼び出されるたびにループの長さを 1 ずつfor()増やします。

于 2012-02-22T07:12:41.133 に答える