Home  /  Products  /  Stata 15  /  Multilevel mixed-effects interval regression

This page announced the new features in Stata 15. Please see our Stata 18 page for the new features in Stata 18.

Multilevel mixed-effects interval regression

Highlights

  • Interval-measured outcomes, including
    • Interval-censored
    • Left-censored
    • Right-censored
  • Random effects
    • Random intercepts
    • Random coefficients
  • Multilevel: two, three, or more levels
  • Support for complex survey data
  • Graph marginal means and marginal effects
  • Intraclass correlation
  • Support for Bayesian estimation

What's this about?

The new meintreg command fits models in which the outcome is interval measured (interval-censored) and the observations are clustered.

Interval measured means that rather than the outcome (y) being observed precisely, it is known only that yl ≤ y ≤ yu in some or all observations. Observations can also be left-censored (y ≤ yl) or right-censored (y ≥ yu). An interval-measured outcome might be income, recorded in income brackets, or weekly minutes of exercise, recorded as less than 30 minutes, 31-59 minutes, 60-89 minutes, etc.

Multilevel mixed effects means that the fitted model accounts for clustering, such as when people live near each other or students attend the same school or students are tested repeatedly.

Let's see it work

We have fictional data on 8,424 people living in 50 cities. Different cities have different policies. One aspect of this is that some cities are more walker-friendly than others. Our data include a binary variable for being walker-friendly.

In these data, respondents were asked how many minutes per week they engage in physical activity outside of their jobs and were offered the following response categories: 0–30, 30–59, 60–89, 90–119, 120–239, and more than 240 minutes. This variable will be our dependent variable. It is both interval- and right-censored. The dataset contains variables exerlo and exerup. The values recorded in the variables are

Walking per week
exerlo exerup Meaning
0 30 0-30 minutes
30 59 30-59 minutes
60 89 60-89 minutes
90 119 90-119 minutes
120 239 120-239 minutes
240 . more than 240 minutes

Our independent variables will be

walkwalker-friendly city
agerespondent's age in years
workhours/week worked
kidschildren under age 12 in household

Even so, we expect that there are unmeasured characteristics about cities that matter—like weather, for example. Unmeasured characteristics would result in people who live in the same city behaving in overly similar ways, and so we will allow for city-level random effects. Variable cid records city ID number.

To fit a random-intercept model using exercise in the range exerlo and exerup, we type

. meintreg exerlo exerup age work kids walk || cid:

Fitting fixed-effects model:

Iteration 0:   log likelihood = -16228.837  
Iteration 1:   log likelihood = -16211.002  
Iteration 2:   log likelihood = -16209.922  
Iteration 3:   log likelihood =  -16209.92  
Iteration 4:   log likelihood =  -16209.92  

Refining starting values:

Grid node 0:   log likelihood =  -16183.23

Fitting full model:

Iteration 0:   log likelihood =  -16183.23  
Iteration 1:   log likelihood = -16021.782  
Iteration 2:   log likelihood = -15902.566  
Iteration 3:   log likelihood = -15856.315  
Iteration 4:   log likelihood = -15836.069  
Iteration 5:   log likelihood =  -15830.93  
Iteration 6:   log likelihood =  -15830.08  
Iteration 7:   log likelihood = -15830.052  
Iteration 8:   log likelihood = -15830.052  

Mixed-effects interval regression               Number of obs     =      8,424
                                                   Uncensored     =          0
                                                   Left-censored  =          0
                                                   Right-censored =         54
                                                   Interval-cens. =      8,370

Group variable:             cid                 Number of groups  =         50
                                                Obs per group:
                                                              min =         81
                                                              avg =      168.5
                                                              max =        249

Integration method: mvaghermite                 Integration pts.  =          7

                                                Wald chi2(4)      =     125.15
Log likelihood = -15830.052                     Prob > chi2       =     0.0000

Coef. Std. Err. z P>|z| [95% Conf. Interval]
age -.1745526 .0695285 -2.51 0.012 -.3108259 -.0382792
work -.3635766 .0794423 -4.58 0.000 -.5192807 -.2078726
kids -10.80659 1.102927 -9.80 0.000 -12.96829 -8.644894
walk 13.7689 6.514321 2.11 0.035 1.001061 26.53673
_cons 78.62265 5.01203 15.69 0.000 68.79925 88.44605
cid
var(_cons) 269.5222 56.84988 178.2592 407.5089
var(e.exerlo) 2440.81 40.2559 2363.172 2520.999
LR test vs. interval model: chibar2(01) = 759.74 Prob >= chibar2 = 0.0000

Just above the coefficient table is a header. It reports the number of censored and uncensored observations, reports the number of observations per city, and provides information about the fitted model.

The coefficient table itself is divided into three or four parts. The first part reports the fixed part of the model, that is, the effects of variables age, work, etc.

The second part reports the variance of the random intercept. We fit a model with a random intercept (_cons) whose estimated mean is 79 and standard deviation is 16 (a variance of 270).

The third part, var(e.exerlo), reports the variance of the model's error. Stata labels errors as e. followed by the dependent variable's name. Even so, e.exerlo is a bit misleading. A better name would have been e.exer, with exer being the imaginary name of the unobserved variable that is bounded by exerlo and exerup.

The fourth part—the part below the table—reports a test of whether the intercept is random. Anyone would reject the null hypothesis that it has variance equal to zero.

To obtain standard deviations instead of variances for the random portion of the model, we can type estat sd:

. estat sd
Coef. Std. Err. z P>|z| [95% Conf. Interval]
cid
sd(_cons) 16.41713 1.731419 13.35137 20.18685
sd(e.exerlo) 49.40456 .4074108 48.61246 50.20955

We can also obtain the residual intraclass correlation.

. estat icc

Residual intraclass correlation

Level ICC Std. Err. [95% Conf. Interval]
cid .0994425 .0189461 .0679834 .1432221

Tell me more

You can also fit Bayesian multilevel interval regression using the bayes prefix.

Learn more about Stata's multilevel mixed-effects models features.

Read more about multilevel interval regression in the Stata Multilevel Mixed-Effects Reference Manual.