Estimating a non-linear model in OpenMx

Posted on
No user picture. pgok Joined: 11/13/2019
Hi, I'm pretty new to OpenMx, having only used it off-and-on over the past few years. That said, I am pretty interested in using it more, especially given the multilevel modeling capabilities (e.g., as discussed in Pritkin et al., 2017). However, my question isn't about the multilevel modeling capabilities, my question is about fitting non-linear models.

I am looking at trying to fit an exponential decay/growth model. It would look something like this: y~b2+(b0-b2)\*exp(b1\*time). Where b1 is constrained to be negative (so that everything after b2 eventually drops out of the model as time goes to infinity). The problem I'm running into is that it isn't clear to me how to include a parameter in the exponential part of the model. I ran across [this thread](https://openmx.ssri.psu.edu/thread/735) on including non-linear constraints, but I can't seem to figure out how to adapt that code to my issue. A (not working) minimal example is below:


time=seq(0, 100, by=0.01)
y<-1+(5-1)*exp(-.1*time)+rnorm(length(time),sd=.2)

data_list<-data.frame(
time=time,
y=y
)

plot(time,y)

modelExp<-mxModel("exp", type="RAM", manifestVars=c(names(data_list)),
latentVars = "b1",
mxPath(from="b1", to="b1", free=F, value=0, arrows=2),
mxPath(from="one", to="b1", value=-.1),
mxAlgebra(exp(b1*data.time), name="expPart"),
mxPath(from="one", to="y", label="b2"),
mxPath(from="time", to="y", label="expPart"),
mxPath(from="one", to="time"),
mxData(data_list, type="raw")
)
fit<-mxRun(modelExp)

Part of the problem is that I can't figure out how to reference the parameter b1 (which I am using a zero variance latent variable for) in the mxAlgebra. The other part of the problem is, based on the previous thread I mentioned, I don't see a way to estimate the (b0-b2) part. I don't think I actually need to estimate those as two separate parameters (I already estimate b2 as a kind of intercept, so I can just estimate b0-b2 as some single parameter b* and then do the algebra separately), but it looks like the model with try to make "expPart" equal to b0-b2.

It's entirely possible that I am either trying to do something entirely unreasonable, or something trivially easy in an overly complicated way, but regardless it escapes me how I would do this and any help would be greatly appreciated.

Replied on Fri, 11/12/2021 - 10:37
Picture of user. mhunter Joined: 07/31/2009

Hi!

Here's how I would do the model you proposed.


modelExp<-mxModel("exp", type="RAM", manifestVars='y',
mxMatrix('Full', nrow=1, ncol=3,
labels=c('b0', 'b1', 'b2'), name='params',
free=TRUE, values=c(.5, -1, 1)),
mxAlgebra(b2+(b0-b2)*exp(b1*data.time), name='predY'),
mxPath(from="one", to="y", labels='predY[1,1]', free=FALSE),
mxPath(from="y", arrows=2, labels='resid', values=1),
mxData(data_list, type="raw")
)
fit<-mxRun(modelExp)

It seems to do okay.


> summary(fit)
Summary of exp

free parameters:
name matrix row col Estimate Std.Error A
1 resid S y y 0.04087680 0.0005780564
2 b0 params 1 1 5.00507642 0.0130146711
3 b1 params 1 2 -0.09956287 0.0005180620
4 b2 params 1 3 0.99726107 0.0026133824

Replied on Mon, 11/15/2021 - 13:20
No user picture. pgok Joined: 11/13/2019

In reply to by mhunter

Follow up question: is it possible to include a latent variable in the mxAlgebra expression? For example, if I had some latent variable "L" that I wanted to include in an arbitrary expression [say, y=b1*log(b2*L)], I am guessing that there is a way to do that, but it doesn't seem to have a clear analog to the manifest variable case [i.e., there isn't a way to do b1*log(b2*latent.L)].

Replied on Thu, 11/18/2021 - 09:19
Picture of user. AdminNeale Joined: 03/01/2013

In reply to by pgok

A latent variable doesn’t have a value per se. it is a variable with a distribution. One possibility would be to create a mixture distribution of say 12 latent trait values weighted by Gaussian quadrature weights. This is currently one way of handling such things, but it doesn’t scale well as with 2 latent moderators it would take 144 evaluations to get the likelihood of a data vector, and generally with m quadrature points and d latent moderators, m^d gets large fast. There may be an alternative in OpenMx 3.0, but it’s not ready yet. Some time in 2022 looks likely though. The other route would be via MCMC, but there appear to be many pitfalls and limitations in that area.
Replied on Thu, 11/18/2021 - 18:44
No user picture. pgok Joined: 11/13/2019

In reply to by AdminNeale

Thanks for the response. I can't say I'm terribly surprised that it gets radically more complex for the latent variable case, but I figured I'd ask anyways. I'll be interested to see what OpenMx 3.0 has in store for this though. Thanks again!