Binary variables

Binary variables have 2 levels, typically 0 (false/failure) and 1 (true/success).

Below are 8 observations of binary variable \(y\):

\[y=(0,1,1,1,1,0,1,1)\]

The mean of a 0/1 variable is the proportion of ones:

\[\bar{y}=\frac{\sum_iy}{n}=.75\]

The mean is also an estimate of the probability of ones in the population, \(p=Pr(y=1)\):

\[\bar{y}=\hat{p}=.75\]

The \(\hat{}\) symbol denotes a sample estimate of a population parameter.

Binomial distribution

With \(n\) observations of a binary variable, we can observe between 0 and \(n\) ones (successes).

Over repeated samples of size \(n\), the number of ones we observe in each sample is determined by the probability of success, \(p\).

For samples of size \(n\) with probability of success \(p\), the binomial distribution describes the probability of observing each number of successes, from 0 to \(n\).

The binomial distribution has 2 parameters:

• \(n\), the sample size or number of trials (typically known so not estimated)
• \(p\), the probability of success

A histogram of the number of successes (ones) in 1000 samples of size \(n=8\) from a population with \(p=.25\) approximates the binomial distribution for \(n=8\), \(p=.25\):

When the \(p\) is low, most samples contain more zeroes than ones.

When \(p\) is higher, say \(p=.75\), most samples contain more ones than zeroes.

Let \(y^*\) be the number of successes in a sample of size \(n\) of binary variable \(y\).

The exact probability of observing \(k\) successes in samples of size \(n\) where the probability of success is \(p\) is given by:

\[Pr(y^*=k|n,p) = {n\choose k}p^k(1-p)^{(n-k)}\]

Mean and variance of binomial distribution

For samples of binomially-distributed \(y\) of size \(n\) with probability of success \(p\), the mean (expected) number of successes is \(np\):

\[E(y|n,p) = np\]

For samples with \(n=8\) and probability \(p=.75\), the mean is \(6\):

\[E(y|n=8, p=.75) = (8)(.75) = 6\]

The variance of \(y\) is:

\[Var(y) = np(1-p)\]

Both mean and variance of \(y\) vary with \(p\), so the binomial distribution is heteroskedastic.

• \(Var(y)\) is maximized when \(p=.5\): half zeroes and half ones.
• \(Var(y)\) is minimized when \(p=0\) of \(p=1\), all zeroes or all ones, respectively.

Representation of binary data

In most data sets, each row will be a single observation of a binary variable (y below):

Above, each row is a single observation of a binomial variable y, so \(n=1\)^† and \(y\) will be either 0 or 1.

When \(n=1\), the mean of binomially-distributed \(y\) is \(p\) itself:

\[E(y|n=1,p) = np = (1)p = p\]

Regression models mean outcomes, so we are modeling \(p\) with \(n=1\).

Stata’s logit command, discussed extensively later, expects data with \(n=1\) per row.

Aggregated data

Regression assumes that observations with the same predictor values (x above) come from the same population with a common \(p=Pr(y=1)\), so we could aggregate the data like so:

Now \(y\) is the count of ones and \(n\) is the sample size.

In either style of data, \(n\) is known so only \(p\) will be modeled.

The 2 datasets above will produce the same logistic regression estimates, but aggregated data must be analyzed in blogit in Stata rather than logit.

^†The binomial distribution with \(n=1\) is also known as the Bernoulli distribution.

Analyzing binomially distributed variables in linear regression

Regression models estimate the relationship between the mean outcome and a set of predictors.

As seen above, for binary outcomes the mean is \(p\) itself.

We could model the relationships between predictors and \(p\) in a linear regression, but some assumptions will be violated:

• homoskedasticity is assumed, but binomial variables are heteroskedastic
• normally-distributed errors are assumed, which cannot be true with binary outcomes
• outcome is assumed unbounded, but probabilities are bound between 0 and 1
  • transforming the outcome can extend its bounds, but the other problems persist

Logistic regression addresses each issue above.

Odds

Regression models linear relationships between the mean outcome and the predictors, which can be problematic for bounded outcomes:

Transforming the outcome can extend its range.

In logistic regression, two transformations on the mean outcome (probability) change its bounds from \([0,1]\) to \((-\infty, \infty)\).

