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.
# ---------------------------------------------------------------------
# Program: UnivariateRaw-OpenMx100214.R
# Author: Steven M. Boker
# Date: Sun Feb 14 13:18:47 EST 2010
#
# This program fits a FIML univariate model to the
# multiData simulated data.
#
#
# ---------------------------------------------------------------------
# Revision History
# -- Sun Feb 14 13:18:51 EST 2010
# Created UnivariateRaw-OpenMx100214.R.
#
# ---------------------------------------------------------------------
# ----------------------------------
# Read libraries and set options.
require(OpenMx)
# ----------------------------------
# Read the data and print descriptive statistics.
data(multiData1)
# ----------------------------------
# Build an OpenMx univariate regression model using y and x1
manifests <- c("x1", "y")
uniRegModelRaw <- mxModel("ur",
type="RAM",
manifestVars=manifests,
mxPath(from="x1", to="y", arrows=1,
free=TRUE, values=.2, labels="b1"),
mxPath(from=manifests,
arrows=2, free=TRUE, values=.8,
labels=c("VarX1", "VarE")),
mxPath(from="one", to=manifests,
arrows=1, free=TRUE, values=.1,
labels=c("MeanX1", "MeanY")),
mxData(observed=multiData1, type="raw")
)
uniRegModelRawOut <- mxRun(uniRegModelRaw)
summary(uniRegModelRawOut)
#---------------------
# check values: uniRegModelRawOut
expectVal <- c(0.669179, 1.13643, 1.647629, 0.984894, 3.189368)
expectSE <-c(0.053849, 0.071873, 0.104204, 0.047674, 0.078154)
expectMin <- 3151.492
omxCheckCloseEnough(expectVal, uniRegModelRawOut$output$estimate, 0.001)
omxCheckCloseEnough(expectSE,
as.vector(uniRegModelRawOut$output[['standardErrors']]), 0.001)
omxCheckCloseEnough(expectMin, uniRegModelRawOut$output$minimum, 0.001)
# -----------------
uniRegModelRaw$fitfunction$rowDiagnostics <- TRUE
urm1 <- mxRun(uniRegModelRaw)
lk <- attr(urm1$fitfunction$result, 'likelihoods')
omxCheckTrue(all(attr(urm1$fitfunction$result, "rowObs") == 2))
omxCheckCloseEnough(urm1$output$fit, expectMin, 1e-3)
omxCheckCloseEnough(-2 * sum(log(lk)), expectMin, 1e-3)
m <- attr(mxEval(fitfunction, urm1), "rowDist")
qqplot(d <- qchisq(ppoints(nrow(urm1$data$observed)), df = 2), m,
main = expression("Q-Q plot for" ~~ {chi^2}[nu == 2]))
qqline(m, distribution = function(p) qchisq(p, df = 2),
prob = c(0.1, 0.6), col = 2)
mc <- hist(m, breaks=c(0:10, 15), plot=FALSE)$counts
dc <- hist(d, breaks=c(0:10, 15), plot=FALSE)$counts
obsChi <- sum((mc - dc)**2/dc)
omxCheckTrue(pchisq(obsChi, df=length(mc)-1, lower.tail=FALSE) > .65)
# The goodness of fit chi-square probability on the residuals compared to the expected distribution of residuals should be high (e.g. more than .65).
algModel <- mxModel("wrap",
mxModel(uniRegModelRaw, mxFitFunctionML(vector=TRUE, rowDiagnostics=TRUE)),
mxAlgebra(-2*sum(log(ur.fitfunction)), 'm2ll'),
mxFitFunctionAlgebra('m2ll'))
urm2 <- mxRun(algModel)
omxCheckCloseEnough(urm2$output$fit, expectMin, 1e-3)
omxCheckEquals(max(abs(urm2$ur$fitfunction$info$likelihoods -
urm2$ur$fitfunction$result)),0)
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.