## time poison treat
## Min. :0.1800 I :16 A:12
## 1st Qu.:0.3000 II :16 B:12
## Median :0.4000 III:16 C:12
## Mean :0.4794 D:12
## 3rd Qu.:0.6225
## Max. :1.2400
## Start: AIC=-129.85
## log(time) ~ (poison + treat)^2
##
## Df Sum of Sq RSS AIC
## - poison:treat 6 0.39575 2.3423 -132.96
## <none> 1.9465 -129.85
##
## Step: AIC=-132.96
## log(time) ~ poison + treat
##
## Df Sum of Sq RSS AIC
## <none> 2.3423 -132.964
## - treat 3 3.5572 5.8994 -94.625
## - poison 2 5.2375 7.5797 -80.595
##
## Call:
## lm(formula = log(time) ~ poison + treat, data = rats)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.41827 -0.13020 -0.01135 0.10405 0.58670
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.88731 0.08349 -10.627 1.77e-13 ***
## poisonII -0.18666 0.08349 -2.236 0.0307 *
## poisonIII -0.77515 0.08349 -9.284 9.82e-12 ***
## treatB 0.70465 0.09641 7.309 5.27e-09 ***
## treatC 0.19671 0.09641 2.040 0.0476 *
## treatD 0.50707 0.09641 5.260 4.57e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2362 on 42 degrees of freedom
## Multiple R-squared: 0.7897, Adjusted R-squared: 0.7646
## F-statistic: 31.54 on 5 and 42 DF, p-value: 3.401e-13
## Start: AIC=-74.11
## time ~ (poison + treat)^2
##
## Df Deviance AIC
## - poison:treat 6 2.4036 -76.840
## <none> 1.9205 -74.113
##
## Step: AIC=-75.26
## time ~ poison + treat
##
## Df Deviance AIC
## <none> 2.4036 -75.262
## - treat 3 6.1046 -20.866
## - poison 2 7.5009 3.922
##
## Call:
## glm(formula = time ~ poison + treat, family = Gamma(link = log),
## data = rats)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.42216 -0.16456 -0.02846 0.10147 0.58524
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.87536 0.08752 -10.002 1.12e-12 ***
## poisonII -0.16096 0.08752 -1.839 0.0730 .
## poisonIII -0.78792 0.08752 -9.003 2.34e-11 ***
## treatB 0.72218 0.10106 7.146 8.99e-09 ***
## treatC 0.19855 0.10106 1.965 0.0561 .
## treatD 0.52281 0.10106 5.173 6.05e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.06127824)
##
## Null deviance: 11.5710 on 47 degrees of freedom
## Residual deviance: 2.4036 on 42 degrees of freedom
## AIC: -75.262
##
## Number of Fisher Scoring iterations: 5
## Gamma Dispersion: 0.04966362
##
## Call:
## glm(formula = time ~ .^2, family = inverse.gaussian(link = "identity"),
## data = rats)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.6987 -0.2145 0.0238 0.1703 0.5229
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.41250 0.04062 10.155 4.12e-12 ***
## poisonII -0.09250 0.04920 -1.880 0.0682 .
## poisonIII -0.20250 0.04322 -4.685 3.91e-05 ***
## treatB 0.46750 0.13293 3.517 0.0012 **
## treatC 0.15500 0.07712 2.010 0.0520 .
## treatD 0.19750 0.08359 2.363 0.0237 *
## poisonII:treatB 0.02750 0.17655 0.156 0.8771
## poisonIII:treatB -0.34250 0.13701 -2.500 0.0171 *
## poisonII:treatC -0.10000 0.08920 -1.121 0.2697
## poisonIII:treatC -0.13000 0.08044 -1.616 0.1148
## poisonII:treatD 0.15000 0.12145 1.235 0.2248
## poisonIII:treatD -0.08250 0.08951 -0.922 0.3628
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for inverse.gaussian family taken to be 0.09403976)
##
## Null deviance: 25.7437 on 47 degrees of freedom
## Residual deviance: 3.6418 on 36 degrees of freedom
## AIC: -84.797
##
## Number of Fisher Scoring iterations: 3
## Gamma Dispersion: 0.03974636
Code used in analysis
knitr::opts_chunk$set(
echo = FALSE,
message = FALSE,
warning = FALSE
)
#ELMR 7.2
library(MASS)
library(faraway)
library(nnet)
data(rats)
summary(rats)
par(mfrow=(c(1,2)))
plot(time~.^2,rats)
llmd<-lm(log(time)~.^2,rats)
rlmdl<-step(llmd)
summary(rlmdl)
gmdl<-glm(time~.^2,family=Gamma(link=log),rats)
rgmdl<-step(gmdl)
summary(rgmdl)
cat("Gamma Dispersion:", gamma.dispersion(rgmdl))
igmdl<-glm(time~.^2,family=inverse.gaussian(link="identity"),rats)
summary(igmdl)
cat("Gamma Dispersion:", gamma.dispersion(igmdl))
par(mfrow=(c(1,2)))
plot(residuals(gmdl)~log(fitted(gmdl)),ylab="Deviance residuals",xlab=expression(log(hat(mu))))
abline(h=0)
plot(residuals(igmdl)~log(fitted(igmdl)),ylab="Inverse Deviance residuals",xlab=expression(log(hat(mu))))
abline(h=0)