----------------------------------------------------------------------------------------------------------------------------------------------- name: log: c:\vhm812-data\l13b-intro_cluster_discrete.txt log type: text opened on: 5 Apr 2021, 09:07:17 . . *Random effects binary data . *Pig respiratory disease data . use pig_adg.dta, clear (Growth Performance and Abattoir Surveillance Data on Pigs, Real Data) . * generate dichotomous atrophic rhinits variable . egen ar_g1=cut(ar), at(0, 1.5, 99) icodes . . * 2x2 table analysis . cc pn ar_g1 Proportion | Exposed Unexposed | Total exposed -----------------+------------------------+------------------------ Cases | 109 77 | 186 0.5860 Controls | 66 89 | 155 0.4258 -----------------+------------------------+------------------------ Total | 175 166 | 341 0.5132 | | | Point estimate | [95% Conf. Interval] |------------------------+------------------------ Odds ratio | 1.908894 | 1.21155 3.009556 (exact) Attr. frac. ex. | .4761365 | .1746111 .6677251 (exact) Attr. frac. pop | .2790262 | +------------------------------------------------- chi2(1) = 8.69 Pr>chi2 = 0.0032 . . * ordinary logistic regression . logit pn i.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 | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- 1.ar_g1 | .6465241 .2203379 2.93 0.003 .2146697 1.078378 _cons | -.1448309 .1556373 -0.93 0.352 -.4498744 .1602125 ------------------------------------------------------------------------------ . logit pn i.ar_g1,or 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 | Odds Ratio Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- 1.ar_g1 | 1.908894 .4206017 2.93 0.003 1.239452 2.939909 _cons | .8651685 .1346525 -0.93 0.352 .6377082 1.17376 ------------------------------------------------------------------------------ Note: _cons estimates baseline odds. . . * random effect logistic regression . 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 | Coef. 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 . melogit pn i.ar_g1 || farm:, or 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 | Odds Ratio Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- 1.ar_g1 | 1.547951 .3995959 1.69 0.091 .9333053 2.567383 _cons | 1.019849 .3068993 0.07 0.948 .5654408 1.839437 -------------+---------------------------------------------------------------- farm | var(_cons)| .8774246 .4325507 .333877 2.305861 ------------------------------------------------------------------------------ Note: Estimates are transformed only in the first equation. Note: _cons estimates baseline odds (conditional on zero random effects). LR test vs. logistic model: chibar2(01) = 34.16 Prob >= chibar2 = 0.0000 . . *ICC // approximation to real ICC by pi^2/3 = 3.29 . estat icc Residual intraclass correlation ------------------------------------------------------------------------------ Level | ICC Std. Err. [95% Conf. Interval] -----------------------------+------------------------------------------------ farm | .2105503 .0819422 .0921359 .4120752 ------------------------------------------------------------------------------ . . * PA logistic regression . xtgee pn ar_g1, fam(binomial) link(logit) i(farm) robust Iteration 1: tolerance = .17657314 Iteration 2: tolerance = .0013339 Iteration 3: tolerance = .00001851 Iteration 4: tolerance = 2.251e-07 GEE population-averaged model Number of obs = 341 Group variable: farm Number of groups = 15 Link: logit Obs per group: Family: binomial min = 14 Correlation: exchangeable avg = 22.7 max = 28 Wald chi2(1) = 2.71 Scale parameter: 1 Prob > chi2 = 0.1000 (Std. Err. adjusted for clustering on farm) ------------------------------------------------------------------------------ | Robust pn | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- ar_g1 | .3539583 .2151855 1.64 0.100 -.0677976 .7757142 _cons | .0183926 .2714184 0.07 0.946 -.5135777 .5503629 ------------------------------------------------------------------------------ . . * Random effects models for count data . * open the tb_real dataset . use tb_real, clear . * Poisson model with no random effects . glm reactors i.type i.sex i.age, exp(tar) link(log) fam(poisson) Iteration 0: log likelihood = -338.47585 Iteration 1: log likelihood = -246.77274 Iteration 2: log likelihood = -239.32758 Iteration 3: log likelihood = -238.67153 Iteration 4: log likelihood = -238.662 Iteration 5: log likelihood = -238.662 Generalized linear models Number of obs = 134 Optimization : ML Residual df = 127 Scale parameter = 1 Deviance = 348.3526568 (1/df) Deviance = 2.742934 Pearson = 1105.702579 (1/df) Pearson = 8.70632 Variance function: V(u) = u [Poisson] Link function : g(u) = ln(u) [Log] AIC = 3.666597 Log likelihood = -238.6620031 BIC = -273.673 ------------------------------------------------------------------------------ | OIM reactors | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | .4422157 .2364116 1.87 0.061 -.0211426 .905574 cervid | 1.066241 .2333805 4.57 0.000 .6088238 1.523658 other | .4383852 .6149186 0.71 0.476 -.7668331 1.643603 | sex | male | -.3618833 .1953974 -1.85 0.064 -.7448552 .0210886 | age | 12-24 mo | 2.673416 .721739 3.70 0.000 1.258834 4.087999 >24 mont | 2.601192 .7136354 3.64 0.000 1.202493 3.999892 | _cons | -11.68987 .7397843 -15.80 0.000 -13.13982 -10.23992 ln(tar) | 1 (exposure) ------------------------------------------------------------------------------ . estimates store pois . nbreg reactors i.type i.sex i.age, exp(tar) Fitting Poisson model: Iteration 0: log likelihood = -244.95075 Iteration 1: log likelihood = -239.03967 Iteration 2: log likelihood = -238.66707 Iteration 3: log likelihood = -238.662 Iteration 4: log likelihood = -238.662 Fitting constant-only model: Iteration 0: log likelihood = -169.94856 Iteration 1: log likelihood = -163.31099 Iteration 2: log likelihood = -163.12383 Iteration 3: log likelihood = -163.12339 Iteration 4: log likelihood = -163.12339 Fitting full model: Iteration 0: log likelihood = -163.12339 Iteration 1: log likelihood = -158.08548 Iteration 2: log likelihood = -157.73682 Iteration 3: log likelihood = -157.73596 Iteration 4: log likelihood = -157.73596 Negative binomial regression Number of obs = 134 LR chi2(6) = 10.77 Dispersion = mean Prob > chi2 = 0.0956 Log likelihood = -157.73596 Pseudo R2 = 0.0330 ------------------------------------------------------------------------------ reactors | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | .6046115 .6749332 0.90 0.370 -.7182333 1.927456 cervid | .6657202 .6838382 0.97 0.330 -.6745781 2.006018 other | .8000348 1.118942 0.71 0.475 -1.393052 2.993121 | sex | male | -.0574791 .4050815 -0.14 0.887 -.8514242 .736466 | age | 12-24 mo | 2.253036 .9033002 2.49 0.013 .4826001 4.023472 >24 mont | 2.480955 .881572 2.81 0.005 .7531053 4.208804 | _cons | -11.18145 1.060641 -10.54 0.000 -13.26026 -9.102626 ln(tar) | 1 (exposure) -------------+---------------------------------------------------------------- /lnalpha | .5541007 .2534715 .0573057 1.050896 -------------+---------------------------------------------------------------- alpha | 1.740375 .4411355 1.05898 2.860212 ------------------------------------------------------------------------------ LR test of alpha=0: chibar2(01) = 161.85 Prob >= chibar2 = 0.000 . . * negative binomial model no random effects . glm reactors i.type i.sex i.age, exp(tar) link(log) fam(nbin ml) Iteration 0: log likelihood = -168.6314 Iteration 1: log likelihood = -157.89374 Iteration 2: log likelihood = -157.73611 Iteration 3: log likelihood = -157.73596 Iteration 4: log likelihood = -157.73596 Generalized linear models Number of obs = 134 Optimization : ML Residual df = 127 Scale parameter = 1 Deviance = 99.35977885 (1/df) Deviance = .7823605 Pearson = 374.8585955 (1/df) Pearson = 2.951642 Variance function: V(u) = u+(1.7404)u^2 [Neg. Binomial] Link function : g(u) = ln(u) [Log] AIC = 2.458746 Log likelihood = -157.735955 BIC = -522.6659 ------------------------------------------------------------------------------ | OIM reactors | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | .6046115 .674193 0.90 0.370 -.7167824 1.926005 cervid | .6657202 .6821121 0.98 0.329 -.671195 2.002635 other | .8000348 1.11801 0.72 0.474 -1.391225 2.991295 | sex | male | -.0574791 .4044774 -0.14 0.887 -.8502403 .735282 | age | 12-24 mo | 2.253036 .8976472 2.51 0.012 .4936798 4.012392 >24 mont | 2.480955 .8779162 2.83 0.005 .7602705 4.201639 | _cons | -11.18145 1.052762 -10.62 0.000 -13.24482 -9.118069 ln(tar) | 1 (exposure) ------------------------------------------------------------------------------ Note: Negative binomial parameter estimated via ML and treated as fixed once estimated. . estimates store nb . * Pearson dispersion parameter still large (2.95) . . * Poisson model with normal distributed random effects . mepoisson reactors i.type i.sex i.age, exp(tar) || farm_id: Fitting fixed-effects model: Iteration 0: log likelihood = -602.90391 Iteration 1: log likelihood = -247.78547 Iteration 2: log likelihood = -239.70423 Iteration 3: log likelihood = -238.67796 Iteration 4: log likelihood = -238.66201 Iteration 5: log likelihood = -238.662 Refining starting values: Grid node 0: log likelihood = -150.38337 Fitting full model: Iteration 0: log likelihood = -150.38337 Iteration 1: log likelihood = -144.07026 Iteration 2: log likelihood = -143.57132 Iteration 3: log likelihood = -143.55978 Iteration 4: log likelihood = -143.55988 Iteration 5: log likelihood = -143.55988 Mixed-effects Poisson regression Number of obs = 134 Group variable: farm_id Number of groups = 30 Obs per group: min = 1 avg = 4.5 max = 13 Integration method: mvaghermite Integration pts. = 7 Wald chi2(6) = 18.42 Log likelihood = -143.55988 Prob > chi2 = 0.0053 ------------------------------------------------------------------------------ reactors | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | -.3941614 .3320059 -1.19 0.235 -1.044881 .2565583 cervid | -.2394465 .4847969 -0.49 0.621 -1.189631 .7107379 other | -.1058558 .7985319 -0.13 0.895 -1.67095 1.459238 | sex | male | -.3392105 .2083118 -1.63 0.103 -.7474941 .0690731 | age | 12-24 mo | 2.716639 .7470885 3.64 0.000 1.252372 4.180906 >24 mont | 2.466631 .7256478 3.40 0.001 1.044388 3.888875 | _cons | -11.05318 .8283607 -13.34 0.000 -12.67674 -9.429625 ln(tar) | 1 (exposure) -------------+---------------------------------------------------------------- farm_id | var(_cons)| 1.697739 .6022648 .8470544 3.402754 ------------------------------------------------------------------------------ LR test vs. Poisson model: chibar2(01) = 190.20 Prob >= chibar2 = 0.0000 . mepoisson reactors i.type i.sex i.age, exp(tar) irr || farm_id: Fitting fixed-effects model: Iteration 0: log likelihood = -602.90391 Iteration 1: log likelihood = -247.78547 Iteration 2: log likelihood = -239.70423 Iteration 3: log likelihood = -238.67796 Iteration 4: log likelihood = -238.66201 Iteration 5: log likelihood = -238.662 Refining starting values: Grid node 0: log likelihood = -150.38337 Fitting full model: Iteration 0: log likelihood = -150.38337 Iteration 1: log likelihood = -144.07026 Iteration 2: log likelihood = -143.57132 Iteration 3: log likelihood = -143.55978 Iteration 4: log likelihood = -143.55988 Iteration 5: log likelihood = -143.55988 Mixed-effects Poisson regression Number of obs = 134 Group variable: farm_id Number of groups = 30 Obs per group: min = 1 avg = 4.5 max = 13 Integration method: mvaghermite Integration pts. = 7 Wald chi2(6) = 18.42 Log likelihood = -143.55988 Prob > chi2 = 0.0053 ------------------------------------------------------------------------------ reactors | IRR Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | .6742453 .2238534 -1.19 0.235 .3517337 1.292474 cervid | .7870634 .3815659 -0.49 0.621 .3043336 2.035493 other | .8995544 .7183229 -0.13 0.895 .1880684 4.30268 | sex | male | .7123325 .1483872 -1.63 0.103 .4735517 1.071514 | age | 12-24 mo | 15.12939 11.30299 3.64 0.000 3.498633 65.42507 >24 mont | 11.78269 8.550082 3.40 0.001 2.841659 48.85589 | _cons | .0000158 .0000131 -13.34 0.000 3.12e-06 .0000803 ln(tar) | 1 (exposure) -------------+---------------------------------------------------------------- farm_id | var(_cons)| 1.697739 .6022648 .8470544 3.402754 ------------------------------------------------------------------------------ Note: Estimates are transformed only in the first equation. Note: _cons estimates baseline incidence rate (conditional on zero random effects). LR test vs. Poisson model: chibar2(01) = 190.20 Prob >= chibar2 = 0.0000 . **same with meglm . meglm reactors i.type i.sex i.age, exp(tar) || farm_id:, fam(poisson) Fitting fixed-effects model: Iteration 0: log likelihood = -602.90391 Iteration 1: log likelihood = -247.78547 Iteration 2: log likelihood = -239.70423 Iteration 3: log likelihood = -238.67796 Iteration 4: log likelihood = -238.66201 Iteration 5: log likelihood = -238.662 Refining starting values: Grid node 0: log likelihood = -150.38337 Fitting full model: Iteration 0: log likelihood = -150.38337 Iteration 1: log likelihood = -144.07026 Iteration 2: log likelihood = -143.57132 Iteration 3: log likelihood = -143.55978 Iteration 4: log likelihood = -143.55988 Iteration 5: log likelihood = -143.55988 Mixed-effects GLM Number of obs = 134 Family: Poisson Link: log Group variable: farm_id Number of groups = 30 Obs per group: min = 1 avg = 4.5 max = 13 Integration method: mvaghermite Integration pts. = 7 Wald chi2(6) = 18.42 Log likelihood = -143.55988 Prob > chi2 = 0.0053 ------------------------------------------------------------------------------ reactors | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | -.3941614 .3320059 -1.19 0.235 -1.044881 .2565583 cervid | -.2394465 .4847969 -0.49 0.621 -1.189631 .7107379 other | -.1058558 .7985319 -0.13 0.895 -1.67095 1.459238 | sex | male | -.3392105 .2083118 -1.63 0.103 -.7474941 .0690731 | age | 12-24 mo | 2.716639 .7470885 3.64 0.000 1.252372 4.180906 >24 mont | 2.466631 .7256478 3.40 0.001 1.044388 3.888875 | _cons | -11.05318 .8283607 -13.34 0.000 -12.67674 -9.429625 ln(tar) | 1 (exposure) -------------+---------------------------------------------------------------- farm_id | var(_cons)| 1.697739 .6022648 .8470544 3.402754 ------------------------------------------------------------------------------ LR test vs. Poisson model: chibar2(01) = 190.20 Prob >= chibar2 = 0.0000 . estimates store pois_norm . . * Negative binomial with Normal random effects . meglm reactors i.type i.sex i.age, exp(tar) || farm_id:, fam(nbin) Fitting fixed-effects model: Iteration 0: log likelihood = -175.65303 Iteration 1: log likelihood = -158.06607 Iteration 2: log likelihood = -157.73916 Iteration 3: log likelihood = -157.73596 Iteration 4: log likelihood = -157.73596 Refining starting values: Grid node 0: log likelihood = -153.88154 Fitting full model: Iteration 0: log likelihood = -153.88154 Iteration 1: log likelihood = -147.66002 Iteration 2: log likelihood = -144.59242 Iteration 3: log likelihood = -143.77217 Iteration 4: log likelihood = -143.59593 Iteration 5: log likelihood = -143.5688 Iteration 6: log likelihood = -143.56194 Iteration 7: log likelihood = -143.56037 Iteration 8: log likelihood = -143.55999 Iteration 9: log likelihood = -143.55991 Iteration 10: log likelihood = -143.55989 Iteration 11: log likelihood = -143.55988 Mixed-effects GLM Number of obs = 134 Family: negative binomial Link: log Overdispersion: mean Group variable: farm_id Number of groups = 30 Obs per group: min = 1 avg = 4.5 max = 13 Integration method: mvaghermite Integration pts. = 7 Wald chi2(6) = 18.42 Log likelihood = -143.55988 Prob > chi2 = 0.0053 ------------------------------------------------------------------------------ reactors | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- type | beef | -.3942466 .3320071 -1.19 0.235 -1.044968 .2564753 cervid | -.2395109 .4847974 -0.49 0.621 -1.189696 .7106746 other | -.1059246 .7985326 -0.13 0.894 -1.67102 1.459171 | sex | male | -.3392058 .2083119 -1.63 0.103 -.7474896 .069078 | age | 12-24 mo | 2.716643 .7470858 3.64 0.000 1.252382 4.180905 >24 mont | 2.466621 .7256448 3.40 0.001 1.044384 3.888859 | _cons | -11.05311 .8283585 -13.34 0.000 -12.67667 -9.429561 ln(tar) | 1 (exposure) -------------+---------------------------------------------------------------- /lnalpha | -16.02436 946.7293 -1871.58 1839.531 -------------+---------------------------------------------------------------- farm_id | var(_cons)| 1.697762 .6022739 .8470651 3.402805 ------------------------------------------------------------------------------ LR test vs. nbinomial model: chibar2(01) = 28.35 Prob >= chibar2 = 0.0000 . estimates store nb_norm . * the overdispersion lnalpha is non signfincant which suggests that the . * overdispersion is due to clustering since this model is exactly . * the same as the Poisson with random effects . . estimates store pois_norm . estimates table pois nb pois_norm nb_norm, se(%4.3f) b(%4.3f) ------------------------------------------------------ Variable | pois nb pois_~m nb_norm -------------+---------------------------------------- reactors | type | beef | 0.442 0.605 -0.394 -0.394 | 0.236 0.674 0.332 0.332 cervid | 1.066 0.666 -0.240 -0.240 | 0.233 0.682 0.485 0.485 other | 0.438 0.800 -0.106 -0.106 | 0.615 1.118 0.799 0.799 | sex | male | -0.362 -0.057 -0.339 -0.339 | 0.195 0.404 0.208 0.208 | age | 12-24 mo | 2.673 2.253 2.717 2.717 | 0.722 0.898 0.747 0.747 >24 mont | 2.601 2.481 2.467 2.467 | 0.714 0.878 0.726 0.726 | _cons | -11.690 -11.181 -11.053 -11.053 | 0.740 1.053 0.828 0.828 -------------+---------------------------------------- /lnalpha | -16.024 -16.024 | 946.729 946.729 var( | _cons[far~d])| 1.698 1.698 | 0.602 0.602 ------------------------------------------------------ legend: b/se . estimates stats pois nb pois_norm nb_norm Akaike's information criterion and Bayesian information criterion ----------------------------------------------------------------------------- Model | N ll(null) ll(model) df AIC BIC -------------+--------------------------------------------------------------- pois | 134 . -238.662 7 491.324 511.6089 nb | 134 . -157.736 7 329.4719 349.7568 pois_norm | 134 . -143.5599 9 305.1198 331.2003 nb_norm | 134 . -143.5599 9 305.1198 331.2003 ----------------------------------------------------------------------------- Note: BIC uses N = number of observations. See [R] BIC note. . . log close name: log: c:\vhm812-data\l13b-intro_cluster_discrete.txt log type: text closed on: 5 Apr 2021, 09:07:19 -----------------------------------------------------------------------------------------------------------------------------------------------