Stata allows you to fit multilevel mixed-effects probit models with meprobit. A multilevel mixed-effects probit model is an example of a multilevel mixed-effects generalized linear model (GLM). You can fit the latter in Stata using meglm.

 

 

Let’s fit a crossed-effects probit model. A crossed-effects model is a multilevel model in which the levels of random effects are not nested. We investigate the extent to which two salamander populations, whiteside and roughbutt, cross-breed. We label whiteside males wsm, whiteside females wsf, roughbutt males rbm, and roughbutt females rbf. Our dependent variable y is coded 1 if there was a successful mating and 0 otherwise. Let’s fit our model:

 

. webuse salamander

. meprobit y wsm##wsf || _all: R.male || female:
note: crossed random effects model specified; option intmethod(laplace) implied
Fitting fixed-effects model:
Iteration 0: log likelihood = -223.01026
Iteration 1: log likelihood = -222.78736
Iteration 2: log likelihood = -222.78735
Refining starting values:
Grid node 0: log likelihood = -216.49485
Fitting full model:
Iteration 0: log likelihood = -216.49485 (not concave)
Iteration 1: log likelihood = -214.34477
Iteration 2: log likelihood = -209.72558 (not concave)
Iteration 3: log likelihood = -209.70603 (not concave)
Iteration 4: log likelihood = -209.70407
Iteration 5: log likelihood = -208.1606
Iteration 6: log likelihood = -208.11432
Iteration 7: log likelihood = -208.1121
Iteration 8: log likelihood = -208.11182
Iteration 9: log likelihood = -208.11182
Mixed-effects probit regression                                                  Number of obs = 360
Group Variable No. of Groups Observaation per Group Minimum Observaation per Group Average Observaation per Group Maximum
_all 1 360 360.0 360
female 60 6 6.0 6
Integration method: laplace Wald chi2(3) = 41.02 Log likelihood = -208.11182 Prob > chi2 = 0.0000
y Coef. Std. Err. z P>|z| [95% Conf. Interval] [95% Conf. Interval]
1.wsm -.4122025 .2736828 -1.51 0.132 -.9486109 .1242058
1.wsf -1.720337 .321778 -5.35 0.000 -2.351011 -1.089664
wsm#wsf

1 1

 2.121122 .361068 5.84 0.000 1.409446 2.832799
_cons .5950999 .2313942 2.57 0.010 .1415755 1.048624
_all<male

var (_cons)

 .3867405 .1792564 .1559132 .9593043
female

var (_cons)

 .4464109 .1982622 .1869368 1.066043

 

LR test vs. probit regression: chi2(2) = 29.35                                                 Prob > chi2 = 0.0000

Note: LR test is conservative and provided only for reference.
Our model has two random-effects equations, separated by ||. We use the _all notation that identifies all the observations as one big group. We use the R. notation to tell Stata to treat male as an indicator variable.

 

The output table includes the fixed-effect portion of our model and the estimated variance components. The estimates of the random intercepts suggest that the heterogeneity among the female salamanders is larger than the heterogeneity among the male salamanders.

 

If we wish, we can constrain the two random intercepts to be equal.

 

. constraint 1 _b[var(_cons[_all>male]):_cons] = _b[var(_cons[female]):_cons]

. meprobit y wsm##wsf || _all: R.male || female:, constraint(1) nolog

 

note: crossed random effects model specified; option intmethod(laplace) implied

 

Fitting fixed-effects model:

Iteration 0: log likelihood = -223.01026

Iteration 1: log likelihood = -222.78736

Iteration 2: log likelihood = -222.78735

 

Refining starting values:

Grid node 0:         log likelihood = -216.49485

 

Fitting full model:

Iteration 0: log likelihood = -216.49485 (not concave)

Iteration 1: log likelihood = -214.33961

Iteration 2: log likelihood = -209.6337

Iteration 3: log likelihood = -208.1889

Iteration 4: log likelihood = -208.14633

Iteration 5: log likelihood = -208.1479

Iteration 6: log likelihood = -208.14476

Iteration 7: log likelihood = -208.14476

 

Mixed-effects probit regression                                          Number of obs = 360

 

Group variable No. of Groups Observations per Group Minimum Observations per Group Average Observations per Group Maximum
_all 1 360 360.0 360
female 60 6 6.0 6

 

Integration method: laplace

 

Wald chi2(3) = 41.01

Log likelihood = -208.14476                                                       Prob > chi2 = 0.0000

( 1)   [var(_cons[_all>male])]_cons – [var(_cons[female])]_cons = 0

 

y Coef. Std. Err. z P>|z| [95% Conf. Interval] [95% Conf. Interval]
1.wsm -.4129849 .2728817 -1.51 0.130 -.9478232 .1193184
1.wsf -1.720523 .3219742 -5.34 0.000 -2.351581 -1.114613
wsm#wsf

1 1

 2.119226 .3605337 5.88 0.000 1.412553 2.825898
_cons .5962792 .2343886 2.54 0.011 .136886 1.055672
_all>male

var (_cons)

 .4155883 .1478019 .2069845 .8344279
female

var (_cons)

 .4155883 .1478019 .2069845 .8344279