Multilevel parallel mediation using gsem - convergence issues

[Joran Veldkamp]

Dear readers,

I run into model convergence issues when estimating a parellel mediation model using the gsem command. I want to predict delinquent behavior (count variable, variety score) in adolescence from a binary peer deviancy variable (which is is a "treatment" variable, created with propensity score matching analysis using 1-1 matching). In order to test which mechanism is producing this effect, I want to specify a multilevel path model with (1) normative attitudes towards delinquency and (2) positive reinforcements for delinquency as mediators, which are both continuous scales, and the matched pair indicator as 2nd level variable.

Below you can see a preview of all these variables. The numbers in the front indicate the waves (e.g. k7 = wave 7), while "treated" is thus the peer deviancy variable (from wave 5 actually).

* Example generated by -dataex-. To install: ssc install dataex
clear
input float(k7_variety treated k6_norms k6_reinforce)
 2 0  3.2  2.833333
 . 0  4.6  2.416667
 3 0  3.4       2.5
 2 0  4.6  2.583333
 1 0  5.8 2.3333333
 2 0  3.4 1.9166666
 3 0  4.2  2.583333
 2 0    3      2.25
 4 0    6  2.666667
 0 0    5  3.083333
 2 0    4      2.25
 2 0    4  3.166667
 3 0  5.4      3.25
 3 0  5.4      3.75
 0 0  4.4 2.0833333
 3 0  3.6 2.3333333
 4 0  4.2       3.5
 4 0  4.2      1.75
 2 0  3.6  3.416667
 0 0    5  2.833333
 3 0  2.8      2.25
 4 0  2.6       2.5
 . 0  4.8         3
 6 0  4.4         2
 6 0  1.8      1.75
 4 0  5.8         3
 2 0  2.6         2
 5 0  2.4      2.25
 7 0  4.4         2
 6 0  3.2 2.3333333
 8 0    3      1.75
 2 0  3.6 1.6666666
 2 0  6.6         3
 4 0  3.6 2.0833333
 2 0    7       3.5
 0 0    5  3.083333
 4 0  4.4 1.6666666
 2 0  5.8 1.8333334
 1 0    6  3.083333
 2 0  4.8  2.583333
 1 0    5 1.8333334
 0 0    6  2.833333
 5 0  4.2       2.5
 5 0  5.6       2.5
 4 0  4.6      2.25
 . 0    5 2.1666667
 4 0    6      2.75
 3 0  4.6  3.083333
 5 0  2.6 2.1666667
 3 0  4.2  3.583333
 1 0  5.2  3.166667
 2 0  5.6 2.1666667
 . 0    6  2.833333
 3 0  3.8 2.1666667
 1 0    6      2.25
 1 0  2.8 1.6666666
 2 0  4.2  2.916667
 3 0  3.8      1.75
 5 0    3 1.9166666
 2 0    4  3.083333
 1 0  2.2 1.9166666
 1 0  4.8  2.416667
 4 0    4         3
 3 0  2.8         3
 2 0  3.8 2.3333333
 1 0  4.8  2.666667
15 0    4         2
 1 0  2.6 2.1666667
 2 0  3.8  2.583333
 3 0  3.8 2.1666667
 2 0  6.2  3.166667
 2 0  3.4 1.8333334
 4 0    4  2.583333
 0 0    6      3.75
 4 0  5.2      2.75
 3 0  6.2  3.083333
 4 0  4.6 2.0833333
 8 0    4         2
 2 0  4.4         3
 4 0  5.8         3
 2 0  4.2 1.9166666
 1 0  3.6         2
 0 0  4.8 1.9166666
 3 0  4.2       2.5
 3 0  4.2  3.083333
 3 0  4.8      2.25
11 0    5 1.8333334
 3 0  3.8  2.666667
 0 0    5  2.416667
 2 0  4.8 2.3333333
 3 0    4  2.583333
 0 0  5.4  3.166667
 2 0  6.6  3.833333
 2 0  4.2      2.25
 6 0 2.25 2.0833333
 2 0  5.6  2.833333
 2 0    4      2.25
 7 0  1.8 1.6666666
 4 0  4.2         3
 2 0    5 2.3333333
end

However, when running the gsem multilevel model as specified below, I get the error message "cannot compute an improvement -- discontinuous region encountered" after 39 iterations. I am not sure what to change about the options, because I do need the covariance structure like that and the two errors of the mediators to be allowed to have covariance (since they are correlated). Note that the "if tots1a==4" ensures we only select the matched pairs.

gsem (k7_variety      <- treated k6_norms k6_reinforce  M1[pair], ///
                                                                        family(nb) link(log)) ///
         (k6_norms       <- treated                                        M2[pair]) ///
         (k6_reinforce   <- treated                                        M3[pair]) ///
                                                                                   if tots1a==4, ///
     covstructure(M1[pair]*M2[pair]@0 M1[pair]*M3[pair]@0 M2[pair]*M3[pair]@0 ///
     e.k6_norms*e.k6_reinforce)

Breaking the code up and running it with one mediator at a time gives some more information about where the issue is coming from, because the model with k6_norms as mediator runs just fine, while the model with k6_reinforce as mediator does give the same convergence issues.I am not sure where the problem with the k6_reinforce comes from, because specifying other models (e.g. mixed k7_variety treated k6_reinforce if tots1a==4 || pair:, var) gives no issues, so its only when the mediation pathways are specified.

Does anyone know if I need to change some options of this model? Or any other suggestion to get successful convergence, which makes intuitively sense as well?

Thank you in advance for an answer!

Kind regards,
Joran Veldkamp