version 12.1

* Stata program for citation impact/margins analysis, updated Feb 2013
* esttab from SSC needs to be installed.
* Also, fre from SSC
* Some tables will need additional formatting by hand

use "https://www3.nd.edu/~rwilliam/margins/bornmann.dta", clear
* Define vars, select sample
encode oecdma, gen(noecdma)
encode doctype, gen(ndoctype)
fre n*
drop if noecdma == 1 | noecdma == 3 | noecdma == 6
* For py2, 1980 = year 31, 2010 = 1
gen py2 = 2011 - py
* Create dummies for labeling purposes in the tables
tab1 univ, gen(University)
tab1 noecdma, gen(subject)
ren subject1 Engineering
ren subject2 MedHealth
ren subject3 NatSciences
tab1 ndoctype, gen(doctype)
rename doctype1 Article
rename doctype2 Note
ren doctype3 Proceedings
ren doctype4 Review
clonevar PPTop10 = pr2
label variable jif "Journal Impact Factor"
label variable py2 "# Years since publication" 
label variable PPTop10 "Highly Cited Publication"
label define univ 1 "Univ 1", modify
label define univ 2 "Univ 2", modify
label define univ 3 "Univ 3", modify
label define univ 4 "Univ 4", modify

* Table 1  - Descriptives
sum PPTop10 University1-University4 Engineering-NatSciences Article-Review jif py2 auth prange

*** Table 2
* Model 1 - univ only
logit PPTop10 i.univ
est store m1

* Model 2 - All vars, no transformations
logit PPTop10 i.univ i.noecdma i.ndoctype jif py2 auth prange  
est store m2

* Model 3 - Final Model with Squared Terms
logit PPTop10 i.univ i.noecdma i.ndoctype jif py2 auth prange c.jif#c.jif c.prange#c.prange
est store m3

* Create Table 2. It will still need hand-editing
esttab m1 m2 m3 using table2, z pr2 aic bic scalars(chi2 df_m) mti("Baseline" "All vars" "Squared trms") rtf

* Table 3 - AAPs and AMES, Categorical Variables
* Table will still need hand-editing
est restore m3
margins univ noecdma ndoctype, post
est store m11
est restore m3
margins, dydx(univ noecdma ndoctype) post
est store m12
esttab m11 m12 using table3, z mti(AAPS AMES) rtf title("Table 3. Average Adjusted Predictions and Average Marginal Effects") replace

*** Figures follow

* Graph 1 AAPs for JIF
est restore m3
margins , at(jif = (0 (1) 35))
marginsplot, scheme(s1mono) yti(`"Pr(Highly Cited Publication)"') title("") recastci(rarea) ci1opts(fintensity(30)) ylabel(#10) xlabel(#10)

* Graph 2 MERs for JIF
est restore m3
margins , dydx(jif) at(jif = (0 (1) 35))
marginsplot, scheme(s1mono) yti(`"Effects on Pr(Highly Cited Publication)"') title("") recastci(rarea) ci1opts(fintensity(30)) ylabel(#10) xlabel(#10)

* Graph 3 APRs for # Pages
est restore m3
margins , at(prange = (1 (1) 120))
marginsplot, scheme(s1mono) yti(`"Pr(Highly Cited Publication)"') title("") recastci(rarea) ci1opts(fintensity(30)) ylabel(#10) xlabel(#20)

* Graph 4 MERs for # Pages
est restore m3
margins , dydx(prange) at(prange = (1 (1) 120))
marginsplot, scheme(s1mono) yti(`"Effects on Pr(Highly Cited Publication)"') title("") recastci(rarea) ci1opts(fintensity(30)) ylabel(#10) xlabel(#20)

******************************************************************************

* Graph 5 APRs for univ, JIF
est restore m3
margins univ, at(jif = (0 (1) 13))
marginsplot, scheme(s1mono) yti(`"Pr(Highly Cited Publication)"') title("") ylabel(#10) xlabel(#20)

* Graph 6 MERs for Univ, JIF
est restore m3
margins r.univ, at(jif = (0 (1) 13))
marginsplot, scheme(s1mono) yti(`"Effects on Pr(Highly Cited Publication)"') title("") ylabel(#10) xlabel(#20)

* Graph 7 APRs for univ, # Pages
est restore m3
margins univ, at(prange = (1 (1) 25))
marginsplot, scheme(s1mono) yti(`"Pr(Highly Cited Publication)"') title("") ylabel(#10) xlabel(#10)

* Graph 8 MERs for univ, # Pages
est restore m3
margins r.univ, at(prange = (1 (1) 25))
marginsplot, scheme(s1mono) yti(`"Effects on Pr(Highly Cited Publication)"') title("") ylabel(#10) xlabel(#20)