We first transform probabilities to odds:

\[odds(p) = \frac{p}{1-p}\]

Odds lie within \([0,\infty)\) and have a positive, monotonic relationship with probability \(p\).

Odds are interpreted as the expected number of successes per failure.

For example, if \(p=.75\), then \(odds(p)=\frac{.75}{.25}=3\), which means we expect 3 successes for every failure.

The next transformation extends the left bound to \(-\infty\).

Logits

If variable \(x\) is defined over \((0,\infty)\), then \(log(x)\) is defined over \((-\infty, \infty)\).

The log of the odds, or logit, is thus defined over \((-\infty, \infty)\).

\[logit(p) = log(odds(p)) = log(\frac{p}{1-p})\]

The logit has an S-shaped (sigmoidal) relationship with the probability:

The logit function is often referred to as a link function, which transforms the mean outcome to allow a linear relationship with the predictors.

The Logistic Regression Model

The logistic regression model asserts a linear relationship between the logit (log-odds) of the outcome variable and the predictors \(x_1\) through \(x_k\) for observation \(i\):

\[logit(p_i) = log(\frac{p_i}{1-p_i}) = \beta_0 + \beta_1x_{i1} + ... + \beta_kx_{ik}\]

The linear predictor, \(\beta_0 + \beta_1x_{i1}+...+\beta_kx_{ik}\), expresses the expected outcome in logits (log odds).

The inverse logit function transforms a logit back into a probability:

\[logit^{-1}(x) = \frac{exp(x)}{1+exp(x)}\]

\[logit^{-1}(logit(p_i)) = p_i\]

Thus, the inverse logit of the linear predictor expresses the expected outcome as a probability:

\[p_i = \frac{exp(\beta_0 + \beta_1x_{i1}+...+\beta_kx_{ik})}{1+exp(\beta_0 + \beta_1x_{i1}+...+\beta_kx_{ik})}\]

In the above expression, we see that \(p_i\) has a non-linear relationship with the predictors.

• logit of the outcome has a linear relationship with the predictors
• logit has an S-shaped relationship with probability
• thus, probability has an S-shaped relationship with the (continuous) predictors

Because of this nonlinear relationship, the amount of change in \(p\) due to a one-unit increase in \(x\) depends on the value of \(x\).

Interpreting logistic regression coefficients

Logistic regression coefficients can be interpreted in at least 2 ways.

The first interpretation expresses additive changes in the logits, or log odds of the outcome.

Imagine a simple logistic regression model with a single predictor \(x\):

\[logit(\frac{p_i}{1-p_i}) = \beta_0 + \beta_1x_i\]

Setting \(x=0\), the expected outcome is \(\beta_0\):

\[\begin{align} logit(\frac{p_{x=0}}{1-p_{x=0}}) &= \beta_0 + \beta_1(0) \\ &= \beta_0 \end{align}\]

If we set \(x=1\), the expected outcome is \(\beta_0 + \beta_1\):

\[\begin{align} log(\frac{p_{x=1}}{1-p_{x=1}}) &= \beta_0 + \beta_1(1) \\ &= \beta_0 + \beta_1 \end{align}\]

Taking the difference between these expectations gives us an interpretation of \(\beta_1\):

\[\begin{align} log(\frac{p_{x=1}}{1-p_{x=1}}) - log(\frac{p_{x=0}}{1-p_{x=0}}) &= (\beta_0 + \beta_1) - \beta_0 \\ &= \beta_1 \end{align}\]

Thus, \(\beta_1\) is the expected change in the logit (log-odds) of the outcome for a one-unit increase in the predictor \(x\).

Interpreting coefficients as change in logits is simple because the relationships are linear, but most audiences will struggle with understanding magnitudes in logits.

Odds ratios

A more common way to interpret logistic regression coefficients is to exponentiate them to express odds ratios.

Recall that the exponential function is the inverse of the (natural) logarithm function:

\[exp(log(x)) = x\] Also recall that the difference of two logs is the log of the ratio:

\[log(a) - log(b) = log(\frac{a}{b})\]

Logistic regression coefficients denote the change in the log of the odds of the outcome for a unit increase in the predictor, which we can express as the log of the ratio of the odds:

