Decomposing the difference between two means

// GSS 1990-2004
use year race sibs reg16 educ maeduc if inrange(year, 1990, 2004) ///
  & year != 2002 ///  // 2002 measures race in a non-standard way
    using "gss7221_r1.dta", clear

// In Table 7.8 Treiman states that he's got 17,090 cases (14,985 non-black + 
// 2,105 black). There seems to be no way I can get there with the data
// because maeduc, mother's education, only has 15,996 valid observations.

// The reason for this seems to be that in his do-file, he starts with
// the 1990 GSS file and then appends all other GSS year files -- including
// the 1990 file, thus including the 1990 file twice. Running the do-file with

// expand 2 if year == 1990

// uncommented (almost) replicates the numbers reported in the book.

drop year // Not needed

// Race variables
generate black    = (race == 2)
generate nonblack = !black
drop race

// Truncate number of siblings at 15
replace sibs = 15 if sibs > 15  & !missing(sibs)
label var sibs "Sibsize"

// South
gen south = (inrange(reg16, 5, 7))
label var south "Southern origin"
drop reg16

// Education
label var educ "Education"
label var maeduc "Mother's education"

// Listwise deletion
mark touse if !missing(educ, maeduc, sibs)
keep if touse



// Table 7.8a
eststo clear
local vlist educ maeduc sibs south 
local upper

local lower `vlist'

foreach v of local vlist {
  estpost correlate `v' `lower' if !black
  local nnonblack = e(N)
  foreach m in b rho p count {
   matrix `m' = e(`m')
  }
  if "`upper'"!="" {
   estpost correlate `v' `upper' if black
   local nblack = e(N)
   foreach m in b rho p count {
    matrix `m' = e(`m'), `m'
   }
  }
  ereturn post b
  foreach m in rho p count {
   quietly estadd matrix `m' = `m'
  }
  eststo `v'
  local lower: list lower - v
  local upper `upper' `v'
}

