Hello All, I have searched many help forums, message boards, etc. and I just can't apply the comments to what I need my program to do. I am running R 2.5.1 on an XP system, and my desire is to produce replicate datasets for a simulation study I am running. Essentially, I have sets of parameters (a's, b's, and c's) that define a function which produces a decimal value. This value is compared to a random uniform value, and is coded a 1 if the function is greater than the uniform value, 0 if it is <= to the uniform value. My code thus far works great, but I just need it to run several times. Here we go:
library(mvtnorm) library(sm) library(ltm) library(irtoys) k<- 5000 set.seed(271828) t <- rmvnorm(n=k,mean=c(-1,0,1),sigma=matrix(c(1,.8,.5,.8,1,.8,.5,.8,1),3,3)) #Using mv here because of the likely association of ability (theta = t) across time. t1<-as.matrix(t[,1]) t2<-as.matrix(t[,2]) t3<-as.matrix(t[,3]) set.seed(271828) # Population item parameters (n=54) from which we will select relevant items # These are the parameters that are used in the function a <- c(1.18120, 0.92613, 0.96886, 0.80503, 1.12384, 0.84073, 0.85544, 0.86801, 1.01054, 0.82278, 1.10353, 0.78865, 0.98421, 1.76071, 0.89603, 0.84671, 0.89737, 0.74775, 0.32190, 0.69730, 0.72059, 1.16762, 1.29257, 1.32902, 0.59540, 0.51022, 0.59259, 0.93951, 0.68568, 0.55649, 0.88084, 0.52940, 0.45735, 0.57560, 1.11779, 0.96984, 1.19692, 0.99102, 1.25847, 1.62555, 0.63049, 1.07807, 1.04897, 1.23138, 1.14014, 1.25230, 1.14844, 0.59287, 0.83143, 0.81723, 0.52141, 0.61980, 0.49945, 1.02749) b <- c(-2.51737, -1.95897, -1.72667, -0.82988, -0.36093, 0.72554, 0.91442, 0.78061, 0.06088, 0.75733, -0.76371, 0.24552, -0.42050, 0.88232, -0.81761, 0.06466, -0.43866, -0.46042, 0.21636, -0.73147, -1.44086, -1.03718, 0.07275, -0.17197, 1.53796, -0.45631, -1.69826, -0.66506, 0.98921, 0.30714, -0.62245, 0.97253, 1.95894, 0.21277, 1.96346, 1.18825, 1.59917, -0.28401, -1.23530, -0.09671, -0.31581, -0.66149, -0.81284, -0.35399, -0.07623, 1.06442, -0.68559, 1.07591, 0.97458, 0.06436, 1.25622, 1.73954, 1.75052, 2.34088) c <- c(0.00000, 0.00000, 0.00000, 0.00000, 0.19648, 0.31302, 0.26454, 0.19714, 0.06813, 0.21344, 0.00000, 0.03371, 0.00000, 0.16581, 0.11054, 0.08756, 0.07115, 0.26892, 0.00000, 0.06883, 0.00000, 0.14815, 0.32389, 0.19616, 0.17597, 0.00000, 0.00000, 0.04337, 0.19949, 0.20377, 0.00000, 0.06243, 0.13639, 0.00000, 0.18166, 0.15996, 0.20184, 0.08331, 0.24453, 0.26114, 0.16434, 0.20750, 0.32658, 0.31870, 0.45227, 0.35039, 0.31178, 0.17999, 0.22774, 0.21675, 0.10153, 0.17764, 0.15205, 0.19858) # Item parameters for generating 3PL data for all five testing occasions: # This selects the relevant parameters for a particular data generation run # Only parameters for the first testing occasion are shown to save space a1 <- as.matrix(a[c(1:5,15:20,22:24,38:44)]) b1 <- as.matrix(b[c(1:5,15:20,22:24,38:44)]) c1 <- as.matrix(c[c(1:5,15:20,22:24,38:44)]) # Here is where I would like to begin my replications, but don't know how to make R do it. # The code below produces a matrix of 0's and 1's (which will be used by another program) # I would like to nest this in a "do loop" such that, say, 30 replicate datasets are produced using the # same parameters. N <- nrow(t1) # number of examinees n <- nrow(a1) # number of items d <- 1.7 theta <- t1 response <- matrix (0,N,n) uni <- matrix (runif(N*n),nrow = N) for (i in 1:N) { for (j in 1:n) { if ( c1[j]+(1-c1[j])/(1+exp(-d*a1[j]*(theta[i]-b1[j]))) > uni[i,j] ) response[i,j] = 1 else response[i,j] = 0 } } write.table(response, file="C:/responses.dat", sep=" ",row.names=FALSE, col.names=FALSE) I tried earlier nesting this in another for loop, but that indexes elements of matrices and vectors, and doesn't seem to apply to a "global" loop methodology. I am attempting to use replicate as we speak, but documentation is sparse (help("replicate") is nested in lapply information). Any guidance is greatly appreciated. Thanks in advance, Jonathan Beard -- View this message in context: http://www.nabble.com/Generating-Replicate-Datasets-%28using-loops-or-other-means%29-tf4418768.html#a12603580 Sent from the R help mailing list archive at Nabble.com. ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.