Compute condition number of information matrix
[openmx:openmx.git] / models / nightly / 3LatentMultiRegWith2LevelModerator-a.R
1 # ---------------------------------------------------------------------
2 # Program: 3LatentMultiRegWithModerator100521.R
3 #  Author: Steven M. Boker
4 #    Date: Sat May 22 11:13:51 EDT 2010
5 #
6 # This program tests variations on a latent variable multiple regression
7 #    using a standard RAM.
8 #
9 # ---------------------------------------------------------------------
10 # Revision History
11 #    -- Sat May 22 11:13:51 EDT 2010
12 #      Created 3LatentMultiRegWithModerator100521.R.
13 #
14 # ---------------------------------------------------------------------
15
16 # ----------------------------------
17 # Read libraries and set options.
18
19 library(OpenMx)
20
21 options(width=100)
22
23 # ---------------------------------------------------------------------
24 # Data for multiple regression of F3 on F1 and F2 with moderator variable Z.
25
26 numberSubjects <- 1000
27 numberIndicators <- 12
28 numberFactors <- 3
29
30 fixedBMatrixF <- matrix(c(.4, .2), 2, 1, byrow=TRUE)
31 randomBMatrixF <- matrix(c(.3, .5), 2, 1, byrow=TRUE)
32 XMatrixF <- matrix(rnorm(numberSubjects*2, mean=0, sd=1), numberSubjects, 2)
33 UMatrixF <- matrix(rnorm(numberSubjects*1, mean=0, sd=1), numberSubjects, 1)
34 Z <- matrix(floor(runif(numberSubjects, min=0, max=1.999)), nrow=numberSubjects, ncol=2)
35
36 XMatrix <- cbind(XMatrixF, XMatrixF %*% fixedBMatrixF + (XMatrixF*Z) %*% randomBMatrixF + UMatrixF)
37
38 BMatrix <- matrix(c( 1, .6, .7, .8,  0,  0,  0,  0,  0,  0,  0,  0,
39                      0,  0,  0,  0,  1, .5, .6, .7,  0,  0,  0,  0,
40                      0,  0,  0,  0,  0,  0,  0,  0,  1, .7, .6, .5), numberFactors, numberIndicators, byrow=TRUE)
41 UMatrix <- matrix(rnorm(numberSubjects*numberIndicators, mean=0, sd=1), numberSubjects, numberIndicators)
42 YMatrix <- XMatrix %*% BMatrix + UMatrix
43
44 cor(cbind(XMatrix,Z[,1]))
45
46 dimnames(YMatrix) <- list(NULL, paste("X", 1:numberIndicators, sep=""))
47
48 latentMultiRegModerated1 <- data.frame(YMatrix,Z=Z[,1])
49
50 round(cor(latentMultiRegModerated1), 3)
51 round(cov(latentMultiRegModerated1), 3)
52
53 latentMultiRegModerated1$Z <- latentMultiRegModerated1$Z - mean(latentMultiRegModerated1$Z)
54
55 numberFactors <- 3
56 numberIndicators <- 12
57 numberModerators <- 1
58 indicators <- paste("X", 1:numberIndicators, sep="")
59 moderators <- c("Z")
60 totalVars <- numberIndicators + numberFactors + numberModerators
61
62 # ----------------------------------
63 # Build an orthogonal simple structure factor model
64
65 latents <- paste("F", 1:numberFactors, sep="")
66
67 uniqueLabels <- paste("U_", indicators, sep="")
68 meanLabels <- paste("M_", latents, sep="")
69 factorVarLabels <- paste("Var_", latents, sep="")
70
71 latents1 <- latents[1]
72 indicators1 <- indicators[1:4]
73 loadingLabels1 <- paste("b_F1", indicators[1:4], sep="") 
74 latents2 <- latents[2]
75 indicators2 <- indicators[5:8]
76 loadingLabels2 <- paste("b_F2", indicators[5:8], sep="") 
77 latents3 <- latents[3]
78 indicators3 <- indicators[9:12]
79 loadingLabels3 <- paste("b_F3", indicators[9:12], sep="") 
80
81 threeLatentOrthogonal <- mxModel("threeLatentOrthogonal",
82     type="RAM",
83     manifestVars=c(indicators),
84     latentVars=c(latents,"dummy1"),
85     mxPath(from=latents1, to=indicators1, 
86            arrows=1, connect="all.pairs", 
87            free=TRUE, values=.2, 
88            labels=loadingLabels1),
89     mxPath(from=latents2, to=indicators2, 
90            arrows=1, connect="all.pairs",
91            free=TRUE, values=.2, 
92            labels=loadingLabels2),
93     mxPath(from=latents3, to=indicators3, 
94            arrows=1, connect="all.pairs",
95            free=TRUE, values=.2, 
96            labels=loadingLabels3),
97     mxPath(from=latents1, to=indicators1[1], 
98            arrows=1, 
99            free=FALSE, values=1),
100     mxPath(from=latents2, to=indicators2[1], 
101            arrows=1, 
102            free=FALSE, values=1),
103     mxPath(from=latents3, to=indicators3[1], 
104            arrows=1, 
105            free=FALSE, values=1),
106     mxPath(from=indicators, 
107            arrows=2, 
108            free=TRUE, values=.2, 
109            labels=uniqueLabels),
110     mxPath(from=latents,
111            arrows=2, 
112            free=TRUE, values=.8, 
113            labels=factorVarLabels),
114     mxPath(from="one", to=indicators, 
115            arrows=1, free=FALSE, values=0),
116     mxPath(from="one", to=c(latents), 
117            arrows=1, free=TRUE, values=.1, 
118            labels=meanLabels),
119     mxData(observed=latentMultiRegModerated1, type="raw")
120     )
121
122 threeLatentOrthogonalOut <- mxRun(threeLatentOrthogonal)
123