*! version 1.1.0 4-23-93 J. Hilbe, W Rogers(as gpredict command) STB-16 sg16.3 * RESIDUAL DIAGNOSTICS FOR GLM program define gpredict version 3.0 local varlist "req new max(1)" local in "opt" local if "opt" local weight "fweight aweight" local options "Pearson DEVianc Likelih ANscomb DBeta DPears DDevian Eta Mu Hat SPear SDev" parse "`*'" parse "`varlist'", parse(" ") cap drop `1' /* varlist created it */ local nopt = ("`pearson'"~="") + ("`devianc'"~="") + ("`likelih'"~="") + /* */ ("`anscomb'"~="") + ("`dbeta'"~="") + ("`dpears'"~="") + ("`ddevian'"~="") + /* */ ("`eta'"~="") + ("`mu'"~="") + ("`hat'"~="") + ("`spear'"~="") + ("`sdev'"~="") if (`nopt'~=1) { di in re "One and only one option required" exit 198 } tempvar myeta mymu m touse wt w yhat nbv tempvar Dpr Ddev Dbeta Presid Dresid Aresid Lresid Hat stdp pp Spear Sdev local y $S_E_depv local m $S_E_m local k $S_E_k local delta $S_E_del local disp $S_E_disp qui { gen byte `touse'=1 `if' `in' replace `touse'=0 if `y'==. | `touse'==. if "`exp'"~="" { gen float `wt' `exp' if `touse' replace `wt'=. if `wt'<=0 replace `touse'=0 if `wt'==. } else gen float `wt' = `touse' if `touse' predict `myeta' if `touse' replace `myeta' = `myeta' $S_E_po if ("`eta'"~="") { rename `myeta' `1' lab var `1' "Eta" exit } _glmilnk `myeta' `mymu' $S_E_pow `k' if ("`mu'" ~="") { if ("$S_E_fam"=="bin" & "$S_E_grp"!="") { replace `mymu' = `mymu'*`m' } rename `mymu' `1' lab var `1' "Mu" exit } local mu `mymu' gen `w' = . _glmwgt $S_E_fam $S_E_link `myeta' `mymu' `w' $S_E_pow `k' `touse' `disp' if ("$S_E_fam"=="bin" & "$S_E_grp"=="") { /* bernoulli */ gen `Presid' = (`y'-`mu')/sqrt(`mu'*(1-`mu')) gen `Dresid' = `y'*log(`y'/`mu') if `y'==1 replace `Dresid' = (1+`y')*log((1-`y')/(1-`mu')) if `y'==0 replace `Dresid' = sqrt(2*`Dresid') replace `Dresid' = -`Dresid' if (`y'-`mu')<0 predict `stdp', stdp gen `Hat' =(`stdp'^2*`mu'*(1-`mu'))/`delta' gen `Dpr' = (`Presid'/sqrt(1-`Hat'))^2 gen `Ddev'= (`Dresid'/sqrt(1-`Hat'))^2 gen `Dbeta' = (`Presid'^2*`Hat')/(1-`Hat')^2 gen `Lresid' = sqrt(`Hat'*`Dpr'+(1-`Hat')*`Ddev') replace `Lresid' = -`Lresid' if (`y'-`mu')<0 gen `Aresid' = (2.0533902*ibeta(2/3,2/3,`y'))-(2.0533902*ibeta(2/3,2/3,`mu')) replace `Aresid' = `Aresid'/(`mu'*(1-`mu'))^(1/6) gen `Spear' = `Presid'/sqrt(`delta'*(1-`Hat')) gen `Sdev' = `Dresid'/sqrt(`delta'*(1-`Hat')) } else if ("$S_E_fam"=="bin") { /* binomial */ gen `yhat' = `mu'*`m' gen `pp' = `yhat'/`m' predict `stdp',stdp gen `Hat'=(`stdp'^2 * `yhat'*(1-`mu'))/`delta' gen `Presid' = (`y'-`yhat')/sqrt(`yhat'*(1-`mu')) gen `Dresid' = `y'*log(`y'/`yhat') + (`m'-`y')*log((`m'-`y')/(`m'-`yhat')) replace `Dresid' = sqrt(2*`Dresid') replace `Dresid' = -`Dresid' if (`y'-`yhat')<0 gen `Aresid' = (2.0533902*ibeta(2/3,2/3,(`y'/`m')))-(2.0533902*ibeta(2/3,2/3,(`pp'))) replace `Aresid' = `Aresid'/(((`pp'*(1-`pp'))^(1/6))*sqrt((1-`Hat')/`m')) gen `Dbeta' = (`Presid'^2*`Hat')/(1-`Hat')^2 gen `Dpr'=`Presid'/sqrt(1-`Hat') gen `Ddev'=`Dresid'/sqrt(1-`Hat') gen `Lresid' = sqrt(`Hat'*`Dpr'^2+(1-`Hat')*`Ddev'^2) replace `Lresid' = -`Lresid' if (`y'-`mu')<0 gen `Spear' = `Presid'/sqrt(`delta'*(1-`Hat')) gen `Sdev' = `Dresid'/sqrt(`delta'*(1-`Hat')) } else if ("$S_E_fam"=="gau") { /* gaussian */ gen `Presid' = (`y'-`mu')/sqrt(`mu') * sqrt(`w') gen `Dresid' = (`y'-`mu')/sqrt(`mu') * sqrt(`w') predict `stdp',stdp gen `Hat' = `stdp'^2*`mu'^2/`delta' gen `Spear' = `Presid'/sqrt(`delta'*(1-`Hat')) gen `Sdev' = `Spear' } else if ("$S_E_fam"=="poi") { /* poisson */ gen `Presid' = (`y'-`mu')/sqrt(`mu')*sqrt(`wt') gen `Aresid' = (3/2)*(`y'^(3/2) - `mu'^(3/2)) / `mu'^(1/6) * sqrt(`wt') gen `Dresid' = sqrt(2*(`y'*log(`y'/`mu')-`y'+`mu')) * sqrt(`wt') replace `Dresid' = -`Dresid' if (`y'-`mu')<0 predict `stdp',stdp gen `Hat' = (`stdp'^2*`mu')/`delta' gen `Spear' = (`y'-`mu')/sqrt(`mu' *`delta'* (1-`Hat')) gen `Sdev' = sqrt(2*(`y'*log(`y'/`mu')-(`y'-`mu'))/(`delta'*(1-`Hat'))) replace `Sdev' = -`Sdev' if (`y'-`mu')<0 } else if ("$S_E_fam"=="gam") { /* gamma */ gen `Presid' = (`y'-`mu')/sqrt(`mu'^2) * sqrt(`wt') gen `Aresid' = (3*(`y'^(1/3) - `mu'^(1/3))/`mu'^(1/3)) * sqrt(`wt') gen `Dresid' = sqrt(2*(-log(`y'/`mu')+(`y'-`mu')/`mu')) * sqrt(`wt') replace `Dresid' = -`Dresid' if (`y'-`mu')<0 predict `stdp',stdp gen `Hat' = (`stdp'^2*`mu'^2)/`delta' gen `Spear' = (`y'-`mu')/sqrt(`mu'^2 * `delta'* (1-`Hat')) gen `Sdev' = sqrt(2*(-log(`y'/`mu')+(`y'-`mu')/`mu')/(`delta'*(1-`Hat'))) replace `Sdev' = -`Sdev' if (`y'-`mu')<0 } else if ("$S_E_fam"=="ivg") { /* inverse gaussian */ gen `Presid' = (`y'-`mu')/(sqrt(`mu'^3)) * sqrt(`wt') gen `Aresid' = (log(`y')-log(`mu'))/`mu'^(1/2) * sqrt(`wt') gen `Dresid' = sqrt((`y'-`mu')^2/(`y'*`mu'^2)) * sqrt(`wt') replace `Dresid' = -`Dresid' if (`y'-`mu')<0 predict `stdp',stdp gen `Hat' = (`stdp'^2*`mu'^3)/`delta' gen `Spear' = (`y'-`mu')/sqrt(`mu'^3 * `delta'* (1-`Hat')) gen `Sdev' = sqrt(((`y'/`mu')^2/(`y'*`mu'^2))/(`delta'*(1-`Hat'))) replace `Sdev' = -`Sdev' if (`y'-`mu')<0 } else if ("$S_E_fam"=="nb") { /* negative binomial */ gen `Presid' = (`y'-`mu')/sqrt(`mu'+`k'*`mu'* `mu') * sqrt(`wt') gen `Dresid' = sqrt(2*(`y'*log(`y'/`mu')-(1+`k'*`y')/`k' * log((1+`k'*`y')/(1+`k'*`mu')))) if `y'>0 replace `Dresid' = sqrt(2*(log(1+`k'*`mu')/`k')) if `y'==0 replace `Dresid' = -`Dresid' if (`y'-`mu')<0 gen `Aresid'=(3/(2*`k')*((1+`k'*`y')^(2/3)-(1+`k'*`mu')^(2/3))+3/2*(`y'^(2/3)-`mu'^(2/3)))/((`mu'+`k'*`mu'^2)^(1/6)) predict `stdp',stdp gen `nbv' = `mu'+`k'*`mu'*`mu' gen `Hat' = (`stdp'^2*`nbv')/`delta' gen `Spear' = (`y'-`mu')/sqrt(`nbv' * `delta'* (1-`Hat')) gen `Sdev' = sqrt(2*(`y'*log(`y'/`mu')-(1+`k'*`y')/`k' * log((1+`k'*`y')/(1+`k'*`mu')))/(`delta'*(1-`Hat'))) if `y'>0 replace `Sdev' = sqrt(2*(log(1+`k'*`mu')/`k')/(`delta'*(1-`Hat'))) if `y'==0 replace `Sdev' = -`Sdev' if (`y'-`mu')<0 } } if ("`pearson'"~="") { rename `Presid' `1' lab var `1' "Pearson residual" } else if ("`devianc'"~="") { rename `Dresid' `1' lab var `1' "Deviance residual" } else if ("`likelih'"~="") { rename `Lresid' `1' lab var `1' "Likelihood residual" } else if ("`anscomb'"~="") { rename `Aresid' `1' lab var `1' "Anscombe residual" } else if ("`hat'"~="") { rename `Hat' `1' lab var `1' "Hat" } else if ("`dbeta'"~="") { rename `Dbeta' `1' lab var `1' "Delta Beta residual" } else if ("`dpears'"~="") { rename `Dpr' `1' lab var `1' "Delta Pearson residual" } else if ("`ddevian'"~="") { rename `Ddev' `1' lab var `1' "Delta deviance residual" } else if ("`spear'"~="") { rename `Spear' `1' lab var `1' "Standardized Pearson residual" } else if ("`sdev'"~="") { rename `Sdev' `1' lab var `1' "Standardized deviance residual" } end