You are here

How to get estimates of latent random slopes and random intercept

3 posts / 0 new
Last post
Victor Feagins's picture
Offline
Joined: 05/02/2024 - 15:29
How to get estimates of latent random slopes and random intercept

Hello I fit a random slope and random intercept hlm and I am wondering how to get the latent values of the random slope and intercept. From what I can tell the issue is that the latent variables exist in second level so the function mxFactorScores does not know what to do. In this example U_00 is my random intercept and U_01 is my random slope.

Here is code to simulate HLM data

library(dplyr)
library(tidyr)
library(OpenMx)
 
# Simulate Data -------------
J <- 100 #number of clusters 
n_j <- 30#number of observations 
 
gamma_00 <- 5
gamma_10 <- 3
mu_x <- 4
 
 
sim.dat <- data.frame(J_ID = 1:J) |> 
  mutate(U_0j = rnorm(J, mean = 0 , sd = 1),
         U_1j = rnorm(J, mean = 0, sd = 1),
         gamma_00 = gamma_00,
         gamma_10 = gamma_10,
         n_j = n_j) |> 
  rowwise() |> 
  mutate(X = list(rnorm(n_j, mu_x, sd = 1)),
         R = list(rnorm(n_j, 0, sd = 1))) |> 
  unnest(cols = c(X, R)) |> 
  mutate(Y = gamma_00 + U_0j + (gamma_10 + U_1j) * X + R)
 
 
sim.dat.obs <- sim.dat|> 
  select(J_ID, X, Y) |> 
  mutate(J_ID = as.factor(J_ID))

Here is code for openmx

## J level (Teacher) --------
J_obs <- sim.dat.obs |> 
  select(J_ID) |> 
  distinct()
 
J_cluster_name <- "J_ID"
 
J_latent <- c("U_0", "U_1")
 
J_data <- mxData(J_obs, type = "raw",
                 primaryKey = J_cluster_name)
 
J_var <- mxPath(from = J_latent,
                arrows = 2, values = 1,
                connect = "single")
J_cov <- mxPath(from = J_latent,
                arrows = 2, values = 0,
                connect = "unique.bivariate",
                free = FALSE)
 
J_model <- mxModel("Level_J", type = "RAM",
                   latentVars = J_latent,
                   J_data,
                   J_var,
                   J_cov)
 
## I level (Student) --------
I_obs <- mxData(sim.dat.obs, type = "raw")
 
I_manifest <- c("X", "Y")
 
 
I_means <- mxPath(from = "one",
                  to = c("X"),
                  labels = c("mu_x"),
                  value = 1)
I_var <- mxPath(from = c("X", "Y"),
                labels = c("sd_x", "R_ysd"),
                arrow = 2,
                value = 1)
 
 
 
I_fixed <- mxPath(from = c("one","X"),
                    to = "Y",
                    value = 1,
                    labels = c("gamma_00", "gamma_10"))
 
I_randomslope <- mxPath(from = c("Level_J.U_1"),
                          to = c("Y"),
                          labels = c("data.X"),
                          free = FALSE,
                          joinKey = J_cluster_name)
 
I_randomint <- mxPath(from = c("Level_J.U_0"),
                        to = c("Y"),
                        value = 1,
                        free = FALSE,
                        joinKey = J_cluster_name)
 
I_model<- mxModel("Level_I", type = "RAM",
                  manifestVars = I_manifest,
                  J_model,
                  I_obs,
                  I_means,
                  I_var,
                  I_fixed,
                  I_randomslope,
                  I_randomint)
I.fit <- mxRun(I_model)
#> Running Level_I with 7 parameters
 
summary(I.fit)
#> Summary of Level_I 
#>  
#> free parameters:
#>             name    matrix row col  Estimate  Std.Error A
#> 1       gamma_10         A   Y   X 3.0702452 0.09300885  
#> 2           sd_x         S   X   X 0.9924763 0.02562636  
#> 3          R_ysd         S   Y   Y 1.0691329 0.02857744  
#> 4           mu_x         M   1   X 4.0054015 0.01818861  
#> 5       gamma_00         M   1   Y 5.0011756 0.13307218  
#> 6 Level_J.S[1,1] Level_J.S U_0 U_0 1.1217698 0.24712431 !
#> 7 Level_J.S[2,2] Level_J.S U_1 U_1 0.8267078 0.12055658  
#> 
#> Model Statistics: 
#>                |  Parameters  |  Degrees of Freedom  |  Fit (-2lnL units)
#>        Model:              7                   5993              17907.01
#>    Saturated:             NA                     NA                    NA
#> Independence:             NA                     NA                    NA
#> Number of observations/statistics: 3100/6000
#> 
#> Information Criteria: 
#>       |  df Penalty  |  Parameters Penalty  |  Sample-Size Adjusted
#> AIC:       5921.013               17921.01                 17921.05
#> BIC:     -30271.658               17963.29                 17941.04
#> CFI: NA 
#> TLI: 1   (also known as NNFI) 
#> RMSEA:  0  [95% CI (NA, NA)]
#> Prob(RMSEA <= 0.05): NA
#> To get additional fit indices, see help(mxRefModels)
#> timestamp: 2024-05-03 20:11:30 
#> Wall clock time: 9.318122 secs 
#> optimizer:  SLSQP 
#> OpenMx version number: 2.21.8 
#> Need help?  See help(mxSummary)
 
mxFactorScores(I.fit, "WeightedML")
#> Running Container with 0 parameters
#> Error in `[.data.frame`(got, , paste0(names(coef(fit)), "SE"), drop = FALSE): undefined columns selected
 
mxFactorScores(I.fit@submodels$Level_J, "WeightedML")
#> Error in mxEvalByName(model$expectation$M, model, compute = TRUE): 'name' argument must be a character argument

One work around I tried was to specify dummy latent variables which are the Level_J.U1 times one at the level 1 but the results were not coherent. Observations within a cluster had different values of U_1.

Victor Feagins's picture
Offline
Joined: 05/02/2024 - 15:29
I also asked on Stackoverflow

I also asked on Stackoverflow https://stackoverflow.com/questions/78427513/how-to-get-estimates-of-random-slope-and-intercept-in-openmx

AdminHunter's picture
Offline
Joined: 03/01/2013 - 11:03
Not Implemented

Hi Victor,

Unfortunately, no such feature is currently implemented. In your specific case, it wouldn't be terribly hard for OpenMx to figure out the factor scores, but I don't think the general case (arbitrarily many levels with multiple groups and different factor structures everywhere) is solved. I don't see this feature being implemented anytime soon.

I've created an Issue on GitHub to give a better error message though: https://github.com/OpenMx/OpenMx/issues/391