Seizure Counts for Epileptics

Thall and Vail (1990) give a data set on two-week seizure counts for 59 epileptics. The number of seizures was recorded for a baseline period of 8 weeks, and then patients were randomly assigned to a treatment group or a control group. Counts were then recorded for four successive two-week periods. The subject's age is the only covariate.

epil

Format

This data frame has 236 rows and the following 9 columns:

y

the count for the 2-week period.

trt

treatment, "placebo" or "progabide".

base

the counts in the baseline 8-week period.

age

subject's age, in years.

V4

0/1 indicator variable of period 4.

subject

subject number, 1 to 59.

period

period, 1 to 4.

lbase

log-counts for the baseline period, centred to have zero mean.

lage

log-ages, centred to have zero mean.

Source

Thall, P. F. and Vail, S. C. (1990) Some covariance models for longitudinal count data with over-dispersion. Biometrics 46, 657--671.

References

Venables, W. N. and Ripley, B. D. (2002) Modern Applied Statistics with S. Fourth Edition. Springer.

Examples

summary(glm(y ~ lbase*trt + lage + V4, family = poisson,
            data = epil), cor = FALSE)
#> 
#> Call:
#> glm(formula = y ~ lbase * trt + lage + V4, family = poisson, 
#>     data = epil)
#> 
#> Deviance Residuals: 
#>     Min       1Q   Median       3Q      Max  
#> -5.0915  -1.4126  -0.2739   0.7580  10.7711  
#> 
#> Coefficients:
#>                    Estimate Std. Error z value Pr(>|z|)    
#> (Intercept)         1.89791    0.04260  44.552  < 2e-16 ***
#> lbase               0.94862    0.04360  21.759  < 2e-16 ***
#> trtprogabide       -0.34588    0.06100  -5.670 1.42e-08 ***
#> lage                0.88760    0.11650   7.619 2.56e-14 ***
#> V4                 -0.15977    0.05458  -2.927  0.00342 ** 
#> lbase:trtprogabide  0.56154    0.06352   8.841  < 2e-16 ***
#> ---
#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#> 
#> (Dispersion parameter for poisson family taken to be 1)
#> 
#>     Null deviance: 2517.83  on 235  degrees of freedom
#> Residual deviance:  869.07  on 230  degrees of freedom
#> AIC: 1647
#> 
#> Number of Fisher Scoring iterations: 5
#> 
epil2 <- epil[epil$period == 1, ]
epil2["period"] <- rep(0, 59); epil2["y"] <- epil2["base"]
epil["time"] <- 1; epil2["time"] <- 4
epil2 <- rbind(epil, epil2)
epil2$pred <- unclass(epil2$trt) * (epil2$period > 0)
epil2$subject <- factor(epil2$subject)
epil3 <- aggregate(epil2, list(epil2$subject, epil2$period > 0),
   function(x) if(is.numeric(x)) sum(x) else x[1])
epil3$pred <- factor(epil3$pred,
   labels = c("base", "placebo", "drug"))

contrasts(epil3$pred) <- structure(contr.sdif(3),
    dimnames = list(NULL, c("placebo-base", "drug-placebo")))
summary(glm(y ~ pred + factor(subject) + offset(log(time)),
            family = poisson, data = epil3), cor = FALSE)
