*! Version 1.1.2 - 27 May 2014 *! By Joaquim J.S. Ramalho * Please email jsr@uevora.pt for help and support * The software is provided as is, without warranty of any kind, express or implied, including * but not limited to the warranties of merchantability, fitness for a particular purpose and * noninfringement. In no event shall the authors be liable for any claim, damages or other * liability, whether in an action of contract, tort or otherwise, arising from, out of or in * connection with the software or the use or other dealings in the software. program define frm_reset, rclass if(c(stata_version) >= 12.1) version 12.1 else version 11.0 syntax [anything] [, Lastpower(integer 3) Wald lm lr /// DIFficult TECHnique(string) ITERate(integer 1000) from(string) /// TOLerance(real 1e-6) LTOLerance(real 1e-7) NRTOLerance(real 1e-5) /// NONRTOLerance TRace GRADient showstep HESSian SHOWTOLerance /// ml irls fisher(integer 1) search] * Preliminaries (defaults, variables, sample and possible errors) local nmod: word count `anything' if(`nmod'>=2) { display as error "Only one model at a time can be tested" exit 198 } if(`lastpower'<2) { display as error "Lastpower cannot be lower than 2" exit 198 } if(`nmod'==1) qui estimates restore `anything' capture if(`nmod'==0) qui estimates store keepresults if("`e(cmd)'"!="frm") { di as error "results for frm not found or both components of a two-part model estimated together" exit 301 } local model=e(model) if("`model'"=="1P" | "`model'"=="2Pfrac") local linkfrac=e(linkfrac) if("`model'"=="2Pbin") local linkbin=e(linkbin) local inflation=e(inflation) if("`model'"!="2Pbin" & "`lr'"!="") { di as error "LR tests are not applicable to fractional models" exit 198 } if("`model'"=="2Pbin" & "`lr'"!="" & ("`irls'"!="" | "`e(opt)'"=="irls")) { di as error "LR tests require ml estimation of both models; irls estimation is not allowed" exit 198 } if("`wald'"=="" & "`lm'"=="" & "`lr'"=="") local lm "lm" if("`lr'"!="") local LLh0=e(ll) if("`wald'"!="") { if("`e(vce)'"!="cluster") local vce `e(vce)' if("`e(vce)'"=="cluster") local vce `e(vce)' `e(clustvar)' } tempvar ybin XB touse qui gen byte `touse'=e(sample) local yobs=e(depvar) if("`model'"=="2Pbin" & `inflation'==0) qui gen `ybin'=(`yobs'>0) if("`model'"=="2Pbin" & `inflation'==1) qui gen `ybin'=(`yobs'==1) local _rhs: colnames(e(b)) local cons "_cons" local _rhs: list _rhs-cons * Tests - all versions if("`model'"=="1P") display _newline(1) "*** Fractional `linkfrac' regression model ***" if("`model'"=="2Pbin") display _newline(1) "*** Binary component of two-part model - `linkbin ***' specification ***" if("`model'"=="2Pfrac") display _newline(1) "*** Fractional component of two-part model - `linkfrac' specification ***" di in text "{hline 11}{c TT}{hline 21}" di in text %10s "Version" _col(12) "{c |} Statistic p-value" di in text "{hline 11}{c +}{hline 21}" qui predict `XB' if (`touse'), xb if("`lm'"!="") { if("`model'"=="1P" | "`model'"=="2Pfrac") local link `linkfrac' if("`model'"=="2Pbin") local link `linkbin' tempvar G g w uw gw gwx ones qui predict `G' if (`touse') if("`link'"=="cauchit") qui gen `g'=1/(_pi*(`XB'^2+1)) if (`touse') if("`link'"=="logit") qui gen `g'=exp(`XB')/((1+exp(`XB'))^2) if (`touse') if("`link'"=="probit") qui gen `g'=normalden(`XB') if (`touse') if("`link'"=="loglog") qui gen `g'=exp(-`XB')*exp(-exp(-`XB')) if (`touse') if("`link'"=="cloglog") qui gen `g'=exp(`XB')*exp(-exp(`XB')) if (`touse') qui gen `w'=sqrt(`G'*(1-`G')) if (`touse') qui replace `w'=1e-10 if `w'<1e-10 if("`model'"=="1P" | "`model'"=="2Pfrac") qui gen `uw'=(`yobs'-`G')/`w' if (`touse') if("`model'"=="2Pbin") qui gen `uw'=(`ybin'-`G')/`w' if (`touse') qui gen `gw'=`g'/`w' if (`touse') local gwx `gw' foreach var of varlist `_rhs' { tempvar XXX`var' qui gen `XXX`var''=`var'*`g'/`w' if (`touse') local gwx `gwx' `XXX`var'' } if("`model'"=="1P" | "`model'"=="2Pfrac") qui gen `ones'=1 if (`touse') } * Tests - each version forvalues i=2/`lastpower' { tempvar XB`i' qui gen `XB`i''=`XB'^`i' if (`touse') if("`wald'"!="" | "`lr'"!="") { if(`i'==2) local resetvars `XB`i'' if(`i'!=2) local resetvars `resetvars' `XB`i'' if("`linkbin'"=="cauchit") local linkbin "frm_cauchit" if("`linkfrac'"=="cauchit") local linkfrac "frm_cauchit" if("`model'"=="1P") qui frm `yobs' `_rhs' `resetvars' if (`touse'), /// model(`model') linkfrac(`linkfrac') vcefrac(`vce') /// `difficult' technique(`technique') iterate (`iterate') from(`from') /// tolerance(`tolerance') ltolerance(`ltolerance') nrtolerance(`nrtolerance') /// `nonrtolerance' `trace' `gradient' `showstep' `hessian' `showtolerance' /// `ml' `irls' fisher(`fisher') `search' if("`model'"=="2Pfrac") qui frm `yobs' `_rhs' `resetvars' if (`touse'), /// model(`model') linkfrac(`linkfrac') vcefrac(`vce') inflation(`inflation') /// `difficult' technique(`technique') iterate (`iterate') from(`from') /// tolerance(`tolerance') ltolerance(`ltolerance') nrtolerance(`nrtolerance') /// `nonrtolerance' `trace' `gradient' `showstep' `hessian' `showtolerance' /// `ml' `irls' fisher(`fisher') `search' if("`model'"=="2Pbin") qui frm `ybin' `_rhs' `resetvars' if (`touse'), /// model(`model') linkbin(`linkbin') vcebin(`vce') inflation(`inflation') /// `difficult' technique(`technique') iterate (`iterate') from(`from') /// tolerance(`tolerance') ltolerance(`ltolerance') nrtolerance(`nrtolerance') /// `nonrtolerance' `trace' `gradient' `showstep' `hessian' `showtolerance' /// `ml' `irls' fisher(`fisher') `search' if("`link'"=="frm_cauchit") local link "cauchit" if(e(converged)==1 | "`irls'"!="") { if("`wald'"!="") { qui test `resetvars' di in text %10s "Wald (`i')" _col(12) "{c |}" as result %10.3f r(chi2) _col(25) %8.4f r(p) return scalar W`i'=r(chi2) return scalar W`i'p=r(p) } if("`lr'"!="") { scalar LR=2*(e(ll)-`LLh0') scalar LRp=chi2tail(`i'-1,LR) di in text %10s "LR (`i')" _col(12) "{c |}" as result %10.3f LR _col(25) %8.4f LRp return scalar LR`i'=LR return scalar LR`i'p=LRp } } else { if("`wald'"!="") di in text %10s "Wald (`i')" _col(12) "{c |}" _col(21) as result "na" _col(31) %8.4f "na" if("`lr'"!="") di in text %10s "LR (`i')" _col(12) "{c |}" _col(21) as result "na" _col(31) %8.4f "na" } } if("`lm'"!="") { tempvar gwz r uwr qui gen `gwz'=`g'*`XB`i''/`w' if (`touse') if("`model'"=="2Pbin") { if(`i'==2) local lmvars `gwz' if(`i'!=2) local lmvars `lmvars' `gwz' qui regress `uw' `gwx' `lmvars' if (`touse'), nocons scalar LM=e(mss) } if("`model'"=="1P" | "`model'"=="2Pfrac") { qui regress `gwz' `gwx' if (`touse'), nocons qui predict `r' if (`touse'), resid qui gen `uwr'=`r'*`uw' if (`touse') if(`i'==2) local lmvars `uwr' if(`i'!=2) local lmvars `lmvars' `uwr' qui regress `ones' `lmvars' if (`touse'), nocons scalar LM=e(N)-e(rss) } scalar LMp=chi2tail(`i'-1,LM) di in text %10s "LM (`i')" _col(12) "{c |}" as result %10.3f LM _col(25) %8.4f LMp return scalar LM`i'=LM return scalar LM`i'p=LMp } } di in text "{hline 11}{c BT}{hline 21}" * Managing memory if(`nmod'==0) qui estimates restore keepresults if(`nmod'==1) qui estimates restore `anything' scalar drop _all end