A2.7 Assessing nonlinearities

twoway scatter job phd, jitter(3) ///
    title("Job Prestige and Prestige of PhD Program")

Plotting across groups of a nominal variable

  • Relationships across groups of a categorical variable can be examined
    • You can overlay multiple plots with || and note different categories via different markers on the plot
    • Alternatively, you can use the by(groupvar) option for side-by-side plots separately for each group
twoway scatter job phd if fem == 0, jitter(3) || ///
    scatter job phd if fem == 1, jitter(3) || , ///
    title("Job Prestige and Prestige of PhD Program") ///
    legend(order(1 "Men" 2 "Women"))

twoway scatter job phd, by(fem) jitter(3) ///
    title("Job Prestige and Prestige of PhD Program")

Lowess

  • To fit a lowess use lowess DV IV
    • Specify the bandwidth in option bw( )
      • I generally find 0.3 or 0.4 works well
    • Sometimes useful to save the lowess line predictions (e.g. for plotting later alongside a model’s predictions)
      • Do this with the gen( ) option
lowess job phd, bw(0.4)

Higher order polynomial models

  • The lowess suggested three turns in the regression line, suggesting we try including phd1, phd2, phd3, and phd4
    • With the c. prefix and the ## Stata makes it easy to add higher order terms to the model
use "https://tdmize.github.io/data/data/regjob3", clear
(Long's data on academic jobs of biochemists \ 2009-03-13)
regress job c.phd##c.phd##c.phd##c.phd
      Source |       SS           df       MS      Number of obs   =       408
-------------+----------------------------------   F(4, 403)       =     17.38
       Model |  56.7624316         4  14.1906079   Prob > F        =    0.0000
    Residual |   329.03396       403   .81646144   R-squared       =    0.1471
-------------+----------------------------------   Adj R-squared   =    0.1387
       Total |  385.796392       407  .947902683   Root MSE        =    .90358

------------------------------------------------------------------------------
         job | Coefficient  Std. err.      t    P>|t|     [95% conf. interval]
-------------+----------------------------------------------------------------
         phd |     -7.696      4.939   -1.558   0.120      -17.405       2.014
             |
 c.phd#c.phd |      5.166      2.801    1.844   0.066       -0.341      10.672
             |
 c.phd#c.phd#|
       c.phd |     -1.333      0.668   -1.996   0.047       -2.646      -0.020
             |
 c.phd#c.phd#|
 c.phd#c.phd |      0.119      0.057    2.100   0.036        0.008       0.231
             |
       _cons |      5.210      3.069    1.697   0.090       -0.824      11.243
------------------------------------------------------------------------------

Local polynomial smoothing

  • lpoly in Stata fits a local polynomial smoothed line
    • By default, a scatterplot is included. I find this distracting and prefer to remove it with the noscatter option
    • I always find confidence intervals helpful. Include them with the ci option
    • The legend isn’t useful here so turn it off with option legend(off)
    • If desired, can try different options for the polynomial with the degree( ) option
      • See other options for smoothing the line in the help file
lpoly depressSS role

* I prefer to omit the scatterplot, add CIs, etc.
lpoly depressSS role, noscatter ci legend(off)

Binned scatterplots

  • The Stata user-written programs binscatter and binscatterhist are both extremely useful
    • Check out the examples in the help files—they nicely display the features of each program
*Binned scatterplot
binscatter PHscale age, ///
    title("Binned scatterplot for physical health by age") ///
    ytitle("Physical health scale") line(none)

*Add histogram
binscatterhist PHscale age, hist(PHscale age) xmin(15)

*Try a quadratic line fit
binscatterhist PHscale age, line(qfit) ///
    title("Binned scatterplot for physical health by age") ///
    subtitle("Including a quadratic line fit") ///
    ytitle("Physical health scale")

Treating continuous & ordinal IVs as nominal

  • Enter the continuous variable as nominal in the model with i.
    • Then, plot predictions as if IV was entered as continuous
reg MHscale i.roles ///
    i.male c.hhincome c.age i.educ i.wave, cluster(ID)

qui margins, at(roles=(0(1)8))

marginsplot, recast(scatter) ciopts(color(gs10)) ///
    title("Predictions from Categorical Form of Social Roles") ///
    ytitle("Mental Health Scale") xtitle("Number of Social Roles")
  • If IV has many categories, bin it first (egen… cut( ) is useful)
    • E.g. To make age categories in 10-year increments (e.g. 20-30, 31-40, 41-50, etc.)
egen age_bins = cut(age), at(20(10)100)

Ordinal independent variables

  • The strategy for assessing the linearity of the effects of ordinal IVs is the same as the strategy for treating a continuous variable as nominal
    • Here, do not want to bin—want to keep all values and enter ordinal IV as c.var in one model and i.var in the next model
    • Then, make predictions with margins and plot as usual with marginsplot
reg stressSS c.health
margins, at(health=(1(1)5))
marginsplot, ///
    ylab(-1(.5)1) ytitle("Stress scale score") ///
    title("Overall stress levels by health self-assessment", size(medlarge)) ///
    subtitle("Linear (continuous) IV specification of health")

reg stressSS i.health
margins, at(health=(1(1)5))
marginsplot, ///
    ylab(-1(.5)1) ytitle("Stress scale score") ///
    title("Overall stress levels by health self-assessment", size(medlarge)) ///
    subtitle("Nominal IV specification of health")
Back to top