Generalized Linear Mixed Models and GEE

---

Back to schedule

---

Read Pairwise Dataset

library(dyadr)
library(dplyr)
library(lme4)

riggsp  <- read.csv("riggsp.csv")

Artificially Dichotomize Satisfaction

riggsp <- riggsp %>%
  mutate(Sat_bin_A = ifelse(Sat_A < 40, 0, 1))

Logistic Multilevel Modeling (Binary variables)

To account for the nonindependence, we can make use of the glmer() function from the lme4 package. Note that we are asking for the variance of intercepts across dyads, that is the random intercept in traditonal multilevel modeling. The gls() function in the nlme package does not have an option for specifying a link function (i.e., there is no family = option). The syntax of glmer() differs a bit from gls() in that the random effects are specified within the formula: + (1/cuplid).

Indistinguishable Dyads

apim_bin <- glmer(Sat_bin_A ~ Anxiety_A + Anxiety_P 
                  + (1|Dyad),
                  data = riggsp,
                  family = binomial,
                  na.action = na.omit)

summary(apim_bin)

## Generalized linear mixed model fit by maximum likelihood (Laplace
##   Approximation) [glmerMod]
##  Family: binomial  ( logit )
## Formula: Sat_bin_A ~ Anxiety_A + Anxiety_P + (1 | Dyad)
##    Data: riggsp
## 
##      AIC      BIC   logLik deviance df.resid 
##    297.3    312.3   -144.7    289.3      306 
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.29838  0.06707  0.18594  0.29717  1.09984 
## 
## Random effects:
##  Groups Name        Variance Std.Dev.
##  Dyad   (Intercept) 4.909    2.216   
## Number of obs: 310, groups:  Dyad, 155
## 
## Fixed effects:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   6.2562     1.4174   4.414 1.02e-05 ***
## Anxiety_A    -0.6669     0.2148  -3.105 0.001902 ** 
## Anxiety_P    -0.7453     0.2218  -3.361 0.000776 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##           (Intr) Anxt_A
## Anxiety_A -0.815       
## Anxiety_P -0.835  0.494

Distinguishable Dyads

Interaction approach.

apim_bin_di <- glmer(Sat_bin_A ~ Anxiety_A + Anxiety_P + Gender_A 
                     + Anxiety_A*Gender_A + Anxiety_P*Gender_A +
                     (1|Dyad),
                     data = riggsp,
                     family = binomial,
                     na.action = na.omit)

summary(apim_bin_di)

## Generalized linear mixed model fit by maximum likelihood (Laplace
##   Approximation) [glmerMod]
##  Family: binomial  ( logit )
## Formula: 
## Sat_bin_A ~ Anxiety_A + Anxiety_P + Gender_A + Anxiety_A * Gender_A +  
##     Anxiety_P * Gender_A + (1 | Dyad)
##    Data: riggsp
## 
##      AIC      BIC   logLik deviance df.resid 
##    301.4    327.5   -143.7    287.4      303 
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -3.02593  0.07182  0.18324  0.29759  1.12100 
## 
## Random effects:
##  Groups Name        Variance Std.Dev.
##  Dyad   (Intercept) 4.946    2.224   
## Number of obs: 310, groups:  Dyad, 155
## 
## Fixed effects:
##                    Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          6.4419     1.4666   4.393 1.12e-05 ***
## Anxiety_A           -0.7056     0.2235  -3.157  0.00159 ** 
## Anxiety_P           -0.7447     0.2257  -3.299  0.00097 ***
## Gender_A             0.6170     0.6902   0.894  0.37132    
## Anxiety_A:Gender_A  -0.2624     0.2052  -1.279  0.20090    
## Anxiety_P:Gender_A   0.0902     0.2008   0.449  0.65332    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Anxt_A Anxt_P Gndr_A A_A:G_
## Anxiety_A   -0.818                            
## Anxiety_P   -0.826  0.481                     
## Gender_A     0.186 -0.217 -0.074              
## Anxty_A:G_A -0.253  0.211  0.139 -0.528       
## Anxty_P:G_A  0.079 -0.029 -0.045 -0.492 -0.438

Generalized Estimating Equations (GEE) with Indistinguishable Dyads

library(gee)

apim_gee <- gee(Sat_bin_A ~ Anxiety_A + Anxiety_P, 
                id = Dyad,
                data = riggsp, 
                na.action = na.omit,
                family = binomial, 
                corstr = "unstructured")

## (Intercept)   Anxiety_A   Anxiety_P 
##   3.7106315  -0.4020834  -0.4521423

summary(apim_gee)

## 
##  GEE:  GENERALIZED LINEAR MODELS FOR DEPENDENT DATA
##  gee S-function, version 4.13 modified 98/01/27 (1998) 
## 
## Model:
##  Link:                      Logit 
##  Variance to Mean Relation: Binomial 
##  Correlation Structure:     Unstructured 
## 
## Call:
## gee(formula = Sat_bin_A ~ Anxiety_A + Anxiety_P, id = Dyad, data = riggsp, 
##     na.action = na.omit, family = binomial, corstr = "unstructured")
## 
## Summary of Residuals:
##         Min          1Q      Median          3Q         Max 
## -0.95045489  0.04035748  0.14130931  0.24507491  0.59890463 
## 
## 
## Coefficients:
##               Estimate Naive S.E.   Naive z Robust S.E.  Robust z
## (Intercept)  3.7107906  0.5985024  6.200126   0.5475257  6.777382
## Anxiety_A   -0.4028614  0.1167490 -3.450663   0.1149101 -3.505884
## Anxiety_P   -0.4515040  0.1177346 -3.834931   0.1165553 -3.873732
## 
## Estimated Scale Parameter:  0.9875912
## Number of Iterations:  2
## 
## Working Correlation
##           [,1]      [,2]
## [1,] 1.0000000 0.3130389
## [2,] 0.3130389 1.0000000

---

Back to schedule

---