Added copyright information to test suite models
[openmx:openmx.git] / demo / BivariateSaturated_MatrixRawCholesky.R
1 #
2 #   Copyright 2007-2010 The OpenMx Project
3 #
4 #   Licensed under the Apache License, Version 2.0 (the "License");
5 #   you may not use this file except in compliance with the License.
6 #   You may obtain a copy of the License at
7
8 #        http://www.apache.org/licenses/LICENSE-2.0
9
10 #   Unless required by applicable law or agreed to in writing, software
11 #   distributed under the License is distributed on an "AS IS" BASIS,
12 #   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 #   See the License for the specific language governing permissions and
14 #   limitations under the License.
15
16
17 # -----------------------------------------------------------------------
18 # Program: BivariateSaturated_MatrixRawCholesky.R  
19 #  Author: Hermine Maes
20 #    Date: 08 01 2009 
21 #
22 # Bivariate Saturated model to estimate means and (co)variances 
23 # using Cholesky Decomposition
24 # Matrix style model input - Raw data input
25 #
26 # Revision History
27 #   Hermine Maes -- 10 08 2009 updated & reformatted
28 # -----------------------------------------------------------------------
29
30 require(OpenMx)
31
32 #Simulate Data
33 # -----------------------------------------------------------------------
34 require(MASS)
35 set.seed(200)
36 rs=.5
37 xy <- mvrnorm (1000, c(0,0), matrix(c(1,rs,rs,1),2,2))
38 testData <- xy
39 selVars <- c("X","Y")
40 dimnames(testData) <- list(NULL, selVars)
41 summary(testData)
42 cov(testData)
43
44 #example 6: Saturated Model with Raw Data and Matrix-Style Input
45 # -----------------------------------------------------------------------
46 bivSatModel6 <- mxModel("bivSat6",
47     mxMatrix(
48                 type="Full", 
49                 nrow=2, 
50                 ncol=2, 
51                 free=c(TRUE,TRUE,FALSE,TRUE), 
52                 values=c(1,.2,0,1), 
53                 name="Chol"
54         ),
55         mxMatrix(
56                 type="Full", 
57                 nrow=1, 
58                 ncol=2, 
59                 free=TRUE, 
60                 values=c(0,0), 
61                 name="expMean"
62         ),
63     mxAlgebra(
64                 expression=Chol %*% t(Chol), 
65                 name="expCov"
66         ),
67         mxData(
68                 observed=testData, 
69                 type="raw"
70         ),
71         mxFIMLObjective(
72                 covariance="expCov",
73                 means="expMean",
74                 dimnames=selVars
75         )
76 )
77
78 bivSatFit6 <- mxRun(bivSatModel6)
79 EM <- mxEval(expMean, bivSatFit6)
80 EC <- mxEval(expCov, bivSatFit6)
81 LL <- mxEval(objective,bivSatFit6)
82
83 #Mx answers hard-coded
84 # -----------------------------------------------------------------------
85 #example Mx..2: Saturated Model with Raw Data
86 Mx.EM <- matrix(c(0.03211188, -0.004889211),1,2)
87 Mx.EC <- matrix(c(1.0092891, 0.4813504, 0.4813504, 0.9935366),2,2)
88 Mx.LL <- 5415.772
89
90 # Compare OpenMx results to Mx results 
91 # -----------------------------------------------------------------------
92 # (LL: likelihood; EC: expected covariance, EM: expected means)
93 #6: RawMat Cholesky
94 omxCheckCloseEnough(LL,Mx.LL,.001)
95 omxCheckCloseEnough(EC,Mx.EC,.001)
96 omxCheckCloseEnough(EM,Mx.EM,.001)