/* ex7_5.do

   example 7.5

   data required: psa_dre.dta
      Smith - Prostate Ca data: paired PSA & DRE cohort data

   last update:  01 Nov 2005
*/
version 9
set more off
cap semt_profile
cap log close
*************************************
cap pro drop doit
pro def doit
    log using ${semt_log}ex7_5, replace
    use ${semt_data}ex7_5, clear
    di
    di "{txt}example 7.5"
    di
    di "           Table 7.5"
    table ya yb d [fw=nobs], m by(site)
    di

    qui {
        gen byte x_s1 = (site==1)
        gen byte x_s2 = (site==2)

        expand = nobs
        gen patid=_n
        expand=2       // generate a separate record for each test

        bys patid: gen byte test = _n-1
        gen u=(d==1)*(ya==1) if test == 1
        replace  u=(d==1)*(yb==1) if test == 0
        gen x_s1t = x_s1*test
        gen x_s2t = x_s2*test
    }

    di "{txt}*************
    di "{txt} Loglinear model of detection probability;  full model
    di "{txt}     with site x test type interaction terms
    di
    di "{txt} log DP = alpha_0 + alpha_1*X_s1 + alpha_2*X_s2 + beta_0*X_test "
    di "{txt}          + beta_1*X_s1*X_test + beta_2*X_s2*X_test : "

    glm u x_s1 x_s2 test x_s1t x_s2t, family(binomial) link(log) cluster(patid) nolog
    di
    di "{txt} exponentiated coefficient estimates: "
    glm, eform noheader
    di
    di "{txt} joint test for contribution of site x test type interaction:"
    testparm x_s1t x_s2t
    di
    di "{txt}*************
    di "{txt} reduced model
    di "{txt}  log DP = alpha_0 + alpha_1*X_s1 + alpha_2*X_s2 + beta_0*X_test: "

    glm u x_s1 x_s2 t  ,family(binomial) link(log) cluster(patid) nolog
    di
    di "{txt} exponentiated coefficient estimates: "
    glm, eform noheader

    qui log close
end
*************************
doit