You are here

Growth Mixture Modeling Question

9 posts / 0 new
Last post
Ryne's picture
Offline
Joined: 07/31/2009 - 15:12
Growth Mixture Modeling Question
AttachmentSize
Binary Data Gompertz Mixture.R8.19 KB

Hi all,

Having some trouble with a growth mixture model. I looked to previous examples in the /models/passing directory (specifically LCAlazarsfeld.R, Acemix.R and Acemix2.R) and thought the code below would work. To summarize, I created two models for the two classes (differing in a small set of parameters a, r and d), and put them both in a parent model that contains the data and an objective function that sums the (vectorized) objective functions from the two child models. When run, the code below fixes the class probabilities at .5 for each class and yields the same -2LL as the single-class model. Parameters stayed fixed at their starting values for the class-dependent parameters, with standard errors at zero.

My conceptual problem with the code is with the class probability variable as a scalar value, when I would think it would be individually varying. Any help you can give would be appreciated.

tbrick's picture
Offline
Joined: 07/31/2009 - 15:10
Your objective function needs

Your objective function needs to be changed.

The code you have:
 mxAlgebra(-2*sum(log(pclass1%x%Gompertz1.objective + pclass2%x%Gompertz1.objective)), name="mixtureObj")
adds the Gompertz1 objective to itself. Change the second one to Gompertz2, and it should work.

Ryne's picture
Offline
Joined: 07/31/2009 - 15:12
PEBKAC error. Everything

PEBKAC error. Everything comes up roses when I try to mix two different classes instead of two of the same class. Pending my permission to post the data, I'll wrap this up into something for the models/passing directory and the documentation.

Mariya's picture
Offline
Joined: 07/21/2010 - 20:29
Hi all, I am following Ryne's

Hi all,
I am following Ryne's code on fitting GMM, but ran into a few error messages. My model is simpler: I fit a Latent Basis 2-class GMM with 9 repeated assessments, constraining time loadings at times 1 & 4. I am starting with a model that has no within-class intercept & slope variances (so, my 'phi' matrix is constrained to 0).

I get several error messages.
First, when I am creating Model Expectations for Mixture (Parent) Model by specifying the sum of class proportions equal to 1, I get the following error:

> classP <- mxMatrix ("Full", 2, 1, free = TRUE, lbound = 0,
+ labels = c ("pclass1", "pclass2"), name = "classProbs")
>
> classS <- mxAlgebra (sum(classProbs), name = "classSum")
> constant <- mxMatrix ("Iden", 1, name = "con")
> classC <- mxConstraint ("classSum", "=", "con", name = "classCon")
Error in mxConstraint("classSum", "=", "con", name = "classCon") :
unused argument(s) ("=", "con")

Of course, subsequently, I can't specify the MxModel, since I don't have the classC object, and I get an error when running the full mixture model.

I am attaching the complete R syntax file and would greatly appreciate any leads on how to modify it.

Mariya

mspiegel's picture
Offline
Joined: 07/31/2009 - 15:24
The mxConstraint() syntax has

The mxConstraint() syntax has changed. Try using:

mxConstraint(classSum == con, name = "classCon")

Note the absence of quotes in the first argument. MxConstraints and now specified much like MxAlgebras are specified.

Mariya's picture
Offline
Joined: 07/21/2010 - 20:29
Thank you! That solved the

Thank you! That solved the issure of constrains.
However, when I run the final model with
> mixedResults <- mxRun (mixedModel)
I am getting the following error message:

Error: The job for model 'Group Based GMM 2cl' exited abnormally with the error message: Expected covariance matrix is not positive-definite in data row 30 at iteration 266.
In addition: Warning message:
In model 'Group Based GMM 2cl' NPSOL returned a non-zero status code 6. The model does not satisfy the first-order optimality conditions to the required accuracy, and no improved point for the merit function could be found during the final linesearch (Mx status RED)

neale's picture
Offline
Joined: 07/31/2009 - 15:14
I would be tempted to keep

I would be tempted to keep theta elements positive, thusly:

replace

manCov <- mxMatrix ("Diag", 9, 9, free = TRUE,
values = rep (25, times = 9),
labels = 'residual',
name = 'theta')
with

manCov <- mxMatrix ("Diag", 9, 9, free = TRUE,
values = rep (25, times = 9),
labels = 'residual',
lbound = 0,
name = 'theta')

And I wonder if this has anything to do with a "mancave" :)

Mariya's picture
Offline
Joined: 07/21/2010 - 20:29
Thank you,

Thank you, Michael.
specifying lower bounds on the theta matrix did help the estimation process. What I realized, however, is that factor loadings for the slope parameter were equal across latent classes, thus, representation for one of classes was 0% (makes sense, since trajectories were specified as being the same).

When specifying class-specific loadings for the slope (so that trajectories would differ across classes), I run into a new problem of R complaining:

#

Error: Unknown reference 'lambdaS1' detected in the entity 'lambda1' in model 'groupGMM1'

#

Here is the modified syntax I'm using to which, I think, the error message is referring:

factor loadings for slope cl1

time1 <- mxMatrix ("Full", 9, 1, free = c(FALSE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, TRUE, TRUE),
values = c(0, .5, .8, 1, 1, 1, 1, 1, 1),
label = paste ('t', 1:9, 'cl1', sep = ""),
name = 'lambdaS1')

factor loadings for slope cl2

time2 <- mxMatrix ("Full", 9, 1, free = c(FALSE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, TRUE, TRUE),
values = c(0, .5, .8, 1, .8, .7, .6, .5, .4),
label = paste ('t', 1:9, 'cl2', sep = ""),
name = 'lambdaS2')

factor loading matrix lambda

loadings1 <- mxAlgebra (cbind (lambdaInt, lambdaS1), name = 'lambda1')
loadings2 <- mxAlgebra (cbind (lambdaInt, lambdaS2), name = 'lambda2')

model expectations:

meanAlg1 <- mxAlgebra (alpha1 %% t(lambda1), dimnames = list(NULL, names),
name = 'mean1')
covAlg1 <- mxAlgebra (lambda1 %
% phi %*% t(lambda1) + theta,
dimnames = list(names, names),
name = 'cov1')

meanAlg2 <- mxAlgebra (alpha2 %% t(lambda2), dimnames = list(NULL, names),
name = 'mean2')
covAlg2 <- mxAlgebra (lambda2 %
% phi %*% t(lambda2) + theta,
dimnames = list(names, names),
name = 'cov2')

model1 <- mxModel ("groupGMM1",
factorMeans1, factorCov, manCov,
unit, time, loadings1,
meanAlg1, covAlg1,
mxFIMLObjective ('cov1', 'mean1'))
etc.

Steve's picture
Offline
Joined: 07/30/2009 - 14:03
time1 is not included in your

time1 is not included in your model. That should do the trick.