gsem not estimating

[Nina Schnyder]

I'm trying to run a SEM model in Stata.
The N of dataset = 1,500.
The outcome (for the final model) is binary (ever_treated), and latent predictors consists of likert scale responses to several questions.
There are some missing values, both on outcome and predictors (not exceeding 2%).


At next step I'd like to estimate the final `gsem` model as follow:

gsem (Psychosocial -> fb2s fb2d fb2l fb2p fb2q, ologit) ///
(Bio -> fb2n fb2c, ologit) (Narcotics -> fb2j fb2m, ologit) ///
(Socialization -> fb2b fb2a fb2f fb2r, ologit) ///
(Own -> fb2k fb2g fb2h, ologit) ///
(Alternative -> new_fb5j new_fb8b new_fb8c new_fb8d new_fb8f, ologit) ///
(Established -> new_fb5d new_fb5e new_fb5h new_fb8a new_b8e, ologit) ///
(Non_prof -> new_fb5a new_fb5c new_fb5f new_fb5g, ologit) ///
(Other_medical -> new_fb5b new_fb8g, ologit) ///
(WSD -> fb15a-fb15g, ologit) ///
(HelpA -> fb16 fb17 fb18 cat_fb12 new_fb19, ologit) ///
(PublicS -> nfb21_1 nfb21_2 nfb21_3 nfb21_4 nfb21_8 nfb21_9, ologit) ///
(PublicHelpA -> fb22_4 fb22_5 fb22_6 fb22_7 nfb22_9 fb22_10, ologit) ///
(Established <- Psychosocial Socialization Bio Own Narcotics) ///
(Alternative <- Psychosocial Socialization Bio Own Narcotics) ///
(Non_prof <- Psychosocial Socialization Bio Own Narcotics) ///
(Other_medical <- Psychosocial Socialization Bio Own Narcotics) ///
(HelpA <- PublicHelpA PublicS Psychosocial Own Socialization Bio Narcotics) ///
(WSD <- PublicHelpA PublicS Psychosocial Own Socialization Bio Narcotics), ///
(ever_treated <- Established Alternative Non_prof Other_medical WSD HelpA, logit) ///
latent(Psychosocial Own Socialization Bio Narcotics Alternative Established Non_prof Other_medical WSD HelpA PublicS PublicHelpA)

In order to speed up the computations I tired Stata's solution of gradually increasing complexity of the path model for each single latent variable and storing estimates in matrices which are then reused (as per Stata manual recomendation). So for example 'Alternative' 'Established' latent variables are constructed (realtively fast) in this way:

* Treatment recommendations
* Alternative and Established (matrix b1)
gsem  (Alternative -> new_fb5j new_fb8b new_fb8c new_fb8d new_fb8f, ologit) ///
(Established -> new_fb5d new_fb5e new_fb5h new_fb8a new_b8e, ologit), ///
latent (Alternative Established)
matrix b1 = e(b)

* Other Medical and Non Professional (matrix b2)
gsem  (Non_prof -> new_fb5a new_fb5c new_fb5f new_fb5g, ologit) ///
(Other_medical -> new_fb5b new_fb8g, ologit), ///
latent (Non_prof Other_medical)
matrix b2 = e(b)

* ALL treatment recommendations
gsem  (Alternative -> new_fb5j new_fb8b new_fb8c new_fb8d new_fb8f, ologit) ///
(Established -> new_fb5d new_fb5e new_fb5h new_fb8a new_b8e, ologit) ///
(Non_prof -> new_fb5a new_fb5c new_fb5f new_fb5g, ologit) ///
(Other_medical -> new_fb5b new_fb8g, ologit), ///
latent (Alternative Established Non_prof Other_medical) from(b1 b2) 

est save "SEM_Treatment Recommendation", replace

matrix b = e(b)

In a similar vein I also constructed matrices a and c containing other latent variables. At next step I tried to run the following model and after two days of work it doesn't give any results whatsoever:

gsem  (Psychosocial -> fb2s fb2d fb2l fb2p fb2q, ologit) (Bio -> fb2n fb2c, ologit) ///
(Narcotics -> fb2j fb2m, ologit) ///
(Socialization -> fb2b fb2a fb2f fb2r, ologit) ///
(Own -> fb2k fb2g fb2h, ologit) ///
(Alternative -> new_fb5j new_fb8b new_fb8c new_fb8d new_fb8f, ologit) ///
(Established -> new_fb5d new_fb5e new_fb5h new_fb8a new_b8e, ologit) ///
(Non_prof -> new_fb5a new_fb5c new_fb5f new_fb5g, ologit) ///
(Other_medical -> new_fb5b new_fb8g, ologit) ///
(WSD -> fb15a-fb15g, ologit) ///
(HelpA -> fb16 fb17 fb18 cat_fb12 new_fb19, ologit) ///
(PublicS -> nfb21_1 nfb21_2 nfb21_3 nfb21_4 nfb21_8 nfb21_9, ologit) ///
(PublicHelpA -> fb22_4 fb22_5 fb22_6 fb22_7 nfb22_9 fb22_10, ologit) ///
(Established <- Psychosocial Socialization Bio Own Narcotics) ///
(Alternative <- Psychosocial Socialization Bio Own Narcotics) ///
(Non_prof <- Psychosocial Socialization Bio Own Narcotics) ///
(Other_medical <- Psychosocial Socialization Bio Own Narcotics) ///
(HelpA <- PublicHelpA PublicS Psychosocial Own Socialization Bio Narcotics) ///
(WSD <- PublicHelpA PublicS Psychosocial Own Socialization Bio Narcotics), ///
latent(Psychosocial Bio Narcotics Socialization Own Alternative Established Non_prof Other_medical WSD HelpA PublicS PublicHelpA) ///
from(a b c)

Of note, I tried running equivalent `sem` model (using `mlmv` method) and it estimates results without any issues (took few hours though).

Is there any mistake in the code causing that?
Or is the complexity of the model too much for Stata?

[Richard Williams]

Are you sure it is not estimating as opposed to taking forever to estimate? Does the iteration log give you any clues?

Sometimes adding the difficult option will work miracles. Or, you might see if it will run with a subsample of, say, 300 cases.

It wouldn't surprise me if this is too difficult for Stata. Of course it might be too difficult for any software. Do you have access to any other programs, like MPlus, that you could try?

[Marcos Almeida]

Hello Nina,

Welcome to the Stata Forum.

I hope you will get more optimistic replies. Personally, I gather the problem related to the estimations it is not due to "complexity of the model too much for Stata", but complexity in modeling itself.

Edited, to make it clear: it seems there are many variables (latent, inclusive) and many paths, yet the N is less than 1500 (due to missing data). Shall we get an output, I wonder how large would that be, and how tough would it be to become "comprehensive" for the audience.

[Weiwen Ng]

I hope Nina doesn't mind if I piggyback on the discussion, as I tried something like this myself, but with only one level of latent variables. I see in the original question there appear to be 2 levels of latent variables, and that moreover, Nina is then using the second level latent variables to go predict an observed one in a logit model). So, this is a relatively complex model.

In particular, I see that the first level latent variables are composed of several Likert items each. Is it that far wrong to just estimate them using SEM, thus treating them as continuous? I am pretty sure I see people do this all the time. Yes, I'm well aware that strictly speaking, Nina has ordinal data (and she is treating the data as such, to her great credit). But, again, people do this all the time in health services research, and moreover, most people would not bother fitting a complicated SEM like this one. They'd just run a factor analysis (confirmatory or exploratory), probably use regular SEM, then just go sum the scores for each domain and go run a logit on the DV.

For a regular SEM, even if a total of 4% of observations get discarded due to missing items, I see nothing wrong with trying a first run without the MLMV option, which (I tried this once) really adds computation time. In health services research, 4% is considered an acceptable level of missingness (and Nina probably has fewer than 4% missing total, going from what she said). Also, Nina could use the asymptotic distribution free estimation method in SEM, which I have been told (probably by either WBuchanan or Clyde Schechter) is a better estimation method for categorical data. Now, any estimation method in SEM, save MLMV, will drop missing cases listwise, but again, 4% is quite acceptable in my field. Am I off base here? I would like to know.

Moreover, Nina is trying to use the final level latent variables to predict an observed variable. Would it be acceptable to leave that last set of paths (i.e. from the final level latent variables to the observed treatment decision) off the GSEM model, and estimate that? Now, it does look like Nina tried to do this, and it didn't help with convergence, but I'm asking in general if this would be acceptable. I can see the problem in just going to predict each latent factor separately from the others, as that may not properly estimate the covariance between the latent factors.

I would normally have asked if the model could be unidentified, but if it converged in SEM, it's identified.

[Weiwen Ng]

I was unable to edit my prior post, so I'm adding two things I just noticed.


Maybe sheer complexity is a genuine problem. I looked through Nina's code again, and I think I count about 30 Likert items and 12 level 1 latent variables (each of which appear to have unique Likert items loading on them). But then, if I am reading right, for 6 of those same level 1 latent variables, Nina is also saying that the other 6 level 1 latent variables are determinants of them. Is that correct? That sounds like a recursive path, or whatever the terminology is, and now I'm way out of my depth. I would ask if that was theoretically justified in the first place. Also, I have no idea how SEM or GSEM deal with that recursivity (although apparently SEM was able to solve the model...).

