------------------------------------------------------------------------------------------------------------------------------- name: log: c:\acdbookrevision\stata_final_programs_2013\racd06p3.txt log type: text opened on: 26 Jan 2013, 09:50:55 . . ********** OVERVIEW OF racd06p3.do ********** . . * STATA Program . * copyright C 2013 by A. Colin Cameron and Pravin K. Trivedi . * used for "Regression Analyis of Count Data" SECOND EDITION . * by A. Colin Cameron and Pravin K. Trivedi (2013) . * Cambridge University Press . . * Chapter 6.5 only . * 6.5 COMPLETED FERTILITY . . * To run you need files . * racd06data3fertilityswiss.dta . * racd06data4fertilitybritish.dta . * and user-written Stata addons . * fmm and hnblogit . * in your directory . . ********** SETUP ********** . . set more off . version 12 . clear all . set mem 10m set memory ignored. Memory no longer needs to be set in modern Statas; memory adjustments are performed on the fly automatically. . * set linesize 82 . set scheme s1mono /* Graphics scheme */ . . ********** DATA DESCRIPTION . . * Two datasets ... . * (1) Swiss Household Panel W1 1999 N = 1878 . * (2) British Household Panel updated to Wave 18 N = 6782 . * See ?? for more detailed discussion . * Also see racd06makedata3fertility.do for further details . . ********** 6.5 COMPLETED FERTILITY . . ********** SWISS DATA SUMMARY . . use racd06data3fertilityswiss.dta (Completed fertility data --- Swiss Household Panel W1 [1999]) . describe Contains data from racd06data3fertilityswiss.dta obs: 1,878 Completed fertility data --- Swiss Household Panel W1 [1999] vars: 31 12 Dec 2011 17:06 size: 84,510 ------------------------------------------------------------------------------------------------------------------------------- storage display value variable name type format label variable label ------------------------------------------------------------------------------------------------------------------------------- idhous99 long %12.0g idhous99 identification number of household idpers_new byte %9.0g idpers long %12.0g idpers identification number of person sex byte %8.0g sex sex birthm byte %8.0g birthm date of birth: month birthy int %8.0g birthy date of birth: year age byte %8.0g age99 age in year of interview marstat byte %21.0g marstat civil status in year of interview intlang byte %13.0g intlang interview language religion byte %8.0g p99r01 religion disrel_fr byte %13.0g disrel_sp discussion about religion with friends: frequency health byte %16.0g health health status work byte %18.0g work working status income double %27.0g income yearly work income, gross cfriends byte %16.0g cfriends contact with close friends: times per month reading byte %22.0g reading leisure: reading: frequency inchild byte %9.0g children byte %9.0g education byte %31.0g edu highest level of education achieved fath15 byte %13.0g fath15 age15: living with father moth15 byte %13.0g moth15 age15: living with mother stepp15 byte %13.0g stepp15 age 15: living with mother/father's new partner sib15 byte %13.0g sib15 age15: living with brothers and sisters inear15 byte %39.0g inear15 age15: main income earner ocfath15 byte %40.0g ocfath15 isco classification: fathers job: 1-digit-position edufath byte %31.0g edufath father's highest level of education achieved mthwk15 byte %13.0g mthwk15 age15: employment: mother ocmth15 byte %40.0g ocmth15 isco classification: mothers job: 1-digit-position edumth byte %31.0g edumth mother's highest level of education achieved fathpol byte %22.0g fathpol political position: left, right: father mothpol byte %22.0g mothpol political position: left, right: mother ------------------------------------------------------------------------------------------------------------------------------- Sorted by: . summarize Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- idhous99 | 1878 72703.23 42171.01 251 146601 idpers_new | 1878 1 0 1 1 idpers | 1878 7270324 4217101 25101 1.47e+07 sex | 1878 2 0 2 2 birthm | 1878 6.381257 3.491143 -3 12 -------------+-------------------------------------------------------- birthy | 1878 1940.257 10.35302 1908 1954 age | 1878 58.74334 10.35302 45 91 marstat | 1878 2.632588 1.244377 1 5 intlang | 1878 1.775825 .5236152 1 3 religion | 1878 2.310969 2.168746 1 8 -------------+-------------------------------------------------------- disrel_fr | 1858 3.085576 2.719462 0 10 health | 1878 2.048988 .7791547 1 5 work | 1878 2.037274 .9995713 1 3 income | 1878 18716.7 42596.07 -8 1170000 cfriends | 1574 6.101652 6.695623 1 60 -------------+-------------------------------------------------------- reading | 1874 1.238527 .700232 1 5 inchild | 1878 .5351438 .9293591 0 6 children | 1878 1.937167 1.452605 0 11 education | 1878 3.707135 2.211864 1 10 fath15 | 1878 1.124068 .3297471 1 2 -------------+-------------------------------------------------------- moth15 | 1878 1.078275 .268675 1 2 stepp15 | 1878 -3 0 -3 -3 sib15 | 1878 1.13951 .3465703 1 2 inear15 | 1873 1.565937 1.051075 1 5 ocfath15 | 1805 5.634349 2.369507 1 9 -------------+-------------------------------------------------------- edufath | 1727 3.451071 3.203382 0 15 mthwk15 | 1854 1.604099 .4891752 1 2 ocmth15 | 1223 5.693377 2.140334 1 9 edumth | 1776 2.328829 2.012911 1 15 fathpol | 1878 4.678914 2.782603 1 10 -------------+-------------------------------------------------------- mothpol | 1878 4.2082 2.412849 1 10 . . *** TABLE 6.15 - FREQUENCIES AND PREDICTED PROBABILITIES (SWISS) . . summarize children, detail children ------------------------------------------------------------- Percentiles Smallest 1% 0 0 5% 0 0 10% 0 0 Obs 1878 25% 1 0 Sum of Wgt. 1878 50% 2 Mean 1.937167 Largest Std. Dev. 1.452605 75% 3 8 90% 4 9 Variance 2.110062 95% 4 10 Skewness .9637311 99% 7 11 Kurtosis 5.743735 . tabulate children children | Freq. Percent Cum. ------------+----------------------------------- 0 | 379 20.18 20.18 1 | 262 13.95 34.13 2 | 684 36.42 70.55 3 | 353 18.80 89.35 4 | 128 6.82 96.17 5 | 35 1.86 98.03 6 | 16 0.85 98.88 7 | 8 0.43 99.31 8 | 10 0.53 99.84 9 | 1 0.05 99.89 10 | 1 0.05 99.95 11 | 1 0.05 100.00 ------------+----------------------------------- Total | 1,878 100.00 . nbreg children Fitting Poisson model: Iteration 0: log likelihood = -3238.338 Iteration 1: log likelihood = -3238.338 Fitting constant-only model: Iteration 0: log likelihood = -3537.6471 Iteration 1: log likelihood = -3235.3125 Iteration 2: log likelihood = -3235.0786 Iteration 3: log likelihood = -3235.0752 Iteration 4: log likelihood = -3235.0752 Fitting full model: Iteration 0: log likelihood = -3235.0752 Iteration 1: log likelihood = -3235.0752 Negative binomial regression Number of obs = 1878 LR chi2(0) = 0.00 Dispersion = mean Prob > chi2 = . Log likelihood = -3235.0752 Pseudo R2 = 0.0000 ------------------------------------------------------------------------------ children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- _cons | .6612267 .0172682 38.29 0.000 .6273817 .6950717 -------------+---------------------------------------------------------------- /lnalpha | -3.128509 .4221316 -3.955871 -2.301146 -------------+---------------------------------------------------------------- alpha | .043783 .0184822 .019142 .100144 ------------------------------------------------------------------------------ Likelihood-ratio test of alpha=0: chibar2(01) = 6.53 Prob>=chibar2 = 0.005 . forvalues i = 0/12 { 2. predict nbfit`i', pr(`i') 3. } . sum nbfit* Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- nbfit0 | 1878 .1557684 0 .1557684 .1557684 nbfit1 | 1878 .2781575 0 .2781575 .2781575 nbfit2 | 1878 .2592283 0 .2592283 .2592283 nbfit3 | 1878 .167814 0 .167814 .167814 nbfit4 | 1878 .0847571 0 .0847571 .0847571 -------------+-------------------------------------------------------- nbfit5 | 1878 .0355717 0 .0355717 .0355717 nbfit6 | 1878 .0129044 0 .0129044 .0129044 nbfit7 | 1878 .0041567 0 .0041567 .0041567 nbfit8 | 1878 .0012122 0 .0012122 .0012122 nbfit9 | 1878 .0003248 0 .0003248 .0003248 -------------+-------------------------------------------------------- nbfit10 | 1878 .0000808 0 .0000808 .0000808 nbfit11 | 1878 .0000189 0 .0000189 .0000189 nbfit12 | 1878 4.16e-06 0 4.16e-06 4.16e-06 . . *** BRITISH DATA SUMMARY . . use racd06data4fertilitybritish.dta, clear (Completed fertility Data --- BHPS, updated up to to wave W18) . . describe Contains data from racd06data4fertilitybritish.dta obs: 6,782 Completed fertility Data --- BHPS, updated up to to wave W18 vars: 29 12 Dec 2011 17:06 size: 400,138 ------------------------------------------------------------------------------------------------------------------------------- storage display value variable name type format label variable label ------------------------------------------------------------------------------------------------------------------------------- pid long %12.0g cross-wave person identifier sex byte %8.0g sex sex dobm byte %8.0g dobm month of birth doby int %8.0g doby year of birth memorig byte %8.0g memorig sample origin plbornc byte %8.0g plbornc country of birth race byte %8.0g race ethnic group membership paju byte %8.0g paju father not working when resp. aged 14 maju byte %8.0g maju mother not working when resp. aged 14 lprnt byte %8.0g lprnt natural parent of children lnprnt byte %8.0g lnprnt no. of children resp. natural parent to ch1bm byte %8.0g ch1bm month first child born ch1by int %8.0g ch1by year first child born scend byte %8.0g scend school leaving age bwtgm int %8.0g bwtgm birth weight: grammes pagold byte %8.0g pagold goldthorpe social class: father's job pargsc byte %8.0g pargsc rg social class: father's job magold byte %8.0g magold goldthorpe social class: mother's job margsc byte %8.0g margsc rg social class: mother's job j1gold byte %8.0g j1gold goldthorpe social class: first job j1rgsc byte %8.0g j1rgsc rg social class: first job fsource float %14.0g fsource source of fertility data children float %9.0g completed fertility ssex2 float %9.0g same-sex siblings at partity 2 ssex3 float %9.0g same-sex siblings at parity 3 twin1 float %9.0g n-plets at parity one (number of) twin2 float %9.0g n-plets at parity two (number of) twin3 float %9.0g n-plets at parity three (number of) age float %9.0g ------------------------------------------------------------------------------------------------------------------------------- Sorted by: pid . summarize Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- pid | 6782 6.22e+07 4.65e+07 1.00e+07 1.89e+08 sex | 6782 2 0 2 2 dobm | 6782 6.41271 3.454293 -9 12 doby | 6782 1950.625 8.777781 1934 1964 memorig | 6782 3.114126 2.478388 1 7 -------------+-------------------------------------------------------- plbornc | 6782 -5.160572 12.79145 -9 92 race | 6782 .432026 2.544695 -9 9 paju | 6782 -7.097611 2.623091 -9 2 maju | 6782 -4.092156 4.024285 -9 2 lprnt | 6782 -4.972427 4.462329 -9 2 -------------+-------------------------------------------------------- lnprnt | 6782 -4.969773 4.984932 -9 18 ch1bm | 6782 -3.906665 6.921281 -9 12 ch1by | 6782 574.9611 903.7228 -9 2001 scend | 6782 14.63801 5.53207 -9 22 bwtgm | 6782 -8 0 -8 -8 -------------+-------------------------------------------------------- pagold | 6782 1.897818 7.572607 -9 11 pargsc | 6782 1.222943 4.970773 -9 7 magold | 6782 -2.715128 7.297597 -9 11 margsc | 6782 -3.165585 5.97136 -9 6 j1gold | 6782 -3.678119 6.500148 -9 11 -------------+-------------------------------------------------------- j1rgsc | 6782 -2.640666 5.827781 -9 6 fsource | 6782 31.33117 65.53704 0 200 children | 6782 1.856532 1.508928 0 11 ssex2 | 1726 1 0 1 1 ssex3 | 806 1 0 1 1 -------------+-------------------------------------------------------- twin1 | 39 2.025641 .1601282 2 3 twin2 | 50 2 0 2 2 twin3 | 23 2 0 2 2 age | 6782 58.37541 8.777781 45 75 . . **** TABLE 6.16 - FREQUENCIES AND PREDICTED PROBABILITIES (BRITISH) . . summarize children, detail completed fertility ------------------------------------------------------------- Percentiles Smallest 1% 0 0 5% 0 0 10% 0 0 Obs 6782 25% 0 0 Sum of Wgt. 6782 50% 2 Mean 1.856532 Largest Std. Dev. 1.508928 75% 3 10 90% 4 10 Variance 2.276863 95% 4 10 Skewness .725737 99% 6 11 Kurtosis 4.13197 . tabulate children completed | fertility | Freq. Percent Cum. ------------+----------------------------------- 0 | 1,764 26.01 26.01 1 | 805 11.87 37.88 2 | 2,160 31.85 69.73 3 | 1,262 18.61 88.34 4 | 487 7.18 95.52 5 | 181 2.67 98.19 6 | 69 1.02 99.20 7 | 33 0.49 99.69 8 | 12 0.18 99.87 9 | 5 0.07 99.94 10 | 3 0.04 99.99 11 | 1 0.01 100.00 ------------+----------------------------------- Total | 6,782 100.00 . nbreg children Fitting Poisson model: Iteration 0: log likelihood = -11962.845 Iteration 1: log likelihood = -11962.845 Fitting constant-only model: Iteration 0: log likelihood = -12543.881 Iteration 1: log likelihood = -11878.561 Iteration 2: log likelihood = -11878.441 Iteration 3: log likelihood = -11878.441 Fitting full model: Iteration 0: log likelihood = -11878.441 Iteration 1: log likelihood = -11878.441 Negative binomial regression Number of obs = 6782 LR chi2(0) = 0.00 Dispersion = mean Prob > chi2 = . Log likelihood = -11878.441 Pseudo R2 = 0.0000 ------------------------------------------------------------------------------ children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- _cons | .6187102 .0100431 61.61 0.000 .5990262 .6383943 -------------+---------------------------------------------------------------- /lnalpha | -1.928155 .0926511 -2.109748 -1.746562 -------------+---------------------------------------------------------------- alpha | .1454162 .013473 .1212685 .1743724 ------------------------------------------------------------------------------ Likelihood-ratio test of alpha=0: chibar2(01) = 168.81 Prob>=chibar2 = 0.000 . forvalues i = 0/12 { 2. predict nbfit`i', pr(`i') 3. } . sum nbfit* Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- nbfit0 | 6782 .1933002 0 .1933002 .1933002 nbfit1 | 6782 .2825799 0 .2825799 .2825799 nbfit2 | 6782 .236583 0 .236583 .236583 nbfit3 | 6782 .1488131 0 .1488131 .1488131 nbfit4 | 6782 .0781124 0 .0781124 .0781124 -------------+-------------------------------------------------------- nbfit5 | 6782 .0361221 0 .0361221 .0361221 nbfit6 | 6782 .0152 0 .0152 .0152 nbfit7 | 6782 .005944 0 .005944 .005944 nbfit8 | 6782 .0021918 0 .0021918 .0021918 nbfit9 | 6782 .0007702 0 .0007702 .0007702 -------------+-------------------------------------------------------- nbfit10 | 6782 .0002599 0 .0002599 .0002599 nbfit11 | 6782 .0000848 0 .0000848 .0000848 nbfit12 | 6782 .0000268 0 .0000268 .0000268 . . *** Histogram of data for the two data sets . . use racd06data3fertilityswiss.dta, clear (Completed fertility data --- Swiss Household Panel W1 [1999]) . summarize children Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- children | 1878 1.937167 1.452605 0 11 . di "Mean: " r(mean) " Variance: " r(Var) Mean: 1.9371672 Variance: 2.1100617 . label variable children "Number of children (Swiss)" . histogram children, discrete frequency barwidth(0.8) saving(racd06graph1, replace) xlabel(#6) (start=0, width=1) (file racd06graph1.gph saved) . . use racd06data4fertilitybritish.dta, clear (Completed fertility Data --- BHPS, updated up to to wave W18) . summarize children Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- children | 6782 1.856532 1.508928 0 11 . di "Mean: " r(mean) " Variance: " r(Var) Mean: 1.856532 Variance: 2.2768627 . label variable children "Number of children (British)" . histogram children, discrete frequency barwidth(0.8) saving(racd06graph2, replace) xlabel(#6) (start=0, width=1) (file racd06graph2.gph saved) . . graph combine racd06graph1.gph racd06graph2.gph, iscale(0.7) ysize(3) xsize(6) xcommon . . ********** VARIOUS INTERCEPT-ONLY COUNT MODELS FOR SWISS FERTILITY DATA . . use racd06data3fertilityswiss.dta, clear (Completed fertility data --- Swiss Household Panel W1 [1999]) . . * Poisson model . poisson children, vce(robust) Iteration 0: log pseudolikelihood = -3238.338 Iteration 1: log pseudolikelihood = -3238.338 Poisson regression Number of obs = 1878 Wald chi2(0) = . Prob > chi2 = . Log pseudolikelihood = -3238.338 Pseudo R2 = 0.0000 ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- _cons | .6612267 .0173034 38.21 0.000 .6273126 .6951408 ------------------------------------------------------------------------------ . estimates store POISSON . . * Negative binomial model . nbreg children, vce(robust) Fitting Poisson model: Iteration 0: log pseudolikelihood = -3238.338 Iteration 1: log pseudolikelihood = -3238.338 Fitting constant-only model: Iteration 0: log pseudolikelihood = -3537.6471 Iteration 1: log pseudolikelihood = -3235.3125 Iteration 2: log pseudolikelihood = -3235.0786 Iteration 3: log pseudolikelihood = -3235.0752 Iteration 4: log pseudolikelihood = -3235.0752 Fitting full model: Iteration 0: log pseudolikelihood = -3235.0752 Iteration 1: log pseudolikelihood = -3235.0752 Negative binomial regression Number of obs = 1878 Dispersion = mean Wald chi2(0) = . Log pseudolikelihood = -3235.0752 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- _cons | .6612267 .0173034 38.21 0.000 .6273126 .6951408 -------------+---------------------------------------------------------------- /lnalpha | -3.128509 .5480945 -4.202754 -2.054263 -------------+---------------------------------------------------------------- alpha | .043783 .0239972 .0149543 .1281872 ------------------------------------------------------------------------------ . estimates store NB . . * Finite mixtures Poisson - 2 components . * Default start values . fmm children, components(2) mixtureof(poisson) vce(robust) Fitting Poisson model: Iteration 0: log likelihood = -3238.338 Iteration 1: log likelihood = -3238.338 Fitting 2 component Poisson model: Iteration 0: log pseudolikelihood = -3238.5983 (not concave) Iteration 1: log pseudolikelihood = -3238.2633 (not concave) Iteration 2: log pseudolikelihood = -3238.2413 (not concave) Iteration 3: log pseudolikelihood = -3236.8122 (not concave) Iteration 4: log pseudolikelihood = -3235.5237 (not concave) Iteration 5: log pseudolikelihood = -3234.948 (not concave) Iteration 6: log pseudolikelihood = -3232.0024 Iteration 7: log pseudolikelihood = -3227.5113 Iteration 8: log pseudolikelihood = -3226.1319 Iteration 9: log pseudolikelihood = -3225.6678 Iteration 10: log pseudolikelihood = -3225.6625 Iteration 11: log pseudolikelihood = -3225.6625 2 component Poisson regression Number of obs = 1878 Wald chi2(0) = . Log pseudolikelihood = -3225.6625 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- component1 | _cons | .6309416 .0176963 35.65 0.000 .5962575 .6656258 -------------+---------------------------------------------------------------- component2 | _cons | 1.883328 .0786408 23.95 0.000 1.729195 2.037461 -------------+---------------------------------------------------------------- /imlogitpi1 | 4.385301 .3626456 12.09 0.000 3.674529 5.096073 ------------------------------------------------------------------------------ pi1 | .9876942 .0044077 .9752659 .9939165 pi2 | .0123058 .0044077 .0060835 .0247341 ------------------------------------------------------------------------------ . if (e(converged) == 0) display " *** FMM DID NOT CONVERGE *** " . estimates store FMP2a . matrix bfmp2a = e(b) . scalar mu1 = exp(bfmp2a[1,1]) . scalar mu2 = exp(bfmp2a[1,2]) . scalar pi = exp(bfmp2a[1,3])/(1+exp(bfmp2a[1,3])) . display "Mixture probability = " pi " and Poisson means = " mu1 " and " mu2 Mixture probability = .98769418 and Poisson means = 1.8793794 and 6.5753525 . . * Finite mixtures Poisson - 2 components . * Different start values leads to a higher log-likleihood . * This is at the boundary . fmm children, components(2) mixtureof(poisson) vce(robust) from(1 1 3) Fitting 2 component Poisson model: Iteration 0: log pseudolikelihood = -3472.814 (not concave) Iteration 1: log pseudolikelihood = -3226.3697 (not concave) Iteration 2: log pseudolikelihood = -3210.4199 Iteration 3: log pseudolikelihood = -3203.3573 Iteration 4: log pseudolikelihood = -3202.8832 Iteration 5: log pseudolikelihood = -3202.7758 Iteration 6: log pseudolikelihood = -3202.7525 Iteration 7: log pseudolikelihood = -3202.7467 Iteration 8: log pseudolikelihood = -3202.7456 Iteration 9: log pseudolikelihood = -3202.7453 Iteration 10: log pseudolikelihood = -3202.7452 Iteration 11: log pseudolikelihood = -3202.7452 2 component Poisson regression Number of obs = 1878 Wald chi2(0) = . Log pseudolikelihood = -3202.7452 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- component1 | _cons | .7617174 .0179492 42.44 0.000 .7265375 .7968972 -------------+---------------------------------------------------------------- component2 | _cons | -16.88518 .1460314 -115.63 0.000 -17.17139 -16.59896 -------------+---------------------------------------------------------------- /imlogitpi1 | 2.246974 .1326573 16.94 0.000 1.98697 2.506977 ------------------------------------------------------------------------------ pi1 | .9043892 .0114708 .8794222 .9246295 pi2 | .0956108 .0114708 .0753705 .1205778 ------------------------------------------------------------------------------ . if (e(converged) == 0) display " *** FMM DID NOT CONVERGE *** " . estimates store FMP2b . matrix bfmp2a = e(b) . scalar mu1 = exp(bfmp2a[1,1]) . scalar mu2 = exp(bfmp2a[1,2]) . scalar pi = exp(bfmp2a[1,3])/(1+exp(bfmp2a[1,3])) . display "Mixture probability = " pi " and Poisson means = " mu1 " and " mu2 Mixture probability = .90438917 and Poisson means = 2.1419516 and 4.644e-08 . . * Following does not converge. Included for completeness. . * Finite mixtures NB - 2 components . fmm children, components(2) mixtureof(negbin2) vce(robust) iter(20) Fitting Negative Binomial-2 model: Iteration 0: log likelihood = -3238.338 Iteration 1: log likelihood = -3238.338 Iteration 0: log likelihood = -3537.6471 Iteration 1: log likelihood = -3235.3125 Iteration 2: log likelihood = -3235.0786 Iteration 3: log likelihood = -3235.0752 Iteration 4: log likelihood = -3235.0752 Iteration 0: log likelihood = -3235.0752 Iteration 1: log likelihood = -3235.0752 Fitting 2 component Negative Binomial-2 model: Iteration 0: log pseudolikelihood = -3235.3868 (not concave) Iteration 1: log pseudolikelihood = -3230.5866 (not concave) Iteration 2: log pseudolikelihood = -3209.5363 (not concave) Iteration 3: log pseudolikelihood = -3205.2835 (not concave) Iteration 4: log pseudolikelihood = -3203.6939 (not concave) Iteration 5: log pseudolikelihood = -3203.0756 (not concave) Iteration 6: log pseudolikelihood = -3202.8565 (not concave) Iteration 7: log pseudolikelihood = -3202.7824 (not concave) Iteration 8: log pseudolikelihood = -3202.7574 (not concave) Iteration 9: log pseudolikelihood = -3202.7494 (not concave) Iteration 10: log pseudolikelihood = -3202.7465 (not concave) Iteration 11: log pseudolikelihood = -3202.7457 (not concave) Iteration 12: log pseudolikelihood = -3202.7454 (not concave) Iteration 13: log pseudolikelihood = -3202.7453 (not concave) Iteration 14: log pseudolikelihood = -3202.7453 (not concave) Iteration 15: log pseudolikelihood = -3202.7452 (not concave) Iteration 16: log pseudolikelihood = -3202.7452 (not concave) Iteration 17: log pseudolikelihood = -3200.3783 (not concave) Iteration 18: log pseudolikelihood = -3198.4506 (not concave) Iteration 19: log pseudolikelihood = -3189.15 (not concave) Iteration 20: log pseudolikelihood = -3176.5235 (not concave) convergence not achieved 2 component Negative Binomial-2 regression Number of obs = 1878 Wald chi2(0) = . Log pseudolikelihood = -3176.5235 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- component1 | _cons | .7617287 .017948 42.44 0.000 .7265513 .7969062 -------------+---------------------------------------------------------------- component2 | _cons | -202.1108 . . . . . -------------+---------------------------------------------------------------- /imlogitpi1 | 2.246952 .1326344 16.94 0.000 1.986993 2.506911 /lnalpha1 | -28.58679 .0550694 -519.10 0.000 -28.69473 -28.47886 /lnalpha2 | -347.3441 . . . . . ------------------------------------------------------------------------------ Warning: convergence not achieved alpha1 | 3.85e-13 2.12e-14 3.45e-13 4.28e-13 alpha2 | 1.4e-151 . . . pi1 | .9043873 .011469 .8794247 .9246249 pi2 | .0956127 .011469 .0753751 .1205753 ------------------------------------------------------------------------------ . if (e(converged) == 0) display " *** FMM DID NOT CONVERGE *** " *** FMM DID NOT CONVERGE *** . estimates store FMNB2 . . * Following may not converge. Included for completeness. . * Finite mixtures Poisson - 3 components . quietly fmm children, components(2) mixtureof(poisson) vce(robust) . fmm children, components(3) mixtureof(poisson) vce(robust) iter(20) Fitting 3 component Poisson model: Iteration 0: log pseudolikelihood = -3225.6819 (not concave) Iteration 1: log pseudolikelihood = -3225.6735 (not concave) Iteration 2: log pseudolikelihood = -3225.6655 (not concave) Iteration 3: log pseudolikelihood = -3225.6627 Iteration 4: log pseudolikelihood = -3225.6625 (not concave) Iteration 5: log pseudolikelihood = -3225.6625 3 component Poisson regression Number of obs = 1878 Wald chi2(0) = . Log pseudolikelihood = -3225.6625 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- component1 | _cons | .6309378 .0176963 35.65 0.000 .5962537 .6656219 -------------+---------------------------------------------------------------- component2 | _cons | 1.88333 .0786404 23.95 0.000 1.729198 2.037462 -------------+---------------------------------------------------------------- component3 | _cons | .6309144 .0177005 35.64 0.000 .596222 .6656069 -------------+---------------------------------------------------------------- /imlogitpi1 | 3.020858 .0259259 116.52 0.000 2.970044 3.071671 /imlogitpi2 | -1.316853 .3636751 -3.62 0.000 -2.029643 -.6040623 ------------------------------------------------------------------------------ pi1 | .9417741 .0043462 .9326477 .9497305 pi2 | .0123056 .0044079 .0060831 .0247347 pi3 | .0459203 .0011549 .0436568 .0481839 ------------------------------------------------------------------------------ . if (e(converged) == 0) display " *** FMM DID NOT CONVERGE *** " . estimates store FMNB2 . . * Hurdle Poisson model . hplogit children, vce(robust) initial: log pseudolikelihood = -4119.0571 alternative: log pseudolikelihood = -3418.0053 rescale: log pseudolikelihood = -3307.5599 rescale eq: log pseudolikelihood = -3305.4942 Iteration 0: log pseudolikelihood = -3305.4942 Iteration 1: log pseudolikelihood = -3206.1791 Iteration 2: log pseudolikelihood = -3202.7486 Iteration 3: log pseudolikelihood = -3202.7452 Iteration 4: log pseudolikelihood = -3202.7452 Poisson-Logit Hurdle Regression Number of obs = 1878 Wald chi2(0) = . Log pseudolikelihood = -3202.7452 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- logit | _cons | 1.375017 .05751 23.91 0.000 1.2623 1.487735 -------------+---------------------------------------------------------------- poisson | _cons | .7617264 .0179489 42.44 0.000 .7265471 .7969057 ------------------------------------------------------------------------------ AIC Statistic = 3.412 . estimates store HP . . * Hurdle negative binomial model . hnblogit children, vce(robust) initial: log pseudolikelihood = -3823.3998 alternative: log pseudolikelihood = -3644.3622 rescale: log pseudolikelihood = -3644.3622 rescale eq: log pseudolikelihood = -3300.5545 Iteration 0: log pseudolikelihood = -3300.5545 Iteration 1: log pseudolikelihood = -3229.719 Iteration 2: log pseudolikelihood = -3203.8036 Iteration 3: log pseudolikelihood = -3203.0173 Iteration 4: log pseudolikelihood = -3202.8057 Iteration 5: log pseudolikelihood = -3202.7587 Iteration 6: log pseudolikelihood = -3202.7482 Iteration 7: log pseudolikelihood = -3202.7455 Iteration 8: log pseudolikelihood = -3202.7453 Iteration 9: log pseudolikelihood = -3202.7451 Iteration 10: log pseudolikelihood = -3202.745 Iteration 11: log pseudolikelihood = -3202.745 Negative Binomial-Logit Hurdle Regression Number of obs = 1878 Wald chi2(0) = . Log pseudolikelihood = -3202.745 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- logit | _cons | 1.375017 .05751 23.91 0.000 1.2623 1.487735 -------------+---------------------------------------------------------------- negbinomial | _cons | .7617274 .0169528 44.93 0.000 .7285006 .7949543 -------------+---------------------------------------------------------------- /lnalpha | -18.45167 .0027756 -6647.81 0.000 -18.45711 -18.44623 ------------------------------------------------------------------------------ AIC Statistic = 3.412 . estimates store HNB . . * The following predicted frequencies are reported in the text . * Zero Inflated Poisson Model . zip children, inflate(_cons) vce(robust) Fitting constant-only model: Iteration 0: log pseudolikelihood = -3754.8579 Iteration 1: log pseudolikelihood = -3231.043 Iteration 2: log pseudolikelihood = -3203.6841 Iteration 3: log pseudolikelihood = -3202.7496 Iteration 4: log pseudolikelihood = -3202.7452 Iteration 5: log pseudolikelihood = -3202.7452 Fitting full model: Iteration 0: log pseudolikelihood = -3202.7452 Iteration 1: log pseudolikelihood = -3202.7452 Zero-inflated Poisson regression Number of obs = 1878 Nonzero obs = 1499 Zero obs = 379 Inflation model = logit Wald chi2(0) = . Log pseudolikelihood = -3202.745 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- children | _cons | .7617264 .0179489 42.44 0.000 .7265471 .7969057 -------------+---------------------------------------------------------------- inflate | _cons | -2.24693 .1326501 -16.94 0.000 -2.506919 -1.98694 ------------------------------------------------------------------------------ . estimates store ZIP . forvalues i = 0/12 { 2. predict zipfit`i', pr(`i') 3. } . sum zipfit* Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- zipfit0 | 1878 .2018104 0 .2018104 .2018104 zipfit1 | 1878 .2274683 0 .2274683 .2274683 zipfit2 | 1878 .2436153 0 .2436153 .2436153 zipfit3 | 1878 .1739389 0 .1739389 .1739389 zipfit4 | 1878 .093143 0 .093143 .093143 -------------+-------------------------------------------------------- zipfit5 | 1878 .0399019 0 .0399019 .0399019 zipfit6 | 1878 .0142448 0 .0142448 .0142448 zipfit7 | 1878 .0043588 0 .0043588 .0043588 zipfit8 | 1878 .0011671 0 .0011671 .0011671 zipfit9 | 1878 .0002778 0 .0002778 .0002778 -------------+-------------------------------------------------------- zipfit10 | 1878 .0000595 0 .0000595 .0000595 zipfit11 | 1878 .0000116 0 .0000116 .0000116 zipfit12 | 1878 2.07e-06 0 2.07e-06 2.07e-06 . . * Zero Inflated NB Model . zip children, inflate(_cons) vce(robust) Fitting constant-only model: Iteration 0: log pseudolikelihood = -3754.8579 Iteration 1: log pseudolikelihood = -3231.043 Iteration 2: log pseudolikelihood = -3203.6841 Iteration 3: log pseudolikelihood = -3202.7496 Iteration 4: log pseudolikelihood = -3202.7452 Iteration 5: log pseudolikelihood = -3202.7452 Fitting full model: Iteration 0: log pseudolikelihood = -3202.7452 Iteration 1: log pseudolikelihood = -3202.7452 Zero-inflated Poisson regression Number of obs = 1878 Nonzero obs = 1499 Zero obs = 379 Inflation model = logit Wald chi2(0) = . Log pseudolikelihood = -3202.745 Prob > chi2 = . ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- children | _cons | .7617264 .0179489 42.44 0.000 .7265471 .7969057 -------------+---------------------------------------------------------------- inflate | _cons | -2.24693 .1326501 -16.94 0.000 -2.506919 -1.98694 ------------------------------------------------------------------------------ . estimates store ZINB . . * Ordered probit . generate childrange = children . replace childrange = 6 if children >= 6 (21 real changes made) . tabulate childrange childrange | Freq. Percent Cum. ------------+----------------------------------- 0 | 379 20.18 20.18 1 | 262 13.95 34.13 2 | 684 36.42 70.55 3 | 353 18.80 89.35 4 | 128 6.82 96.17 5 | 35 1.86 98.03 6 | 37 1.97 100.00 ------------+----------------------------------- Total | 1,878 100.00 . oprobit childrange, vce(robust) Iteration 0: log pseudolikelihood = -3031.9739 Iteration 1: log pseudolikelihood = -3031.9739 (backed up) Ordered probit regression Number of obs = 1878 Wald chi2(0) = . Prob > chi2 = . Log pseudolikelihood = -3031.9739 Pseudo R2 = 0.0000 ------------------------------------------------------------------------------ | Robust childrange | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- /cut1 | -.835172 .0329113 -.899677 -.770667 /cut2 | -.4088618 .0298247 -.4673171 -.3504064 /cut3 | .5403954 .0305173 .4805826 .6002081 /cut4 | 1.245379 .038758 1.169414 1.321343 /cut5 | 1.770299 .0532375 1.665956 1.874643 /cut6 | 2.059947 .0671009 1.928431 2.191462 ------------------------------------------------------------------------------ . estimates store OPROBIT . predict pop0 pop1 pop2 pop3 pop4 pop5 pop6 (option pr assumed; predicted probabilities) . summarize pop* Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- pop0 | 1878 .2018104 0 .2018104 .2018104 pop1 | 1878 .1395101 0 .1395101 .1395101 pop2 | 1878 .3642173 0 .3642173 .3642173 pop3 | 1878 .1879659 0 .1879659 .1879659 pop4 | 1878 .0681576 0 .0681576 .0681576 -------------+-------------------------------------------------------- pop5 | 1878 .0186368 0 .0186368 .0186368 pop6 | 1878 .0197018 0 .0197018 .0197018 . . *** TABLE 6.17 FIT (and estimates) FOR VARIOUS MODELS (FMP2b preferred to FMP2a) . . estimates table POISSON NB FMP2a FMP2b HP ZIP, b(%10.3f) t(%10.2f) stats(ll aic bic N k) -------------------------------------------------------------------------------------------- Variable | POISSON NB FMP2a FMP2b HP ZIP -------------+------------------------------------------------------------------------------ children | _cons | 0.661 0.661 0.762 | 38.21 38.21 42.44 -------------+------------------------------------------------------------------------------ lnalpha | _cons | -3.129 | -5.71 -------------+------------------------------------------------------------------------------ component1 | _cons | 0.631 0.762 | 35.65 42.44 -------------+------------------------------------------------------------------------------ component2 | _cons | 1.883 -16.885 | 23.95 -115.63 -------------+------------------------------------------------------------------------------ imlogitpi1 | _cons | 4.385 2.247 | 12.09 16.94 -------------+------------------------------------------------------------------------------ logit | _cons | 1.375 | 23.91 -------------+------------------------------------------------------------------------------ poisson | _cons | 0.762 | 42.44 -------------+------------------------------------------------------------------------------ inflate | _cons | -2.247 | -16.94 -------------+------------------------------------------------------------------------------ Statistics | ll | -3238.338 -3235.075 -3225.663 -3202.745 -3202.745 -3202.745 aic | 6478.676 6474.150 6457.325 6411.490 6409.490 6409.490 bic | 6484.214 6485.226 6473.939 6428.104 6420.566 6420.566 N | 1878 1878 1878 1878 1878 1878 k | 1.000 2.000 3.000 3.000 2.000 2.000 -------------------------------------------------------------------------------------------- legend: b/t . estimates table OPROBIT, b(%10.3f) t(%10.2f) stats(ll aic bic N k) --------------------------- Variable | OPROBIT -------------+------------- cut1 | _cons | -0.835 | -25.38 -------------+------------- cut2 | _cons | -0.409 | -13.71 -------------+------------- cut3 | _cons | 0.540 | 17.71 -------------+------------- cut4 | _cons | 1.245 | 32.13 -------------+------------- cut5 | _cons | 1.770 | 33.25 -------------+------------- cut6 | _cons | 2.060 | 30.70 -------------+------------- Statistics | ll | -3031.974 aic | 6075.948 bic | 6109.176 N | 1878 k | 6.000 --------------------------- legend: b/t . . ********** MODEL WITH REGRESSORS - results do not change much . . * Do ZIP with regressors: . * educational dummies, religion dummies, language dunnies, age, age-squared . tabulate education, generate(deduc) highest level of education | achieved | Freq. Percent Cum. --------------------------------+----------------------------------- compulsory or less | 441 23.48 23.48 domestic science course, 1yr | 197 10.49 33.97 general training/apprenticeship | 643 34.24 68.21 high school | 483 25.72 93.93 university/postgrad | 114 6.07 100.00 --------------------------------+----------------------------------- Total | 1,878 100.00 . tabulate intlang, generate(dlang) interview | language | Freq. Percent Cum. --------------+----------------------------------- french | 515 27.42 27.42 german | 1,269 67.57 94.99 italian | 94 5.01 100.00 --------------+----------------------------------- Total | 1,878 100.00 . tabulate religion, generate(drelig) religion | Freq. Percent Cum. ------------------------------+----------------------------------- protestant or reformed church | 858 45.69 45.69 roman catholic | 763 40.63 86.32 other | 100 5.32 91.64 no denomination or religion | 157 8.36 100.00 ------------------------------+----------------------------------- Total | 1,878 100.00 . generate agesq = age^2 . drop deduc1 dlang1 drelig1 // drop one dummy in each category . global XLIST age agesq deduc* dlang* drelig* . poisson children $XLIST, vce(robust) Iteration 0: log pseudolikelihood = -3199.9228 Iteration 1: log pseudolikelihood = -3199.9228 Poisson regression Number of obs = 1878 Wald chi2(11) = 66.62 Prob > chi2 = 0.0000 Log pseudolikelihood = -3199.9228 Pseudo R2 = 0.0119 ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- age | .0244859 .0212412 1.15 0.249 -.0171461 .0661179 agesq | -.0001866 .000173 -1.08 0.281 -.0005256 .0001524 deduc2 | -.0496454 .0641372 -0.77 0.439 -.1753519 .0760611 deduc3 | -.1882767 .0469727 -4.01 0.000 -.2803414 -.096212 deduc4 | -.170919 .0508149 -3.36 0.001 -.2705144 -.0713237 deduc5 | -.2835221 .0839582 -3.38 0.001 -.4480772 -.118967 dlang2 | .1621923 .0393305 4.12 0.000 .085106 .2392787 dlang3 | -.0889937 .0731753 -1.22 0.224 -.2324146 .0544273 drelig2 | .089554 .0374536 2.39 0.017 .0161462 .1629618 drelig3 | .0881417 .0865674 1.02 0.309 -.0815273 .2578106 drelig4 | -.1865608 .0756966 -2.46 0.014 -.3349235 -.0381981 _cons | -.1237865 .6407632 -0.19 0.847 -1.379659 1.132086 ------------------------------------------------------------------------------ . zip children $XLIST, inflate($XLIST) vce(robust) Fitting constant-only model: Iteration 0: log pseudolikelihood = -3754.8579 Iteration 1: log pseudolikelihood = -3220.1107 Iteration 2: log pseudolikelihood = -3187.9915 Iteration 3: log pseudolikelihood = -3186.0322 Iteration 4: log pseudolikelihood = -3185.9825 Iteration 5: log pseudolikelihood = -3185.9822 Fitting full model: Iteration 0: log pseudolikelihood = -3185.9822 Iteration 1: log pseudolikelihood = -3147.4632 Iteration 2: log pseudolikelihood = -3144.4006 Iteration 3: log pseudolikelihood = -3144.38 Iteration 4: log pseudolikelihood = -3144.38 Zero-inflated Poisson regression Number of obs = 1878 Nonzero obs = 1499 Zero obs = 379 Inflation model = logit Wald chi2(11) = 97.62 Log pseudolikelihood = -3144.38 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ | Robust children | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- children | age | .0082131 .0206712 0.40 0.691 -.0323017 .0487278 agesq | -6.93e-06 .0001686 -0.04 0.967 -.0003374 .0003235 deduc2 | -.0175648 .06108 -0.29 0.774 -.1372795 .1021499 deduc3 | -.1865751 .0466976 -4.00 0.000 -.2781008 -.0950495 deduc4 | -.1210487 .0506028 -2.39 0.017 -.2202284 -.0218691 deduc5 | -.1817327 .0831892 -2.18 0.029 -.3447804 -.0186849 dlang2 | .2210996 .0409973 5.39 0.000 .1407463 .3014529 dlang3 | -.1250416 .0734466 -1.70 0.089 -.2689943 .018911 drelig2 | .1060754 .0363455 2.92 0.004 .0348394 .1773113 drelig3 | .0786934 .0971528 0.81 0.418 -.1117225 .2691093 drelig4 | -.1261827 .0740568 -1.70 0.088 -.2713314 .018966 _cons | .2200338 .6228478 0.35 0.724 -1.000725 1.440793 -------------+---------------------------------------------------------------- inflate | age | .0154545 .1300876 0.12 0.905 -.2395125 .2704215 agesq | .0003402 .0009798 0.35 0.728 -.0015801 .0022606 deduc2 | .2786851 .4415598 0.63 0.528 -.5867562 1.144127 deduc3 | -.0360749 .3892149 -0.09 0.926 -.7989222 .7267724 deduc4 | .4910239 .3983761 1.23 0.218 -.2897788 1.271827 deduc5 | .9742744 .5710739 1.71 0.088 -.14501 2.093559 dlang2 | .7392459 .4067751 1.82 0.069 -.0580186 1.53651 dlang3 | -.730013 1.35412 -0.54 0.590 -3.384039 1.924013 drelig2 | .1649629 .2759427 0.60 0.550 -.3758749 .7058006 drelig3 | -.0962057 .9194979 -0.10 0.917 -1.898389 1.705977 drelig4 | .6351458 .4818118 1.32 0.187 -.3091879 1.57948 _cons | -5.419593 4.281145 -1.27 0.206 -13.81048 2.971298 ------------------------------------------------------------------------------ . forvalues i = 0/12 { 2. predict zipregfit`i', pr(`i') 3. } . sum zipregfit* Variable | Obs Mean Std. Dev. Min Max -------------+-------------------------------------------------------- zipregfit0 | 1878 .2065877 .0555878 .1200961 .6293052 zipregfit1 | 1878 .2302628 .0525122 .0789659 .3405415 zipregfit2 | 1878 .2379544 .0238273 .1078057 .2659954 zipregfit3 | 1878 .1685309 .0222554 .0845832 .2057762 zipregfit4 | 1878 .0920025 .0249215 .0302283 .1422742 -------------+-------------------------------------------------------- zipregfit5 | 1878 .0412731 .0175869 .0079013 .0857975 zipregfit6 | 1878 .0158386 .0093647 .0017211 .0453241 zipregfit7 | 1878 .0053433 .0040979 .0003213 .0212788 zipregfit8 | 1878 .001616 .0015421 .0000525 .0088514 zipregfit9 | 1878 .0004446 .0005131 7.62e-06 .0032729 -------------+-------------------------------------------------------- zipregfit10 | 1878 .0001125 .0001538 9.96e-07 .0010891 zipregfit11 | 1878 .0000264 .0000421 1.18e-07 .0003295 zipregfit12 | 1878 5.79e-06 .0000106 1.29e-08 .0000914 . . ********** CLOSE OUTPUT . . * log close . * clear . * exit . . . end of do-file . exit, clear