program define ztpflex_ll args ll xb lnsigma lntheta quietly { tempname sigma theta sca `sigma' = exp(`lnsigma') sca `theta' = exp(`lntheta') /* Get points and weights for Gaussian-Hermite quadrature. */ tempvar x w gen `x'=. gen `w'=. _GetQuad, avar(`x') wvar(`w') quad($R) replace `x' = sqrt(2)*`x' tempvar li lnpri pri lambda gen double `lambda'=. gen double `lnpri'=. gen double `pri'=. gen double `li'=0 forvalues r=1/$R{ replace `lambda' = exp(`xb'+`sigma'* sinh(`theta'*`x'[`r'])/`theta' ) replace `lnpri' = -`lambda' + $ML_y1*ln(`lambda') - ln(1-exp(-`lambda')) - lngamma($ML_y1 + 1) replace `pri' = exp(`lnpri') replace `li' = `li' + `pri'*`w'[`r'] } replace $f=`li' /*for vuong test*/ replace `ll' = ln((1/sqrt(_pi))*`li') if $ML_samp } end