\[\begin{align} \beta_1 &= log(\frac{p_{x=1}}{1-p_{x=1}}) - log(\frac{p_{x=0}}{1-p_{x=0}}) \\ &= log(odds_{x=1}) - log(odds_{x=0}) \\ &= log(\frac{odds_{x=1}}{odds_{x=0}}) \\ \end{align}\]

Now, if we exponentiate both sides, we see how the exponentiated coefficient is an odds ratio:

\[\begin{align} exp(\beta_1) &= exp(log(\frac{odds_{x=1}}{odds_{x=0}})) \\ &= \frac{odds_{x=1}}{odds_{x=0}} \\ &= \frac{\frac{p_{x=1}}{1-p_{x=1}}}{\frac{p_{x=0}}{1-p_{x=0}}} \\ \end{align}\]

Odds ratios express the expected multiplicative change in the odds of the outcome for a unit increase in the predictor.

Workshop data set

We use a data set of birth weights of 180 infants to demonstrate logistic regression in Stata, where infants weighing less than 2500g are classified as low birth weight (lbw).

Outcome

• low: 1=low birth weight, 0=normal weight

Predictors

• age: mother’s age in years (14-45)
• ptl: count of previous premature labors (0-3)
• ht: mother’s history of hypertension, 1=hypertensive, 0=non-hypertensive
• smoke: mother’s history of smoking, 1=smoker, 0=non-smoker
• race: 1=white, 2=black, 3=other

Logistic regression will estimate linear relationships of each predictor above with the log odds of lbw.

Loading and exploring the low birth weight data set in Stata

The low birth weight data set comes with Stata and can be loaded with the command:

/* load low birth weight data */
webuse lbw, clear

We will explore the variables with histograms

/* histograms of predictors */
hist low, discrete freq

hist age, freq

hist ptl, freq discrete

hist ht, freq discrete

hist smoke, freq discrete

hist race, freq discrete

The logit command

Stat’s logit command performs logistic regression using the following syntax:

logit depvar indepvarlist

where depvar is the dependent variable (outcome) and indepvarlist is a list of one or more independent variables (predictors).

• depvar is generally coded 0/1, and Stata will model the probability that depvar=1
• Precede categorical predictors with i. in indepvarlist to request that Stata enter dummy variables representing each category but the first (by default) into the regression
  • Stata refers to variables preceded by i. as factors
• Add the option or after a , after indepvarlist to produce odds ratios; otherwise coefficients are in the logit metric by default

Stata’s logistic command estimates the same model as logit but can only output odds ratios.

Logistic regression of lbw data

We will estimate the relationships of mother’s age, previous premature labors, history of hypertension, smoking, and race with the log odds of having a low birth weight (lbw) baby.

• linear slopes (in logits) are estimated for age and ptl
• i. precedes each of ht, smoke and race to request that dummy variables for each category but the first be entered into the model

/* logistic regression with low as outcome,
   age and ptl as continuous predictors,
   ht, smoke and race as categorical predictors */
logit low age ptl i.ht i.smoke i.race

Iteration 0:  Log likelihood =   -117.336  
Iteration 1:  Log likelihood =   -105.286  
Iteration 2:  Log likelihood = -104.90511  
Iteration 3:  Log likelihood = -104.90441  
Iteration 4:  Log likelihood = -104.90441  

Logistic regression                                     Number of obs =    189
                                                        LR chi2(6)    =  24.86
                                                        Prob > chi2   = 0.0004
Log likelihood = -104.90441                             Pseudo R2     = 0.1059

------------------------------------------------------------------------------
         low | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         age |  -.0469041   .0353279    -1.33   0.184    -.1161454    .0223372
         ptl |   .7501393   .3298036     2.27   0.023     .1037362    1.396542
        1.ht |   1.227762   .6246189     1.97   0.049     .0035312    2.451992
             |
       smoke |
     Smoker  |   .9728681   .3863386     2.52   0.012     .2156584    1.730078
             |
        race |
      Black  |    .973844   .5013314     1.94   0.052    -.0087474    1.956435
      Other  |   1.021987   .4191365     2.44   0.015     .2004942    1.843479
             |
       _cons |  -.9237957    .891843    -1.04   0.300    -2.671776    .8241845
------------------------------------------------------------------------------