Back to problem solving on the original proposed GSEM model. I see the SEM manual says that in some cases, changing the integration method to mode and curvature adaptive Gauss Hermite, i.e.

..., method(mcaghermite)

can work in multilevel models. Specifically, "It is accurate but not as fast. If you are fitting a multilevel model, there are cases where method 1 [that is, the default mvaghermite method] will not work but this method will work."

If not for that, the other integration methods (these being ghermite and laplace) may work; they are less accurate, but they can then be used as starting values. The manual also says that using fewer integration points (default is 7) to find a starting solution is a possible solution. I wonder if these solutions might help! I will note, though, that the manual says these are potential solutions to convergence problems, e.g. initial values not feasible or an infinite iteration log with a lot of (not concave) notations. It's not clear that this is the problem; it may be that the model is simply too complex to converge as a GSEM.

http://www.stata.com/manuals13/semin...pdf#semintro12

[Nina Schnyder]

Are you sure it is not estimating as opposed to taking forever to estimate? Does the iteration log give you any clues?

Sometimes adding the difficult option will work miracles. Or, you might see if it will run with a subsample of, say, 300 cases.

It wouldn't surprise me if this is too difficult for Stata. Of course it might be too difficult for any software. Do you have access to any other programs, like MPlus, that you could try?

Dear Richard

thank you so much for your time and your suggestions. I think there's not really a problem with estimating per se (as the model run with the sem command). But the iteration log doesn't give me any clues unfortunately, it didn't move even though we run the model for days now.

I run parts of the model in R with the lavaan package and it worked much faster than in stata. Maybe it will be possible to estimate the full model in R.

[Nina Schnyder]

Hello Nina,

Welcome to the Stata Forum.

I hope you will get more optimistic replies. Personally, I gather the problem related to the estimations it is not due to "complexity of the model too much for Stata", but complexity in modeling itself.

Edited, to make it clear: it seems there are many variables (latent, inclusive) and many paths, yet the N is less than 1500 (due to missing data). Shall we get an output, I wonder how large would that be, and how tough would it be to become "comprehensive" for the audience.

Dear Marcos

thank you for your comments and you nice welcome!

I know that this model is quite complex. Conceptually it makes perfect sense since all the separate pieces of the model were already shown to be related before. But so far not all together. You are right, the output (of the sem command model) is quite big but I'm quite positive that it would be comprehensive for the audience.

[Nina Schnyder]

I was unable to edit my prior post, so I'm adding two things I just noticed.


Maybe sheer complexity is a genuine problem. I looked through Nina's code again, and I think I count about 30 Likert items and 12 level 1 latent variables (each of which appear to have unique Likert items loading on them). But then, if I am reading right, for 6 of those same level 1 latent variables, Nina is also saying that the other 6 level 1 latent variables are determinants of them. Is that correct? That sounds like a recursive path, or whatever the terminology is, and now I'm way out of my depth. I would ask if that was theoretically justified in the first place. Also, I have no idea how SEM or GSEM deal with that recursivity (although apparently SEM was able to solve the model...).

Back to problem solving on the original proposed GSEM model. I see the SEM manual says that in some cases, changing the integration method to mode and curvature adaptive Gauss Hermite, i.e.

..., method(mcaghermite)

can work in multilevel models. Specifically, "It is accurate but not as fast. If you are fitting a multilevel model, there are cases where method 1 [that is, the default mvaghermite method] will not work but this method will work."

If not for that, the other integration methods (these being ghermite and laplace) may work; they are less accurate, but they can then be used as starting values. The manual also says that using fewer integration points (default is 7) to find a starting solution is a possible solution. I wonder if these solutions might help! I will note, though, that the manual says these are potential solutions to convergence problems, e.g. initial values not feasible or an infinite iteration log with a lot of (not concave) notations. It's not clear that this is the problem; it may be that the model is simply too complex to converge as a GSEM.

http://www.stata.com/manuals13/semin...pdf#semintro12

Dear Weiwen

thank you for your ideas. i will try them out.

[Roman Mostazir]

Defining this option may also help to converge, you can give it a go:

 technique(bhhh) iterate(1000) diff

[Nina Schnyder]

Dear all

thank you so much for your help. Unfortunately all suggestions didn't work. I was able to estimate the same model in R (lavaan package) within seconds.

Best
Nina

[Weiwen Ng]

Dear all

thank you so much for your help. Unfortunately all suggestions didn't work. I was able to estimate the same model in R (lavaan package) within seconds.

Best
Nina

I am disappointed that I was not able to help more! It looks like I should go familiarize myself with R.