. * do-file for lecture 12a of VHM 812/802, Winter 2022 . version 17 /* works also with versions 14-16 */ . set more off . cd "r:\" r:\ . . * unadjusted and mixed model results (for reference) . use simcont_clustherd.dta, clear . regress milk X Source | SS df MS Number of obs = 11,626 -------------+---------------------------------- F(1, 11624) = 317.72 Model | 36598.5078 1 36598.5078 Prob > F = 0.0000 Residual | 1338999 11,624 115.192618 R-squared = 0.0266 -------------+---------------------------------- Adj R-squared = 0.0265 Total | 1375597.51 11,625 118.330968 Root MSE = 10.733 ------------------------------------------------------------------------------ milk | Coefficient Std. err. t P>|t| [95% conf. interval] -------------+---------------------------------------------------------------- X | 3.55661 .199534 17.82 0.000 3.16549 3.94773 _cons | 30.0215 .1457715 205.95 0.000 29.73576 30.30723 ------------------------------------------------------------------------------ . mixed milk X || herd:, reml Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -40902.479 Iteration 1: log restricted-likelihood = -40902.479 (backed up) Computing standard errors ... Mixed-effects REML regression Number of obs = 11,626 Group variable: herd Number of groups = 100 Obs per group: min = 20 avg = 116.3 max = 311 Wald chi2(1) = 6.44 Log restricted-likelihood = -40902.479 Prob > chi2 = 0.0112 ------------------------------------------------------------------------------ milk | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- X | 3.796004 1.495942 2.54 0.011 .8640104 6.727997 _cons | 31.13696 1.058717 29.41 0.000 29.06191 33.21201 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herd: Identity | var(_cons) | 54.91494 7.998609 41.2771 73.05868 -----------------------------+------------------------------------------------ var(Residual) | 64.20087 .8457062 62.56453 65.88001 ------------------------------------------------------------------------------ LR test vs. linear model: chibar2(01) = 6374.40 Prob >= chibar2 = 0.0000 . collapse milk X, by(herd) . regress milk X Source | SS df MS Number of obs = 100 -------------+---------------------------------- F(1, 98) = 6.37 Model | 356.97798 1 356.97798 Prob > F = 0.0132 Residual | 5493.56229 98 56.056758 R-squared = 0.0610 -------------+---------------------------------- Adj R-squared = 0.0514 Total | 5850.54027 99 59.0963664 Root MSE = 7.4871 ------------------------------------------------------------------------------ milk | Coefficient Std. err. t P>|t| [95% conf. interval] -------------+---------------------------------------------------------------- X | 3.778772 1.497421 2.52 0.013 .8071885 6.750356 _cons | 31.16586 1.058837 29.43 0.000 29.06463 33.26708 ------------------------------------------------------------------------------ . use simcont_clustcow.dta, clear . regress milk X Source | SS df MS Number of obs = 11,626 -------------+---------------------------------- F(1, 11624) = 624.90 Model | 72138.7619 1 72138.7619 Prob > F = 0.0000 Residual | 1341880.62 11,624 115.440522 R-squared = 0.0510 -------------+---------------------------------- Adj R-squared = 0.0509 Total | 1414019.39 11,625 121.636076 Root MSE = 10.744 ------------------------------------------------------------------------------ milk | Coefficient Std. err. t P>|t| [95% conf. interval] -------------+---------------------------------------------------------------- X | 4.982006 .1992962 25.00 0.000 4.591352 5.37266 _cons | 29.25664 .1412627 207.11 0.000 28.97974 29.53354 ------------------------------------------------------------------------------ . mixed milk X || herd:, reml Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -40947.175 Iteration 1: log restricted-likelihood = -40947.175 (backed up) Computing standard errors ... Mixed-effects REML regression Number of obs = 11,626 Group variable: herd Number of groups = 100 Obs per group: min = 20 avg = 116.3 max = 311 Wald chi2(1) = 1108.56 Log restricted-likelihood = -40947.175 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ milk | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- X | 4.968194 .1492174 33.30 0.000 4.675733 5.260655 _cons | 30.64647 .7281274 42.09 0.000 29.21936 32.07357 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herd: Identity | var(_cons) | 51.41187 7.459585 38.68643 68.32319 -----------------------------+------------------------------------------------ var(Residual) | 64.71069 .8524578 63.06129 66.40324 ------------------------------------------------------------------------------ LR test vs. linear model: chibar2(01) = 6310.00 Prob >= chibar2 = 0.0000 . use simbin_clustherd.dta, clear . logit Y X Iteration 0: log likelihood = -6894.3552 Iteration 1: log likelihood = -6815.0583 Iteration 2: log likelihood = -6814.7785 Iteration 3: log likelihood = -6814.7785 Logistic regression Number of obs = 11,626 LR chi2(1) = 159.15 Prob > chi2 = 0.0000 Log likelihood = -6814.7785 Pseudo R2 = 0.0115 ------------------------------------------------------------------------------ Y | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- X | .5287317 .0423191 12.49 0.000 .4457877 .6116757 _cons | -1.241768 .0325699 -38.13 0.000 -1.305604 -1.177932 ------------------------------------------------------------------------------ . melogit Y X || herd: Fitting fixed-effects model: Iteration 0: log likelihood = -6828.9777 Iteration 1: log likelihood = -6814.7876 Iteration 2: log likelihood = -6814.7785 Iteration 3: log likelihood = -6814.7785 Refining starting values: Grid node 0: log likelihood = -6065.269 Fitting full model: Iteration 0: log likelihood = -6065.269 Iteration 1: log likelihood = -6065.089 Iteration 2: log likelihood = -6065.0864 Iteration 3: log likelihood = -6065.0864 Mixed-effects logistic regression Number of obs = 11,626 Group variable: herd Number of groups = 100 Obs per group: min = 20 avg = 116.3 max = 311 Integration method: mvaghermite Integration pts. = 7 Wald chi2(1) = 9.26 Log likelihood = -6065.0864 Prob > chi2 = 0.0023 ------------------------------------------------------------------------------ Y | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- X | .6199967 .2037578 3.04 0.002 .2206389 1.019355 _cons | -1.305448 .1454551 -8.97 0.000 -1.590534 -1.020361 -------------+---------------------------------------------------------------- herd | var(_cons)| .9417563 .1493109 .6902154 1.284968 ------------------------------------------------------------------------------ LR test vs. logistic model: chibar2(01) = 1499.38 Prob >= chibar2 = 0.0000 . use simbin_clustcow.dta, clear . logit Y X Iteration 0: log likelihood = -6910.3442 Iteration 1: log likelihood = -6811.48 Iteration 2: log likelihood = -6811.0741 Iteration 3: log likelihood = -6811.0741 Logistic regression Number of obs = 11,626 LR chi2(1) = 198.54 Prob > chi2 = 0.0000 Log likelihood = -6811.0741 Pseudo R2 = 0.0144 ------------------------------------------------------------------------------ Y | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- X | .5863084 .0419748 13.97 0.000 .5040393 .6685775 _cons | -1.25032 .0316033 -39.56 0.000 -1.312261 -1.188379 ------------------------------------------------------------------------------ . melogit Y X || herd: Fitting fixed-effects model: Iteration 0: log likelihood = -6824.417 Iteration 1: log likelihood = -6811.0819 Iteration 2: log likelihood = -6811.0741 Iteration 3: log likelihood = -6811.0741 Refining starting values: Grid node 0: log likelihood = -5999.0535 Fitting full model: Iteration 0: log likelihood = -5999.0535 Iteration 1: log likelihood = -5995.9716 Iteration 2: log likelihood = -5995.9694 Iteration 3: log likelihood = -5995.9694 Mixed-effects logistic regression Number of obs = 11,626 Group variable: herd Number of groups = 100 Obs per group: min = 20 avg = 116.3 max = 311 Integration method: mvaghermite Integration pts. = 7 Wald chi2(1) = 229.28 Log likelihood = -5995.9694 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ Y | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- X | .6974798 .046063 15.14 0.000 .6071979 .7877616 _cons | -1.361196 .1111563 -12.25 0.000 -1.579059 -1.143334 -------------+---------------------------------------------------------------- herd | var(_cons)| 1.068314 .1682536 .7845836 1.45465 ------------------------------------------------------------------------------ LR test vs. logistic model: chibar2(01) = 1630.21 Prob >= chibar2 = 0.0000 . . * 2-level somatic cell count data . use scc40_2level, clear . mixed t_lnscc h_size c_heifer i.t_season t_dim || herdid:, reml Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -3624.9622 Iteration 1: log restricted-likelihood = -3624.9622 Computing standard errors ... Mixed-effects REML regression Number of obs = 2,178 Group variable: herdid Number of groups = 40 Obs per group: min = 12 avg = 54.5 max = 105 Wald chi2(6) = 244.36 Log restricted-likelihood = -3624.9622 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ t_lnscc | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- h_size | .0040837 .0037726 1.08 0.279 -.0033105 .0114778 c_heifer | -.7367168 .0554447 -13.29 0.000 -.8453863 -.6280472 | t_season | apr-jun | .1609431 .0906574 1.78 0.076 -.0167422 .3386285 jul-sep | .0015031 .0863774 0.02 0.986 -.1677935 .1707997 oct-dec | .0014582 .0918936 0.02 0.987 -.1786499 .1815663 | t_dim | .0027731 .0004991 5.56 0.000 .0017949 .0037513 _cons | 4.641202 .1974215 23.51 0.000 4.254263 5.028141 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herdid: Identity | var(_cons) | .1491533 .0436191 .0840821 .2645832 -----------------------------+------------------------------------------------ var(Residual) | 1.557228 .0477206 1.466451 1.653625 ------------------------------------------------------------------------------ LR test vs. linear model: chibar2(01) = 97.01 Prob >= chibar2 = 0.0000 . mixed t_lnscc h_size c_heifer i.t_season t_dim || herdid:, reml stddev Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -3624.9622 Iteration 1: log restricted-likelihood = -3624.9622 Computing standard errors ... Mixed-effects REML regression Number of obs = 2,178 Group variable: herdid Number of groups = 40 Obs per group: min = 12 avg = 54.5 max = 105 Wald chi2(6) = 244.36 Log restricted-likelihood = -3624.9622 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ t_lnscc | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- h_size | .0040837 .0037726 1.08 0.279 -.0033105 .0114778 c_heifer | -.7367168 .0554447 -13.29 0.000 -.8453863 -.6280472 | t_season | apr-jun | .1609431 .0906574 1.78 0.076 -.0167422 .3386285 jul-sep | .0015031 .0863774 0.02 0.986 -.1677935 .1707997 oct-dec | .0014582 .0918936 0.02 0.987 -.1786499 .1815663 | t_dim | .0027731 .0004991 5.56 0.000 .0017949 .0037513 _cons | 4.641202 .1974215 23.51 0.000 4.254263 5.028141 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herdid: Identity | sd(_cons) | .3862037 .0564716 .2899691 .5143765 -----------------------------+------------------------------------------------ sd(Residual) | 1.247889 .0191205 1.210971 1.285933 ------------------------------------------------------------------------------ LR test vs. linear model: chibar2(01) = 97.01 Prob >= chibar2 = 0.0000 . di 1.96*sqrt(.1491533) .75695926 . di .1491533/(.1491533+1.557228) /* ICC and VPC */ .08740913 . testparm i.t_season /* Wald test for season */ ( 1) [t_lnscc]2.t_season = 0 ( 2) [t_lnscc]3.t_season = 0 ( 3) [t_lnscc]4.t_season = 0 chi2( 3) = 6.21 Prob > chi2 = 0.1017 . estimates store full . mixed t_lnscc h_size c_heifer i.t_season t_dim, reml Mixed-effects REML regression Number of obs = 2,178 Wald chi2(6) = 250.46 Log restricted-likelihood = -3673.4664 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ t_lnscc | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- h_size | .0073914 .001605 4.61 0.000 .0042456 .0105371 c_heifer | -.7258892 .0567071 -12.80 0.000 -.8370332 -.6147453 | t_season | apr-jun | .1058737 .0903724 1.17 0.241 -.0712529 .2830004 jul-sep | -.0454035 .0876933 -0.52 0.605 -.2172792 .1264723 oct-dec | -.0385968 .093818 -0.41 0.681 -.2224768 .1452832 | t_dim | .0030497 .0005049 6.04 0.000 .0020601 .0040393 _cons | 4.492766 .1067432 42.09 0.000 4.283553 4.701979 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ var(Residual) | 1.679409 .0508912 1.582569 1.782176 ------------------------------------------------------------------------------ . lrtest full Likelihood-ratio test Assumption: . nested within full LR chi2(1) = 97.01 Prob > chi2 = 0.0000 Note: The reported degrees of freedom assumes the null hypothesis is not on the boundary of the parameter space. If this is not true, then the reported test is conservative. Note: LR tests based on REML are valid only when the fixed-effects specification is identical for both models. . mixed t_lnscc || herdid:, reml /* "null" model ~ no predictors */ Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -3722.5104 Iteration 1: log restricted-likelihood = -3722.5104 Computing standard errors ... Mixed-effects REML regression Number of obs = 2,178 Group variable: herdid Number of groups = 40 Obs per group: min = 12 avg = 54.5 max = 105 Wald chi2(0) = . Log restricted-likelihood = -3722.5104 Prob > chi2 = . ------------------------------------------------------------------------------ t_lnscc | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- _cons | 4.746742 .0679798 69.83 0.000 4.613505 4.87998 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herdid: Identity | var(_cons) | .1482566 .0426236 .0843907 .2604553 -----------------------------+------------------------------------------------ var(Residual) | 1.730413 .0529421 1.629699 1.837352 ------------------------------------------------------------------------------ LR test vs. linear model: chibar2(01) = 101.05 Prob >= chibar2 = 0.0000 . . * Reunion Island data, 3-level model for (natural) log cfs . use reu_cfs, clear . mixed lncfs i.heifer || herd: || cow:, reml Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -1500.2224 Iteration 1: log restricted-likelihood = -1495.8851 Iteration 2: log restricted-likelihood = -1495.8804 Iteration 3: log restricted-likelihood = -1495.8804 Computing standard errors ... Mixed-effects REML regression Number of obs = 3,027 Grouping information ------------------------------------------------------------- | No. of Observations per group Group variable | groups Minimum Average Maximum ----------------+-------------------------------------------- herd | 50 13 60.5 226 cow | 1,575 1 1.9 5 ------------------------------------------------------------- Wald chi2(1) = 0.12 Log restricted-likelihood = -1495.8804 Prob > chi2 = 0.7341 ------------------------------------------------------------------------------ lncfs | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- heifer | primiparous | -.0058409 .017194 -0.34 0.734 -.0395405 .0278587 _cons | 4.21903 .0195052 216.30 0.000 4.180801 4.25726 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herd: Identity | var(_cons) | .0145193 .0036025 .0089279 .0236125 -----------------------------+------------------------------------------------ cow: Identity | var(_cons) | .0200572 .0040814 .0134605 .0298868 -----------------------------+------------------------------------------------ var(Residual) | .1341146 .0048118 .1250077 .143885 ------------------------------------------------------------------------------ LR test vs. linear model: chi2(2) = 221.39 Prob > chi2 = 0.0000 Note: LR test is conservative and provided only for reference. . * NOT: mixed lncfs i.heifer || cow: || herd:, reml . * ICC: same herd . di .0145193/(.0145193+.0200572+.1341146) .08607034 . * ICC: same cow . di (.0145193+.0200572)/(.0145193+.0200572+.1341146) .20496932 . * ICC calculation by Stata . estat icc Residual intraclass correlation ------------------------------------------------------------------------------ Level | ICC Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herd | .0860701 .0196752 .0545356 .1332689 cow|herd | .2049693 .0285332 .1546347 .2665226 ------------------------------------------------------------------------------ . * 4-level model with region random effects . mixed lncfs i.heifer || region: || herd: || cow:, reml Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -1500.237 Iteration 1: log restricted-likelihood = -1495.8809 Iteration 2: log restricted-likelihood = -1495.8694 Iteration 3: log restricted-likelihood = -1495.8694 Computing standard errors ... Mixed-effects REML regression Number of obs = 3,027 Grouping information ------------------------------------------------------------- | No. of Observations per group Group variable | groups Minimum Average Maximum ----------------+-------------------------------------------- region | 5 193 605.4 960 herd | 50 13 60.5 226 cow | 1,575 1 1.9 5 ------------------------------------------------------------- Wald chi2(1) = 0.12 Log restricted-likelihood = -1495.8694 Prob > chi2 = 0.7308 ------------------------------------------------------------------------------ lncfs | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- heifer | primiparous | -.0059168 .0171944 -0.34 0.731 -.0396173 .0277836 _cons | 4.218036 .0206238 204.52 0.000 4.177614 4.258458 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ region: Identity | var(_cons) | .0001844 .0013431 1.16e-10 292.2396 -----------------------------+------------------------------------------------ herd: Identity | var(_cons) | .0144022 .0036567 .008756 .0236892 -----------------------------+------------------------------------------------ cow: Identity | var(_cons) | .0200442 .0040818 .0134478 .0298764 -----------------------------+------------------------------------------------ var(Residual) | .1341222 .0048124 .1250142 .1438939 ------------------------------------------------------------------------------ LR test vs. linear model: chi2(3) = 221.41 Prob > chi2 = 0.0000 Note: LR test is conservative and provided only for reference. . * 4-level model with region fixed effects . mixed lncfs i.heifer i.region || herd: || cow:, reml Performing EM optimization ... Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -1505.8762 Iteration 1: log restricted-likelihood = -1501.496 Iteration 2: log restricted-likelihood = -1501.4927 Iteration 3: log restricted-likelihood = -1501.4927 Computing standard errors ... Mixed-effects REML regression Number of obs = 3,027 Grouping information ------------------------------------------------------------- | No. of Observations per group Group variable | groups Minimum Average Maximum ----------------+-------------------------------------------- herd | 50 13 60.5 226 cow | 1,575 1 1.9 5 ------------------------------------------------------------- Wald chi2(5) = 4.10 Log restricted-likelihood = -1501.4927 Prob > chi2 = 0.5347 ------------------------------------------------------------------------------ lncfs | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- heifer | primiparous | -.0065108 .0172027 -0.38 0.705 -.0402274 .0272058 | region | 2 | .0673308 .0854623 0.79 0.431 -.1001721 .2348338 3 | .0578045 .1019399 0.57 0.571 -.141994 .2576029 4 | .1382983 .0832458 1.66 0.097 -.0248604 .301457 5 | .0780979 .082739 0.94 0.345 -.0840675 .2402634 | _cons | 4.131972 .075647 54.62 0.000 3.983706 4.280237 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects parameters | Estimate Std. err. [95% conf. interval] -----------------------------+------------------------------------------------ herd: Identity | var(_cons) | .0146174 .0037663 .0088217 .0242209 -----------------------------+------------------------------------------------ cow: Identity | var(_cons) | .0199372 .0040785 .0133517 .0297709 -----------------------------+------------------------------------------------ var(Residual) | .134201 .0048163 .1250856 .1439806 ------------------------------------------------------------------------------ LR test vs. linear model: chi2(2) = 203.76 Prob > chi2 = 0.0000 Note: LR test is conservative and provided only for reference. . . * pig data, binary outcome . use pig_adg, clear (Growth Performance and Abattoir Surveillance Data on Pigs, Real Data) . generate ar_g1=ar>1 if ar~=. . * ordinary logistic regression . logit pn ar_g1 Iteration 0: log likelihood = -234.95215 Iteration 1: log likelihood = -230.59287 Iteration 2: log likelihood = -230.59173 Iteration 3: log likelihood = -230.59173 Logistic regression Number of obs = 341 LR chi2(1) = 8.72 Prob > chi2 = 0.0031 Log likelihood = -230.59173 Pseudo R2 = 0.0186 ------------------------------------------------------------------------------ pn | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- ar_g1 | .6465241 .2203379 2.93 0.003 .2146697 1.078378 _cons | -.1448309 .1556373 -0.93 0.352 -.4498744 .1602125 ------------------------------------------------------------------------------ . * random effects model . melogit pn ar_g1 || farm: Fitting fixed-effects model: Iteration 0: log likelihood = -230.79889 Iteration 1: log likelihood = -230.59177 Iteration 2: log likelihood = -230.59173 Iteration 3: log likelihood = -230.59173 Refining starting values: Grid node 0: log likelihood = -213.91338 Fitting full model: Iteration 0: log likelihood = -213.91338 Iteration 1: log likelihood = -213.52 Iteration 2: log likelihood = -213.51177 Iteration 3: log likelihood = -213.51176 Mixed-effects logistic regression Number of obs = 341 Group variable: farm Number of groups = 15 Obs per group: min = 14 avg = 22.7 max = 28 Integration method: mvaghermite Integration pts. = 7 Wald chi2(1) = 2.86 Log likelihood = -213.51176 Prob > chi2 = 0.0905 ------------------------------------------------------------------------------ pn | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- ar_g1 | .4369322 .2581451 1.69 0.091 -.0690229 .9428872 _cons | .0196548 .3009262 0.07 0.948 -.5701497 .6094593 -------------+---------------------------------------------------------------- farm | var(_cons)| .8774246 .4325507 .333877 2.305861 ------------------------------------------------------------------------------ LR test vs. logistic model: chibar2(01) = 34.16 Prob >= chibar2 = 0.0000 . . * Reunion Island data, 3-level model for binary outcome fscr . use reu_cfs, clear . melogit fscr i.heifer i.ai lncfs || herd: || cow: Fitting fixed-effects model: Iteration 0: log likelihood = -2011.8253 Iteration 1: log likelihood = -2009.1421 Iteration 2: log likelihood = -2009.1412 Iteration 3: log likelihood = -2009.1412 Refining starting values: Grid node 0: log likelihood = -2028.472 Fitting full model: Iteration 0: log likelihood = -2028.472 (not concave) Iteration 1: log likelihood = -2016.2289 (not concave) Iteration 2: log likelihood = -2008.7841 Iteration 3: log likelihood = -2000.2338 Iteration 4: log likelihood = -1998.6784 Iteration 5: log likelihood = -1998.6736 Iteration 6: log likelihood = -1998.6736 Mixed-effects logistic regression Number of obs = 3,027 Grouping information ------------------------------------------------------------- | No. of Observations per group Group variable | groups Minimum Average Maximum ----------------+-------------------------------------------- herd | 50 13 60.5 226 cow | 1,575 1 1.9 5 ------------------------------------------------------------- Integration method: mvaghermite Integration pts. = 7 Wald chi2(3) = 82.69 Log likelihood = -1998.6736 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ fscr | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- heifer | primiparous | -.0620822 .0974796 -0.64 0.524 -.2531386 .1289743 | ai | ai | -1.023486 .1309939 -7.81 0.000 -1.280229 -.7667427 lncfs | .4971351 .1019534 4.88 0.000 .2973101 .6969601 _cons | -1.515241 .4460314 -3.40 0.001 -2.389447 -.6410359 -------------+---------------------------------------------------------------- herd | var(_cons)| .0732867 .0359684 .0280068 .1917725 -------------+---------------------------------------------------------------- herd>cow | var(_cons)| .2799605 .1219134 .1192414 .6573041 ------------------------------------------------------------------------------ LR test vs. logistic model: chi2(2) = 20.94 Prob > chi2 = 0.0000 Note: LR test is conservative and provided only for reference. .