Starting towards the top of the output (after the likelihood iterations we see):

• 189 observations were used
• The likelihood ratio \(\chi^2\) statistic equals 24.86, which assesses the null hypothesis that all model coefficients (except the intercept _cons)
  • tests if current model improves fit over model with no predictors
  • p-value=0.004 suggests that at least one coefficient is not zero
• The \(Pseudo R^2\) statistic, a goodness-of-fit measure similar to \(R^2\) for linear regression, but not interpreted as proportion of outcome variance explained
• A table of regression coefficients and their associated standard errors, test statistics (z) and p-values, and 95% confidence intervals
  • z and p-values assess null hypothesis that coefficient is 0

Interpreting logistic regression coefficients

logit low age ptl i.ht i.smoke i.race

Iteration 0:  Log likelihood =   -117.336  
Iteration 1:  Log likelihood =   -105.286  
Iteration 2:  Log likelihood = -104.90511  
Iteration 3:  Log likelihood = -104.90441  
Iteration 4:  Log likelihood = -104.90441  

Logistic regression                                     Number of obs =    189
                                                        LR chi2(6)    =  24.86
                                                        Prob > chi2   = 0.0004
Log likelihood = -104.90441                             Pseudo R2     = 0.1059

------------------------------------------------------------------------------
         low | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         age |  -.0469041   .0353279    -1.33   0.184    -.1161454    .0223372
         ptl |   .7501393   .3298036     2.27   0.023     .1037362    1.396542
        1.ht |   1.227762   .6246189     1.97   0.049     .0035312    2.451992
             |
       smoke |
     Smoker  |   .9728681   .3863386     2.52   0.012     .2156584    1.730078
             |
        race |
      Black  |    .973844   .5013314     1.94   0.052    -.0087474    1.956435
      Other  |   1.021987   .4191365     2.44   0.015     .2004942    1.843479
             |
       _cons |  -.9237957    .891843    -1.04   0.300    -2.671776    .8241845
------------------------------------------------------------------------------

The estimated coefficients above are interpreted as change in logits (log odds of being lbw).

Continuous predictors

• for each additional year of age, the log odds of lbw decrease by .047
• for each additional premature labor, the log odds of lbw increase by .75

Categorical predictors

• mothers with a history of hypertension are expected to have an increase of 1.23 in the log odds of lbw compared to mothers without a history
• mothers who smoked during pregnancy are expected to have an increase of .973 in the log odds of lbw compared to non-smokers
• black mothers are expected to have an increase of .973 in the log odds of lbw compared to white mothers
• other race mothers are expected to have an increase of 1.02 in the log odds of lbw compared to white mothers

Interpreting odds ratios

Odds ratios are interpreted in a couple ways (let \(OR\) be some odds ratio):

• \(OR\) is the factor or multiplicative change in the odds of the outcome
• \((OR-1)\times100\%\) is the percent change in the odds of the outcome

An odds ratio above one denotes an increase in the odds, below one denotes a decrease in the odds, and equal to one denotes no change in the odds.

We repeat our logistic regression of low, now with odds ratios by adding the or option:

/* rerunning model to express odds ratios */
logit low age ptl i.ht i.smoke i.race, or

Iteration 0:  Log likelihood =   -117.336  
Iteration 1:  Log likelihood =   -105.286  
Iteration 2:  Log likelihood = -104.90511  
Iteration 3:  Log likelihood = -104.90441  
Iteration 4:  Log likelihood = -104.90441  

Logistic regression                                     Number of obs =    189
                                                        LR chi2(6)    =  24.86
                                                        Prob > chi2   = 0.0004
Log likelihood = -104.90441                             Pseudo R2     = 0.1059

------------------------------------------------------------------------------
         low | Odds ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         age |   .9541789   .0337091    -1.33   0.184     .8903458    1.022589
         ptl |   2.117295   .6982914     2.27   0.023     1.109308    4.041203
        1.ht |   3.413581   2.132187     1.97   0.049     1.003537    11.61146
             |
       smoke |
     Smoker  |   2.645521   1.022067     2.52   0.012     1.240679    5.641093
             |
        race |
      Black  |   2.648104   1.327578     1.94   0.052     .9912907    7.074066
      Other  |   2.778709   1.164658     2.44   0.015     1.222007    6.318482
             |
       _cons |   .3970093   .3540699    -1.04   0.300     .0691294    2.280021
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.

