R Under development (unstable) (2022-12-20 r83481) -- "Unsuffered Consequences" Copyright (C) 2022 The R Foundation for Statistical Computing Platform: x86_64-pc-linux-gnu (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > ###-- 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 depression 1 -2.6692211 2 -2.8898179 3 -2.0074308 4 -2.0074308 5 1.3015211 6 1.3015211 7 1.9633114 8 -0.9044468 9 3.5074889 10 2.4045050 attr(,"constant") [1] 6.743646 > > stopifnot(exprs = { + all.equal(residuals(roller.glm0, type = "partial"), + residuals(roller.lm0, type = "partial"), tol = 1e-14) # 4.0e-16 + all.equal(deviance(roller.lm), + deviance(roller.glm), tol = 1e-14) # 2.4e-16 + all.equal(weighted.residuals(roller.lm), + residuals (roller.glm), tol = 2e-14) # 9.17e-16 + all.equal(deviance(roller.lm0), + deviance(roller.glm0), tol = 1e-14) # 2.78e-16 + all.equal(weighted.residuals(roller.lm0, drop=FALSE), + residuals (roller.glm0), tol = 2e-14) # 6.378e-16 + }) > > (im.lm0 <- influence.measures(roller.lm0)) Influence measures of lm(formula = weight ~ depression, data = roller, weights = 0:9) : dfb.1_ dfb.dprs dffit cov.r cook.d hat inf 2 -0.0530 0.0482 -0.0530 1.551 0.00164 0.1277 3 -0.1769 0.1538 -0.1782 1.569 0.01806 0.1742 4 0.0130 -0.0113 0.0131 1.842 0.00010 0.2613 5 -0.1174 -0.0185 -0.3370 1.049 0.05552 0.0892 6 -0.1031 -0.0162 -0.2960 1.229 0.04577 0.1114 7 -0.0118 -0.1891 -0.5010 1.070 0.11932 0.1555 8 0.7572 -0.5948 0.7972 1.467 0.30998 0.3510 9 0.1225 -0.2067 -0.2664 2.391 0.04079 0.4470 * 10 -0.3348 1.1321 2.0937 0.233 0.89584 0.2826 * > > stopifnot(exprs = { + 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) + }) Analysis of Variance Table Response: weight Df Sum Sq Mean Sq F value Pr(>F) depression 1 158.41 158.408 5.4302 0.05259 Residuals 7 204.20 29.172 > > > ### 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)) Influence measures of lm(formula = sr ~ pop15 + pop75 + dpi + ddpi, data = LifeCycleSavings) : dfb.1_ dfb.pp15 dfb.pp75 dfb.dpi dfb.ddpi dffit cov.r Australia 0.01232 -0.01044 -0.02653 0.04534 -0.000159 0.0627 1.193 Austria -0.01005 0.00594 0.04084 -0.03672 -0.008182 0.0632 1.268 Belgium -0.06416 0.05150 0.12070 -0.03472 -0.007265 0.1878 1.176 Bolivia 0.00578 -0.01270 -0.02253 0.03185 0.040642 -0.0597 1.224 Brazil 0.08973 -0.06163 -0.17907 0.11997 0.068457 0.2646 1.082 Canada 0.00541 -0.00675 0.01021 -0.03531 -0.002649 -0.0390 1.328 Chile -0.19941 0.13265 0.21979 -0.01998 0.120007 -0.4554 0.655 China 0.02112 -0.00573 -0.08311 0.05180 0.110627 0.2008 1.150 Colombia 0.03910 -0.05226 -0.02464 0.00168 0.009084 -0.0960 1.167 Costa Rica -0.23367 0.28428 0.14243 0.05638 -0.032824 0.4049 0.968 Denmark -0.04051 0.02093 0.04653 0.15220 0.048854 0.3845 0.934 Ecuador 0.07176 -0.09524 -0.06067 0.01950 0.047786 -0.1695 1.139 Finland -0.11350 0.11133 0.11695 -0.04364 -0.017132 -0.1464 1.203 France -0.16600 0.14705 0.21900 -0.02942 0.023952 0.2765 1.226 Germany -0.00802 0.00822 0.00835 -0.00697 -0.000293 -0.0152 1.226 Greece -0.14820 0.16394 0.02861 0.15713 -0.059599 -0.2811 1.140 Guatamala 0.01552 -0.05485 0.00614 0.00585 0.097217 -0.2305 1.085 Honduras -0.00226 0.00984 -0.01020 0.00812 -0.001887 0.0482 1.186 Iceland 0.24789 -0.27355 -0.23265 -0.12555 0.184698 -0.4768 0.866 India 0.02105 -0.01577 -0.01439 -0.01374 -0.018958 0.0381 1.202 Ireland -0.31001 0.29624 0.48156 -0.25733 -0.093317 0.5216 1.268 Italy 0.06619 -0.07097 0.00307 -0.06999 -0.028648 0.1388 1.162 Japan 0.63987 -0.65614 -0.67390 0.14610 0.388603 0.8597 1.085 Korea -0.16897 0.13509 0.21895 0.00511 -0.169492 -0.4303 0.870 Luxembourg -0.06827 0.06888 0.04380 -0.02797 0.049134 -0.1401 1.196 Malta 0.03652 -0.04876 0.00791 -0.08659 0.153014 0.2386 1.128 Norway 0.00222 -0.00035 -0.00611 -0.01594 -0.001462 -0.0522 1.168 Netherlands 0.01395 -0.01674 -0.01186 0.00433 0.022591 0.0366 1.229 New Zealand -0.06002 0.06510 0.09412 -0.02638 -0.064740 0.1469 1.134 Nicaragua -0.01209 0.01790 0.00972 -0.00474 -0.010467 0.0397 1.174 Panama 0.02828 -0.05334 0.01446 -0.03467 -0.007889 -0.1775 1.067 Paraguay -0.23227 0.16416 0.15826 0.14361 0.270478 -0.4655 0.873 Peru -0.07182 0.14669 0.09148 -0.08585 -0.287184 0.4811 0.831 Philippines -0.15707 0.22681 0.15743 -0.11140 -0.170674 0.4884 0.818 Portugal -0.02140 0.02551 -0.00380 0.03991 -0.028011 -0.0690 1.233 South Africa 0.02218 -0.02030 -0.00672 -0.02049 -0.016326 0.0343 1.195 South Rhodesia 0.14390 -0.13472 -0.09245 -0.06956 -0.057920 0.1607 1.313 Spain -0.03035 0.03131 0.00394 0.03512 0.005340 -0.0526 1.208 Sweden 0.10098 -0.08162 -0.06166 -0.25528 -0.013316 -0.4526 1.086 Switzerland 0.04323 -0.04649 -0.04364 0.09093 -0.018828 0.1903 1.147 Turkey -0.01092 -0.01198 0.02645 0.00161 0.025138 -0.1445 1.100 Tunisia 0.07377 -0.10500 -0.07727 0.04439 0.103058 -0.2177 1.131 United Kingdom 0.04671 -0.03584 -0.17129 0.12554 0.100314 -0.2722 1.189 United States 0.06910 -0.07289 0.03745 -0.23312 -0.032729 -0.2510 1.655 Venezuela -0.05083 0.10080 -0.03366 0.11366 -0.124486 0.3071 1.095 Zambia 0.16361 -0.07917 -0.33899 0.09406 0.228232 0.7482 0.512 Jamaica 0.10958 -0.10022 -0.05722 -0.00703 -0.295461 -0.3456 1.200 Uruguay -0.13403 0.12880 0.02953 0.13132 0.099591 -0.2051 1.187 Libya 0.55074 -0.48324 -0.37974 -0.01937 -1.024477 -1.1601 2.091 Malaysia 0.03684 -0.06113 0.03235 -0.04956 -0.072294 -0.2126 1.113 cook.d hat inf Australia 8.04e-04 0.0677 Austria 8.18e-04 0.1204 Belgium 7.15e-03 0.0875 Bolivia 7.28e-04 0.0895 Brazil 1.40e-02 0.0696 Canada 3.11e-04 0.1584 Chile 3.78e-02 0.0373 * China 8.16e-03 0.0780 Colombia 1.88e-03 0.0573 Costa Rica 3.21e-02 0.0755 Denmark 2.88e-02 0.0627 Ecuador 5.82e-03 0.0637 Finland 4.36e-03 0.0920 France 1.55e-02 0.1362 Germany 4.74e-05 0.0874 Greece 1.59e-02 0.0966 Guatamala 1.07e-02 0.0605 Honduras 4.74e-04 0.0601 Iceland 4.35e-02 0.0705 India 2.97e-04 0.0715 Ireland 5.44e-02 0.2122 Italy 3.92e-03 0.0665 Japan 1.43e-01 0.2233 Korea 3.56e-02 0.0608 Luxembourg 3.99e-03 0.0863 Malta 1.15e-02 0.0794 Norway 5.56e-04 0.0479 Netherlands 2.74e-04 0.0906 New Zealand 4.38e-03 0.0542 Nicaragua 3.23e-04 0.0504 Panama 6.33e-03 0.0390 Paraguay 4.16e-02 0.0694 Peru 4.40e-02 0.0650 Philippines 4.52e-02 0.0643 Portugal 9.73e-04 0.0971 South Africa 2.41e-04 0.0651 South Rhodesia 5.27e-03 0.1608 Spain 5.66e-04 0.0773 Sweden 4.06e-02 0.1240 Switzerland 7.33e-03 0.0736 Turkey 4.22e-03 0.0396 Tunisia 9.56e-03 0.0746 United Kingdom 1.50e-02 0.1165 United States 1.28e-02 0.3337 * Venezuela 1.89e-02 0.0863 Zambia 9.66e-02 0.0643 * Jamaica 2.40e-02 0.1408 Uruguay 8.53e-03 0.0979 Libya 2.68e-01 0.5315 * Malaysia 9.11e-03 0.0652 > summary(IM) Potentially influential observations of lm(formula = sr ~ pop15 + pop75 + dpi + ddpi, data = LifeCycleSavings) : dfb.1_ dfb.pp15 dfb.pp75 dfb.dpi dfb.ddpi dffit cov.r cook.d Chile -0.20 0.13 0.22 -0.02 0.12 -0.46 0.65_* 0.04 United States 0.07 -0.07 0.04 -0.23 -0.03 -0.25 1.66_* 0.01 Zambia 0.16 -0.08 -0.34 0.09 0.23 0.75 0.51_* 0.10 Libya 0.55 -0.48 -0.38 -0.02 -1.02_* -1.16_* 2.09_* 0.27 hat Chile 0.04 United States 0.33_* Zambia 0.06 Libya 0.53_* > ## colnames will differ in the next line > stopifnot( + all.equal(dfbetas(lm.SR), IM$infmat[, 1:5], check.attributes = FALSE, + tolerance = 1e-12) + ) > signif(dfbeta(lm.SR), 3) (Intercept) pop15 pop75 dpi ddpi Australia 0.0916 -1.53e-03 -0.02910 4.27e-05 -3.16e-05 Austria -0.0747 8.69e-04 0.04470 -3.46e-05 -1.62e-03 Belgium -0.4750 7.50e-03 0.13200 -3.26e-05 -1.44e-03 Bolivia 0.0429 -1.86e-03 -0.02470 3.00e-05 8.06e-03 Brazil 0.6600 -8.92e-03 -0.19400 1.12e-04 1.34e-02 Canada 0.0402 -9.87e-04 0.01120 -3.32e-05 -5.25e-04 Chile -1.4000 1.83e-02 0.22700 -1.78e-05 2.25e-02 China 0.1560 -8.33e-04 -0.09060 4.85e-05 2.18e-02 Colombia 0.2900 -7.63e-03 -0.02700 1.58e-06 1.80e-03 Costa Rica -1.7000 4.07e-02 0.15300 5.19e-05 -6.37e-03 Denmark -0.2940 2.99e-03 0.04980 1.40e-04 9.46e-03 Ecuador 0.5310 -1.39e-02 -0.06620 1.83e-05 9.44e-03 Finland -0.8420 1.62e-02 0.12800 -4.10e-05 -3.39e-03 France -1.2300 2.14e-02 0.23900 -2.76e-05 4.73e-03 Germany -0.0597 1.20e-03 0.00915 -6.56e-06 -5.82e-05 Greece -1.0900 2.38e-02 0.03110 1.47e-04 -1.17e-02 Guatamala 0.1140 -7.95e-03 0.00667 5.46e-06 1.91e-02 Honduras -0.0168 1.44e-03 -0.01120 7.64e-06 -3.74e-04 Iceland 1.7800 -3.87e-02 -0.24700 -1.14e-04 3.55e-02 India 0.1560 -2.31e-03 -0.01580 -1.29e-05 -3.76e-03 Ireland -2.2800 4.28e-02 0.52200 -2.40e-04 -1.83e-02 Italy 0.4910 -1.04e-02 0.00335 -6.57e-05 -5.67e-03 Japan 4.6300 -9.33e-02 -0.71800 1.34e-04 7.49e-02 Korea -1.2200 1.91e-02 0.23300 4.66e-06 -3.26e-02 Luxembourg -0.5070 1.01e-02 0.04790 -2.63e-05 9.73e-03 Malta 0.2700 -7.08e-03 0.00861 -8.09e-05 3.01e-02 Norway 0.0165 -5.12e-05 -0.00670 -1.50e-05 -2.90e-04 Netherlands 0.1040 -2.45e-03 -0.01300 4.07e-06 4.48e-03 New Zealand -0.4440 9.48e-03 0.10300 -2.47e-05 -1.28e-02 Nicaragua -0.0899 2.62e-03 0.01060 -4.46e-06 -2.08e-03 Panama 0.2080 -7.73e-03 0.01570 -3.24e-05 -1.55e-03 Paraguay -1.6700 2.33e-02 0.16800 1.31e-04 5.20e-02 Peru -0.5150 2.07e-02 0.09670 -7.79e-05 -5.49e-02 Philippines -1.1200 3.19e-02 0.16600 -1.01e-04 -3.26e-02 Portugal -0.1590 3.73e-03 -0.00416 3.76e-05 -5.55e-03 South Africa 0.1650 -2.97e-03 -0.00737 -1.93e-05 -3.24e-03 South Rhodesia 1.0700 -1.97e-02 -0.10100 -6.54e-05 -1.15e-02 Spain -0.2260 4.58e-03 0.00432 3.31e-05 1.06e-03 Sweden 0.7390 -1.17e-02 -0.06650 -2.37e-04 -2.60e-03 Switzerland 0.3200 -6.77e-03 -0.04760 8.52e-05 -3.72e-03 Turkey -0.0807 -1.74e-03 0.02880 1.51e-06 4.96e-03 Tunisia 0.5450 -1.53e-02 -0.08410 4.15e-05 2.03e-02 United Kingdom 0.3450 -5.21e-03 -0.18700 1.17e-04 1.98e-02 United States 0.5130 -1.06e-02 0.04100 -2.19e-04 -6.48e-03 Venezuela -0.3740 1.46e-02 -0.03650 1.06e-04 -2.44e-02 Zambia 1.1200 -1.06e-02 -0.34100 8.14e-05 4.16e-02 Jamaica 0.8080 -1.45e-02 -0.06220 -6.57e-06 -5.81e-02 Uruguay -0.9920 1.88e-02 0.03220 1.23e-04 1.97e-02 Libya 4.0400 -6.98e-02 -0.41100 -1.80e-05 -2.01e-01 Malaysia 0.2720 -8.88e-03 0.03520 -4.63e-05 -1.42e-02 > covratio (lm.SR) Australia Austria Belgium Bolivia Brazil 1.1928303 1.2678392 1.1761879 1.2238199 1.0823332 Canada Chile China Colombia Costa Rica 1.3283009 0.6547098 1.1498637 1.1666845 0.9681384 Denmark Ecuador Finland France Germany 0.9344047 1.1393880 1.2031561 1.2262654 1.2256855 Greece Guatamala Honduras Iceland India 1.1396174 1.0852720 1.1855450 0.8658808 1.2024438 Ireland Italy Japan Korea Luxembourg 1.2680432 1.1624611 1.0845999 0.8695843 1.1961844 Malta Norway Netherlands New Zealand Nicaragua 1.1282611 1.1680616 1.2285315 1.1336998 1.1742677 Panama Paraguay Peru Philippines Portugal 1.0667255 0.8732040 0.8312741 0.8177726 1.2331038 South Africa South Rhodesia Spain Sweden Switzerland 1.1945449 1.3130954 1.2081541 1.0864869 1.1471125 Turkey Tunisia United Kingdom United States Venezuela 1.1003557 1.1314365 1.1886236 1.6554816 1.0945955 Zambia Jamaica Uruguay Libya Malaysia 0.5116454 1.1995171 1.1872025 2.0905736 1.1126445 > > ## 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" Potentially influential observations of lm(formula = reacttime ~ 1) : , , dfb.1_ deg0NA deg4NA deg8NA deg0NP deg4NP deg8NP 1 -0.25 -0.37 -0.21 -0.04 -0.18 0.05 2 -0.25 -0.12 -0.21 -0.58 -0.67 -0.53 3 0.11 -0.12 0.05 0.84 0.39 0.05 4 -0.25 0.12 0.05 -0.04 0.39 0.43 5 0.52 0.73 0.05 -0.04 0.00 -0.12 6 -0.76 -0.37 -1.06_* -0.58 -0.67 -0.85 7 0.11 -0.12 0.32 0.18 0.18 0.23 8 0.11 0.37 0.68 0.18 0.18 0.43 9 0.52 0.37 0.05 -0.04 0.18 0.05 10 0.11 -0.37 0.05 0.18 0.00 0.05 , , dffit deg0NA deg4NA deg8NA deg0NP deg4NP deg8NP 1 -0.25 -0.37 -0.21 -0.04 -0.18 0.05 2 -0.25 -0.12 -0.21 -0.58 -0.67 -0.53 3 0.11 -0.12 0.05 0.84 0.39 0.05 4 -0.25 0.12 0.05 -0.04 0.39 0.43 5 0.52 0.73 0.05 -0.04 0.00 -0.12 6 -0.76 -0.37 -1.06_* -0.58 -0.67 -0.85 7 0.11 -0.12 0.32 0.18 0.18 0.23 8 0.11 0.37 0.68 0.18 0.18 0.43 9 0.52 0.37 0.05 -0.04 0.18 0.05 10 0.11 -0.37 0.05 0.18 0.00 0.05 , , cov.r deg0NA deg4NA deg8NA deg0NP deg4NP deg8NP 1 0.08_* 0.16_* 0.15_* 0.19_* 0.29_* 0.34_* 2 0.08_* 0.18_* 0.15_* 0.14_* 0.20_* 0.26_* 3 0.08_* 0.18_* 0.15_* 0.11_* 0.25_* 0.34_* 4 0.08_* 0.18_* 0.15_* 0.19_* 0.25_* 0.28_* 5 0.06_* 0.11_* 0.15_* 0.19_* 0.30_* 0.33_* 6 0.05_* 0.16_* 0.07_* 0.14_* 0.20_* 0.19_* 7 0.08_* 0.18_* 0.14_* 0.19_* 0.29_* 0.32_* 8 0.08_* 0.16_* 0.10_* 0.19_* 0.29_* 0.28_* 9 0.06_* 0.16_* 0.15_* 0.19_* 0.29_* 0.34_* 10 0.08_* 0.16_* 0.15_* 0.19_* 0.30_* 0.34_* , , cook.d deg0NA deg4NA deg8NA deg0NP deg4NP deg8NP 1 0.00 0.02 0.01 0.00 0.01 0.00 2 0.00 0.00 0.01 0.04 0.08 0.06 3 0.00 0.00 0.00 0.07 0.04 0.00 4 0.00 0.00 0.00 0.00 0.04 0.05 5 0.01 0.06 0.00 0.00 0.00 0.00 6 0.03 0.02 0.07 0.04 0.08 0.12 7 0.00 0.00 0.01 0.01 0.01 0.01 8 0.00 0.02 0.04 0.01 0.01 0.05 9 0.01 0.02 0.00 0.00 0.01 0.00 10 0.00 0.02 0.00 0.01 0.00 0.00 , , hat deg0NA deg4NA deg8NA deg0NP deg4NP deg8NP 1 0.10 0.10 0.10 0.10 0.10 0.10 2 0.10 0.10 0.10 0.10 0.10 0.10 3 0.10 0.10 0.10 0.10 0.10 0.10 4 0.10 0.10 0.10 0.10 0.10 0.10 5 0.10 0.10 0.10 0.10 0.10 0.10 6 0.10 0.10 0.10 0.10 0.10 0.10 7 0.10 0.10 0.10 0.10 0.10 0.10 8 0.10 0.10 0.10 0.10 0.10 0.10 9 0.10 0.10 0.10 0.10 0.10 0.10 10 0.10 0.10 0.10 0.10 0.10 0.10 > > > > ## 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) [1] TRUE >