A2.7 Assessing nonlinearities
- To examine the bivariate relationship between two variables, scatterplots can be helpful
- Note this is most helpful when the relationship between variables is strong and when you have a relatively small dataset (n.b. can plot random subset of data when N is large)
- In Stata: twoway scatter var1 var2
- var1 will be on the y-axis; var2 on the x-axis
- Often, you will want to “jitter” the plot so that observations do not overlap
- This adds random noise to where exactly the observation is plotted on the scatterplot (see the following two examples)
- This is helpful to avoid falsely concluding there are few observations at places with lots of overlapping data
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
- Specify the bandwidth in option bw( )
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")