capture log close
log using gcpaper04-model3-tab4-fig910, replace text

//  program:    gcpaper04-model3-tab4-fig910.do
//  task:       model 3: table 4 & figures 9 and 10
//              full model
//  project:    Group comparisons in logit and probit
//              using predicted probabilities
//  author:     jsl \ 2009-05-26

//  #0
//  program setup

version 10
clear all 
set linesize 80
matrix drop _all
set scheme s2manual

local pgm "gcpaper04"
local dte "2009-05-26"

//  #1
//  load the tenure data

use gc_tenure01, clear
datasignature confirm
keep if allisonsmr // use allison sample

//  #2
//  Table 4 - Model 4: all variables by gender

//  #2a - matrix for table 4

local nrows = 8
matrix table3 = J(`nrows',6,.)
matrix colnames table3 = Wb Wexpb Wz Mb Mexpb Mz
matrix rownames table3 = Constant year yearsq select articles prestige Log-likelihood N

//  #2b - create table 4

* for women
logit tenure year yearsq select articles prestige if female==1, nolog
logit, or
    mat Mb = e(b)
    mat MV = e(V)
    mat table3[1,1] = Mb[1,6]                   // constant
    mat table3[2,1] = Mb[1,1]                   // b-year
    mat table3[3,1] = Mb[1,2]                   // b-yearsq
    mat table3[4,1] = Mb[1,3]                   // b-select
    mat table3[5,1] = Mb[1,4]                   // b-articles
    mat table3[6,1] = Mb[1,5]                   // b-prestige
    mat table3[2,2] = exp(Mb[1,1])              // odds ratios
    mat table3[3,2] = exp(Mb[1,2])
    mat table3[4,2] = exp(Mb[1,3])
    mat table3[5,2] = exp(Mb[1,4])
    mat table3[6,2] = exp(Mb[1,5])
    mat table3[1,3] = Mb[1,6]/sqrt(MV[6,6])     // z-test
    mat table3[2,3] = Mb[1,1]/sqrt(MV[1,1])
    mat table3[3,3] = Mb[1,2]/sqrt(MV[2,2])
    mat table3[4,3] = Mb[1,3]/sqrt(MV[3,3])
    mat table3[5,3] = Mb[1,4]/sqrt(MV[4,4])
    mat table3[6,3] = Mb[1,5]/sqrt(MV[5,5])
    mat table3[7,1] = e(ll)                     // LL
    mat table3[8,1] = e(N)                      // N
* for men
logit tenure year yearsq select articles prestige if female==0, nolog
logit, or
    mat Mb = e(b)
    mat MV = e(V)
    mat table3[1,4] = Mb[1,6]                   // constant
    mat table3[2,4] = Mb[1,1]                   // b-year
    mat table3[3,4] = Mb[1,2]                   // b-yearsq
    mat table3[4,4] = Mb[1,3]                   // b-select
    mat table3[5,4] = Mb[1,4]                   // b-articles
    mat table3[6,4] = Mb[1,5]                   // b-prestige
    mat table3[2,5] = exp(Mb[1,1])              // odds ratios
    mat table3[3,5] = exp(Mb[1,2])
    mat table3[4,5] = exp(Mb[1,3])
    mat table3[5,5] = exp(Mb[1,4])
    mat table3[6,5] = exp(Mb[1,5])
    mat table3[1,6] = Mb[1,6]/sqrt(MV[6,6])     // z-test
    mat table3[2,6] = Mb[1,1]/sqrt(MV[1,1])
    mat table3[3,6] = Mb[1,2]/sqrt(MV[2,2])
    mat table3[4,6] = Mb[1,3]/sqrt(MV[3,3])
    mat table3[5,6] = Mb[1,4]/sqrt(MV[4,4])
    mat table3[6,6] = Mb[1,5]/sqrt(MV[5,5])
    mat table3[7,4] = e(ll)                     // LL
    mat table3[8,4] = e(N)                      // N

//  #2c - print table 4

mat list table3, format(%9.3f)

//  #3
//  Model 3 using all predictors

//  #3a - gcsetup & gcmodel

gcsetup year yearsq select articles prestige, details grpvar(male) ///
    gr0var(f0) gr0label(women) gr1var(m1) gr1label(men)

gcmodel year yearsq select articles prestige, label(M3_full) details  

//  #3b - estimate model using global rhsinteract from -gcmodel-

logit tenure $rhsinteract, nocon nolog

//  #4
//  Figure 9 - Model 3

//  #4a - predictions over articles given prestige with year = 7

foreach prst in 1 2 3 4 5 { // loop over prestige

    gcprgen articles, from(0) to(50) gap(2) ///
        gen(M3prst`prst')  x(year=7 yearsq=47 prestige=`prst')
    label var M3prst`prst'Value  "Number of Articles"
    label var M3prst`prst'men    "Men"
    label var M3prst`prst'women  "Women"
    label var M3prst`prst'Diff   "Male-Female difference"
    label var M3prst`prst'Lower  "95% confidence interval"
    label var M3prst`prst'Upper  "95% confidence interval"

    * create data for significant and non-signifcant cases
    gen D`prst'nonsig = M3prst`prst'Diff 
    gen D`prst'sig = M3prst`prst'Diff if M3prst`prst'Sig_95==1
    label var D`prst'sig "Significant difference"
    label var D`prst'nonsig "Non-significant difference"

} // loop over prestige

