mxComputeEM: Fit a model using DLR's (1977) Expectation-Maximization (EM)...

In OpenMx: Extended Structural Equation Modelling

Fit a model using DLR's (1977) Expectation-Maximization (EM) algorithm

Description

The EM algorithm constitutes the following steps: Start with an initial parameter vector. Predict the missing data to form a completed data model. Optimize the completed data model to obtain a new parameter vector. Repeat these steps until convergence criteria are met.

Usage

mxComputeEM(
  expectation = NULL,
  predict = NA_character_,
  mstep,
  observedFit = "fitfunction",
  ...,
  maxIter = 500L,
  tolerance = 1e-09,
  verbose = 0L,
  freeSet = NA_character_,
  accel = "varadhan2008",
  information = NA_character_,
  infoArgs = list(),
  estep = NULL
)

Arguments

expectation	a vector of expectation names \lifecycledeprecated
predict	what to predict from the observed data \lifecycledeprecated
mstep	a compute plan to optimize the completed data model
observedFit	the name of the observed data fit function (defaults to "fitfunction")
...	Not used. Forces remaining arguments to be specified by name.
maxIter	maximum number of iterations
tolerance	optimization is considered converged when the maximum relative change in fit is less than tolerance
verbose	integer. Level of run-time diagnostic output. Set to zero to disable
freeSet	names of matrices containing free variables
accel	name of acceleration method ("varadhan2008" or "ramsay1975")
information	name of information matrix approximation method
infoArgs	arguments to control the information matrix method
estep	a compute plan to perform the expectation step

Details

The arguments to this function have evolved. The old style mxComputeEM(e,p,mstep=m) is equivalent to the new style mxComputeEM(estep=mxComputeOnce(e,p), mstep=m). This change allows the API to more closely match the literature on the E-M method. You might use mxAlgebra(..., recompute='onDemand') to contain the results of the E-step and then cause this algebra to be recomputed using mxComputeOnce.

This compute plan does not work with any and all expectations. It requires a special kind of expectation that can predict its missing data to create a completed data model.

The EM algorithm does not produce a parameter covariance matrix for standard errors. The Oakes (1999) direct method and S-EM, an implementation of Meng & Rubin (1991), are included.

Ramsay (1975) was recommended in Bock, Gibbons, & Muraki (1988).

References

Bock, R. D., Gibbons, R., & Muraki, E. (1988). Full-information item factor analysis. Applied Psychological Measurement, 6(4), 431-444.

Dempster, A. P., Laird, N. M., & Rubin, D. B. (1977). Maximum likelihood from incomplete data via the EM algorithm. Journal of the Royal Statistical Society. Series B (Methodological), 1-38.

Meng, X.-L. & Rubin, D. B. (1991). Using EM to obtain asymptotic variance-covariance matrices: The SEM algorithm. Journal of the American Statistical Association, 86 (416), 899-909.

Oakes, D. (1999). Direct calculation of the information matrix via the EM algorithm. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 61(2), 479-482.

Ramsay, J. O. (1975). Solving implicit equations in psychometric data analysis. Psychometrika, 40 (3), 337-360.

Varadhan, R. & Roland, C. (2008). Simple and globally convergent methods for accelerating the convergence of any EM algorithm. Scandinavian Journal of Statistics, 35, 335-353.

See Also

MxAlgebra, mxComputeOnce

Examples

library(OpenMx)
set.seed(190127)

N <- 200
x <- matrix(c(rnorm(N/2,0,1),
              rnorm(N/2,3,1)),ncol=1,dimnames=list(NULL,"x"))
data4mx <- mxData(observed=x,type="raw")

class1 <- mxModel("Class1",
	mxMatrix(type="Full",nrow=1,ncol=1,free=TRUE,values=0,name="Mu"),
	mxMatrix(type="Full",nrow=1,ncol=1,free=TRUE,values=4,name="Sigma"),
	mxExpectationNormal(covariance="Sigma",means="Mu",dimnames="x"),
	mxFitFunctionML(vector=TRUE))

class2 <- mxRename(class1, "Class2")

mm <- mxModel(
	"Mixture", data4mx, class1, class2,
	mxAlgebra((1-Posteriors) * Class1.fitfunction, name="PL1"),
	mxAlgebra(Posteriors * Class2.fitfunction, name="PL2"),
	mxAlgebra(PL1 + PL2, name="PL"),
	mxAlgebra(PL2 / PL,  recompute='onDemand',
	          initial=matrix(runif(N,.4,.6), nrow=N, ncol = 1), name="Posteriors"),
	mxAlgebra(-2*sum(log(PL)), name="FF"),
	mxFitFunctionAlgebra(algebra="FF"),
	mxComputeEM(
	  estep=mxComputeOnce("Mixture.Posteriors"),
	  mstep=mxComputeGradientDescent(fitfunction="Mixture.fitfunction")))

mm <- mxOption(mm, "Max minutes", 1/20)  # remove this line to find optimum
mmfit <- mxRun(mm)
summary(mmfit)

OpenMx documentation built on Oct. 19, 2024, 9:06 a.m.