You are here

Example Script for LDE estimation

8 posts / 0 new
Last post
Steve's picture
Offline
Joined: 07/30/2009 - 14:03
Example Script for LDE estimation
AttachmentSize
Plain text icon demoOscillator.txt1.03 KB
Binary Data LDE_2ndOrder_OneIndicator.R3.52 KB

Here is a very basic script for a univariate second order differential equation that Pascal Deboeck and I wrote yesterday. This is the infamous damped linear oscillator.

In order to give you some sample data to run it on, I have attached a data file and used a read.table file that accesses the data file in this post using its URL.

If you have comments, I'd like to hear them!

amullings09's picture
Offline
Joined: 12/21/2014 - 17:45
Help with LDE function

Hey,

I am a new graduate student and am interested in using this model for my thesis. I am trying to practice with your demo oscillator data, and then plan on testing it with some of my own data. So far the function provided has worked, except it seems that one part of the function is outdated:

Warning message:
In mxMLObjective("R") :
Objective functions like mxMLObjective() have been deprecated in favor of expectation and fit functions.
Please use mxExpectationNormal(covariance= , ...) instead, and add a call to mxFitFunctionML()
See examples at help(mxExpectationNormal)

I was wondering whether someone could direct/guide me on what fit function would be best to use for the LDE estimator? I am new to SEM (e.g., First year MA candidate), so I am at a bit of a learning curve working with these functions in R.

mhunter's picture
Offline
Joined: 07/31/2009 - 15:26
Change script or ignore warning

You can ignore the warning if you want. The model still works fine. If the warning message bothers you, then you can replace the line

    mxMLObjective("R"),

with the pair of lines

    mxExpectationNormal("R"),
    mxFitFunctionML(),

The single line above and the pair of lines are exactly equivalent in terms of the modeling procedures performed by OpenMx. You can verify this for yourself by running the model both ways and comparing the results.

We made the expectation-fit split in 2.0 because we realized that things like the mxMLObjective were really making two decisions: one, the model is for a covariance matrix "R" with a normal distribution (expectation normal), and two, we're going to fit the normal distribution with maximum likelihood (fit function ML).

Hayafh's picture
Offline
Joined: 12/01/2018 - 12:49
Warnings and no Standard errors

Hi,
I tried running the LDE model using the code and the dataset given. However, I get the following warnings.
Warning message:
In mxMLObjective("R") :
Objective functions like mxMLObjective() have been deprecated in favor of expectation and fit functions.
Please use mxExpectationNormal(covariance= , ...) instead, and add a call to mxFitFunctionML()
See examples at help(mxExpectationNormal)
> ldeModel1Fit <- mxRun(ldeModel1)
Running LDE_Model_1 with 11 parameters
Warning message:
In model 'LDE_Model_1' Optimizer returned a non-zero status code 5. The Hessian at the solution does not appear to be convex. See ?mxCheckIdentification for possible diagnosis (Mx status RED).

With these warnings, I get extreme values for eta and zeta (I’m assuming that is what the estimates for matrix A are): -8 and -70, and there are NA’s for standard errors.

I have tried mxCheckIdentification(ldeModel1Fit) and got the following:
Model is locally identified
$status
[1] TRUE

$non_identified_parameters
[1] "LDE_Model_1.A[3,1]" "LDE_Model_1.A[3,2]" "LDE_Model_1.S[1,1]" "LDE_Model_1.S[1,2]" "LDE_Model_1.S[2,2]"
[6] "LDE_Model_1.S[3,3]" "LDE_Model_1.U[3,3]" "LDE_Model_1.U[4,4]" "LDE_Model_1.U[5,5]"

AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Not sure what your version of

Not sure what your version of OpenMx is (?), but the example works fine for me once I changed the url for the data from

demoOscillator <- read.table("http://openmx.psyc.virginia.edu/sites/default/files/demoOscillator.txt", col.names="x")

to

demoOscillator <- read.table("http://openmx.ssri.psu.edu/sites/default/files/demoOscillator.txt", col.names="x")

:

> mxCheckIdentification(ldeModel1Fit)
Model is locally identified
$status
[1] TRUE

$jacobian
LDE_Model_1.A[3,1] LDE_Model_1.A[3,2] LDE_Model_1.S[1,1] LDE_Model_1.S[1,2] LDE_Model_1.S[2,2] LDE_Model_1.S[3,3] LDE_Model_1.U[1,1] LDE_Model_1.U[2,2]
cov1_1 1.054973e-01 -0.0083329529 1.322023938 45.5491149 392.338181 0.164025 1 0
...snip...

$non_identified_parameters
[1] "None"

mhunter's picture
Offline
Joined: 07/31/2009 - 15:26
Bounds

For the warning about mxMLObjective see the previous solution provided.

For the Mx status RED, I'd recommend setting a lower bound of about -10 or so on eta and zeta.

tbates's picture
Offline
Joined: 07/31/2009 - 14:25
Error in eigen

Updated for modern fit function and new URL, and adding mxTryhard for code Red, I get this error

ldeModel1Fit <- mxTryHard(ldeModel1)
Fit attempt 0, fit=-1521.67613260047, new current best! (was 201.317056439741)
Error in eigen(fit$output$calculatedHessian, symmetric = T, only.values = T) : 
  infinite or missing values in 'x'

Guessing that the error should be suppressed in mxTryHard?

OpenMx version: 2.11.5 [GIT v2.11.5]
R version: R version 3.5.1 (2018-07-02)
MacOS: 10.14.2

File attachments: 
AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
I reproduce the behavior you

I reproduce the behavior you report.

Guessing that the error should be suppressed in mxTryHard?

Arguably, yes. But historically, I don't believe it ever has been. It probably only looks strange when you're running in an interactive session (in which case, nowadays, argument silent to mxTryHard() defaults to TRUE).