clear // Install rcall command *github install haghish/rcall, stable // Clear memory of R session rcall clear //#install.packages("essurvey"); rcall: library(essurvey) // Download data rcall: download_rounds(c(9), /// ess_email = "REGISTERED E-MAIL ADDRESS HERE", /// output_dir = "data", /// format = 'stata') // Open data set local vars cntry gndr polintr dweight use `vars' using "data/ESS9/ESS9e01_2.dta", clear // Variables generate interest = (inlist(polintr, 1, 2)) replace interest = . if missing(polintr) generate female = (gndr == 2) // Calculate gender gap per country tempname foo postfile `foo' str50 cntry men women gap gap_se using deleteme.dta, replace levelsof cntry, local(cntry) foreach i of local cntry { qui regress interest i.female [pw = dweight] if cntry == `"`i'"' qui margins i.female, post local men = _b[0.female] * 100 local women = _b[1.female] * 100 qui regress interest i.female [pw = dweight] if cntry == `"`i'"' qui margins r.female, post local gap = _b[r1vs0.female] * -100 local gap_se = _se[r1vs0.female] * -100 post `foo' ("`i'") (`men') (`women') (`gap') (`gap_se') } postclose `foo' // Open data set with gender gap per country use deleteme.dta, clear kountry cntry, from(iso2c) egen country = rank(gap), unique labmask country, value(NAMES_STD) generate ub = gap + 1.96 * gap_se generate lb = gap - 1.96 * gap_se generate y1 = -women generate y2 = 0 generate y3 = -women + men generate x1 = country - 0.2 generate x2 = country + 0.2 mylabels 0(20)60, myscale(-@) local(show) local spaces = 30 * " " twoway (scatter y1 country, msymbol(p) yaxis(1 2)) /// (rbar y1 y2 x1, blcolor(gs10%40) bfcolor(gs10%40) barwidth(0.4)) /// (rbar y1 y3 x2, blcolor(gs5%40) bfcolor(gs5%40) barwidth(0.4)) /// (rbar gap y2 x2, blcolor(gs3%40) bfcolor(gs3%40) barwidth(0.4)) /// (rspike ub lb x2) /// , ylabel(`show', grid ang(h)) /// ylabel(0(10)40, grid ang(h) axis(2)) /// ytitle("% interested in politics") /// ytitle("`spaces' Gender gap", axis(2)) /// xlabel(1/19, noticks ang(45) valuelabel labsize(small)) /// legend(order(2 "Men" 3 "Women" 4 "Gender gap" 5 "Gender gap 95% CI") /// col(4) pos(12) ring(0)) name(gendergap, replace)
Showing posts with label postfile. Show all posts
Showing posts with label postfile. Show all posts
May 19, 2020
Random graphs (146): Bar graph
Oct 26, 2019
Specification curve analysis
use "ZA4612_v1-0-1.dta", clear
do "ZA4612_patch_v1-0-1.do" // Some patch from data provider
keep v399-v404 v933 v827 v301 v298
// Outcomes
mvdecode v399-v404, mv(9 = .a )
generate stress = 5 - v399
generate depress = 5 - v400
generate calm = v401 - 1
generate energy = v402 - 1
generate pain = 5 - v403
generate lonely = 5 - v404
// Predictors
// Topbot
mvdecode v933 v827, mv(96 = .a \ 99 = .b)
generate topbot = v933
replace topbot = v827 if missing(topbot)
// Age
mvdecode v301, mv(999 = .a)
rename v301 age
// Female
generate female = (v298 == 2)
drop v399-v404 v933 v827 v298 // Drop unnecessary variables
// Post definition
tempname foo
postfile `foo' str50 spec k1 k2 k3 b se using deleteme.dta, replace
// Loop
forvalues k1 = 1/6 {
if `k1' == 1 local y "stress"
if `k1' == 2 local y "depress"
if `k1' == 3 local y "pain"
if `k1' == 4 local y "calm"
if `k1' == 5 local y "energy"
if `k1' == 6 local y "lonely"
forvalues k2 = 1/4 {
if `k2' == 1 local agecontrol "age"
if `k2' == 2 local agecontrol "c.age##c.age"
if `k2' == 3 local agecontrol "c.age##c.age##c.age"
if `k2' == 4 local agecontrol "c.age##c.age##c.age##c.age"
forvalues k3 = 1/3 {
if `k3' == 1 local ifs " "
if `k3' == 2 local ifs "if female == 1"
if `k3' == 3 local ifs "if female == 0"
local spec regress `y' topbot `agecontrol' female `ifs'
qui `spec'
local b = _b[topbot]
local se = _se[topbot]
post `foo' ("`spec'") (`k1') (`k2') (`k3') (`b') (`se')
}
}
}
postclose `foo'
// Plot specification curve
use deleteme, clear
// Generate ranked analytical choice variable
sort b
generate sk = _n
label var sk "Specification (sorted by coefficient size)"
// Calculate CI's
generate ub = b + 1.96 * se
generate lb = b - 1.96 * se
// Remind yourself what the variables mean
label var k1 "Outcome"
label var k2 "Age control"
label var k3 "Subsample"
// Stack indicators
generate k3c = k3
generate k2c = k2 + 1 + 3 // 3 because K3 has 3 categories, 1 for title
generate k1c = k1 + 1 + 3 + 4 + 1 // K2 has 4 categories
// Calculate some things for size of second axis
qui summarize b
global brange = r(max) - r(min)
global bmin = r(min)
global bmax = r(max)
global from_y = $bmin - (4.5 * $brange)
// Plot
twoway (scatter k1c k2c k3c sk, msymbol(o o o) msize(vsmall vsmall vsmall) ///
yscale(range(1 22)) ///
ylabel( 1 "Full sample" 2 "Females only" 3 "Males only" 4 "{bf:Subsample}" ///
5 "Age" 6 "Age squared" 7 "Age cubed" 8 "Age quartic" 9 "{bf:Age control}" ///
10 "Stress" 11 "Depression" 12 "Pain" ///
13 "Calmness" 14 "Energy" 15 "Loneliness" 16 "{bf:Outcome}", tstyle(notick) axis(1))) ///
(rcap b b sk, yaxis(2) yscale(range($from_y $bmax) axis(2)) ///
ylab(, format(%2.1g) axis(2)) ///
ytitle("{bf:Coefficient size}", axis(2) placement(north)) ///
yline(0, axis(2))) ///
(rspike ub lb sk, yaxis(2)), ///
xlabel(none) ///
legend(off) name(curve, replace)
Reference
Simonsohn, Uri, Joseph P. Simmons, and Leif D. Nelson. 2015. Specification Curve. Descriptive and Inferential Statistics on All Reasonable Specifications. University of Pennsylvania. doi: 10.2139/ssrn.2694998
Nov 17, 2017
Random graphs (119): Line plots and dot plots
// Prepare ESS round 8
use cntry essround wrkctra dweight using "ESS8e01.dta", clear
recode wrkctra (6 = .a) (7 = .b) (8 = .c) (9 = .d)
// Add ESS rounds 1-7
append using "ESS1-7e01.dta", keep(cntry essround wrkctr wrkctra dweight)
// Prepare variables
generate nocontract = (wrkctra == 3) if !missing(wrkctra)
kountry cntry, from(iso2c)
rename NAMES_STD country
// Prepare files for post commands
tempname foo1
tempname foo2
postfile `foo1' str20 country essround ll nocontract ul using `foo2', replace
qui levelsof country, local(country)
// Loop
foreach x of local country {
foreach i of numlist 1/8 {
capture logit nocontract [pw = dweight] if country == "`x'" & essround == `i'
if _rc == 0 {
qui margins [pw = dweight]
matrix prevs = r(table)
local prev = prevs[1,1] * 100
local ll = prevs[5,1] * 100
local ul = prevs[6,1] * 100
*di "`x'" _skip(5) `i' _skip(5) `ll' _skip(5) `prev' _skip(5) `ul'
post `foo1' ("`x'") (`i') (`ll') (`prev') (`ul')
}
}
}
postclose `foo1'
// Plot estimates
use `foo2', clear
// Fix value label
label define essround 1 "2002" 2 "2004" 3 "2006" 4 "2008" ///
5 "2010" 6 "2012" 7 "2014" 8 "2016", modify
label val essround essround
// First plot
sort country essround
twoway (rarea ll ul essround, lcolor(white)) ///
(connected nocontract essround), ///
by(country, ///
note("") ///
legend(off)) ///
xlabel(1/8, val ang(90)) ///
xtitle("") ytitle("% without job contract") ///
name(figure2a, replace) ysize(9)
// Second plot
reshape wide ll ul nocontract, i(country) j(essround)
scores average = mean(nocontract*)
egen order_ = rank(average), unique
labmask order_, value(country)
twoway (dot nocontract1 order_, horizontal) ///
(dot nocontract2 order_, horizontal) ///
(dot nocontract3 order_, horizontal) ///
(dot nocontract4 order_, horizontal) ///
(dot nocontract5 order_, horizontal) ///
(dot nocontract6 order_, horizontal) ///
(dot nocontract7 order_, horizontal) ///
(dot nocontract8 order_, horizontal) ///
(rspike ul1 ll1 order_, horizontal) ///
(rspike ul2 ll2 order_, horizontal) ///
(rspike ul3 ll3 order_, horizontal) ///
(rspike ul4 ll4 order_, horizontal) ///
(rspike ul5 ll5 order_, horizontal) ///
(rspike ul6 ll6 order_, horizontal) ///
(rspike ul7 ll7 order_, horizontal) ///
(rspike ul8 ll8 order_, horizontal), ///
ylabel(1/32, val) ytitle("") xscale(alt) ///
xtitle("% without job contract") ///
legend(order(1 "2002" 2 "2004" 3 "2006" 4 "2008" ///
5 "2010" 6 "2012" 7 "2014" 8 "2016") ///
pos(5) ring(0)) ///
name(figure2b, replace) ysize(9)
graph combine figure2a figure2b, ///
col(2) name(figure2, replace) ///
note(" " "{it:Source:} European Social Survey 2002-16, weighted data. {it:Note:} Error bands/spikes denote 95% confidence intervals", span size(*.7))
Aug 12, 2017
Random graphs (109): Dropping one country at a time
// Open ESS 6 data
use stflife eduyrs agea gndr cntry using "ESS6e02_2.dta", clear
// Life satisfaction
recode stflife (77 88 99 = .), gen(lifesat)
// Education
sum eduyrs if eduyrs < 77, detail
generate education = eduyrs if eduyrs < 77
replace education = r(p99) if education >= r(p99) & !missing(education)
// Age
recode agea (999 = .), gen(age)
// Gender
generate female = (gndr == 2) if gndr != 9
// Fit model across all countries and save estimates
regress lifesat education c.age##c.age i.female, cluster(cntry)
local opointest = _b[education]
local olb = _b[education] - 1.96 * _se[education]
local oub = _b[education] + 1.96 * _se[education]
// Define temporary objects
tempname foo
tempname foox
postfile `foo' str3 cntry pointest lb ub using `foox', replace
// Drop one country at a time
levelsof cntry, local(country)
foreach i of local country {
qui regress lifesat education c.age##c.age i.female if cntry != "`i'"
local pointest = _b[education]
local lb = _b[education] - 1.96 * _se[education]
local ub = _b[education] + 1.96 * _se[education]
post `foo' ("`i'") (`pointest') (`lb') (`ub')
}
postclose `foo'
// Open estimates
use `foox', clear
// Creats country variable
kountry cntry, from(iso2c)
rename NAMES_STD geo
encode geo, generate(country)
label define country 29 "Kosovo", modify
// Plot estimates
sort country
twoway (rarea lb ub country, horizontal color(gs14)) ///
(function y = `opointest', horizontal range(country) lpattern(solid)) ///
(function y = `oub', horizontal range(country) lpattern(dash)) ///
(function y = `olb', horizontal range(country) lpattern(dash)) ///
(dot pointest country, horizontal) ///
, ylabel(1/29, val) ytitle("") xscale(alt) ///
legend(order(5 "Point estimate after excluding country" ///
2 "Point estimate from complete sample" ///
1 "95% CI after excluding country" ///
3 "95% CI (cluster robust) from complete sample") pos(6) span) ///
ysize(7) name(robustness, replace)
Apr 4, 2016
Random graphs (71): Proportions with confidence intervals
use sharew5_rel1-0-0_it.dta, clear
// Prepare variables
recode it004_ (1 = 1 "Yes") (5 = 0 "No") (-1 -2 . = .), generate(internetuser)
label var internetuser "Internet use in last 7 days"
decode country, gen(cntry)
// Define stuff
tempname foo
postfile `foo' str25 cntry perc perc_ll perc_ul using "C:\Windows\Temp\test.dta", replace
levelsof cntry, local(levels)
// Calculate proportions by country
foreach country of local levels {
capture proportion internetuser if cntry == "`country'"
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) `fperc_ll' _skip(2) `fperc' _skip(2) `fperc_ul'
if "`fperc'" != "" { // Make sure that the loop doesn't break if empty
post `foo' ("`country'") (`fperc') (`fperc_ll') (`fperc_ul')
}
}
postclose `foo'
// Use posted data set
use "C:\Windows\Temp\test.dta", clear
// Sort countries by size
egen order_ = rank(perc), unique
labmask order_, value(cntry)
// Plot
twoway (dot perc order_, horizontal) ///
(rspike perc_ll perc_ul order_, horizontal), ///
ylabel(1/15, val) legend(off) xscale(range(20 80) alt) ///
ytitle("") xtitle("% Internet users (in last 7 days)") ///
xlabel(20(10)80, grid) ///
note(" " "{it:Source:} SHARE wave 5, doi:10.6103/SHARE.w5.100", span)
Mar 14, 2016
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
Feb 25, 2016
Random graphs (61): Dot plot for two groups
tempname foo
tempname foox
postfile `foo' str2 cntry perc_f perc_f_ll perc_f_ul ///
perc_m perc_m_ll perc_m_ul using `foox' , replace
levelsof cntry, local(levels1)
// Calculate proportions by country
foreach country of local levels1 {
capture proportion nonesep_f if cntry == "`country'"
*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
capture proportion nonesep_m if cntry == "`country'"
*matrix list r(table)
matrix mcoefs = r(table)
local mperc = mcoefs[1,2] * 100
local mperc_ll = mcoefs[5,2] * 100
local mperc_ul = mcoefs[6,2] * 100
di "`country'" _skip(2) `fperc_ll' _skip(2) `fperc' _skip(2) `fperc_ul'
di "`country'" _skip(2) `mperc_ll' _skip(2) `mperc' _skip(2) `mperc_ul'
if "`fperc'" != "" { // Make sure that the loop doesn't break if empty
post `foo' ("`country'") (`fperc') (`fperc_ll') (`fperc_ul') (`mperc') (`mperc_ll') (`mperc_ul')
}
}
postclose `foo'
// Use posted data set
use `foox', clear
// Reverse sign of one variable
replace perc_m = -perc_m
replace perc_m_ll = -perc_m_ll
replace perc_m_ul = -perc_m_ul
// Sort countries by size
egen order_ = rank(-perc_f), unique
labmask order_, value(cntry)
// Plot
twoway (dot perc_m order_, horizontal) ///
(dot perc_f order_, horizontal) ///
, ylabel(1/31, val) legend(order(1 "Men" 2 "Women") pos(6) row(1)) ///
ytitle("") xtitle("% of men/women pooling all income") ///
xline(0) ///
xlabel(-100 "100" -50 "50" 0 "0" 50 "50" 100 "100") ///
name(somefigure, replace) ysize(8)
Jan 30, 2016
Random graphs (56): Intraclass Correlation Coefficients with confidence intervals
webuse productivity, clear
tempname foo
postfile `foo' str100 commandline_str str100 description_str icc icclb iccub N N_groups using "C:\Windows\Temp\test.dta", replace
qui mixed gsp || region: if year <= 1973
qui estat icc
matrix groups = e(N_g)
scalar n_g = groups[1,1]
matrix cis = r(ci2)
scalar icclb = cis[1,1]
scalar iccub = cis[1,2]
local desc "ICC 1970 to 1973"
post `foo' (e(cmdline)) ("`desc'") (r(icc2)) (icclb) (iccub) (e(N)) (n_g)
qui mixed gsp || region: if year >= 1974 | year <= 1977
qui estat icc
matrix groups = e(N_g)
scalar n_g = groups[1,1]
matrix cis = r(ci2)
scalar icclb = cis[1,1]
scalar iccub = cis[1,2]
local desc "ICC 1974 to 1977"
post `foo' (e(cmdline)) ("`desc'") (r(icc2)) (icclb) (iccub) (e(N)) (n_g)
qui mixed gsp || region: if year >= 1978 | year <= 1981
qui estat icc
matrix groups = e(N_g)
scalar n_g = groups[1,1]
matrix cis = r(ci2)
scalar icclb = cis[1,1]
scalar iccub = cis[1,2]
local desc "ICC 1978 to 1981"
post `foo' (e(cmdline)) ("`desc'") (r(icc2)) (icclb) (iccub) (e(N)) (n_g)
qui mixed gsp || region: if year >= 1981 | year <= 1984
qui estat icc
matrix groups = e(N_g)
scalar n_g = groups[1,1]
matrix cis = r(ci2)
scalar icclb = cis[1,1]
scalar iccub = cis[1,2]
local desc "ICC 1981 to 1984"
post `foo' (e(cmdline)) ("`desc'") (r(icc2)) (icclb) (iccub) (e(N)) (n_g)
qui mixed gsp || region: if year >= 1985 | year <= 1986
qui estat icc
matrix groups = e(N_g)
scalar n_g = groups[1,1]
matrix cis = r(ci2)
scalar icclb = cis[1,1]
scalar iccub = cis[1,2]
local desc "ICC 1985 to 1986"
post `foo' (e(cmdline)) ("`desc'") (r(icc2)) (icclb) (iccub) (e(N)) (n_g)
postclose `foo'
use "C:\Windows\Temp\test.dta", clear
list
encode description_str, gen(description)
twoway (dot icc description, horizontal) ///
(rcap icclb iccub description, horizontal) ///
, ylabel(1/5, val) ytitle("") legend(off) ///
xlabel(0 (.1) 1) xtitle("ICC") name(figure, replace)
erase "C:\Windows\Temp\test.dta"
Aug 18, 2015
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
Jun 20, 2015
Random graphs (49): Means with confidence intervals
tempname foo
tempname cogscore
postfile `foo' time ART cog coglb cogub using `cogscore', replace
levelsof(time), local(timepoint)
levelsof(ART), local(art)
foreach y of local art {
foreach x of local timepoint {
qui reg COG if time == `x' & ART == `y'
local cog = _b[_cons]
local coglb = _b[_cons] - (1.96 * _se[_cons])
local cogub = _b[_cons] + (1.96 * _se[_cons])
di `timepoint' _skip(5) `art' _skip(5) `cog' _skip(5) `coglb' _skip(5) `cogub'
post `foo' (`x') (`y') (`cog') (`coglb') (`cogub')
}
}
postclose `foo'
use `cogscore', clear
label define time 1 "Age 3" 2 "Age 5" 3 "Age 7" 4 "Age 11"
label val time time
twoway (rarea coglb cogub time if ART == 1) ///
(rarea coglb cogub time if ART == 0) ///
(line cog time if ART == 0) ///
(line cog time if ART == 1) ///
, legend(label(3 "Planned natural conception") ///
label(4 "ART conception") col(1) ring(0) pos(11) order(4 3)) ///
ytitle("Average cognitive development") ///
xtitle(" ") xlabel(1(1)4, valuelabels) ///
note("{it:Note:} Shaded areas denote 95 % confidence intervals." ///
"{it:Source:} Millennium Cohort Study.", span) name(bands, replace)
Apr 20, 2015
Random graphs (45): Re-creating -graph dot- using -twoway dot-
This plot recreates this one, but adds confidence intervals.
preserve
// Calculate country-specific means to be plotted
tempname foo
tempname idealage
postfile `foo' quintile incocomp incocomp_lb incocomp_ub using `foo2', replace
levelsof(quintile), local(quintile)
foreach x of local quintile {
qui reg incocomp if quintile == `x', cluster(cntry)
local incocomp = _b[_cons]
local incocomplb = _b[_cons] - (1.96 * _se[_cons])
local incocompub = _b[_cons] + (1.96 * _se[_cons])
post `foo' (`x') (`incocomp') (`incocomplb') (`incocompub')
}
postclose `foo'
// Plot country-specific means
use `foo2', clear
twoway (dot incocomp quintile) /// -twoway scatter- doesn't have the lines
(rcap incocomp_lb incocomp_ub quintile) ///
, yscale(range(1.5 3.0)) ///
ylabel(1.5 (.25) 3.0, format(%6.2f)) ///
legend(off)
title("Income quintiles")
xlabel(1 `""1" "(poorest)""' /// Double-compound quotes allow for the line break
2 3 4 5 `""5" "(richest)"') ///
fxsize(40) /// Forces x-size to be 40 percent of its original size
xtitle("") ///
xscale(range(0.2 5.8)) /// Making this a bit bigger makes room for the labels
name(byquintile_hor, replace)
restore
// Calculate country-specific means to be plotted 2
preserve
tempname foo
tempname idealage
postfile `foo' str2 cntry incocomp incocomp_lb incocomp_ub using `foo2', replace
levelsof(cntry), local(cntry)
foreach x of local cntry {
qui reg incocomp if cntry == "`x'"
local incocomp = _b[_cons]
local incocomplb = _b[_cons] - (1.96 * _se[_cons])
local incocompub = _b[_cons] + (1.96 * _se[_cons])
post `foo' ("`x'") (`incocomp') (`incocomplb') (`incocompub')
}
postclose `foo'
use `foo2', clear
// Create sorted-by-size variable
egen order_ = rank(-incocomp), unique
labmask order_, value(cntry)
// Plot country-specific means 2
twoway (dot incocomp order_) ///
(rcap incocomp_lb incocomp_ub order_) ///
, yscale(off) /// Not necessary, as we're using the one from the other graph
xlabel(1/23, val alt) ///
title("Countries") ///
xtitle("") ///
legend(off) ///
/* graphregion(margin(b=16)) */ /// Could help as -ycommon- isn't accounting for the different label sizes
name(bycntry_hor, replace)
// Combine plots
graph combine byquintile_hor bycntry_hor, ///
l1title("Mean income comparison orientation") ///
ycommon row(1) ///
imargin(zero) // Sets margin between both plots to 0
restore
Jan 9, 2015
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))
Nov 30, 2014
Random graphs (37): Proportions with confidence intervals
use acountry asex aage ayear a602 a612 aweight using GGS_Wave1_V.4.2.dta, clear // Create outcome variable recode a612 (1 2 1601 = 1 "(Likely) infertile") /// (3 4 1501 1602 = 0 "(Likely) fertile") /// (2001 . .a .b .c = .x "Don't know"), /// gen(infertility) replace infertility = 0 if a602 == 1 // Also fertile when currently pregnant // Create country identifier gen cntry = "" replace cntry = "BG" if acountry == 11 replace cntry = "RU" if acountry == 12 replace cntry = "GE" if acountry == 13 replace cntry = "DE" if acountry == 14 replace cntry = "FR" if acountry == 15 replace cntry = "HU" if acountry == 16 replace cntry = "IT" if acountry == 17 replace cntry = "NL" if acountry == 18 replace cntry = "RO" if acountry == 19 replace cntry = "NO" if acountry == 20 replace cntry = "AT" if acountry == 21 replace cntry = "EE" if acountry == 22 replace cntry = "BE" if acountry == 23 replace cntry = "AU" if acountry == 24 replace cntry = "LT" if acountry == 25 replace cntry = "PL" if acountry == 26 replace cntry = "CZ" if acountry == 28 // Recode gender variable gen sex = "" replace sex = "Men" if asex == 1 replace sex = "Women" if asex == 2 // Drop unnecessary countries drop if inlist(cntry, "IT", "AU", "NL") // Prevalence plot by country preserve // Drop unnecesary age groups keep if aage <= 39 keep if aage >= 34 svyset [pweight = aweight] // Declare weight // Create temporary files and postfile tempname foo tempname infert postfile `foo' str2 cntry str5 sex perc_infert perc_infert_ll perc_infert_ul using `infert' , replace levelsof cntry, local(levels1) levelsof sex, local(sex1) // Calculate proportions by country and sex foreach sex of local sex1 { foreach country of local levels1 { // Calculate proportions capture proportion infertility if cntry == "`country'" & sex == "`sex'" *matrix list r(table) matrix coefs = r(table) local perc = coefs[1,2] * 100 local perc_ll = coefs[5,2] * 100 local perc_ul = coefs[6,2] * 100 di "`country'" _skip(2) `perc_ll' _skip(2) `perc' _skip(2) `perc_ul' if "`perc'" != "" { // Make sure that the loop doesn't break if empty post `foo' ("`country'") ("`sex'") (`perc') (`perc_ll') (`perc_ul') } } } postclose `foo' // Use posted data set use `infert', clear // Sort countries by size for women egen order_ = rank(-perc_infert) if sex == "Women", unique labmask order_, value(cntry) // Plot for women twoway (dot perc_infert order_ if sex == "Women", horizontal) /// (rcap perc_infert_ll perc_infert_ul order_ if sex == "Women", horizontal) /// , ylabel(1/14, val) legend(off) /// ytitle("") xtitle("% of {bf:women} aged 35{c 0150}39 reporting (suspected) infertility") /// /*caption("{it:Note:} Error bars denote 95 % confidence intervals", size(small) span)*/ /// name(infert3539women, replace) // Drop women's order drop order_ // Sort countries by size for men egen order_ = rank(-perc_infert) if sex == "Men", unique labmask order_, value(cntry) // Plot for men twoway (dot perc_infert order_ if sex == "Men", horizontal) /// (rcap perc_infert_ll perc_infert_ul order_ if sex == "Men", horizontal) /// , ylabel(1/13, val) legend(off) /// ytitle("") xtitle("% of {bf:men} aged 35{c 0150}39 reporting (suspected) infertility") /// /*caption("{it:Note:} Error bars denote 95 % confidence intervals", size(small) span)*/ /// name(infert3539men, replace) // Combine graphs graph combine infert3539women infert3539men /// , ysize(8) col(1) xcommon /// caption("{it:Source:} Generations and Gender Survey Wave 1 (2005{c 0150}2011}" /// "{it:Note:} Error bars denote 95 % confidence intervals", size(vsmall) span) restore// Plot infertility by age // Recode age into groups recode aage (20/24 = 1 "20–24") /// (25/29 = 2 "25–29") /// (30/34 = 3 "30–34") /// (35/39 = 4 "35–39") /// (40/45 = 5 "40–45") /// (1/19 = .a "Below 20 y") /// (46/85 = .b "46+ y") /// , gen(age) // Save valuelabels in temporary do-file tempfile valuelabels label save age using `valuelabels' svyset [pweight = aweight] // Declare weight // Create temporary files and postfile tempname foo tempname infertage postfile `foo' str2 cntry str5 sex age perc_infert perc_infert_ll perc_infert_ul using `infertage' , replace levelsof cntry, local(levels1) levelsof sex, local(sex1) // Calculate proportion by age group, country, and sex foreach sex of local sex1 { foreach country of local levels1 { forvalues x = 1/5 { // Calculate proportions capture proportion infertility if age == `x' & cntry == "`country'" & sex == "`sex'" *matrix list r(table) matrix coefs = r(table) local perc = coefs[1,2] * 100 local perc_ll = coefs[5,2] * 100 local perc_ul = coefs[6,2] * 100 di "`country'" _skip(2) `x' _skip(2) `perc_ll' _skip(2) `perc' _skip(2) `perc_ul' if "`perc'" != "" { // Make sure that the loop doesn't break if post `foo' ("`country'") ("`sex'") (`x') (`perc') (`perc_ll') (`perc_ul') } } } } postclose `foo' // Use posted values use `infertage', clear drop if perc_infert == . // Drop empty rows // Define and set saved value labels again do `valuelabels' label val age age // Plot for women by country and age twoway (rarea perc_infert_ll perc_infert_ul age if sex == "Women", fcolor(gs15) cmissing(no)) /// (line perc_infert age if sex == "Women", cmissing(no)) /// , by(cntry, cols(3) legend(off) note("") /// caption("{it:Source:} Generations and Gender Survey Wave 1 (2005{c 0150}2011)" /// "{it:Note:} Gray areas denote 95 % confidence intervals", size(vsmall) span)) /// ytitle("% of {bf:women} reporting (suspected) infertility") xtitle("Age group") /// xlabel(1/5, val alternate) /// ysize(8) name(infertage, replace)
Sep 14, 2014
Random graphs (29): Line plots with confidence bands
use mainstat age country refyear coeff if refyear == 2011 using eulfs.dta, replace
drop refyear // Drop unnecessary variable
generate retired = (mainstat == 4) if mainstat != . // Create dummy variable for retired
svyset [pweight = coeff] // Declare weight
// Create temporary files and postfile
tempname foo
tempname retage
postfile `foo' cntry age perc_retired perc_retired_ll perc_retired_ul using `retage' , replace
levelsof country, local(levels1)
foreach country of local levels1 {
forvalues x = 57 (5) 77 {
// Calculate proportions
capture proportion retired if age == `x' & country == `country'
*matrix list r(table)
matrix coefs = r(table)
local perc = coefs[1,2] * 100
local perc_ll = coefs[5,2] * 100
local perc_ul = coefs[6,2] * 100
qui di `country' _skip(2) `x' _skip(2) `perc_ll' _skip(2) `perc_' _skip(2) `perc_ul'
if "`perc'" != "" { // Make sure that the loop doesn't break if
post `foo' (`country') (`x') (`perc') (`perc_ll') (`perc_ul')
}
}
}
postclose `foo'
use `retage', clear
drop if perc_retired == . // Drop empty rows
// Label variables
label define age 57 "55{c 150}59 y" 62 "60{c 150}64 y" 67 "65{c 150}69 y" ///
72 "70{c 150}74 y" 77 "75{c 150}79 y", modify
label value age age
label define cntry 1 "AT" 2 "BE" 3 "BG" 4 "CH" 5 "CY" ///
6 "CZ" 7 "DE" 8 "DK" 9 "EE" 10 "ES" ///
11 "FI" 12 "FR" 13 "GR" 14 "HU" 15 "IE" ///
16 "IS" 17 "IT" 18 "LT" 19 "LU" 20 "LV" ///
21 "MT" 22 "NL" 23 "NO" 24 "PL" 25 "PT" ///
26 "RO" 27 "SE" 28 "SI" 29 "SK" 30 "UK" ///
, modify
label value cntry cntry
// Plot plots
twoway (rarea perc_retired_ll perc_retired_ul age, fcolor(gs14)) ///
(line perc_retired age) ///
, by(cntry, cols(3) legend(off) note("") ///
caption("{it:Source:} EU-LFS 2011" ///
"{it:Notes:} Gray areas denote 95 % confidence intervals"))) ///
ytitle("% retired") xtitle("Age group") ///
xlabel(57 (5) 77, val ang(v)) ///
ysize(10)
/* ,by(, note("")) suppresses the default "Graphs by" note) */
/* legend(off) needs to be in by(, legend(off)) to work */
Subscribe to:
Posts (Atom)