//  #4b - plot dc over articles for levels of prestige

local title "In year 7 with average selectivity"

local plttag "`pgm'.do jsl `dte'"
local pltfil "`pgm'_F9_M3_dc_byprst_overart.png"
local caption "Source: `pltfil' from `plttag', size(tiny)"

label var D1sig "Weak (prestige=1)"
label var D2sig "Adequate (prestige=2)"
label var D3sig "Good (prestige=3)"
label var D4sig "Strong (prestige=4)"
label var D5sig "Distinguished (prestige=5)"
local optlegend "pos(11) order(1 3 5 7 9) ring(0) cols(1) region(ls(none))"

//  colors, lines, and symbols

* vvthin vthin thin medthin medium medthick thick vthick vvthick vvvthick
local w1 thin
local w2 medthin
local w3 medium
local w4 medthick
local w5 thick

* color for prestige levels
    * loc c1 red loc c2 orange loc c3 green loc c4 blue loc c5 purple loc c5 black

* line patterns
local lin1non "clpat(dash)  clwid(`w1') clcol(black)" // dash if not sig
local lin1sig "clpat(solid) clwid(`w1') clcol(black)" // solid if sig
local lin2non "clpat(dash)  clwid(`w2') clcol(black)"
local lin2sig "clpat(solid) clwid(`w2') clcol(black)"
local lin3non "clpat(dash)  clwid(`w3') clcol(black)"
local lin3sig "clpat(solid) clwid(`w3') clcol(black)"
local lin4non "clpat(dash)  clwid(`w4') clcol(black)"
local lin4sig "clpat(solid) clwid(`w4') clcol(black)"
local lin5non "clpat(dash)  clwid(`w5') clcol(black)"
local lin5sig "clpat(solid) clwid(`w5') clcol(black)"

* symbols - invisible
local optnosym "msym(i) mcol(none)"