Note: test (z) statistics and p-values are the same whether regression coefficients or odds ratios are reported

Continuous predictors

• For each additional year of age, the odds of lbw decrease by 4.6% (\((.954-1)\times100\%=-4.6\%\))
• For each additional premature labor, the odds of lbw increase by a factor of 2.12

Categorical predictors

• mothers with a history of hypertension are expected to have 3.41 times the odds of lbw compared to mothers without a history
• mothers who smoked during pregnancy are expected to have a 2.65 factor increase in the odds of lbw compared to non-smokers
• black mothers are expected to have 165% (\((2.65-1)\times100\%=165\%\)) increase in the odds of lbw compared to white mothers
• other race mothers are expected to have 178% (\((2.78-1)\times100\%=178\%\)) increase in the odds of lbw compared to white mothers

Exponentiated interaction coefficients

Interactions in regression models allow the effect of one predictor to vary with another predictor.

For example, we might believe that the effect of age is different for smoking vs nonsmoking mothers.

Interaction terms are generated by taking the product of the component variables.

Stata provides the # and ## operator for interactions, where the former requests only the interaction and the latter requests both the lower order terms and the interaction.

Below, we add the interaction of smoke and age to the model:

/* adding interaction of age and smoke */
logit low c.age##i.smoke ptl i.ht i.race, or

Iteration 0:  Log likelihood =   -117.336  
Iteration 1:  Log likelihood = -105.18462  
Iteration 2:  Log likelihood = -104.70994  
Iteration 3:  Log likelihood =  -104.7077  
Iteration 4:  Log likelihood =  -104.7077  

Logistic regression                                     Number of obs =    189
                                                        LR chi2(7)    =  25.26
                                                        Prob > chi2   = 0.0007
Log likelihood = -104.7077                              Pseudo R2     = 0.1076

------------------------------------------------------------------------------
         low | Odds ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         age |   .9336296   .0471014    -1.36   0.173       .84573    1.030665
             |
       smoke |
     Smoker  |   .9512806   1.593868    -0.03   0.976     .0356552    25.38017
             |
 smoke#c.age |
     Smoker  |   1.045964   .0751283     0.63   0.532     .9086104    1.204082
             |
         ptl |    2.10686    .694485     2.26   0.024     1.104215    4.019921
        1.ht |   3.469858   2.170839     1.99   0.047     1.018067    11.82625
             |
        race |
      Black  |    2.48018   1.269646     1.77   0.076     .9093663    6.764375
      Other  |   2.704271   1.141077     2.36   0.018     1.182722    6.183265
             |
       _cons |    .660881   .8026954    -0.34   0.733     .0611319    7.144616
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.

In the logit metric, interaction coefficients are interpreted as additive changes in a predictor’s effects due to a one-unit increase in the interacting predictor.

Exponentiated interaction coefficients are not odds ratios, but instead are the ratio of one odds ratio vs another.

The exponentiated interaction coefficient above, 1.046, can be interpreted as:

• the ratio comparing the odds ratio of age for smokers to nonsmokers
  • the odds ratio of age for smokers is 1.046 times larger (or 4.6% larger) than the odds ratio of age for nonsmokers (which is .934m the OR for smoke in the table above)
• the ratio comparing the odds ratio of smoke between two mothers whose ages differ by one year (older mother in numerator)
  • a mother who is one year older than another mother will have an odds ratio for smoking that is 4.6% larger (holding other predictors constant)

Exercise 1

To practice logistic regression in Stata, we will use the nhanes2 dataset, coming from the National Health and Nutrition Examination Survey, that can be loaded with:

webuse nhanes2, clear

The model variables are:

• highbp: high blood pressure, 1=high blood pressure, 0=normal
• region: region of US, 1=NE, 2=MW, 3=S, 4=W
• rural: rural residence, 0=urban, 1= rural
• weight: weight in kg
• age: age in years

1. Perform a logistic regression of outcome highbp with categorical predictors (factors) region and rural, and linear relationships for continuous predictors weight and age, with coefficients expressed in logits (log odds). Interpret the coefficients.
2. Rerun the model with coefficients expressed as odds ratios and interpret the odds ratios.

