rm(list=ls())
## load dataset
data(breslow.dat, package = "robust")
names(breslow.dat)
 [1] "ID"    "Y1"    "Y2"    "Y3"    "Y4"    "Base"  "Age"   "Trt"  
 [9] "Ysum"  "sumY"  "Age10" "Base4"
summary(breslow.dat[c(6,7,8,10)])
      Base             Age               Trt          sumY       
 Min.   :  6.00   Min.   :18.00   placebo  :28   Min.   :  0.00  
 1st Qu.: 12.00   1st Qu.:23.00   progabide:31   1st Qu.: 11.50  
 Median : 22.00   Median :28.00                  Median : 16.00  
 Mean   : 31.22   Mean   :28.34                  Mean   : 33.05  
 3rd Qu.: 41.00   3rd Qu.:32.00                  3rd Qu.: 36.00  
 Max.   :151.00   Max.   :42.00                  Max.   :302.00  
## look at the response variable in more detail
opar <- par(no.readonly = TRUE)
par(mfrow = c(1,2))
attach(breslow.dat)
The following objects are masked from breslow.dat (pos = 3):

    Age, Age10, Base, Base4, ID, sumY, Trt, Y1, Y2, Y3, Y4,
    Ysum

The following objects are masked from breslow.dat (pos = 4):

    Age, Age10, Base, Base4, ID, sumY, Trt, Y1, Y2, Y3, Y4,
    Ysum
hist(sumY, breaks = 20, xlab = "Seizure Count", main = "Distribution of Seizures")
boxplot(sumY ~ Trt, xlab = "Treatment", main = "Group Comparisons")
par(opar)

## from above figure, we can see the skewed nature of the dependent variable and the possible presence of outliers.
## fit the poisson regression
fit <- glm(sumY ~ Base + Age + Trt, data = breslow.dat, family = poisson())
summary(fit)

Call:
glm(formula = sumY ~ Base + Age + Trt, family = poisson(), data = breslow.dat)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-6.0569  -2.0433  -0.9397   0.7929  11.0061  

Coefficients:
               Estimate Std. Error z value Pr(>|z|)    
(Intercept)   1.9488259  0.1356191  14.370  < 2e-16 ***
Base          0.0226517  0.0005093  44.476  < 2e-16 ***
Age           0.0227401  0.0040240   5.651 1.59e-08 ***
Trtprogabide -0.1527009  0.0478051  -3.194   0.0014 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 2122.73  on 58  degrees of freedom
Residual deviance:  559.44  on 55  degrees of freedom
AIC: 850.71

Number of Fisher Scoring iterations: 5
detach(breslow.dat)
LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQpgYGB7cn0Kcm0obGlzdD1scygpKQojIyBsb2FkIGRhdGFzZXQKZGF0YShicmVzbG93LmRhdCwgcGFja2FnZSA9ICJyb2J1c3QiKQpuYW1lcyhicmVzbG93LmRhdCkKc3VtbWFyeShicmVzbG93LmRhdFtjKDYsNyw4LDEwKV0pCgojIyBsb29rIGF0IHRoZSByZXNwb25zZSB2YXJpYWJsZSBpbiBtb3JlIGRldGFpbApvcGFyIDwtIHBhcihuby5yZWFkb25seSA9IFRSVUUpCnBhcihtZnJvdyA9IGMoMSwyKSkKYXR0YWNoKGJyZXNsb3cuZGF0KQpoaXN0KHN1bVksIGJyZWFrcyA9IDIwLCB4bGFiID0gIlNlaXp1cmUgQ291bnQiLCBtYWluID0gIkRpc3RyaWJ1dGlvbiBvZiBTZWl6dXJlcyIpCmJveHBsb3Qoc3VtWSB+IFRydCwgeGxhYiA9ICJUcmVhdG1lbnQiLCBtYWluID0gIkdyb3VwIENvbXBhcmlzb25zIikKcGFyKG9wYXIpCgojIyBmcm9tIGFib3ZlIGZpZ3VyZSwgd2UgY2FuIHNlZSB0aGUgc2tld2VkIG5hdHVyZSBvZiB0aGUgZGVwZW5kZW50IHZhcmlhYmxlIGFuZCB0aGUgcG9zc2libGUgcHJlc2VuY2Ugb2Ygb3V0bGllcnMuCgojIyBmaXQgdGhlIHBvaXNzb24gcmVncmVzc2lvbgpmaXQgPC0gZ2xtKHN1bVkgfiBCYXNlICsgQWdlICsgVHJ0LCBkYXRhID0gYnJlc2xvdy5kYXQsIGZhbWlseSA9IHBvaXNzb24oKSkKc3VtbWFyeShmaXQpCmRldGFjaChicmVzbG93LmRhdCkKYGBg