ifa: Store latent mean & cov in regular MxMatrix
[openmx:openmx.git] / models / failing / irt-sat12.R
1 # echo 0 > /proc/self/coredump_filter  # normally 023
2 # R --vanilla --no-save -f models/failing/bock-aitkin-1981.R
3 # R -d gdb --vanilla --no-save -f models/failing/bock-aitkin-1981.R
4
5 #options(error = browser)
6 require(mirt)
7 require(OpenMx)
8 require(rpf)
9 library(stringr)
10
11 data(SAT12)
12 data <- key2binary(SAT12,
13                    key = c(1,4,5,2,3,1,2,1,3,1,2,4,2,1,5,3,4,4,1,4,3,3,4,1,3,5,1,3,1,5,4,5))
14 df <- list()
15 for (cx in 1:dim(data)[2]) {
16         df[[cx]] <- factor(data[,cx], levels=c(0,1))
17 }
18 df <- as.data.frame(df)
19 colnames(df) <- str_replace_all(colnames(data), "\\.", "")
20
21 numItems <- 32
22 numPersons <- dim(data)[1]
23 maxDim <- 2
24
25 items <- vector("list", numItems)
26 for (ix in 1:numItems) {
27         items[[ix]] <- rpf.drm(dimensions=maxDim)
28 }
29
30 maxParam <- max(vapply(items, function(i) i@numParam, 0))
31 maxOutcomes <- max(vapply(items, function(i) i@numOutcomes, 0))
32
33 design <- matrix(c(rep(1,numItems),
34                    1+c(2,3,2,3,3,2,1,2,1,1,1,3,1,3,1,2,1,1,3,3,1,1,3,1,3,3,1,3,2,3,1,2)), byrow=TRUE, nrow=2)
35
36 spec <- mxMatrix(name="ItemSpec", nrow=3, ncol=numItems,
37          values=c(rep(mxLookupIRTItemModelID("drm"), numItems),
38                  rep(2, numItems),
39                  rep(2, numItems)),
40          free=FALSE, byrow=TRUE)
41
42 design <- mxMatrix(name="Design", nrow=maxDim, ncol=numItems,
43                    values=design)
44
45 ip.mat <- mxMatrix(name="ItemParam", nrow=maxParam, ncol=numItems,
46                    values=c(
47                      rep(1.414,numItems),
48                      rep(1,numItems),
49                      rep(0,numItems),
50                      rep(0,numItems)),
51          free=c(rep(TRUE, numItems*3),
52                 rep(FALSE, numItems*1)),
53          byrow=TRUE)
54
55 #ip.mat@values[2,1] <- correct.mat[2,1]
56 #ip.mat@free[2,1] <- FALSE
57
58 m1 <- mxModel(model="sat12",
59           spec, design,
60           ip.mat,
61           mxData(observed=df, type="raw"),
62           mxExpectationBA81(
63              ItemSpec="ItemSpec",
64              Design="Design",
65              ItemParam="ItemParam"),
66           mxFitFunctionBA81()
67 )
68
69 m1 <- mxOption(m1, "Analytic Gradients", 'no')
70 if (1) {
71         m1 <- mxOption(m1, "Analytic Gradients", 'yes')
72         m1 <- mxOption(m1, "Verify level", '-1')
73 }
74 m1 <- mxOption(m1, "Function precision", '1.0E-5')
75 m1 <- mxOption(m1, "Calculate Hessian", "No")
76 m1 <- mxOption(m1, "Standard Errors", "No")
77
78 m1 <- mxRun(m1, silent=TRUE)
79 print(m1@matrices$ItemParam@values)
80 q()
81
82 -2LL 17405.6316