Random graphs (144): Automatically labeling figures with letters

// Open ESS round 5
use agea stflife cntry if inlist(cntry, "DE", "GB", "NL", "SE") using "ESS5e03_4.dta", clear

// Country variable
kountry cntry, from(iso2c)
rename NAMES_STD country

// Restrict to complete cases
keep if !missing(agea, stflife, cntry)

// Center age
qui sum agea, detail
replace agea = `r(p99)' if agea > `r(p99)' & !missing(agea)
bys country: center agea, gen(age)
drop agea
label var age "Age (centered)"

// Letters
tokenize "`c(ALPHA)'"
local i 1

levelsof country, local(country)

foreach z of local country {
  qui sum age if country == `"`z'"', detail 
  foreach x of numlist 10 25 50 75 90 {
    local x`x': di %9.1f r(p`x')
    di `x`x''
  }
  qui regress stflife c.age##c.age if country == `"`z'"'
  
  qui margins, at(age = (`x10' `x25' `x50' `x75' `x90'))
  marginsplot, title(`"{bf:``i''} `z'"') ///
               ytitle("Life satisfaction" "(predicted)") ///
               xtitle("Age (percentiles)") ///
               xlabel(`x10' "10th" ///
                      `x25' "25th" ///
                      `x50' "50th" ///
                      `x75' "75th" ///
                      `x90' "90th") ///
               ylabel(, format(%6.1f)) ///
               recastci(rarea) ciopts(color(gs10)) ///
               name("``i''", replace) nodraw
  local ++i
}

graph combine A B C D, col(2) ycommon

Random graphs (65): Dot plots with confidence intervals

preserve
// Keep those countries for which three waves are available
keep if inlist(cntry, "AT", "AU", "BG", "CZ", "DEE", "DEW", "ES", "GB", "HU") | ///
        inlist(cntry, "IE", "IL", "JP", "NL", "NO", "PH", "PL", "RU", "SE")   | ///
  inlist(cntry, "SI", "US")

// Define temporary objects
tempname foo
tempname foox
postfile `foo' str3 cntry year perc perc_ll perc_ul using `foox', replace

levelsof cntry, local(levels1)
levelsof year, local(levels2)

// Calculate proportions by country and year
foreach year of local levels2 {
 foreach country of local levels1 {

    capture proportion separate if cntry == "`country'" & year == `year'
    *matrix list r(table)
    matrix fcoefs  = r(table)
    local fperc    = fcoefs[1,2] * 100
    local fperc_ll = fcoefs[5,2] * 100
    local fperc_ul = fcoefs[6,2] * 100

    di "`country'"  _skip(2) `year' _skip(2)  `fperc_ll' _skip(2) `fperc' _skip(2) `fperc_ul'

    if "`fperc'" != "" {    // Make sure that the loop doesn't break if empty
    post `foo' ("`country'") (`year') (`fperc') (`fperc_ll') (`fperc_ul')
    }
 }
}
postclose `foo'

// Use posted data set
use `foox', clear

// Generate country name variable using -kountry-
kountry cntry, from(iso2c)
rename NAMES_STD country
replace country = "Germany (West)" if country == "dew"
replace country = "Germany (East)" if country == "dee"

// Plot average change over time
twoway (scatter perc year, connect(direct) msymbol(o)) ///
       (rcap perc_ll perc_ul year) ///
      , by(country, ///
           note("{it:Source:} ISSP 'Family and Changing Gender Roles II{c 150}IV. {it:Note:} Error bars denote 95% confidence intervals.", span) ///
           legend(off)) ///
        xlabel(1994 2002 2012) xtitle("") ///
        ytitle("% of couples keeping incomes separate") ///
        ylabel(0 10 20 30 40 50, gstyle(minor)) yscale(r(0))
        // yscale(r(0)) adds missing gridline according to
        // http://www.stata.com/statalist/archive/2013-04/msg00904.html    
restore

Random Graphs (50): Dot graphs with confidence intervals

// Create temporary object

tempname bqlrer

// Define postfile

postfile `bqlrer' count str100 commandline str100 descri str10 entity ///
                  prev prev_loci prev_hici differ differ_loci differ_hici ///
                  loed mide hied ratio_ ratio_loci ratio_hici n ///
                  using "results\results", replace

   // Count variable
local count = 0

  // Levels for loop
levelsof entity, local(entity_levels)
*di `entity_levels'

  // Loop: 1 round per entity
foreach Y of local entity_levels  {
     local entity_level = "`Y'"
     di "`entity_level'"
  *qui sum poorhealth if entity == "`Y'"
  *scalar proportion = r(mean)
  
  // Model 1: OV: PSRH lowest2, adjusted -- Prevalence rate only, education not in model
  qui logit poorhealth c_age c_female if entity == "`Y'"
  di _rc
  if _rc == 0 {
  *prev prev_loci prev_hici
  qui margins if entity == "`Y'", post
     scalar prev = 100 * _b[_cons]
     scalar prev_loci = 100 * (_b[_cons] - 1.96 * _se[_cons])
     scalar prev_hici = 100 * (_b[_cons] + 1.96 * _se[_cons])
  di `prev' 
  di `prev_loci' 
  di `prev_hici'
  }
  
  // Model 2: OV: PSRH lowest2, adjusted
     local count = `count' + 1
     local desc "DV: Poor health (bottom 2), IV: Education, Age, Sex"
     qui logit poorhealth i.education c_age c_female if entity == "`Y'"
  local commandline = e(cmdline)
  di _rc
  if _rc == 0 {
  
  // Differences Model 2
  qui margins r.education if entity == "`Y'"
  scalar differ      = -10 * _b[r3vs1.education]
  scalar differ_loci = -10 * (_b[r3vs1.education] + (1.96 * _se[r3vs1.education]))
  scalar differ_hici = -10 * (_b[r3vs1.education] - (1.96 * _se[r3vs1.education]))

  // Predicted probabilities Model 2
  qui margins i.education if entity == "`Y'", post
  matrix preds = r(table)
  scalar loed = preds[1,1]
     scalar mide = preds[1,2]
     scalar hied = preds[1,3]
     matrix drop preds
  
  // Ratios Model 2
  qui nlcom _b[1.education]/_b[3.education], post
  scalar ratio_     = _b[_nl_1]
  scalar ratio_loci = _b[_nl_1] - 1.96 * _se[_nl_1]
  scalar ratio_hici = _b[_nl_1] + 1.96 * _se[_nl_1]
  
  }
  else {
  scalar differ = -99
  scalar differ_loci = -99
  scalar differ_hici = -99
  scalar ratio_ = -99
   scalar ratio_loci = -99
     scalar ratio_hici = -99
  scalar loed = -99
  scalar mide = -99
  scalar hied = -99  
  }
     post `bqlrer' (`count') (`"`commandline'"') ("`desc'") ("`Y'") ///
                   (prev) (prev_loci) (prev_hici) ///
                   (differ) (differ_loci) (differ_hici) ///
                   (loed) (mide) (hied) ///
                   (ratio_) (ratio_loci) (ratio_hici) (e(N))
}

// Fit one model for entire EVS
// Model 1: OV: PSRH lowest2, adjusted -- Prevalence rate only, education not in model
qui logit poorhealth c_age c_female [pw = eu_weight] if survey == "EVS"
qui margins if survey == "EVS", post
scalar prev = 100 * _b[_cons]
scalar prev_loci = 100 * (_b[_cons] - 1.96 * _se[_cons])
scalar prev_hici = 100 * (_b[_cons] + 1.96 * _se[_cons])
di `prev' 
di `prev_loci' 
di `prev_hici'

// Model 2: OV: PSRH lowest2, adjusted
local count = `count' + 1
local desc "DV: Poor health (bottom 2), IV: Education, Age, Sex"
qui logit poorhealth i.education c_age c_female /*[pw = eu_weight]*/ if survey == "EVS"
local commandline = e(cmdline)

// Differences Model 2
qui margins r.education if survey == "EVS"
scalar differ      = -10 * _b[r3vs1.education]
scalar differ_loci = -10 * (_b[r3vs1.education] + (1.96 * _se[r3vs1.education]))
scalar differ_hici = -10 * (_b[r3vs1.education] - (1.96 * _se[r3vs1.education]))

  // Predicted probabilities Model 2
qui margins i.education if survey == "EVS", post
matrix preds = r(table)
scalar loed = preds[1,1]
scalar mide = preds[1,2]
scalar hied = preds[1,3]
matrix drop preds
  
  // Ratios Model 2
qui nlcom _b[1.education]/_b[3.education], post
scalar ratio_     = _b[_nl_1]
scalar ratio_loci = _b[_nl_1] - 1.96 * _se[_nl_1]
scalar ratio_hici = _b[_nl_1] + 1.96 * _se[_nl_1]

local Y = "EVS"

post `bqlrer' (`count') (`"`commandline'"') ("`desc'") ("`Y'") ///
              (prev) (prev_loci) (prev_hici) ///
              (differ) (differ_loci) (differ_hici) ///
              (loed) (mide) (hied) ///
              (ratio_) (ratio_loci) (ratio_hici) (e(N))
  
postclose `bqlrer'

// Create data set with FIPS codes and state names
preserve
tempfile codescheme
egen pickone = tag(entity)
keep if pickone
keep entity entity_num
save `codescheme'
restore

// Merge FIPS codes to results
use results\results, clear
merge 1:1 entity using `codescheme', keepusing(entity entity_num)
drop _merge

label var count "Generic counter"
label var commandline "Command line"
label var descri "Description of model fitted"
label var entity "Entity (string)"
label var prev   "Proportion poor health"
label var prev_loci "Proportion 95 % CI (lower)"
label var prev_hici "Proportion 95 % CI (higher)"
label var loed   "Proportion poor health lower educ."
label var mide   "Proportion poor health mid educ."
label var hied   "Proportion poor health higher educ."
label var differ "Difference lower educated - higher"
label var differ_loci  "Difference 95 % CI (lower)"
label var differ_hici  "Difference 95 % CI (higher)"
label var ratio_      "Ratio lower educated / higher"
label var ratio_loci  "Ratio 95 % CI (lower)"
label var ratio_hici  "Ratio 95 % CI (higher)"
label var n           "Size of entity sample"

// Save results, also as csv-file
save results\results, replace
outsheet using "results\results.csv", comma replace nolabel


use results\results, replace

// Create US identifier
gen usa = 0 
replace usa = 1 if regexm(entity, "US-")

replace entity = "{bf:Europe}" if entity == "EVS"

egen order_ratio = rank(-ratio_), unique 
labmask order_ratio, value(entity)

egen order_differ = rank(-differ), unique 
labmask order_differ, value(entity)

egen order_prev = rank(-prev), unique
labmask order_prev, value(entity)

twoway (dot prev order_prev if usa == 0, horizontal ndots(20)) ///
       (dot prev order_prev if usa == 1, horizontal ndots(20)) ///
    (rspike prev_loci prev_hici order_prev, horizontal) ///
    , ///
    legend(label(1 "Europe") label(2 "US") label(3 "95% CI") pos(1) ring(0)) ///
    ylabel(1/97, valuelabels ang(h) labsize(*.65)) ///
    ytitle("") ///
    xtitle("Prevalence poor health") ///
    name(prevalence, replace) ///
    xsize(2.7) ysize(10)

generate ratio_loci_p = ratio_loci               // Shorten CI's to make graph fit better
replace  ratio_loci_p =  0 if ratio_loci_p <  0  
generate ratio_hici_p = ratio_hici
replace  ratio_hici_p = 6 if ratio_hici_p > 6
 
twoway (dot ratio_ order_ratio if usa == 0, horizontal ndots(20) msymbol(smplus)) ///
       (dot ratio_ order_ratio if usa == 1, horizontal ndots(20) msymbol(smx)) ///
    (rspike ratio_loci_p ratio_hici_p order_ratio, horizontal) ///
    , ///
    legend(label(1 "Europe") label(2 "US") label(3 "95 % CI") pos(1) ring(0)) ///
    ylabel(1/97, valuelabels ang(h) labsize(*.65)) ///
    xscale(range(0 6)) xlabel(0 1 3 5)  /// 
    ytitle("") ///
    xtitle("Relative" "inequalities") ///
    xline(1) ///
    xline(1.81019, lpattern(dash)) ///
    name(relative, replace) ///
    xsize(2.7) ysize(10)

generate differ_loci_p = differ_loci              // Shorten CI's to make graph fit better
replace  differ_loci_p = -1 if differ_loci_p < -1
    
twoway (dot differ order_differ if usa == 0, horizontal ndots(20) msymbol(smplus)) ///
       (dot differ order_differ if usa == 1, horizontal ndots(20) msymbol(smx)) ///
       (rspike differ_loci_p differ_hici order_differ, horizontal) ///
    , ///
    legend(label(1 "Europe") label(2 "US") label(3 "95 % CI") pos(1) ring(0)) ///
    xscale(range(0 30)) xlabel(0(10)30)  /// 
    ylabel(1/97, valuelabels ang(h) labsize(*.65)) ///
    ytitle("") ///
    xtitle("Absolute" "inequalities") ///
    xline(6.8754, lpattern(dash)) ///
    name(absolute, replace) ///
    xsize(2.7) ysize(10)

twoway (dot prev order_prev if usa == 0, horizontal ndots(20) msymbol(smplus)) ///
       (dot prev order_prev if usa == 1, horizontal ndots(20) msymbol(smx)) ///
    (rspike prev_loci prev_hici order_prev, horizontal) ///
    , ///
    legend(label(1 "Europe") label(2 "US") label(3 "95% CI") row(1) pos(12) size(vsmall) region(lwidth(vthin) lcolor(black)) bmargin(tiny) colgap(*.3)) ///
    ylabel(1/97, valuelabels ang(h) labsize(*.65)) ///
    ytitle("") ///
    xtitle("Prevalence" "poor health") ///
    xline(11.26938, lpattern(dash)) ///
    name(prevalence_leg, replace) ///
    xsize(2.7) ysize(10)    
    
// 8.27 × 11.69
*graph combine prevalence absolute relative, xsize(8.27) ysize(10) col(3) name(comb, replace)

grc1leg prevalence_leg absolute relative, xsize(8.27) ysize(10) col(3) imargin(small) legendfrom(prevalence_leg) name(oneleg, replace) pos(6) span // ignores size-commands

graph display oneleg, xsize(6.27) ysize(9.69)  // Thus redraw so that size commands take effect

Random graphs (43): Means with confidence intervals

use "ESS3e03_5.dta", clear

// Restrict to respondents 25-42 y
keep if agea >= 25 & agea <= 42 

// Generate variables of interest

// Split ballot identifier
*fre icsbfm

// Sex
generate female = (gndr == 2) if gndr != .a
drop if female == .

// Country
replace cntry = "UK" if cntry == "GB"

// iagpnt "In your opinion, what is the ideal age for a XXX 
//         to become a mother/father?"
recode iagpnt (  0 = .a "No ideal age") ///
              (777 = .b "Refusal") ///
              (888 = .c "Don't know") ///
              (999 = .d "No answer") ///
              (998 = .e "Split ballot") ///
              , gen(idealageparent)

clonevar idealagefather = idealageparent
replace  idealagefather = .e if icsbfm == 1

clonevar idealagemother = idealageparent
replace  idealagemother = .e if icsbfm == 2

// tochld "After what age would you say a woman/man is generally too old to
//        "consider having any more children?"
recode tochld (  0 = .a "Never too old") ///
              (777 = .b "Refusal") ///
              (888 = .c "Don't know") ///
              (999 = .d "No answer") ///
              (998 = .e "Split ballot") ///
              (997 = .f "Wrong age group") ///
              , gen(toooldforchild)

clonevar  toooldforchildf = toooldforchild
replace   toooldforchildf = .e if icsbfm == 1
label var toooldforchildf "Man too old for a(nother) child"

clonevar  toooldforchildm = toooldforchild
replace   toooldforchildm = .e if icsbfm == 2
label var toooldforchildm "Woman too old for a(nother) child"

// Set outliers and "never too old" to country-specific 99th percentile
levelsof(cntry), local(country)
foreach x of varlist toooldforchildf toooldforchildm {
  foreach y of local country {
   qui sum `x' if cntry == "`y'", detail
  *di `x' _skip(2) "`y'" _skip(2) r(p95) _skip(2) r(p99) 
   replace `x' = r(p99) if `x'   == .a ///
                         & cntry == "`y'"
   replace `x' = r(p99) if `x'    > r(p99) ///
                         & !missing(`x') ///
                         & cntry == "`y'"
 }
}

// Create temporary files and postfile
tempname foo
tempname idealage
postfile `foo' str2 cntry idealagem idealagemlb idealagemub ///
                          idealagef idealageflb idealagefub ///
                          toooldforchildfm toooldforchildflb toooldforchildfub  ///
                          toooldforchildmm toooldforchildmlb toooldforchildmub ///
                          using `idealage', replace

levelsof(cntry), local(country)
foreach x of local country {
 qui reg idealagemother if cntry == "`x'"
 local idealagem   = _b[_cons]
 local idealagemlb = _b[_cons] - (1.96 * _se[_cons])
 local idealagemub = _b[_cons] + (1.96 * _se[_cons])

 qui reg idealagefather if cntry == "`x'"
 local idealagef = _b[_cons]
 local idealageflb = _b[_cons] - (1.96 * _se[_cons])
 local idealagefub = _b[_cons] + (1.96 * _se[_cons])

 qui reg toooldforchildf if cntry == "`x'"
 local toooldforchildfm   = _b[_cons]
 local toooldforchildflb = _b[_cons] - (1.96 * _se[_cons])
 local toooldforchildfub = _b[_cons] + (1.96 * _se[_cons])

 qui reg toooldforchildm if cntry == "`x'"
 local toooldforchildmm = _b[_cons]
 local toooldforchildmlb = _b[_cons] - (1.96 * _se[_cons])
 local toooldforchildmub = _b[_cons] + (1.96 * _se[_cons])

post `foo' ("`x'") (`idealagem') (`idealagemlb') (`idealagemub') ///
                   (`idealagef') (`idealageflb') (`idealagefub') ///
                   (`toooldforchildfm') (`toooldforchildflb') (`toooldforchildfub') ///
                   (`toooldforchildmm') (`toooldforchildmlb') (`toooldforchildmub') 
}
postclose `foo'

use `idealage', clear

egen order_ = rank(-toooldforchildmm), unique
labmask order_, value(cntry)

twoway (scatter toooldforchildmm order_) ///
       (rcap toooldforchildmub toooldforchildmlb order_) ///
       (scatter toooldforchildfm order_) ///
       (rcap toooldforchildfub toooldforchildflb order_) ///       
       , legend(label(1 "... women") ///
                label(3 "... men") ///
                order(3 1) pos(1) ring(0)) ///
         xlabel(1/23, val alt) ///
         ylabel(40(5)60) ///
         xtitle(" ") ytitle("Age in years") ///
         title("Age when one is too old to have a(nother) child for ...") ///
         name(tooold, replace)
   
drop order_

egen order_ = rank(-idealagem), unique
labmask order_, value(cntry)

twoway (scatter idealagem order_) ///
       (rcap idealagemub idealagemlb order_) ///
       (scatter idealagef order_) ///       
       (rcap idealagefub idealageflb order_) ///
       , legend(label(1 "... mother") ///
                label(3 "... father") ///
                order(3 1) pos(1) ring(0)) ///
         xlabel(1/23, val alt) ///
         ylabel(20(5)40) ///
         xtitle(" ") ytitle("Age in years") ///
         title("Ideal age to become a ...") ///
         name(idealage, replace)

graph combine idealage tooold, ///
          col(1) ysize(8) /// 
          note("{it:Source:} European Social Survey Round 3, own calculations" ///
               "{it:Notes:} Respondents age 25{c 150}42 y only. Error bars denote 95 % CI's", span size(small))