twoway ///
    (connected D1sig    M3prst1Value, `lin1sig' `optnosym') ///
    (connected D1nonsig M3prst1Value, `lin1non' `optnosym') ///
    (connected D2sig    M3prst1Value, `lin2sig' `optnosym') ///
    (connected D2nonsig M3prst1Value, `lin2non' `optnosym') ///
    (connected D3sig    M3prst1Value, `lin3sig' `optnosym') ///
    (connected D3nonsig M3prst1Value, `lin3non' `optnosym') ///
    (connected D4sig    M3prst1Value, `lin4sig' `optnosym') ///
    (connected D4nonsig M3prst1Value, `lin4non' `optnosym') ///
    (connected D5sig    M3prst1Value, `lin5sig' `optnosym') ///
    (connected D5nonsig M3prst1Value, `lin5non' `optnosym') ///
    , ///
     title("`title'", size(medium)) ///
     xlabel(0(10)50) ///
     ylabel(-.1(.1).6) yline(0) ytitle("`ytitledc'") ///
     legend(`optlegend') caption(`caption')
graph export `pltfil', replace width(1200)

//  #5
//  Figure 10 - Model 3

//  #5a - predictions over prestige given article with year = 7

foreach art in 0 10 20 30 40 { // loop over articles

    gcprgen prestige, from(1) to(5) gap(.25) ///
        gen(M3art`art')  x(year=7 yearsq=47 articles=`art')
    label var M3art`art'Value  "Departmental Prestige"
    label var M3art`art'men    "Men"
    label var M3art`art'women  "Women"
    label var M3art`art'Diff   "Male-Female difference"
    label var M3art`art'Lower  "95% confidence interval"
    label var M3art`art'Upper  "95% confidence interval"

    * create data for significant and non-signifcant cases
    gen D`art'nonsig = M3art`art'Diff 
    gen D`art'sig = M3art`art'Diff if M3art`art'Sig_95==1
    label var D`art'sig "Significant difference"
    label var D`art'nonsig "Non-significant difference"

} // loop over articles

//  #5b - plot dc over prestige for levels of articles

local title "In year 7 with average selectivity"

local plttag "`pgm'.do jsl `dte'"
local pltfil "`pgm'_F10_M3_dc_byart_overprst.png"
local caption "Source: `pltfil' from `plttag', size(tiny)"

label var D0sig "0 articles"
label var D10sig "10 articles"
label var D20sig "20 articles"
label var D30sig "30 articles"
label var D40sig "40 articles"
local optlegend "pos(11) order(1 3 5 7 9) ring(0) cols(1) region(ls(none))"

//  colors, lines, and symbols

* vvthin vthin thin medthin medium medthick thick vthick vvthick vvvthick
local w1 thin
local w2 medthin
local w3 medium
local w4 medthick
local w5 thick

* color for prestige levels
* loc c1 red loc c2 orange loc c3 green loc c4 blue loc c5 purple loc c5 black

* line patterns
local lin1non "clpat(dash)  clwid(`w1') clcol(black)" // dash if not sig
local lin1sig "clpat(solid) clwid(`w1') clcol(black)" // solid if sig
local lin2non "clpat(dash)  clwid(`w2') clcol(black)"
local lin2sig "clpat(solid) clwid(`w2') clcol(black)"
local lin3non "clpat(dash)  clwid(`w3') clcol(black)"
local lin3sig "clpat(solid) clwid(`w3') clcol(black)"
local lin4non "clpat(dash)  clwid(`w4') clcol(black)"
local lin4sig "clpat(solid) clwid(`w4') clcol(black)"
local lin5non "clpat(dash)  clwid(`w5') clcol(black)"
local lin5sig "clpat(solid) clwid(`w5') clcol(black)"

* symbols - invisible
local optnosym "msym(i) mcol(none)"

twoway ///
    (connected D0sig     M3art0Value, `lin1sig' `optnosym') ///
    (connected D0nonsig  M3art0Value, `lin1non' `optnosym') ///
    (connected D10sig    M3art0Value, `lin2sig' `optnosym') ///
    (connected D10nonsig M3art0Value, `lin2non' `optnosym') ///
    (connected D20sig    M3art0Value, `lin3sig' `optnosym') ///
    (connected D20nonsig M3art0Value, `lin3non' `optnosym') ///
    (connected D30sig    M3art0Value, `lin4sig' `optnosym') ///
    (connected D30nonsig M3art0Value, `lin4non' `optnosym') ///
    (connected D40sig    M3art0Value, `lin5sig' `optnosym') ///
    (connected D40nonsig M3art0Value, `lin5non' `optnosym') ///
    , ///
     title("`title'", size(medium)) ///
     xlabel(1(1)5) ///
     ylabel(-.1(.1).6) yline(0) ytitle("`ytitledc'") ///
     legend(`optlegend') caption(`caption')
graph export `pltfil', replace width(1200)

log close
exit