#> 
#> Call:
#> glm(formula = y ~ pred + factor(subject) + offset(log(time)), 
#>     family = poisson, data = epil3)
#> 
#> Deviance Residuals: 
#>     Min       1Q   Median       3Q      Max  
#> -5.2928  -0.7350   0.0000   0.6997   4.7145  
#> 
#> Coefficients:
#>                     Estimate Std. Error z value Pr(>|z|)    
#> (Intercept)        1.122e+00  2.008e-01   5.590 2.28e-08 ***
#> predplacebo-base   1.087e-01  4.691e-02   2.318 0.020474 *  
#> preddrug-placebo  -1.016e-01  6.507e-02  -1.561 0.118431    
#> factor(subject)2  -2.300e-15  2.828e-01   0.000 1.000000    
#> factor(subject)3  -3.857e-01  3.144e-01  -1.227 0.219894    
#> factor(subject)4  -1.744e-01  2.960e-01  -0.589 0.555847    
#> factor(subject)5   1.577e+00  2.197e-01   7.178 7.08e-13 ***
#> factor(subject)6   6.729e-01  2.458e-01   2.738 0.006182 ** 
#> factor(subject)7  -4.082e-02  2.858e-01  -0.143 0.886411    
#> factor(subject)8   1.758e+00  2.166e-01   8.117 4.77e-16 ***
#> factor(subject)9   5.878e-01  2.494e-01   2.356 0.018454 *  
#> factor(subject)10  5.423e-01  2.515e-01   2.156 0.031060 *  
#> factor(subject)11  1.552e+00  2.202e-01   7.048 1.81e-12 ***
#> factor(subject)12  9.243e-01  2.364e-01   3.910 9.22e-05 ***
#> factor(subject)13  3.075e-01  2.635e-01   1.167 0.243171    
#> factor(subject)14  1.212e+00  2.278e-01   5.320 1.04e-07 ***
#> factor(subject)15  1.765e+00  2.164e-01   8.153 3.54e-16 ***
#> factor(subject)16  9.708e-01  2.348e-01   4.134 3.57e-05 ***
#> factor(subject)17 -4.082e-02  2.858e-01  -0.143 0.886411    
#> factor(subject)18  2.236e+00  2.104e-01  10.629  < 2e-16 ***
#> factor(subject)19  2.776e-01  2.651e-01   1.047 0.295060    
#> factor(subject)20  3.646e-01  2.603e-01   1.401 0.161324    
#> factor(subject)21  3.922e-02  2.801e-01   0.140 0.888645    
#> factor(subject)22 -8.338e-02  2.889e-01  -0.289 0.772894    
#> factor(subject)23  1.823e-01  2.708e-01   0.673 0.500777    
#> factor(subject)24  8.416e-01  2.393e-01   3.517 0.000436 ***
#> factor(subject)25  2.069e+00  2.123e-01   9.750  < 2e-16 ***
#> factor(subject)26 -5.108e-01  3.266e-01  -1.564 0.117799    
#> factor(subject)27 -2.231e-01  3.000e-01  -0.744 0.456990    
#> factor(subject)28  1.386e+00  2.236e-01   6.200 5.66e-10 ***
#> factor(subject)29  1.604e+00  2.227e-01   7.203 5.90e-13 ***
#> factor(subject)30  1.023e+00  2.372e-01   4.313 1.61e-05 ***
#> factor(subject)31  9.149e-02  2.821e-01   0.324 0.745700    
#> factor(subject)32 -3.111e-02  2.909e-01  -0.107 0.914822    
#> factor(subject)33  4.710e-01  2.597e-01   1.814 0.069736 .  
#> factor(subject)34  3.887e-01  2.640e-01   1.473 0.140879    
#> factor(subject)35  1.487e+00  2.250e-01   6.609 3.87e-11 ***
#> factor(subject)36  3.598e-01  2.656e-01   1.355 0.175551    
#> factor(subject)37 -1.221e-01  2.979e-01  -0.410 0.681943    
#> factor(subject)38  1.344e+00  2.283e-01   5.889 3.90e-09 ***
#> factor(subject)39  1.082e+00  2.354e-01   4.596 4.30e-06 ***
#> factor(subject)40 -7.687e-01  3.634e-01  -2.116 0.034384 *  
#> factor(subject)41  1.656e-01  2.772e-01   0.597 0.550234    
#> factor(subject)42  5.227e-02  2.848e-01   0.184 0.854388    
#> factor(subject)43  1.543e+00  2.239e-01   6.891 5.54e-12 ***
#> factor(subject)44  9.605e-01  2.393e-01   4.014 5.96e-05 ***
#> factor(subject)45  1.177e+00  2.326e-01   5.061 4.18e-07 ***
#> factor(subject)46 -5.275e-01  3.355e-01  -1.572 0.115840    
#> factor(subject)47  1.053e+00  2.363e-01   4.456 8.35e-06 ***
#> factor(subject)48 -5.275e-01  3.355e-01  -1.572 0.115840    
#> factor(subject)49  2.949e+00  2.082e-01  14.168  < 2e-16 ***
#> factor(subject)50  3.887e-01  2.640e-01   1.473 0.140879    
#> factor(subject)51  1.038e+00  2.367e-01   4.385 1.16e-05 ***
#> factor(subject)52  5.711e-01  2.548e-01   2.241 0.025023 *  
#> factor(subject)53  1.670e+00  2.215e-01   7.538 4.76e-14 ***
#> factor(subject)54  4.443e-01  2.611e-01   1.702 0.088759 .  
#> factor(subject)55  2.674e-01  2.709e-01   0.987 0.323618    
#> factor(subject)56  1.124e+00  2.341e-01   4.800 1.59e-06 ***
#> factor(subject)57  2.674e-01  2.709e-01   0.987 0.323618    
#> factor(subject)58 -6.017e-01  3.436e-01  -1.751 0.079911 .  
#> factor(subject)59 -7.556e-02  2.942e-01  -0.257 0.797331    
#> ---
#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#> 
#> (Dispersion parameter for poisson family taken to be 1)
#> 
#>     Null deviance: 3180.82  on 117  degrees of freedom
#> Residual deviance:  303.16  on  57  degrees of freedom
#> AIC: 1003.5
#> 
#> Number of Fisher Scoring iterations: 5
#> 