Predicted outcomes

Understanding the magnitude of predictor effects as changes in logits or as odds ratios can be difficult.

Expressing model effects in probabilities can make magnitudes clearer.

In this section we show how Stata’s predict and margins commands predict model-based probabilities of the outcome and how to use these to interpret model results.

Because of the nonlinear relationship between the predictors and the probability of the outcome:

• a one-unit change in a predictor will not always lead to the same change in probability of the outcome
• the change in probability due to a change in a predictor depends on the value of the predictor itself as well as values of the other predictors

At times, it can also be useful to predict outcomes in logits, for example, to assess assumptions of linear relationships.

predict after logistic regression

Stata’s predict command is primarily used to make model-based predictions for the current data.

The basic syntax is:

predict newvar

where newvar is the name of new variable that will hold the predictions.

By default, after logit, predict will estimate predicted probabilities, but many other types of quantities can be predicted with various options (see logit postestimation).

After running predict, newvar will be added to the data set with the predictions.

/* predp will be predicted probabilities */
predict predp

(option pr assumed; Pr(low))

Here are the model variables and predicted probabilities from the model for the first 10 observations:

/* display model variables and predictions for first 10 obs */
list low age ptl ht smoke race predp in 1/10

     | low   age   ptl   ht       smoke    race      predp |
     |-----------------------------------------------------|
  1. |   0    19     0    0   Nonsmoker   Black   .3012971 |
  2. |   0    33     0    0   Nonsmoker   Other   .1900565 |
  3. |   0    20     0    0      Smoker   White   .2913144 |
  4. |   0    21     0    0      Smoker   White   .2817266 |
  5. |   0    18     0    0      Smoker   White    .311053 |
     |-----------------------------------------------------|
  6. |   0    21     0    0   Nonsmoker   Other   .2917718 |
  7. |   0    22     0    0   Nonsmoker   White   .1239348 |
  8. |   0    17     0    0   Nonsmoker   Other   .3319944 |
  9. |   0    29     0    0      Smoker   White   .2122952 |
 10. |   0    26     0    0      Smoker   White   .2367767 |
     +-----------------------------------------------------+

A model with high predictive power will predict probabilities near zero for observations where low=0, and near one for observations where low=1.

Let’s examine the distribution of predicted probabilities by observed outcome low:

/* compare distribution of predp between observations with low=1 vs low=0 */
graph box predp, over(low) b1title("observed low") ytitle("predicted Pr(low=1)")

The distributions are not very separated, indicating that the model has low predictive power.

margins and marginal estimates

Stata’s margins command is one its most powerful and versatile commands, which estimates many kinds of model predictions, including:

• predictive margins: predicted outcomes at fixed values of some predictors and averaging over the distribution of other predictors
  • also known as marginal means
  • represents the average predicted outcome in the population or some subpopulation
• adjusted predictions: predicted outcomes at fixed values of all predictors
  • also known as conditional margins
  • represent the predicted outcome for a specific type of person (or unit)
• marginal effects: effects of predictors, holding some predictors fixed and averaging over the distributions of others
  • typically estimated as a difference on the original scale of the outcome (difference in probabilities for logistic regression)
  • are equal to partial derivatives for continuous predictors

Because of its versatility, margins has an extensive syntax that can be difficult to master (the manual entry for margins is 78 pages long!).

Nevertheless, margins functionality is one of the primary reasons to use Stata.

margins can only be used after a model has been estimated.

Predictive margins

Predictive margins are predicted outcomes averaged over the distribution of some (or all) predictors and at fixed levels of the others.

Generally, when specifying margins, whichever predictors are not fixed at specific levels will be averaged over when predicting the outcome.

Thus, simply specifying margins after our logistic regression of lbw will result in a predicted probability of lbw, averaging over the distributions of all predictors in the model:

/* probability of lbw, averaging over all predictors */
margins

Predictive margins                                         Number of obs = 189
Model VCE: OIM

Expression: Pr(low), predict()

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
       _cons |   .3121693   .0314715     9.92   0.000     .2504864    .3738523
------------------------------------------------------------------------------

Note: z-statistics and p-values assessing null hypotheses testing predictive margins against zero are generally ignored (we generally are not interested in testing whether a predicted probability is different from zero).

