*! version 4.0 May2023

Title

binscatterhist -- Binned scatterplots with variables distribution

Syntax

binscatterhist
varlist [if] [in]
[weight]
[, options]

where varlist is

y_1 [y_2 [...]] x

options

Main
binsreg (recommended) creates the scatterplot points by running the program binsreg in the background (binsreg), for a correct data visualization in the presence of covariates

by(varname) plot separate series for each group (see important notes below)

medians plot within-bin medians instead of means

legend works normally when a histogram is not included. With histogram, has to take the form legend(order(1 ".." 2 ".."))

Bins
nquantiles(#) number of equal-sized bins to be created; default is 20

genxq(varname) generate quantile variable containing the bins

discrete each x-value to be used as a separate bin

xq(varname) variable which already contains bins; bins therefore not recomputed

Residuals Computation
regtype(regtype) can be reghdfe or areg, requires absorb() to be specified; default is reghdfe

SE/Robust
cluster(varname) clustered s.e. affect both the sample for residualization and the computation of s.e. for slope reporting

vce(robust) robust s.e. affect both the sample for residualization and the computation of s.e. for slope reporting

Controls
controls(varlist) residualize the x & y variables on controls before plotting

absorb(varlist) residualize the x & y variables on a categorical variable. Can include only one variable, if regtype(areg) is specified

noaddmean do not add the mean of each variable back to its residuals

Fit Line
linetype(binscatterhist##linetype:linetype) type of fit line; default is lfit, may also be qfit, connect, or none

rd(numlist) create regression discontinuity at x-values

reportreg display the regressions used to estimate the fit lines

Coefficient and Sample Reporting
coefficient(#) reports the slope of the fitted line with its standard error, rounded at # using round(coef,#). See round

ci(#) reports the #% confidence interval, rounded as the coefficient()

pvalue reports the pvalue of the regression on residualized variables

sample reports the sample size of the regression on residualized variables

stars(stars) signals p-value using stars; stars can be: nostars, 1 (*5% **1%), 2 (+10% *5% **1%), 3 (+10% *5% **1% ***0.1%), 4 (*5% **1% ***0.1%); default is stars(1)

boxposition(x y) sets x and y position for the text box

Graph Style
colors(colorstylelist) ordered list of colors

mcolors(colorstylelist) overriding ordered list of colors for the markers

lcolors(colorstylelist) overriding ordered list of colors for the lines

msymbols(symbolstylelist) ordered list of symbols

msize(markersizestylelist) size of markers

twoway_options titles, legends, axes, added lines and text,
regions, name, aspect ratio, etc.

Background scatter
scatter(bins) if bins="N", it creates a scatterplot of the graphed variables, otherwise, a binned scatterplot with the specified number of bin.

scattermsymbols(symbolstylelist) ordered list of symbols for background scatter

scattermcolors(colorstylelist) overriding ordered list of colors for the markers of the background scatter

scattermsize(markersizestylelist) size of markers for background scatter

Histogram
histogram(varlist) plots a histogram for earliest vars (max=2). Selected variables have to be the scattered ones

xmin(value) sets the base position of the y histogram in terms of xaxis

ymin(value) sets the base position of the x histogram in terms of yaxis

xhistbarheight(value) sets the height of the x histogram as a percentage; default is xhistbarheight(10)

yhistbarheight(value) sets the height of the y histogram as a percentage; default is yhistbarheight(10)

xhistbarwidth(value) sets the width of the x histogram as a percentage; default is xhistbarwidth(100)

yhistbarwidth(value) sets the width of the y histogram as a percentage; default is yhistbarwidth(100)

xhistbins(#) number of bins to be created; default is 20

yhistbins(#) number of bins to be created; default is 20

Following options require axis to be specified: x or y, e.g. xcolor() ycolor()
axiscolor(colorstyle) outline and fill color and opacity; default is teal%50 for xcolor() and maroon%50 for ycolor()

axisfcolor(colorstyle) fill color and opacity

axisfintensity(intensitystyle) fill intensity

axislcolor(colorstyle) outline color and opacity

axislwidth(linewidthstyle) thickness of outline

axislpattern(linepatternstyle) outline pattern (solid, dashed, etc.)

axislalign(linealignmentstyle) outline alignment (inside, outside, center)

axislstyle(linestyle) overall look of outline

axisbstyle(areastyle) overall look of bar, all settings above

axispstyle(pstyle) overall plot style, including areastyle

Save Output
savegraph(filename) save graph to file; format automatically detected from extension [ex: .gph .jpg .png]

savedata(filename) save filename.csv containg scatterpoint data, and filename.do to process data into graph

replace overwrite existing files

fastxtile options
nofastxtile use xtile instead of fastxtile

randvar(varname) use varname to sample observations when computing quantile boundaries

randcut(#) upper bound on randvar() used to cut the sample; default is randcut(1)

randn(#) number of observations to sample when computing quantile boundaries

aweights and fweights are allowed;
see weight.

Description

binscatterhist generates binned scatterplots, with the option to plot the variables underlying distribution and report estimation results.

Binned scatterplots provide a non-parametric way of visualizing the relationship between two variables,
binscatterhist adds several features to the popular program used to produce such scatterplots, binscatter.
Binscatterhist uses by default reghdfe to calculate residuals, therefore allowing for multiple fixed effects, but keeps areg as an alternative option.
Binscatterhist allows, in addition, to create histograms of the scattered variables, and include them in the
graph for a more complete representation of the data. Finally, it allows for an automatic reporting of the slope and standard error of the fitted
line, with options for robust and clustered standard errors.
As binscatter, binscatterhist solves the problem of scatterplots
with a large number of observations: it groups the x-axis variable into equal-sized bins, computes the
mean of the x-axis and y-axis variables within each bin, then creates a scatterplot of these data points.

Binscatterhist keeps as a base, the same options as binscatter: it provides built-in options to control for covariates before plotting the relationship
(see Controls), and will plot fit lines based on the underlying data, and can automatically
handle regression discontinuities (see Fit Line).

Options

Main

binsreg (recommended) creates the scatterplot points by running the program binsreg in the background (binsreg), for a correct data visualization in the presence of covariates.

by(varname) plots a separate series for each by-value. Both numeric and string by-variables
are supported, but numeric by-variables will have faster run times.

Users should be aware of the two ways in which binscatterhist does not condition on by-values:

1) When combined with controls() or absorb(), the program residualizes using the restricted model in which each covariate
has the same coefficient in each by-value sample. It does not run separate regressions for each by-value. If you wish to control for
covariates using a different model, you can residualize your x- and y-variables beforehand using your desired model then run binscatterhist
on the residuals you constructed.

2) When not combined with discrete or xq(), the program constructs a single set of bins
using the unconditional quantiles of the x-variable. It does not bin the x-variable separately for each by-value.
If you wish to use a different binning procedure (such as constructing equal-sized bins separately for each
by-value), you can construct a variable containing your desired bins beforehand, then run binscatterhist with xq().

medians creates the binned scatterplot using the median x- and y-value within each bin, rather than the mean.
This option only affects the scatter points; it does not, for instance, cause linetype(lfit)
to use quantile regression instead of OLS when drawing a fit line.

legend works normally when a histogram is not included. When a histogram is included, it has to take the form legend(order(1 ".." 2 "..")). This is is due to a compatibility issue in the auxiliary package addplot.

Bins

nquantiles(#) specifies the number of equal-sized bins to be created. This is equivalent to the number of
points in each series. The default is 20. If the x-variable has fewer
unique values than the number of bins specified, then discrete will be automatically invoked, and no
binning will be performed.
This option cannot be combined with discrete or xq().

Binning is performed after residualization when combined with controls() or absorb().
Note that the binning procedure is equivalent to running xtile, which in certain cases will generate
fewer quantile categories than specified. (e.g. sysuse auto; xtile temp=mpg, nq(20); tab temp)

genxq(varname) creates a categorical variable containing the computed bins.
This option cannot be combined with discrete or xq().

discrete specifies that the x-variable is discrete and that each x-value is to be treated as
a separate bin. binscatterhist will therefore plot the mean y-value associated with each x-value.
This option cannot be combined with nquantiles(), genxq() or xq().

In most cases, discrete should not be combined with controls() or absorb(), since residualization occurs before binning,
and in general the residual of a discrete variable will not be discrete.

xq(varname) specifies a categorical variable that contains the bins to be used, instead of binscatterhist generating them.
This option is typically used to avoid recomputing the bins needlessly when binscatterhist is being run repeatedly on the same sample
and with the same x-variable.
It may be convenient to use genxq(binvar) in the first iteration, and specify xq(binvar) in subsequent iterations.
Computing quantiles is computationally intensive in large datasets, so avoiding repetition can reduce run times considerably.
This option cannot be combined with nquantiles(), genxq() or discrete.

Care should be taken when combining xq() with controls() or absorb(). Binning takes place after residualization,
so if the sample changes or the control variables change, the bins ought to be recomputed as well.

Residuals Computation

regtype(regtype) can be reghdfe or areg, requires absorb() to be specified; default is reghdfe. When reghdfe is specified, absorb()
allows for more than one carname, however, interactions are not allowed, including tricks like e.g. one##control, to include controls in the absorb. Such controls
must be included in the controls() option. reghdfe drops singleton observations w.r.t. the included fixed effects, therefore sample size might differ between
reghdfe and areg.

SE/Robust

cluster(varname) clustered s.e. affect both the sample for residualization and the computation of s.e. for slope reporting.

vce(robust) robust s.e. affect both the sample for residualization and the computation of s.e. for slope reporting.

Controls

controls(varlist) residualizes the x-variable and y-variables on the specified controls before binning and plotting.
To do so, binscatterhist runs a regression of each variable on the controls, generates the residuals, and adds the sample mean of
each variable back to its residuals.

absorb(varname) absorbs fixed effects in the categorical variable from the x-variable and y-variables before binning and plotting,
To do so, binscatterhist runs an areg of each variable with absorb(varname) and any controls() specified. It then generates the
residuals and adds the sample mean of each variable back to its residuals.

noaddmean prevents the sample mean of each variable from being added back to its residuals, when combined with controls() or absorb().

Fit Line

linetype(binscatterhist##linetype:linetype) specifies the type of line plotted on each series.
The default is lfit, which plots a linear fit line. Other options are qfit for a quadratic fit line,
connect for connected points, and none for no line.

Linear or quadratic fit lines are estimated using the underlying data, not the binned scatter points. When combined with
controls() or absorb(), the fit line is estimated after the variables have been residualized.

rd(numlist) draws a dashed vertical line at the specified x-values and generates regression discontinuities when combined with line(lfit|qfit).
Separate fit lines will be estimated below and above each discontinuity. These estimations are performed using the underlying data, not the binned scatter points.

The regression discontinuities do not affect the binned scatter points in any way.
Specifically, a bin may contain a discontinuity within its range, and therefore include data from both sides of the discontinuity.

reportreg displays the regressions used to estimate the fit lines in the results window.

Coefficient and Sample Reporting

coefficient reports the slope of the fitted line with its standard error.

ci(#) reports the #% confidence interval, rounded as the coefficient().

pvalue reports the pvalue of the regression on residualized variables.

sample reports the sample size of the regression on residualized variables.

stars(stars) signals p-value using stars; stars can be: nostars, 1 (*5% **1%), 2 (+10% *5% **1%), 3 (+10% *5% **1% ***0.1%), 4 (*5% **1% ***0.1%); default is stars(1).

boxposition(x y) sets x and y position for the text box (both have to be set)

Graph Style

colors(colorstylelist) specifies an ordered list of colors for each series

mcolors(colorstylelist) specifies an ordered list of colors for the markers of each series, which overrides any list provided in colors()

lcolors(colorstylelist) specifies an ordered list of colors for the line of each series, which overrides any list provided in colors()

msymbols(symbolstylelist) specifies an ordered list of symbols for each series

msize(markersizestylelist) size of markers for background scatter

twoway_options:

Any unrecognized options added to binscatterhist are appended to the end of the twoway command which generates the
binned scatter plot.

These can be used to control the graph titles,
legends, axes, added lines and text,
regions, name, aspect ratio, etc.

Background scatter

scatter(bins) if bins="N", it creates a scatterplot of the graphed variables, otherwise, a binned scatterplot with the specified number of bins

scattermsymbols(symbolstylelist)ordered list of symbols for background scatter

scattermcolors(colorstylelist) overriding ordered list of colors for the markers of the background scatter

scattermsize(markersizestylelist) size of markers for background scatter

Histogram

histogram(varnames) plots a histogram for each of the selected variables (max=2). Selected variables have to be the scattered ones. Requires xmin() if yvar is plotted and ymin() if xvar() is plotted

xmin(value) set as minimum value of the xaxis in the binscatterhist without histogram, it correctly positions a ygraph

ymin(value) set as minimum value of the yaxis in the binscatterhist without histogram, it correctly positions a xgraph

xhistbarheight(value) sets the height of the x histogram as a percentage of the binscatterhist main area height; default is xhistbarheight(10)

yhistbarheight(value) sets the height of the y histogram as a percentage of the binscatterhist main area height; default is yhistbarheight(10)

xhistbarwidth(value) sets the width of the x histogram as a percentage; default is xhistbarwidth(100)

yhistbarwidth(value) sets the width of the y histogram as a percentage; default is xhistbarwidth(100)

xhistbins(#) number of bins to be created; default is 20

yhistbins(#) number of bins to be created; default is 20

Following options require axis to be specified: x or y, e.g. xcolor() ycolor()

axiscolor(colorstyle) outline and fill color and opacity; default is teal%50 for xcolor and maroon%50 for ycolor

axisfcolor(colorstyle) fill color and opacity

axisfintensity(intensitystyle) fill intensity

axislcolor(colorstyle) outline color and opacity

axislwidth(linewidthstyle) thickness of outline

axislpattern(linepatternstyle) outline pattern (solid, dashed, etc.)

axislalign(linealignmentstyle) outline alignment (inside, outside, center)

axislstyle(linestyle) overall look of outline

axisbstyle(areastyle) overall look of bar, all settings above

axispstyle(pstyle) overall plot style, including areastyle

Save Output

savegraph(filename) saves the graph to a file. The format is automatically detected from the extension specified [ex: .gph .jpg .png],
and either graph save or graph export is run. If no file extension is specified .gph is assumed.

savedata(filename) saves filename.csv containing the binned scatterpoint data, and filename.do which
loads the scatterpoint data, labels the variables, and plots the binscatterhist graph.

Note that the saved result e(cmd) provides an alternative way of capturing the binscatterhist graph and editing it.

replace specifies that files be overwritten if they alredy exist

fastxtile options

nofastxtile forces the use of xtile instead of fastxtile to compute bins. There is no situation where this should
be necessary or useful. The fastxile program generates identical results to xtile, but runs faster on large datasets, and has
additional options for random sampling which may be useful to increase speed.

fastxtile is built into the binscatterhist code, but may also be installed
separately from SSC (click here to install) for use outside of binscatterhist.

randvar(varname) requests that varname be used to select a
sample of observations when computing the quantile boundaries. Sampling increases
the speed of the binning procedure, but generates bins which are only approximately equal-sized
due to sampling error. It is possible to omit this option and still perform random sampling from U[0,1]
as described below in randcut() and randn().

randcut(#) specifies the upper bound on the variable contained
in randvar(varname). Quantile boundaries are approximated using observations for which
randvar() <= #. If no variable is specified in randvar(),
a standard uniform random variable is generated. The default is randcut(1).
This option cannot be combined with randn().

randn(#) specifies an approximate number of observations to sample when
computing the quantile boundaries. Quantile boundaries are approximated using observations
for which a uniform random variable is <= #/N. The exact number of observations
sampled may therefore differ from #, but it equals # in expectation. When this option is
combined with randvar(varname), varname ought to be distributed U[0,1].
Otherwise, a standard uniform random variable is generated. This option cannot be combined
with randcut().

Examples of main features (Stepner, 2013)

Load the 1988 extract of the National Longitudinal Survey of Young Women and Mature Women.

. sysuse nlsw88

. keep if inrange(age,35,44) & inrange(race,1,2)

What is the relationship between job tenure and wages?

. scatter wage tenure

. binscatterhist wage tenure

The scatter was too crowded to be easily interpetable. The binscatterhist is cleaner, but a linear fit looks unreasonable.

Try a quadratic fit.

. binscatterhist wage tenure, line(qfit)

We can also plot a linear regression discontinuity.

. binscatterhist wage tenure, rd(2.5)

What is the relationship between age and wages?

. scatter wage age

. binscatterhist wage age

The binscatterhist is again much easier to interpret. (Note that binscatterhist automatically
used each age as a discrete bin, since there are fewer than 20 unique values.)

How does the relationship vary by race?

. binscatterhist wage age, by(race)

The relationship between age and wages is very different for whites and blacks. But what if we control for occupation?

. binscatterhist wage age, by(race) absorb(occupation)

A very different picture emerges. Let's label this graph nicely.

. binscatterhist wage age, by(race) absorb(occupation) msymbols(O T) xtitle(Age) ytitle(Hourly Wage) legend(lab(1 White) lab(2 Black))

Examples of histogram features

Load the National Longitudinal Survey of Women 1988

. webuse nlsw88, clear

The basic binscatterhist works exactly as binscatter

. binscatterhist wage tenure

Let's add the distribution of the x variable tenure

. binscatterhist wage tenure, histogram(tenure)

The default position of the x graph is not pleasant, let's fix that and add both variables distribution this time

. binscatterhist wage tenure, histogram(wage tenure) ymin(4)

We want to experiment with the look, let's try a simpler look and a smaller width

. binscatterhist wage tenure, histogram(wage tenure) ymin(4) yhistbarwidth(50) xhistbarwidth(50) ybstyle(outline) xbstyle(outline)

Let's now try some further options: increasing number of bins and height of the x and y distribution

. binscatterhist wage tenure, histogram(wage tenure) ymin(4) xhistbarheight(15) yhistbarheight(15) xhistbins(40) yhistbins(40)

Examples of further features

Let's report the estimation results, using robust standard errors and with grade fixed effects

. binscatterhist wage tenure, absorb(grade) vce(robust) coef(0.01) sample xmin(-2.2) ymin(5) histogram(wage tenure) xhistbarheight(15) yhistbarheight(15) xhistbins(40) yhistbins(40)

Let's use now areg - therefore keeping singleton fixed effects. With a negative slope, the reported coefficient and sample adjust automatically their position We further report 95% CI and pvalue

. replace tenure=-abs(tenure)

. binscatterhist wage tenure, regtype(areg) absorb(grade) vce(robust) coef(0.01) ci(95) pvalue sample xmin(-22) ymin(5) histogram(wage tenure) xhistbarheight(15) yhistbarheight(15) xhistbins(40) yhistbins(40)

In the presence of covariates, binscatterhist may not correctly visualise the data. Scatterpoints can be correctly created by running the package binsreg in the background

. replace tenure=abs(tenure)

. binscatterhist wage tenure, controls(age) histogram(wage tenure) ymin(8) coef(0.01) pvalue sample binsreg

Adding the option scatter(N) we can add a background binned scatterplot with number of bins equal to the sample size, or directly time the desired number

. binscatterhist wage tenure, absorb(grade) vce(robust) coef(0.01) ci(95) pvalue sample xmin(-10) ymin(-10) histogram(wage tenure) xhistbarheight(15) yhistbarheight(15) xhistbins(40) yhistbins(40) scatter(N)

Saved Results

binscatterhist keeps in e() the info stored by reg/areg and in r(), the following info:

Scalars

r(N) number of observations

Macros

r(graphcmd) twoway command used to generate graph, which does not depend on loaded data

Note: it is often important to reference this result using `"`e(graphcmd)'"'
rather than r(graphcmd) in order to avoid truncation due to Stata's character limit for strings.

Matrices

r(byvalues) ordered list of by-values (if numeric by-variable specified)

r(rdintervals) ordered list of rd intervals (if rd specified)

r(y#_coefs) fit line coefficients for #th y-variable (if lfit or qfit specified)

Author

Matteo Pinna

matteo.pinna@gess.ethz.ch

Acknowledgements

The author would like to thank one anonymous referee, Elliott Ash, Chistopher Baum, Suresh Naidu, Sergio Galletta, Malka Guillot, Susanna B. Berkouwer, Guohui Jiang, Jack Glaser,Fabián Enrique Gonzáles and Guillermo Cruces for useful feedback and bug finding.

The present version of binscatterhist is based on a program in Michael Stepner (2013): "BINSCATTER: Stata module to generate binned scatterplots" - https://EconPapers.repec.org/RePEc:boc:bocode:s457709; a program by Ben Jann (2014): "ADDPLOT: Stata module to add twoway plot objects to an existing twoway graph," Statistical Software Components S457917, Boston College Department of Economics, revised 28 Jan 2015 <https://ideas.repec.org/c/boc/bocode/s457917.html>; and uses the BINSREG program: Cattaneo, M. D., R. K. Crump, M. H. Farrell, and Y. Feng. 2019b. Binscatter Regressions. arXiv:1902.09615.