national_pride_replication

* Stata .do file to replicate the analysis in http://kashin.guru/2014/07/01/pride/ based on WVS question on national pride
* ssc inst coefplot
* Data is available at http://www.worldvaluessurvey.org/WVSDocumentationWVL.jsp (file WVS_Longitudinal_1981-2014_stata_dta_v_2014_06_17_Beta)
cd "~/Downloads/"
use WVS_Longitudinal_1981-2014_stata_dta_v_2014_06_17_Beta.dta, clear

* Generate regressors
gen male =.
replace male = 0 if X001 > 0
replace male = 1 if X001 == 1

gen age =.
replace age = X003 if X003 > 13

gen maritalstatus =.
replace maritalstatus = X007 if X007 > 0

gen children =.
replace children = X011 if X011 >= 0

gen education =.
replace education = X025 if X025 > 0

gen employment =.
replace employment = X028 if X028 > 0

gen profession =.
replace profession = X036 if X036 > 0

gen incomescale =.
replace incomescale = X047 if X047 > 0

gen nationalpride =.
replace nationalpride = G006 if G006 > 0

gen region =.
replace region = X048 if X048 > 0

gen year = S020

gen country = S003

label values age X003
label values maritalstatus X007
label values children X011
label values education X025
label values employment X028
label values profession X036
label values incomescale X047
label values nationalpride G006
label values region X048
label values country S003

* Rescale nationalpride so that greater values indicate higher patriotism
gen nationalpride_res =.
forval i = 1/4 {
	qui replace nationalpride_res = `i' if nationalpride == 5-`i'
}

* Generate dummies
tab maritalstatus, gen(maritalstatus_d)
tab education, gen(education_d)
tab employment, gen(employment_d)
tab profession, gen(profession_d)
tab incomescale, gen(incomescale_d)
drop maritalstatus_d1 education_d1 employment_d1 profession_d1 incomescale_d1


* OLS regression
reg nationalpride_res age maritalstatus_d* education_d* employment_d* profession_d* incomescale_d* i.region if inlist(country,724,804,840,643,380,276,356) & year > 1996, r
est store pride_by_region
coefplot pride_by_region, drop(_cons age maritalstatus_d* education_d* employment_d* profession_d* incomescale_d*) ylab(, labs(tiny)) xlab(,labs(tiny)) xtitle("region fixed effect",size(tiny)) note("Dots are point estimates for region fixed effects in the regression of national pride on age, education, employment, profession, marital status, number" "of children, income scale for World Values Survey waves 4-6 (1999-2014) carried out in Ukraine, Germany, Italy, India, Russia, Spain and United States." "Grey lines are 95% CI (Huber-Eicker-White standard errors).", size(tiny) span) ysize(10)

* Robustness: ordered logit regression
ologit nationalpride_res age maritalstatus_d* education_d* employment_d* profession_d* incomescale_d* i.region i.year if inlist(country,724,804,840,643,380,276,356) & year > 1996, r
*margins, dydx(region) post
est store pride_by_region_logit
coefplot pride_by_region_logit, drop(_cons age maritalstatus_d* education_d* employment_d* profession_d* incomescale_d*) ylab(, labs(tiny)) if(\`b'<=15) eform