*! v2.23 3/09/2021. (FRA) Correcting Bug for PIP *! v2.22 9/21/2020. (FRA) Corrected names when requesting sig_s(i) or sig_r(i). and adds sig_r(i) and sig_s(i) for models 1 and 2 *! v2.21 7/24/2020. (FRA) Fix typo from Rho_s to Rho r. *! v2.2 7/14/2020. (FRA) adds model 8 *! v2.12 7/13/2020. (FRA) corrected for problem predicting values. *! v2.11 7/13/2020. (FRA) adds Model 7. pip1-pip6 ready *! v2.1 7/12/2020. (FRA) adds Model 7. PIP1-PIP6 are yet to be added *! v2 KY predictor. (FRA) Baseline * Improved Predictor. It dropped unnecessary code, and now can be used for "Score" program define ky_p syntax anything(id="newvarname") [if] [in] , [star * ] version 15 /*[pi_s pi_r pi_w pi_v /// Probabilities of each event pi_r1 pi_r2 pi_r3 pi_s1 pi_s2 pi_s3 /// Probabilities by category pip_r1 pip_r2 pip_r3 pip_s1 pip_s2 pip_s3 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 pi_5 pi_6 pi_7 pi_8 pi_9 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 pip_5 pip_6 pip_7 pip_8 pip_9 /// Probabilities by group classification lf_1 lf_2 lf_3 lf_4 lf_5 lf_6 lf_7 lf_8 lf_9 /// lf_s1 lf_s2 lf_s3 lf_r1 lf_r2 lf_r3 /// rho_s rho_r /// RHO's (Regression to the mean factors) mean_e mean_n mean_w mean_t mean_v /// Conditional Means sig_e sig_n sig_w sig_t sig_v /// conditional Standard deviations mean_s1 mean_s2 mean_s3 sig_s1 sig_s2 sig_s3 /// intra group mean and sigma mean_r1 mean_r2 mean_r3 sig_r1 sig_r2 sig_r3 /// intra group mean and sigma bclass * ]*/ *lf_r1 lf_r2 lf_r3 lf_s1 lf_s2 lf_s3 /// These will contain the Likelihood function by S/R type *lf_1 lf_2 lf_3 lf_4 lf_5 lf_6 lf_7 lf_8 lf_9 /// These will contain the Likelihood function by combination type marksample touse, novarlist *parse "`0'", parse(",") *local opts `3' ** ky_p_1 to ky_p_4 still missing. if "`star'"!="" { ky_star, `options' prefix(`anything') } else if e(method_c)==1 { ky_p_1 `0' } else if e(method_c)==2 { ky_p_2 `0' } else if e(method_c)==3 { ky_p_3 `0' } else if e(method_c)==4 { ky_p_4 `0' } else if e(method_c)==5 { ky_p_5 `0' } else if e(method_c)==6 { ky_p_6 `0' } else if e(method_c)==7 { ky_p_7 `0' } else if e(method_c)==8 { ky_p_8 `0' } else { display in red "error, method not allowed" exit 1 } end ** Note: I think it is now far more efficient! program define ky_p_1 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s /// Probabilities of each event pi_s1 pi_s2 /// Probabilities by category pip_s1 pip_s2 /// Probabilities by category, on S only. rho_s /// RHO's (Regression to the mean factors) mean_e mean_n /// Conditional Means sig_r1 sig_s1 sig_s2 /// Conditional Sig for r1 s1 s2 sig_e sig_n bclass_s * ] // conditional Standard deviations marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** This are all Basic variables. Everything else is constructed from this ones tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r qui { ** First Options that are not included above. if "`options'"!="" { ml_p `0' exit } ** then options that only depend on a single parameter else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } *** Sigmas else if "`sig_e'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=exp(`_sig_e') if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { _predict double `_sig_n' if `touse', xb eq(ln_sig_n) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=exp(`_sig_n') if `touse' label var `varlist' "Conditional stdev for component n" } else if "`rho_s'" != "" { _predict double `_rho_s' if `touse', xb eq(arho_s) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=tanh(`_rho_s') if `touse' label var `varlist' "Rho s: RTM Survey data" } *** pi_s else if "`pi_s'" != "" { _predict double `_pi_s' if `touse', xb eq(lpi_s) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=invlogit(`_pi_s') if `touse' label var `varlist' "A priori Pi_s" } ** From here we have the ones that derive from everything else else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pi_s replace `varlist'=1-`varlist' label var `varlist' "Latent Class pi s type 2" } ** means by type else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean mean_r1" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean mean_s1" } else if "`mean_s2'"!="" { _predict double `_mean_e' if `touse', xb eq(mu_e) _predict double `_mean_n' if `touse', xb eq(mu_n) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean mean_s2" } else if "`sig_r1'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=exp(`_sig_e') if `touse' label var `varlist' "Conditional sigma sig_r1" } else if "`sig_s1'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=exp(`_sig_e') if `touse' label var `varlist' "Conditional sigma sig_s1" } else if "`sig_s2'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) _predict double `_sig_n' if `touse', xb eq(ln_sig_n) _predict double `_rho_s' if `touse', xb eq(arho_s) replace `_sig_e'=exp(`_sig_e') replace `_sig_n'=exp(`_sig_n') replace `_rho_s'=tanh(`_rho_s') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2) if `touse' label var `varlist' "Conditional sigma sig_s2" } else if "`pip_s1'" != "" { ky_p double `_pi_s' if `touse', pi_s ky_p double `_mean_e' if `touse', mean_e ky_p double `_mean_n' if `touse', mean_n ky_p double `_sig_e' if `touse', sig_e ky_p double `_sig_n' if `touse', sig_n ky_p double `_rho_s' if `touse', rho_s ** LL tempvar lls1 lls2 gen double `lls1'=normalden(`_s',`_mean_e',`_sig_e') gen double `lls2'=normalden(`_s',`_mean_e'+`_mean_n',sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(`_pi_s'*`lls1')/(`_pi_s'*`lls1'+(1-`_pi_s')*`lls2') label var `varlist' "Posterior Latent Class pi S-type 1" } else if "`pip_s2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pip_s1 replace `varlist'=1-`varlist' label var `varlist' "Posterior Latent Class pi s type 2" } else if "`bclass_s'"!="" { tempvar _pip_s1 ky_p double `_pip_s1' if `touse', pip_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1+`_pip_s1'<.5 if `touse' label var `varlist' "two-step classification based on S only" } } // End qui end *********************** *********************** *********************** program define ky_p_2 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_w /// Probabilities of each event pi_s1 pi_s2 pi_s3 /// Probabilities by category pip_s1 pip_s2 pip_s3 /// Probabilities by category. Posterior Prob rho_s /// RHO's (Regression to the mean factors) mean_e mean_n mean_w /// Conditional Means sig_e sig_n sig_w /// conditional Standard deviations sig_r1 sig_s1 sig_s2 sig_s3 /// Conditional Sig for r1 s1 s2 s3 bclass_s * ] // two step classification if only S is observed marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r qui { if "`options'"!="" { ml_p `0' exit } else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=(invlogit(xb(lpi_w))) if `touse' label var `varlist' "A priori Pi_w" } else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_w'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_w) label var `varlist' "Conditional mean for component w" } else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] _predict double `_sig_e', xb eq(ln_sig_e) predictnl `typlist' `varlist'=exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] _predict double `_sig_n', xb eq(ln_sig_n) gen `typlist' `varlist'=exp(`_sig_n') if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_w'"!="" { syntax newvarname [if] [in] [, * ] _predict double `_sig_w', xb eq(ln_sig_w) gen `typlist' `varlist'=exp(`_sig_w') if `touse' label var `varlist' "Conditional stdev for component w" } else if "`rho_s'"!="" { syntax newvarname [if] [in] [, * ] _predict double `_rho_s', xb eq(arho_s) gen `typlist' `varlist'=tanh(`_rho_s') if `touse' label var `varlist' "Rho s: RTM Survey data" } else if "`sig_r1'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=exp(`_sig_e') if `touse' label var `varlist' "Conditional sigma sig_r1" } else if "`sig_s1'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=exp(`_sig_e') if `touse' label var `varlist' "Conditional sigma sig_s1" } else if "`sig_s2'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) _predict double `_sig_n' if `touse', xb eq(ln_sig_n) _predict double `_rho_s' if `touse', xb eq(arho_s) replace `_sig_e'=exp(`_sig_e') replace `_sig_n'=exp(`_sig_n') replace `_rho_s'=tanh(`_rho_s') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2) if `touse' label var `varlist' "Conditional sigma sig_s2" } else if "`sig_s3'"!="" { _predict double `_sig_e' if `touse', xb eq(ln_sig_e) _predict double `_sig_n' if `touse', xb eq(ln_sig_n) _predict double `_sig_w' if `touse', xb eq(ln_sig_w) _predict double `_rho_s' if `touse', xb eq(arho_s) replace `_sig_e'=exp(`_sig_e') replace `_sig_n'=exp(`_sig_n') replace `_sig_w'=exp(`_sig_w') replace `_rho_s'=tanh(`_rho_s') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2) if `touse' label var `varlist' "Conditional sigma sig_s3" } *** A priori Prob *** Prior by Group *** This may depend on the model. Assume we start with model 5. And go backward else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { syntax newvarname [if] [in] [, * ] ky_p double `_pi_s' , pi_s ky_p double `_pi_w' , pi_w gen `typlist' `varlist'=(1-`_pi_s')*(1-`_pi_w') if `touse' label var `varlist' "Latent Class pi s type 2" } else if "`pi_s3'" != "" { syntax newvarname [if] [in] [, * ] ky_p double `_pi_s' , pi_s ky_p double `_pi_w' , pi_w gen `typlist' `varlist'=(1-`_pi_s')*( `_pi_w') if `touse' label var `varlist' "Latent Class pi s type 2" } else if "`pi_1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s1 label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s2 label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s2 label var `varlist' "Latent Class pi_3" } *** posterior over Stype. else if "`pip_s1'`pip_s2'`pip_s3'" != "" { ky_p double `_pi_s' , pi_s ky_p double `_pi_w' , pi_w local _pi_s1 ( `_pi_s') local _pi_s2 ((1-`_pi_s')*(1-`_pi_w')) local _pi_s3 ((1-`_pi_s')*( `_pi_w')) ky_p double `_mean_e' , mean_e ky_p double `_mean_n' , mean_n ky_p double `_mean_w' , mean_w ky_p double `_sig_e' , sig_e ky_p double `_sig_n' , sig_n ky_p double `_sig_w' , sig_w ky_p double `_rho_s' , rho_s tempvar _mean_s1 _mean_s2 _mean_s3 tempvar _sig_s1 _sig_s2 _sig_s3 gen double `_mean_s1'= (`_mean_e') gen double `_mean_s2'= (`_mean_e'+`_mean_n') gen double `_mean_s3'= (`_mean_e'+`_mean_n'+`_mean_w') gen double `_sig_s1' =(`_sig_e') gen double `_sig_s2' =(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) gen double `_sig_s3' =(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2)) tempvar _ls1 _ls2 _ls3 gen double `_ls1'= (normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_ls2'= (normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_ls3'= (normalden(`_s',`_mean_s3',`_sig_s3')) syntax newvarname [if] [in] [, * ] if "`pip_s1'"!="" { gen `typlist' `varlist'=(`_pi_s1'*`_ls1')/(`_pi_s1'*`_ls1'+`_pi_s2'*`_ls2'+`_pi_s3'*`_ls3') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'"!="" { gen `typlist' `varlist'=(`_pi_s2'*`_ls2')/(`_pi_s1'*`_ls1'+`_pi_s2'*`_ls2'+`_pi_s3'*`_ls3') if `touse' label var `varlist' "Posterior pi s type 2" } else if "`pip_s3'"!="" { gen `typlist' `varlist'=(`_pi_s3'*`_ls3')/(`_pi_s1'*`_ls1'+`_pi_s2'*`_ls2'+`_pi_s3'*`_ls3') if `touse' label var `varlist' "Posterior pi s type 3" } } ********************* else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'" != "" { ky_p double `_pi_s' , pi_s ky_p double `_pi_w' , pi_w local _pi_s2 ((1-`_pi_s')*(1-`_pi_w')) local _pi_s3 ((1-`_pi_s')*( `_pi_w')) ky_p double `_mean_e' , mean_e ky_p double `_mean_n' , mean_n ky_p double `_mean_w' , mean_w ky_p double `_sig_e' , sig_e ky_p double `_sig_n' , sig_n ky_p double `_sig_w' , sig_w ky_p double `_rho_s' , rho_s tempvar _mean_s1 _mean_s2 _mean_s3 gen double `_mean_s1'= (`_mean_e') gen double `_mean_s2'= (`_mean_e'+`_mean_n') gen double `_mean_s3'= (`_mean_e'+`_mean_n'+`_mean_w') gen double `_sig_s1' =(`_sig_e') gen double `_sig_s2' =(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) gen double `_sig_s3' =(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2)) tempvar _rho_r1s2 _rho_r1s3 gen double `_rho_r1s2' = ((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2') gen double `_rho_r1s3' = ((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s3') tempvar _std_r1 _std_s1 _std_s2 _std_s3 gen double `_std_r1' =((`_r'-`_mean_s1')/`_sig_s1') *gen double `_std_s1' =((`_s'-`_mean_s1')/`_sig_s1') gen double `_std_s2' =((`_s'-`_mean_s2')/`_sig_s2') gen double `_std_s3' =((`_s'-`_mean_s3')/`_sig_s3') tempvar _lf_2 _lf_3 gen double `_lf_2' = (1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_lf_3' = (1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(`_pi_s2'*`_lf_2')/(`_pi_s2'*`_lf_2'+`_pi_s3'*`_lf_3') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } else if "`pip_3'" != "" { tempvar _pip_2 ky_p double `_pip_2' , pip_2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1-`_pip_2' if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } else if "`bclass_s'"!="" { tempvar _pip_s1 _pip_s2 ky_p double `_pip_s1' , pip_s1 ky_p double `_pip_s2' , pip_s2 local _pip_s3 (1-`_pip_s1'-`_pip_s2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_s`i''>`_high' replace `_high' =`_pip_s`i'' if `_pip_s`i''>`_high' } label var `varlist' "two-step classification based on S only" } } end *********************** *********************** *********************** program define ky_p_3 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_r /// Probabilities of each event pi_r1 pi_r2 pi_s1 pi_s2 /// Probabilities by category pip_r1 pip_r2 pip_s1 pip_s2 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 /// Probabilities by group classification rho_s /// RHO's (Regression to the mean factors) mean_e mean_n mean_t /// Conditional Means sig_e sig_n sig_t /// conditional Standard deviations mean_r1 mean_r2 mean_s1 mean_s2 /// mean by group sig_r1 sig_r2 sig_s1 sig_s2 /// sigma by group bclass_r bclass_s bclass * ] marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r qui { if "`options'"!="" { ml_p `0' exit } *** A priori Prob else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_r))) if `touse' label var `varlist' "A priori Pi_r" } **************** else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_t'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_t) label var `varlist' "Conditional mean for component t" } *** Sigmas else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_n)) if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_t'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_t)) if `touse' label var `varlist' "Conditional stdev for component t" } else if "`rho_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_s)) if `touse' label var `varlist' "Rho s: RTM Survey data" } else if "`pi_r1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r label var `varlist' "Latent Class pi r type 1" } else if "`pi_r2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r replace `varlist'=1-`varlist' label var `varlist' "Latent Class pi r type 2" } else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s replace `varlist'=1-`varlist' label var `varlist' "Latent Class pi s type 2" } *** dependent from above else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean R1" } else if "`mean_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_t label var `varlist' "Conditional mean R2" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean S1" } else if "`mean_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean S2" } *** else if "`sig_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma R1" } else if "`sig_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_t label var `varlist' "Conditional sigma R2" } else if "`sig_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma S1" } else if "`sig_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) if `touse' label var `varlist' "Conditional sigma S2" } else if "`pi_1'" != "" { ky_p double `_pi_r', pi_r ky_p double `_pi_s', pi_s syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { ky_p double `_pi_r', pi_r ky_p double `_pi_s', pi_s syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*(1-`_pi_s') if `touse' label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { ky_p double `_pi_r', pi_r ky_p double `_pi_s', pi_s syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(1-`_pi_r')*`_pi_s' if `touse' label var `varlist' "Latent Class pi_3" } else if "`pi_4'" != "" { ky_p double `_pi_r', pi_r ky_p double `_pi_s', pi_s syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(1-`_pi_r')*(1-`_pi_s') if `touse' label var `varlist' "Latent Class pi_4" } *** posterior. assuming we see only 1 thing. else if "`pip_r1'" != "" { tempvar _mean_r1 _mean_r2 _sig_r1 _sig_r2 _pi_r1 ky_p double `_mean_r1', mean_r1 ky_p double `_sig_r1', sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2', sig_r2 ky_p double `_pi_r1', pi_r1 tempvar _pip_1 _pip_2 gen double `_pip_1' = `_pi_r1'*normalden(`_r',`_mean_r1',`_sig_r1') gen double `_pip_2' =(1-`_pi_r1')*normalden(`_r',`_mean_r2',`_sig_r2') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2') if `touse' label var `varlist' "Posterior pi r type 1" } else if "`pip_r2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pip_r1 replace `varlist'=1-`varlist' label var `varlist' "Posterior pi r type 2" } else if "`pip_s1'" != "" { syntax newvarname [if] [in] [, * ] tempvar _mean_s1 _mean_s2 _sig_s1 _sig_s2 _pi_s1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_pi_s1' , pi_s1 tempvar _pip_1 _pip_2 gen double `_pip_1' = `_pi_s1' *normalden(`_s',`_mean_s1',`_sig_s1') gen double `_pip_2' =(1- `_pi_s1')*normalden(`_s',`_mean_s2',`_sig_s2') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pip_s1 replace `varlist'=1-`varlist' label var `varlist' "Posterior pi s type 2" } ************************************************************ else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'`pip_3'`pip_4'" != "" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_r2 _sig_r2 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_rho_s' , rho_s tempvar _pi_2 _pi_3 _pi_4 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_sig_e' , sig_e local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r3s1 0 local _rho_r3s2 0 tempvar _pip2 _pip3 _pip4 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip4'=`_pi_4'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s2',`_sig_s2')) syntax newvarname [if] [in] [, * ] if "`pip_2'" != "" { gen `typlist' `varlist'=`_pip2'/(`_pip2'+`_pip3'+`_pip4') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } else if "`pip_3'" != "" { gen `typlist' `varlist'=`_pip3'/(`_pip2'+`_pip3'+`_pip4') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } else if "`pip_4'" != "" { gen `typlist' `varlist'=`_pip4'/(`_pip2'+`_pip3'+`_pip4') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_4" } } **** mean and sigmas for 6 elements? else if "`bclass_s'"!="" { tempvar _pip_s1 ky_p `typlist' `_pip_s1' , pip_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1+`_pip_s1'<.5 if `touse' } else if "`bclass_r'"!="" { tempvar _pip_r1 ky_p `typlist' `_pip_r1' , pip_r1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1+`_pip_r1'<.5 if `touse' } else if "`bclass'"!="" { tempvar _high gen double `_high'=. replace `_high'=0 if `_l'==0 tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_r2 _sig_r2 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_rho_s' , rho_s tempvar _pi_2 _pi_3 _pi_4 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_sig_e' , sig_e local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r3s1 0 local _rho_r3s2 0 tempvar _pip2 _pip3 _pip4 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip4'=`_pi_4'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s2',`_sig_s2')) tempvar _pip_2 _pip_3 _pip_4 gen double `_pip_2'=`_pip2'/(`_pip2'+`_pip3'+`_pip4') if `touse' gen double `_pip_3'=`_pip3'/(`_pip2'+`_pip3'+`_pip4') if `touse' gen double `_pip_4'=`_pip4'/(`_pip2'+`_pip3'+`_pip4') if `touse' syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1 if `_l'==1 forvalues i=2/4 { replace `varlist'=`i' if `_pip_`i''>`_high' & `_l'==0 replace `_high' =`_pip_`i'' if `_pip_`i''>`_high' & `_l'==0 } } } end *********************** *********************** *********************** program define ky_p_4 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_r pi_w /// Probabilities of each event rho_s /// RHO's (Regression to the mean factors) mean_e mean_n mean_w mean_t /// Conditional Means sig_e sig_n sig_w sig_t /// conditional Standard deviations mean_r1 mean_r2 mean_s1 mean_s2 mean_s3 /// mean by group sig_r1 sig_r2 sig_s1 sig_s2 sig_s3 /// sigma by group pi_r1 pi_r2 pi_s1 pi_s2 pi_s3 /// Probabilities by category pip_r1 pip_r2 pip_s1 pip_s2 pip_s3 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 pi_5 pi_6 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 pip_5 pip_6 /// Probabilities by group classification bclass_r bclass_s bclass * ] marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r *** Moments based on parameters: qui { if "`options'"!="" { ml_p `0' exit } else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_r))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_w))) if `touse' label var `varlist' "A priori Pi_r" } **************** else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_w'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_w) label var `varlist' "Conditional mean for component w" } else if "`mean_t'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_t) label var `varlist' "Conditional mean for component t" } *** Sigmas else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_n)) if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_w'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_w)) if `touse' label var `varlist' "Conditional stdev for component w" } else if "`sig_t'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_t)) if `touse' label var `varlist' "Conditional stdev for component t" } else if "`rho_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_s)) if `touse' label var `varlist' "Rho s: RTM Survey data" } *********************************************** else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean R1" } else if "`mean_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_t label var `varlist' "Conditional mean R2" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean S1" } else if "`mean_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean S2" } else if "`mean_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n ky_p double `_mean_w', mean_w gen `typlist' `varlist'=`_mean_e'+`_mean_n'+`_mean_w' if `touse' label var `varlist' "Conditional mean S3" } else if "`sig_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma R1" } else if "`sig_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_t label var `varlist' "Conditional sigma R2" } else if "`sig_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma S1" } else if "`sig_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) if `touse' label var `varlist' "Conditional sigma S2" } else if "`sig_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_sig_w', sig_w ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2)) if `touse' label var `varlist' "Conditional sigma S3" } *************************************** else if "`pi_r1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r label var `varlist' "Latent Class pi r type 1" } else if "`pi_r2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r replace `varlist'=1-`varlist' if `touse' label var `varlist' "Latent Class pi r type 2" } else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*(1-`_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 2" } else if "`pi_s3'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*( `_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 3" } *************************************** else if "`pi_1'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_3" } else if "`pi_4'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_4" } else if "`pi_5'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_5" } else if "`pi_6'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_6" } *** posterior. assuming we see only 1 thing. else if "`pip_r1'" != "" { tempvar _mean_r1 _mean_r2 _sig_r1 _sig_r2 _pi_r1 ky_p double `_mean_r1', mean_r1 ky_p double `_sig_r1', sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2', sig_r2 ky_p double `_pi_r1', pi_r1 tempvar _pip_1 _pip_2 gen double `_pip_1' = `_pi_r1'*normalden(`_r',`_mean_r1',`_sig_r1') gen double `_pip_2' =(1-`_pi_r1')*normalden(`_r',`_mean_r2',`_sig_r2') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2') if `touse' label var `varlist' "Posterior pi r type 1" } else if "`pip_r2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pip_r1 replace `varlist'=1-`varlist' label var `varlist' "Posterior pi r type 2" } else if "`pip_s1'`pip_s2'`pip_s3'" != "" { tempvar _mean_s1 _mean_s2 _mean_s3 /// _sig_s1 _sig_s2 _sig_s3 /// _pi_s1 _pi_s2 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 ky_p double `_pi_s1' , pi_s1 ky_p double `_pi_s2' , pi_s2 local _pi_s3 (1-`_pi_s1'-`_pi_s2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_s1' *normalden(`_s',`_mean_s1',`_sig_s1') gen double `_pip_2' = `_pi_s2' *normalden(`_s',`_mean_s2',`_sig_s2') gen double `_pip_3' = `_pi_s3' *normalden(`_s',`_mean_s3',`_sig_s3') if "`pip_s1'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 2" } else if "`pip_s3'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 3" } } ******************************************** else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'`pip_3'`pip_4'`pip_5'`pip_6'" != "" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_sig_e' , sig_e ky_p double `_rho_s' , rho_s tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r1s3 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s3')) local _rho_r2s1 0 local _rho_r2s2 0 local _rho_r2s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip5'=`_pi_5'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip6'=`_pi_6'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s3',`_sig_s3')) syntax newvarname [if] [in] [, * ] if "`pip_2'"!="" { gen `typlist' `varlist'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } if "`pip_3'"!="" { gen `typlist' `varlist'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } if "`pip_4'"!="" { gen `typlist' `varlist'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_4" } if "`pip_5'"!="" { gen `typlist' `varlist'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_5" } if "`pip_6'"!="" { gen `typlist' `varlist'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_6" } } ************************************************************************* else if "`bclass_s'"!="" { tempvar _pip_s1 _pip_s2 ky_p double `_pip_s1' , pip_s1 ky_p double `_pip_s2' , pip_s2 local _pip_s3 (1-`_pip_s1'-`_pip_s2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_s`i''>`_high' replace `_high' =`_pip_s`i'' if `_pip_s`i''>`_high' } label var `varlist' "two-step classification based on S only" } else if "`bclass_r'"!="" { tempvar _pip_r1 ky_p `typlist' `_pip_r1' , pip_r1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1+`_pip_r1'<.5 if `touse' } else if "`bclass'"!="" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_rho_s' , rho_s ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r1s3 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s3')) local _rho_r2s1 0 local _rho_r2s2 0 local _rho_r2s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip5'=`_pi_5'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip6'=`_pi_6'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s3',`_sig_s3')) tempvar _pip_2 _pip_3 _pip_4 _pip_5 _pip_6 gen double `_pip_2'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_3'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_4'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_5'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_6'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' tempvar _high gen double `_high'=. replace `_high'=0 if `_l'==0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1 if `_l'==1 forvalues i=2/6 { replace `varlist'=`i' if `_pip_`i''>`_high' & `_l'==0 replace `_high' =`_pip_`i'' if `_pip_`i''>`_high' & `_l'==0 } } } end *********************** *********************** *********************** program define ky_p_5 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_r pi_w pi_v /// Probabilities of each event rho_s rho_r /// RHO's (Regression to the mean factors) mean_e mean_n mean_w mean_t mean_v /// Conditional Means sig_e sig_n sig_w sig_t sig_v /// conditional Standard deviations mean_r1 mean_r2 mean_r3 mean_s1 mean_s2 mean_s3 /// mean by group sig_r1 sig_r2 sig_r3 sig_s1 sig_s2 sig_s3 /// sigma by group pi_r1 pi_r2 pi_r3 pi_s1 pi_s2 pi_s3 /// Probabilities by category pip_r1 pip_r2 pip_r3 pip_s1 pip_s2 pip_s3 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 pi_5 pi_6 pi_7 pi_8 pi_9 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 pip_5 pip_6 pip_7 pip_8 pip_9 /// Probabilities by group classification bclass_r bclass_s bclass * ] marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: *** Creating this variables first may make this program more efficient. tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r *************************************************************** qui { if "`options'"!="" { ml_p `0' exit } else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_r))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_w))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_v'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_v))) if `touse' label var `varlist' "A priori Pi_v" } **************** else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_w'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_w) label var `varlist' "Conditional mean for component w" } else if "`mean_t'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_t) label var `varlist' "Conditional mean for component t" } else if "`mean_v'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_v) label var `varlist' "Conditional mean for component v" } *** Sigmas else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_n)) if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_w'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_w)) if `touse' label var `varlist' "Conditional stdev for component w" } else if "`sig_t'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_t)) if `touse' label var `varlist' "Conditional stdev for component t" } else if "`sig_v'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_v)) if `touse' label var `varlist' "Conditional stdev for component v" } else if "`rho_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_r)) if `touse' label var `varlist' "Rho r: RTM Survey data" } else if "`rho_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_s)) if `touse' label var `varlist' "Rho s: RTM Survey data" } *********************************************** else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean R1" } else if "`mean_r2'"!="" { ky_p double `_mean_e', mean_e ky_p double `_mean_v', mean_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_mean_e'+`_mean_v' if `touse' label var `varlist' "Conditional mean R2" } else if "`mean_r3'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_t label var `varlist' "Conditional mean R3" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean S1" } else if "`mean_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean S2" } else if "`mean_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n ky_p double `_mean_w', mean_w gen `typlist' `varlist'=`_mean_e'+`_mean_n'+`_mean_w' if `touse' label var `varlist' "Conditional mean S3" } else if "`sig_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma R1" } else if "`sig_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_v', sig_v ky_p double `_rho_r', rho_r gen `typlist' `varlist'=(sqrt((1+`_rho_r')^2*`_sig_e'^2+`_sig_v'^2)) if `touse' label var `varlist' "Conditional sigma R2" } else if "`sig_r3'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_t label var `varlist' "Conditional sigma R3" } else if "`sig_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma S1" } else if "`sig_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) if `touse' label var `varlist' "Conditional sigma S2" } else if "`sig_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_sig_w', sig_w ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2)) if `touse' label var `varlist' "Conditional sigma S3" } *************************************** else if "`pi_r1'" != "" { ky_p double `_pi_r' if `touse' , pi_r ky_p double `_pi_v' if `touse' , pi_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(( `_pi_r')*( `_pi_v')) if `touse' label var `varlist' "Latent Class pi r type 1" } else if "`pi_r2'" != "" { ky_p double `_pi_r' if `touse' , pi_r ky_p double `_pi_v' if `touse' , pi_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(( `_pi_r')*(1-`_pi_v')) if `touse' label var `varlist' "Latent Class pi r type 2" } else if "`pi_r3'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r replace `varlist'=1-`varlist' if `touse' label var `varlist' "Latent Class pi r type 3" } else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*(1-`_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 2" } else if "`pi_s3'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*( `_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 3" } *************************************** else if "`pi_1'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_3" } else if "`pi_4'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_4" } else if "`pi_5'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_5" } else if "`pi_6'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_6" } else if "`pi_7'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_7" } else if "`pi_8'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_8" } else if "`pi_9'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_9" } *** posterior. assuming we see only 1 thing. else if "`pip_r1'`pip_r2'`pip_r3'" != "" { tempvar _mean_r1 _mean_r2 _mean_r3 /// _sig_r1 _sig_r2 _sig_r3 /// _pi_r1 _pi_r2 ky_p double `_mean_r1', mean_r1 ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_pi_r1' , pi_r1 ky_p double `_pi_r2' , pi_r2 local _pi_r3 (1-`_pi_r1'-`_pi_r2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_r1' *normalden(`_r',`_mean_r1',`_sig_r1') gen double `_pip_2' = `_pi_r2' *normalden(`_r',`_mean_r2',`_sig_r2') gen double `_pip_3' = `_pi_r3' *normalden(`_r',`_mean_r3',`_sig_r3') if "`pip_r1'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 1" } else if "`pip_r2'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 2" } else if "`pip_r3'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 3" } } else if "`pip_s1'`pip_s2'`pip_s3'" != "" { tempvar _mean_s1 _mean_s2 _mean_s3 /// _sig_s1 _sig_s2 _sig_s3 /// _pi_s1 _pi_s2 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 ky_p double `_pi_s1' , pi_s1 ky_p double `_pi_s2' , pi_s2 local _pi_s3 (1-`_pi_s1'-`_pi_s2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_s1' *normalden(`_s',`_mean_s1',`_sig_s1') gen double `_pip_2' = `_pi_s2' *normalden(`_s',`_mean_s2',`_sig_s2') gen double `_pip_3' = `_pi_s3' *normalden(`_s',`_mean_s3',`_sig_s3') if "`pip_s1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 2" } else if "`pip_s3'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 3" } } ******************************************** else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'`pip_3'`pip_4'`pip_5'`pip_6'`pip_7'`pip_8'`pip_9'" != "" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 /// _mean_r3 _sig_r3 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_rho_s' , rho_s ky_p double `_rho_r' , rho_r ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 _pi_7 _pi_8 _pi_9 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 ky_p double `_pi_7' , pi_7 ky_p double `_pi_8' , pi_8 ky_p double `_pi_9' , pi_9 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r1s3 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s3')) local _rho_r2s1 (( (1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s1')) local _rho_r2s2 (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s2')) local _rho_r2s3 (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s3')) local _rho_r3s1 0 local _rho_r3s2 0 local _rho_r3s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 _pip7 _pip8 _pip9 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_r3 ((`_r'-`_mean_r3')/`_sig_r3') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(1/(2*_pi*`_sig_r2'*`_sig_s1'*sqrt(1-`_rho_r2s1'^2))*exp(-1/(2*(1-`_rho_r2s1'^2))*(`_std_r2'^2+`_std_s1'^2-2*`_rho_r2s1'*`_std_r2'*`_std_s1'))) gen double `_pip5'=`_pi_5'*(1/(2*_pi*`_sig_r2'*`_sig_s2'*sqrt(1-`_rho_r2s2'^2))*exp(-1/(2*(1-`_rho_r2s2'^2))*(`_std_r2'^2+`_std_s2'^2-2*`_rho_r2s2'*`_std_r2'*`_std_s2'))) gen double `_pip6'=`_pi_6'*(1/(2*_pi*`_sig_r2'*`_sig_s3'*sqrt(1-`_rho_r2s3'^2))*exp(-1/(2*(1-`_rho_r2s3'^2))*(`_std_r2'^2+`_std_s3'^2-2*`_rho_r2s3'*`_std_r2'*`_std_s3'))) gen double `_pip7'=`_pi_7'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip8'=`_pi_8'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip9'=`_pi_9'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s3',`_sig_s3')) syntax newvarname [if] [in] [, * ] if "`pip_2'"!="" { gen `typlist' `varlist'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } if "`pip_3'"!="" { gen `typlist' `varlist'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } if "`pip_4'"!="" { gen `typlist' `varlist'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_4" } if "`pip_5'"!="" { gen `typlist' `varlist'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_5" } if "`pip_6'"!="" { gen `typlist' `varlist'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_6" } if "`pip_7'"!="" { gen `typlist' `varlist'=`_pip7'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_7" } if "`pip_8'"!="" { gen `typlist' `varlist'=`_pip8'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_8" } if "`pip_9'"!="" { gen `typlist' `varlist'=`_pip9'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_9" } } ************************************************************************* else if "`bclass_s'"!="" { tempvar _pip_s1 _pip_s2 ky_p double `_pip_s1' , pip_s1 ky_p double `_pip_s2' , pip_s2 local _pip_s3 (1-`_pip_s1'-`_pip_s2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_s`i''>`_high' replace `_high' =`_pip_s`i'' if `_pip_s`i''>`_high' } label var `varlist' "two-step classification based on S only" } else if "`bclass_r'"!="" { tempvar _pip_r1 _pip_r2 ky_p double `_pip_r1' , pip_r1 ky_p double `_pip_r2' , pip_r2 local _pip_r3 (1-`_pip_r1'-`_pip_r2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_r`i''>`_high' replace `_high' =`_pip_r`i'' if `_pip_r`i''>`_high' } label var `varlist' "two-step classification based on R only" } else if "`bclass'"!="" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 /// _mean_r3 _sig_r3 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_rho_s' , rho_s ky_p double `_rho_r' , rho_r ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 _pi_7 _pi_8 _pi_9 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 ky_p double `_pi_7' , pi_7 ky_p double `_pi_8' , pi_8 ky_p double `_pi_9' , pi_9 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r1s3 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s3')) local _rho_r2s1 (( (1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s1')) local _rho_r2s2 (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s2')) local _rho_r2s3 (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s3')) local _rho_r3s1 0 local _rho_r3s2 0 local _rho_r3s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 _pip7 _pip8 _pip9 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_r3 ((`_r'-`_mean_r3')/`_sig_r3') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(1/(2*_pi*`_sig_r2'*`_sig_s1'*sqrt(1-`_rho_r2s1'^2))*exp(-1/(2*(1-`_rho_r2s1'^2))*(`_std_r2'^2+`_std_s1'^2-2*`_rho_r2s1'*`_std_r2'*`_std_s1'))) gen double `_pip5'=`_pi_5'*(1/(2*_pi*`_sig_r2'*`_sig_s2'*sqrt(1-`_rho_r2s2'^2))*exp(-1/(2*(1-`_rho_r2s2'^2))*(`_std_r2'^2+`_std_s2'^2-2*`_rho_r2s2'*`_std_r2'*`_std_s2'))) gen double `_pip6'=`_pi_6'*(1/(2*_pi*`_sig_r2'*`_sig_s3'*sqrt(1-`_rho_r2s3'^2))*exp(-1/(2*(1-`_rho_r2s3'^2))*(`_std_r2'^2+`_std_s3'^2-2*`_rho_r2s3'*`_std_r2'*`_std_s3'))) gen double `_pip7'=`_pi_7'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip8'=`_pi_8'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip9'=`_pi_9'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s3',`_sig_s3')) tempvar _pip_2 _pip_3 _pip_4 _pip_5 _pip_6 _pip_7 _pip_8 _pip_9 gen double `_pip_2'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_3'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_4'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_5'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_6'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_7'=`_pip7'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_8'=`_pip8'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_9'=`_pip9'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' tempvar _high gen double `_high'=. replace `_high'=0 if `_l'==0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1 if `_l'==1 forvalues i=2/9 { replace `varlist'=`i' if `_pip_`i''>`_high' & `_l'==0 replace `_high' =`_pip_`i'' if `_pip_`i''>`_high' & `_l'==0 } } } end *********************** *********************** *********************** program define ky_p_6 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_r pi_v /// Probabilities of each event rho_s rho_r /// RHO's (Regression to the mean factors) mean_e mean_n mean_t mean_v /// Conditional Means sig_e sig_n sig_t sig_v /// conditional Standard deviations mean_r1 mean_r2 mean_r3 mean_s1 mean_s2 /// mean by group sig_r1 sig_r2 sig_r3 sig_s1 sig_s2 /// sigma by group pi_r1 pi_r2 pi_r3 pi_s1 pi_s2 /// Probabilities by category pip_r1 pip_r2 pip_r3 pip_s1 pip_s2 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 pi_5 pi_6 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 pip_5 pip_6 /// Probabilities by group classification bclass_r bclass_s bclass * ] marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: *** Creating this variables first may make this program more efficient. tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r *************************************************************** qui { if "`options'"!="" { ml_p `0' exit } else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_r))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_v'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_v))) if `touse' label var `varlist' "A priori Pi_v" } **************** else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_t'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_t) label var `varlist' "Conditional mean for component t" } else if "`mean_v'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_v) label var `varlist' "Conditional mean for component v" } *** Sigmas else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_n)) if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_t'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_t)) if `touse' label var `varlist' "Conditional stdev for component t" } else if "`sig_v'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_v)) if `touse' label var `varlist' "Conditional stdev for component v" } else if "`rho_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_r)) if `touse' label var `varlist' "Rho r: RTM Survey data" } else if "`rho_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_s)) if `touse' label var `varlist' "Rho s: RTM Survey data" } *********************************************** else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean R1" } else if "`mean_r2'"!="" { ky_p double `_mean_e', mean_e ky_p double `_mean_v', mean_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_mean_e'+`_mean_v' if `touse' label var `varlist' "Conditional mean R2" } else if "`mean_r3'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_t label var `varlist' "Conditional mean R3" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean S1" } else if "`mean_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean S2" } else if "`sig_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma R1" } else if "`sig_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_v', sig_v ky_p double `_rho_r', rho_r gen `typlist' `varlist'=(sqrt((1+`_rho_r')^2*`_sig_e'^2+`_sig_v'^2)) if `touse' label var `varlist' "Conditional sigma R2" } else if "`sig_r3'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_t label var `varlist' "Conditional sigma R3" } else if "`sig_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma S1" } else if "`sig_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) if `touse' label var `varlist' "Conditional sigma S2" } *************************************** else if "`pi_r1'" != "" { ky_p double `_pi_r' if `touse' , pi_r ky_p double `_pi_v' if `touse' , pi_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(( `_pi_r')*( `_pi_v')) if `touse' label var `varlist' "Latent Class pi r type 1" } else if "`pi_r2'" != "" { ky_p double `_pi_r' if `touse' , pi_r ky_p double `_pi_v' if `touse' , pi_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(( `_pi_r')*(1-`_pi_v')) if `touse' label var `varlist' "Latent Class pi r type 2" } else if "`pi_r3'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r replace `varlist'=1-`varlist' if `touse' label var `varlist' "Latent Class pi r type 3" } else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { ky_p double `_pi_s' if `touse' , pi_s syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(1-`_pi_s') if `touse' label var `varlist' "Latent Class pi s type 2" } *************************************** else if "`pi_1'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_3" } else if "`pi_4'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_4" } else if "`pi_5'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_5" } else if "`pi_6'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_6" } *** posterior. assuming we see only 1 thing. else if "`pip_r1'`pip_r2'`pip_r3'" != "" { tempvar _mean_r1 _mean_r2 _mean_r3 /// _sig_r1 _sig_r2 _sig_r3 /// _pi_r1 _pi_r2 ky_p double `_mean_r1', mean_r1 ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_pi_r1' , pi_r1 ky_p double `_pi_r2' , pi_r2 local _pi_r3 (1-`_pi_r1'-`_pi_r2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_r1' *normalden(`_r',`_mean_r1',`_sig_r1') gen double `_pip_2' = `_pi_r2' *normalden(`_r',`_mean_r2',`_sig_r2') gen double `_pip_3' = `_pi_r3' *normalden(`_r',`_mean_r3',`_sig_r3') if "`pip_r1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 1" } else if "`pip_r2'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 2" } else if "`pip_r3'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 3" } } else if "`pip_s1'`pip_s2'" != "" { tempvar _mean_s1 _mean_s2 _mean_s3 /// _sig_s1 _sig_s2 _sig_s3 /// _pi_s1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_pi_s1' , pi_s1 local _pi_s2 (1-`_pi_s1') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_s1' *normalden(`_s',`_mean_s1',`_sig_s1') gen double `_pip_2' = `_pi_s2' *normalden(`_s',`_mean_s2',`_sig_s2') if "`pip_s1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2') if `touse' label var `varlist' "Posterior pi s type 2" } } ******************************************* else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'`pip_3'`pip_4'`pip_5'`pip_6'" != "" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_r2 _sig_r2 /// _mean_r3 _sig_r3 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_rho_s' , rho_s ky_p double `_rho_r' , rho_r ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 _pi_7 _pi_8 _pi_9 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r2s1 (( (1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s1')) local _rho_r2s2 (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s2')) local _rho_r3s1 0 local _rho_r3s2 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 _pip7 _pip8 _pip9 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_r3 ((`_r'-`_mean_r3')/`_sig_r3') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r2'*`_sig_s1'*sqrt(1-`_rho_r2s1'^2))*exp(-1/(2*(1-`_rho_r2s1'^2))*(`_std_r2'^2+`_std_s1'^2-2*`_rho_r2s1'*`_std_r2'*`_std_s1'))) gen double `_pip4'=`_pi_4'*(1/(2*_pi*`_sig_r2'*`_sig_s2'*sqrt(1-`_rho_r2s2'^2))*exp(-1/(2*(1-`_rho_r2s2'^2))*(`_std_r2'^2+`_std_s2'^2-2*`_rho_r2s2'*`_std_r2'*`_std_s2'))) gen double `_pip5'=`_pi_5'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip6'=`_pi_6'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s2',`_sig_s2')) syntax newvarname [if] [in] [, * ] if "`pip_2'"!="" { gen `typlist' `varlist'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } if "`pip_3'"!="" { gen `typlist' `varlist'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } if "`pip_4'"!="" { gen `typlist' `varlist'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_4" } if "`pip_5'"!="" { gen `typlist' `varlist'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_5" } if "`pip_6'"!="" { gen `typlist' `varlist'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_6" } } ************************************************************************* else if "`bclass_s'"!="" { tempvar _pip_s1 ky_p `typlist' `_pip_s1' , pip_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1+`_pip_s1'<.5 if `touse' label var `varlist' "two-step classification based on S only" } else if "`bclass_r'"!="" { tempvar _pip_r1 _pip_r2 ky_p double `_pip_r1' , pip_r1 ky_p double `_pip_r2' , pip_r2 local _pip_r3 (1-`_pip_r1'-`_pip_r2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_r`i''>`_high' replace `_high' =`_pip_r`i'' if `_pip_r`i''>`_high' } label var `varlist' "two-step classification based on R only" } else if "`bclass'"!="" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_r2 _sig_r2 /// _mean_r3 _sig_r3 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_rho_s' , rho_s ky_p double `_rho_r' , rho_r ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 _pi_7 _pi_8 _pi_9 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r2s1 (( (1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s1')) local _rho_r2s2 (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s2')) local _rho_r3s1 0 local _rho_r3s2 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 _pip7 _pip8 _pip9 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_r3 ((`_r'-`_mean_r3')/`_sig_r3') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r2'*`_sig_s1'*sqrt(1-`_rho_r2s1'^2))*exp(-1/(2*(1-`_rho_r2s1'^2))*(`_std_r2'^2+`_std_s1'^2-2*`_rho_r2s1'*`_std_r2'*`_std_s1'))) gen double `_pip4'=`_pi_4'*(1/(2*_pi*`_sig_r2'*`_sig_s2'*sqrt(1-`_rho_r2s2'^2))*exp(-1/(2*(1-`_rho_r2s2'^2))*(`_std_r2'^2+`_std_s2'^2-2*`_rho_r2s2'*`_std_r2'*`_std_s2'))) gen double `_pip5'=`_pi_5'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip6'=`_pi_6'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s2',`_sig_s2')) tempvar _pip_2 _pip_3 _pip_4 _pip_5 _pip_6 gen double `_pip_2'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_3'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_4'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_5'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_6'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' tempvar _high gen double `_high'=. replace `_high'=0 if `_l'==0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1 if `_l'==1 forvalues i=2/6 { replace `varlist'=`i' if `_pip_`i''>`_high' & `_l'==0 replace `_high' =`_pip_`i'' if `_pip_`i''>`_high' & `_l'==0 } } } end program define ky_p_7 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_r pi_w /// Probabilities of each event rho_s rho_w /// RHO's (Regression to the mean factors) mean_e mean_n mean_w mean_t /// Conditional Means sig_e sig_n sig_w sig_t /// conditional Standard deviations mean_r1 mean_r2 mean_s1 mean_s2 mean_s3 /// mean by group sig_r1 sig_r2 sig_s1 sig_s2 sig_s3 /// sigma by group pi_r1 pi_r2 pi_s1 pi_s2 pi_s3 /// Probabilities by category pip_r1 pip_r2 pip_s1 pip_s2 pip_s3 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 pi_5 pi_6 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 pip_5 pip_6 /// Probabilities by group classification bclass_r bclass_s bclass * ] marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r _rho_w *** Moments based on parameters: qui { if "`options'"!="" { ml_p `0' exit } else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_r))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_w))) if `touse' label var `varlist' "A priori Pi_r" } **************** else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_w'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_w) label var `varlist' "Conditional mean for component w" } else if "`mean_t'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_t) label var `varlist' "Conditional mean for component t" } *** Sigmas else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_n)) if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_w'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_w)) if `touse' label var `varlist' "Conditional stdev for component w" } else if "`sig_t'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_t)) if `touse' label var `varlist' "Conditional stdev for component t" } else if "`rho_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_s)) if `touse' label var `varlist' "Rho s: RTM Survey data" } else if "`rho_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_w)) if `touse' label var `varlist' "Rho w: corr(e_i,w_i)" } *********************************************** else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean R1" } else if "`mean_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_t label var `varlist' "Conditional mean R2" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean S1" } else if "`mean_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean S2" } else if "`mean_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n ky_p double `_mean_w', mean_w gen `typlist' `varlist'=`_mean_e'+`_mean_n'+`_mean_w' if `touse' label var `varlist' "Conditional mean S3" } else if "`sig_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma R1" } else if "`sig_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_t label var `varlist' "Conditional sigma R2" } else if "`sig_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma S1" } else if "`sig_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) if `touse' label var `varlist' "Conditional sigma S2" } else if "`sig_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_sig_w', sig_w ky_p double `_rho_s', rho_s ky_p double `_rho_w', rho_w gen `typlist' `varlist'=sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2+2*(1+`_rho_s')*`_rho_w'*`_sig_w'*`_sig_e') if `touse' label var `varlist' "Conditional sigma S3" } *************************************** else if "`pi_r1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r label var `varlist' "Latent Class pi r type 1" } else if "`pi_r2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r replace `varlist'=1-`varlist' if `touse' label var `varlist' "Latent Class pi r type 2" } else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*(1-`_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 2" } else if "`pi_s3'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*( `_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 3" } *************************************** else if "`pi_1'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_3" } else if "`pi_4'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_4" } else if "`pi_5'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_5" } else if "`pi_6'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_6" } *** posterior. assuming we see only 1 thing. else if "`pip_r1'" != "" { tempvar _mean_r1 _mean_r2 _sig_r1 _sig_r2 _pi_r1 ky_p double `_mean_r1', mean_r1 ky_p double `_sig_r1', sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2', sig_r2 ky_p double `_pi_r1', pi_r1 tempvar _pip_1 _pip_2 gen double `_pip_1' = `_pi_r1'*normalden(`_r',`_mean_r1',`_sig_r1') gen double `_pip_2' =(1-`_pi_r1')*normalden(`_r',`_mean_r2',`_sig_r2') syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2') if `touse' label var `varlist' "Posterior pi r type 1" } else if "`pip_r2'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse', pip_r1 replace `varlist'=1-`varlist' label var `varlist' "Posterior pi r type 2" } else if "`pip_s1'`pip_s2'`pip_s3'" != "" { tempvar _mean_s1 _mean_s2 _mean_s3 /// _sig_s1 _sig_s2 _sig_s3 /// _pi_s1 _pi_s2 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 ky_p double `_pi_s1' , pi_s1 ky_p double `_pi_s2' , pi_s2 local _pi_s3 (1-`_pi_s1'-`_pi_s2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_s1' *normalden(`_s',`_mean_s1',`_sig_s1') gen double `_pip_2' = `_pi_s2' *normalden(`_s',`_mean_s2',`_sig_s2') gen double `_pip_3' = `_pi_s3' *normalden(`_s',`_mean_s3',`_sig_s3') if "`pip_s1'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 2" } else if "`pip_s3'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 3" } } ******************************************** else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'`pip_3'`pip_4'`pip_5'`pip_6'" != "" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_rho_s' , rho_s ky_p double `_rho_w' , rho_w * ky_p double `_sig_w' , sig_w ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r1s3 (((1+`_rho_s')*`_sig_e'+`_rho_w'*`_sig_w')/`_sig_s3') local _rho_r2s1 0 local _rho_r2s2 0 local _rho_r2s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip5'=`_pi_5'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip6'=`_pi_6'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s3',`_sig_s3')) syntax newvarname [if] [in] [, * ] if "`pip_2'"!="" { gen `typlist' `varlist'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } if "`pip_3'"!="" { gen `typlist' `varlist'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } if "`pip_4'"!="" { gen `typlist' `varlist'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_4" } if "`pip_5'"!="" { gen `typlist' `varlist'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_5" } if "`pip_6'"!="" { gen `typlist' `varlist'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_6" } } ************************************************************************* else if "`bclass_s'"!="" { tempvar _pip_s1 _pip_s2 ky_p double `_pip_s1' , pip_s1 ky_p double `_pip_s2' , pip_s2 local _pip_s3 (1-`_pip_s1'-`_pip_s2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_s`i''>`_high' replace `_high' =`_pip_s`i'' if `_pip_s`i''>`_high' } label var `varlist' "two-step classification based on S only" } else if "`bclass_r'"!="" { tempvar _pip_r1 ky_p `typlist' `_pip_r1' , pip_r1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1+`_pip_r1'<.5 if `touse' } else if "`bclass'"!="" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_rho_s' , rho_s ky_p double `_rho_w' , rho_w ky_p double `_sig_w' , sig_w ky_p double `_sig_e' , sig_e tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') local _rho_r1s1 1 local _rho_r1s2 (((1+`_rho_s')*`_sig_e'^2)/(`_sig_r1'*`_sig_s2')) local _rho_r1s3 (((1+`_rho_s')*`_sig_e'+`_rho_w'*`_sig_w')/`_sig_s3') local _rho_r2s1 0 local _rho_r2s2 0 local _rho_r2s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip5'=`_pi_5'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip6'=`_pi_6'*(normalden(`_r',`_mean_r2',`_sig_r2')*normalden(`_s',`_mean_s3',`_sig_s3')) tempvar _pip_2 _pip_3 _pip_4 _pip_5 _pip_6 gen double `_pip_2'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_3'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_4'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_5'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' gen double `_pip_6'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6') if `touse' tempvar _high gen double `_high'=. replace `_high'=0 if `_l'==0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1 if `_l'==1 forvalues i=2/6 { replace `varlist'=`i' if `_pip_`i''>`_high' & `_l'==0 replace `_high' =`_pip_`i'' if `_pip_`i''>`_high' & `_l'==0 } } } end program define ky_p_8 version 15 syntax anything(id="newvarname") [if] [in] , [pi_s pi_r pi_w pi_v /// Probabilities of each event rho_s rho_r rho_w /// RHO's (Regression to the mean factors) mean_e mean_n mean_w mean_t mean_v /// Conditional Means sig_e sig_n sig_w sig_t sig_v /// conditional Standard deviations mean_r1 mean_r2 mean_r3 mean_s1 mean_s2 mean_s3 /// mean by group sig_r1 sig_r2 sig_r3 sig_s1 sig_s2 sig_s3 /// sigma by group pi_r1 pi_r2 pi_r3 pi_s1 pi_s2 pi_s3 /// Probabilities by category pip_r1 pip_r2 pip_r3 pip_s1 pip_s2 pip_s3 /// Probabilities by category. Posterior Prob pi_1 pi_2 pi_3 pi_4 pi_5 pi_6 pi_7 pi_8 pi_9 /// Probabilities by group classification pip_1 pip_2 pip_3 pip_4 pip_5 pip_6 pip_7 pip_8 pip_9 /// Probabilities by group classification bclass_r bclass_s bclass * ] marksample touse, novarlist *** getting dep variables tokenize "`e(depvar)'" local _r `1' local _s `2' local _l `3' *** NOW ALL cases *** defining all elements: *** Creating this variables first may make this program more efficient. tempvar _pi_s _pi_r _pi_w _pi_v tempvar _mean_e _mean_n _mean_t _mean_w _mean_v tempvar _sig_e _sig_n _sig_t _sig_w _sig_v tempvar _rho_s _rho_r _rho_w *************************************************************** qui { if "`options'"!="" { ml_p `0' exit } else if "`pi_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_s))) if `touse' label var `varlist' "A priori Pi_s" } else if "`pi_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_r))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_w))) if `touse' label var `varlist' "A priori Pi_r" } else if "`pi_v'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =(invlogit(xb(lpi_v))) if `touse' label var `varlist' "A priori Pi_v" } **************** else if "`mean_e'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_e) label var `varlist' "Conditional mean for component e" } else if "`mean_n'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_n) label var `varlist' "Conditional mean for component n" } else if "`mean_w'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_w) label var `varlist' "Conditional mean for component w" } else if "`mean_t'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_t) label var `varlist' "Conditional mean for component t" } else if "`mean_v'"!="" { syntax newvarname [if] [in] [, * ] _predict `typlist' `varlist' if `touse', xb eq(mu_v) label var `varlist' "Conditional mean for component v" } *** Sigmas else if "`sig_e'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_e)) if `touse' label var `varlist' "Conditional stdev for component e" } else if "`sig_n'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_n)) if `touse' label var `varlist' "Conditional stdev for component n" } else if "`sig_w'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_w)) if `touse' label var `varlist' "Conditional stdev for component w" } else if "`sig_t'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_t)) if `touse' label var `varlist' "Conditional stdev for component t" } else if "`sig_v'"!="" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist' =exp(xb(ln_sig_v)) if `touse' label var `varlist' "Conditional stdev for component v" } else if "`rho_r'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_r)) if `touse' label var `varlist' "Rho r: RTM Survey data" } else if "`rho_s'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_s)) if `touse' label var `varlist' "Rho s: RTM Survey data" } else if "`rho_w'" != "" { syntax newvarname [if] [in] [, * ] predictnl `typlist' `varlist'=tanh(xb(arho_w)) if `touse' label var `varlist' "Rho w: corr(e_i,w_i)" } *********************************************** else if "`mean_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean R1" } else if "`mean_r2'"!="" { ky_p double `_mean_e', mean_e ky_p double `_mean_v', mean_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_mean_e'+`_mean_v' if `touse' label var `varlist' "Conditional mean R2" } else if "`mean_r3'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_t label var `varlist' "Conditional mean R3" } else if "`mean_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , mean_e label var `varlist' "Conditional mean S1" } else if "`mean_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n gen `typlist' `varlist'=`_mean_e'+`_mean_n' if `touse' label var `varlist' "Conditional mean S2" } else if "`mean_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_mean_e', mean_e ky_p double `_mean_n', mean_n ky_p double `_mean_w', mean_w gen `typlist' `varlist'=`_mean_e'+`_mean_n'+`_mean_w' if `touse' label var `varlist' "Conditional mean S3" } else if "`sig_r1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional mean R1" } else if "`sig_r2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_v', sig_v ky_p double `_rho_r', rho_r gen `typlist' `varlist'=(sqrt((1+`_rho_r')^2*`_sig_e'^2+`_sig_v'^2)) if `touse' label var `varlist' "Conditional sigma R2" } else if "`sig_r3'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_t label var `varlist' "Conditional sigma R3" } else if "`sig_s1'"!="" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , sig_e label var `varlist' "Conditional sigma S1" } else if "`sig_s2'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_rho_s', rho_s gen `typlist' `varlist'=(sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2)) if `touse' label var `varlist' "Conditional sigma S2" } else if "`sig_s3'"!="" { syntax newvarname [if] [in] [, * ] ky_p double `_sig_e', sig_e ky_p double `_sig_n', sig_n ky_p double `_sig_w', sig_w ky_p double `_rho_s', rho_s ky_p double `_rho_w', rho_w gen `typlist' `varlist'=sqrt((1+`_rho_s')^2*`_sig_e'^2+`_sig_n'^2+`_sig_w'^2+2*(1+`_rho_s')*`_rho_w'*`_sig_w'*`_sig_e') if `touse' label var `varlist' "Conditional sigma S3" } *************************************** else if "`pi_r1'" != "" { ky_p double `_pi_r' if `touse' , pi_r ky_p double `_pi_v' if `touse' , pi_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(( `_pi_r')*( `_pi_v')) if `touse' label var `varlist' "Latent Class pi r type 1" } else if "`pi_r2'" != "" { ky_p double `_pi_r' if `touse' , pi_r ky_p double `_pi_v' if `touse' , pi_v syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=(( `_pi_r')*(1-`_pi_v')) if `touse' label var `varlist' "Latent Class pi r type 2" } else if "`pi_r3'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_r replace `varlist'=1-`varlist' if `touse' label var `varlist' "Latent Class pi r type 3" } else if "`pi_s1'" != "" { syntax newvarname [if] [in] [, * ] ky_p `typlist' `varlist' if `touse' , pi_s label var `varlist' "Latent Class pi s type 1" } else if "`pi_s2'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*(1-`_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 2" } else if "`pi_s3'" != "" { ky_p double `_pi_s' if `touse' , pi_s ky_p double `_pi_w' if `touse' , pi_w syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=((1-`_pi_s')*( `_pi_w')) if `touse' label var `varlist' "Latent Class pi s type 3" } *************************************** else if "`pi_1'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_1" } else if "`pi_2'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_2" } else if "`pi_3'" != "" { ky_p double `_pi_r', pi_r1 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_3" } else if "`pi_4'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_4" } else if "`pi_5'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_5" } else if "`pi_6'" != "" { ky_p double `_pi_r', pi_r2 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_6" } else if "`pi_7'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s1 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_7" } else if "`pi_8'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s2 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_8" } else if "`pi_9'" != "" { ky_p double `_pi_r', pi_r3 ky_p double `_pi_s', pi_s3 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pi_r'*`_pi_s' if `touse' label var `varlist' "Latent Class pi_9" } *** posterior. assuming we see only 1 thing. else if "`pip_r1'`pip_r2'`pip_r3'" != "" { tempvar _mean_r1 _mean_r2 _mean_r3 /// _sig_r1 _sig_r2 _sig_r3 /// _pi_r1 _pi_r2 ky_p double `_mean_r1', mean_r1 ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_pi_r1' , pi_r1 ky_p double `_pi_r2' , pi_r2 local _pi_r3 (1-`_pi_r1'-`_pi_r2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_r1' *normalden(`_r',`_mean_r1',`_sig_r1') gen double `_pip_2' = `_pi_r2' *normalden(`_r',`_mean_r2',`_sig_r2') gen double `_pip_3' = `_pi_r3' *normalden(`_r',`_mean_r3',`_sig_r3') if "`pip_r1'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 1" } else if "`pip_r2'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 2" } else if "`pip_r3'"!="" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi r type 3" } } else if "`pip_s1'`pip_s2'`pip_s3'" != "" { tempvar _mean_s1 _mean_s2 _mean_s3 /// _sig_s1 _sig_s2 _sig_s3 /// _pi_s1 _pi_s2 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 ky_p double `_pi_s1' , pi_s1 ky_p double `_pi_s2' , pi_s2 local _pi_s3 (1-`_pi_s1'-`_pi_s2') tempvar _pip_1 _pip_2 _pip_3 gen double `_pip_1' = `_pi_s1' *normalden(`_s',`_mean_s1',`_sig_s1') gen double `_pip_2' = `_pi_s2' *normalden(`_s',`_mean_s2',`_sig_s2') gen double `_pip_3' = `_pi_s3' *normalden(`_s',`_mean_s3',`_sig_s3') if "`pip_s1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_1'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 1" } else if "`pip_s2'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_2'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 2" } else if "`pip_s3'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_pip_3'/(`_pip_1'+`_pip_2'+`_pip_3') if `touse' label var `varlist' "Posterior pi s type 3" } } ******************************************** else if "`pip_1'" != "" { syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=`_l' if `touse' label var `varlist' "Posterior prob pi_1 Constant and equal to 1" } else if "`pip_2'`pip_3'`pip_4'`pip_5'`pip_6'`pip_7'`pip_8'`pip_9'" != "" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 /// _mean_r3 _sig_r3 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_rho_s' , rho_s ky_p double `_rho_r' , rho_r ky_p double `_rho_w' , rho_w ky_p double `_sig_e' , sig_e ky_p double `_sig_w' , sig_w tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 _pi_7 _pi_8 _pi_9 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 ky_p double `_pi_7' , pi_7 ky_p double `_pi_8' , pi_8 ky_p double `_pi_9' , pi_9 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') tempvar _rho_r1s2 _rho_r1s3 _rho_r2s1 _rho_r2s2 _rho_r2s3 local _rho_r1s1 1 gen double `_rho_r1s2' = (((1+`_rho_s')* `_sig_e'^2)/(`_sig_r1'*`_sig_s2')) gen double `_rho_r1s3' = (((1+`_rho_s')* `_sig_e'^2+`_rho_w'*`_sig_e'*`_sig_w')/(`_sig_r1'*`_sig_s3')) gen double `_rho_r2s1' = (( (1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s1')) gen double `_rho_r2s2' = (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s2')) gen double `_rho_r2s3' = (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2+`_rho_w'*(1+`_rho_r')*`_sig_e'*`_sig_w')/(`_sig_r2'*`_sig_s3')) local _rho_r3s1 0 local _rho_r3s2 0 local _rho_r3s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 _pip7 _pip8 _pip9 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_r3 ((`_r'-`_mean_r3')/`_sig_r3') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(1/(2*_pi*`_sig_r2'*`_sig_s1'*sqrt(1-`_rho_r2s1'^2))*exp(-1/(2*(1-`_rho_r2s1'^2))*(`_std_r2'^2+`_std_s1'^2-2*`_rho_r2s1'*`_std_r2'*`_std_s1'))) gen double `_pip5'=`_pi_5'*(1/(2*_pi*`_sig_r2'*`_sig_s2'*sqrt(1-`_rho_r2s2'^2))*exp(-1/(2*(1-`_rho_r2s2'^2))*(`_std_r2'^2+`_std_s2'^2-2*`_rho_r2s2'*`_std_r2'*`_std_s2'))) gen double `_pip6'=`_pi_6'*(1/(2*_pi*`_sig_r2'*`_sig_s3'*sqrt(1-`_rho_r2s3'^2))*exp(-1/(2*(1-`_rho_r2s3'^2))*(`_std_r2'^2+`_std_s3'^2-2*`_rho_r2s3'*`_std_r2'*`_std_s3'))) gen double `_pip7'=`_pi_7'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip8'=`_pi_8'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip9'=`_pi_9'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s3',`_sig_s3')) syntax newvarname [if] [in] [, * ] if "`pip_2'"!="" { gen `typlist' `varlist'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_2" } if "`pip_3'"!="" { gen `typlist' `varlist'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_3" } if "`pip_4'"!="" { gen `typlist' `varlist'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_4" } if "`pip_5'"!="" { gen `typlist' `varlist'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_5" } if "`pip_6'"!="" { gen `typlist' `varlist'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_6" } if "`pip_7'"!="" { gen `typlist' `varlist'=`_pip7'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_7" } if "`pip_8'"!="" { gen `typlist' `varlist'=`_pip8'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_8" } if "`pip_9'"!="" { gen `typlist' `varlist'=`_pip9'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' replace `varlist'=0 if `touse' & `_l'==1 label var `varlist' "Posterior prob pi_9" } } ************************************************************************* else if "`bclass_s'"!="" { tempvar _pip_s1 _pip_s2 ky_p double `_pip_s1' , pip_s1 ky_p double `_pip_s2' , pip_s2 local _pip_s3 (1-`_pip_s1'-`_pip_s2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_s`i''>`_high' replace `_high' =`_pip_s`i'' if `_pip_s`i''>`_high' } label var `varlist' "two-step classification based on S only" } else if "`bclass_r'"!="" { tempvar _pip_r1 _pip_r2 ky_p double `_pip_r1' , pip_r1 ky_p double `_pip_r2' , pip_r2 local _pip_r3 (1-`_pip_r1'-`_pip_r2') tempvar _high gen double `_high'=0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=0 forvalues i=1/3 { replace `varlist'=`i' if `_pip_r`i''>`_high' replace `_high' =`_pip_r`i'' if `_pip_r`i''>`_high' } label var `varlist' "two-step classification based on R only" } else if "`bclass'"!="" { tempvar _mean_s1 _sig_s1 /// _mean_s2 _sig_s2 /// _mean_s3 _sig_s3 /// _mean_r2 _sig_r2 /// _mean_r3 _sig_r3 *_mean_r1 sig_r1 ky_p double `_mean_s1', mean_s1 ky_p double `_sig_s1' , sig_s1 ky_p double `_mean_s2', mean_s2 ky_p double `_sig_s2' , sig_s2 ky_p double `_mean_s3', mean_s3 ky_p double `_sig_s3' , sig_s3 *ky_p double `_mean_r1', mean_r1 *ky_p double `_sig_r1' , sig_r1 ky_p double `_mean_r2', mean_r2 ky_p double `_sig_r2' , sig_r2 ky_p double `_mean_r3', mean_r3 ky_p double `_sig_r3' , sig_r3 ky_p double `_rho_s' , rho_s ky_p double `_rho_r' , rho_r ky_p double `_rho_w' , rho_w ky_p double `_sig_e' , sig_e ky_p double `_sig_w' , sig_w tempvar _pi_2 _pi_3 _pi_4 _pi_5 _pi_6 _pi_7 _pi_8 _pi_9 ky_p double `_pi_2' , pi_2 ky_p double `_pi_3' , pi_3 ky_p double `_pi_4' , pi_4 ky_p double `_pi_5' , pi_5 ky_p double `_pi_6' , pi_6 ky_p double `_pi_7' , pi_7 ky_p double `_pi_8' , pi_8 ky_p double `_pi_9' , pi_9 local _mean_r1 (`_mean_s1') local _sig_r1 (`_sig_s1') tempvar _rho_r1s2 _rho_r1s3 _rho_r2s1 _rho_r2s2 _rho_r2s3 local _rho_r1s1 1 gen double `_rho_r1s2' = (((1+`_rho_s')* `_sig_e'^2)/(`_sig_r1'*`_sig_s2')) gen double `_rho_r1s3' = (((1+`_rho_s')* `_sig_e'^2+`_rho_w'*`_sig_e'*`_sig_w')/(`_sig_r1'*`_sig_s3')) gen double `_rho_r2s1' = (( (1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s1')) gen double `_rho_r2s2' = (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2)/(`_sig_r2'*`_sig_s2')) gen double `_rho_r2s3' = (((1+`_rho_s')*(1+`_rho_r')*`_sig_e'^2+`_rho_w'*(1+`_rho_r')*`_sig_e'*`_sig_w')/(`_sig_r2'*`_sig_s3')) local _rho_r3s1 0 local _rho_r3s2 0 local _rho_r3s3 0 tempvar _pip2 _pip3 _pip4 _pip5 _pip6 _pip7 _pip8 _pip9 local _std_r1 ((`_r'-`_mean_r1')/`_sig_r1') local _std_r2 ((`_r'-`_mean_r2')/`_sig_r2') local _std_r3 ((`_r'-`_mean_r3')/`_sig_r3') local _std_s1 ((`_s'-`_mean_s1')/`_sig_s1') local _std_s2 ((`_s'-`_mean_s2')/`_sig_s2') local _std_s3 ((`_s'-`_mean_s3')/`_sig_s3') gen double `_pip2'=`_pi_2'*(1/(2*_pi*`_sig_r1'*`_sig_s2'*sqrt(1-`_rho_r1s2'^2))*exp(-1/(2*(1-`_rho_r1s2'^2))*(`_std_r1'^2+`_std_s2'^2-2*`_rho_r1s2'*`_std_r1'*`_std_s2'))) gen double `_pip3'=`_pi_3'*(1/(2*_pi*`_sig_r1'*`_sig_s3'*sqrt(1-`_rho_r1s3'^2))*exp(-1/(2*(1-`_rho_r1s3'^2))*(`_std_r1'^2+`_std_s3'^2-2*`_rho_r1s3'*`_std_r1'*`_std_s3'))) gen double `_pip4'=`_pi_4'*(1/(2*_pi*`_sig_r2'*`_sig_s1'*sqrt(1-`_rho_r2s1'^2))*exp(-1/(2*(1-`_rho_r2s1'^2))*(`_std_r2'^2+`_std_s1'^2-2*`_rho_r2s1'*`_std_r2'*`_std_s1'))) gen double `_pip5'=`_pi_5'*(1/(2*_pi*`_sig_r2'*`_sig_s2'*sqrt(1-`_rho_r2s2'^2))*exp(-1/(2*(1-`_rho_r2s2'^2))*(`_std_r2'^2+`_std_s2'^2-2*`_rho_r2s2'*`_std_r2'*`_std_s2'))) gen double `_pip6'=`_pi_6'*(1/(2*_pi*`_sig_r2'*`_sig_s3'*sqrt(1-`_rho_r2s3'^2))*exp(-1/(2*(1-`_rho_r2s3'^2))*(`_std_r2'^2+`_std_s3'^2-2*`_rho_r2s3'*`_std_r2'*`_std_s3'))) gen double `_pip7'=`_pi_7'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s1',`_sig_s1')) gen double `_pip8'=`_pi_8'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s2',`_sig_s2')) gen double `_pip9'=`_pi_9'*(normalden(`_r',`_mean_r3',`_sig_r3')*normalden(`_s',`_mean_s3',`_sig_s3')) tempvar _pip_2 _pip_3 _pip_4 _pip_5 _pip_6 _pip_7 _pip_8 _pip_9 gen double `_pip_2'=`_pip2'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_3'=`_pip3'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_4'=`_pip4'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_5'=`_pip5'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_6'=`_pip6'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_7'=`_pip7'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_8'=`_pip8'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' gen double `_pip_9'=`_pip9'/(`_pip2'+`_pip3'+`_pip4'+`_pip5'+`_pip6'+`_pip7'+`_pip8'+`_pip9') if `touse' tempvar _high gen double `_high'=. replace `_high'=0 if `_l'==0 syntax newvarname [if] [in] [, * ] gen `typlist' `varlist'=1 if `_l'==1 forvalues i=2/9 { replace `varlist'=`i' if `_pip_`i''>`_high' & `_l'==0 replace `_high' =`_pip_`i'' if `_pip_`i''>`_high' & `_l'==0 } } } end