. do "C:\Users\DEFAUL~1.SID\AppData\Local\Temp\STD00000000.tmp"

. * do-file for part of lecture 11b of VHM 802/812, Winter 2015
. version 13

. set more off

. cd "h:\vhm\vhm802\data_stata"
h:\vhm\vhm802\data_stata

.
. * unadjusted and mixed model results (for reference)
. use simcont_clustherd.dta, clear

. regress milk X

      Source |       SS       df       MS              Number of obs =   11626
-------------+------------------------------           F(  1, 11624) =  317.72
       Model |  36598.5078     1  36598.5078           Prob > F      =  0.0000
    Residual |     1338999 11624  115.192618           R-squared     =  0.0266
-------------+------------------------------           Adj R-squared =  0.0265
       Total |  1375597.51 11625  118.330968           Root MSE      =  10.733

-------------+----------------------------------------------------------------
         milk |      Coef.    Std. Err.       t    P>|t|      [95% Con f. 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

Computing standard errors:

Mixed-effects REML regression                   Number of obs      =     11626
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 |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    3.796004    1.495943      2.54    0.011       .864009     6.727999
       _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.91499   7.998621      41.27713    73.05877
-----------------------------+------------------------------------------------
               var(Residual) |   64.20087   .8457062      62.56453    65.88001
-----------------------------+------------------------------------------------
LR test vs. linear regression: chibar2(01) =  6374.40 Prob >= chibar2 =  0.0000

. use simcont_clustcow.dta, clear

. regress milk X

      Source |       SS       df       MS              Number of obs =   11626
-------------+------------------------------           F(  1, 11624) =  624.90
       Model |  72138.7619     1  72138.7619           Prob > F      =  0.0000
    Residual |  1341880.62 11624  115.440522           R-squared     =  0.0510
-------------+------------------------------           Adj R-squared =  0.0509
       Total |  1414019.39 11625  121.636076           Root MSE      =  10.744

-------------+----------------------------------------------------------------
         milk |      Coef.    Std. Err.       t    P>|t|      [95% Con f. 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

Computing standard errors:

Mixed-effects REML regression                   Number of obs      =     11626
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 |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    4.968194    .1492174     33.30    0.000      4.675733     5.260655
       _cons |    30.64647    .7281276     42.09    0.000      29.21936     32.07357
-------------+----------------------------------------------------------------

-----------------------------+------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
herd: Identity               |
                  var(_cons) |   51.41189   7.459591      38.68644    68.32323
-----------------------------+------------------------------------------------
               var(Residual) |   64.71069   .8524578      63.06129    66.40324
-----------------------------+------------------------------------------------
LR test vs. linear regression: 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   =      11626
                                                  LR chi2(1)      =     159.15
                                                  Prob > chi2     =     0.0000
Log likelihood = -6814.7785                       Pseudo R2       =     0.0115

-------------+----------------------------------------------------------------
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5287317    .0423191     12.49    0.000      .4457877     .6116757
       _cons |   -1.241768    .0325699    -38.13    0.000     -1.305604    -1.177932
-------------+----------------------------------------------------------------

. meqrlogit Y X || herd:

Refining starting values:

Iteration 0:   log likelihood = -6065.2694
Iteration 1:   log likelihood = -6065.1778
Iteration 2:   log likelihood = -6065.0871

Performing gradient-based optimization:

Iteration 0:   log likelihood = -6065.0871
Iteration 1:   log likelihood = -6065.0867

Mixed-effects logistic regression               Number of obs      =     11626
Group variable: herd                            Number of groups   =       100

                                                Obs per group: min =        20
                                                               avg =     116.3
                                                               max =       311

Integration points =   7                        Wald chi2(1)       =      9.25
Log likelihood = -6065.0867                     Prob > chi2        =    0.0024

-------------+----------------------------------------------------------------
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |     .619974    .2038516      3.04    0.002      .2204322     1.019516
       _cons |   -1.305417    .1455518     -8.97    0.000     -1.590693    -1.020141
-------------+----------------------------------------------------------------

-----------------------------+------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
herd: Identity               |
                  var(_cons) |   .9416185   .1492971      .6901029    1.284802
-----------------------------+------------------------------------------------
LR test vs. logistic regression: chibar2(01) =  1499.38 Prob>=chibar2 =  0.0000

. di .619974/(sqrt(1+.346*0.9416185)) /* PA estimate, for later use */
.53843647

. 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   =      11626
                                                  LR chi2(1)      =     198.54
                                                  Prob > chi2     =     0.0000
Log likelihood = -6811.0741                       Pseudo R2       =     0.0144

-------------+----------------------------------------------------------------
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5863084    .0419748     13.97    0.000      .5040393     .6685775
       _cons |    -1.25032    .0316033    -39.56    0.000     -1.312261    -1.188379
-------------+----------------------------------------------------------------

. meqrlogit Y X || herd:

Refining starting values:

Iteration 0:   log likelihood = -5999.0538
Iteration 1:   log likelihood = -5995.9727
Iteration 2:   log likelihood = -5995.9698

Performing gradient-based optimization:

Iteration 0:   log likelihood = -5995.9698
Iteration 1:   log likelihood = -5995.9698

Mixed-effects logistic regression               Number of obs      =     11626
Group variable: herd                            Number of groups   =       100

                                                Obs per group: min =        20
                                                               avg =     116.3
                                                               max =       311

Integration points =   7                        Wald chi2(1)       =    229.28
Log likelihood = -5995.9698                     Prob > chi2        =    0.0000

-------------+----------------------------------------------------------------
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |     .697479     .046063     15.14    0.000      .6071972     .7877609
       _cons |   -1.361173    .1112103    -12.24    0.000     -1.579141    -1.143204
-------------+----------------------------------------------------------------

-----------------------------+------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
herd: Identity               |
                  var(_cons) |   1.068172   .1682549      .7844459     1.45452
-----------------------------+------------------------------------------------
LR test vs. logistic regression: chibar2(01) =  1630.21 Prob>=chibar2 =  0.0000

. di .697479/(sqrt(1+.346*1.068172)) /* PA estimate, for later use */
.59598625

.
. * fixed effects modelling
. use simcont_clustcow.dta, clear

. regress milk X i.herd

      Source |       SS       df       MS              Number of obs =   11626
-------------+------------------------------           F(100, 11525) =  103.27
       Model |   668233.58   100   6682.3358           Prob > F      =  0.0000
    Residual |  745785.805 11525  64.7102651           R-squared     =  0.4726
-------------+------------------------------           Adj R-squared =  0.4680
       Total |  1414019.39 11625  121.636076           Root MSE      =  8.0443

-------------+----------------------------------------------------------------
         milk |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    4.967911     .149217     33.29    0.000       4.67542     5.260401
             |
        herd |
          2  |    21.11033    2.513358      8.40    0.000      16.18373     26.03694
          3  |    7.236815    2.413283      3.00    0.003      2.506371     11.96726
          4  |    10.38913    2.392565      4.34    0.000      5.699296     15.07897
          5  |    1.856204    2.373225      0.78    0.434      -2.79572     6.508128
          6  |     22.3249    2.355121      9.48    0.000      17.70846     26.94134
          7  |    19.67363    2.355121      8.35    0.000      15.05719     24.29006
          8  |    20.58936    2.338142      8.81    0.000       16.0062     25.17252
          9  |    4.836125    2.338142      2.07    0.039      .2529695     9.419281
         10  |   -4.982617    2.338142     -2.13    0.033     -9.565773    -.3994617
         11  |    10.95127     2.32218      4.72    0.000      6.399398     15.50313
         12  |     3.61817     2.32218      1.56    0.119     -.9336981     8.170038
         13  |   -.9090294     2.32218     -0.39    0.695     -5.460898     3.642839
         14  |    6.437472    2.307151      2.79    0.005      1.915064     10.95988
         15  |    1.814185    2.292969      0.79    0.429     -2.680425     6.308795
         16  |    3.944011    2.254858      1.75    0.080     -.4758925     8.363915
         17  |     11.9281    2.254858      5.29    0.000      7.508197       16.348
         18  |     3.93309    2.254858      1.74    0.081     -.4868133     8.352994
         19  |    20.00963    2.232587      8.96    0.000      15.63338     24.38588
         20  |     14.6525    2.232587      6.56    0.000      10.27625     19.02875
         21  |     10.9608    2.222254      4.93    0.000      6.604801     15.31679
         22  |    7.190073    2.222254      3.24    0.001      2.834076     11.54607
         23  |    5.612255    2.212409      2.54    0.011      1.275557     9.948953
         24  |    3.203415    2.203014      1.45    0.146     -1.114866     7.521696
         25  |    2.956442    2.194041      1.35    0.178     -1.344251     7.257135
         26  |    11.80699     2.18546      5.40    0.000      7.523121     16.09086
         27  |    .0918443     2.18546      0.04    0.966     -4.192028     4.375716
         28  |   -1.178179     2.18546     -0.54    0.590     -5.462051     3.105693
         29  |    4.578893    2.161833      2.12    0.034      .3413331     8.816453
         30  |    6.631271    2.161833      3.07    0.002      2.393711     10.86883
         31  |    5.551421     2.15459      2.58    0.010      1.328058     9.774783
         32  |   -2.206988     2.15459     -1.02    0.306     -6.430351     2.016374
         33  |    5.821662     2.15459      2.70    0.007      1.598299     10.04502
         34  |   -1.102839    2.147634     -0.51    0.608     -5.312565     3.106888
         35  |    8.300868    2.147634      3.87    0.000      4.091142     12.51059
         36  |    9.826285    2.147634      4.58    0.000      5.616558     14.03601
         37  |    18.43057    2.147634      8.58    0.000      14.22084      22.6403
         38  |   -.1458487    2.147634     -0.07    0.946     -4.355575     4.063878
         39  |    16.57678    2.140945      7.74    0.000      12.38016     20.77339
         40  |    7.425361    2.140945      3.47    0.001      3.228746     11.62198
         41  |    14.93193     2.13451      7.00    0.000      10.74793     19.11593
         42  |    11.21274     2.13451      5.25    0.000      7.028737     15.39674
         43  |    19.93934    2.128314      9.37    0.000      15.76748      24.1112
         44  |    13.77498    2.122344      6.49    0.000      9.614824     17.93513
         45  |    1.333881    2.122344      0.63    0.530     -2.826274     5.494036
         46  |     11.4609    2.116587      5.41    0.000      7.312026     15.60977
         47  |    5.715232    2.111034      2.71    0.007      1.577247     9.853217
         48  |    9.789586    2.111034      4.64    0.000      5.651601     13.92757
         49  |    12.30603    2.100492      5.86    0.000      8.188707     16.42335
         50  |    2.699361    2.100492      1.29    0.199     -1.417959     6.816682
         51  |    7.424841    1.986609      3.74    0.000      3.530749     11.31893
         52  |    8.593715    1.959664      4.39    0.000       4.75244     12.43499
         53  |    2.399339    1.958235      1.23    0.221     -1.439133     6.237812
         54  |    .0153907    1.944042      0.01    0.994     -3.795261     3.826042
         55  |    9.623738    1.944042      4.95    0.000      5.813087     13.43439
         56  |    16.41427    1.936256      8.48    0.000      12.61887     20.20966
         57  |    1.699915    1.931187      0.88    0.379      -2.08554     5.485371
         58  |    14.31867    1.922948      7.45    0.000      10.54936     18.08797
         59  |    7.084435    1.922096      3.69    0.000      3.316801     10.85207
         60  |    1.817823    1.919607      0.95    0.344     -1.944933     5.580578
         61  |     9.88462    1.915675      5.16    0.000      6.129572     13.63967
         62  |    6.613898    1.911991      3.46    0.001      2.866071     10.36173
         63  |    12.70136    1.911991      6.64    0.000      8.953534     16.44919
         64  |    8.013387    1.910582      4.19    0.000      4.268322     11.75845
         65  |   -13.06771     1.90989     -6.84    0.000     -16.81142    -9.323999
         66  |    7.752948    1.909207      4.06    0.000      4.010578     11.49532
         67  |    -5.62285    1.902812     -2.96    0.003     -9.352685    -1.893016
         68  |     7.35862    1.900457      3.87    0.000      3.633402     11.08384
         69  |    10.45819    1.899319      5.51    0.000      6.735205     14.18118
         70  |    9.985111    1.897659      5.26    0.000      6.265378     13.70484
         71  |    8.920595    1.897117      4.70    0.000      5.201923     12.63927
         72  |    12.01857    1.896582      6.34    0.000      8.300952      15.7362
         73  |    7.489065    1.896052      3.95    0.000       3.77248     11.20565
         74  |    20.08057    1.894498     10.60    0.000      16.36703     23.79411
         75  |    -3.11481    1.892014     -1.65    0.100     -6.823479     .5938582
         76  |   -11.08463    1.888746     -5.87    0.000      -14.7869     -7.38237
         77  |    14.00118    1.886979      7.42    0.000      10.30238     17.69998
         78  |    -8.43247    1.883247     -4.48    0.000     -12.12395    -4.740986
         79  |   -3.531187    1.882852     -1.88    0.061     -7.221896     .1595233
         80  |   -.1802846     1.88246     -0.10    0.924     -3.870227     3.509658
         81  |    8.714384    1.882073      4.63    0.000      5.025202     12.40357
         82  |   -3.744767    1.880556     -1.99    0.046     -7.430977    -.0585577
         83  |    5.663561    1.880556      3.01    0.003      1.977352      9.34977
         84  |    .1282536    1.880556      0.07    0.946     -3.557956     3.814463
         85  |     5.25524    1.879818      2.80    0.005      1.570477     8.940003
         86  |     .601706    1.879455      0.32    0.749     -3.082344     4.285756
         87  |    3.362022    1.879094      1.79    0.074     -.3213211     7.045365
         88  |    1.979885    1.879094      1.05    0.292     -1.703458     5.663228
         89  |    3.037238    1.876321      1.62    0.106     -.6406698     6.715146
         90  |    6.268529    1.875005      3.34    0.001        2.5932     9.943857
         91  |    11.79651    1.871607      6.30    0.000      8.127848     15.46518
         92  |    2.798214    1.870732      1.50    0.135     -.8687376     6.465166
         93  |    8.956046    1.870445      4.79    0.000      5.289657     12.62244
         94  |    -9.05919     1.86877     -4.85    0.000      -12.7223    -5.396084
         95  |    .1418222     1.86877      0.08    0.940     -3.521284     3.804928
         96  |    11.83867    1.867171      6.34    0.000        8.1787     15.49864
         97  |    .5145903    1.862333      0.28    0.782     -3.135899     4.165079
         98  |    .2560854    1.861666      0.14    0.891     -3.393095     3.905266
         99  |    14.01656     1.85642      7.55    0.000      10.37766     17.65546
        100  |   -.3351337     1.85569     -0.18    0.857     -3.972601     3.302334
             |
       _cons |    24.32444      1.8003     13.51    0.000      20.79555     27.85333
-------------+----------------------------------------------------------------

. use simbin_clustcow.dta, clear

. logit Y X i.herd, asis

Iteration 0:   log likelihood = -6910.3442
Iteration 1:   log likelihood = -5871.9702
Iteration 2:   log likelihood = -5810.5761
Iteration 3:   log likelihood =  -5808.037
Iteration 4:   log likelihood =  -5807.961
Iteration 5:   log likelihood = -5807.9609

Logistic regression                               Number of obs   =      11626
                                                  LR chi2(100)    =    2204.77
                                                  Prob > chi2     =     0.0000
Log likelihood = -5807.9609                       Pseudo R2       =     0.1595

-------------+----------------------------------------------------------------
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .7037681    .0462786     15.21    0.000      .6130636     .7944726
             |
        herd |
          2  |    1.464052    .7736212      1.89    0.058     -.0522174     2.980322
          3  |    1.170623    .7593863      1.54    0.123     -.3177464     2.658993
          4  |    1.618762    .7468126      2.17    0.030      .1550363     3.082488
          5  |   -.3622129    .8817161     -0.41    0.681     -2.090345     1.365919
          6  |    2.547188    .7525714      3.38    0.001      1.072175     4.022201
          7  |    2.907366     .770109      3.78    0.000       1.39798     4.416752
          8  |    2.272912    .7409149      3.07    0.002      .8207457     3.725079
          9  |    .3856383    .7825535      0.49    0.622     -1.148138     1.919415
         10  |   -1.622873    1.199061     -1.35    0.176     -3.972989     .7272438
         11  |    1.915425    .7323036      2.62    0.009       .480136     3.350713
         12  |    .1271594    .8013326      0.16    0.874     -1.443424     1.697742
         13  |    .1271594    .8013326      0.16    0.874     -1.443424     1.697742
         14  |    .5002686    .7658501      0.65    0.514      -1.00077     2.001307
         15  |   -.5407558    .8772847     -0.62    0.538     -2.260202     1.178691
         16  |   -.3271922    .8268742     -0.40    0.692     -1.947836     1.293451
         17  |    1.349539    .7216192      1.87    0.061     -.0648087     2.763887
         18  |    .3443184     .761653      0.45    0.651     -1.148494     1.837131
         19  |    2.523936    .7244937      3.48    0.000      1.103955     3.943918
         20  |    1.600434    .7143851      2.24    0.025      .2002652     3.000603
         21  |    .8535915    .7275448      1.17    0.241     -.5723701     2.279553
         22  |    1.559832     .712662      2.19    0.029      .1630404     2.956624
         23  |    .2093641    .7584782      0.28    0.783     -1.277226     1.695954
         24  |    .3543194    .7465973      0.47    0.635     -1.108984     1.817623
         25  |   -.2512472    .7929482     -0.32    0.751     -1.805397     1.302903
         26  |     1.38016    .7070654      1.95    0.051     -.0056627     2.765983
         27  |   -.2705637    .7924538     -0.34    0.733     -1.823745     1.282617
         28  |   -.2705637    .7924538     -0.34    0.733     -1.823745     1.282617
         29  |    .9513724    .7103156      1.34    0.180     -.4408206     2.343565
         30  |    1.351809    .7024931      1.92    0.054      -.025052      2.72867
         31  |    .5879134    .7212949      0.82    0.415     -.8257986     2.001625
         32  |   -.6242464     .821278     -0.76    0.447     -2.233922     .9854288
         33  |    1.322313    .7015163      1.88    0.059     -.0526338      2.69726
         34  |   -.4061465     .790207     -0.51    0.607     -1.954924     1.142631
         35  |    1.370376    .6993817      1.96    0.050     -.0003874     2.741138
         36  |    1.088843    .7039312      1.55    0.122     -.2908367     2.468523
         37  |    2.451419    .7039883      3.48    0.000      1.071627      3.83121
         38  |    .2917802    .7335553      0.40    0.691     -1.145962     1.729522
         39  |    2.396693    .7013308      3.42    0.001       1.02211     3.771276
         40  |    1.063799    .7031859      1.51    0.130     -.3144201     2.442018
         41  |     2.71341    .7076245      3.83    0.000      1.326492     4.100329
         42  |    1.473795    .6957181      2.12    0.034      .1102122     2.837377
         43  |    2.649916    .7042126      3.76    0.000      1.269685     4.030147
         44  |    1.892143    .6923309      2.73    0.006      .5351993     3.249086
         45  |    .3222867    .7244805      0.44    0.656     -1.097669     1.742242
         46  |     2.01674    .6917662      2.92    0.004      .6609029     3.372577
         47  |    .6155321    .7087622      0.87    0.385     -.7736163     2.004681
         48  |    1.254889     .693386      1.81    0.070     -.1041221     2.613901
         49  |    1.507312    .6887245      2.19    0.029      .1574369     2.857187
         50  |    .2279885    .7228601      0.32    0.752     -1.188791     1.644768
         51  |    1.093577    .6696441      1.63    0.102     -.2189013     2.406056
         52  |     1.75509    .6608457      2.66    0.008      .4598563     3.050324
         53  |   -.3494651    .7021245     -0.50    0.619     -1.725604     1.026674
         54  |    -.462057    .7012843     -0.66    0.510     -1.836549     .9124349
         55  |      1.2621    .6594879      1.91    0.056     -.0304724     2.554673
         56  |    3.018589    .6663273      4.53    0.000      1.712611     4.324567
         57  |   -.3178017    .6887097     -0.46    0.644     -1.667648     1.032044
         58  |    1.836669    .6539541      2.81    0.005      .5549426     3.118396
         59  |    .5157728    .6628867      0.78    0.437     -.7834612     1.815007
         60  |    .1771838    .6688583      0.26    0.791     -1.133754     1.488122
         61  |    1.650583    .6526415      2.53    0.011      .3714294     2.929737
         62  |    .4430613    .6614298      0.67    0.503     -.8533174      1.73944
         63  |    1.481111    .6522875      2.27    0.023      .2026513     2.759571
         64  |    .9428472    .6549982      1.44    0.150     -.3409256      2.22662
         65  |   -1.203757    .7284121     -1.65    0.098     -2.631419     .2239041
         66  |    1.462062    .6518104      2.24    0.025      .1845372     2.739587
         67  |    -1.03609    .7107443     -1.46    0.145     -2.429124      .356943
         68  |    .5836485    .6564051      0.89    0.374     -.7028818     1.870179
         69  |     1.63526    .6495967      2.52    0.012      .3620738     2.908446
         70  |    1.009096     .651701      1.55    0.122     -.2682148     2.286406
         71  |    1.106806    .6509832      1.70    0.089     -.1690977     2.382709
         72  |    1.764194    .6490061      2.72    0.007      .4921657     3.036223
         73  |    1.454343    .6493218      2.24    0.025      .1816959      2.72699
         74  |    2.645304    .6518276      4.06    0.000      1.367745     3.922863
         75  |   -.4566702    .6762469     -0.68    0.499      -1.78209     .8687494
         76  |   -2.441802    .8586407     -2.84    0.004     -4.124707    -.7588972
         77  |    1.703338    .6472359      2.63    0.008      .4347787     2.971897
         78  |    -1.37405     .717516     -1.92    0.055     -2.780356     .0322552
         79  |   -.6994571    .6802396     -1.03    0.304     -2.032702      .633788
         80  |   -.4161735    .6704249     -0.62    0.535     -1.730182     .8978352
         81  |    1.267766    .6472003      1.96    0.050     -.0007232     2.536255
         82  |   -1.106402    .6978069     -1.59    0.113     -2.474079     .2612743
         83  |    .8230967    .6493033      1.27    0.205     -.4495144     2.095708
         84  |   -.1682734    .6631341     -0.25    0.800     -1.467992     1.131446
         85  |    .8790128    .6487414      1.35    0.175     -.3924968     2.150523
         86  |    .3069123    .6543537      0.47    0.639     -.9755974     1.589422
         87  |   -.1488085    .6621338     -0.22    0.822     -1.446567      1.14895
         88  |    .3585251    .6535491      0.55    0.583     -.9224076     1.639458
         89  |    .3130271    .6533424      0.48    0.632     -.9675004     1.593555
         90  |    .0773685    .6563783      0.12    0.906     -1.209109     1.363846
         91  |    1.641164    .6444201      2.55    0.011      .3781236     2.904204
         92  |    .3000771    .6518837      0.46    0.645     -.9775916     1.577746
         93  |    1.006892    .6459791      1.56    0.119     -.2592041     2.272988
         94  |   -.9668609    .6821546     -1.42    0.156     -2.303859     .3701376
         95  |   -.3907563    .6633882     -0.59    0.556     -1.690973     .9094607
         96  |    1.841699    .6435616      2.86    0.004      .5803409     3.103056
         97  |    .2936752    .6495229      0.45    0.651     -.9793664     1.566717
         98  |    .2063121    .6502998      0.32    0.751     -1.068252     1.480876
         99  |     2.02631      .64172      3.16    0.002      .7685615     3.284058
        100  |    .1299544    .6493153      0.20    0.841      -1.14268     1.402589
             |
       _cons |   -2.129608     .631634     -3.37    0.001     -3.367588    -.8916277
-------------+----------------------------------------------------------------

.
. * Mantel-Haenszel procedure
. use simbin_clustcow.dta, clear

. cc Y X, by(herd)

            herd |       OR       [95% Conf. Interval]   M-H Weight
-----------------+-------------------------------------------------
               1 |          .      .9048609          .            0 (exact)
               2 |       1.25      .1626352   9.892487     1.142857 (exact)
               3 |   5.833333      .7077892   70.40121          .48 (exact)
               4 |   .5357143      .0853896   3.264111     2.153846 (exact)
               5 |          .      .7951185          .            0 (exact)
               6 |       2.75      .4092367    21.5737     .8571429 (exact)
               7 |   3.333333      .4045859    40.9002     .6428571 (exact)
               8 |       2.75      .4622328   17.56168     .9655172 (exact)
               9 |        6.5      .5552238    330.517     .3448276 (exact)
              10 |          .             0          .            0 (exact)
              11 |   1.714286      .3243022   9.247996          1.4 (exact)
              12 |      1.625      .1542296   22.28922           .8 (exact)
              13 |   .6153846      .0448629   6.483884          1.3 (exact)
              14 |   2.954545      .3729396   35.64326     .7096774 (exact)
              15 |   .4666667      .0074234   10.13893        .9375 (exact)
              16 |        3.2      .2206377   178.2772     .4285714 (exact)
              17 |   .9090909      .1914922   4.326952          2.2 (exact)
              18 |   1.333333      .1849703   10.74942          1.2 (exact)
              19 |   1.781818      .3609598   9.173558     1.486486 (exact)
              20 |       2.75      .6028725   12.98858     1.297297 (exact)
              21 |   3.878788      .6954999    26.8677     .8684211 (exact)
              22 |     1.2375      .2872619   5.364267     2.105263 (exact)
              23 |   7.714286      .7574556   374.0486     .3589744 (exact)
              24 |   10.23077       1.04612    484.685         .325 (exact)
              25 |   4.470588      .3785522     231.49     .4146341 (exact)
              26 |   1.818182      .4426808   7.646529     1.833333 (exact)
              27 |   4.705882      .3995133     243.07     .4047619 (exact)
              28 |   4.705882      .3995133     243.07     .4047619 (exact)
              29 |   2.185714      .5033592   10.19969     1.555556 (exact)
              30 |   1.346154      .3473786   5.289449     2.311111 (exact)
              31 |       6.75      1.105351   70.43244     .6086957 (exact)
              32 |        3.3      .2363841   180.9459     .4347826 (exact)
              33 |   1.442308      .3753545   5.614642      2.26087 (exact)
              34 |        4.4      .3808243   226.3113     .4255319 (exact)
              35 |   1.586538      .4215339   6.072896     2.212766 (exact)
              36 |   3.046154      .7310814   13.76903     1.382979 (exact)
              37 |   2.307692      .5715921   9.726768     1.659574 (exact)
              38 |   2.222222      .3953127    15.5077     1.148936 (exact)
              39 |        3.8      .9207688   17.01576         1.25 (exact)
              40 |   2.142857      .5366468    8.95608         1.75 (exact)
              41 |   1.058824      .2550643   4.389841     2.428571 (exact)
              42 |   .7878788       .218652   2.832178     3.367347 (exact)
              43 |      4.125      .9439785   20.85567         1.12 (exact)
              44 |   1.477273      .4268405   5.139578     2.588235 (exact)
              45 |      1.575      .3149334   8.701329     1.568627 (exact)
              46 |   6.296296      1.615629   25.78857     1.038462 (exact)
              47 |   .3285024      .0642865    1.45527      3.90566 (exact)
              48 |   3.589744      .9603942    14.2219     1.471698 (exact)
              49 |        1.7      .5097952    5.73953     2.545455 (exact)
              50 |   .9565217       .190967   4.796742     2.090909 (exact)
              51 |   1.935185      .7341281   5.179395      3.56044 (exact)
              52 |   1.899209      .8232199   4.399525     4.728972 (exact)
              53 |   1.459574       .367796   6.240886     2.175926 (exact)
              54 |   .9814815      .2452764   3.928224     2.722689 (exact)
              55 |   2.933333      1.267148   6.883933     3.781513 (exact)
              56 |   1.890909      .7493979   4.910543     3.928571 (exact)
              57 |   1.142857      .3369456   3.962716     3.099237 (exact)
              58 |          2       .967757   4.145337            6 (exact)
              59 |   1.895425      .7884736   4.681192     4.340426 (exact)
              60 |   1.631579      .6253286   4.400487     3.958333 (exact)
              61 |    3.72205      1.790853   7.781826     4.322148 (exact)
              62 |        2.4      1.004589   5.961365     3.928571 (exact)
              63 |   2.368889      1.172023   4.809167     5.844156 (exact)
              64 |   2.695652       1.25491   5.894686     4.717949 (exact)
              65 |   1.689189      .3145984    11.2302     1.414013 (exact)
              66 |   2.737269      1.357276   5.549848     5.468354 (exact)
              67 |   1.538462      .3485597   7.686526     1.857143 (exact)
              68 |   1.793033      .8293041   3.940985     5.674419 (exact)
              69 |   2.216283      1.155521   4.261833     6.988506 (exact)
              70 |   1.274536      .6457513   2.522527     8.519774 (exact)
              71 |   2.643343       1.32713   5.315913     5.780899 (exact)
              72 |   1.755061      .9313773   3.312375      8.27933 (exact)
              73 |   1.903226      1.001306   3.627241         7.75 (exact)
              74 |   1.463889      .7390993   2.912365     7.868852 (exact)
              75 |   1.818815      .6222934   5.716085     3.053191 (exact)
              76 |          2      .1021853   119.1905     .4923077 (exact)
              77 |   1.728997      .9497046    3.15195     9.271357 (exact)
              78 |   .2717087      .0270621   1.479965     3.432692 (exact)
              79 |      1.125      .3677098   3.503059     3.674641 (exact)
              80 |   3.863636      1.291432   13.87083     2.095238 (exact)
              81 |   2.888889      1.556936   5.393682     6.781991 (exact)
              82 |   .6934813      .1680196   2.638256     3.353488 (exact)
              83 |   1.241429       .656701   2.353241     9.767442 (exact)
              84 |   3.448276      1.324761   10.02351     2.832558 (exact)
              85 |   2.450475      1.283746   4.730079     6.792627 (exact)
              86 |   2.676923      1.244167    5.97691     4.770642 (exact)
              87 |   4.538314      1.686463   14.15979     2.383562 (exact)
              88 |   1.393189      .6801907   2.883222     7.374429 (exact)
              89 |   1.559343      .7615576   3.238469     6.977974 (exact)
              90 |    5.38206      2.165676   15.13253     2.606061 (exact)
              91 |   1.889362      1.097655   3.256022     10.68182 (exact)
              92 |   3.203111      1.529742   7.003147     4.722449 (exact)
              93 |   2.220123      1.237445   4.008197     8.605691 (exact)
              94 |   3.210526      .9342177   13.99336     1.809524 (exact)
              95 |   2.166667      .8767537    5.70404     3.857143 (exact)
              96 |   2.498798      1.468905   4.257137     9.674419 (exact)
              97 |   2.045113      1.057989   4.028883     7.176259 (exact)
              98 |   1.794737      .9199104   3.562613     7.437722 (exact)
              99 |   1.354032      .8402866   2.182873     16.15635 (exact)
             100 |   1.088372      .5740342   2.067865     10.36977 (exact)
-----------------+-------------------------------------------------
           Crude |   1.797341       1.65396   1.953169              (exact)
    M-H combined |   2.008957       1.83547   2.198841
-----------------+-------------------------------------------------
Test of homogeneity (Tarone)   chi2(99) =   101.38  Pr>chi2 = 0.4148

                   Test that combined OR = 1:
                                Mantel-Haenszel chi2(1) =    233.10
                                                Pr>chi2 =    0.0000

. di "approx SE:  " (ln(2.008957)-ln(1.83547))/1.96
approx SE:  .04607913

.
. * robust variance
. use simcont_clustherd.dta, clear

. regress milk X, vce(cluster herd)

Linear regression                                      Number of obs =   11626
                                                       F(  1,    99) =    4.32
                                                       Prob > F      =  0.0403
                                                       R-squared     =  0.0266
                                                       Root MSE      =  10.733

(Std. Err. adjusted for 100 clusters in herd)
-------------+----------------------------------------------------------------
             |               Robust
         milk |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |     3.55661    1.711824      2.08    0.040      .1599804      6.95324
       _cons |     30.0215    1.092129     27.49    0.000      27.85448     32.18852
-------------+----------------------------------------------------------------

. use simcont_clustcow.dta, clear

. regress milk X, vce(cluster herd)

Linear regression                                      Number of obs =   11626
                                                       F(  1,    99) = 1232.30
                                                       Prob > F      =  0.0000
                                                       R-squared     =  0.0510
                                                       Root MSE      =  10.744

(Std. Err. adjusted for 100 clusters in herd)
-------------+----------------------------------------------------------------
             |               Robust
         milk |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    4.982006    .1419208     35.10    0.000      4.700404     5.263607
       _cons |    29.25664    .8738122     33.48    0.000       27.5228     30.99047
-------------+----------------------------------------------------------------

. use simbin_clustherd.dta, clear

. logit Y X, vce(cluster herd)

Iteration 0:   log pseudolikelihood = -6894.3552
Iteration 1:   log pseudolikelihood = -6815.0583
Iteration 2:   log pseudolikelihood = -6814.7785
Iteration 3:   log pseudolikelihood = -6814.7785

Logistic regression                               Number of obs   =      11626
                                                  Wald chi2(1)    =       6.31
                                                  Prob > chi2     =     0.0120
Log pseudolikelihood = -6814.7785                 Pseudo R2       =     0.0115

(Std. Err. adjusted for 100 clusters in herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5287317    .2105488      2.51    0.012      .1160637     .9413998
       _cons |   -1.241768     .145627     -8.53    0.000     -1.527192    -.9563443
-------------+----------------------------------------------------------------

. use simbin_clustcow.dta, clear

. logit Y X, vce(cluster herd)

Iteration 0:   log pseudolikelihood = -6910.3442
Iteration 1:   log pseudolikelihood =   -6811.48
Iteration 2:   log pseudolikelihood = -6811.0741
Iteration 3:   log pseudolikelihood = -6811.0741

Logistic regression                               Number of obs   =      11626
                                                  Wald chi2(1)    =     179.23
                                                  Prob > chi2     =     0.0000
Log pseudolikelihood = -6811.0741                 Pseudo R2       =     0.0144

(Std. Err. adjusted for 100 clusters in herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5863084    .0437941     13.39    0.000      .5004736     .6721432
       _cons |    -1.25032    .1144647    -10.92    0.000     -1.474667    -1.025973
-------------+----------------------------------------------------------------

.
. * GEE for binary outcome
. use simbin_clustherd.dta, clear

. xtgee Y X, fam(bin) link(logit) robust i(herd) corr(indep) /* same as robust variance */

Iteration 1: tolerance = 2.181e-09

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                                logit      Obs per group: min =        20
Family:                           binomial                     avg =     116.3
Correlation:                   independent                     max =       311
                                                Wald chi2(1)       =      6.31
Scale parameter:                         1      Prob > chi2        =    0.0120

Pearson chi2(11626):              11626.00      Deviance           =  13629.56
Dispersion (Pearson):                    1      Dispersion         =  1.172334

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5287317    .2105488      2.51    0.012      .1160637     .9413998
       _cons |   -1.241768     .145627     -8.53    0.000     -1.527192    -.9563443
-------------+----------------------------------------------------------------

. xtgee Y X, fam(bin) robust i(herd) corr(indep) /* logit is default link for fam(bin) */

Iteration 1: tolerance = 2.181e-09

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                                logit      Obs per group: min =        20
Family:                           binomial                     avg =     116.3
Correlation:                   independent                     max =       311
                                                Wald chi2(1)       =      6.31
Scale parameter:                         1      Prob > chi2        =    0.0120

Pearson chi2(11626):              11626.00      Deviance           =  13629.56
Dispersion (Pearson):                    1      Dispersion         =  1.172334

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5287317    .2105488      2.51    0.012      .1160637     .9413998
       _cons |   -1.241768     .145627     -8.53    0.000     -1.527192    -.9563443
-------------+----------------------------------------------------------------

. xtgee Y X, fam(bin) robust i(herd) corr(exch)

Iteration 1: tolerance = .06052491
Iteration 2: tolerance = .00166018
Iteration 3: tolerance = .00001941
Iteration 4: tolerance = 1.325e-07

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                                logit      Obs per group: min =        20
Family:                           binomial                     avg =     116.3
Correlation:                  exchangeable                     max =       311
                                                Wald chi2(1)       =      9.98
Scale parameter:                         1      Prob > chi2        =    0.0016

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5588542    .1768743      3.16    0.002      .2121869     .9055215
       _cons |   -1.109623    .1251571     -8.87    0.000     -1.354926    -.8643196
-------------+----------------------------------------------------------------

. estat wcor, compact /* displays working correlation(s) */

Error structure: exchangeable
Estimated within-herd correlation: .15134071

. use simbin_clustcow.dta, clear

. xtgee Y X, fam(bin) robust i(herd) corr(indep)

Iteration 1: tolerance = 3.102e-09

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                                logit      Obs per group: min =        20
Family:                           binomial                     avg =     116.3
Correlation:                   independent                     max =       311
                                                Wald chi2(1)       =    179.23
Scale parameter:                         1      Prob > chi2        =    0.0000

Pearson chi2(11626):              11626.00      Deviance           =  13622.15
Dispersion (Pearson):                    1      Dispersion         =  1.171697

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5863084    .0437941     13.39    0.000      .5004736     .6721432
       _cons |    -1.25032    .1144647    -10.92    0.000     -1.474667    -1.025973
-------------+----------------------------------------------------------------

. xtgee Y X, fam(bin) robust i(herd) corr(exch)

Iteration 1: tolerance = .06275994
Iteration 2: tolerance = .00152416
Iteration 3: tolerance = .00003718
Iteration 4: tolerance = 1.113e-06
Iteration 5: tolerance = 2.966e-08

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                                logit      Obs per group: min =        20
Family:                           binomial                     avg =     116.3
Correlation:                  exchangeable                     max =       311
                                                Wald chi2(1)       =    184.02
Scale parameter:                         1      Prob > chi2        =    0.0000

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    .5687287    .0419247     13.57    0.000      .4865578     .6508997
       _cons |   -1.112355    .0948099    -11.73    0.000     -1.298179    -.9265314
-------------+----------------------------------------------------------------

. estat wcor, compact

Error structure: exchangeable
Estimated within-herd correlation: .15631633

.
. * GEE for continuous outcome
. use simcont_clustherd.dta, clear

. xtgee milk X, fam(gauss) link(iden) robust i(herd) corr(exch)

Iteration 1: tolerance = .05318109
Iteration 2: tolerance = .00034544
Iteration 3: tolerance = 6.094e-07

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                             identity      Obs per group: min =        20
Family:                           Gaussian                     avg =     116.3
Correlation:                  exchangeable                     max =       311
                                                Wald chi2(1)       =      6.51
Scale parameter:                  116.7291      Prob > chi2        =    0.0107

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
         milk |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    3.797275    1.488196      2.55    0.011      .8804645     6.714085
       _cons |    31.13476    1.040243     29.93    0.000      29.09592      33.1736
-------------+----------------------------------------------------------------

. use simcont_clustcow.dta, clear

. xtgee milk X, fam(gauss) link(iden) robust i(herd) corr(exch)

Iteration 1: tolerance = .04589737
Iteration 2: tolerance = .00008274
Iteration 3: tolerance = 1.942e-07

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                             identity      Obs per group: min =        20
Family:                           Gaussian                     avg =     116.3
Correlation:                  exchangeable                     max =       311
                                                Wald chi2(1)       =   1234.47
Scale parameter:                  117.3372      Prob > chi2        =    0.0000

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
         milk |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |    4.968177    .1414023     35.14    0.000      4.691034     5.245321
       _cons |    30.64796    .7217876     42.46    0.000      29.23328     32.06264
-------------+----------------------------------------------------------------

.
. * analyses for pig data
. 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 |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
       ar_g1 |    .6465241    .2203379      2.93    0.003      .2146697     1.078378
       _cons |   -.1448309    .1556373     -0.93    0.352     -.4498744     .1602125
-------------+----------------------------------------------------------------

. * fixed effects logistic regression
. logit pn ar_g1 i.farm, asis

Iteration 0:   log likelihood = -234.95215
Iteration 1:   log likelihood =  -194.7543
Iteration 2:   log likelihood = -193.88065
Iteration 3:   log likelihood = -193.86711
Iteration 4:   log likelihood =  -193.8671

Logistic regression                               Number of obs   =        341
                                                  LR chi2(15)     =      82.17
                                                  Prob > chi2     =     0.0000
Log likelihood =  -193.8671                       Pseudo R2       =     0.1749

-------------+----------------------------------------------------------------
           pn |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
       ar_g1 |    .3651127    .2678889      1.36    0.173     -.1599398     .8901652
             |
        farm |
          2  |   -.6608706    .7345805     -0.90    0.368     -2.100622     .7788806
          3  |   -3.181252    .7360392     -4.32    0.000     -4.623862    -1.738642
          4  |   -1.658611    .6726559     -2.47    0.014     -2.976992    -.3402299
          5  |   -.8748085    .7372227     -1.19    0.235     -2.319738     .5701213
          6  |    1.411727    1.158223      1.22    0.223     -.8583474     3.681801
          7  |   -1.241977    .7120953     -1.74    0.081     -2.637658     .1537041
          8  |   -.6862302    .7012406     -0.98    0.328     -2.060636      .688176
          9  |   -1.644387    .6834869     -2.41    0.016     -2.983997    -.3047778
         10  |   -1.536526    .7318713     -2.10    0.036     -2.970967    -.1020843
         11  |   -1.310881    .7055684     -1.86    0.063     -2.693769      .072008
         12  |    -1.40703    .7684634     -1.83    0.067      -2.91319     .0991308
         13  |   -2.973781    .7437175     -4.00    0.000      -4.43144    -1.516121
         14  |   -1.904232    .6830497     -2.79    0.005     -3.242985    -.5654794
         15  |   -2.792443    .7512611     -3.72    0.000     -4.264887    -1.319998
             |
       _cons |    1.461934    .5707978      2.56    0.010      .3431912     2.580678
-------------+----------------------------------------------------------------

. * Mantel-Haenszel procedure
. cc pn ar_g1, by(farm)

Farm identificat |       OR       [95% Conf. Interval]   M-H Weight
-----------------+-------------------------------------------------
               1 |          0             0   1.987843     1.230769 (exact)
               2 |   3.333333       .249723   183.4806     .4285714 (exact)
               3 |   1.636364        .15285    22.6999     .7857143 (exact)
               4 |          1      .1657092   6.035157     1.615385 (exact)
               5 |       .125      .0023556   1.791506     2.105263 (exact)
               6 |          .             0          .            0 (exact)
               7 |   2.666667      .1641335   157.1844     .4285714 (exact)
               8 |          8      .3194679    494.706     .1785714 (exact)
               9 |   3.928571      .4498158    49.3709          .56 (exact)
              10 |       1.25      .1282001     12.298     .9411765 (exact)
              11 |   1.428571      .1745708   13.12955            1 (exact)
              12 |   .8333333       .047384    11.9645     .8571429 (exact)
              13 |        4.5      .3690823   64.38143           .4 (exact)
              14 |         .8      .0833334   7.784078         1.25 (exact)
              15 |   2.166667      .1305404   26.50935     .5454545 (exact)
-----------------+-------------------------------------------------
           Crude |   1.908894       1.21155   3.009556              (exact)
    M-H combined |   1.412829      .8470242   2.356586
-----------------+-------------------------------------------------
Test of homogeneity (Tarone)   chi2(14) =    16.06  Pr>chi2 = 0.3096

                   Test that combined OR = 1:
                                Mantel-Haenszel chi2(1) =      1.78
                                                Pr>chi2 =    0.1821

. di "ln OR:  " ln(1.412829)
ln OR:  .34559408

. di "SE:  " (ln(1.412829)-ln(.8470242))/1.96
SE:  .26103066

. * robust variance
. logit pn ar_g1, vce(cluster farm)

Iteration 0:   log pseudolikelihood = -234.95215
Iteration 1:   log pseudolikelihood = -230.59287
Iteration 2:   log pseudolikelihood = -230.59173
Iteration 3:   log pseudolikelihood = -230.59173

Logistic regression                               Number of obs   =        341
                                                  Wald chi2(1)    =       5.48
                                                  Prob > chi2     =     0.0192
Log pseudolikelihood = -230.59173                 Pseudo R2       =     0.0186

(Std. Err. adjusted for 15 clusters in farm)
-------------+----------------------------------------------------------------
             |               Robust
           pn |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
       ar_g1 |    .6465241    .2760887      2.34    0.019      .1054002     1.187648
       _cons |   -.1448309     .278696     -0.52    0.603     -.6910651     .4014032
-------------+----------------------------------------------------------------

. * random effects model
. meqrlogit pn ar_g1 || farm:

Refining starting values:

Iteration 0:   log likelihood = -213.91344
Iteration 1:   log likelihood = -213.51368
Iteration 2:   log likelihood =  -213.5118

Performing gradient-based optimization:

Iteration 0:   log likelihood =  -213.5118
Iteration 1:   log likelihood =  -213.5118

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 points =   7                        Wald chi2(1)       =      2.86
Log likelihood =  -213.5118                     Prob > chi2        =    0.0905

-------------+----------------------------------------------------------------
           pn |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
       ar_g1 |    .4369353    .2581461      1.69    0.091     -.0690217     .9428923
       _cons |    .0196485    .3009417      0.07    0.948     -.5701864     .6094834
-------------+----------------------------------------------------------------

-----------------------------+------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
farm: Identity               |
                  var(_cons) |   .8773258   .4324896      .3338486    2.305537
-----------------------------+------------------------------------------------
LR test vs. logistic regression: chibar2(01) =    34.16 Prob>=chibar2 =  0.0000

. di .4369353/(sqrt(1+.346*.8773258))
.38269474

. * GEE with exchangeable correlation structure
. xtgee pn ar_g1, fam(bin) robust i(farm) corr(exch)

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: min =        14
Family:                           binomial                     avg =      22.7
Correlation:                  exchangeable                     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% Con f. Interval]
-------------+----------------------------------------------------------------
       ar_g1 |    .3539583    .2151855      1.64    0.100     -.0677976     .7757142
       _cons |    .0183926    .2714184      0.07    0.946     -.5135777     .5503629
-------------+----------------------------------------------------------------

. estat wcor, compact

Error structure: exchangeable
Estimated within-farm correlation: .17940234

.
. * margins for linear mixed models
. use simcont_clustcow.dta, clear

. mixed milk i.X || herd:, reml

Performing EM optimization:

Performing gradient-based optimization:

Iteration 0:   log restricted-likelihood = -40947.175
Iteration 1:   log restricted-likelihood = -40947.175

Computing standard errors:

Mixed-effects REML regression                   Number of obs      =     11626
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 |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
         1.X |    4.968194    .1492174     33.30    0.000      4.675733     5.260655
       _cons |    30.64647    .7281276     42.09    0.000      29.21936     32.07357
-------------+----------------------------------------------------------------

-----------------------------+------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
herd: Identity               |
                  var(_cons) |   51.41189   7.459591      38.68644    68.32323
-----------------------------+------------------------------------------------
               var(Residual) |   64.71069   .8524578      63.06129    66.40324
-----------------------------+------------------------------------------------
LR test vs. linear regression: chibar2(01) =  6310.00 Prob >= chibar2 =  0.0000

. margins X /* broad inference space predictions ~ zero random effects */

Adjusted predictions                              Number of obs   =      11626


Expression  : Linear prediction, fixed portion, predict()


-------------+----------------------------------------------------------------
             |            Delta-method
             |     Margin    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |
          0  |    30.64647    .7281276     42.09    0.000      29.21936     32.07357
          1  |    35.61466    .7279903     48.92    0.000      34.18783     37.04149
-------------+----------------------------------------------------------------

.
. * margins for discrete models
. use simbin_clustcow.dta, clear

. * GLMM
. meqrlogit Y i.X || herd:

Refining starting values:

Iteration 0:   log likelihood = -5999.0538
Iteration 1:   log likelihood = -5995.9727
Iteration 2:   log likelihood = -5995.9698

Performing gradient-based optimization:

Iteration 0:   log likelihood = -5995.9698
Iteration 1:   log likelihood = -5995.9698

Mixed-effects logistic regression               Number of obs      =     11626
Group variable: herd                            Number of groups   =       100

                                                Obs per group: min =        20
                                                               avg =     116.3
                                                               max =       311

Integration points =   7                        Wald chi2(1)       =    229.28
Log likelihood = -5995.9698                     Prob > chi2        =    0.0000

-------------+----------------------------------------------------------------
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
         1.X |     .697479     .046063     15.14    0.000      .6071972     .7877609
       _cons |   -1.361173    .1112103    -12.24    0.000     -1.579141    -1.143204
-------------+----------------------------------------------------------------

-----------------------------+------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
herd: Identity               |
                  var(_cons) |   1.068172   .1682549      .7844459     1.45452
-----------------------------+------------------------------------------------
LR test vs. logistic regression: chibar2(01) =  1630.21 Prob>=chibar2 =  0.0000

. * gives an error: margins X
. * medians (broad inference space)
. margins X, expression(1/(1+exp(-predict(xb))))

Adjusted predictions                              Number of obs   =      11626


Expression  : 1/(1+exp(-predict(xb)))


-------------+----------------------------------------------------------------
             |            Delta-method
             |     Margin    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |
          0  |    .2040498     .018062     11.30    0.000      .1686488     .2394508
          1  |    .3399104    .0246436     13.79    0.000      .2916098      .388211
-------------+----------------------------------------------------------------

. * means (broad inference space, PA predictions)
. margins X, expression(1/(1+exp(-predict(xb)/(sqrt(1+.346*1.068314)))))

Adjusted predictions                              Number of obs   =      11626


Expression  : 1/(1+exp(-predict(xb)/(sqrt(1+.346*1.068314))))


-------------+----------------------------------------------------------------
             |            Delta-method
             |     Margin    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |
          0  |    .2381076    .0172389     13.81    0.000      .2043201     .2718952
          1  |    .3619047    .0216728     16.70    0.000      .3194269     .4043825
-------------+----------------------------------------------------------------

. * GEE with PA predictions
. xtgee Y i.X, fam(bin) link(logit) robust i(herd) corr(exch)

Iteration 1: tolerance = .06275994
Iteration 2: tolerance = .00152416
Iteration 3: tolerance = .00003718
Iteration 4: tolerance = 1.113e-06
Iteration 5: tolerance = 2.966e-08

GEE population-averaged model                   Number of obs      =     11626
Group variable:                       herd      Number of groups   =       100
Link:                                logit      Obs per group: min =        20
Family:                           binomial                     avg =     116.3
Correlation:                  exchangeable                     max =       311
                                                Wald chi2(1)       =    184.02
Scale parameter:                         1      Prob > chi2        =    0.0000

(Std. Err. adjusted for clustering on herd)
-------------+----------------------------------------------------------------
             |               Robust
            Y |      Coef.    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
         1.X |    .5687287    .0419247     13.57    0.000      .4865578     .6508997
       _cons |   -1.112355    .0948099    -11.73    0.000     -1.298179    -.9265314
-------------+----------------------------------------------------------------

. margins X

Adjusted predictions                              Number of obs   =      11626
Model VCE    : Robust


Expression  : Pr(Y != 0), predict()


-------------+----------------------------------------------------------------
             |            Delta-method
             |     Margin    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |
          0  |     .247432    .0176545     14.02    0.000      .2128298     .2820342
          1  |    .3673443    .0217279     16.91    0.000      .3247585     .4099302
-------------+----------------------------------------------------------------

. margins X, expression(1/(1+exp(-predict(xb)))) /* same as above */

Adjusted predictions                              Number of obs   =      11626
Model VCE    : Robust


Expression  : 1/(1+exp(-predict(xb)))


-------------+----------------------------------------------------------------
             |            Delta-method
             |     Margin    Std. Err.       z    P>|z|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
           X |
          0  |     .247432    .0176545     14.02    0.000      .2128298     .2820342
          1  |    .3673443    .0217279     16.91    0.000      .3247585     .4099302
-------------+----------------------------------------------------------------

.
end of do-file

.