ANLY 545 HW6

library(vcdExtra)

## Loading required package: vcd

## Loading required package: grid

## Loading required package: gnm

data("Caesar",package="vcdExtra")
Caesar.df <- as.data.frame(Caesar)
Caesar.df$Infect <- as.numeric(Caesar.df$Infection %in% c("Type 1", "Type 2"))

1. Fit the main-effects logit model for the binary response Infect. Note that with the data in the form of a frequency data frame you will need to use weights=Freq in the call to glm(). (It might also be convenient to reorder the levels of the factors so that “No” is the baseline level for each.)

Caesar.df$Risk <- factor(Caesar.df$Risk, levels(Caesar.df$Risk)[c(2,1)])
Caesar.df$Antibiotics <- factor(Caesar.df$Antibiotics, levels(Caesar.df$Antibiotics)[c(2,1)])
Caesar.df$Planned <- factor(Caesar.df$Planned, levels(Caesar.df$Planned)[c(2,1)])
model <- glm(Infect ~ Risk + Antibiotics + Planned, family = binomial, data = Caesar.df, weights = Freq)

2. Use summary () or car::Anova () to test the terms in this model.

summary(model)

## 
## Call:
## glm(formula = Infect ~ Risk + Antibiotics + Planned, family = binomial, 
##     data = Caesar.df, weights = Freq)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -6.7471  -0.4426   0.0000   3.2338   5.4201  
## 
## Coefficients:
##                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)     -0.7935     0.4785  -1.658   0.0972 .  
## RiskYes          1.8270     0.4364   4.186 2.84e-05 ***
## AntibioticsYes  -3.0011     0.4593  -6.535 6.37e-11 ***
## PlannedYes      -0.9064     0.4084  -2.219   0.0265 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 300.85  on 16  degrees of freedom
## Residual deviance: 236.36  on 13  degrees of freedom
## AIC: 244.36
## 
## Number of Fisher Scoring iterations: 6

anova(model)

## Analysis of Deviance Table
## 
## Model: binomial, link: logit
## 
## Response: Infect
## 
## Terms added sequentially (first to last)
## 
## 
##             Df Deviance Resid. Df Resid. Dev
## NULL                           16     300.85
## Risk         1    4.104        15     296.75
## Antibiotics  1   55.163        14     241.59
## Planned      1    5.230        13     236.36

3. Interpret the coefficients in the fitted model in terms of their effect on the odds of infection.

effect <- exp(coef(model)) - 1
effect.perc <- paste(round(100*effect, 2), "%", sep="")
effect.perc.abs <- paste(round(100*abs(effect), 2), "%", sep="")
effect.perc

## [1] "-54.77%" "521.52%" "-95.03%" "-59.6%"

The odds of infection decreases by 95% when antibiotics is used and decreases by 60% when C-section is planned.

4. Make one or more effects plots for this model, showing separate terms, or their combinations.

library(effects)

## Loading required package: carData

## 
## Attaching package: 'carData'

## The following object is masked from 'package:vcdExtra':
## 
##     Burt

## lattice theme set by effectsTheme()
## See ?effectsTheme for details.

Caesar.effect <- allEffects(model)
plot(Caesar.effect)

model2 <- update(model, . ~ . + Antibiotics:Planned)
plot(allEffects(model2))