The above predictive margin:

• is interpreted as the estimated average probability of lbw in the population (averaging over or regardless the predictors)
• equals the observed probability of the outcome
• also equals the mean of our model-predicted probabilities for each observation, predp

/* means of low and predp equals the predictive margin above */
summ low predp

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
         low |        189    .3121693    .4646093          0          1
       predp |        189    .3121693    .1663442   .0458923   .7944289

Predictive margins for subpopulations

We can also estimate predictive margins for specific subpopulations.

The syntax:

margins factorvar

estimates predictive margins for each level of factorvar (predictor preceded by i.), averaging over the distribution of the other predictors.

Thus, the following margins estimates predicted probabilities of lbw for smoking and non-smoking mothers, averaging over the distribution of all other predictors:

/* average probabilities of lbw for smokers and nonsmokers */
margins smoke

Predictive margins                                         Number of obs = 189
Model VCE: OIM

Expression: Pr(low), predict()

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
       smoke |
  Nonsmoker  |   .2400186   .0387952     6.19   0.000     .1639814    .3160559
     Smoker  |   .4283874   .0580954     7.37   0.000     .3145226    .5422522
------------------------------------------------------------------------------

The predictive margin for Nonsmoker above can be estimated manually with the following procedure:

1. Convert everyone in the sample to smoke=0
2. Leave all other predictors at their observed values
3. Predict probabilities for everyone based on the model estimates
4. The mean of these predicted probabilities is the predictive margin

By leaving the other predictors at their observed values when predicting and then averaging the predictions, the resulting mean is averaged over the distribution of those predictors.

The predictive margin for Smoker is obtained using the same procedure except everyone is converted to smoke=1 in the first step.

Use # to obtain predictive margins at each combination of the factors’ levels:

/* predictive margins for each combination of smoke and race */
margins smoke#race

Predictive margins                                         Number of obs = 189
Model VCE: OIM

Expression: Pr(low), predict()

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
  smoke#race |
  Nonsmoker #|
      White  |   .1551203   .0467924     3.32   0.001     .0634088    .2468317
  Nonsmoker #|
      Black  |    .314714   .0899992     3.50   0.000     .1383188    .4911092
  Nonsmoker #|
      Other  |   .3245289    .056837     5.71   0.000     .2131304    .4359273
     Smoker #|
      White  |   .3145167   .0568889     5.53   0.000     .2030165     .426017
     Smoker #|
      Black  |   .5340904   .1121103     4.76   0.000     .3143583    .7538226
     Smoker #|
      Other  |   .5454064   .0971575     5.61   0.000     .3549812    .7358315
------------------------------------------------------------------------------

The predictive margins above are the predicted probabilities of lbw for each combination of smoking status and race, averaging over age, premature labors, and history of hypertension.

Fixing predictor levels in margins with at()

We can fix the levels of any predictor in margins, including continuous predictors, with the at() option:

• specify at() after a ,
• within at()
  • specify the names of predictors whose levels are to be fixed
  • after each predictor’s name, specify = and then the value at which to fix the predictor
  • to fix the predictor at multiple values, after =, put the values inside () using one of the following notations
    • a list of numbers
    • a/b requests every integer from a to b
    • a(c)b requests all values from a to b in increments of c

Remember that any predictor not specified in margins will be averaged over when predicting outcomes.

First, we estimate the predictive margin for mothers who are age 20:

/* average probability of lbw for 20-year old mothers */
margins, at(age=20)

Predictive margins                                         Number of obs = 189
Model VCE: OIM

Expression: Pr(low), predict()
At: age = 20

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
       _cons |   .3365013   .0372488     9.03   0.000     .2634951    .4095075
------------------------------------------------------------------------------

Next are predictive margins for smoking and non-smoking mothers who are aged 20 or 30:

/* average probabilities of lbw for smoking and non-smoking moms age 20 or 30 */
margins smoke, at(age=(20 30))

Predictive margins                                         Number of obs = 189
Model VCE: OIM

Expression: Pr(low), predict()
1._at: age = 20
2._at: age = 30

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
   _at#smoke |
1#Nonsmoker  |   .2612321   .0444322     5.88   0.000     .1741465    .3483176
   1#Smoker  |   .4590947   .0625156     7.34   0.000     .3365665     .581623
