# Spline Model with Estimated Knot Point

7 posts / 0 new Offline
Joined: 08/04/2009 - 15:59
Spline Model with Estimated Knot Point

Hi everyone. I'm trying to fit a bi-linear spline model with an estimated knot point (y_nt = g_0n + g_1n(min(0,time-lambda)) + g_2n(max(0,time-lambda). I've pasted my script below, which runs (with warnings), but provides estimates (below as well) that are not reasonable given the min and max functions and the estimate of the knot point (e.g., min of (0 and 1-9.14) is -3.44). Can the min and max functions be used in this manner?

###### ########### Spline Growth Model with Estimated Knot Point

spl.hght.omx <- mxModel("Spline Growth Model, Path Specification",
type="RAM",
mxData(
observed=hght_wide,
type="raw"),
manifestVars=c('hght01','hght03','hght06','hght09','hght12','hght15','hght18','hght24','hght36'),
latentVars=c('g_0n','g_1n','g_2n'),

# Residual Variances

mxPath(
from=c('hght01','hght03','hght06','hght09','hght12','hght15','hght18','hght24','hght36'),
arrows=2,
free=TRUE,
values =10,
labels='v_u'
),

# Latent Variable Covariances

mxPath(
from=c('g_0n','g_1n','g_2n'),
connect='unique.pairs',
arrows=2,
free=TRUE,
values = c(1,0,0,
1,0,
1),
labels=c('v_0','c_10','c_20','v_1','c_21','v_2')
),

mxPath(
from='g_0n',
to=c('hght01','hght03','hght06','hght09','hght12','hght15','hght18','hght24','hght36'),
arrows=1,
free=FALSE,
values=1
),

mxPath(
from='g_1n',
to=c('hght01','hght03','hght06','hght09','hght12','hght15','hght18','hght24','hght36'),
arrows=1,
free=TRUE,
values=0,
labels=c('L11', 'L12', 'L13', 'L14', 'L15', 'L16', 'L17', 'L18', 'L19')
),

mxPath(
from='g_2n',
to=c('hght01','hght03','hght06','hght09','hght12','hght15','hght18','hght24','hght36'),
arrows=1,
free=TRUE,
values=0,
labels=c('L21', 'L22', 'L23', 'L24', 'L25', 'L26', 'L27', 'L28', 'L29')
),

# Latent Variable Means

mxPath(
from='one',
to=c('g_0n','g_1n','g_2n'),
arrows=1,
free=TRUE,
values=c(70, 3, 1),
labels=c('mu_0','mu_1','mu_2')
),

mxMatrix("Full", 1, 1, free = TRUE, values = 7, labels = 'l', name = 'lambda'),

mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L11', name = 'load11'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L12', name = 'load12'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L13', name = 'load13'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L14', name = 'load14'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L15', name = 'load15'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L16', name = 'load16'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L17', name = 'load17'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L18', name = 'load18'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L19', name = 'load19'),

mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L21', name = 'load21'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L22', name = 'load22'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L23', name = 'load23'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L24', name = 'load24'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L25', name = 'load25'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L26', name = 'load26'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L27', name = 'load27'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L28', name = 'load28'),
mxMatrix("Full", 1, 1, free = TRUE, values = 0, labels = 'L29', name = 'load29'),

) # Close Model

#### Output

free parameters:
name matrix row col Estimate Std.Error lbound ubound
1 L11 A hght01 g_1n -3.438915e+00 NA
2 L12 A hght03 g_1n -2.983294e+00 NA
3 L13 A hght06 g_1n -2.299862e+00 NA
4 L14 A hght09 g_1n 5.426485e-02 NA
5 L15 A hght12 g_1n 7.721895e-02 NA
6 L16 A hght15 g_1n 7.721895e-02 NA
7 L17 A hght18 g_1n 7.721895e-02 NA
8 L18 A hght24 g_1n 7.721895e-02 NA
9 L19 A hght36 g_1n 7.721895e-02 NA
10 L21 A hght01 g_2n 7.721895e-02 NA
11 L22 A hght03 g_2n 7.721895e-02 NA
12 L23 A hght06 g_2n 7.721895e-02 NA
13 L24 A hght09 g_2n -1.542440e+00 NA
14 L25 A hght12 g_2n -9.329994e-01 NA
15 L26 A hght15 g_2n -2.495680e-01 NA
16 L27 A hght18 g_2n 4.338635e-01 NA
17 L28 A hght24 g_2n 1.800726e+00 NA
18 L29 A hght36 g_2n 4.534452e+00 NA
19 v_u S hght01 hght01 3.071683e-06 NA
20 v_0 S g_0n g_0n 4.076779e+01 NA
21 c_10 S g_0n g_1n 1.234159e+01 NA
22 v_1 S g_1n g_1n 3.498987e+01 NA
23 c_20 S g_0n g_2n -2.020484e+01 NA
24 c_21 S g_1n g_2n -1.337786e+01 NA
25 v_2 S g_2n g_2n 2.858856e+01 NA
26 mu_0 M 1 g_0n 8.581006e+01 NA
27 mu_1 M 1 g_1n 2.145233e+01 NA
28 mu_2 M 1 g_2n 4.749725e+01 NA
29 l lambda 1 1 9.149271e+00 NA Offline
Joined: 07/31/2009 - 15:24
I think you want to define