// Table 7.8b
esttab using table2.tex, replace nonumbers noobs mtitles not booktabs label nostar ///
       title(Correlations between study variables: Blacks (\emph{N} = `nblack') above, non-blacks (\emph{N}  = `nnonblack') below the diagonal) ///
mtitle("Education" "Mother's education" "Sibsize" "Southern origin")

eststo clear
eststo: estpost sum educ maeduc sibs south if black
eststo: estpost sum educ maeduc sibs south if !black

esttab using table1.tex, booktabs replace ///
       cells("mean(label(Mean) fmt(2))" "sd(label(SD) fmt(2) par)") ///
       label mtitle("Blacks" "Non-Blacks") ///
       title("Means and standard deviations of study variables")



// Table 7.9 
eststo clear
eststo: regress educ maeduc sibs south if touse & black
eststo: regress educ maeduc sibs south if touse & !black

esttab using table3.tex, booktabs replace ///
       b(2) se(2) r2(2) ///
       label mtitle("Blacks" "Non-Blacks") ///
       title("Coefficients of a model of years of schooling, for blacks and non-blacks, US adults, 1990--2004")



// Table 7.10
eststo clear
eststo: oaxaca educ maeduc sibs south if touse, by(black) detail noisily
eststo: oaxaca educ maeduc sibs south if touse, by(nonblack) detail noisily
esttab using table4.tex, booktabs replace ///
       b(2) nose not label mtitle("Blacks as reference" "Non-Blacks as reference") ///
       eqlabel("\emph{Overall}" ///
               "\emph{Differences in assets}" ///
               "\emph{Differences in returns to assets}" ///
               "\emph{Interactions}") ///
       coeflabels(overall:difference "Difference in years of schooling" ///
                 overall:endowments "Total due to difference in assets" ///
                 overall:coefficients "Total due to difference in returns" ///
                 overall:interaction "Total due to interactions") ///
       drop(overall:group*) ///
       varwidth(30) alignment(D{.}{.}{-1}) ///
       title("Decomposition of the difference in the mean years of schooling by blacks and non-blacks, US adults, 1990--2004")

Expectation Maximization (EM) for missing values using Stata

The code below allows replicating the analyses from Allison (2002, pp. 21-3).

use "https://statisticalhorizons.com/wp-content/uploads/college.dta", clear

// Table 4.1


eststo clear
estpost summarize gradrat csat lenroll private stufac rmbrd act

esttab using test.tex, cells("count(label(Nonmissing cases)) mean(label(Mean) fmt(2)) sd(label(SD) fmt(2))") ///
        nomtitle nonumber ///
        title(Descriptive Statistics for College Data Based on Available Cases) ///
        booktabs replace

// Tabe 4.2

eststo clear
eststo: regress gradrat csat lenroll private stufac rmbrd

#delimit ;
esttab using test.tex, cells("b(fmt(3) label(Coefficient))
                              se(fmt(3) label(Standard Error)) 
                              t(fmt(2) label(t Statistic)) 
                              p(fmt(4) label(p Value))")
                       order(_cons) coeflabel(_cons "Intercept")
         nomtitle nonumber
                       title(Regression that predicts GRADRAT Using Listwise Deletion) 
         booktabs append ;
#delimit cr

// EM imputation
mi set mlong
mi register imputed gradrat csat lenroll private stufac rmbrd act 
mi impute mvn gradrat csat lenroll private stufac rmbrd act, emonly
matrix m = r(Beta_em)' // Transpose matrix of imputed means
matrix C = corr(r(Sigma_em)) // Matrix of correlations
matrix variances = diag((vecdiag(r(Sigma_em)))) // Matrix of variances
matrix sds = vecdiag(cholesky(variances))' // Vector of standard deviations
matrix descriptives = m, sds // Matrix needed for Table 4.3

// Table 4.3 

esttab matrix(descriptives, fmt(2 2)) using test.tex, ///
       nomtitle title("Means and Standard Deviations from the EM Algorithm") ///
       booktabs append

// Table 4.4

esttab matrix(C, fmt(3 3)) using test.tex, ///
       nomtitle title("Correlations from the EM Algorithm") ///
       booktabs append

// Table 4.5
drop *                                         // Get rid of data but not matrices
ssd init gradrat csat lenroll private stufac rmbrd act 
ssd set observations 1302
ssd set means (stata) m
ssd set sd (stata) sds
ssd set corr (stata) C

eststo clear

eststo: sem (gradrat <- csat lenroll private stufac rmbrd) 

#delimit ;
esttab using test.tex, cells("b(fmt(3) label(Coefficient))
               se(fmt(3) label(Standard Error)) 
               t(fmt(2) label(t Statistic)) 
      p(fmt(4) label(p Value))")
     order(_cons) coeflabel(_cons "Intercept")
  nomtitle nonumber title(Regression that predicts GRADRAT Based on the EM Algorithm)
  keep(gradrat:) eqlabels("", none) // Removes equation label
  booktabs append 
  ;
#delimit cr

Reference

Allison, Paul D. 2002. Missing Data. Sage. doi: 10.4135/9781412985079

Random graphs (123): Descriptive statistics

// Open Allbus 2010
use v11 v13 v15 v17 v19 using ZA4612_v1-0-1.dta, clear

// Prepare variables
recode v11 v13 v15 v17 v19 (99 = .)
label define x 1 `""1" "Unimportant""' 7 `""7 Very" "important"'

label var v11 "v11: Job security"
label var v13 "v13: High income"
label var v15 "v15: Opportunities for promotion"
label var v17 "v17: Respectable occupation"
label var v19 "v19: Lots of leisure time"

// Plot histograms
foreach x of varlist v11 v13 v15 v17 v19 {
  label val `x' x
  histogram `x', frequency ///
                 xlabel(1 (1) 7, val) ///
                 ylabel(0 (500) 1500) ///
                 ytick(0 (250) 1750) ///
                 name(`x', replace) nodraw
}

// Calculate correlation matrix
cor v11-v19
matrix correlations = r(C)

// Plot correlation matrix
plotmatrix, mat(correlations) split(-1(0.10)1) color(sienna) ///
            freq formatcells(%6.2f) nodiag legend(title(Pearson's {it:r}) col(1)) ///
            name(cormatrix, replace) nodraw 

// Combine plots
graph combine v11 v13 v15 v17 v19 cormatrix, col(2) ysize(8) xsize(6) altshrink

Random graphs (110): Correlation heatmap

// Open Allbus 2010 data
use "ZA4612_v1-0-1.dta", clear

// Men 25-64 only
keep if v298 == 1
keep if inrange(v301, 20, 64)

// Combine top bottom variables into one 
// (Seems pointless, though, as education in years
//  was only measured in one of the splits)
generate topbot = v933 if v933 <= 10
replace  topbot = v827 if missing(topbot)
replace  topbot =    . if topbot > 10

// Education, truncated to 5th-95th percentile
qui sum v917 if v917 < 94, detail 
generate education = v917 if v917 < 94
replace  education = r(p95) if education >=  r(p95) & !missing(education)
replace  education = r(p95) if education >=  r(p95) & !missing(education)
replace  education =  r(p5) if education <=   r(p5) & !missing(education)

// ISEI
generate isei88 = v358 if !inlist(v358, 0, 99)
replace  isei88 = v392 if !inlist(v392, 0, 99)

// Household income, logged
generate income = v674 if !inlist(v674, 99997, 99999)
generate lnincome = log(income)

// Generate correlation matrix
cor topbot education income isei
matrix correlations = r(C)

// Plot
plotmatrix, mat(correlations) split(-1(0.10)1) color(vermillion) ///
            freq formatcells(%6.2f) nodiag legend(title(Pearson's {it:r}) col(1))

Replicating simple analyses with -corr2data-

A correlation table is a pivotal feature of academic papers in many disciplines such as psychology. The reason for this is that many basic statistical methods rely on correlation matrices only. Thus, the vector of means and standard deviations plus the correlation matrix and the sample size are sufficient to replicate many standard analyses such as OLS regression or principal component analysis. -corr2data- allows to do so by creating an artficial data set based on the correlation matrix (or, covariance matrix), means, standard deviations, and sample size, which can then be used for replicating or modifying the analyses.

Table 1 of Huff-Corzine et al. (1986)

In order to try out -corr2data-, we will replicate a randomly chosen research paper, namely Huff-Corzine et al. (1986). Huff-Corzine and colleagues deal with the North–South differences in homicide in the US, and, more importantly for our purposes, they provide a full correlation table and conduct OLS regression analyses.


First, we need to read in the correlation matrix as reported in Figure 1. (N.B. Stata 12 now comes with the -ssd- command which might allow for a more convenient way of reading in summary statistics.)

version 12
clear
capture which estout
if _rc ssc install estout

// Read in correlation table
#delimit ;
input str30 varnames double (a1-a9);
"1 Homicide rate"
    1.00  0.91  0.83  0.88  0.43 -0.06 -0.30  0.65  0.83;
"2 Structural poverty index"
    0.91  1.00  0.78  0.90  0.40 -0.04 -0.17  0.65  0.84;
"3 Southernness index"
    0.83  0.78  1.00  0.72  0.34  0.14 -0.44  0.71  0.93;
"4 Perc nonwhite"
    0.88  0.90  0.72  1.00  0.32  0.08 -0.14  0.68  0.76;
"5 Perc ages 20-34"
    0.43  0.40  0.34  0.32  1.00 -0.14 -0.21 -0.09  0.32;
"6 Perc rural"
   -0.06 -0.04  0.14  0.08 -0.14  1.00  0.14  0.28  0.10;
"7 Hospital beds/100K"
   -0.30 -0.17 -0.44 -0.14 -0.21  0.14  1.00 -0.20 -0.30;
"8 Gini index"
    0.65  0.65  0.71  0.68 -0.09  0.28 -0.20  1.00  0.73;
"9 Perc born in South"
    0.83  0.84  0.93  0.76  0.32  0.10 -0.30  0.73  1.00;
end;
#delimit cr

// Convert to a matrix called M
mkmat a1-a9, matrix(M) 
matrix list M

// Read in variable names
levelsof(varnames), local(names)

local ednames = ""
  // Create local "ednames" that includes all variable names,
  // use function -strtoname()- to bring variable names into 
  // proper Stata format
foreach x of local names {
  *di strtoname("`x'", 1)
  local ed_name1 = strtoname("`x'", 1)
  *di "`ed_name1 '"
  local ed_name2 = substr("`ed_name1'", 4, .)
  *di "`ed_name2 '"
  local ednames = "`ednames'" + " " + "`ed_name2'"
}

In this step, we create the artificial data set based on the correlation matrix and means and s.d.'s:

// Create data based on correlation matrix M, vectors of
// means and sd's, specify no. of cases and variable names
corr2data `ednames', n(48) clear ///
     corr(M) ///
    means(7.12 7.12 17.71 10.42 19.75 34.22 7.81 .38 7.18) ///
   sds(4.25 3.88  9.22  8.85  1.27 14.37 1.57 .02 1.44)

One shouldn't underestimate the importance of fine-looking labels in a data set:

// Create nice-looking value labels based on the variable names
local varlabels = ""

foreach var of varlist _all {
  local varname = "`var'"
  *di "`varname'"
  local varname = subinstr("`varname'", "_"   , " ", .)
  *di "`varname'"
  local varname = subinstr("`varname'", "Perc", "%", .)
  *di "`varname'"
  lab variable `var'  "`varname'"
  local varlabels = "`varlabels'" + " " + "`varname'"
}

A brief check shows that every thing seems to be alright:

corr _all

tabstat _all, statistics(mean sd)

However, replicating the analyses yields different results:

eststo clear

// Run regression

eststo: regress Homicide_rate ///
          Structural_poverty_index Southernness_index ///
    Perc_nonwhite Perc_ages_20_34 Perc_rural ///
    Hospital_beds_100K Gini_index, beta  
    
eststo: regress Homicide_rate ///
          Structural_poverty_index ///
    Perc_nonwhite Perc_ages_20_34 Perc_rural ///
    Hospital_beds_100K Gini_index Perc_born_in_South, beta
    
estadd beta: est1 est2

esttab, r2 nonumbers mlabels("South index" "% born South") ///
        cells ((b(fmt(3) star) beta(fmt(3)))) label not)

Download the do-file here.

These results differ substantially from those reported by Huff-Corzine and others (see below), however I have a feeling that I should rather be trusting my own result.

Table 2 of Huff-Corzine et al. (1986)

A robustness check offered by -corr2data- is the -seed(#)- option, which allows for generating different artificial data sets. Specifying different seeds and comparing the results can serve as a check whether the summary statistics are really sufficient for replicating an analysis. Using different seeds in our example always leads to the same result.

Reference

Huff-Corzine, Lin, Jay Corzine, and David C. Moore. 1986. "Southern Exposure. Deciphering the South's Influence on Homicide Rates." Social Forces 64(4):906-924. doi: 10.1093/sf/64.4.906

Publication-style correlation tables in Stata

Quite a number of Stata users have engaged in programming commands to create readymade correlation tables to use in publications. Stata's own -correlate- and -pwcorr- lack many desirable features; for instance, means and standard deviations cannot be included in the table automatically and columns cannot be numbered. Furthermore, variable labels cannot be displayed, only variable names can be shown in the table.

sysuse auto
correlate price-foreign, means

-estpost-, -mkcorr-, -corrtab-, and -makematrix- are all user-written commands which aim to improve the Stata default table.

 

-estpost-

-estpost- stems from Ben Jann's mighty -estout- package. The general use of the command is as follows:

sysuse auto
capture which estout
if _rc ssc install estout
estpost correlate price-trunk, matrix
esttab, unstack not nonum compress noobs

In order to display both the upper and the lower triangle of the matrix, the -nohalf- option of -estpost- needs to be used:

estpost correlate price-trunk, matrix nohalf

However, getting -esttab- to use variable labels or to include descriptive statistics doesn't seem to be a trivial task, although the latter can be done. Finally, it should be noted that the default of -estpost correlate- is pairwise deletion, not listwise deletion; so that it is -pwcorr- rather than -correlate-. In order to obtain a correlation matrix based on listwise deletion, the -listwise- option needs to be specified.

 

-mkcorr-

-mkcorr- by Glenn Hoetker comes with many nice features and only few drawbacks.

sysuse auto, clear
capture which mkcorr
if _rc ssc install mkcorr
mkcorr price-trunk, log(corr_table)

It allows including summary statistics (mean, SD, minimum, and maximum) into the table:

mkcorr price-trunk, log(corr_table) replace means

The use of variable labels:

mkcorr price-trunk, log(corr_table) replace lab

The use of numbers in the column headers:

mkcorr price-trunk, log(corr_table) replace num

The inclusion of p-values:

mkcorr price-trunk, log(corr_table) replace means sig 

Furthermore, the number of decimal places for correlation can be manipulated via -cdec(#)- (-mdec(#)- for summary statistics).
The main drawback is that the table is not being displayed in the results window; instead it can only be written into a tab-separated text-only file that then needs additional formatting in Word or Excel. Also, the p-values are rather ugly; an option to include stars to denote statistical significance would be a great improvement. Again, the default setting is to provide correlations based on pairwise deletion; the option -casewise- will yield results based on listwise deletion.

-corrtab-

-corrtab- by Fred Wolfe is a tool that has a somewhat more limited functionality, or, more functions that I don't find to be of great use.

sysuse auto, clear
capture which corrtab
if _rc ssc install corrtab
corrtab price-trunk 

-corrtab- also allows surpressing correlations based on their p-values when using the option -print(#)-:

corrtab price-trunk, print (.10)

Or, based on the absolute value of the correlation coefficient when specifying -above(#)-:

corrtab price-trunk, above(.4)

Correlations can also be surpressed based on the position in the list of variables. The option -var(#)- only lists the first # variables in the columns, later variables are only included in the rows of the table:

corrtab price-trunk, var(3)

Correlations can be sorted by size for a single variable specified in -vsort()-:

corrtab price-trunk, vsort(price)

In order to use variable labels, some complicated additional commands need to be used; and the default for treating missing values is also pairwise deletion.
In sum, there seems to be little reason to ever use -corrtab-.

 

-makematrix-

-makematrix- by Nick Cox is a very flexible tool for all sort of applications, one of them being correlation tables. He explains it a bit in Cox (2003), but it seems to require too much tweaking to get readymade publication-style correlation tables in a short amount of time. Furthermore, there's a bug in the command: The option -listwise- will yield a pairwise correlation matrix instead of a correlation based on listwise deletion.

In sum, the commands reviewed here all leave a lot to be desired; -mkcorr- is the command that seems to be the most useful one to me.

References

Cox, Nicholas J. 2003. "Speaking Stata. Problems with Tables, Part II." Stata Journal 3(4):420-439.