2#Nonsmoker  |    .187449   .0493688     3.80   0.000     .0906879      .28421
   2#Smoker  |   .3577227   .0789483     4.53   0.000     .2029869    .5124585
------------------------------------------------------------------------------

And average predicted probabilities for mothers aged between 20, 25, 30, 35 or 40 and who have had 0, 1, or 2 premature labors:

/* average probabilities of lbw for mothers aged 20 to 40 (in increments of 5) and 
   either 0 to 2 premature labors */
margins, at(age=(20(5)40) ptl=(0/2))

Predictive margins                                         Number of obs = 189
Model VCE: OIM

Expression: Pr(low), predict()
1._at:  age = 20
        ptl =  0
2._at:  age = 20
        ptl =  1
3._at:  age = 20
        ptl =  2
4._at:  age = 25
        ptl =  0
5._at:  age = 25
        ptl =  1
6._at:  age = 25
        ptl =  2
7._at:  age = 30
        ptl =  0
8._at:  age = 30
        ptl =  1
9._at:  age = 30
        ptl =  2
10._at: age = 35
        ptl =  0
11._at: age = 35
        ptl =  1
12._at: age = 35
        ptl =  2
13._at: age = 40
        ptl =  0
14._at: age = 40
        ptl =  1
15._at: age = 40
        ptl =  2

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         _at |
          1  |    .303371   .0387435     7.83   0.000      .227435    .3793069
          2  |   .4632617   .0746991     6.20   0.000     .3168543    .6096692
          3  |   .6311228    .135043     4.67   0.000     .3664434    .8958023
          4  |   .2598657   .0371378     7.00   0.000     .1870771    .3326544
          5  |    .410806   .0679355     6.05   0.000     .2776549     .543957
          6  |   .5799304   .1349761     4.30   0.000     .3153822    .8444786
          7  |   .2204237    .052976     4.16   0.000     .1165927    .3242547
          8  |   .3602354    .079823     4.51   0.000     .2037853    .5166855
          9  |    .527131   .1434264     3.68   0.000     .2460204    .8082417
         10  |   .1852631   .0689992     2.69   0.007     .0500271    .3204991
         11  |   .3124608   .0995509     3.14   0.002     .1173446    .5075769
         12  |    .473814   .1592928     2.97   0.003     .1616059    .7860222
         13  |   .1544039   .0806383     1.91   0.056    -.0036444    .3124521
         14  |   .2681937   .1185453     2.26   0.024     .0358491    .5005383
         15  |   .4211062   .1789255     2.35   0.019     .0704186    .7717937
------------------------------------------------------------------------------

Above we see that Stata uses a shorthand _at variable to label the estimates with a key above the margins table.

Graphing margins estimates with marginsplot

Another compelling reason to use margins is marginsplot, which graphs the estimates of the last margins command.

A plot of predicted outcomes across a range of predictor values can be helpful to interpret the predictor’s effects.

marginsplot will plot the predicted outcome at the predictor values specified in the previous margins and will generally connect them with a line.

Let’s examine a plot of the predictive margins for smoking and non-smoking mothers:

/* estimate marginal means for smoking and non-smoking mothers */
margins smoke
/* plot marginal means */
marginsplot

If more than one predictor is specified in margins, marginsplot will place one predictor on the x-axis and will plot separate lines by the other predictors.

/* estimate marginal means for smoking and non-smoking mothers at ages 20, 30 and 40 */
margins smoke, at(age=(20 30 40))
* plot marginal means
marginsplot

Use the x() (shorthand for xdimension()) option to specify which variable appears on the x-axis:

/* estimate marginal means for smoking and non-smoking mothers at ages 20, 30 and 40 */
margins smoke, at(age=(20 30 40))
/* plot marginal means */
marginsplot, x(smoke)

If many predictors are specified in margins, the subsequent marginsplot may be cluttered with too many lines:

/* margins with many predictors fixed at many levels */
margins smoke#race, at(age=(20 30 40) ptl=(0 1))
/* too many lines in marginsplot */
marginsplot

Specify some of the predictors in a by() option to plot separate graphs by those variables:

/* separate graphs by race and ptl */
marginsplot, by(race ptl)