(a) Fit the main-effects logit model
library(vcdExtra)
## Loading required package: vcd
## Loading required package: grid
## Loading required package: gnm
data ("Caesar", package = "vcdExtra")
summary(Caesar)
## Number of cases in table: 251
## Number of factors: 4
## Test for independence of all factors:
## Chisq = 233.12, df = 18, p-value = 2.172e-39
## Chi-squared approximation may be incorrect
Caesar.df <- as.data.frame(Caesar)
Caesar.df$Infect <- as.numeric(Caesar.df$Infection %in% c("Type1", "Type2"))
str(Caesar.df)
## 'data.frame': 24 obs. of 6 variables:
## $ Infection : Factor w/ 3 levels "Type 1","Type 2",..: 1 2 3 1 2 3 1 2 3 1 ...
## $ Risk : Factor w/ 2 levels "Yes","No": 1 1 1 2 2 2 1 1 1 2 ...
## $ Antibiotics: Factor w/ 2 levels "Yes","No": 1 1 1 1 1 1 2 2 2 2 ...
## $ Planned : Factor w/ 2 levels "Yes","No": 1 1 1 1 1 1 1 1 1 1 ...
## $ Freq : num 0 1 17 0 1 1 11 17 30 4 ...
## $ Infect : num 0 0 0 0 0 0 0 0 0 0 ...
Based on the summary table above, three variabes : Risk, Antibistics and Planned can be reordered.
Caesar.df$Risk <- factor(Caesar.df$Risk, levels(Caesar.df$Risk)[c(2,1)])
Caesar.df$Planned <- factor(Caesar.df$Planned, levels(Caesar.df$Planned)[c(2,1)])
Caesar.df$Antibiotics <- factor(Caesar.df$Antibiotics, levels(Caesar.df$Antibiotics)[c(2,1)])
str(Caesar.df)
## 'data.frame': 24 obs. of 6 variables:
## $ Infection : Factor w/ 3 levels "Type 1","Type 2",..: 1 2 3 1 2 3 1 2 3 1 ...
## $ Risk : Factor w/ 2 levels "No","Yes": 2 2 2 1 1 1 2 2 2 1 ...
## $ Antibiotics: Factor w/ 2 levels "No","Yes": 2 2 2 2 2 2 1 1 1 1 ...
## $ Planned : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 2 2 2 2 2 ...
## $ Freq : num 0 1 17 0 1 1 11 17 30 4 ...
## $ Infect : num 0 0 0 0 0 0 0 0 0 0 ...
regression <- glm(Infect ~ Risk + Antibiotics + Planned, data = Caesar.df, family = binomial, weights=Freq)
(b) test the terms in this model
car::Anova (regression)
## Analysis of Deviance Table (Type II tests)
##
## Response: Infect
## LR Chisq Df Pr(>Chisq)
## Risk 4.3272e-12 1 1
## Antibiotics -7.3439e-12 1 1
## Planned -5.0940e-13 1 1
summary(regression)
##
## Call:
## glm(formula = Infect ~ Risk + Antibiotics + Planned, family = binomial,
## data = Caesar.df, weights = Freq)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.008e-05 -3.481e-06 -2.356e-06 0.000e+00 0.000e+00
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.781e+01 1.411e+05 0 1
## RiskYes -6.504e-01 1.345e+05 0 1
## AntibioticsYes 2.874e-01 1.339e+05 0 1
## PlannedYes -3.153e-02 1.250e+05 0 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 0.000e+00 on 16 degrees of freedom
## Residual deviance: 2.913e-10 on 13 degrees of freedom
## AIC: 8
##
## Number of Fisher Scoring iterations: 25
c. Interpret the coefficients
coef(regression)
## (Intercept) RiskYes AntibioticsYes PlannedYes
## -27.80625338 -0.65042629 0.28744709 -0.03152661
exp(coef(regression))
## (Intercept) RiskYes AntibioticsYes PlannedYes
## 8.392621e-13 5.218233e-01 1.333020e+00 9.689652e-01
- when risk, odds of infection will increase by 52%
- when antibiotic is taken, odds of infection will increase by 133%
- when planned, odds of infection will increase by 96%
d.
install.packages("effects",repos = "http://cran.us.r-project.org")
##
## There is a binary version available but the source version is
## later:
## binary source needs_compilation
## effects 4.0-0 4.0-1 FALSE
## installing the source package 'effects'
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.
plot(allEffects(regression), rows=1, cols=3)
