--------------------------------------------------------------------------------------------- name: log: U:\projects\misc\stata_lf_programs\tobit\two_tobits.txt log type: text opened on: 31 Oct 2014, 15:40:29 . . * Simulate mixture of two Tobits . set seed 12345 . set obs 10000 obs was 0, now 10000 . capture drop g1 g2 y1 y2 mix x . gen g1 = rnormal(0, 10) . gen g2 = rnormal(0, 20) . gen x = rnormal(100, 20) . gen y1 = -300 + 5*x + g1 . replace y1 = 0 if y1 < 0 (252 real changes made) . gen y2 = 10 + 10*x + g2 . replace y2 =0 if y2 <0 (0 real changes made) . gen mix = y1 . replace mix = y2 if uniform() > 0.7 (2979 real changes made) . hist mix, scale(0.7) kdensity (bin=40, start=0, width=44.686481) . . * Set up likelihood estimator . capture program drop twotobits_lf . program twotobits_lf 1. version 12 2. args todo b lnf 3. tempvar xb1 xb2 lnsigma1 lnsigma2 sigma1 sigma2 lj lpr1 lpr2 pr1 pr2 4. . mleval `xb1'=`b', eq(1) 5. mleval `xb2'=`b', eq(2) 6. mleval `lpr1'=`b', eq(3) 7. mleval `lnsigma1'=`b', scalar eq(4) 8. mleval `lnsigma2'=`b', scalar eq(5) 9. . quietly { 10. gen double `sigma1' = exp(`lnsigma1') 11. gen double `sigma2' = exp(`lnsigma2') 12. qui gen double `pr1' = exp(`lpr1') / (1+exp(`lpr1')) 13. qui gen double `pr2' = 1 - `pr1' 14. . gen double `lj'=`pr1'*normal(-`xb1'/`sigma1')+`pr2'*normal(-`xb2'/`sigma2') /* > */ if $ML_y1 ==0 15. replace `lj'=`pr1'*normalden($ML_y1,`xb1',`sigma1')+`pr2'* normalden($ML_y1,` > xb2',`sigma2') /* > */ if $ML_y1 > 0 16. mlsum `lnf' = ln(`lj') 17. } 18. . if (`todo'==0 | `lnf' >=.) exit 19. end . . * Estimate . ml model d0 twotobits_lf (component1: mix = x) (component2: mix = x) /imlogitpi1 /lnsigma1 > /lnsigma2 . ml search initial: log likelihood = - (could not be evaluated) feasible: log likelihood = -107383 rescale: log likelihood = -107383 rescale eq: log likelihood = -73526.422 . ml maximize, diff initial: log likelihood = -73526.422 rescale: log likelihood = -73526.422 rescale eq: log likelihood = -72786.141 Iteration 0: log likelihood = -72786.141 (not concave) Iteration 1: log likelihood = -69800.188 (not concave) Iteration 2: log likelihood = -59368.989 (not concave) Iteration 3: log likelihood = -56796.94 (not concave) Iteration 4: log likelihood = -56426.548 (not concave) Iteration 5: log likelihood = -53821.508 (not concave) Iteration 6: log likelihood = -48553.615 (not concave) Iteration 7: log likelihood = -46558.118 (not concave) Iteration 8: log likelihood = -45843.681 Iteration 9: log likelihood = -45105.228 Iteration 10: log likelihood = -44797.086 Iteration 11: log likelihood = -44795.711 Iteration 12: log likelihood = -44795.71 Number of obs = 10000 Wald chi2(1) = 608540.07 Log likelihood = -44795.71 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- component1 | x | 5.005779 .0064169 780.09 0.000 4.993202 5.018356 _cons | -300.7093 .6568419 -457.81 0.000 -301.9967 -299.4219 -------------+---------------------------------------------------------------- component2 | x | 9.991051 .017897 558.25 0.000 9.955973 10.02613 _cons | 10.91177 1.822194 5.99 0.000 7.340334 14.4832 -------------+---------------------------------------------------------------- imlogitpi1 | _cons | .8573181 .0218658 39.21 0.000 .8144619 .9001743 -------------+---------------------------------------------------------------- lnsigma1 | _cons | 2.310929 .0085393 270.62 0.000 2.294192 2.327666 -------------+---------------------------------------------------------------- lnsigma2 | _cons | 2.993702 .0129554 231.08 0.000 2.96831 3.019094 ------------------------------------------------------------------------------ . . * check mixture probability (should be .7) . di exp(.8573181) / (1 +exp(.8573181)) .70210002 . . log close name: log: U:\projects\misc\stata_lf_programs\tobit\two_tobits.txt log type: text closed on: 31 Oct 2014, 15:40:42 ---------------------------------------------------------------------------------------------