Nothing
#
# Copyright 2007-2019 by the individuals mentioned in the source code history
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
require(MASS)
require(OpenMx)
# Simulate univariate data
set.seed(100)
x <- rnorm(1000, 0, 1)
univData <- as.matrix(x)
dimnames(univData) <- list(NULL, "X")
summary(univData)
mean(univData)
var(univData)
# Simulate bivariate data
set.seed(200)
rs = .5
xy <- mvtnorm::rmvnorm(1000, c(0, 0), matrix(c(1, rs, rs, 1), 2, 2))
bivData <- xy
dimnames(bivData) <- list(NULL, c('X','Y'))
summary(bivData)
colMeans(bivData)
cov(bivData)
# hist(univData)
# Create ordinal and binary data from continuous data
univDataOrd <- data.frame(X=cut(univData[,1], breaks=5, ordered_result=TRUE, labels=c(0,1,2,3,4)) )
table(univDataOrd)
univDataBin <- data.frame(X=ifelse(univData[,1] >.5,1,0))
table(univDataBin)
bivDataOrd <- data.frame(bivData)
for (i in 1:2) {
bivDataOrd[,i] <- cut(bivData[,i], breaks=5, ordered_result=TRUE, labels=c(0,1,2,3,4))
}
table(bivDataOrd[,1],bivDataOrd[,2])
bivDataBin <- data.frame(bivData)
for (i in 1:2) {
bivDataBin[,i] <- ifelse(bivData[,i] >.5,1,0)
}
table(bivDataBin[,1],bivDataBin[,2])
# Create data objects for alternative data types
obsRawData <- mxData(observed=univData, type="raw")
# obsRawData
obsBivData <- mxData(observed=bivData, type="raw")
univDataOrdF <- mxFactor( x=univDataOrd, levels=c(0:4) )
univDataBinF <- mxFactor( x=univDataBin, levels=c(0,1) )
bivDataOrdF <- mxFactor( x=bivDataOrd, levels=c(0:4) )
bivDataBinF <- mxFactor( x=bivDataBin, levels=c(0,1) )
obsRawDataOrd <- mxData(observed = univDataOrdF, type="raw")
obsRawDataBin <- mxData(observed = univDataBinF, type="raw")
obsBivDataOrd <- mxData(observed = bivDataOrdF , type="raw")
obsBivDataBin <- mxData(observed = bivDataBinF , type="raw")
univDataCov <- var(univData)
obsCovData <- mxData(observed=univDataCov, type="cov", numObs=1000)
obsCovData <- mxData(observed=var(univData), type="cov", numObs=1000)
obsCovData
# Fit saturated model to covariance matrices using the path method
selVars <- c("X")
manifestVars = selVars
expVariance <- mxPath(from=c("X"), arrows=2, free=TRUE, values=1, lbound=.01, labels="vX" )
expVariance
univSatModel1 <- mxModel("univSat1", manifestVars=selVars, obsCovData, expVariance, type="RAM" )
univSatModel1
univSatModel1$matrices$S
univSatFit1 <- mxRun(univSatModel1)
univSatFit1$matrices$S$values
univSatFit1[['S']]$values
EC1 <- mxEval(S, univSatFit1)
SL1 <- univSatFit1$output$SaturatedLikelihood
LL1 <- mxEval(objective, univSatFit1)
Chi1 <- LL1-SL1
EC1
SL1
LL1
Chi1
summary(univSatFit1)
# Fit saturated model to covariance matrices and means using the path method
expMean <- mxPath(from="one", to="X", arrows=1, free = TRUE, values=0, labels="mX")
expMean
obsCovMeanData <- mxData( observed=var(univData), type="cov", numObs=1000, means=colMeans(univData) )
univSatModel1M <- mxModel(univSatModel1, name="univSat1M", expMean, obsCovMeanData )
univSatModel1M
univSatFit1M <- mxRun(univSatModel1M)
EM1m <- mxEval(M, univSatFit1M)
summary(univSatFit1M)
# Fit saturated model to raw data using the path method
univSatModel2 <- mxModel(univSatModel1M, obsRawData )
univSatFit2 <- mxRun(univSatModel2)
summary(univSatFit2)
EM2 <- mxEval(M, univSatFit2)
EC2 <- mxEval(S, univSatFit2)
LL2 <- mxEval(objective, univSatFit2)
EM2
EC2
LL2
# Fit saturated model to covariance matrices using the matrix method
expCovMat <- mxMatrix( type="Symm", nrow=1, ncol=1, free=TRUE, values=1, name="expCov" )
expCovMat
expectation <- mxExpectationNormal( covariance="expCov", dimnames=selVars )
MLobjective <- mxFitFunctionML( )
MLobjective
univSatModel3 <- mxModel("univSat3", obsCovData, expCovMat, expectation, MLobjective)
univSatFit3 <- mxRun(univSatModel3)
univSatSumm3 <- summary(univSatFit3)
univSatSumm3
# Fit saturated model to covariance matrices and means using the matrix method
expMeanMat <- mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE, values=0, name="expMean" )
expMeanMat
expectationM <- mxExpectationNormal( covariance="expCov", means="expMean", dimnames=selVars )
MLobjectiveM <- mxFitFunctionML( )
MLobjectiveM
univSatModel3M <- mxModel(univSatModel3, name="univSat3M", obsCovMeanData, expMeanMat, expectationM, MLobjectiveM)
univSatFit3M <- mxRun(univSatModel3M)
univSatSumm3M <- summary(univSatFit3M)
# Fit saturated model to raw data using the matrix method
NORMexpectation <- mxExpectationNormal( covariance="expCov", means="expMean", dimnames=selVars)
FIMLobjective <- mxFitFunctionML()
FIMLobjective
univSatModel4 <- mxModel("univSat4", obsRawData, expCovMat, expMeanMat, NORMexpectation, FIMLobjective)
univSatFit4 <- mxRun(univSatModel4)
summary(univSatFit4)
# Fit saturated model to raw binary data using the matrix method
expCovMatBin <- mxMatrix( type="Stand", nrow=1, ncol=1, free=TRUE, values=.5, name="expCov" )
expMeanMatBin <- mxMatrix( type="Zero", nrow=1, ncol=1, name="expMean" )
expThreMatBin <- mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE, values=0, name="expThre" )
expThreMatBin
expectationBin <- mxExpectationNormal( covariance="expCov", means="expMean", threshold="expThre", dimnames=selVars )
FIMLobjectiveBin <- mxFitFunctionML( )
univSatModel5 <- mxModel("univSat5", obsRawDataBin, expCovMatBin, expMeanMatBin, expThreMatBin, expectationBin, FIMLobjectiveBin )
univSatFit5 <- mxRun(univSatModel5)
summary(univSatFit5)
# Fit saturated model to raw ordinal data using the matrix method
nv <- 1
nth <- 4
expCovMatOrd <- mxMatrix( type="Stand", nrow=nv, ncol=nv, free=TRUE, values=.5, name="expCov" )
expMeanMatOrd <- mxMatrix( type="Zero", nrow=1, ncol=nv, name="expMean" )
expThreMatOrd <- mxMatrix( type="Full", nrow=nth, ncol=nv, free=TRUE, values=c(-1.5,-.5,.5,1.5), name="expThre" )
expThreMatOrd
expectationOrd <- mxExpectationNormal( covariance="expCov", means="expMean", threshold="expThre", dimnames=selVars )
FIMLobjectiveOrd <- mxFitFunctionML()
univSatModel6 <- mxModel("univSat6", obsRawDataOrd, expCovMatOrd, expMeanMatOrd, expThreMatOrd, expectationOrd, FIMLobjectiveOrd )
univSatFit6 <- mxRun(univSatModel6)
summary(univSatFit6)
# Fit saturated model to raw ordinal data using the matrix method with deviations/increments
lth <- -1.5; ith <- 1
thValues <- matrix(rep(c(lth,(rep(ith,nth-1)))),nrow=nth,ncol=nv)
thLBound <- matrix(rep(c(-3,(rep(0.001,nth-1))),nv),nrow=nth,ncol=nv)
Inc <- mxMatrix( type="Lower", nrow=nth, ncol=nth, free=FALSE, values=1, name="Inc" )
Thre <- mxMatrix( type="Full", nrow=nth, ncol=1, free=TRUE, values=thValues, lbound=thLBound, name="Thre" )
expThreMatOrd <- mxAlgebra( expression= Inc %*% Thre, name="expThre" )
expectationOrd <- mxExpectationNormal( covariance="expCov", means="expMean", threshold="expThre", dimnames=selVars )
FIMLobjectiveOrd <- mxFitFunctionML( )
univSatModel6I <- mxModel("univSat6", obsRawDataOrd, expCovMatOrd, expMeanMatOrd, Inc, Thre, expThreMatOrd, expectationOrd, FIMLobjectiveOrd )
univSatFit6I <- mxRun(univSatModel6I)
summary(univSatFit6I)
# Fit bivariate saturated model to raw data using the path method
nv <- 2
selVars <- c('X','Y')
obsBivData <- mxData(observed=bivData, type="raw" )
expVars <- mxPath(from=c("X", "Y"), arrows=2, free=TRUE, values=1, lbound=.01, labels=c("varX","varY") )
expCovs <- mxPath(from="X", to="Y", arrows=2, free=TRUE, values=.2, lbound=.01, labels="covXY" )
expMeans <- mxPath(from="one", to=c("X", "Y"), arrows=1, free=TRUE, values=.01, labels=c("meanX","meanY") )
bivSatModel1 <- mxModel("bivSat1", manifestVars=selVars, obsBivData, expVars, expCovs, expMeans, type="RAM" )
bivSatFit1 <- mxRun(bivSatModel1)
summary(bivSatFit1)
# Fit bivariate saturated model to raw data using the matrix method in the piecemeal style
expCovM <- mxMatrix(type="Symm", nrow=nv, ncol=nv, free=TRUE, values=c(1,0.5,1), labels=c('V1','Cov','V2'), name="expCov" )
expMeanM <- mxMatrix(type="Full", nrow=1, ncol=nv, free=TRUE, values=c(0,0), labels=c('M1','M2'), name="expMean" )
lowerTriM <- mxMatrix(type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.5, name="Chol" )
expCovMA <- mxAlgebra(expression=Chol %*% t(Chol), name="expCov" )
expCovM
expMeanM
expCovMA
expectationBiv <- mxExpectationNormal( covariance="expCov", means="expMean", dimnames=selVars )
FIMLobjective <- mxFitFunctionML()
bivSatModel2 <- mxModel("bivSat2", obsBivData, lowerTriM, expCovMA, expMeanM, expectationBiv, FIMLobjective )
bivSatFit2 <- mxRun(bivSatModel2)
summary(bivSatFit2)
# Fit bivariate saturated model to raw data using the matrix method in the recursive style
bivSatModel3 <- mxModel("bivSat3", mxData( observed=bivData, type="raw" ) )
bivSatModel3 <- mxModel(bivSatModel3, mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.5, name="Chol" ) )
bivSatModel3 <- mxModel(bivSatModel3, mxAlgebra( expression=Chol %*% t(Chol), name="expCov" ) )
bivSatModel3 <- mxModel(bivSatModel3, mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=c(0,0), labels=c('M1','M2'), name="expMean" ) )
bivSatModel3 <- mxModel(bivSatModel3, mxExpectationNormal( covariance="expCov", means="expMean", dimnames=selVars ), mxFitFunctionML() )
bivSatFit3 <- mxRun(bivSatModel3)
summary(bivSatFit3)
# Fit bivariate saturated model to raw data using the matrix method in the classic style
bivSatModel4 <- mxModel("bivSat4",
mxData(observed=bivData, type="raw"),
mxMatrix(type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.5, name="Chol"),
mxAlgebra(expression=Chol %*% t(Chol), name="expCov"),
mxMatrix(type="Full", nrow=1, ncol=nv, free=TRUE, values=c(0,0), name="expMean"),
mxExpectationNormal(covariance="expCov", means="expMean", dimnames=selVars),
mxFitFunctionML()
)
bivSatFit4 <- mxRun(bivSatModel4)
summary(bivSatFit4)
# Fit bivariate saturated model to raw binary data using the matrix method
expCovMatBin <- mxMatrix(type="Stand", nrow=nv, ncol=nv, free=TRUE, values=.5, name="expCov" )
expMeanMatBin <- mxMatrix(type="Zero", nrow=1, ncol=nv, name="expMean" )
expThreMatBin <- mxMatrix(type="Full", nrow=1, ncol=nv, free=TRUE, values=0, name="expThre" )
expectationBin <- mxExpectationNormal(covariance="expCov", means="expMean", threshold="expThre", dimnames=selVars )
FIMLobjectiveBin <- mxFitFunctionML()
bivSatModel5 <- mxModel("bivSat5", obsBivDataBin, expCovMatBin, expMeanMatBin, expThreMatBin, expectationBin, FIMLobjectiveBin )
bivSatFit5 <- mxRun(bivSatModel5)
summary(bivSatFit5)
# Fit bivariate saturated model to raw ordinal data using the matrix method
expCovMatOrd <- mxMatrix(type="Stand", nrow=nv, ncol=nv, free=TRUE, values=.5, name="expCov" )
expMeanMatOrd <- mxMatrix(type="Zero", nrow=1, ncol=nv, name="expMean" )
expThreMatOrd <- mxMatrix(type="Full", nrow=nth, ncol=nv, free=TRUE, values=c(-1.5,-.5,.5,1.5), name="expThre" )
expectationOrd <- mxExpectationNormal(covariance="expCov", means="expMean", threshold="expThre", dimnames=selVars )
FIMLobjectiveOrd <- mxFitFunctionML()
bivSatModel6 <- mxModel("bivSat6", obsBivDataOrd, expCovMatOrd, expMeanMatOrd, expThreMatOrd, expectationOrd, FIMLobjectiveOrd )
bivSatFit6 <- mxRun(bivSatModel6)
summary(bivSatFit6)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.