source("http://192.38.117.59/~linearpredictors/datafiles/readFever.R")

## Intercept corresponds to mage=30

model <- glm(episodes ~ parity + I(mage - 30),
             family = poisson, data = fever)

data <- subset(fever, !is.na(mage) & !is.na(parity) & !is.na(episodes))
## reschi?

## qroups according to the deciles of the fitted probabilities (= expected no. of episodes)
#data$episodegroup <- cut(model$fitted,
#                         breaks =quantile(model$fitted, probs = seq(0,1,length.out=11)),
#                         labels=1:10)

## number of women:
#table(data$predictedgroup)
#O <- tapply(data$episodes, data$predictedgroup, sum)
#E <- tapply(model$fitted, data$predictedgroup, sum)

#(O-E)/sqrt(E)

(deviation <- sum(O-E)/sqrt(sum(E)))
deviation^2

## calculation of poisson probabilities for expected values

expected <- predict(model, type = "response")
E <- tapply(expected, data$fever, sum)
tapply(expected, data$fever, length)