| |
| R version 3.6.2 Patched (2020-02-12 r77795) -- "Dark and Stormy Night" |
| Copyright (C) 2020 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. |
| |
| Natural language support but running in an English locale |
| |
| 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. |
| |
| > pkgname <- "datasets" |
| > source(file.path(R.home("share"), "R", "examples-header.R")) |
| > options(warn = 1) |
| > library('datasets') |
| > |
| > base::assign(".oldSearch", base::search(), pos = 'CheckExEnv') |
| > base::assign(".old_wd", base::getwd(), pos = 'CheckExEnv') |
| > cleanEx() |
| > nameEx("AirPassengers") |
| > ### * AirPassengers |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: AirPassengers |
| > ### Title: Monthly Airline Passenger Numbers 1949-1960 |
| > ### Aliases: AirPassengers |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > ## Not run: |
| > ##D ## These are quite slow and so not run by example(AirPassengers) |
| > ##D |
| > ##D ## The classic 'airline model', by full ML |
| > ##D (fit <- arima(log10(AirPassengers), c(0, 1, 1), |
| > ##D seasonal = list(order = c(0, 1, 1), period = 12))) |
| > ##D update(fit, method = "CSS") |
| > ##D update(fit, x = window(log10(AirPassengers), start = 1954)) |
| > ##D pred <- predict(fit, n.ahead = 24) |
| > ##D tl <- pred$pred - 1.96 * pred$se |
| > ##D tu <- pred$pred + 1.96 * pred$se |
| > ##D ts.plot(AirPassengers, 10^tl, 10^tu, log = "y", lty = c(1, 2, 2)) |
| > ##D |
| > ##D ## full ML fit is the same if the series is reversed, CSS fit is not |
| > ##D ap0 <- rev(log10(AirPassengers)) |
| > ##D attributes(ap0) <- attributes(AirPassengers) |
| > ##D arima(ap0, c(0, 1, 1), seasonal = list(order = c(0, 1, 1), period = 12)) |
| > ##D arima(ap0, c(0, 1, 1), seasonal = list(order = c(0, 1, 1), period = 12), |
| > ##D method = "CSS") |
| > ##D |
| > ##D ## Structural Time Series |
| > ##D ap <- log10(AirPassengers) - 2 |
| > ##D (fit <- StructTS(ap, type = "BSM")) |
| > ##D par(mfrow = c(1, 2)) |
| > ##D plot(cbind(ap, fitted(fit)), plot.type = "single") |
| > ##D plot(cbind(ap, tsSmooth(fit)), plot.type = "single") |
| > ## End(Not run) |
| > |
| > |
| > cleanEx() |
| > nameEx("BOD") |
| > ### * BOD |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: BOD |
| > ### Title: Biochemical Oxygen Demand |
| > ### Aliases: BOD |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > ## Don't show: |
| > options(show.nls.convergence=FALSE) |
| > old <- options(digits = 5) |
| > ## End(Don't show) |
| > require(stats) |
| > # simplest form of fitting a first-order model to these data |
| > fm1 <- nls(demand ~ A*(1-exp(-exp(lrc)*Time)), data = BOD, |
| + start = c(A = 20, lrc = log(.35))) |
| > coef(fm1) |
| A lrc |
| 19.14258 -0.63282 |
| > fm1 |
| Nonlinear regression model |
| model: demand ~ A * (1 - exp(-exp(lrc) * Time)) |
| data: BOD |
| A lrc |
| 19.143 -0.633 |
| residual sum-of-squares: 26 |
| > # using the plinear algorithm |
| > fm2 <- nls(demand ~ (1-exp(-exp(lrc)*Time)), data = BOD, |
| + start = c(lrc = log(.35)), algorithm = "plinear", trace = TRUE) |
| 32.946 : -1.0498 22.1260 |
| 25.992 : -0.62572 19.10319 |
| 25.99 : -0.6327 19.1419 |
| 25.99 : -0.63282 19.14256 |
| > # using a self-starting model |
| > fm3 <- nls(demand ~ SSasympOrig(Time, A, lrc), data = BOD) |
| > summary(fm3) |
| |
| Formula: demand ~ SSasympOrig(Time, A, lrc) |
| |
| Parameters: |
| Estimate Std. Error t value Pr(>|t|) |
| A 19.143 2.496 7.67 0.0016 ** |
| lrc -0.633 0.382 -1.65 0.1733 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 2.55 on 4 degrees of freedom |
| |
| > ## Don't show: |
| > options(old) |
| > ## End(Don't show) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("ChickWeight") |
| > ### * ChickWeight |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: ChickWeight |
| > ### Title: Weight versus age of chicks on different diets |
| > ### Aliases: ChickWeight |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > |
| > cleanEx() |
| > nameEx("DNase") |
| > ### * DNase |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: DNase |
| > ### Title: Elisa assay of DNase |
| > ### Aliases: DNase |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## Don't show: |
| > options(show.nls.convergence=FALSE) |
| > ## End(Don't show) |
| > coplot(density ~ conc | Run, data = DNase, |
| + show.given = FALSE, type = "b") |
| > coplot(density ~ log(conc) | Run, data = DNase, |
| + show.given = FALSE, type = "b") |
| > ## fit a representative run |
| > fm1 <- nls(density ~ SSlogis( log(conc), Asym, xmid, scal ), |
| + data = DNase, subset = Run == 1) |
| > ## compare with a four-parameter logistic |
| > fm2 <- nls(density ~ SSfpl( log(conc), A, B, xmid, scal ), |
| + data = DNase, subset = Run == 1) |
| > summary(fm2) |
| |
| Formula: density ~ SSfpl(log(conc), A, B, xmid, scal) |
| |
| Parameters: |
| Estimate Std. Error t value Pr(>|t|) |
| A -0.007897 0.017200 -0.459 0.654 |
| B 2.377239 0.109516 21.707 5.35e-11 *** |
| xmid 1.507403 0.102080 14.767 4.65e-09 *** |
| scal 1.062579 0.056996 18.643 3.16e-10 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.01981 on 12 degrees of freedom |
| |
| > anova(fm1, fm2) |
| Analysis of Variance Table |
| |
| Model 1: density ~ SSlogis(log(conc), Asym, xmid, scal) |
| Model 2: density ~ SSfpl(log(conc), A, B, xmid, scal) |
| Res.Df Res.Sum Sq Df Sum Sq F value Pr(>F) |
| 1 13 0.0047896 |
| 2 12 0.0047073 1 8.2314e-05 0.2098 0.6551 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Formaldehyde") |
| > ### * Formaldehyde |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Formaldehyde |
| > ### Title: Determination of Formaldehyde |
| > ### Aliases: Formaldehyde |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > plot(optden ~ carb, data = Formaldehyde, |
| + xlab = "Carbohydrate (ml)", ylab = "Optical Density", |
| + main = "Formaldehyde data", col = 4, las = 1) |
| > abline(fm1 <- lm(optden ~ carb, data = Formaldehyde)) |
| > summary(fm1) |
| |
| Call: |
| lm(formula = optden ~ carb, data = Formaldehyde) |
| |
| Residuals: |
| 1 2 3 4 5 6 |
| -0.006714 0.001029 0.002771 0.007143 0.007514 -0.011743 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 0.005086 0.007834 0.649 0.552 |
| carb 0.876286 0.013535 64.744 3.41e-07 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.008649 on 4 degrees of freedom |
| Multiple R-squared: 0.999, Adjusted R-squared: 0.9988 |
| F-statistic: 4192 on 1 and 4 DF, p-value: 3.409e-07 |
| |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0)) |
| > plot(fm1) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("HairEyeColor") |
| > ### * HairEyeColor |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: HairEyeColor |
| > ### Title: Hair and Eye Color of Statistics Students |
| > ### Aliases: HairEyeColor |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > ## Full mosaic |
| > mosaicplot(HairEyeColor) |
| > ## Aggregate over sex (as in Snee's original data) |
| > x <- apply(HairEyeColor, c(1, 2), sum) |
| > x |
| Eye |
| Hair Brown Blue Hazel Green |
| Black 68 20 15 5 |
| Brown 119 84 54 29 |
| Red 26 17 14 14 |
| Blond 7 94 10 16 |
| > mosaicplot(x, main = "Relation between hair and eye color") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Harman23.cor") |
| > ### * Harman23.cor |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Harman23.cor |
| > ### Title: Harman Example 2.3 |
| > ### Aliases: Harman23.cor |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > (Harman23.FA <- factanal(factors = 1, covmat = Harman23.cor)) |
| |
| Call: |
| factanal(factors = 1, covmat = Harman23.cor) |
| |
| Uniquenesses: |
| height arm.span forearm lower.leg weight |
| 0.158 0.135 0.190 0.187 0.760 |
| bitro.diameter chest.girth chest.width |
| 0.829 0.877 0.801 |
| |
| Loadings: |
| Factor1 |
| height 0.918 |
| arm.span 0.930 |
| forearm 0.900 |
| lower.leg 0.902 |
| weight 0.490 |
| bitro.diameter 0.413 |
| chest.girth 0.351 |
| chest.width 0.446 |
| |
| Factor1 |
| SS loadings 4.064 |
| Proportion Var 0.508 |
| |
| Test of the hypothesis that 1 factor is sufficient. |
| The chi square statistic is 611.44 on 20 degrees of freedom. |
| The p-value is 1.12e-116 |
| > for(factors in 2:4) print(update(Harman23.FA, factors = factors)) |
| |
| Call: |
| factanal(factors = factors, covmat = Harman23.cor) |
| |
| Uniquenesses: |
| height arm.span forearm lower.leg weight |
| 0.170 0.107 0.166 0.199 0.089 |
| bitro.diameter chest.girth chest.width |
| 0.364 0.416 0.537 |
| |
| Loadings: |
| Factor1 Factor2 |
| height 0.865 0.287 |
| arm.span 0.927 0.181 |
| forearm 0.895 0.179 |
| lower.leg 0.859 0.252 |
| weight 0.233 0.925 |
| bitro.diameter 0.194 0.774 |
| chest.girth 0.134 0.752 |
| chest.width 0.278 0.621 |
| |
| Factor1 Factor2 |
| SS loadings 3.335 2.617 |
| Proportion Var 0.417 0.327 |
| Cumulative Var 0.417 0.744 |
| |
| Test of the hypothesis that 2 factors are sufficient. |
| The chi square statistic is 75.74 on 13 degrees of freedom. |
| The p-value is 6.94e-11 |
| |
| Call: |
| factanal(factors = factors, covmat = Harman23.cor) |
| |
| Uniquenesses: |
| height arm.span forearm lower.leg weight |
| 0.127 0.005 0.193 0.157 0.090 |
| bitro.diameter chest.girth chest.width |
| 0.359 0.411 0.490 |
| |
| Loadings: |
| Factor1 Factor2 Factor3 |
| height 0.886 0.267 -0.130 |
| arm.span 0.937 0.195 0.280 |
| forearm 0.874 0.188 |
| lower.leg 0.877 0.230 -0.145 |
| weight 0.242 0.916 -0.106 |
| bitro.diameter 0.193 0.777 |
| chest.girth 0.137 0.755 |
| chest.width 0.261 0.646 0.159 |
| |
| Factor1 Factor2 Factor3 |
| SS loadings 3.379 2.628 0.162 |
| Proportion Var 0.422 0.329 0.020 |
| Cumulative Var 0.422 0.751 0.771 |
| |
| Test of the hypothesis that 3 factors are sufficient. |
| The chi square statistic is 22.81 on 7 degrees of freedom. |
| The p-value is 0.00184 |
| |
| Call: |
| factanal(factors = factors, covmat = Harman23.cor) |
| |
| Uniquenesses: |
| height arm.span forearm lower.leg weight |
| 0.137 0.005 0.191 0.116 0.138 |
| bitro.diameter chest.girth chest.width |
| 0.283 0.178 0.488 |
| |
| Loadings: |
| Factor1 Factor2 Factor3 Factor4 |
| height 0.879 0.277 -0.115 |
| arm.span 0.937 0.194 0.277 |
| forearm 0.875 0.191 |
| lower.leg 0.887 0.209 0.135 -0.188 |
| weight 0.246 0.882 0.111 -0.109 |
| bitro.diameter 0.187 0.822 |
| chest.girth 0.117 0.729 0.526 |
| chest.width 0.263 0.644 0.141 |
| |
| Factor1 Factor2 Factor3 Factor4 |
| SS loadings 3.382 2.595 0.323 0.165 |
| Proportion Var 0.423 0.324 0.040 0.021 |
| Cumulative Var 0.423 0.747 0.787 0.808 |
| |
| Test of the hypothesis that 4 factors are sufficient. |
| The chi square statistic is 4.63 on 2 degrees of freedom. |
| The p-value is 0.0988 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Harman74.cor") |
| > ### * Harman74.cor |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Harman74.cor |
| > ### Title: Harman Example 7.4 |
| > ### Aliases: Harman74.cor |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > (Harman74.FA <- factanal(factors = 1, covmat = Harman74.cor)) |
| |
| Call: |
| factanal(factors = 1, covmat = Harman74.cor) |
| |
| Uniquenesses: |
| VisualPerception Cubes PaperFormBoard |
| 0.677 0.866 0.830 |
| Flags GeneralInformation PargraphComprehension |
| 0.768 0.487 0.491 |
| SentenceCompletion WordClassification WordMeaning |
| 0.500 0.514 0.474 |
| Addition Code CountingDots |
| 0.818 0.731 0.824 |
| StraightCurvedCapitals WordRecognition NumberRecognition |
| 0.681 0.833 0.863 |
| FigureRecognition ObjectNumber NumberFigure |
| 0.775 0.812 0.778 |
| FigureWord Deduction NumericalPuzzles |
| 0.816 0.612 0.676 |
| ProblemReasoning SeriesCompletion ArithmeticProblems |
| 0.619 0.524 0.593 |
| |
| Loadings: |
| Factor1 |
| VisualPerception 0.569 |
| Cubes 0.366 |
| PaperFormBoard 0.412 |
| Flags 0.482 |
| GeneralInformation 0.716 |
| PargraphComprehension 0.713 |
| SentenceCompletion 0.707 |
| WordClassification 0.697 |
| WordMeaning 0.725 |
| Addition 0.426 |
| Code 0.519 |
| CountingDots 0.419 |
| StraightCurvedCapitals 0.565 |
| WordRecognition 0.408 |
| NumberRecognition 0.370 |
| FigureRecognition 0.474 |
| ObjectNumber 0.434 |
| NumberFigure 0.471 |
| FigureWord 0.429 |
| Deduction 0.623 |
| NumericalPuzzles 0.569 |
| ProblemReasoning 0.617 |
| SeriesCompletion 0.690 |
| ArithmeticProblems 0.638 |
| |
| Factor1 |
| SS loadings 7.438 |
| Proportion Var 0.310 |
| |
| Test of the hypothesis that 1 factor is sufficient. |
| The chi square statistic is 622.91 on 252 degrees of freedom. |
| The p-value is 2.28e-33 |
| > for(factors in 2:5) print(update(Harman74.FA, factors = factors)) |
| |
| Call: |
| factanal(factors = factors, covmat = Harman74.cor) |
| |
| Uniquenesses: |
| VisualPerception Cubes PaperFormBoard |
| 0.650 0.864 0.844 |
| Flags GeneralInformation PargraphComprehension |
| 0.778 0.375 0.316 |
| SentenceCompletion WordClassification WordMeaning |
| 0.319 0.503 0.258 |
| Addition Code CountingDots |
| 0.670 0.608 0.581 |
| StraightCurvedCapitals WordRecognition NumberRecognition |
| 0.567 0.832 0.850 |
| FigureRecognition ObjectNumber NumberFigure |
| 0.743 0.770 0.625 |
| FigureWord Deduction NumericalPuzzles |
| 0.792 0.629 0.579 |
| ProblemReasoning SeriesCompletion ArithmeticProblems |
| 0.634 0.539 0.553 |
| |
| Loadings: |
| Factor1 Factor2 |
| VisualPerception 0.506 0.306 |
| Cubes 0.304 0.209 |
| PaperFormBoard 0.297 0.260 |
| Flags 0.327 0.339 |
| GeneralInformation 0.240 0.753 |
| PargraphComprehension 0.171 0.809 |
| SentenceCompletion 0.163 0.809 |
| WordClassification 0.344 0.615 |
| WordMeaning 0.148 0.849 |
| Addition 0.563 0.115 |
| Code 0.591 0.207 |
| CountingDots 0.647 |
| StraightCurvedCapitals 0.612 0.241 |
| WordRecognition 0.315 0.263 |
| NumberRecognition 0.328 0.205 |
| FigureRecognition 0.457 0.218 |
| ObjectNumber 0.431 0.209 |
| NumberFigure 0.601 0.116 |
| FigureWord 0.399 0.222 |
| Deduction 0.379 0.477 |
| NumericalPuzzles 0.604 0.237 |
| ProblemReasoning 0.390 0.462 |
| SeriesCompletion 0.486 0.474 |
| ArithmeticProblems 0.544 0.389 |
| |
| Factor1 Factor2 |
| SS loadings 4.573 4.548 |
| Proportion Var 0.191 0.190 |
| Cumulative Var 0.191 0.380 |
| |
| Test of the hypothesis that 2 factors are sufficient. |
| The chi square statistic is 420.24 on 229 degrees of freedom. |
| The p-value is 2.01e-13 |
| |
| Call: |
| factanal(factors = factors, covmat = Harman74.cor) |
| |
| Uniquenesses: |
| VisualPerception Cubes PaperFormBoard |
| 0.500 0.793 0.662 |
| Flags GeneralInformation PargraphComprehension |
| 0.694 0.352 0.316 |
| SentenceCompletion WordClassification WordMeaning |
| 0.300 0.502 0.256 |
| Addition Code CountingDots |
| 0.200 0.586 0.494 |
| StraightCurvedCapitals WordRecognition NumberRecognition |
| 0.569 0.838 0.848 |
| FigureRecognition ObjectNumber NumberFigure |
| 0.643 0.780 0.635 |
| FigureWord Deduction NumericalPuzzles |
| 0.788 0.590 0.580 |
| ProblemReasoning SeriesCompletion ArithmeticProblems |
| 0.597 0.498 0.500 |
| |
| Loadings: |
| Factor1 Factor2 Factor3 |
| VisualPerception 0.176 0.656 0.198 |
| Cubes 0.122 0.428 |
| PaperFormBoard 0.145 0.563 |
| Flags 0.239 0.487 0.107 |
| GeneralInformation 0.745 0.191 0.237 |
| PargraphComprehension 0.780 0.249 0.118 |
| SentenceCompletion 0.802 0.175 0.160 |
| WordClassification 0.571 0.327 0.256 |
| WordMeaning 0.821 0.248 |
| Addition 0.162 -0.118 0.871 |
| Code 0.198 0.219 0.572 |
| CountingDots 0.179 0.688 |
| StraightCurvedCapitals 0.190 0.381 0.499 |
| WordRecognition 0.231 0.253 0.210 |
| NumberRecognition 0.158 0.299 0.195 |
| FigureRecognition 0.108 0.557 0.186 |
| ObjectNumber 0.178 0.267 0.342 |
| NumberFigure 0.427 0.424 |
| FigureWord 0.167 0.355 0.240 |
| Deduction 0.392 0.472 0.181 |
| NumericalPuzzles 0.178 0.406 0.473 |
| ProblemReasoning 0.382 0.473 0.182 |
| SeriesCompletion 0.379 0.528 0.283 |
| ArithmeticProblems 0.377 0.226 0.554 |
| |
| Factor1 Factor2 Factor3 |
| SS loadings 3.802 3.488 3.186 |
| Proportion Var 0.158 0.145 0.133 |
| Cumulative Var 0.158 0.304 0.436 |
| |
| Test of the hypothesis that 3 factors are sufficient. |
| The chi square statistic is 295.59 on 207 degrees of freedom. |
| The p-value is 5.12e-05 |
| |
| Call: |
| factanal(factors = factors, covmat = Harman74.cor) |
| |
| Uniquenesses: |
| VisualPerception Cubes PaperFormBoard |
| 0.438 0.780 0.644 |
| Flags GeneralInformation PargraphComprehension |
| 0.651 0.352 0.312 |
| SentenceCompletion WordClassification WordMeaning |
| 0.283 0.485 0.257 |
| Addition Code CountingDots |
| 0.240 0.551 0.435 |
| StraightCurvedCapitals WordRecognition NumberRecognition |
| 0.491 0.646 0.696 |
| FigureRecognition ObjectNumber NumberFigure |
| 0.549 0.598 0.593 |
| FigureWord Deduction NumericalPuzzles |
| 0.762 0.592 0.583 |
| ProblemReasoning SeriesCompletion ArithmeticProblems |
| 0.601 0.497 0.500 |
| |
| Loadings: |
| Factor1 Factor2 Factor3 Factor4 |
| VisualPerception 0.160 0.689 0.187 0.160 |
| Cubes 0.117 0.436 |
| PaperFormBoard 0.137 0.570 0.110 |
| Flags 0.233 0.527 |
| GeneralInformation 0.739 0.185 0.213 0.150 |
| PargraphComprehension 0.767 0.205 0.233 |
| SentenceCompletion 0.806 0.197 0.153 |
| WordClassification 0.569 0.339 0.242 0.132 |
| WordMeaning 0.806 0.201 0.227 |
| Addition 0.167 -0.118 0.831 0.166 |
| Code 0.180 0.120 0.512 0.374 |
| CountingDots 0.210 0.716 |
| StraightCurvedCapitals 0.188 0.438 0.525 |
| WordRecognition 0.197 0.553 |
| NumberRecognition 0.122 0.116 0.520 |
| FigureRecognition 0.408 0.525 |
| ObjectNumber 0.142 0.219 0.574 |
| NumberFigure 0.293 0.336 0.456 |
| FigureWord 0.148 0.239 0.161 0.365 |
| Deduction 0.378 0.402 0.118 0.301 |
| NumericalPuzzles 0.175 0.381 0.438 0.223 |
| ProblemReasoning 0.366 0.399 0.123 0.301 |
| SeriesCompletion 0.369 0.500 0.244 0.239 |
| ArithmeticProblems 0.370 0.158 0.496 0.304 |
| |
| Factor1 Factor2 Factor3 Factor4 |
| SS loadings 3.647 2.872 2.657 2.290 |
| Proportion Var 0.152 0.120 0.111 0.095 |
| Cumulative Var 0.152 0.272 0.382 0.478 |
| |
| Test of the hypothesis that 4 factors are sufficient. |
| The chi square statistic is 226.68 on 186 degrees of freedom. |
| The p-value is 0.0224 |
| |
| Call: |
| factanal(factors = factors, covmat = Harman74.cor) |
| |
| Uniquenesses: |
| VisualPerception Cubes PaperFormBoard |
| 0.450 0.781 0.639 |
| Flags GeneralInformation PargraphComprehension |
| 0.649 0.357 0.288 |
| SentenceCompletion WordClassification WordMeaning |
| 0.277 0.485 0.262 |
| Addition Code CountingDots |
| 0.215 0.386 0.444 |
| StraightCurvedCapitals WordRecognition NumberRecognition |
| 0.256 0.639 0.706 |
| FigureRecognition ObjectNumber NumberFigure |
| 0.550 0.614 0.596 |
| FigureWord Deduction NumericalPuzzles |
| 0.764 0.521 0.564 |
| ProblemReasoning SeriesCompletion ArithmeticProblems |
| 0.580 0.442 0.478 |
| |
| Loadings: |
| Factor1 Factor2 Factor3 Factor4 Factor5 |
| VisualPerception 0.161 0.658 0.136 0.182 0.199 |
| Cubes 0.113 0.435 0.107 |
| PaperFormBoard 0.135 0.562 0.107 0.116 |
| Flags 0.231 0.533 |
| GeneralInformation 0.736 0.188 0.192 0.162 |
| PargraphComprehension 0.775 0.187 0.251 0.113 |
| SentenceCompletion 0.809 0.208 0.136 |
| WordClassification 0.568 0.348 0.223 0.131 |
| WordMeaning 0.800 0.215 0.224 |
| Addition 0.175 -0.100 0.844 0.176 |
| Code 0.185 0.438 0.451 0.426 |
| CountingDots 0.222 0.690 0.101 0.140 |
| StraightCurvedCapitals 0.186 0.425 0.458 0.559 |
| WordRecognition 0.197 0.557 |
| NumberRecognition 0.121 0.130 0.508 |
| FigureRecognition 0.400 0.529 |
| ObjectNumber 0.145 0.208 0.562 |
| NumberFigure 0.306 0.325 0.452 |
| FigureWord 0.147 0.242 0.145 0.364 |
| Deduction 0.370 0.452 0.139 0.287 -0.190 |
| NumericalPuzzles 0.170 0.402 0.439 0.230 |
| ProblemReasoning 0.358 0.423 0.126 0.302 |
| SeriesCompletion 0.360 0.549 0.256 0.223 -0.107 |
| ArithmeticProblems 0.371 0.185 0.502 0.307 |
| |
| Factor1 Factor2 Factor3 Factor4 Factor5 |
| SS loadings 3.632 2.964 2.456 2.345 0.663 |
| Proportion Var 0.151 0.124 0.102 0.098 0.028 |
| Cumulative Var 0.151 0.275 0.377 0.475 0.503 |
| |
| Test of the hypothesis that 5 factors are sufficient. |
| The chi square statistic is 186.82 on 166 degrees of freedom. |
| The p-value is 0.128 |
| > Harman74.FA <- factanal(factors = 5, covmat = Harman74.cor, |
| + rotation = "promax") |
| > print(Harman74.FA$loadings, sort = TRUE) |
| |
| Loadings: |
| Factor1 Factor2 Factor3 Factor4 Factor5 |
| VisualPerception 0.831 -0.127 0.230 |
| Cubes 0.534 |
| PaperFormBoard 0.736 -0.290 0.136 |
| Flags 0.647 -0.104 |
| SeriesCompletion 0.555 0.126 0.127 |
| GeneralInformation 0.764 |
| PargraphComprehension 0.845 -0.140 0.140 |
| SentenceCompletion 0.872 -0.140 |
| WordClassification 0.277 0.505 0.104 |
| WordMeaning 0.846 -0.108 |
| Addition -0.334 1.012 |
| CountingDots 0.206 -0.200 0.722 0.185 |
| ArithmeticProblems 0.197 0.500 0.139 |
| WordRecognition -0.126 0.127 -0.103 0.657 |
| NumberRecognition 0.568 |
| FigureRecognition 0.399 -0.142 -0.207 0.562 |
| ObjectNumber -0.108 0.107 0.613 |
| StraightCurvedCapitals 0.542 0.247 0.618 |
| Code 0.112 0.288 0.486 0.424 |
| NumberFigure 0.255 -0.230 0.211 0.413 |
| FigureWord 0.187 0.347 |
| Deduction 0.404 0.169 0.117 -0.203 |
| NumericalPuzzles 0.393 0.368 |
| ProblemReasoning 0.381 0.188 0.169 |
| |
| Factor1 Factor2 Factor3 Factor4 Factor5 |
| SS loadings 3.529 3.311 2.367 2.109 0.762 |
| Proportion Var 0.147 0.138 0.099 0.088 0.032 |
| Cumulative Var 0.147 0.285 0.384 0.471 0.503 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("InsectSprays") |
| > ### * InsectSprays |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: InsectSprays |
| > ### Title: Effectiveness of Insect Sprays |
| > ### Aliases: InsectSprays |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > boxplot(count ~ spray, data = InsectSprays, |
| + xlab = "Type of spray", ylab = "Insect count", |
| + main = "InsectSprays data", varwidth = TRUE, col = "lightgray") |
| > fm1 <- aov(count ~ spray, data = InsectSprays) |
| > summary(fm1) |
| Df Sum Sq Mean Sq F value Pr(>F) |
| spray 5 2669 533.8 34.7 <2e-16 *** |
| Residuals 66 1015 15.4 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0)) |
| > plot(fm1) |
| > fm2 <- aov(sqrt(count) ~ spray, data = InsectSprays) |
| > summary(fm2) |
| Df Sum Sq Mean Sq F value Pr(>F) |
| spray 5 88.44 17.688 44.8 <2e-16 *** |
| Residuals 66 26.06 0.395 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > plot(fm2) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("JohnsonJohnson") |
| > ### * JohnsonJohnson |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: JohnsonJohnson |
| > ### Title: Quarterly Earnings per Johnson & Johnson Share |
| > ### Aliases: JohnsonJohnson |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > |
| > cleanEx() |
| > nameEx("LifeCycleSavings") |
| > ### * LifeCycleSavings |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: LifeCycleSavings |
| > ### Title: Intercountry Life-Cycle Savings Data |
| > ### Aliases: LifeCycleSavings |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > pairs(LifeCycleSavings, panel = panel.smooth, |
| + main = "LifeCycleSavings data") |
| > fm1 <- lm(sr ~ pop15 + pop75 + dpi + ddpi, data = LifeCycleSavings) |
| > summary(fm1) |
| |
| Call: |
| lm(formula = sr ~ pop15 + pop75 + dpi + ddpi, data = LifeCycleSavings) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -8.2422 -2.6857 -0.2488 2.4280 9.7509 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 28.5660865 7.3545161 3.884 0.000334 *** |
| pop15 -0.4611931 0.1446422 -3.189 0.002603 ** |
| pop75 -1.6914977 1.0835989 -1.561 0.125530 |
| dpi -0.0003369 0.0009311 -0.362 0.719173 |
| ddpi 0.4096949 0.1961971 2.088 0.042471 * |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 3.803 on 45 degrees of freedom |
| Multiple R-squared: 0.3385, Adjusted R-squared: 0.2797 |
| F-statistic: 5.756 on 4 and 45 DF, p-value: 0.0007904 |
| |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Loblolly") |
| > ### * Loblolly |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Loblolly |
| > ### Title: Growth of Loblolly pine trees |
| > ### Aliases: Loblolly |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > plot(height ~ age, data = Loblolly, subset = Seed == 329, |
| + xlab = "Tree age (yr)", las = 1, |
| + ylab = "Tree height (ft)", |
| + main = "Loblolly data and fitted curve (Seed 329 only)") |
| > fm1 <- nls(height ~ SSasymp(age, Asym, R0, lrc), |
| + data = Loblolly, subset = Seed == 329) |
| > age <- seq(0, 30, length.out = 101) |
| > lines(age, predict(fm1, list(age = age))) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Nile") |
| > ### * Nile |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Nile |
| > ### Title: Flow of the River Nile |
| > ### Aliases: Nile |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > par(mfrow = c(2, 2)) |
| > plot(Nile) |
| > acf(Nile) |
| > pacf(Nile) |
| > ar(Nile) # selects order 2 |
| |
| Call: |
| ar(x = Nile) |
| |
| Coefficients: |
| 1 2 |
| 0.4081 0.1812 |
| |
| Order selected 2 sigma^2 estimated as 21247 |
| > cpgram(ar(Nile)$resid) |
| > par(mfrow = c(1, 1)) |
| > arima(Nile, c(2, 0, 0)) |
| |
| Call: |
| arima(x = Nile, order = c(2, 0, 0)) |
| |
| Coefficients: |
| ar1 ar2 intercept |
| 0.4096 0.1987 919.8397 |
| s.e. 0.0974 0.0990 35.6410 |
| |
| sigma^2 estimated as 20291: log likelihood = -637.98, aic = 1283.96 |
| > |
| > ## Now consider missing values, following Durbin & Koopman |
| > NileNA <- Nile |
| > NileNA[c(21:40, 61:80)] <- NA |
| > arima(NileNA, c(2, 0, 0)) |
| |
| Call: |
| arima(x = NileNA, order = c(2, 0, 0)) |
| |
| Coefficients: |
| ar1 ar2 intercept |
| 0.3622 0.1678 918.3103 |
| s.e. 0.1273 0.1323 39.5037 |
| |
| sigma^2 estimated as 23676: log likelihood = -387.7, aic = 783.41 |
| > plot(NileNA) |
| > pred <- |
| + predict(arima(window(NileNA, 1871, 1890), c(2, 0, 0)), n.ahead = 20) |
| > lines(pred$pred, lty = 3, col = "red") |
| > lines(pred$pred + 2*pred$se, lty = 2, col = "blue") |
| > lines(pred$pred - 2*pred$se, lty = 2, col = "blue") |
| > pred <- |
| + predict(arima(window(NileNA, 1871, 1930), c(2, 0, 0)), n.ahead = 20) |
| > lines(pred$pred, lty = 3, col = "red") |
| > lines(pred$pred + 2*pred$se, lty = 2, col = "blue") |
| > lines(pred$pred - 2*pred$se, lty = 2, col = "blue") |
| > |
| > ## Structural time series models |
| > par(mfrow = c(3, 1)) |
| > plot(Nile) |
| > ## local level model |
| > (fit <- StructTS(Nile, type = "level")) |
| |
| Call: |
| StructTS(x = Nile, type = "level") |
| |
| Variances: |
| level epsilon |
| 1469 15099 |
| > lines(fitted(fit), lty = 2) # contemporaneous smoothing |
| > lines(tsSmooth(fit), lty = 2, col = 4) # fixed-interval smoothing |
| > plot(residuals(fit)); abline(h = 0, lty = 3) |
| > ## local trend model |
| > (fit2 <- StructTS(Nile, type = "trend")) ## constant trend fitted |
| |
| Call: |
| StructTS(x = Nile, type = "trend") |
| |
| Variances: |
| level slope epsilon |
| 1427 0 15047 |
| > pred <- predict(fit, n.ahead = 30) |
| > ## with 50% confidence interval |
| > ts.plot(Nile, pred$pred, |
| + pred$pred + 0.67*pred$se, pred$pred -0.67*pred$se) |
| > |
| > ## Now consider missing values |
| > plot(NileNA) |
| > (fit3 <- StructTS(NileNA, type = "level")) |
| |
| Call: |
| StructTS(x = NileNA, type = "level") |
| |
| Variances: |
| level epsilon |
| 685.8 17899.8 |
| > lines(fitted(fit3), lty = 2) |
| > lines(tsSmooth(fit3), lty = 3) |
| > plot(residuals(fit3)); abline(h = 0, lty = 3) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("Orange") |
| > ### * Orange |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Orange |
| > ### Title: Growth of Orange Trees |
| > ### Aliases: Orange |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > coplot(circumference ~ age | Tree, data = Orange, show.given = FALSE) |
| > fm1 <- nls(circumference ~ SSlogis(age, Asym, xmid, scal), |
| + data = Orange, subset = Tree == 3) |
| > plot(circumference ~ age, data = Orange, subset = Tree == 3, |
| + xlab = "Tree age (days since 1968/12/31)", |
| + ylab = "Tree circumference (mm)", las = 1, |
| + main = "Orange tree data and fitted model (Tree 3 only)") |
| > age <- seq(0, 1600, length.out = 101) |
| > lines(age, predict(fm1, list(age = age))) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("OrchardSprays") |
| > ### * OrchardSprays |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: OrchardSprays |
| > ### Title: Potency of Orchard Sprays |
| > ### Aliases: OrchardSprays |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > pairs(OrchardSprays, main = "OrchardSprays data") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("PlantGrowth") |
| > ### * PlantGrowth |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: PlantGrowth |
| > ### Title: Results from an Experiment on Plant Growth |
| > ### Aliases: PlantGrowth |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > ## One factor ANOVA example from Dobson's book, cf. Table 7.4: |
| > require(stats); require(graphics) |
| > boxplot(weight ~ group, data = PlantGrowth, main = "PlantGrowth data", |
| + ylab = "Dried weight of plants", col = "lightgray", |
| + notch = TRUE, varwidth = TRUE) |
| Warning in bxp(list(stats = c(4.17, 4.53, 5.155, 5.33, 6.11, 3.59, 4.17, : |
| some notches went outside hinges ('box'): maybe set notch=FALSE |
| > anova(lm(weight ~ group, data = PlantGrowth)) |
| Analysis of Variance Table |
| |
| Response: weight |
| Df Sum Sq Mean Sq F value Pr(>F) |
| group 2 3.7663 1.8832 4.8461 0.01591 * |
| Residuals 27 10.4921 0.3886 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Puromycin") |
| > ### * Puromycin |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Puromycin |
| > ### Title: Reaction Velocity of an Enzymatic Reaction |
| > ### Aliases: Puromycin |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## Don't show: |
| > options(show.nls.convergence=FALSE) |
| > ## End(Don't show) |
| > plot(rate ~ conc, data = Puromycin, las = 1, |
| + xlab = "Substrate concentration (ppm)", |
| + ylab = "Reaction velocity (counts/min/min)", |
| + pch = as.integer(Puromycin$state), |
| + col = as.integer(Puromycin$state), |
| + main = "Puromycin data and fitted Michaelis-Menten curves") |
| > ## simplest form of fitting the Michaelis-Menten model to these data |
| > fm1 <- nls(rate ~ Vm * conc/(K + conc), data = Puromycin, |
| + subset = state == "treated", |
| + start = c(Vm = 200, K = 0.05)) |
| > fm2 <- nls(rate ~ Vm * conc/(K + conc), data = Puromycin, |
| + subset = state == "untreated", |
| + start = c(Vm = 160, K = 0.05)) |
| > summary(fm1) |
| |
| Formula: rate ~ Vm * conc/(K + conc) |
| |
| Parameters: |
| Estimate Std. Error t value Pr(>|t|) |
| Vm 2.127e+02 6.947e+00 30.615 3.24e-11 *** |
| K 6.412e-02 8.281e-03 7.743 1.57e-05 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 10.93 on 10 degrees of freedom |
| |
| > summary(fm2) |
| |
| Formula: rate ~ Vm * conc/(K + conc) |
| |
| Parameters: |
| Estimate Std. Error t value Pr(>|t|) |
| Vm 1.603e+02 6.480e+00 24.734 1.38e-09 *** |
| K 4.771e-02 7.782e-03 6.131 0.000173 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 9.773 on 9 degrees of freedom |
| |
| > ## add fitted lines to the plot |
| > conc <- seq(0, 1.2, length.out = 101) |
| > lines(conc, predict(fm1, list(conc = conc)), lty = 1, col = 1) |
| > lines(conc, predict(fm2, list(conc = conc)), lty = 2, col = 2) |
| > legend(0.8, 120, levels(Puromycin$state), |
| + col = 1:2, lty = 1:2, pch = 1:2) |
| > |
| > ## using partial linearity |
| > fm3 <- nls(rate ~ conc/(K + conc), data = Puromycin, |
| + subset = state == "treated", start = c(K = 0.05), |
| + algorithm = "plinear") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("Theoph") |
| > ### * Theoph |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Theoph |
| > ### Title: Pharmacokinetics of Theophylline |
| > ### Aliases: Theoph |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## Don't show: |
| > options(show.nls.convergence=FALSE) |
| > ## End(Don't show) |
| > coplot(conc ~ Time | Subject, data = Theoph, show.given = FALSE) |
| > Theoph.4 <- subset(Theoph, Subject == 4) |
| > fm1 <- nls(conc ~ SSfol(Dose, Time, lKe, lKa, lCl), |
| + data = Theoph.4) |
| > summary(fm1) |
| |
| Formula: conc ~ SSfol(Dose, Time, lKe, lKa, lCl) |
| |
| Parameters: |
| Estimate Std. Error t value Pr(>|t|) |
| lKe -2.4365 0.2257 -10.797 4.77e-06 *** |
| lKa 0.1583 0.2297 0.689 0.51 |
| lCl -3.2861 0.1448 -22.695 1.51e-08 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.8465 on 8 degrees of freedom |
| |
| > plot(conc ~ Time, data = Theoph.4, |
| + xlab = "Time since drug administration (hr)", |
| + ylab = "Theophylline concentration (mg/L)", |
| + main = "Observed concentrations and fitted model", |
| + sub = "Theophylline data - Subject 4 only", |
| + las = 1, col = 4) |
| > xvals <- seq(0, par("usr")[2], length.out = 55) |
| > lines(xvals, predict(fm1, newdata = list(Time = xvals)), |
| + col = 4) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("Titanic") |
| > ### * Titanic |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: Titanic |
| > ### Title: Survival of passengers on the Titanic |
| > ### Aliases: Titanic |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > mosaicplot(Titanic, main = "Survival on the Titanic") |
| > ## Higher survival rates in children? |
| > apply(Titanic, c(3, 4), sum) |
| Survived |
| Age No Yes |
| Child 52 57 |
| Adult 1438 654 |
| > ## Higher survival rates in females? |
| > apply(Titanic, c(2, 4), sum) |
| Survived |
| Sex No Yes |
| Male 1364 367 |
| Female 126 344 |
| > ## Use loglm() in package 'MASS' for further analysis ... |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("ToothGrowth") |
| > ### * ToothGrowth |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: ToothGrowth |
| > ### Title: The Effect of Vitamin C on Tooth Growth in Guinea Pigs |
| > ### Aliases: ToothGrowth |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > coplot(len ~ dose | supp, data = ToothGrowth, panel = panel.smooth, |
| + xlab = "ToothGrowth data: length vs dose, given type of supplement") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("UCBAdmissions") |
| > ### * UCBAdmissions |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: UCBAdmissions |
| > ### Title: Student Admissions at UC Berkeley |
| > ### Aliases: UCBAdmissions |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > ## Data aggregated over departments |
| > apply(UCBAdmissions, c(1, 2), sum) |
| Gender |
| Admit Male Female |
| Admitted 1198 557 |
| Rejected 1493 1278 |
| > mosaicplot(apply(UCBAdmissions, c(1, 2), sum), |
| + main = "Student admissions at UC Berkeley") |
| > ## Data for individual departments |
| > opar <- par(mfrow = c(2, 3), oma = c(0, 0, 2, 0)) |
| > for(i in 1:6) |
| + mosaicplot(UCBAdmissions[,,i], |
| + xlab = "Admit", ylab = "Sex", |
| + main = paste("Department", LETTERS[i])) |
| > mtext(expression(bold("Student admissions at UC Berkeley")), |
| + outer = TRUE, cex = 1.5) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("UKDriverDeaths") |
| > ### * UKDriverDeaths |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: UKDriverDeaths |
| > ### Title: Road Casualties in Great Britain 1969-84 |
| > ### Aliases: UKDriverDeaths Seatbelts |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## work with pre-seatbelt period to identify a model, use logs |
| > work <- window(log10(UKDriverDeaths), end = 1982+11/12) |
| > par(mfrow = c(3, 1)) |
| > plot(work); acf(work); pacf(work) |
| > par(mfrow = c(1, 1)) |
| > (fit <- arima(work, c(1, 0, 0), seasonal = list(order = c(1, 0, 0)))) |
| |
| Call: |
| arima(x = work, order = c(1, 0, 0), seasonal = list(order = c(1, 0, 0))) |
| |
| Coefficients: |
| ar1 sar1 intercept |
| 0.4378 0.6281 3.2274 |
| s.e. 0.0764 0.0637 0.0131 |
| |
| sigma^2 estimated as 0.00157: log likelihood = 300.85, aic = -593.7 |
| > z <- predict(fit, n.ahead = 24) |
| > ts.plot(log10(UKDriverDeaths), z$pred, z$pred+2*z$se, z$pred-2*z$se, |
| + lty = c(1, 3, 2, 2), col = c("black", "red", "blue", "blue")) |
| > |
| > ## now see the effect of the explanatory variables |
| > X <- Seatbelts[, c("kms", "PetrolPrice", "law")] |
| > X[, 1] <- log10(X[, 1]) - 4 |
| > arima(log10(Seatbelts[, "drivers"]), c(1, 0, 0), |
| + seasonal = list(order = c(1, 0, 0)), xreg = X) |
| |
| Call: |
| arima(x = log10(Seatbelts[, "drivers"]), order = c(1, 0, 0), seasonal = list(order = c(1, |
| 0, 0)), xreg = X) |
| |
| Coefficients: |
| ar1 sar1 intercept kms PetrolPrice law |
| 0.3348 0.6672 3.3539 0.0082 -1.2224 -0.0963 |
| s.e. 0.0775 0.0612 0.0441 0.0902 0.3839 0.0166 |
| |
| sigma^2 estimated as 0.001476: log likelihood = 349.73, aic = -685.46 |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("UKLungDeaths") |
| > ### * UKLungDeaths |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: UKLungDeaths |
| > ### Title: Monthly Deaths from Lung Diseases in the UK |
| > ### Aliases: UKLungDeaths ldeaths fdeaths mdeaths |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) # for time |
| > plot(ldeaths) |
| > plot(mdeaths, fdeaths) |
| > ## Better labels: |
| > yr <- floor(tt <- time(mdeaths)) |
| > plot(mdeaths, fdeaths, |
| + xy.labels = paste(month.abb[12*(tt - yr)], yr-1900, sep = "'")) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("UKgas") |
| > ### * UKgas |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: UKgas |
| > ### Title: UK Quarterly Gas Consumption |
| > ### Aliases: UKgas |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > ## maybe str(UKgas) ; plot(UKgas) ... |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("USArrests") |
| > ### * USArrests |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: USArrests |
| > ### Title: Violent Crime Rates by US State |
| > ### Aliases: USArrests |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > summary(USArrests) |
| Murder Assault UrbanPop Rape |
| Min. : 0.800 Min. : 45.0 Min. :32.00 Min. : 7.30 |
| 1st Qu.: 4.075 1st Qu.:109.0 1st Qu.:54.50 1st Qu.:15.07 |
| Median : 7.250 Median :159.0 Median :66.00 Median :20.10 |
| Mean : 7.788 Mean :170.8 Mean :65.54 Mean :21.23 |
| 3rd Qu.:11.250 3rd Qu.:249.0 3rd Qu.:77.75 3rd Qu.:26.18 |
| Max. :17.400 Max. :337.0 Max. :91.00 Max. :46.00 |
| > |
| > require(graphics) |
| > pairs(USArrests, panel = panel.smooth, main = "USArrests data") |
| > |
| > ## Difference between 'USArrests' and its correction |
| > USArrests["Maryland", "UrbanPop"] # 67 -- the transcription error |
| [1] 67 |
| > UA.C <- USArrests |
| > UA.C["Maryland", "UrbanPop"] <- 76.6 |
| > |
| > ## also +/- 0.5 to restore the original <n>.5 percentages |
| > s5u <- c("Colorado", "Florida", "Mississippi", "Wyoming") |
| > s5d <- c("Nebraska", "Pennsylvania") |
| > UA.C[s5u, "UrbanPop"] <- UA.C[s5u, "UrbanPop"] + 0.5 |
| > UA.C[s5d, "UrbanPop"] <- UA.C[s5d, "UrbanPop"] - 0.5 |
| > |
| > ## ==> UA.C is now a *C*orrected version of USArrests |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("USJudgeRatings") |
| > ### * USJudgeRatings |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: USJudgeRatings |
| > ### Title: Lawyers' Ratings of State Judges in the US Superior Court |
| > ### Aliases: USJudgeRatings |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > pairs(USJudgeRatings, main = "USJudgeRatings data") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("USPersonalExpenditure") |
| > ### * USPersonalExpenditure |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: USPersonalExpenditure |
| > ### Title: Personal Expenditure Data |
| > ### Aliases: USPersonalExpenditure |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) # for medpolish |
| > USPersonalExpenditure |
| 1940 1945 1950 1955 1960 |
| Food and Tobacco 22.200 44.500 59.60 73.2 86.80 |
| Household Operation 10.500 15.500 29.00 36.5 46.20 |
| Medical and Health 3.530 5.760 9.71 14.0 21.10 |
| Personal Care 1.040 1.980 2.45 3.4 5.40 |
| Private Education 0.341 0.974 1.80 2.6 3.64 |
| > medpolish(log10(USPersonalExpenditure)) |
| 1: 1.126317 |
| 2: 1.032421 |
| Final: 1.032421 |
| |
| Median Polish Results (Dataset: "log10(USPersonalExpenditure)") |
| |
| Overall: 0.9872192 |
| |
| Row Effects: |
| Food and Tobacco Household Operation Medical and Health Personal Care |
| 0.7880270 0.4327608 0.0000000 -0.5606543 |
| Private Education |
| -0.7319467 |
| |
| Column Effects: |
| 1940 1945 1950 1955 1960 |
| -0.4288933 -0.2267967 0.0000000 0.1423128 0.3058289 |
| |
| Residuals: |
| 1940 1945 1950 1955 1960 |
| Food and Tobacco 0.000000 0.0999105 0.000000 -0.053048 -0.142555 |
| Household Operation 0.030103 -0.0028516 0.042418 0.000000 -0.061167 |
| Medical and Health -0.010551 0.0000000 0.000000 0.016596 0.031234 |
| Personal Care 0.019362 0.0968971 -0.037399 -0.037399 0.000000 |
| Private Education -0.293625 -0.0399168 0.000000 0.017388 0.000000 |
| |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("VADeaths") |
| > ### * VADeaths |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: VADeaths |
| > ### Title: Death Rates in Virginia (1940) |
| > ### Aliases: VADeaths |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > n <- length(dr <- c(VADeaths)) |
| > nam <- names(VADeaths) |
| > d.VAD <- data.frame( |
| + Drate = dr, |
| + age = rep(ordered(rownames(VADeaths)), length.out = n), |
| + gender = gl(2, 5, n, labels = c("M", "F")), |
| + site = gl(2, 10, labels = c("rural", "urban"))) |
| > coplot(Drate ~ as.numeric(age) | gender * site, data = d.VAD, |
| + panel = panel.smooth, xlab = "VADeaths data - Given: gender") |
| > summary(aov.VAD <- aov(Drate ~ .^2, data = d.VAD)) |
| Df Sum Sq Mean Sq F value Pr(>F) |
| age 4 6288 1572.1 590.858 8.55e-06 *** |
| gender 1 648 647.5 243.361 9.86e-05 *** |
| site 1 77 76.8 28.876 0.00579 ** |
| age:gender 4 86 21.6 8.100 0.03358 * |
| age:site 4 43 10.6 3.996 0.10414 |
| gender:site 1 73 73.0 27.422 0.00636 ** |
| Residuals 4 11 2.7 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0)) |
| > plot(aov.VAD) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("WWWusage") |
| > ### * WWWusage |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: WWWusage |
| > ### Title: Internet Usage per Minute |
| > ### Aliases: WWWusage |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > work <- diff(WWWusage) |
| > par(mfrow = c(2, 1)); plot(WWWusage); plot(work) |
| > ## Not run: |
| > ##D require(stats) |
| > ##D aics <- matrix(, 6, 6, dimnames = list(p = 0:5, q = 0:5)) |
| > ##D for(q in 1:5) aics[1, 1+q] <- arima(WWWusage, c(0, 1, q), |
| > ##D optim.control = list(maxit = 500))$aic |
| > ##D for(p in 1:5) |
| > ##D for(q in 0:5) aics[1+p, 1+q] <- arima(WWWusage, c(p, 1, q), |
| > ##D optim.control = list(maxit = 500))$aic |
| > ##D round(aics - min(aics, na.rm = TRUE), 2) |
| > ## End(Not run) |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("WorldPhones") |
| > ### * WorldPhones |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: WorldPhones |
| > ### Title: The World's Telephones |
| > ### Aliases: WorldPhones |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > matplot(rownames(WorldPhones), WorldPhones, type = "b", log = "y", |
| + xlab = "Year", ylab = "Number of telephones (1000's)") |
| > legend(1951.5, 80000, colnames(WorldPhones), col = 1:6, lty = 1:5, |
| + pch = rep(21, 7)) |
| > title(main = "World phones data: log scale for response") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("ability.cov") |
| > ### * ability.cov |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: ability.cov |
| > ### Title: Ability and Intelligence Tests |
| > ### Aliases: ability.cov |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > |
| > cleanEx() |
| > nameEx("airmiles") |
| > ### * airmiles |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: airmiles |
| > ### Title: Passenger Miles on Commercial US Airlines, 1937-1960 |
| > ### Aliases: airmiles |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(airmiles, main = "airmiles data", |
| + xlab = "Passenger-miles flown by U.S. commercial airlines", col = 4) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("airquality") |
| > ### * airquality |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: airquality |
| > ### Title: New York Air Quality Measurements |
| > ### Aliases: airquality |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > pairs(airquality, panel = panel.smooth, main = "airquality data") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("anscombe") |
| > ### * anscombe |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: anscombe |
| > ### Title: Anscombe's Quartet of 'Identical' Simple Linear Regressions |
| > ### Aliases: anscombe |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > summary(anscombe) |
| x1 x2 x3 x4 y1 |
| Min. : 4.0 Min. : 4.0 Min. : 4.0 Min. : 8 Min. : 4.260 |
| 1st Qu.: 6.5 1st Qu.: 6.5 1st Qu.: 6.5 1st Qu.: 8 1st Qu.: 6.315 |
| Median : 9.0 Median : 9.0 Median : 9.0 Median : 8 Median : 7.580 |
| Mean : 9.0 Mean : 9.0 Mean : 9.0 Mean : 9 Mean : 7.501 |
| 3rd Qu.:11.5 3rd Qu.:11.5 3rd Qu.:11.5 3rd Qu.: 8 3rd Qu.: 8.570 |
| Max. :14.0 Max. :14.0 Max. :14.0 Max. :19 Max. :10.840 |
| y2 y3 y4 |
| Min. :3.100 Min. : 5.39 Min. : 5.250 |
| 1st Qu.:6.695 1st Qu.: 6.25 1st Qu.: 6.170 |
| Median :8.140 Median : 7.11 Median : 7.040 |
| Mean :7.501 Mean : 7.50 Mean : 7.501 |
| 3rd Qu.:8.950 3rd Qu.: 7.98 3rd Qu.: 8.190 |
| Max. :9.260 Max. :12.74 Max. :12.500 |
| > |
| > ##-- now some "magic" to do the 4 regressions in a loop: |
| > ff <- y ~ x |
| > mods <- setNames(as.list(1:4), paste0("lm", 1:4)) |
| > for(i in 1:4) { |
| + ff[2:3] <- lapply(paste0(c("y","x"), i), as.name) |
| + ## or ff[[2]] <- as.name(paste0("y", i)) |
| + ## ff[[3]] <- as.name(paste0("x", i)) |
| + mods[[i]] <- lmi <- lm(ff, data = anscombe) |
| + print(anova(lmi)) |
| + } |
| Analysis of Variance Table |
| |
| Response: y1 |
| Df Sum Sq Mean Sq F value Pr(>F) |
| x1 1 27.510 27.5100 17.99 0.00217 ** |
| Residuals 9 13.763 1.5292 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| Analysis of Variance Table |
| |
| Response: y2 |
| Df Sum Sq Mean Sq F value Pr(>F) |
| x2 1 27.500 27.5000 17.966 0.002179 ** |
| Residuals 9 13.776 1.5307 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| Analysis of Variance Table |
| |
| Response: y3 |
| Df Sum Sq Mean Sq F value Pr(>F) |
| x3 1 27.470 27.4700 17.972 0.002176 ** |
| Residuals 9 13.756 1.5285 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| Analysis of Variance Table |
| |
| Response: y4 |
| Df Sum Sq Mean Sq F value Pr(>F) |
| x4 1 27.490 27.4900 18.003 0.002165 ** |
| Residuals 9 13.742 1.5269 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > |
| > ## See how close they are (numerically!) |
| > sapply(mods, coef) |
| lm1 lm2 lm3 lm4 |
| (Intercept) 3.0000909 3.000909 3.0024545 3.0017273 |
| x1 0.5000909 0.500000 0.4997273 0.4999091 |
| > lapply(mods, function(fm) coef(summary(fm))) |
| $lm1 |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 3.0000909 1.1247468 2.667348 0.025734051 |
| x1 0.5000909 0.1179055 4.241455 0.002169629 |
| |
| $lm2 |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 3.000909 1.1253024 2.666758 0.025758941 |
| x2 0.500000 0.1179637 4.238590 0.002178816 |
| |
| $lm3 |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 3.0024545 1.1244812 2.670080 0.025619109 |
| x3 0.4997273 0.1178777 4.239372 0.002176305 |
| |
| $lm4 |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 3.0017273 1.1239211 2.670763 0.025590425 |
| x4 0.4999091 0.1178189 4.243028 0.002164602 |
| |
| > |
| > ## Now, do what you should have done in the first place: PLOTS |
| > op <- par(mfrow = c(2, 2), mar = 0.1+c(4,4,1,1), oma = c(0, 0, 2, 0)) |
| > for(i in 1:4) { |
| + ff[2:3] <- lapply(paste0(c("y","x"), i), as.name) |
| + plot(ff, data = anscombe, col = "red", pch = 21, bg = "orange", cex = 1.2, |
| + xlim = c(3, 19), ylim = c(3, 13)) |
| + abline(mods[[i]], col = "blue") |
| + } |
| > mtext("Anscombe's 4 Regression data sets", outer = TRUE, cex = 1.5) |
| > par(op) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("attenu") |
| > ### * attenu |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: attenu |
| > ### Title: The Joyner-Boore Attenuation Data |
| > ### Aliases: attenu |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > ## check the data class of the variables |
| > sapply(attenu, data.class) |
| event mag station dist accel |
| "numeric" "numeric" "factor" "numeric" "numeric" |
| > summary(attenu) |
| event mag station dist |
| Min. : 1.00 Min. :5.000 117 : 5 Min. : 0.50 |
| 1st Qu.: 9.00 1st Qu.:5.300 1028 : 4 1st Qu.: 11.32 |
| Median :18.00 Median :6.100 113 : 4 Median : 23.40 |
| Mean :14.74 Mean :6.084 112 : 3 Mean : 45.60 |
| 3rd Qu.:20.00 3rd Qu.:6.600 135 : 3 3rd Qu.: 47.55 |
| Max. :23.00 Max. :7.700 (Other):147 Max. :370.00 |
| NA's : 16 |
| accel |
| Min. :0.00300 |
| 1st Qu.:0.04425 |
| Median :0.11300 |
| Mean :0.15422 |
| 3rd Qu.:0.21925 |
| Max. :0.81000 |
| |
| > pairs(attenu, main = "attenu data") |
| > coplot(accel ~ dist | as.factor(event), data = attenu, show.given = FALSE) |
| > coplot(log(accel) ~ log(dist) | as.factor(event), |
| + data = attenu, panel = panel.smooth, show.given = FALSE) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("attitude") |
| > ### * attitude |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: attitude |
| > ### Title: The Chatterjee-Price Attitude Data |
| > ### Aliases: attitude |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > pairs(attitude, main = "attitude data") |
| > summary(attitude) |
| rating complaints privileges learning raises |
| Min. :40.00 Min. :37.0 Min. :30.00 Min. :34.00 Min. :43.00 |
| 1st Qu.:58.75 1st Qu.:58.5 1st Qu.:45.00 1st Qu.:47.00 1st Qu.:58.25 |
| Median :65.50 Median :65.0 Median :51.50 Median :56.50 Median :63.50 |
| Mean :64.63 Mean :66.6 Mean :53.13 Mean :56.37 Mean :64.63 |
| 3rd Qu.:71.75 3rd Qu.:77.0 3rd Qu.:62.50 3rd Qu.:66.75 3rd Qu.:71.00 |
| Max. :85.00 Max. :90.0 Max. :83.00 Max. :75.00 Max. :88.00 |
| critical advance |
| Min. :49.00 Min. :25.00 |
| 1st Qu.:69.25 1st Qu.:35.00 |
| Median :77.50 Median :41.00 |
| Mean :74.77 Mean :42.93 |
| 3rd Qu.:80.00 3rd Qu.:47.75 |
| Max. :92.00 Max. :72.00 |
| > summary(fm1 <- lm(rating ~ ., data = attitude)) |
| |
| Call: |
| lm(formula = rating ~ ., data = attitude) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -10.9418 -4.3555 0.3158 5.5425 11.5990 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 10.78708 11.58926 0.931 0.361634 |
| complaints 0.61319 0.16098 3.809 0.000903 *** |
| privileges -0.07305 0.13572 -0.538 0.595594 |
| learning 0.32033 0.16852 1.901 0.069925 . |
| raises 0.08173 0.22148 0.369 0.715480 |
| critical 0.03838 0.14700 0.261 0.796334 |
| advance -0.21706 0.17821 -1.218 0.235577 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 7.068 on 23 degrees of freedom |
| Multiple R-squared: 0.7326, Adjusted R-squared: 0.6628 |
| F-statistic: 10.5 on 6 and 23 DF, p-value: 1.24e-05 |
| |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0), |
| + mar = c(4.1, 4.1, 2.1, 1.1)) |
| > plot(fm1) |
| > summary(fm2 <- lm(rating ~ complaints, data = attitude)) |
| |
| Call: |
| lm(formula = rating ~ complaints, data = attitude) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -12.8799 -5.9905 0.1783 6.2978 9.6294 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 14.37632 6.61999 2.172 0.0385 * |
| complaints 0.75461 0.09753 7.737 1.99e-08 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 6.993 on 28 degrees of freedom |
| Multiple R-squared: 0.6813, Adjusted R-squared: 0.6699 |
| F-statistic: 59.86 on 1 and 28 DF, p-value: 1.988e-08 |
| |
| > plot(fm2) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("beavers") |
| > ### * beavers |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: beavers |
| > ### Title: Body Temperature Series of Two Beavers |
| > ### Aliases: beavers beaver1 beaver2 |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > (yl <- range(beaver1$temp, beaver2$temp)) |
| [1] 36.33 38.35 |
| > |
| > beaver.plot <- function(bdat, ...) { |
| + nam <- deparse(substitute(bdat)) |
| + with(bdat, { |
| + # Hours since start of day: |
| + hours <- time %/% 100 + 24*(day - day[1]) + (time %% 100)/60 |
| + plot (hours, temp, type = "l", ..., |
| + main = paste(nam, "body temperature")) |
| + abline(h = 37.5, col = "gray", lty = 2) |
| + is.act <- activ == 1 |
| + points(hours[is.act], temp[is.act], col = 2, cex = .8) |
| + }) |
| + } |
| > op <- par(mfrow = c(2, 1), mar = c(3, 3, 4, 2), mgp = 0.9 * 2:0) |
| > beaver.plot(beaver1, ylim = yl) |
| > beaver.plot(beaver2, ylim = yl) |
| > par(op) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("cars") |
| > ### * cars |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: cars |
| > ### Title: Speed and Stopping Distances of Cars |
| > ### Aliases: cars |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > plot(cars, xlab = "Speed (mph)", ylab = "Stopping distance (ft)", |
| + las = 1) |
| > lines(lowess(cars$speed, cars$dist, f = 2/3, iter = 3), col = "red") |
| > title(main = "cars data") |
| > plot(cars, xlab = "Speed (mph)", ylab = "Stopping distance (ft)", |
| + las = 1, log = "xy") |
| > title(main = "cars data (logarithmic scales)") |
| > lines(lowess(cars$speed, cars$dist, f = 2/3, iter = 3), col = "red") |
| > summary(fm1 <- lm(log(dist) ~ log(speed), data = cars)) |
| |
| Call: |
| lm(formula = log(dist) ~ log(speed), data = cars) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -1.00215 -0.24578 -0.02898 0.20717 0.88289 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) -0.7297 0.3758 -1.941 0.0581 . |
| log(speed) 1.6024 0.1395 11.484 2.26e-15 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.4053 on 48 degrees of freedom |
| Multiple R-squared: 0.7331, Adjusted R-squared: 0.7276 |
| F-statistic: 131.9 on 1 and 48 DF, p-value: 2.259e-15 |
| |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0), |
| + mar = c(4.1, 4.1, 2.1, 1.1)) |
| > plot(fm1) |
| > par(opar) |
| > |
| > ## An example of polynomial regression |
| > plot(cars, xlab = "Speed (mph)", ylab = "Stopping distance (ft)", |
| + las = 1, xlim = c(0, 25)) |
| > d <- seq(0, 25, length.out = 200) |
| > for(degree in 1:4) { |
| + fm <- lm(dist ~ poly(speed, degree), data = cars) |
| + assign(paste("cars", degree, sep = "."), fm) |
| + lines(d, predict(fm, data.frame(speed = d)), col = degree) |
| + } |
| > anova(cars.1, cars.2, cars.3, cars.4) |
| Analysis of Variance Table |
| |
| Model 1: dist ~ poly(speed, degree) |
| Model 2: dist ~ poly(speed, degree) |
| Model 3: dist ~ poly(speed, degree) |
| Model 4: dist ~ poly(speed, degree) |
| Res.Df RSS Df Sum of Sq F Pr(>F) |
| 1 48 11354 |
| 2 47 10825 1 528.81 2.3108 0.1355 |
| 3 46 10634 1 190.35 0.8318 0.3666 |
| 4 45 10298 1 336.55 1.4707 0.2316 |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("chickwts") |
| > ### * chickwts |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: chickwts |
| > ### Title: Chicken Weights by Feed Type |
| > ### Aliases: chickwts |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > boxplot(weight ~ feed, data = chickwts, col = "lightgray", |
| + varwidth = TRUE, notch = TRUE, main = "chickwt data", |
| + ylab = "Weight at six weeks (gm)") |
| Warning in bxp(list(stats = c(216, 271.5, 342, 373.5, 404, 108, 136, 151.5, : |
| some notches went outside hinges ('box'): maybe set notch=FALSE |
| > anova(fm1 <- lm(weight ~ feed, data = chickwts)) |
| Analysis of Variance Table |
| |
| Response: weight |
| Df Sum Sq Mean Sq F value Pr(>F) |
| feed 5 231129 46226 15.365 5.936e-10 *** |
| Residuals 65 195556 3009 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0), |
| + mar = c(4.1, 4.1, 2.1, 1.1)) |
| > plot(fm1) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("co2") |
| > ### * co2 |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: co2 |
| > ### Title: Mauna Loa Atmospheric CO2 Concentration |
| > ### Aliases: co2 |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(co2, ylab = expression("Atmospheric concentration of CO"[2]), |
| + las = 1) |
| > title(main = "co2 data set") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("crimtab") |
| > ### * crimtab |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: crimtab |
| > ### Title: Student's 3000 Criminals Data |
| > ### Aliases: crimtab |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > dim(crimtab) |
| [1] 42 22 |
| > utils::str(crimtab) |
| 'table' int [1:42, 1:22] 0 0 0 0 0 0 1 0 0 0 ... |
| - attr(*, "dimnames")=List of 2 |
| ..$ : chr [1:42] "9.4" "9.5" "9.6" "9.7" ... |
| ..$ : chr [1:22] "142.24" "144.78" "147.32" "149.86" ... |
| > ## for nicer printing: |
| > local({cT <- crimtab |
| + colnames(cT) <- substring(colnames(cT), 2, 3) |
| + print(cT, zero.print = " ") |
| + }) |
| 42 44 47 49 52 54 57 60 62 65 67 70 72 75 77 80 82 85 87 90 93 95 |
| 9.4 |
| 9.5 1 |
| 9.6 |
| 9.7 |
| 9.8 1 |
| 9.9 1 1 1 |
| 10 1 1 2 2 1 |
| 10.1 1 3 1 1 1 |
| 10.2 2 2 2 1 2 1 |
| 10.3 1 1 3 2 2 3 5 |
| 10.4 1 1 2 3 3 4 3 3 |
| 10.5 1 3 7 6 4 3 1 3 1 1 |
| 10.6 1 4 5 9 14 6 3 1 1 |
| 10.7 1 2 4 9 14 16 15 7 3 1 2 |
| 10.8 2 5 6 14 27 10 7 1 2 1 |
| 10.9 2 6 14 24 27 14 10 4 1 |
| 11 2 6 12 15 31 37 27 17 10 6 |
| 11.1 3 3 12 22 26 24 26 24 7 4 1 |
| 11.2 3 2 7 21 30 38 29 27 20 4 1 1 |
| 11.3 1 5 10 24 26 39 26 24 7 2 |
| 11.4 3 4 9 29 56 58 26 22 10 11 |
| 11.5 5 11 17 33 57 38 34 25 11 2 |
| 11.6 2 1 4 13 37 39 48 38 27 12 2 2 1 |
| 11.7 2 9 17 30 37 48 45 24 9 9 2 |
| 11.8 1 2 11 15 35 41 34 29 10 5 1 |
| 11.9 1 1 2 12 10 27 32 35 19 10 9 3 1 |
| 12 1 4 8 19 42 39 22 16 8 2 2 |
| 12.1 2 4 13 22 28 15 27 10 4 1 |
| 12.2 1 2 5 6 23 17 16 11 8 1 1 |
| 12.3 4 8 10 13 20 23 6 5 |
| 12.4 1 1 1 2 7 12 4 7 7 1 1 |
| 12.5 1 1 3 12 11 8 6 8 2 |
| 12.6 1 3 5 7 8 6 3 1 1 |
| 12.7 1 1 7 5 5 8 2 2 |
| 12.8 1 2 3 1 8 5 3 1 1 |
| 12.9 1 2 2 1 1 |
| 13 3 1 1 2 1 |
| 13.1 1 1 |
| 13.2 1 1 1 3 |
| 13.3 1 1 |
| 13.4 |
| 13.5 1 |
| > |
| > ## Repeat Student's experiment: |
| > |
| > # 1) Reconstitute 3000 raw data for heights in inches and rounded to |
| > # nearest integer as in Student's paper: |
| > |
| > (heIn <- round(as.numeric(colnames(crimtab)) / 2.54)) |
| [1] 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 |
| > d.hei <- data.frame(height = rep(heIn, colSums(crimtab))) |
| > |
| > # 2) shuffle the data: |
| > |
| > set.seed(1) |
| > d.hei <- d.hei[sample(1:3000), , drop = FALSE] |
| > |
| > # 3) Make 750 samples each of size 4: |
| > |
| > d.hei$sample <- as.factor(rep(1:750, each = 4)) |
| > |
| > # 4) Compute the means and standard deviations (n) for the 750 samples: |
| > |
| > h.mean <- with(d.hei, tapply(height, sample, FUN = mean)) |
| > h.sd <- with(d.hei, tapply(height, sample, FUN = sd)) * sqrt(3/4) |
| > |
| > # 5) Compute the difference between the mean of each sample and |
| > # the mean of the population and then divide by the |
| > # standard deviation of the sample: |
| > |
| > zobs <- (h.mean - mean(d.hei[,"height"]))/h.sd |
| > |
| > # 6) Replace infinite values by +/- 6 as in Student's paper: |
| > |
| > zobs[infZ <- is.infinite(zobs)] # none of them |
| named numeric(0) |
| > zobs[infZ] <- 6 * sign(zobs[infZ]) |
| > |
| > # 7) Plot the distribution: |
| > |
| > require(grDevices); require(graphics) |
| > hist(x = zobs, probability = TRUE, xlab = "Student's z", |
| + col = grey(0.8), border = grey(0.5), |
| + main = "Distribution of Student's z score for 'crimtab' data") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("discoveries") |
| > ### * discoveries |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: discoveries |
| > ### Title: Yearly Numbers of Important Discoveries |
| > ### Aliases: discoveries |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(discoveries, ylab = "Number of important discoveries", |
| + las = 1) |
| > title(main = "discoveries data set") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("esoph") |
| > ### * esoph |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: esoph |
| > ### Title: Smoking, Alcohol and (O)esophageal Cancer |
| > ### Aliases: esoph |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > require(graphics) # for mosaicplot |
| > summary(esoph) |
| agegp alcgp tobgp ncases ncontrols |
| 25-34:15 0-39g/day:23 0-9g/day:24 Min. : 0.000 Min. : 1.00 |
| 35-44:15 40-79 :23 10-19 :24 1st Qu.: 0.000 1st Qu.: 3.00 |
| 45-54:16 80-119 :21 20-29 :20 Median : 1.000 Median : 6.00 |
| 55-64:16 120+ :21 30+ :20 Mean : 2.273 Mean :11.08 |
| 65-74:15 3rd Qu.: 4.000 3rd Qu.:14.00 |
| 75+ :11 Max. :17.000 Max. :60.00 |
| > ## effects of alcohol, tobacco and interaction, age-adjusted |
| > model1 <- glm(cbind(ncases, ncontrols) ~ agegp + tobgp * alcgp, |
| + data = esoph, family = binomial()) |
| > anova(model1) |
| Analysis of Deviance Table |
| |
| Model: binomial, link: logit |
| |
| Response: cbind(ncases, ncontrols) |
| |
| Terms added sequentially (first to last) |
| |
| |
| Df Deviance Resid. Df Resid. Dev |
| NULL 87 227.241 |
| agegp 5 88.128 82 139.112 |
| tobgp 3 19.085 79 120.028 |
| alcgp 3 66.054 76 53.973 |
| tobgp:alcgp 9 6.489 67 47.484 |
| > ## Try a linear effect of alcohol and tobacco |
| > model2 <- glm(cbind(ncases, ncontrols) ~ agegp + unclass(tobgp) |
| + + unclass(alcgp), |
| + data = esoph, family = binomial()) |
| > summary(model2) |
| |
| Call: |
| glm(formula = cbind(ncases, ncontrols) ~ agegp + unclass(tobgp) + |
| unclass(alcgp), family = binomial(), data = esoph) |
| |
| Deviance Residuals: |
| Min 1Q Median 3Q Max |
| -1.7628 -0.6426 -0.2709 0.3043 2.0421 |
| |
| Coefficients: |
| Estimate Std. Error z value Pr(>|z|) |
| (Intercept) -4.01097 0.31224 -12.846 < 2e-16 *** |
| agegp.L 2.96113 0.65092 4.549 5.39e-06 *** |
| agegp.Q -1.33735 0.58918 -2.270 0.02322 * |
| agegp.C 0.15292 0.44792 0.341 0.73281 |
| agegp^4 0.06668 0.30776 0.217 0.82848 |
| agegp^5 -0.20288 0.19523 -1.039 0.29872 |
| unclass(tobgp) 0.26162 0.08198 3.191 0.00142 ** |
| unclass(alcgp) 0.65308 0.08452 7.727 1.10e-14 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| (Dispersion parameter for binomial family taken to be 1) |
| |
| Null deviance: 227.241 on 87 degrees of freedom |
| Residual deviance: 59.277 on 80 degrees of freedom |
| AIC: 222.76 |
| |
| Number of Fisher Scoring iterations: 6 |
| |
| > ## Re-arrange data for a mosaic plot |
| > ttt <- table(esoph$agegp, esoph$alcgp, esoph$tobgp) |
| > o <- with(esoph, order(tobgp, alcgp, agegp)) |
| > ttt[ttt == 1] <- esoph$ncases[o] |
| > tt1 <- table(esoph$agegp, esoph$alcgp, esoph$tobgp) |
| > tt1[tt1 == 1] <- esoph$ncontrols[o] |
| > tt <- array(c(ttt, tt1), c(dim(ttt),2), |
| + c(dimnames(ttt), list(c("Cancer", "control")))) |
| > mosaicplot(tt, main = "esoph data set", color = TRUE) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("euro") |
| > ### * euro |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: euro |
| > ### Title: Conversion Rates of Euro Currencies |
| > ### Aliases: euro euro.cross |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > cbind(euro) |
| euro |
| ATS 13.760300 |
| BEF 40.339900 |
| DEM 1.955830 |
| ESP 166.386000 |
| FIM 5.945730 |
| FRF 6.559570 |
| IEP 0.787564 |
| ITL 1936.270000 |
| LUF 40.339900 |
| NLG 2.203710 |
| PTE 200.482000 |
| > |
| > ## These relations hold: |
| > euro == signif(euro, 6) # [6 digit precision in Euro's definition] |
| ATS BEF DEM ESP FIM FRF IEP ITL LUF NLG PTE |
| TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE |
| > all(euro.cross == outer(1/euro, euro)) |
| [1] TRUE |
| > |
| > ## Convert 20 Euro to Belgian Franc |
| > 20 * euro["BEF"] |
| BEF |
| 806.798 |
| > ## Convert 20 Austrian Schilling to Euro |
| > 20 / euro["ATS"] |
| ATS |
| 1.453457 |
| > ## Convert 20 Spanish Pesetas to Italian Lira |
| > 20 * euro.cross["ESP", "ITL"] |
| [1] 232.7443 |
| > |
| > require(graphics) |
| > dotchart(euro, |
| + main = "euro data: 1 Euro in currency unit") |
| > dotchart(1/euro, |
| + main = "euro data: 1 currency unit in Euros") |
| > dotchart(log(euro, 10), |
| + main = "euro data: log10(1 Euro in currency unit)") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("faithful") |
| > ### * faithful |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: faithful |
| > ### Title: Old Faithful Geyser Data |
| > ### Aliases: faithful |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > f.tit <- "faithful data: Eruptions of Old Faithful" |
| > |
| > ne60 <- round(e60 <- 60 * faithful$eruptions) |
| > all.equal(e60, ne60) # relative diff. ~ 1/10000 |
| [1] "Mean relative difference: 9.515332e-05" |
| > table(zapsmall(abs(e60 - ne60))) # 0, 0.02 or 0.04 |
| |
| 0 0.02 0.04 |
| 106 163 3 |
| > faithful$better.eruptions <- ne60 / 60 |
| > te <- table(ne60) |
| > te[te >= 4] # (too) many multiples of 5 ! |
| ne60 |
| 105 108 110 112 113 120 216 230 240 245 249 250 255 260 261 262 265 270 272 275 |
| 6 4 7 8 4 4 4 5 6 5 4 4 4 5 4 4 4 8 5 4 |
| 276 282 288 |
| 4 6 6 |
| > plot(names(te), te, type = "h", main = f.tit, xlab = "Eruption time (sec)") |
| > |
| > plot(faithful[, -3], main = f.tit, |
| + xlab = "Eruption time (min)", |
| + ylab = "Waiting time to next eruption (min)") |
| > lines(lowess(faithful$eruptions, faithful$waiting, f = 2/3, iter = 3), |
| + col = "red") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("freeny") |
| > ### * freeny |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: freeny |
| > ### Title: Freeny's Revenue Data |
| > ### Aliases: freeny freeny.x freeny.y |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > summary(freeny) |
| y lag.quarterly.revenue price.index income.level |
| Min. :8.791 Min. :8.791 Min. :4.278 Min. :5.821 |
| 1st Qu.:9.045 1st Qu.:9.020 1st Qu.:4.392 1st Qu.:5.948 |
| Median :9.314 Median :9.284 Median :4.510 Median :6.061 |
| Mean :9.306 Mean :9.281 Mean :4.496 Mean :6.039 |
| 3rd Qu.:9.591 3rd Qu.:9.561 3rd Qu.:4.605 3rd Qu.:6.139 |
| Max. :9.794 Max. :9.775 Max. :4.710 Max. :6.200 |
| market.potential |
| Min. :12.97 |
| 1st Qu.:13.01 |
| Median :13.07 |
| Mean :13.07 |
| 3rd Qu.:13.12 |
| Max. :13.17 |
| > pairs(freeny, main = "freeny data") |
| > # gives warning: freeny$y has class "ts" |
| > |
| > summary(fm1 <- lm(y ~ ., data = freeny)) |
| |
| Call: |
| lm(formula = y ~ ., data = freeny) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -0.0259426 -0.0101033 0.0003824 0.0103236 0.0267124 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) -10.4726 6.0217 -1.739 0.0911 . |
| lag.quarterly.revenue 0.1239 0.1424 0.870 0.3904 |
| price.index -0.7542 0.1607 -4.693 4.28e-05 *** |
| income.level 0.7675 0.1339 5.730 1.93e-06 *** |
| market.potential 1.3306 0.5093 2.613 0.0133 * |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.01473 on 34 degrees of freedom |
| Multiple R-squared: 0.9981, Adjusted R-squared: 0.9978 |
| F-statistic: 4354 on 4 and 34 DF, p-value: < 2.2e-16 |
| |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0), |
| + mar = c(4.1, 4.1, 2.1, 1.1)) |
| > plot(fm1) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("infert") |
| > ### * infert |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: infert |
| > ### Title: Infertility after Spontaneous and Induced Abortion |
| > ### Aliases: infert |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > model1 <- glm(case ~ spontaneous+induced, data = infert, family = binomial()) |
| > summary(model1) |
| |
| Call: |
| glm(formula = case ~ spontaneous + induced, family = binomial(), |
| data = infert) |
| |
| Deviance Residuals: |
| Min 1Q Median 3Q Max |
| -1.6678 -0.8360 -0.5772 0.9030 1.9362 |
| |
| Coefficients: |
| Estimate Std. Error z value Pr(>|z|) |
| (Intercept) -1.7079 0.2677 -6.380 1.78e-10 *** |
| spontaneous 1.1972 0.2116 5.657 1.54e-08 *** |
| induced 0.4181 0.2056 2.033 0.042 * |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| (Dispersion parameter for binomial family taken to be 1) |
| |
| Null deviance: 316.17 on 247 degrees of freedom |
| Residual deviance: 279.61 on 245 degrees of freedom |
| AIC: 285.61 |
| |
| Number of Fisher Scoring iterations: 4 |
| |
| > ## adjusted for other potential confounders: |
| > summary(model2 <- glm(case ~ age+parity+education+spontaneous+induced, |
| + data = infert, family = binomial())) |
| |
| Call: |
| glm(formula = case ~ age + parity + education + spontaneous + |
| induced, family = binomial(), data = infert) |
| |
| Deviance Residuals: |
| Min 1Q Median 3Q Max |
| -1.7603 -0.8162 -0.4956 0.8349 2.6536 |
| |
| Coefficients: |
| Estimate Std. Error z value Pr(>|z|) |
| (Intercept) -1.14924 1.41220 -0.814 0.4158 |
| age 0.03958 0.03120 1.269 0.2046 |
| parity -0.82828 0.19649 -4.215 2.49e-05 *** |
| education6-11yrs -1.04424 0.79255 -1.318 0.1876 |
| education12+ yrs -1.40321 0.83416 -1.682 0.0925 . |
| spontaneous 2.04591 0.31016 6.596 4.21e-11 *** |
| induced 1.28876 0.30146 4.275 1.91e-05 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| (Dispersion parameter for binomial family taken to be 1) |
| |
| Null deviance: 316.17 on 247 degrees of freedom |
| Residual deviance: 257.80 on 241 degrees of freedom |
| AIC: 271.8 |
| |
| Number of Fisher Scoring iterations: 4 |
| |
| > ## Really should be analysed by conditional logistic regression |
| > ## which is in the survival package |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("iris") |
| > ### * iris |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: iris |
| > ### Title: Edgar Anderson's Iris Data |
| > ### Aliases: iris iris3 |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > dni3 <- dimnames(iris3) |
| > ii <- data.frame(matrix(aperm(iris3, c(1,3,2)), ncol = 4, |
| + dimnames = list(NULL, sub(" L.",".Length", |
| + sub(" W.",".Width", dni3[[2]])))), |
| + Species = gl(3, 50, labels = sub("S", "s", sub("V", "v", dni3[[3]])))) |
| > all.equal(ii, iris) # TRUE |
| [1] TRUE |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("islands") |
| > ### * islands |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: islands |
| > ### Title: Areas of the World's Major Landmasses |
| > ### Aliases: islands |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > dotchart(log(islands, 10), |
| + main = "islands data: log10(area) (log10(sq. miles))") |
| > dotchart(log(islands[order(islands)], 10), |
| + main = "islands data: log10(area) (log10(sq. miles))") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("longley") |
| > ### * longley |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: longley |
| > ### Title: Longley's Economic Regression Data |
| > ### Aliases: longley |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## give the data set in the form it is used in S-PLUS: |
| > longley.x <- data.matrix(longley[, 1:6]) |
| > longley.y <- longley[, "Employed"] |
| > pairs(longley, main = "longley data") |
| > summary(fm1 <- lm(Employed ~ ., data = longley)) |
| |
| Call: |
| lm(formula = Employed ~ ., data = longley) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -0.41011 -0.15767 -0.02816 0.10155 0.45539 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) -3.482e+03 8.904e+02 -3.911 0.003560 ** |
| GNP.deflator 1.506e-02 8.492e-02 0.177 0.863141 |
| GNP -3.582e-02 3.349e-02 -1.070 0.312681 |
| Unemployed -2.020e-02 4.884e-03 -4.136 0.002535 ** |
| Armed.Forces -1.033e-02 2.143e-03 -4.822 0.000944 *** |
| Population -5.110e-02 2.261e-01 -0.226 0.826212 |
| Year 1.829e+00 4.555e-01 4.016 0.003037 ** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.3049 on 9 degrees of freedom |
| Multiple R-squared: 0.9955, Adjusted R-squared: 0.9925 |
| F-statistic: 330.3 on 6 and 9 DF, p-value: 4.984e-10 |
| |
| > opar <- par(mfrow = c(2, 2), oma = c(0, 0, 1.1, 0), |
| + mar = c(4.1, 4.1, 2.1, 1.1)) |
| > plot(fm1) |
| > par(opar) |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("morley") |
| > ### * morley |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: morley |
| > ### Title: Michelson Speed of Light Data |
| > ### Aliases: morley |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > michelson <- transform(morley, |
| + Expt = factor(Expt), Run = factor(Run)) |
| > xtabs(~ Expt + Run, data = michelson) # 5 x 20 balanced (two-way) |
| Run |
| Expt 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 |
| 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 |
| 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 |
| 3 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 |
| 4 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 |
| 5 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 |
| > plot(Speed ~ Expt, data = michelson, |
| + main = "Speed of Light Data", xlab = "Experiment No.") |
| > fm <- aov(Speed ~ Run + Expt, data = michelson) |
| > summary(fm) |
| Df Sum Sq Mean Sq F value Pr(>F) |
| Run 19 113344 5965 1.105 0.36321 |
| Expt 4 94514 23629 4.378 0.00307 ** |
| Residuals 76 410166 5397 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > fm0 <- update(fm, . ~ . - Run) |
| > anova(fm0, fm) |
| Analysis of Variance Table |
| |
| Model 1: Speed ~ Expt |
| Model 2: Speed ~ Run + Expt |
| Res.Df RSS Df Sum of Sq F Pr(>F) |
| 1 95 523510 |
| 2 76 410166 19 113344 1.1053 0.3632 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("mtcars") |
| > ### * mtcars |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: mtcars |
| > ### Title: Motor Trend Car Road Tests |
| > ### Aliases: mtcars |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > pairs(mtcars, main = "mtcars data", gap = 1/4) |
| > coplot(mpg ~ disp | as.factor(cyl), data = mtcars, |
| + panel = panel.smooth, rows = 1) |
| > ## possibly more meaningful, e.g., for summary() or bivariate plots: |
| > mtcars2 <- within(mtcars, { |
| + vs <- factor(vs, labels = c("V", "S")) |
| + am <- factor(am, labels = c("automatic", "manual")) |
| + cyl <- ordered(cyl) |
| + gear <- ordered(gear) |
| + carb <- ordered(carb) |
| + }) |
| > summary(mtcars2) |
| mpg cyl disp hp drat |
| Min. :10.40 4:11 Min. : 71.1 Min. : 52.0 Min. :2.760 |
| 1st Qu.:15.43 6: 7 1st Qu.:120.8 1st Qu.: 96.5 1st Qu.:3.080 |
| Median :19.20 8:14 Median :196.3 Median :123.0 Median :3.695 |
| Mean :20.09 Mean :230.7 Mean :146.7 Mean :3.597 |
| 3rd Qu.:22.80 3rd Qu.:326.0 3rd Qu.:180.0 3rd Qu.:3.920 |
| Max. :33.90 Max. :472.0 Max. :335.0 Max. :4.930 |
| wt qsec vs am gear carb |
| Min. :1.513 Min. :14.50 V:18 automatic:19 3:15 1: 7 |
| 1st Qu.:2.581 1st Qu.:16.89 S:14 manual :13 4:12 2:10 |
| Median :3.325 Median :17.71 5: 5 3: 3 |
| Mean :3.217 Mean :17.85 4:10 |
| 3rd Qu.:3.610 3rd Qu.:18.90 6: 1 |
| Max. :5.424 Max. :22.90 8: 1 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("nhtemp") |
| > ### * nhtemp |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: nhtemp |
| > ### Title: Average Yearly Temperatures in New Haven |
| > ### Aliases: nhtemp |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > plot(nhtemp, main = "nhtemp data", |
| + ylab = "Mean annual temperature in New Haven, CT (deg. F)") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("nottem") |
| > ### * nottem |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: nottem |
| > ### Title: Average Monthly Temperatures at Nottingham, 1920-1939 |
| > ### Aliases: nottem |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("npk") |
| > ### * npk |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: npk |
| > ### Title: Classical N, P, K Factorial Experiment |
| > ### Aliases: npk |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > |
| > cleanEx() |
| > nameEx("occupationalStatus") |
| > ### * occupationalStatus |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: occupationalStatus |
| > ### Title: Occupational Status of Fathers and their Sons |
| > ### Aliases: occupationalStatus |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > |
| > plot(occupationalStatus) |
| > |
| > ## Fit a uniform association model separating diagonal effects |
| > Diag <- as.factor(diag(1:8)) |
| > Rscore <- scale(as.numeric(row(occupationalStatus)), scale = FALSE) |
| > Cscore <- scale(as.numeric(col(occupationalStatus)), scale = FALSE) |
| > modUnif <- glm(Freq ~ origin + destination + Diag + Rscore:Cscore, |
| + family = poisson, data = occupationalStatus) |
| > |
| > summary(modUnif) |
| |
| Call: |
| glm(formula = Freq ~ origin + destination + Diag + Rscore:Cscore, |
| family = poisson, data = occupationalStatus) |
| |
| Deviance Residuals: |
| Min 1Q Median 3Q Max |
| -2.6521 -0.6267 0.0000 0.5913 2.0964 |
| |
| Coefficients: |
| Estimate Std. Error z value Pr(>|z|) |
| (Intercept) 0.568592 0.183358 3.101 0.001929 ** |
| origin2 0.431314 0.149415 2.887 0.003893 ** |
| origin3 1.461862 0.131141 11.147 < 2e-16 *** |
| origin4 1.788731 0.126588 14.130 < 2e-16 *** |
| origin5 0.441011 0.144844 3.045 0.002329 ** |
| origin6 2.491735 0.121219 20.556 < 2e-16 *** |
| origin7 1.127536 0.129032 8.738 < 2e-16 *** |
| origin8 0.796445 0.131863 6.040 1.54e-09 *** |
| destination2 0.873202 0.166902 5.232 1.68e-07 *** |
| destination3 1.813992 0.153861 11.790 < 2e-16 *** |
| destination4 2.082515 0.150887 13.802 < 2e-16 *** |
| destination5 1.366383 0.155590 8.782 < 2e-16 *** |
| destination6 2.816369 0.146054 19.283 < 2e-16 *** |
| destination7 1.903918 0.147810 12.881 < 2e-16 *** |
| destination8 1.398585 0.151658 9.222 < 2e-16 *** |
| Diag1 1.665495 0.237383 7.016 2.28e-12 *** |
| Diag2 0.959681 0.212122 4.524 6.06e-06 *** |
| Diag3 0.021750 0.156530 0.139 0.889490 |
| Diag4 0.226399 0.124364 1.820 0.068689 . |
| Diag5 0.808646 0.229754 3.520 0.000432 *** |
| Diag6 0.132277 0.077191 1.714 0.086597 . |
| Diag7 0.506709 0.115936 4.371 1.24e-05 *** |
| Diag8 0.221880 0.134803 1.646 0.099771 . |
| Rscore:Cscore 0.136974 0.007489 18.289 < 2e-16 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| (Dispersion parameter for poisson family taken to be 1) |
| |
| Null deviance: 4679.004 on 63 degrees of freedom |
| Residual deviance: 58.436 on 40 degrees of freedom |
| AIC: 428.78 |
| |
| Number of Fisher Scoring iterations: 4 |
| |
| > plot(modUnif) # 4 plots, with warning about h_ii ~= 1 |
| Warning: not plotting observations with leverage one: |
| 1, 10, 19, 28, 37, 46, 55, 64 |
| Warning: not plotting observations with leverage one: |
| 1, 10, 19, 28, 37, 46, 55, 64 |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("precip") |
| > ### * precip |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: precip |
| > ### Title: Annual Precipitation in US Cities |
| > ### Aliases: precip |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > dotchart(precip[order(precip)], main = "precip data") |
| > title(sub = "Average annual precipitation (in.)") |
| > |
| > ## Old ("wrong") version of dataset (just name change): |
| > precip.O <- local({ |
| + p <- precip; names(p)[names(p) == "Cincinnati"] <- "Cincinati" ; p }) |
| > stopifnot(all(precip == precip.O), |
| + match("Cincinnati", names(precip)) == 46, |
| + identical(names(precip)[-46], names(precip.O)[-46])) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("presidents") |
| > ### * presidents |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: presidents |
| > ### Title: Quarterly Approval Ratings of US Presidents |
| > ### Aliases: presidents |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > plot(presidents, las = 1, ylab = "Approval rating (%)", |
| + main = "presidents data") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("pressure") |
| > ### * pressure |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: pressure |
| > ### Title: Vapor Pressure of Mercury as a Function of Temperature |
| > ### Aliases: pressure |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(pressure, xlab = "Temperature (deg C)", |
| + ylab = "Pressure (mm of Hg)", |
| + main = "pressure data: Vapor Pressure of Mercury") |
| > plot(pressure, xlab = "Temperature (deg C)", log = "y", |
| + ylab = "Pressure (mm of Hg)", |
| + main = "pressure data: Vapor Pressure of Mercury") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("quakes") |
| > ### * quakes |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: quakes |
| > ### Title: Locations of Earthquakes off Fiji |
| > ### Aliases: quakes |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > pairs(quakes, main = "Fiji Earthquakes, N = 1000", cex.main = 1.2, pch = ".") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("randu") |
| > ### * randu |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: randu |
| > ### Title: Random Numbers from Congruential Generator RANDU |
| > ### Aliases: randu |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("sleep") |
| > ### * sleep |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: sleep |
| > ### Title: Student's Sleep Data |
| > ### Aliases: sleep |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > ## Student's paired t-test |
| > with(sleep, |
| + t.test(extra[group == 1], |
| + extra[group == 2], paired = TRUE)) |
| |
| Paired t-test |
| |
| data: extra[group == 1] and extra[group == 2] |
| t = -4.0621, df = 9, p-value = 0.002833 |
| alternative hypothesis: true difference in means is not equal to 0 |
| 95 percent confidence interval: |
| -2.4598858 -0.7001142 |
| sample estimates: |
| mean of the differences |
| -1.58 |
| |
| > |
| > ## The sleep *prolongations* |
| > sleep1 <- with(sleep, extra[group == 2] - extra[group == 1]) |
| > summary(sleep1) |
| Min. 1st Qu. Median Mean 3rd Qu. Max. |
| 0.00 1.05 1.30 1.58 1.70 4.60 |
| > stripchart(sleep1, method = "stack", xlab = "hours", |
| + main = "Sleep prolongation (n = 10)") |
| > boxplot(sleep1, horizontal = TRUE, add = TRUE, |
| + at = .6, pars = list(boxwex = 0.5, staplewex = 0.25)) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("stackloss") |
| > ### * stackloss |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: stackloss |
| > ### Title: Brownlee's Stack Loss Plant Data |
| > ### Aliases: stackloss stack.loss stack.x |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats) |
| > summary(lm.stack <- lm(stack.loss ~ stack.x)) |
| |
| Call: |
| lm(formula = stack.loss ~ stack.x) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -7.2377 -1.7117 -0.4551 2.3614 5.6978 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) -39.9197 11.8960 -3.356 0.00375 ** |
| stack.xAir.Flow 0.7156 0.1349 5.307 5.8e-05 *** |
| stack.xWater.Temp 1.2953 0.3680 3.520 0.00263 ** |
| stack.xAcid.Conc. -0.1521 0.1563 -0.973 0.34405 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 3.243 on 17 degrees of freedom |
| Multiple R-squared: 0.9136, Adjusted R-squared: 0.8983 |
| F-statistic: 59.9 on 3 and 17 DF, p-value: 3.016e-09 |
| |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("sunspot.month") |
| > ### * sunspot.month |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: sunspot.month |
| > ### Title: Monthly Sunspot Data, from 1749 to "Present" |
| > ### Aliases: sunspot.month |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## Compare the monthly series |
| > plot (sunspot.month, |
| + main="sunspot.month & sunspots [package'datasets']", col=2) |
| > lines(sunspots) # -> faint differences where they overlap |
| > |
| > ## Now look at the difference : |
| > all(tsp(sunspots) [c(1,3)] == |
| + tsp(sunspot.month)[c(1,3)]) ## Start & Periodicity are the same |
| [1] TRUE |
| > n1 <- length(sunspots) |
| > table(eq <- sunspots == sunspot.month[1:n1]) #> 132 are different ! |
| |
| FALSE TRUE |
| 143 2677 |
| > i <- which(!eq) |
| > rug(time(eq)[i]) |
| > s1 <- sunspots[i] ; s2 <- sunspot.month[i] |
| > cbind(i = i, time = time(sunspots)[i], sunspots = s1, ss.month = s2, |
| + perc.diff = round(100*2*abs(s1-s2)/(s1+s2), 1)) |
| i time sunspots ss.month perc.diff |
| [1,] 55 1753.500 22.2 22.0 0.9 |
| [2,] 838 1818.750 31.7 31.6 0.3 |
| [3,] 841 1819.000 32.5 32.8 0.9 |
| [4,] 862 1820.750 9.0 8.9 1.1 |
| [5,] 864 1820.917 9.7 9.1 6.4 |
| [6,] 866 1821.083 4.3 4.2 2.4 |
| [7,] 876 1821.917 0.0 0.2 200.0 |
| [8,] 901 1824.000 21.6 21.7 0.5 |
| [9,] 917 1825.333 15.4 15.5 0.6 |
| [10,] 920 1825.583 25.4 25.7 1.2 |
| [11,] 943 1827.500 42.9 42.3 1.4 |
| [12,] 946 1827.750 57.2 56.1 1.9 |
| [13,] 955 1828.500 54.3 54.2 0.2 |
| [14,] 960 1828.917 46.6 46.9 0.6 |
| [15,] 965 1829.333 67.5 67.4 0.1 |
| [16,] 968 1829.583 78.3 77.6 0.9 |
| [17,] 976 1830.250 107.1 106.3 0.7 |
| [18,] 988 1831.250 54.6 54.5 0.2 |
| [19,] 992 1831.583 54.9 55.0 0.2 |
| [20,] 994 1831.750 46.2 46.3 0.2 |
| [21,] 998 1832.083 55.5 55.6 0.2 |
| [22,] 1003 1832.500 13.9 14.0 0.7 |
| [23,] 1047 1836.167 98.1 98.2 0.1 |
| [24,] 1061 1837.333 111.3 111.7 0.4 |
| [25,] 1081 1839.000 107.6 105.6 1.9 |
| [26,] 1087 1839.500 84.7 84.8 0.1 |
| [27,] 1090 1839.750 90.8 90.9 0.1 |
| [28,] 1092 1839.917 63.6 63.7 0.2 |
| [29,] 1095 1840.167 55.5 67.8 20.0 |
| [30,] 1102 1840.750 49.8 55.0 9.9 |
| [31,] 1105 1841.000 24.0 24.1 0.4 |
| [32,] 1108 1841.250 42.6 40.2 5.8 |
| [33,] 1109 1841.333 67.4 67.5 0.1 |
| [34,] 1113 1841.667 35.1 36.5 3.9 |
| [35,] 1124 1842.583 26.5 26.6 0.4 |
| [36,] 1125 1842.667 18.5 18.4 0.5 |
| [37,] 1132 1843.250 8.8 9.5 7.7 |
| [38,] 1145 1844.333 12.0 11.6 3.4 |
| [39,] 1149 1844.667 6.9 7.0 1.4 |
| [40,] 1156 1845.250 56.9 57.0 0.2 |
| [41,] 1168 1846.250 69.2 69.3 0.1 |
| [42,] 1185 1847.667 161.2 160.9 0.2 |
| [43,] 1191 1848.167 108.9 108.6 0.3 |
| [44,] 1194 1848.417 123.8 129.0 4.1 |
| [45,] 1196 1848.583 132.5 132.6 0.1 |
| [46,] 1200 1848.917 159.9 159.5 0.3 |
| [47,] 1201 1849.000 156.7 157.0 0.2 |
| [48,] 1202 1849.083 131.7 131.8 0.1 |
| [49,] 1203 1849.167 96.5 96.2 0.3 |
| [50,] 1206 1849.417 81.2 81.1 0.1 |
| [51,] 1208 1849.583 61.3 67.7 9.9 |
| [52,] 1211 1849.833 99.7 99.0 0.7 |
| [53,] 1224 1850.917 60.0 61.0 1.7 |
| [54,] 1235 1851.833 50.9 51.0 0.2 |
| [55,] 1238 1852.083 67.5 66.4 1.6 |
| [56,] 1243 1852.500 42.0 42.1 0.2 |
| [57,] 1256 1853.583 50.4 50.5 0.2 |
| [58,] 1258 1853.750 42.3 42.4 0.2 |
| [59,] 1264 1854.250 26.4 26.5 0.4 |
| [60,] 1270 1854.750 12.7 12.6 0.8 |
| [61,] 1272 1854.917 21.4 21.6 0.9 |
| [62,] 1282 1855.750 9.7 9.6 1.0 |
| [63,] 1283 1855.833 4.3 4.2 2.4 |
| [64,] 1290 1856.417 5.0 5.2 3.9 |
| [65,] 1301 1857.333 29.2 28.5 2.4 |
| [66,] 1333 1860.000 81.5 82.4 1.1 |
| [67,] 1334 1860.083 88.0 88.3 0.3 |
| [68,] 1346 1861.083 77.8 77.7 0.1 |
| [69,] 1350 1861.417 87.8 88.1 0.3 |
| [70,] 1366 1862.750 42.0 41.9 0.2 |
| [71,] 1407 1866.167 24.6 24.5 0.4 |
| [72,] 1424 1867.583 4.9 4.8 2.1 |
| [73,] 1427 1867.833 9.3 9.6 3.2 |
| [74,] 1429 1868.000 15.6 15.5 0.6 |
| [75,] 1430 1868.083 15.8 15.7 0.6 |
| [76,] 1435 1868.500 28.6 29.0 1.4 |
| [77,] 1437 1868.667 43.8 47.2 7.5 |
| [78,] 1438 1868.750 61.7 61.6 0.2 |
| [79,] 1442 1869.083 59.3 59.9 1.0 |
| [80,] 1445 1869.333 104.0 103.9 0.1 |
| [81,] 1450 1869.750 59.4 59.3 0.2 |
| [82,] 1451 1869.833 77.4 78.1 0.9 |
| [83,] 1452 1869.917 104.3 104.4 0.1 |
| [84,] 1455 1870.167 159.4 157.5 1.2 |
| [85,] 1472 1871.583 110.0 110.1 0.1 |
| [86,] 1476 1871.917 90.3 90.4 0.1 |
| [87,] 1486 1872.750 103.5 102.6 0.9 |
| [88,] 1497 1873.667 47.5 47.1 0.8 |
| [89,] 1498 1873.750 47.4 47.1 0.6 |
| [90,] 1514 1875.083 22.2 21.5 3.2 |
| [91,] 1527 1876.167 31.2 30.6 1.9 |
| [92,] 1539 1877.167 11.7 11.9 1.7 |
| [93,] 1541 1877.333 21.2 21.6 1.9 |
| [94,] 1542 1877.417 13.4 14.2 5.8 |
| [95,] 1543 1877.500 5.9 6.0 1.7 |
| [96,] 1545 1877.667 16.4 16.9 3.0 |
| [97,] 1547 1877.833 14.5 14.2 2.1 |
| [98,] 1548 1877.917 2.3 2.2 4.4 |
| [99,] 1550 1878.083 6.0 6.6 9.5 |
| [100,] 1553 1878.333 5.8 5.9 1.7 |
| [101,] 1561 1879.000 0.8 1.0 22.2 |
| [102,] 1571 1879.833 12.9 13.1 1.5 |
| [103,] 1572 1879.917 7.2 7.3 1.4 |
| [104,] 1574 1880.083 27.5 27.2 1.1 |
| [105,] 1575 1880.167 19.5 19.3 1.0 |
| [106,] 1576 1880.250 19.3 19.5 1.0 |
| [107,] 1588 1881.250 51.7 51.6 0.2 |
| [108,] 1592 1881.583 58.0 58.4 0.7 |
| [109,] 1594 1881.750 64.0 64.4 0.6 |
| [110,] 1598 1882.083 69.3 69.5 0.3 |
| [111,] 1599 1882.167 67.5 66.8 1.0 |
| [112,] 1613 1883.333 32.1 31.5 1.9 |
| [113,] 1614 1883.417 76.5 76.3 0.3 |
| [114,] 1623 1884.167 86.8 87.5 0.8 |
| [115,] 1643 1885.833 33.3 30.9 7.5 |
| [116,] 1656 1886.917 12.4 13.0 4.7 |
| [117,] 1663 1887.500 23.3 23.4 0.4 |
| [118,] 1683 1889.167 7.0 6.7 4.4 |
| [119,] 1687 1889.500 9.7 9.4 3.1 |
| [120,] 1712 1891.583 33.2 33.0 0.6 |
| [121,] 1716 1891.917 32.3 32.5 0.6 |
| [122,] 1723 1892.500 76.8 76.5 0.4 |
| [123,] 1734 1893.417 88.2 89.9 1.9 |
| [124,] 1735 1893.500 88.8 88.6 0.2 |
| [125,] 1738 1893.750 79.7 80.0 0.4 |
| [126,] 1774 1896.750 28.4 28.7 1.1 |
| [127,] 1837 1902.000 5.2 5.5 5.6 |
| [128,] 2126 1926.083 70.0 69.9 0.1 |
| [129,] 2151 1928.167 85.4 85.5 0.1 |
| [130,] 2153 1928.333 76.9 77.0 0.1 |
| [131,] 2162 1929.083 64.1 62.8 2.0 |
| [132,] 2174 1930.083 49.2 49.9 1.4 |
| [133,] 2233 1935.000 18.9 18.6 1.6 |
| [134,] 2315 1941.833 38.3 38.4 0.3 |
| [135,] 2329 1943.000 12.4 12.5 0.8 |
| [136,] 2378 1947.083 113.4 133.4 16.2 |
| [137,] 2427 1951.167 59.9 55.9 6.9 |
| [138,] 2498 1957.083 130.2 130.3 0.1 |
| [139,] 2552 1961.583 55.9 55.8 0.2 |
| [140,] 2556 1961.917 40.0 39.9 0.3 |
| [141,] 2594 1965.083 14.2 14.3 0.7 |
| [142,] 2790 1981.417 90.0 90.9 1.0 |
| [143,] 2819 1983.833 33.3 33.4 0.3 |
| > |
| > ## How to recreate the "old" sunspot.month (R <= 3.0.3): |
| > .sunspot.diff <- cbind( |
| + i = c(1202L, 1256L, 1258L, 1301L, 1407L, 1429L, 1452L, 1455L, |
| + 1663L, 2151L, 2329L, 2498L, 2594L, 2694L, 2819L), |
| + res10 = c(1L, 1L, 1L, -1L, -1L, -1L, 1L, -1L, |
| + 1L, 1L, 1L, 1L, 1L, 20L, 1L)) |
| > ssm0 <- sunspot.month[1:2988] |
| > with(as.data.frame(.sunspot.diff), ssm0[i] <<- ssm0[i] - res10/10) |
| > sunspot.month.0 <- ts(ssm0, start = 1749, frequency = 12) |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("sunspot.year") |
| > ### * sunspot.year |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: sunspot.year |
| > ### Title: Yearly Sunspot Data, 1700-1988 |
| > ### Aliases: sunspot.year |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > utils::str(sm <- sunspots)# the monthly version we keep unchanged |
| Time-Series [1:2820] from 1749 to 1984: 58 62.6 70 55.7 85 83.5 94.8 66.3 75.9 75.5 ... |
| > utils::str(sy <- sunspot.year) |
| Time-Series [1:289] from 1700 to 1988: 5 11 16 23 36 58 29 20 10 8 ... |
| > ## The common time interval |
| > (t1 <- c(max(start(sm), start(sy)), 1)) # Jan 1749 |
| [1] 1749 1 |
| > (t2 <- c(min( end(sm)[1],end(sy)[1]), 12)) # Dec 1983 |
| [1] 1983 12 |
| > s.m <- window(sm, start=t1, end=t2) |
| > s.y <- window(sy, start=t1, end=t2[1]) # {irrelevant warning} |
| > stopifnot(length(s.y) * 12 == length(s.m), |
| + ## The yearly series *is* close to the averages of the monthly one: |
| + all.equal(s.y, aggregate(s.m, FUN = mean), tol = 0.0020)) |
| > ## NOTE: Strangely, correctly weighting the number of days per month |
| > ## (using 28.25 for February) is *not* closer than the simple mean: |
| > ndays <- c(31, 28.25, rep(c(31,30, 31,30, 31), 2)) |
| > all.equal(s.y, aggregate(s.m, FUN = mean)) # 0.0013 |
| [1] "Mean relative difference: 0.001312539" |
| > all.equal(s.y, aggregate(s.m, FUN = weighted.mean, w = ndays)) # 0.0017 |
| [1] "Mean relative difference: 0.001692215" |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("sunspots") |
| > ### * sunspots |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: sunspots |
| > ### Title: Monthly Sunspot Numbers, 1749-1983 |
| > ### Aliases: sunspots |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(sunspots, main = "sunspots data", xlab = "Year", |
| + ylab = "Monthly sunspot numbers") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("swiss") |
| > ### * swiss |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: swiss |
| > ### Title: Swiss Fertility and Socioeconomic Indicators (1888) Data |
| > ### Aliases: swiss |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > pairs(swiss, panel = panel.smooth, main = "swiss data", |
| + col = 3 + (swiss$Catholic > 50)) |
| > summary(lm(Fertility ~ . , data = swiss)) |
| |
| Call: |
| lm(formula = Fertility ~ ., data = swiss) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -15.2743 -5.2617 0.5032 4.1198 15.3213 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 66.91518 10.70604 6.250 1.91e-07 *** |
| Agriculture -0.17211 0.07030 -2.448 0.01873 * |
| Examination -0.25801 0.25388 -1.016 0.31546 |
| Education -0.87094 0.18303 -4.758 2.43e-05 *** |
| Catholic 0.10412 0.03526 2.953 0.00519 ** |
| Infant.Mortality 1.07705 0.38172 2.822 0.00734 ** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 7.165 on 41 degrees of freedom |
| Multiple R-squared: 0.7067, Adjusted R-squared: 0.671 |
| F-statistic: 19.76 on 5 and 41 DF, p-value: 5.594e-10 |
| |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("trees") |
| > ### * trees |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: trees |
| > ### Title: Diameter, Height and Volume for Black Cherry Trees |
| > ### Aliases: trees |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > pairs(trees, panel = panel.smooth, main = "trees data") |
| > plot(Volume ~ Girth, data = trees, log = "xy") |
| > coplot(log(Volume) ~ log(Girth) | Height, data = trees, |
| + panel = panel.smooth) |
| > summary(fm1 <- lm(log(Volume) ~ log(Girth), data = trees)) |
| |
| Call: |
| lm(formula = log(Volume) ~ log(Girth), data = trees) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -0.205999 -0.068702 0.001011 0.072585 0.247963 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) -2.35332 0.23066 -10.20 4.18e-11 *** |
| log(Girth) 2.19997 0.08983 24.49 < 2e-16 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.115 on 29 degrees of freedom |
| Multiple R-squared: 0.9539, Adjusted R-squared: 0.9523 |
| F-statistic: 599.7 on 1 and 29 DF, p-value: < 2.2e-16 |
| |
| > summary(fm2 <- update(fm1, ~ . + log(Height), data = trees)) |
| |
| Call: |
| lm(formula = log(Volume) ~ log(Girth) + log(Height), data = trees) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -0.168561 -0.048488 0.002431 0.063637 0.129223 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) -6.63162 0.79979 -8.292 5.06e-09 *** |
| log(Girth) 1.98265 0.07501 26.432 < 2e-16 *** |
| log(Height) 1.11712 0.20444 5.464 7.81e-06 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 0.08139 on 28 degrees of freedom |
| Multiple R-squared: 0.9777, Adjusted R-squared: 0.9761 |
| F-statistic: 613.2 on 2 and 28 DF, p-value: < 2.2e-16 |
| |
| > step(fm2) |
| Start: AIC=-152.69 |
| log(Volume) ~ log(Girth) + log(Height) |
| |
| Df Sum of Sq RSS AIC |
| <none> 0.1855 -152.685 |
| - log(Height) 1 0.1978 0.3832 -132.185 |
| - log(Girth) 1 4.6275 4.8130 -53.743 |
| |
| Call: |
| lm(formula = log(Volume) ~ log(Girth) + log(Height), data = trees) |
| |
| Coefficients: |
| (Intercept) log(Girth) log(Height) |
| -6.632 1.983 1.117 |
| |
| > ## i.e., Volume ~= c * Height * Girth^2 seems reasonable |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("uspop") |
| > ### * uspop |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: uspop |
| > ### Title: Populations Recorded by the US Census |
| > ### Aliases: uspop |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(uspop, log = "y", main = "uspop data", xlab = "Year", |
| + ylab = "U.S. Population (millions)") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("volcano") |
| > ### * volcano |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: volcano |
| > ### Title: Topographic Information on Auckland's Maunga Whau Volcano |
| > ### Aliases: volcano |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(grDevices); require(graphics) |
| > filled.contour(volcano, color.palette = terrain.colors, asp = 1) |
| > title(main = "volcano data: filled contour map") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("warpbreaks") |
| > ### * warpbreaks |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: warpbreaks |
| > ### Title: The Number of Breaks in Yarn during Weaving |
| > ### Aliases: warpbreaks |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > summary(warpbreaks) |
| breaks wool tension |
| Min. :10.00 A:27 L:18 |
| 1st Qu.:18.25 B:27 M:18 |
| Median :26.00 H:18 |
| Mean :28.15 |
| 3rd Qu.:34.00 |
| Max. :70.00 |
| > opar <- par(mfrow = c(1, 2), oma = c(0, 0, 1.1, 0)) |
| > plot(breaks ~ tension, data = warpbreaks, col = "lightgray", |
| + varwidth = TRUE, subset = wool == "A", main = "Wool A") |
| > plot(breaks ~ tension, data = warpbreaks, col = "lightgray", |
| + varwidth = TRUE, subset = wool == "B", main = "Wool B") |
| > mtext("warpbreaks data", side = 3, outer = TRUE) |
| > par(opar) |
| > summary(fm1 <- lm(breaks ~ wool*tension, data = warpbreaks)) |
| |
| Call: |
| lm(formula = breaks ~ wool * tension, data = warpbreaks) |
| |
| Residuals: |
| Min 1Q Median 3Q Max |
| -19.5556 -6.8889 -0.6667 7.1944 25.4444 |
| |
| Coefficients: |
| Estimate Std. Error t value Pr(>|t|) |
| (Intercept) 44.556 3.647 12.218 2.43e-16 *** |
| woolB -16.333 5.157 -3.167 0.002677 ** |
| tensionM -20.556 5.157 -3.986 0.000228 *** |
| tensionH -20.000 5.157 -3.878 0.000320 *** |
| woolB:tensionM 21.111 7.294 2.895 0.005698 ** |
| woolB:tensionH 10.556 7.294 1.447 0.154327 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 10.94 on 48 degrees of freedom |
| Multiple R-squared: 0.3778, Adjusted R-squared: 0.3129 |
| F-statistic: 5.828 on 5 and 48 DF, p-value: 0.0002772 |
| |
| > anova(fm1) |
| Analysis of Variance Table |
| |
| Response: breaks |
| Df Sum Sq Mean Sq F value Pr(>F) |
| wool 1 450.7 450.67 3.7653 0.0582130 . |
| tension 2 2034.3 1017.13 8.4980 0.0006926 *** |
| wool:tension 2 1002.8 501.39 4.1891 0.0210442 * |
| Residuals 48 5745.1 119.69 |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| > |
| > |
| > |
| > graphics::par(get("par.postscript", pos = 'CheckExEnv')) |
| > cleanEx() |
| > nameEx("women") |
| > ### * women |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: women |
| > ### Title: Average Heights and Weights for American Women |
| > ### Aliases: women |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(graphics) |
| > plot(women, xlab = "Height (in)", ylab = "Weight (lb)", |
| + main = "women data: American women aged 30-39") |
| > |
| > |
| > |
| > cleanEx() |
| > nameEx("zCO2") |
| > ### * zCO2 |
| > |
| > flush(stderr()); flush(stdout()) |
| > |
| > ### Name: CO2 |
| > ### Title: Carbon Dioxide Uptake in Grass Plants |
| > ### Aliases: CO2 |
| > ### Keywords: datasets |
| > |
| > ### ** Examples |
| > |
| > require(stats); require(graphics) |
| > ## Don't show: |
| > options(show.nls.convergence=FALSE) |
| > ## End(Don't show) |
| > coplot(uptake ~ conc | Plant, data = CO2, show.given = FALSE, type = "b") |
| > ## fit the data for the first plant |
| > fm1 <- nls(uptake ~ SSasymp(conc, Asym, lrc, c0), |
| + data = CO2, subset = Plant == "Qn1") |
| > summary(fm1) |
| |
| Formula: uptake ~ SSasymp(conc, Asym, lrc, c0) |
| |
| Parameters: |
| Estimate Std. Error t value Pr(>|t|) |
| Asym 38.1398 0.9164 41.620 1.99e-06 *** |
| lrc -34.2766 18.9661 -1.807 0.145 |
| c0 -4.3806 0.2042 -21.457 2.79e-05 *** |
| --- |
| Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
| |
| Residual standard error: 1.663 on 4 degrees of freedom |
| |
| > ## fit each plant separately |
| > fmlist <- list() |
| > for (pp in levels(CO2$Plant)) { |
| + fmlist[[pp]] <- nls(uptake ~ SSasymp(conc, Asym, lrc, c0), |
| + data = CO2, subset = Plant == pp) |
| + } |
| > ## check the coefficients by plant |
| > print(sapply(fmlist, coef), digits = 3) |
| Qn1 Qn2 Qn3 Qc1 Qc3 Qc2 Mn3 Mn2 Mn1 Mc2 Mc3 |
| Asym 38.14 42.87 44.23 36.43 40.68 39.82 28.48 32.13 34.08 13.56 18.54 |
| lrc -34.28 -29.66 -37.63 -9.90 -11.54 -51.53 -17.37 -29.04 -8.81 -1.98 -136.11 |
| c0 -4.38 -4.67 -4.49 -4.86 -4.95 -4.46 -4.59 -4.47 -5.06 -4.56 -3.47 |
| Mc1 |
| Asym 21.79 |
| lrc 2.45 |
| c0 -5.14 |
| > |
| > |
| > |
| > ### * <FOOTER> |
| > ### |
| > cleanEx() |
| > options(digits = 7L) |
| > base::cat("Time elapsed: ", proc.time() - base::get("ptime", pos = 'CheckExEnv'),"\n") |
| Time elapsed: 1.974 0.053 2.037 0 0 |
| > grDevices::dev.off() |
| null device |
| 1 |
| > ### |
| > ### Local variables: *** |
| > ### mode: outline-minor *** |
| > ### outline-regexp: "\\(> \\)?### [*]+" *** |
| > ### End: *** |
| > quit('no') |
