*! NJGW 1.1.0 17 Dec 1999 *! Program to calculate Small-Hsiao test of IIA *! Nicholas Winter (nwinter@umich.edu) program define smhsiao version 6.0 syntax varlist [if] [in] [fw pw iw], /* */ Elim(real) [ Samp(varname) Basecategory(string) Quietly Detail * ] marksample touse local y : word 1 of `varlist' /* get DV & dof for test*/ local nvar : word count of `varlist' local dof = `nvar'-1 if "`quietly'"=="quietly" { local shh2 "quietly" } if "`detail'"=="" { local shh "quietly" } if "`weight'"!="" { local weight "[`weight'`exp']" } tempname Vals Counts /* get values of DV */ qui tab `y' if `touse', matrow(`Vals') matcell(`Counts') local nvals = r(r) local i = 1 while `i' <= `nvals' { local Yval`i' = `Vals'[`i',1] local Yvals "`Yvals' `Yval`i''" if `elim'==`Yval`i'' { local ielim `i' } local i = `i' + 1 } if "`ielim'"=="" { di in r "`y' does not take on value `elim' in estimation sample" error 198 } if "`basecat'"=="" { local xm 0 local i 1 while `i'<=`nvals' { local x = `Counts'[`i',1] if `x'>`xm' { local xm `x' local Yxm = `Vals'[`i',1] /* value of y for most freq */ } local i=`i'+1 } local basecat `Yxm' } local bcat "b(`basecat')" qui tab `y' if `touse' & `y'!=`elim', matrow(`Vals') /* get values of DV, minus elim */ local nEvals = r(r) local i = 1 while `i' <= `nEvals' { local EYval`i' = `Vals'[`i',1] local EYvals "`EYvals' `EYval`i''" local i = `i' + 1 } local Ylab : val lab `y' if "`Ylab'"=="" { tempname junk local Ylab "`junk'" } local i 1 while `i'<=`nvals' { local Ylab`i' : lab `Ylab' `Yval`i'' local i=`i'+1 } tempvar lnL denom tempname b0a b0b b0ab b1b if "`samp'"=="" { tempvar samp qui gen `samp'=round(uniform(),1)+1 if `touse' local samp1 1 local samp2 2 } else { tempname SVal qui ta `samp' if `touse', matrow(`SVal') if `r(r)' != 2 { di in r "samp() variable must take exactly two values in the estimation sample" error 198 } local samp1 = `SVal'[1,1] local samp2 = `SVal'[2,1] } qui ta `y' if `touse' & `samp'==`samp1' local cat1 `r(r)' if `cat1'!=`nvals' { di in r "`y' is missing categories in the first half sample" error 148 } qui ta `y' if `touse' & `samp'==`samp2' local cat2 `r(r)' if `cat2'!=`nvals' { di in r "`y' is missing categories in the second half sample" error 148 } `shh2' mlogit `varlist' `if' `in' `weight', b(`basecat') `options' `shh' di "Estimating first-half sample model:" `shh' mlogit `varlist' if `touse' & `samp'==`samp1' `weight', b(`basecat') `options' mat `b0a' = e(b) `shh' di `shh' di "Estimating second-half sample model:" `shh' mlogit `varlist' if `touse' & `samp'==`samp2' `weight', b(`basecat') `options' mat `b0b' = e(b) mat `b0ab' = (0.70710678)*(`b0a') + (0.29289322)*(`b0b') /* Zhang & Hoffman eq. 9 */ ************************************* *Get LnL for amalgamated coefficients ************************************* *get XBs & assemble denominator qui gen double `denom' = 0 if `touse' local i 1 while `i' <= `nEvals' { /* cycle through values (w/o eliminated one) */ local cury : word `i' of `EYvals' tempvar xb`cury' if `cury' != `basecat' { matrix score double `xb`cury'' = `b0ab' if `touse', eq(`Ylab`cury'') } else { qui gen double `xb`cury''=0 /* because (exp(0)=1) */ } qui replace `denom'=`denom' + exp(`xb`cury'') if `touse' local i=`i'+1 } *create Log likelihood using amalgamated coeff. qui gen double `lnL' = . if `touse' local i 1 while `i'<=`nEvals' { local cury : word `i' of `EYvals' qui replace `lnL' = ln(exp(`xb`cury'')/(`denom')) if `y'==`cury' & `touse' local i=`i'+1 } * GET LnL, only for observations in the 2d sample, and without eliminated observations: sum `lnL' if `touse' & `samp'==`samp2' & `y'!=`elim', meanonly local lnL_0 = r(sum) * Now get b1b: coefficients from restricted choice set `shh' di `shh' di "Estimating coefficients in second sample without `y'==`elim'" `shh' mlogit `varlist' if (`touse' & `samp'==`samp2' & `y'!=`elim') `weight', `options' local lnL_1 = e(ll) local SH = -2 * (`lnL_0' - `lnL_1') local p = chiprob(`dof',`SH') di di in g "Small-Hsiao Test of IIA" di in g "-----------------------" di di in g " LnL for weighted average of coeffients in two half-samples: " in y %6.5f `lnL_0' di in g " LnL based on restricted model (`y'==`elim' omitted): " _col(63) in y %6.5f `lnL_1' di in g " Test chi2(`dof') = " in y %6.3f `SH' in g ", Pr = " in y %5.4f `p' end *end&