> stats:::add1.glm
function (object, scope, scale = 0, test = c("none", "Rao", "LRT",
"Chisq", "F"), x = NULL, k = 2, ...)
{
Fstat <- function(table, rdf) {
dev <- table$Deviance
df <- table$Df
diff <- pmax(0, (dev[1L] - dev)/df)
Fs <- (diff/df)/(dev/(rdf - df))
Is this where the double division is happening?
diff has df in the denominator, then Fs assignment sees
diff divided by df again. if df is 1, the double division
will go unnoticed.
Steven McKinney
Statistician
Molecular Oncology and Breast Cancer Program
British Columbia Cancer Research Centre
________________________________________
From: [email protected] [[email protected]] On Behalf
Of William Dunlap [[email protected]]
Sent: June 25, 2013 12:01 PM
To: [email protected]
Subject: [Rd] F statistic in add1.lm vs add1.glm
Should the F statistic be the same when using add1() on models created by lm
and glm(family=gaussian)?
They are in the single-degree-of-freedom case but not in the
multiple-degree-of-freedom case.
MASS:addterm shows the same discrepancy. It looks like the deviance
(==residual sum of squares) gets
divided by the number of degrees of freedom for the term twice in add1.glm.
Using anova() on the output
of lm and glm(family=gaussian) gives the save F statistic as add1.lm gives.
> # factor(carb) consumes 5 degrees of freedom, am 1, compare their F values.
> fit <- lm(mpg ~ hp, data=mtcars) ; add1(fit, ~ hp + factor(carb) + am,
> test="F")
Single term additions
Model:
mpg ~ hp
Df Sum of Sq RSS AIC F value Pr(>F)
<none> 447.67 88.427
factor(carb) 5 103.54 344.13 90.009 1.5044 0.2241
am 1 202.24 245.44 71.194 23.8952 3.46e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> fit <- glm(mpg ~ hp, data=mtcars) ; add1(fit, ~ hp + factor(carb) + am,
> test="F")
Single term additions
Model:
mpg ~ hp
Df Deviance AIC F value Pr(>F)
<none> 447.67 181.24
factor(carb) 5 344.13 182.82 0.3009 0.9077
am 1 245.44 164.01 23.8952 3.46e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