summary(glmmPQL(y ~ lbase*trt + lage + V4,
                random = ~ 1 | subject,
                family = poisson, data = epil))
#> iteration 1
#> iteration 2
#> iteration 3
#> iteration 4
#> iteration 5
#> Linear mixed-effects model fit by maximum likelihood
#>  Data: epil 
#>   AIC BIC logLik
#>    NA  NA     NA
#> 
#> Random effects:
#>  Formula: ~1 | subject
#>         (Intercept) Residual
#> StdDev:   0.4442704 1.400807
#> 
#> Variance function:
#>  Structure: fixed weights
#>  Formula: ~invwt 
#> Fixed effects: y ~ lbase * trt + lage + V4 
#>                         Value Std.Error  DF   t-value p-value
#> (Intercept)         1.8696677 0.1055620 176 17.711554  0.0000
#> lbase               0.8818228 0.1292834  54  6.820849  0.0000
#> trtprogabide       -0.3095253 0.1490438  54 -2.076740  0.0426
#> lage                0.5335460 0.3463119  54  1.540652  0.1292
#> V4                 -0.1597696 0.0774521 176 -2.062819  0.0406
#> lbase:trtprogabide  0.3415425 0.2033325  54  1.679725  0.0988
#>  Correlation: 
#>                    (Intr) lbase  trtprg lage   V4    
#> lbase              -0.126                            
#> trtprogabide       -0.691  0.089                     
#> lage               -0.103 -0.038  0.088              
#> V4                 -0.162  0.000  0.000  0.000       
#> lbase:trtprogabide  0.055 -0.645 -0.184  0.267  0.000
#> 
#> Standardized Within-Group Residuals:
#>         Min          Q1         Med          Q3         Max 
#> -2.13240534 -0.63871136 -0.08486339  0.41960195  4.97872138 
#> 
#> Number of Observations: 236
#> Number of Groups: 59 
summary(glmmPQL(y ~ pred, random = ~1 | subject,
                family = poisson, data = epil3))
#> iteration 1
#> iteration 2
#> iteration 3
#> iteration 4
#> iteration 5
#> iteration 6
#> iteration 7
#> iteration 8
#> Linear mixed-effects model fit by maximum likelihood
#>  Data: epil3 
#>   AIC BIC logLik
#>    NA  NA     NA
#> 
#> Random effects:
#>  Formula: ~1 | subject
#>         (Intercept) Residual
#> StdDev:   0.7257895  2.16629
#> 
#> Variance function:
#>  Structure: fixed weights
#>  Formula: ~invwt 
#> Fixed effects: y ~ pred 
#>                      Value  Std.Error DF   t-value p-value
#> (Intercept)       3.213631 0.10569117 58 30.405865  0.0000
#> predplacebo-base  0.110855 0.09989089 57  1.109763  0.2718
#> preddrug-placebo -0.105613 0.13480483 57 -0.783450  0.4366
#>  Correlation: 
#>                  (Intr) prdpl-
#> predplacebo-base  0.081       
#> preddrug-placebo -0.010 -0.700
#> 
#> Standardized Within-Group Residuals:
#>        Min         Q1        Med         Q3        Max 
#> -2.0446864 -0.4765135 -0.1975651  0.3145761  2.6532834 
#> 
#> Number of Observations: 118
#> Number of Groups: 59