I think you want to define load11 to 19 and 21 to 29 as MxAlgebra objects instead of MxMatrix objects. Use the right hand side of the constraints as the expressions for the MxAlgebra objects, and then remove the constraints. Offline
Joined: 08/04/2009 - 15:59
I tried to set it up a couple

I tried to set it up a couple of ways, but have failed to figure out an appropriate solution. I think the trouble is that mxAlgebra doesn't accept labels or equations with an '=' sign. Thus, I cannot set the mxAlgebra equal to the actual loading (or maybe I'm missing something).

That is, this is an acceptable mxAlgebra

but it has a problem with

mxAlgebra(L11 = min(0, 1-lambda), name = 'Load11')

or

mxAlgebra(min(0, 1-lambda), labels = 'L11', name = 'Load11')

Additionally, using mxMatrix to create matrices for the loadings and using common names causes issues as well. In other models with nonlinear constraints I've used the approach initially outlined, but for some reason it seems to be having trouble with this model. Offline
Joined: 07/31/2009 - 15:10
Try: mxAlgebra(min(0,

Try:

mxAlgebra(min(0, 1-lambda), name='Load11'),
...

mxPath(
from='g_1n',
to=c('hght01', 'hght02', ...),
arrows=1,
free=FALSE,
values=0,
),

etc.

That way each algebra is calculated, and the paths are populated with the results. It'll speed things up by removing extra free parameters and extra constraints, as well.

On the other hand, if you have sharable data for this script and wouldn't mind posting it, I'd like to look into this further--this script seems like it should have worked the way you wrote it, and I'd like to see if I can figure out why it didn't. Offline
Joined: 08/04/2009 - 15:59
Thanks Tim. I'll give it a

Thanks Tim. I'll give it a shot today. I've also tried to write more succinct code using vectors.

mxMatrix('Full', 9, 1, free = TRUE, values = 0, labels = c('L11', 'L12', 'L13', 'L14', 'L15', 'L16', 'L17', 'L18', 'L19'), name = 'loadings1'),
mxMatrix('Full', 9, 1, free = TRUE, values = 0, labels = c('L21', 'L22', 'L23', 'L24', 'L25', 'L26', 'L27', 'L28', 'L29'), name = 'loadings2'),
mxMatrix('Full', 9, 1, free = FALSE, values = c(1,3,6,9,12,15,18,24,36), name = 'time'),

but got an error message that I couldn't figure out.

In the meantime, I'm attaching the data. Variables are:

s_id
hght01 hght03 hght06 hght09 hght12
hght15 hght18 hght24 hght36

Any thoughts are greatly appreciated. Offline
Joined: 08/04/2009 - 15:59
Hi everyone. The attached

Hi everyone. The attached code ran properly. Thanks for all of your help. If you are able to figure out what went wrong with the initial code I submitted that would be great. In working on a book for these things, I'm trying to be consistent as possible with model specification. I can easily see how I can program the other nonlinear models in this fashion. I just want to know if I should change previous model specifications to this general approach.
Thanks,
Kevin Offline
Joined: 07/31/2009 - 15:24
Whenever possible it is

Whenever possible it is advisable to eliminate equality constraints and replace them with algebraic computations. If 'x' and 'y' are two free parameters in the model, and 'y' can be written as some function of 'x', then the numerical optimizer will explore a smaller search space if 'y' is not a free parameter but is expressed as an algebra. Numerical optimizers have the additional problem that parameters are expressed as floating point approximations of real numbers, and it may be the case that the built-in epsilon value is set too small to be able to successfully converge on the minimum values in the first script.