| ###-- Linear Models, basic functionality -- weights included. |
| |
| ## From John Maindonald : |
| roller <- data.frame( |
| weight = c(1.9, 3.1, 3.3, 4.8, 5.3, 6.1, 6.4, 7.6, 9.8, 12.4), |
| depression = c( 2, 1, 5, 5, 20, 20, 23, 10, 30, 25)) |
| |
| roller.lmu <- lm(weight~depression, data=roller) |
| roller.lsfu <- lsfit(roller$depression, roller$weight) |
| |
| roller.lsf <- lsfit(roller$depression, roller$weight, wt = 1:10) |
| roller.lsf0 <- lsfit(roller$depression, roller$weight, wt = 0:9) |
| roller.lm <- lm(weight~depression, data=roller, weights= 1:10) |
| roller.lm0 <- lm(weight~depression, data=roller, weights= 0:9) |
| roller.lm9 <- lm(weight~depression, data=roller[-1,],weights= 1:9) |
| roller.glm <- glm(weight~depression, data=roller, weights= 1:10) |
| roller.glm0<- glm(weight~depression, data=roller, weights= 0:9) |
| |
| predict(roller.glm0, type="terms")# failed till 2003-03-31 |
| |
| ## FIXME : glm()$residual [1] is NA, lm()'s is ok. |
| ## all.equal(residuals(roller.glm0, type = "partial"), |
| ## residuals(roller.lm0, type = "partial") ) |
| |
| |
| all.equal(deviance(roller.lm), |
| deviance(roller.glm)) |
| all.equal(weighted.residuals(roller.lm), |
| residuals (roller.glm)) |
| |
| all.equal(deviance(roller.lm0), |
| deviance(roller.glm0)) |
| all.equal(weighted.residuals(roller.lm0, drop=FALSE), |
| residuals (roller.glm0)) |
| |
| (im.lm0 <- influence.measures(roller.lm0)) |
| |
| all.equal(unname(im.lm0 $ infmat), |
| unname(cbind( dfbetas (roller.lm0) |
| , dffits (roller.lm0) |
| , covratio (roller.lm0) |
| ,cooks.distance(roller.lm0) |
| ,lm.influence (roller.lm0)$hat) |
| )) |
| |
| all.equal(rstandard(roller.lm9), |
| rstandard(roller.lm0),tolerance = 1e-14) |
| all.equal(rstudent(roller.lm9), |
| rstudent(roller.lm0),tolerance = 1e-14) |
| all.equal(rstudent(roller.lm), |
| rstudent(roller.glm)) |
| all.equal(cooks.distance(roller.lm), |
| cooks.distance(roller.glm)) |
| |
| |
| all.equal(summary(roller.lm0)$coefficients, |
| summary(roller.lm9)$coefficients, tolerance = 1e-14) |
| all.equal(print(anova(roller.lm0), signif.st=FALSE), |
| anova(roller.lm9), tolerance = 1e-14) |
| |
| |
| ### more regression tests for lm(), glm(), etc : |
| |
| ## moved from ?influence.measures: |
| lm.SR <- lm(sr ~ pop15 + pop75 + dpi + ddpi, data = LifeCycleSavings) |
| (IM <- influence.measures(lm.SR)) |
| summary(IM) |
| ## colnames will differ in the next line |
| all.equal(dfbetas(lm.SR), IM$infmat[, 1:5], check.attributes = FALSE, |
| tolerance = 1e-12) |
| signif(dfbeta(lm.SR), 3) |
| covratio (lm.SR) |
| |
| ## Multivariate lm ("mlm") --- Example from ?SSD |
| reacttime <- matrix(c( |
| 420, 420, 480, 480, 600, 780, |
| 420, 480, 480, 360, 480, 600, |
| 480, 480, 540, 660, 780, 780, |
| 420, 540, 540, 480, 780, 900, |
| 540, 660, 540, 480, 660, 720, |
| 360, 420, 360, 360, 480, 540, |
| 480, 480, 600, 540, 720, 840, |
| 480, 600, 660, 540, 720, 900, |
| 540, 600, 540, 480, 720, 780, |
| 480, 420, 540, 540, 660, 780), |
| ncol = 6, byrow = TRUE, |
| dimnames = list(subj = 1:10, |
| cond = c("deg0NA", "deg4NA", "deg8NA", |
| "deg0NP", "deg4NP", "deg8NP"))) |
| mlmfit <- lm(reacttime ~ 1) |
| ImMLM <- influence.measures(mlmfit)## fails in R <= 3.5.1 |
| ## and the print() and summary() methods had failed additionally: |
| oo <- capture.output(ImMLM) # now ok |
| summary(ImMLM) # "ok" |
| |
| |
| |
| ## predict.lm(.) |
| |
| all.equal(predict(roller.lm, se.fit=TRUE)$se.fit, |
| predict(roller.lm, newdata=roller, se.fit=TRUE)$se.fit, tolerance = 1e-14) |
