inst/examples/Ex_glmb.wfit.R

In glmbayes: Bayesian Generalized Linear Models (IID Samples)

set.seed(333)
## Dobson (1990) Page 93: Randomized Controlled Trial :
counts <- c(18, 17, 15, 20, 10, 20, 25, 13, 12)
outcome <- gl(3, 1, 9)
treatment <- gl(3, 3)

ps <- Prior_Setup(counts ~ outcome + treatment, family = poisson())
mu <- ps$mu
V0 <- ps$Sigma
glmb.D93 <- glmb(
  counts ~ outcome + treatment,
  family = poisson(),
  pfamily = dNormal(mu = mu, Sigma = V0)
)

## Start setup here [First output from earlier optim optimization]
betastar <- glmb.D93$coef.mode # Posterior mode from optim
x <- glmb.D93$x
y <- glmb.D93$y
# The fitted object does not currently store an offset, so use zeros here.
offset2 <- 0 * y # Should return this from lower level functions
weights2 <- glmb.D93$prior.weights

## Check influence measures for original model
fit <- glmb.wfit(x, y, weights2, offset2, family = poisson(), Bbar = mu, P = solve(V0), betastar)
influence.measures(fit)

print(fit)
print(glmb.D93$coef.mode)

### Now try a strong prior with poorly chosen intercept
mu1 <- 0 * mu
V1 <- 0.1 * V0
glmb2.D93 <- glmb(
  counts ~ outcome + treatment,
  family = poisson(),
  pfamily = dNormal(mu = mu1, Sigma = V1)
)

Bbar2 <- mu1 # Prior mean
betastar2 <- glmb2.D93$coef.mode # Posterior mode from optim
fit2 <- glmb.wfit(x, y, weights2, offset2, family = poisson(), Bbar2, P = solve(V1), betastar2)

influence.measures(fit2)

print(fit2)
print(glmb2.D93$coef.mode)



influence(glmb2.D93)
influence(glmb2.D93, do.coef = TRUE)
influence(glmb2.D93, do.coef = FALSE)

glmb.influence.measures(glmb2.D93, influence(glmb2.D93))


## Do not need method functions
dfbeta(glmb2.D93, influence(glmb2.D93))
dfbeta(glmb2.D93$fit, influence(glmb2.D93))
hatvalues(glmb2.D93, influence(glmb2.D93))
hatvalues(glmb2.D93$fit, influence(glmb2.D93))

# Implemented methods

rstandard(glmb2.D93, infl = influence(glmb2.D93))
rstandard(glmb2.D93)


# Methods requiring dedicated handling

## Needs a method function
dfbetas(glmb2.D93)
dfbetas(glmb2.D93, influence(glmb2.D93, do.coef = TRUE))
dfbetas(glmb2.D93$fit, influence(glmb2.D93))



# Needs a method function
cooks.distance(glmb2.D93)
cooks.distance(glmb2.D93$fit, influence(glmb2.D93))



# The rstandard method now works directly on glmb objects.


# Needs a method function
rstudent(glmb2.D93)
rstudent(glmb2.D93$fit, influence(glmb2.D93))


# Not a method, separate function
glmb.dffits(glmb2.D93)
dffits(glmb2.D93$fit, influence(glmb2.D93))

# Not a method - separate function
glmb.covratio(glmb2.D93)
covratio(glmb2.D93$fit, influence(glmb2.D93))



# Needs a method function but requires a different approach
# The influence function actually stores this measure
## This seems to require more work to create a method function

infl <- influence(glmb2.D93)
hat2 <- infl$hat
hat2