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 twoway scatter. Show all posts
Showing posts with label twoway scatter. 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
Sep 18, 2019
Random graphs (142): Boxplot using twoway commands with overlaid data
use v358 v392 v585 v829 if inrange(v829, 1, 5) using "data allbus\ZA4612_v1-0-1.dta", clear
// Subjective social mobility
recode v829 (5 = 1 "Much lower status") ///
(4 = 2 "Lower status") ///
(3 = 3 "About equal") ///
(2 = 4 "Higher status") ///
(1 = 5 "Much higher status") ///
(95 96 98 99 = .), gen(ssm)
label var ssm "Status of own job compared to father's"
// Respondent's social status ISEI
clonevar isei = v358
replace isei = . if inlist(isei, 0, 99)
replace isei = v392 if missing(isei)
replace isei = . if inlist(isei, 0, 99)
label var isei "Respondent's ISEI"
// Father's social status ISEI
clonevar fisei = v585
replace fisei = . if inlist(fisei, 0, 99)
label var fisei "Father's ISEI"
// Objective status mobility
gen osm = isei - fisei
label var osm "Objective status mobility"
// Boxplot using -graph- command
graph box osm, over(ssm, descending) horizontal ///
yline(0) l1title("Status of own job compared to father's" "(self-assessment)") ///
name(figure1, replace)
// Boxplot using -twoway- command
// Calculate stuff for boxplot
qui bys ssm: sum osm, detail
qui bys ssm: egen min = pctile(osm), p(25) // Interquartile range
qui bys ssm: egen max = pctile(osm), p(75) // Interquartile range
qui bys ssm: egen med = median(osm)
// Plot
twoway (scatter ssm osm, msymb(o) jitter(2) mcolor(red%10)) ///
(scatter ssm med, msymb(oh) mcolor(black)) ///
(rspike min max ssm, horizontal lcolor(black)), ///
ylabel(1/5, valuelabel) xline(0) ///
ytitle("Status of own job compared to father's" "(self-assessment)") ///
legend(order(2 "Median" 3 "Interquartile range") ring(0) pos(11)) ///
xtitle(Objective social mobility) name(figure2, replace)
drop min max med // Drop stuff for boxplot
Jul 24, 2018
Random graphs (137): Logistic regression
clear
// Generate data
set seed 1
set obs 50
gen hours = rnormal(3, 1) // Number of hours studied
gen e = rnormal(1,1)
gen questions = 2 + 2*hours + 1*e // Questions answered correctly
qui sum questions, detail
generate pass = (questions >= r(p75)) // Passing the exam
// 1) Histogram of outcome variable
twoway (histogram pass, discrete percent), ///
xlabel(0 "[0] Failed" 1 "[1] Passed") xtitle("") ///
ytitle("Percent of students") xsize(4) ysize(4) name(figure5, replace)
// 2) Scatterplot
twoway (scatter pass hours), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
ytitle("Exam success") ///
ylabel(0 "[0] Failed" 1 "[1] Passed") legend(off) ///
xsize(4) ysize(4) name(figure6, replace)
// 3) Scatterplot with regression line
regress pass hours
local intercept = round(_b[_cons], .01)
local x = round(_b[hours], .01)
twoway (scatter pass hours) ///
(lfit pass hours, lpattern(solid) range(1 5)), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
text(.8 2 "y = `intercept' + `x' x + e", size(large)) ///
ytitle("Exam success") ///
ylabel(0 "[0] Failed" 1 "[1] Passed") legend(off) ///
xsize(4) ysize(4) name(figure7, replace)
// 4) Logit curve
logit pass hours
predict yhat
twoway (scatter pass hours) ///
(line yhat hours, lpattern(solid) sort), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
ytitle("Exam success") ///
ylabel(0 "[0] Failed" 1 "[1] Passed") legend(off) ///
xsize(4) ysize(4) name(figure8, replace)
graph combine figure5 figure6 figure7 figure8, ///
col(2) xsize(8) ysize(8) altshrink name(figures58, replace)
Random graphs (136): Linear probability model
clear
// Generate data
set seed 1
set obs 50
gen hours = rnormal(3, 1) // Number of hours studied
gen e = rnormal(1,1)
gen questions = 2 + 2*hours + 1*e // Questions answered correctly
qui sum questions, detail
generate pass = (questions >= r(p75)) // Passing the exam
// 1) Histogram of outcome variable
twoway (histogram pass, discrete percent), ///
xlabel(0 "[0] Failed" 1 "[1] Passed") xtitle("") ///
ytitle("Percent of students") xsize(4) ysize(4) name(figure5, replace)
// 2) Scatterplot
twoway (scatter pass hours), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
ytitle("Exam success") ///
ylabel(0 "[0] Failed" 1 "[1] Passed") legend(off) ///
xsize(4) ysize(4) name(figure6, replace)
// 3) Scatterplot with regression line
regress pass hours
local intercept = round(_b[_cons], .01)
local x = round(_b[hours], .01)
twoway (scatter pass hours) ///
(lfit pass hours, lpattern(solid) range(1 5)), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
text(.8 2 "y = `intercept' + `x' x + e", size(large)) ///
ytitle("Exam success") ///
ylabel(0 "[0] Failed" 1 "[1] Passed") legend(off) ///
xsize(4) ysize(4) name(figure7, replace)
// 4) Heteroskedastic residuals
regress pass hours
predict resid, resid
twoway (scatter resid hours), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
yline(0) name(figure7a, replace)
graph combine figure5 figure6 figure7 figure7a, ///
col(2) xsize(8) ysize(8) altshrink name(figures57, replace)
Random graphs (135): Scatterplot with OLS regression
clear
// Generate data
set seed 1
set obs 50
gen hours = rnormal(3, 1) // Number of hours studied
gen e = rnormal(1,1)
gen questions = 2 + 2*hours + 1*e // Questions answered correctly
// 1) Basic scatterplot
twoway (scatter questions hours), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
ytitle("Number of questions answered correctly") ///
name(figure1, replace) ylabel(0 (5) 15) xsize(4) ysize(4)
// 2) Scatterplot with regression line
twoway (scatter questions hours) ///
(lfit questions hours), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
ytitle("Number of questions answered correctly") ///
ylabel(0 (5) 15) legend(off) name(figure2, replace) xsize(4) ysize(4)
// 3) Scatterplot with regression line and equation
regress questions hours
local intercept = round(_b[_cons], .1)
local x = round(_b[hours], .1)
twoway (scatter questions hours) ///
(lfit questions hours), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
text(14 2 "y = `intercept' + `x' + e", size(large)) ///
ytitle("Number of questions answered correctly") ///
ylabel(0 (5) 15) legend(off) name(figure3, replace) xsize(4) ysize(4)
// 4) Scatterplot with regression line and equation and labels for components
twoway (scatter questions hours) ///
(lfit questions hours, lpattern(solid) range(0 5)) ///
(function y = 1.6 + 2.5, range(1 2)) ///
(function y = 2, range(4.1 6.6) horizontal lpattern(solid) lcolor(red)) ///
(scatteri 1.6 0 (3) "Intercept", msymbol(o) mlabcolor(red)) ///
(scatteri 5 2 "Slope", msymbol(i) mlabcolor(red)), ///
xlabel(0 (1) 5) xtitle("Hours studied for exam") ///
text(14 2 "y = `intercept' + `x' x + e", size(large)) ///
ytitle("Number of questions answered correctly") ///
ylabel(0 1.6 5 10 15) legend(off) name(figure4, replace) xsize(4) ysize(4)
graph combine figure1 figure2 figure3 figure4, col(2) xsize(8) ysize(8) altshrink name(figures14, replace)
Apr 19, 2018
Random graphs (131): Scatterplot
sdmxuse data ESTAT, clear dataset(crim_pris_cap) dimensions(A.PRIS_ACT_CAP.P_HTHAB.)
drop unit indic_cr freq
rename value prison
tempfile x
replace geo = "UK" if geo == "UKC-L"
save `x', replace
sdmxuse data ESTAT, clear dataset(spr_exp_sum) start(2008) end(2015) dimensions(A.TOTALNOREROUTE.PC_GDP..)
drop spdeps unit freq
rename value spending
merge 1:1 geo time using `x', keep(match) nogenerate
replace geo = "GR" if geo == "EL"
kountry geo, from(iso2c)
rename NAMES_STD country
regress prison spending if time == "2015"
local r2 = round(`e(r2)', .01) * 100
twoway (scatter prison spending if time == "2015", mlabel(country)) ///
(lfit prison spending if time == "2015"), ///
ytitle("Prison population per 100,000 inhabitants") ///
xtitle("Total social spending as % of GDP") ///
legend(order(2 "Linear fit, explained variance = `r2'%") pos(1) ring(0)) ///
name(f2015, replace)
Nov 13, 2017
Random graphs (117): Scatterplots
cls
clear
input str25 country cf grades attainment gini
"Germany" .46 . .51 31.4
"Norway" .45 .48 .41 26.8
"Sweden" .48 .52 .44 27.2
"Finland" . . .36 26.8
"United Kingdom" . . .42 34.1
"USA (Add Health)" .57 .42 . 41.0
"USA (PSID)" .55 . .51 41.0
end
twoway scatter cf gini, mlabel(country) ytitle("Cognitive skills" "sibling correlation") ///
xtitle("Income inequality" "(Gini coefficient, World Bank)") ///
xlabel(25 (5) 45) ylabel(, format(%6.2f)) name(A, replace)
twoway scatter grades gini, mlabel(country) ytitle("School grades" "sibling correlation") ///
xtitle("Income inequality" "(Gini coefficient, World Bank)") ///
xlabel(25 (5) 45) ylabel(, format(%6.2f)) name(B, replace)
twoway scatter attainment gini, mlabel(country) ytitle("Final educational attainment" "sibling correlation") ///
xtitle("Income inequality" "(Gini coefficient, World Bank)") ///
xlabel(25 (5) 45) ylabel(, format(%6.2f)) name(C, replace)
graph combine A B C, col(1) ysize(6) xsize(4) xcommon ycommon
Aug 1, 2017
Measurement error in the outcome vs. in the predictor
clear
set seed 1
set obs 5000
// Create variables
gen e = rnormal()
gen x_true = rnormal()
gen y_true = 1 + 2 * x_true + 3 * e
gen y_observed = y_true + e
gen x_observed = x_true + e
label var y_true "Y"
label var x_true "X"
label var y_observed "Y w/ error"
label var x_observed "X w/ error"
// Loop over variables
foreach y of varlist y_observed y_true {
foreach x of varlist x_observed x_true {
regress `y' `x'
local b = round( _b[`x'], .01)
local se = round(_se[`x'], .01)
twoway (scatter `y' `x', msymbol(p)) ///
(lfit `y' `x') ///
, name(`y'`x', replace) ///
ytitle(`: variable label `y'') ///
xtitle(`: variable label `x'') ///
title("{it:B} = `b', {it:SE} = `se'") ///
legend(off) nodraw
}
}
// Combine plots
graph combine y_truex_true ///
y_truex_observed ///
y_observedx_true ///
y_observedx_observed, ///
xcommon ycommon col(2) row(2) ///
title("Measurement error in the outcome vs. in the predictor")
Jul 27, 2017
Random graphs (105): Scatterplot
// Open ISSP data
use C_ALPHAN SEX V51 WEIGHT using ZA4350_v2-0-0.dta, clear
// Fix country variable
replace C_ALPHAN = "BE" if C_ALPHAN == "BE-FLA"
replace C_ALPHAN = "DE" if C_ALPHAN == "DE-E"
replace C_ALPHAN = "DE" if C_ALPHAN == "DE-W"
replace C_ALPHAN = "UK" if C_ALPHAN == "GB-GBN"
// Recode job satisfaction
recode V51 (7 = 0 "Completely dissatisfied") ///
(6 = 1 "Very dissatisfied") ///
(5 = 2 "Fairly dissatisfied") ///
(4 = 3 "Neither satisfied nor dissatisfied") ///
(3 = 4 "Fairly satisfied") ///
(2 = 5 "Very satisfied") ///
(1 = 6 "Completely satisfied"), gen(jobsat)
label var jobsat "Job satisfaction"
// Recode sex
gen female = (SEX == 2) if !missing(SEX)
label define female 0 "Male" 1 "Female"
label val female female
label var female "Sex"
drop if missing(female)
// Get rid of old variables
drop SEX V51
// Collapse and reshape data
collapse jobsat [pw = WEIGHT], by(C_ALPHAN female)
reshape wide jobsat, i(C_ALPHAN) j(female)
label var jobsat0 "Job satisfaction: men"
label var jobsat1 "Job satisfaction: women"
// Plot
twoway (scatter jobsat0 jobsat1 if jobsat0 > jobsat1, mlabel(C_ALPHAN) mlabpos(9)) ///
(scatter jobsat0 jobsat1 if jobsat0 < jobsat1, mlabel(C_ALPHAN) mlabpos(6)) ///
(function y = x, range(3.5 5.1)) ///
, xlabel(3.5 (.5) 5.0, format(%6.2f)) ///
ylabel(3.5 (.5) 5.0, format(%6.2f)) ///
ytitle("Job satisfaction: men") ///
xtitle("Job satisfaction: women") ///
title("Cross-national variation in the gender gap in job satisfaction", span) ///
note(" " "{it:Source:} International Social Survey Program (ISSP) 2005, doi:10.4232/1.11648", span) ///
legend(off) name(figure1, replace)
Jul 10, 2017
Using -sdmxuse- to download Eurostat data
// 1) Download data:
sdmxuse data ESTAT, dataset(ilc_di12) clear attributes
destring time, gen(year)
// Generate country name variable
replace geo = "GR" if geo == "EL"
kountry geo, from(iso2c) marker
rename NAMES_STD country
replace country = "Croatia" if country == "Yugoslavia"
drop if MARKER == 0 // Drop EU-28 etc. entries
drop if country == "European Union"
drop MARKER
// Identify EU-15
generate eu15 = (inlist(geo, "AT", "BE", "DK", "FI", "FR", "DE", "GR", "IE") ///
| inlist(geo, "IT", "LU", "NL", "PT", "ES", "SE", "UK"))
sparkline value year if eu15, over(country) ysize(10) ///
note(" " "{it: Source:} Eurostat, ilc_di12, date of extraction: 2017-07-10", span) ///
xlabel(1995 (10) 2015) xmtick(1995 (5) 2015) ///
xtick(1995 (1) 2016) ///
ytitle(" ") xtitle("") ///
title("Income inequality") ///
subtitle("(Gini of equivalised disposable income)") ///
name(regular, replace)
// 2) Download data as time series:
sdmxuse data ESTAT, dataset(ilc_di12) clear timeseries
destring time, gen(year)
twoway (scatter gini_hnd_uk_a year, connect(L)) ///
(scatter gini_hnd_de_a year, connect(L)) ///
, xtitle("") ///
ytitle("Income inequality" ///
"(Gini of equivalised disposable income)") ///
xlabel(1995 (5) 2015) xtick(1995 (1) 2015) ///
ylabel(25 (5) 35) ytick(25 (1) 35) ///
note(" " "{it: Source:} Eurostat, ilc_di12, date of extraction: 2017-07-10", span) ///
legend(order(1 "United Kingdom" 2 "Germany") pos(5) ring(0)) ///
name(timeseries, replace)
Jun 20, 2017
Random graphs (100): Labeling lines directly
import delimited une_rt_a_1_Data.csv, clear
/*
DATASET: Unemployment by sex and age - annual average [une_rt_a]
LAST UPDATE: 14.06.17 13:10:19
EXTRACTION DATE: 18.06.17 11:35:33
SOURCE OF DATA: Eurostat
*/
drop sex age unit flagandfootnotes
// Generate country variable
kountryadd "Germany (until 1990 former territory of the FRG)" to "Germany" add
kountry geo, from(other) stuck marker
ren _ISO3N_ country
kountry country, from(iso3n) to(iso2c)
ren _ISO2C_ country_str
replace country_str = "UK" if country_str == "GB"
list country_str geo
// Fix unemployment rate
replace value = "" if value == ":" // Fix missing data indicator
destring value, replace // Convert to numeric
twoway (line value time if country_str == "UK") ///
(line value time if country_str == "US") ///
(scatteri 4 2016 "United States" 3.5 2016 "United Kingdom", msymbol(none)) ///
, legend(off) ///
xtitle("") ytitle("Male unemployment rate, 25-74 y.") ///
xlabel(1985(5)2015) ///
xscale(range(1983 2023)) ///
note(" " "{it: Source:} Eurostat, une_rt_a, date of extraction: 2017-06-18", span)
Jun 19, 2017
Random graphs (99): Scatterplots
regress treatments index
local r2 = round(e(r2), .01) // Round coeff.
twoway (scatter treatments index) ///
(lfit treatments index) ///
(scatter treatments index if inlist(cntry, "IT", "AT", "UA", "BE", "DK", "RU"), msymbol(o)) ///
, xtitle(ART comprehensive availability) ytitle("ART treatments" "per million women 15-44 y.") ///
title("{bf:B}", justification(left) bexpand span) ///
xlabel(0 (1) 9) ///
text(1500 6.6 "Denmark") ///
text(1500 8.6 "Belgium") ///
text(470 0.4 "Italy") ///
text(400 1.6 "Austria") ///
text(120 7.5 "Ukraine") ///
text(210 8.0 "Kazakhstan", placement(west)) ///
text(150 8.0 "Russia", placement(east)) ///
legend(order(2 "Linear fit R{char 178} = `r2'") ring(0) pos(7)) name(policy, replace)
Feb 26, 2017
Chapter 2 of Singer and Willett's (2003) book on longitudinal data analysis
There is a different take on the chapter here, but I like mine better.
// Open data set use http://www.ats.ucla.edu/stat/stata/examples/alda/data/tolerance, clear // Figure 2.1 list, sep(0) noobs reshape long tol, i(id) j(age) list if inlist(id, 9, 45, 1653), sepby(id) noobs // Table 2.1 reshape wide quietly estpost cor tol*, matrix listwise esttab, unstack not noobs compress nostar b(%6.2f) // Figure 2.2
reshape long twoway (scatter tol age), by(id) ytitle(Tolerance) xtitle(Age) name(Figure22, replace) // Figure 2.3
twoway (scatter tol age) ///
(lowess tol age) ///
, by(id, note("Graphs by id, lowess curves") legend(off)) ///
ytitle(Tolerance) xtitle(Age) name(Figure23, replace)
// Table 2.2
preserve
tempfile table22
generate time = age - 11
statsby initial = _b[_cons] ///
initial_se = _se[_cons] ///
change = _b[time] ///
change_se = _se[time] ///
residual = (e(rmse)^2) ///
explained = e(r2) ///
, by(id) saving(`table22', replace): regress tol time
drop time
reshape wide
merge 1:1 id using `table22', nogenerate
format initial initial_se change change_se residual explained exposure %6.2f
keep id initial initial_se change change_se residual explained male exposure
order id initial initial_se change change_se residual explained male exposure
list, sep(0) noobs abbrev(15)
tempfile figure28
save `figure28', replace
// Figure 2.4
stem initial, round(.01)
stem change, round(.01)
stem residual, round(.01)
stem explained, round(.01)
restore
// Figure 2.5
twoway (scatter tol age) ///
(lfit tol age) ///
, by(id, note("Graphs by id, OLS curves") legend(off)) ///
ytitle(Tolerance) xtitle(Age) name(Figure25, replace)
// Figure 2.6
xtset id age
xtline tol, overlay t(age) i(id) legend(off) ///
ytitle(Tolerance) ylabel(1 (1) 4) ///
addplot(lowess tol age, lwidth(thick) lpattern(solid)) ///
xtitle(Age) title("Observed data and lowess smoother") ///
xsize(3) name(Figure26A, replace) nodraw
quietly regress tol i.id##c.age
predict tolhat
xtline tolhat, overlay t(age) i(id) legend(off) ///
ytitle(Tolerance) ylabel(1 (1) 4) ///
addplot(lfit tol age, lwidth(thick) lpattern(solid)) ///
xtitle(Age) title("OLS trajectories") ///
xsize(3) name(Figure26B, replace) nodraw
graph combine Figure26A Figure26B, col(2) name(Figure26, replace)
// Table 2.2
preserve
use `table22', clear
tabstat initial change, stat(mean sd) format(%6.2f)
cor initial change
restore
// Figure 2.7
quietly sum exposure, detail
generate highexposure = (exposure >= r(p50))
generate lowexposure = (exposure < r(p50))
generate female = !male
label var male "Males"
label var female "Females"
label var lowexposure "Low exposure"
label var highexposure "High exposure"
foreach z of varlist female male lowexposure highexposure {
capture drop tolhat
quietly regress tol i.id##c.age if `z'
predict tolhat if `z'
xtline tolhat if `z', ///
overlay t(age) i(id) legend(off) ///
ytitle("Predicted tolerance") ylabel(1 (1) 4) ///
addplot(lfit tol age if `z', lwidth(thick) lpattern(solid)) ///
xtitle(Age) title(`: variable label `z'') ///
xsize(3) ysize(2) name(Figure27_`z', replace) nodraw
}
graph combine Figure27_female Figure27_male Figure27_lowexposure Figure27_highexposure, col(2) name(Figure27, replace)
// Figure 2.8
use `figure28', clear
quietly correlate initial male
local pmca = round(r(rho), .01)
quietly correlate initial exposure
local pmcb = round(r(rho), .01)
quietly correlate change male
local pmcc = round(r(rho), .01)
quietly correlate change exposure
local pmcd = round(r(rho), .01)
label define male -1 " " 1 "Male" 0 "Female" 2 " ", modify
label val male male
twoway dot initial male, ///
ytitle("Predicted intercept") ///
xtitle("Gender") ///
xlabel(-1 0 1 2, val) ///
ysize(2) xsize(3) ///
note("{it:r} = `pmca'") ///
name(Figure28A, replace) nodraw
twoway scatter initial exposure, ///
ytitle("Predicted intercept") ///
xtitle("Exposure") ///
ysize(2) xsize(3) ///
note("{it:r} = `pmcb'") ///
name(Figure28B, replace) nodraw
twoway dot change male, ///
ytitle("Predicted change") ///
xtitle("Gender") ///
ysize(2) xsize(3) ///
xlabel(-1 0 1 2, val) ///
note("{it:r} = `pmcc'") ///
name(Figure28C, replace) nodraw
twoway scatter change exposure, ///
ytitle("Predicted change") ///
xtitle("Exposure") ///
ysize(2) xsize(3) ///
note("{it:r} = `pmcd'") ///
name(Figure28D, replace) nodraw
graph combine Figure28A Figure28B Figure28C Figure28D, col(2) row(2) name(Figure28, replace)
Reference
Singer, Judith D., and John B. Willett. 2003. Applied Longitudinal Data Analysis. Modeling Change and Event Occurrence. Oxford University Press. doi: 10.1093/acprof:oso/9780195152968.001.0001
May 28, 2016
Random graphs (87): Line plots
unzipfile ZA4582_v1-0-0.dta.zip
use ZA4582_v1-0-0.dta, clear
// Prepare some variables
recode v442 (0 = .a "Not asked") ///
(1 = 3 "Daily") ///
(2 = 2 "At least once per week") ///
(3 = 1 "At least once per month") ///
(4 5 = 0 "Less often/never") ///
(9 = .b "No answer"), generate(internet_inc1998)
recode v462 (0 = .a "Not asked") ///
(1 = 3 "Daily") ///
(2 = 2 "At least once per week") ///
(3 = 1 "At least once per month") ///
(4 5 = 0 "Less often/never") ///
(9 = .b "No answer"), generate(internet_inc2014)
recode v1910 (0 = .a "Not asked") ///
(1 = 3 "Daily") ///
(2 = 8 "Several times a week") ///
(3 = 2 "At least once per week") ///
(4 = 1 "At least once per month") ///
(5 = 0 "Less than once per week") ///
(9 = .b "No answer"), generate(internet_inc2006)
recode v1909 (0 = .a "Missing") ///
(1 = 1 "Yes") ///
(2 = 0 "No"), gen(internet_prev2006)
// Daily and weekly usage
clonevar internet_inc = internet_inc1998
replace internet_inc = internet_inc2014 if year == 2014
preserve
qui tab internet_inc, gen(usage)
rename usage4 dailyusage
rename usage3 weeklyusage
collapse dailyusage weeklyusage, by(year)
replace dailyusage = dailyusage * 100
replace weeklyusage = weeklyusage * 100
drop if missing(dailyusage)
twoway (scatter weeklyusage year, connect(L)) ///
(scatter dailyusage year, connect(L)) ///
, legend(order(1 "At least once per week" ///
2 "Daily") title("Internet usage", size(*.8)) ring(0) pos(11)) ///
xlabel(1998 2004 2014) xtitle("") ytitle(%) ylabel(0 (10) 60) ///
note(" " "{it:Source:} German General Social Survey (Allbus), doi: 10.4232/1.12439" ///
"{it:Note:} Based on recodings of v442 and v462.", span) name(internet1, replace)
restore
// Daily internet usage based on different response formats
replace internet_inc = internet_inc2006 if year == 2006 | year == 2008
preserve
qui tab internet_inc, gen(usage)
rename usage4 dailyusage
collapse dailyusage, by(year)
replace dailyusage = dailyusage * 100
drop if missing(dailyusage)
twoway (scatter dailyusage year, connect(L)) ///
, xlabel(1998 2004 2006 2008 2014) xtitle("") ytitle("% of respondents reporting" "{bf:daily} internet usage") ylabel(0 (10) 60) ///
note(" " "{it:Source:} German General Social Survey (Allbus), doi: 10.4232/1.12439" ///
"{it:Note:} Based on recodings of v442, v462, and v1910.", span) name(internet2, replace)
restore
// Yes/no question
preserve
collapse internet_prev, by(year)
drop if missing(internet_prev)
replace internet_prev = internet_prev * 100
twoway (scatter internet_prev year, connect(L)) ///
, xlabel(2006 2008 2010 2012 2014) xtitle("") ytitle("% of respondents using the internet") ylabel(0 (10) 80) ///
note(" " "{it:Source:} German General Social Survey (Allbus), doi: 10.4232/1.12439. {it:Note:} Based on v1909.", span) ///
name(internet3, replace)
restore
erase ZA4582_v1-0-0.dta
May 24, 2016
Random graphs (84): Line plot
unzipfile ZA4582_v1-0-0.dta.zip
use ZA4582_v1-0-0.dta, clear
recode v126 (0 = .a "Not asked") ///
(8 = .b "Don't know") ///
(9 = .c "No answer") ///
(1 = 3 "Completely agree") ///
(2 = 2 "Mostly agree") ///
(3 = 1 "Mostly disagree") ///
(4 = 0 "Completely disagree") ///
, gen(goodlife)
recode goodlife (3 = 100) (2 1 0 = 0), gen(completelygood)
generate east = (ost_west == 2)
collapse completelygood, by(year east)
twoway (scatter completelygood year if east == 1, connect(L)) ///
(scatter completelygood year if east == 0, connect(L)) ///
, legend(order(2 "Western Germany" 1 "Eastern Germany") pos(7) ring(0)) ///
title(`" 'All in all, one can live very well in a country like Germany' "') ///
ytitle("% completely agree") xlabel(1985 (5) 2010 2014) ///
xmtick(1984 (1) 2014) xtitle("") ylabel(0 (10) 60)
erase ZA4582_v1-0-0.dta
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
Sep 3, 2015
Calculate inflection point of curvilinear relationships
// Generate data clear set seed 1 set obs 500 generate e = rnormal(0,10) generate x = rnormal(0,2) generate y = (-2*x + 2*(x*x) + e) regress y c.x##c.x matrix list e(b) // Find out names of coefficients
// Plot regression coefficients
coefplot, xline(0) ///
xtitle(" " "Regression coefficients and 95% CI's") ///
drop(_cons) scheme(s1mono) ///
coeflabels(x = "x" c.x#c.x = "x{char 178}") ///
ciopts(recast(rcap)) ///
headings(x = "U-shaped effect of x") ///
name(coeffs, replace)
// Calculate bend
local bend = round(-_b[x]/(2*_b[c.x#c.x]), .01)
display "b1 = " round(_b[x], .01) _newline ///
"b2 = " round(_b[c.x#c.x], .01) _newline ///
"-b1/(2*b2) = " `bend'
predict yhat // Generate predicted values
// Plot results
twoway (line yhat x, sort(x)) ///
(scatter y x) ///
, legend(off) ///
xline(`bend') ///
ytitle("y") ///
xtitle("x") ///
name(scatter, replace) ///
title("Inflection point at `bend'")
May 21, 2015
Random graphs (47): Complex scatterplot
import excel "C:\table+graph field match countries.xlsx", sheet("Sheet2") cellrange(A1:I27) clear firstrow
renvars, lower // Switch variable names to lower case
list // Check data
// Bring variable names into systematic shape
foreach x of varlist healthwelfare education socialsciencesbusinesslaw ///
artshumanities services sciencemathcomputing ///
agricultureveterinary engineeringmanufacturing {
ren `x' var_`x'
}
list // Check data
// Reshape data
reshape long var_, i(country) j(field_str) string
// Convert proportion into percentage
replace var_ = var_ * 100
// Get field_str variable into numerical format
encode field_str, gen(field1)
label define field1 1 "Agriculture and Veterinary" ///
2 "Arts and Humanities" ///
3 "Education" ///
4 "Engineering and Manufacturing" ///
5 "Health and Welfare" ///
6 "Science, Math, and Computing" ///
7 "Services" ///
8 "Social Sciences, Business, and Law", modify
// Sort by average segregation per field
list field1 var if country == "EU average"
recode field1 (5 = 1 "Health and Welfare") ///
(3 = 2 "Education") ///
(8 = 3 "Social Sciences, Business, and Law") ///
(2 = 4 "Arts and Humanities") ///
(7 = 5 "Services") ///
(6 = 6 "Science, Math, and Computing") ///
(1 = 7 "Agriculture and Veterinary") ///
(4 = 8 "Engineering and Manufacturing") ///
, gen(field)
// Plot
sort field var_
twoway (scatter field var_ if country != "EU average") ///
(scatter field var_ if country == "EU average", connect(l) lwidth(thick)) ///
, ylabel(1/8, val) xscale(alt) ///
xtitle("Percentage women by field of study" "in European countries") ///
ytitle("") note(" " "{it:Source:} European Labor Force Survey 2011, own calculations.", span) ///
legend(order(2) label(2 "Average" "across" "countries") ring(0) pos(1))
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))
Subscribe to:
Posts (Atom)