| Title: | Functions for University of Auckland Course STATS 201/208 Data Analysis |
|---|---|
| Description: | A set of functions used in teaching STATS 201/208 Data Analysis at the University of Auckland. The functions are designed to make parts of R more accessible to a large undergraduate population who are mostly not statistics majors. |
| Authors: | Brant Deppa [aut] (Wrote the original R scripts this package is derived from), James Curran [aut, cre] (Wrote the original R package. Current maintainer.), Hannah Yun [ctb], Rachel Fewster [ctb], Russell Millar [ctb], Ben Stevenson [ctb], Andrew Balemi [ctb], Chris Wild [ctb], Sophie Jones [ctb], Dineika Chandra [ctr], Brendan McArdle [ctr] |
| Maintainer: | James Curran <[email protected]> |
| License: | GPL-2 | file LICENSE |
| Version: | 3.2.277 |
| Built: | 2026-05-10 23:49:15 UTC |
| Source: | https://github.com/stats-uoa/s20x |
Number of international airline passengers (in thousands) recorded monthly from January 1949 to December 1960.
A data frame with 144 rows and 4 variables:
Monthly total number of international airline passengers (in thousands).
Integer time index from 1 to 144.
Month of observation as a factor with levels Jan to Dec.
Year of observation as a factor with levels 1949 to 1960.
Produces analysis-of-variance-style tables for 'tslm' objects.
## S3 method for class 'tslm' anova(object, ..., verbose = FALSE)## S3 method for class 'tslm' anova(object, ..., verbose = FALSE)
object |
a fitted 'tslm' object. |
... |
optional additional fitted model objects for model comparisons. |
verbose |
logical. For AR-error models, use 'TRUE' to return the raw underlying [nlme::anova.gls()] output. |
For ordinary 'tslm()' fits without autoregressive error terms, 'anova()' returns the usual analysis of variance table from [stats::anova.lm()].
For AR-error models fitted through [nlme::gls()], the reported tests are Wald-style tests of model terms. These test whether each term contributes to the fitted mean model after allowing for the estimated autocorrelation structure. Because these models do not use the ordinary independent-error sum-of-squares decomposition, the compact table reports 'Df', 'F value', and 'Pr(>F)', but does not report 'Sum Sq' or 'Mean Sq'. Compare nested AR-error models with care: 'verbose = TRUE' exposes the underlying 'nlme' comparison output rather than recreating an ordinary 'lm' ANOVA table.
Use 'verbose = TRUE' to see the underlying [nlme::anova.gls()] output.
An analysis-of-variance-style table.
data(beer.df) fit = tslm(beer ~ t + ar(1), data = beer.df, time = t) anova(fit)data(beer.df) fit = tslm(beer ~ t + ar(1), data = beer.df, time = t) anova(fit)
These data come from a classic long-term experiment conducted at the East Malling Research Station, Kent, which is the centre for research into apple growing in the U.K. Commercial apple trees consist of two parts grafted together. The lowest part, the rootstock, largely determines the size of the tree, while the upper part (the scion) determines the fruit characteristics. Rootstocks propagated by cuttings (i.e. asexually produced) were once thought to result in smaller trees than those propagated from seeds (i.e. sexually produced). This hypothesis was re-examined in an experiment begun in 1918. Several trees of each type of 16 types of rootstock were planted, all trees having the same scion. Rootstocks I-IX were asexually produced, while X-XVI were sexually produced. In the winter of 1933-4 a number of trees were removed to make room for more, and the data presented here consists of the above-ground weights of 104 trees felled in this period. No trees of types VIII, XI or XIV were felled. The description is adapted from Lee (1994). The data are from Andrews and Herzberg (1985).
The data consist of a data frame with 104 observations on 4 variables.
Factor giving the rootstock type (I, II, III, IV, IX, V, VI, VII, X, XII, XIII, XV, XVI).
Integer Above-ground weight of tree (pounds, lb).
Numeric Above-ground weight of tree (kilograms, kg); Weight_kg = Weight * 0.45359237.
Factor giving the propagation method (cutting, seed).
Andrews, D. F. and Herzberg, A. M. (1985). Data: A Collection of Problems from Many Fields for the Student and Research Worker. New York: Springer.
Lee, A. J. (1994). Data Analysis: An Introduction Based on R. University of Auckland.
Data from an experiment to measure the effect of different images on emotional arousal, by measuring changes in pupil diameter. The experiment used 20 males and 20 females. Images included a nude man, nude woman, infant, and a landscape.
A data frame with 160 observations on 3 variables.
Numeric Change in the subject's pupil size.
Factor Subject's gender (female, male)
Factor Picture shown to subject (infant, landscape, nude female, nude male)
Plots current vs lagged residuals along with quadrants dividing these residuals about the value zero.
autocorPlot(fit, main = "Current vs Lagged residuals", ...)autocorPlot(fit, main = "Current vs Lagged residuals", ...)
fit |
output from the function 'lm()'. |
main |
the plot title. |
... |
extra parameters to be passed to the |
Plots current vs lagged residuals along with quadrants dividing these residuals about the value zero.
autocor.plot is deprecated and no longer exported.
Use autocorPlot() in new code.
data(airpass.df) time = 1:144 airpass.fit = lm(passengers ~ time, data = airpass.df) autocorPlot(airpass.fit)data(airpass.df) time = 1:144 airpass.fit = lm(passengers ~ time, data = airpass.df) autocorPlot(airpass.fit)
Monthly United States beer production figures (in millions of 31-gallon barrels) for the period July 1970 to June 1978.
A data frame with 96 rows and 4 variables:
Monthly beer production, expressed in megalitres (converted from millions of 31-US-gallon barrels; 1 million 31-gallon barrels is approximately equal to 117.35 megalitres).
Integer time index from 1 to 96.
Month of observation as a factor with levels
Jul, Aug, Sep, Oct, Nov, Dec,
Jan, Feb, Mar, Apr, May, Jun.
Year of observation as a factor with levels
1970 to 1978.
The original primary source for this monthly beer-production series is not identified in the available package materials.
This dataset originates from a study conducted at the University of Auckland in the early 1990s by Dr. R.A. Marshall and colleagues from the Department of Psychology. The research explored how cultural background and ethnic identity influence body image perceptions within the specific context of Aotearoa New Zealand.
A data frame with 246 observations on 8 variables.
Factor Subject's ethnicity (Asian, Europn, Maori, Pacific)
Factor Whether the subject is married (no, yes)
Factor Subject's rating of themself (slight.uw, right, slight.ow, mod.ow, very.ow)
Factor Whether the subject has ever smoked (no, yes)
Numeric Weight in kg.
Numeric Height in cm.
Numeric Age in years.
Factor Stress level group (low, medium, high)
The study specifically focused on a cohort of women who were generally "thin" (slightly underweight for their body size). This was designed to investigate whether body dissatisfaction and varying self-perceptions persisted even among individuals who already met or approached Western "thin" ideals, and how these perceptions differed across Asian, European, Māori, and Pacific ethnic groups.
Marshall, R.A., Department of Psychology, University of Auckland.
Lee, A. J. (1994). Data Analysis: An Introduction Based on R. University of Auckland.
This data consists of 50 sentence lengths from each of 8 books. The books “Disclosure” and “Rising Sun” were written by Michael Crichton, whilst the others “Four Past Midnight”, “The Dark Half”, “ Eye of the Dragon”, “The Shining”, “The Stand” and “The Tommy-Knockers” were written by Stephen King. The pages and sentences were chosen using a multistage design where the pages were selected at random, and then sentences within each page were selected at random. These data were collected by James Curran.
The data frame consists of 400 observations on 2 variables.
Integer sentence length, measured as the number of words in the sentence.
Factor giving the book from which the sentence was sampled (4.Past.Mid, Dark.Half, Disclosure, Eye.Drag, Rising.Sun, Shining, Stand, T.Knock).
'boxqq()' is deprecated and is no longer exported. It draws boxplots and normal quantile-quantile plots of x for each value of the grouping variable g
boxqq(formula, ...) ## S3 method for class 'formula' boxqq(formula, data = NULL, ...)boxqq(formula, ...) ## S3 method for class 'formula' boxqq(formula, data = NULL, ...)
formula |
A symbolic specification of the form |
... |
Arguments to be passed to methods, such as graphical parameters
(see |
data |
An optional data frame in which to evaluate the formula. |
Returns the plot.
boxqq(formula): Box plots and normal quantile-quantile plots
This is a legacy teaching helper retained for compatibility with older course material. New teaching material should prefer current diagnostic plotting workflows.
## Zoo data data(zoo.df) s20x:::boxqq(attendance ~ day.type, data = zoo.df)## Zoo data data(zoo.df) s20x:::boxqq(attendance ~ day.type, data = zoo.df)
Data for the 2001 Bursary results for 75 secondary schools in the Auckland area. For each school the decile rating of the school is recorded along with the percentage of eligible students who gain a B Bursary or better.
A data frame with 75 observations on 2 variables.
Numeric Decile rating of the school.
Numeric percentage of eligible students who gained a B Bursary or better.
This data gives the mean percentage of butterfat produced by different Canadian pure-bred diary cattle. There are five different breeds and two age groups, two years old and greater than five years old. For each combination of breed and age, there are measurements for 10 cows.
A data frame with 100 observations on 3 variables.
Numeric mean percentage of butterfat per cow.
Factor giving the cattle breed (ayrshire, canadian, guernesy, holst.fres, jersey).
Factor giving the age group (2yo, mature).
A Handbook of Small Data Sets
Hand, D.J., Daly, F., Lunn, A.D., McConway, K.J. and Ostrowski, E. (1994). A Handbook of Small Data Sets. Boca Raton, Florida: Chapman and Hall/CRC.
Sokal, R.R. and Rohlf, F.J. (1981). Biometry, 2nd edition. San Francisco: W.H. Freeman, 368.
66 bluegills were captured from Camp Lake, Minnesota. For each bluegill we have the length of the fish, its age in years and its age in scale radius.
A data frame with 66 observations on 3 variables.
Numeric age of the fish, in years.
Numeric radius of the key scale, in hundredths of a millimetre.
Numeric length at capture, in millimetres.
Renders a specified case study R Markdown file shipped with the package to HTML and optionally opens it in a web browser.
casestudy( id, output_dir = tempfile("s20x_case_study_"), open = interactive(), quiet = TRUE, ... ) cs(...)casestudy( id, output_dir = tempfile("s20x_case_study_"), open = interactive(), quiet = TRUE, ... ) cs(...)
id |
A case study identifier. Flexible formats are accepted, including
|
output_dir |
Directory where the rendered HTML file should be written.
Defaults to a temporary directory. This legacy argument is retained for
compatibility; new code may use the camelCase |
open |
Logical; if |
quiet |
Logical; passed to |
... |
Additional arguments passed to |
Case studies are expected to live in inst/case_studies and to be named
using the pattern CS<chapter>_<number>.Rmd (for example,
CS9_2.Rmd).
The case study is rendered on demand using rmarkdown::render().
Figures and other outputs are generated at render time; users therefore need
any required packages installed for the selected case study.
The rendered HTML file is returned invisibly.
Invisibly returns the path to the rendered HTML file.
## Not run: casestudy("CS9_2") casestudy("9.2") casestudy("9_2", outputDir = tempdir()) cs("9_2") ## End(Not run)## Not run: casestudy("CS9_2") casestudy("9.2") casestudy("9_2", outputDir = tempdir()) cs("9_2") ## End(Not run)
These data involve 11 laboratories and 2 brands of chalk. The laboratories tested the density of the chalk. The main interest was whether the different laboratories yielded the same density for the two different types of chalk.
A data frame with 66 observations on 3 variables.
Numeric density of the chalk.
Integer laboratory identifier.
Factor giving the chalk brand tested (A, B).
Calculates and prints the confidence intervals for the fitted model.
ciReg(fit, conf.level = 0.95, print.out = TRUE)ciReg(fit, conf.level = 0.95, print.out = TRUE)
fit |
an object of class |
conf.level |
confidence level of the intervals. |
print.out |
if |
The function returns a two-column matrix containing the upper and lower endpoints of the intervals.
##Peruvian Indians data data(peru.df) fit = lm(BP ~ age + years + weight + height, data = peru.df) ciReg(fit)##Peruvian Indians data data(peru.df) fit = lm(BP ~ age + years + weight + height, data = peru.df) ciReg(fit)
Data from a test to see if a questionnaire was properly designed. The questionnaire measures managers' technical knowledge of computers. The test has 19 managers complete the questionnaire as well as rate their own technical expertise.
A data frame with 19 observations on 2 variables.
Numeric questionnaire score.
Ordered factor giving the self-assessed level of expertise (1 = low, 2 = medium, 3 = high).
Draws a Cook's distance plot.
cooks20x( x, main = "Cook's Distance plot", xlab = "observation number", ylab = "Cook's distance", line = c(0.5, 1.2, 2), cex.labels = 1, axisOpts = list(xAxis = TRUE, yAxisTight = FALSE), ... )cooks20x( x, main = "Cook's Distance plot", xlab = "observation number", ylab = "Cook's distance", line = c(0.5, 1.2, 2), cex.labels = 1, axisOpts = list(xAxis = TRUE, yAxisTight = FALSE), ... )
x |
an object of class |
main |
the plot title |
xlab |
the x-axis title. |
ylab |
the y-axis title. |
line |
a vector of length 3 controlling the distances of the plot title, the x-axis title and the y-axis title from the axis in line units. |
cex.labels |
a factor controlling the font size of the labels on suspected high influence points. |
axisOpts |
a list of additional arguments that can be used to control the axes. At this point
this list only contains one element |
... |
additional arguments are passed to |
Returns the plot and identifies the three highest Cook's values
# Peruvian Indians data data(peru.df) peru.fit = lm(BP ~ age + years + I(years^2) + weight + height, data = peru.df) cooks20x(peru.fit)# Peruvian Indians data data(peru.df) peru.fit = lm(BP ~ age + years + I(years^2) + weight + height, data = peru.df) cooks20x(peru.fit)
Data from a summer school Stats 20x course. Each observation represents a single student.
A data frame with 146 observations on 15 variables.
Factor Final grade for the course (A, B, C, D)
Factor Passed the course (No, Yes)
Numeric Mark in the final exam.
Factor Degree enrolled in (BA, BCom, BSc, Other)
Factor Gender (Female, Male)
Factor Regularly attended class (No, Yes)
Numeric Assignment mark.
Numeric Test mark.
Numeric Mark for the short answer section of the exam.
Numeric Mark for the long answer section of the exam.
Numeric Mark for the multiple choice section of the exam.
Factor Colour of the exam booklet (Blue, Green, Pink, Yellow)
Factor Stage one grade (A, B, C)
Numeric Number of years since doing Stage 1.
Factor Repeating the paper (No, Yes)
Data from a summer school Stats 20x course. Each observation represents a single student. It is of interest to see if there is a relationship between a student's final examination mark and both their gender and whether they regularly attend lectures.
A data frame with 40 observations on 3 variables.
Numeric Final exam mark (out of 100)
Factor Gender (Female, Male)
Factor Regularly attended or not (No, Yes)
Computes a factor that has a level for each combination of the factors 'fac1' and 'fac2'.
crossFactors(x, fac2 = NULL, ...) ## Default S3 method: crossFactors(x, fac2 = NULL, ...) ## S3 method for class 'formula' crossFactors(formula, fac2 = NULL, data = NULL, ...)crossFactors(x, fac2 = NULL, ...) ## Default S3 method: crossFactors(x, fac2 = NULL, ...) ## S3 method for class 'formula' crossFactors(formula, fac2 = NULL, data = NULL, ...)
x |
the name of the first factor or a formula in the form |
fac2 |
the name of the second factor - ignored if |
... |
Optional arguments |
formula |
a formula in the form |
data |
an optional data frame in which to evaluate the formula |
Returns a vector containing the factor which represents the interaction of the given factors.
crossFactors(default): Crossed Factors
crossFactors(formula): Crossed Factors
This function actually returns a factor now instead of a character string, so coercion into a factor is no longer necessary.
## arousal data: data(arousal.df) gender.picture = crossFactors(arousal.df$gender, arousal.df$picture) gender.picture ## arousal data: data(arousal.df) gender.picture = crossFactors(~ gender * picture, data = arousal.df) gender.picture## arousal data: data(arousal.df) gender.picture = crossFactors(arousal.df$gender, arousal.df$picture) gender.picture ## arousal data: data(arousal.df) gender.picture = crossFactors(~ gender * picture, data = arousal.df) gender.picture
Produces a 2-way table of counts and the corresponding chi-square test of independence or homogeneity.
crosstabs(formula, data)crosstabs(formula, data)
formula |
a symbolic description of the model to be fit: ~ fac1 + fac2; where fac1 and fac2 are vectors to be crosstabulated and treated internally as factors. |
data |
an optional data frame containing the variables in the model. |
Invisibly returns an object of class ct.20x, which is a
list containing the following components:
row.props |
a matrix of row proportions, i.e. cell counts divided by row marginals. |
col.props |
a matrix of column proportions, i.e. cell counts divided by column marginals. |
whole.props |
a matrix of whole-table proportions. |
Totals |
a matrix containing the cell counts and the marginal totals. |
exp |
a matrix of expected counts from the chi-square calculation. |
chi |
a matrix of cell contributions to the chi-square statistic. |
This is a legacy teaching helper retained for compatibility with
older course material. New code should usually prefer table() and
chisq.test() directly, or a purpose-built teaching wrapper.
##body image data: data(body.df) crosstabs(~ ethnicity + married, body.df)##body image data: data(body.df) crosstabs(~ ethnicity + married, body.df)
Prices of ladies' diamond rings from a Singaporean retailer and the weight of their diamond stones.
A data frame with 48 observations on 2 variables.
Numeric Price of ring (Singapore dollars)
Numeric Weight of Diamond (carats)
Displays within-level pairwise comparisons from a two-way ANOVA with
interactions. Note that this is just a display function: it ignores any
cross-level pairs included in allpairs, even though these will have
contributed to the computations for the Tukey adjustments. The purpose is
just to organise the output from emmeans into a more convenient
format.
displayPairs(allpairs, levels1, levels2, brief = TRUE, asDF = FALSE)displayPairs(allpairs, levels1, levels2, brief = TRUE, asDF = FALSE)
allpairs |
pairwise output from a command like |
levels1 |
a character string specifying which within-level comparisons
from |
levels2 |
a character string specifying which within-level comparisons
from |
brief |
either |
asDF |
either |
allpairs is a pairwise output from a command like
pairs(emmeans(fit, ~factor1 * factor2)). If allpairs is not
already a data.frame it will be converted to a data.frame
within this function. It must contain a column called contrast with
text descriptions like 'lev1 lev2 - lev3 lev4' etc. levels1
and levels2 are character strings specifying which within-level
comparisons are wanted, and in which order. They must match the order
specified in emmeans, so if using emmeans(fit, ~factor1 *
factor2) then levels1 must belong to factor1 and
levels2 must belong to factor2. All this function does is to
pick out the rows of allpairs with the requested contrasts, so if
there are no contrasts of the requested format (e.g. because levels1
and levels2 have been switched) it will output a blank list. If
brief = TRUE, columns labelled df, SE, and
t.ratio or z.ratio will be removed for a more succinct
display. If asDF = TRUE, the output is returned as a data-frame
suitable for further manipulation, whereas if asDF = FALSE it is
returned as a list for display only.
Rachel Fewster
## Fit a two-way ANOVA to the arousal data in arousal.df. ## The factors are gender (female, male) and picture shown to ## subject (infant, landscape, nude.f, nude.m): data(arousal.df) arousal.fit = lm(arousal ~ gender * picture, data = arousal.df) ## Create all pairwise comparisons using emmeans, if available. if (requireNamespace("emmeans", quietly = TRUE)) { emmeansFun = getExportedValue("emmeans", "emmeans") arousal.allpairs = pairs( emmeansFun(arousal.fit, ~ gender * picture), infer = TRUE ) ## Display only the within-level comparisons: displayPairs( arousal.allpairs, levels1 = c("female", "male"), levels2 = c("infant", "landscape", "nude.f", "nude.m") ) }## Fit a two-way ANOVA to the arousal data in arousal.df. ## The factors are gender (female, male) and picture shown to ## subject (infant, landscape, nude.f, nude.m): data(arousal.df) arousal.fit = lm(arousal ~ gender * picture, data = arousal.df) ## Create all pairwise comparisons using emmeans, if available. if (requireNamespace("emmeans", quietly = TRUE)) { emmeansFun = getExportedValue("emmeans", "emmeans") arousal.allpairs = pairs( emmeansFun(arousal.fit, ~ gender * picture), infer = TRUE ) ## Display only the within-level comparisons: displayPairs( arousal.allpairs, levels1 = c("female", "male"), levels2 = c("infant", "landscape", "nude.f", "nude.m") ) }
Plots the residuals versus the fitted (or predicted) values from a linear
model. A horizontal line is drawn at y = 0, reflecting the fact that we
expect the residuals to have a mean of zero. An optional lowess line is
drawn if smoother is set to TRUE. This can be useful in determining whether
a trend still exists in the residuals. An optional pair of lines is drawn at
+/- 2 times the standard deviation of the residuals - which is estimated
from the Residual Mean Sqare (Within group mean square = WGMS). This can be
useful in highlighting potential outliers. If the model has one or two
factors and no continous variables, i.e. if it is a oneway or twoway ANOVA
model, and levene = TRUE then the P-value from Levene's test for
equality variance is displayed in the top left hand corner,as long as the
number of observations per group exceeds two.
eovcheck(x, ...) ## S3 method for class 'formula' eovcheck( x, data = NULL, xlab = "Fitted values", ylab = "Residuals", col = NULL, smoother = FALSE, twosd = FALSE, levene = FALSE, ... ) ## S3 method for class 'lm' eovcheck(x, smoother = FALSE, twosd = FALSE, levene = FALSE, ...)eovcheck(x, ...) ## S3 method for class 'formula' eovcheck( x, data = NULL, xlab = "Fitted values", ylab = "Residuals", col = NULL, smoother = FALSE, twosd = FALSE, levene = FALSE, ... ) ## S3 method for class 'lm' eovcheck(x, smoother = FALSE, twosd = FALSE, levene = FALSE, ...)
x |
A linear model formula. Alternatively, a fitted lm object from a linear model. |
... |
Optional arguments |
data |
A data frame in which to evaluate the formula. |
xlab |
a title for the x axis: see |
ylab |
a title for the y axis: see |
col |
a colour for the lowess smoother line. |
smoother |
if TRUE then a smoothed lowess line will be added to the plot |
twosd |
if |
levene |
if |
Draws the residual-versus-fitted diagnostic plot. The function is called for its plotting side effects and does not provide a stable data return object.
eovcheck(formula): Testing for equality of variance plot
eovcheck(lm): Testing for equality of variance plot
# one way ANOVA - oysters data(oysters.df) oyster.fit = lm(Oysters ~ Site, data = oysters.df) eovcheck(oyster.fit) # Same model as the previous example, but using eovcheck.formula data(oysters.df) eovcheck(Oysters ~ Site, data = oysters.df) # A two-way model without interaction data(soyabean.df) soya.fit=lm(yield ~ planttime + cultivar, data = soyabean.df) eovcheck(soya.fit) # A two-way model with interaction data(arousal.df) arousal.fit = lm(arousal ~ gender * picture, data = arousal.df) eovcheck(arousal.fit) # A regression model data(peru.df) peru.fit = lm(BP ~ height + weight + age + years, data = peru.df) eovcheck(peru.fit) # A time series model data(airpass.df) t = 1:144 month = factor(rep(1:12, 12)) airpass.df = data.frame(passengers = airpass.df$passengers, t = t, month = month) airpass.fit = lm(log(passengers)[-1] ~ t[-1] + month[-1] + log(passengers)[-144], data = airpass.df) eovcheck(airpass.fit)# one way ANOVA - oysters data(oysters.df) oyster.fit = lm(Oysters ~ Site, data = oysters.df) eovcheck(oyster.fit) # Same model as the previous example, but using eovcheck.formula data(oysters.df) eovcheck(Oysters ~ Site, data = oysters.df) # A two-way model without interaction data(soyabean.df) soya.fit=lm(yield ~ planttime + cultivar, data = soyabean.df) eovcheck(soya.fit) # A two-way model with interaction data(arousal.df) arousal.fit = lm(arousal ~ gender * picture, data = arousal.df) eovcheck(arousal.fit) # A regression model data(peru.df) peru.fit = lm(BP ~ height + weight + age + years, data = peru.df) eovcheck(peru.fit) # A time series model data(airpass.df) t = 1:144 month = factor(rep(1:12, 12)) airpass.df = data.frame(passengers = airpass.df$passengers, t = t, month = month) airpass.fit = lm(log(passengers)[-1] ~ t[-1] + month[-1] + log(passengers)[-144], data = airpass.df) eovcheck(airpass.fit)
Calculates and prints Tukey multiple confidence intervals for contrasts in one or two-way ANOVA.
estimateContrasts( contrast.matrix, fit, row = TRUE, alpha = 0.05, L = NULL, FUN = identity )estimateContrasts( contrast.matrix, fit, row = TRUE, alpha = 0.05, L = NULL, FUN = identity )
contrast.matrix |
A matrix of contrast coefficients. Separate rows of the matrix contain the contrast coefficients for that particular contrast, and a column for each level of the factor. |
fit |
Output from the [lm()] function. |
row |
If 'TRUE', and the ANOVA is two-way, then contrasts in the row effects are printed, otherwise contrasts in the column effects are printed. Ignored if the ANOVA is one-way. |
alpha |
The nominal error rate for the multiple confidence intervals. |
L |
Number of contrasts. If 'NULL', 'L' will be set to the number of rows in the contrast matrix, otherwise 'L' will be as specified. |
FUN |
Optional function to be applied to estimates and confidence intervals. Typically used for back-transformation operations. |
Returns a matrix whose rows correspond to the different contrasts being estimated and whose columns correspond to the point estimate of the contrast, the Tukey lower and upper limits of the confidence interval, the unadjusted p-value, and the Tukey and Bonferroni p-values.
[summary1way()], [summary2way()], [multipleComp()]
## computer data: data(computer.df) computer.df = within(computer.df, {selfassess = factor(selfassess)}) computer.fit = lm(score ~ selfassess, data = computer.df) contrast.matrix = matrix(c(-1 / 2, -1 / 2, 1), byrow = TRUE, nrow = 1, ncol = 3) contrast.matrix estimateContrasts(contrast.matrix, computer.fit)## computer data: data(computer.df) computer.df = within(computer.df, {selfassess = factor(selfassess)}) computer.fit = lm(score ~ selfassess, data = computer.df) contrast.matrix = matrix(c(-1 / 2, -1 / 2, 1), byrow = TRUE, nrow = 1, ncol = 3) contrast.matrix estimateContrasts(contrast.matrix, computer.fit)
House damage and distance from the fire station, of 15 house fires. Data collected by an insurance company for homes in a particular area.
A data frame with 15 observations on 3 variables.
Numeric Damage (1000s of dollars)
Numeric Distance from the fire station (miles)
Numeric Distance from the fire station (kilometres); distance_km = distance * 1.60934.
If hypothprob is absent: prints confidence intervals for the true proportions, a Chi-square test for uniformity, confidence intervals for differences in proportions (with no corrections for multiple comparisons), and plots the proportions.
freq1way( counts, hypothprob, conf.level = 0.95, addCIs = TRUE, digits = 4, arrowwid = 0.1, estimated = 0 )freq1way( counts, hypothprob, conf.level = 0.95, addCIs = TRUE, digits = 4, arrowwid = 0.1, estimated = 0 )
counts |
A 1-way frequency table as produced by |
hypothprob |
If present, a set of probabilities to test the cell counts against. |
conf.level |
confidence level for the confidence interval, expressed as a decimal. |
addCIs |
If true, adds confidence limits to plot of sample proportions. |
digits |
used to control rounding of printout. |
arrowwid |
controls width of arrowheads. |
estimated |
default is |
If hypothprob is present: prints confidence intervals for the true proportions, a Chi-square test for the hypothesised probabilities, and plots the sample proportions (with attached confidence limits) alongside the corresponding hypothesised probabilities.
An invisible list containing the following components:
CIs |
a matrix containing the confidence intervals. |
exp |
a vector of the expected counts. |
chi |
a vector of the components of Chi-square. |
These confidence intervals have been Bonferroni adjusted for multiple comparisons. This is a legacy teaching helper retained for compatibility with older course material.
##Body image data: data(body.df) eth.table = with(body.df, table(ethnicity)) freq1way(eth.table) freq1way(eth.table,hypothprob=c(0.2,0.4,0.3,0.1))##Body image data: data(body.df) eth.table = with(body.df, table(ethnicity)) freq1way(eth.table) freq1way(eth.table,hypothprob=c(0.2,0.4,0.3,0.1))
This data gives fecundity for female fruitflies, Drosophila melanogaster. The fecundity is the number of eggs laid, per day, for the fruitfly's first 14 days of life. There are three strains: A control group, NS, Nonselected Strain, as well as RS, a strain bred for resistance to DDT and SS, a strain bred for susceptibility to DDT. Each strain contains 25 measurements. It is of interest to compare the level of fecundity across strains.
A data frame with 75 observations on 2 variables.
Numeric Number of eggs laid, per day, per fruitfly.
Factor Strain of fruitfly (NS, RS, SS)
A Handbook of Small Data Sets
Hand, D.J., Daly, F., Lunn, A.D., McConway, K.J. and Ostrowski, E. (1994). A Handbook of Small Data Sets. Boca Raton, Florida: Chapman and Hall/CRC.
Sokal, R.R. and Rohlf, F.J. (1981). Biometry, 2nd edition. San Francisco: W.H. Freeman, 239.
Returns the version number of the s20x package. This is useful if a student has problems running commands and the maintainer needs to check the version number.
getVersion()getVersion()
getVersion()getVersion()
A random sample of 100 houses recently sold in Mt Eden, Auckland. For each house we have the advertised price and the actual sale price.
A data frame with 100 observations on 2 variables.
Numeric Advertised price (dollars)
Numeric Final sale price (dollars)
Random sample of 152 families giving their mean income (1000s of dollars). The sample was taken by an advertising agency over their area of operations.
A data frame with 152 observations on 1 variable.
Numeric mean family income, in thousands of dollars.
Displays data with intervals for each combination of the two factors and shows the mean differences between levels of the first factor for each level of the second factor. Note that there should be more than one observation for each combination of factors.
interactionPlots(y, ...) ## Default S3 method: interactionPlots( y, fac1 = NULL, fac2 = NULL, xlab = NULL, xlab2 = NULL, ylab = NULL, data.order = TRUE, exlim = 0.1, jitter = 0.02, conf.level = 0.95, interval.type = c("tukey", "hsd", "lsd", "ci"), pooled = TRUE, tick.length = 0.1, interval.distance = 0.2, col.width = 2/3, xlab.distance = 0.1, xlen = 1.5, ylen = 1, ... ) ## S3 method for class 'formula' interactionPlots( y, data = NULL, xlab = NULL, xlab2 = NULL, ylab = NULL, data.order = TRUE, exlim = 0.1, jitter = 0.02, conf.level = 0.95, interval.type = c("tukey", "hsd", "lsd", "ci"), pooled = TRUE, tick.length = 0.1, interval.distance = 0.2, col.width = 2/3, xlab.distance = 0.1, xlen = 1.5, ylen = 1, ... )interactionPlots(y, ...) ## Default S3 method: interactionPlots( y, fac1 = NULL, fac2 = NULL, xlab = NULL, xlab2 = NULL, ylab = NULL, data.order = TRUE, exlim = 0.1, jitter = 0.02, conf.level = 0.95, interval.type = c("tukey", "hsd", "lsd", "ci"), pooled = TRUE, tick.length = 0.1, interval.distance = 0.2, col.width = 2/3, xlab.distance = 0.1, xlen = 1.5, ylen = 1, ... ) ## S3 method for class 'formula' interactionPlots( y, data = NULL, xlab = NULL, xlab2 = NULL, ylab = NULL, data.order = TRUE, exlim = 0.1, jitter = 0.02, conf.level = 0.95, interval.type = c("tukey", "hsd", "lsd", "ci"), pooled = TRUE, tick.length = 0.1, interval.distance = 0.2, col.width = 2/3, xlab.distance = 0.1, xlen = 1.5, ylen = 1, ... )
y |
either a formula of the form: y~fac1+fac2 where y is the response and fac1 and fac2 are the two explanatory variables used as factors, or a single response vector |
... |
optional arguments. |
fac1 |
if 'y' is a vector, then fac1 contains the levels of factor 1 which correspond to the y value |
fac2 |
if 'y' is a vector, then fac2 contains the levels of factor 2 which correspond to the y value |
xlab |
an optional label for the x-axis. If not specified the name of fac1 will be used. |
xlab2 |
an optional label for the lines. If not specified the name of fac2 will be used. |
ylab |
An optional label for the y-axis. If not specified the name of y will be used. |
data.order |
if TRUE the levels of fac1 and fac2 will be set to unique(fac1) and unique(fac2) respectively. |
exlim |
provide extra limits. |
jitter |
the amount of horizontal jitter to show in the plot. The actual jitter is determined as the function is called, and will likely be different each time the function is used. |
conf.level |
confidence level of the intervals. |
interval.type |
four options for intervals appearing on plot: 'tukey', 'hsd', 'lsd' or 'ci'. |
pooled |
two options: pooled or unpooled standard deviation used for plotted intervals. |
tick.length |
size of tick, in inches. |
interval.distance |
distance, as a fraction of the column width, between the points and interval. This is in addition to the extra space allocated for the jitter. |
col.width |
width of a factor ‘column’, as a fraction of the space between the centres of two columns. |
xlab.distance |
distance of x-axis labels from bottom of plot, as a fraction of the overall height of the plot. |
xlen, ylen
|
character interspacing factor for horizontal (x) and vertical (y) spacing of the legend. |
data |
an optional data frame containing the variables in the model. |
interactionPlots(default): Interactions Plot for Two-way Analysis of Variance
interactionPlots(formula): Interactions Plot for Two-way Analysis of Variance
data(arousal.df) interactionPlots(arousal ~ gender + picture, data = arousal.df) ## This usage is deprecated. with(arousal.df, interactionPlots(arousal, gender, picture))data(arousal.df) interactionPlots(arousal ~ gender + picture, data = arousal.df) ## This usage is deprecated. with(arousal.df, interactionPlots(arousal, gender, picture))
The ages and lengths of 78 bluegills captured from Lake Mary, Minnesota.
A data frame with 78 observations on 2 variables.
Numeric Age of the fish (years)
Numeric Length at capture (mm)
Annual rainfall (in inches) for Los Angeles from 1908 to 1973.
A data frame with 66 rows and 4 variables:
Annual rainfall in Los Angeles, measured in inches.
Annual rainfall in Los Angeles, measured in millimetres (mm); rain_mm = LA.Rain * 25.4.
Integer time index from 1 to 66.
Year of observation as an integer from 1908 to 1973.
Allows a 'numRows' by 'numCols' matrix of plots to be displayed in a single plot. If the function is called with no arguments, then the plotting device layout will be reset to a single plot.
layout20x(numRows = 1, numCols = 1)layout20x(numRows = 1, numCols = 1)
numRows |
Number of rows in the plot array. |
numCols |
Number of columns in the plot array. |
No return value.
This is a legacy convenience wrapper retained for compatibility with
older teaching material. New code can use par(mfrow = ...)
directly.
data(course.df) layout20x(1, 2) stripchart(course.df$Exam) boxplot(course.df$Exam)data(course.df) layout20x(1, 2) stripchart(course.df$Exam) boxplot(course.df$Exam)
Perform a Levene test for equal group variances in both one-way and two-way ANOVA. A table with the results is (normally) displayed.
levene.test(formula, data, digit = 5, show.table = TRUE)levene.test(formula, data, digit = 5, show.table = TRUE)
formula |
a symbolic description of the model to be fitted: response ~ fac1 + fac2. |
data |
an optional data frame containing the variables in the model. |
digit |
the number of decimal places to display. |
show.table |
If this argument is FALSE then the output will be suppressed |
A list with the following elements:
df |
degrees of freedom. |
ss |
sum squares. |
ms |
mean squares. |
f.value |
F-statistic value. |
p.value |
P-value. |
## data(computer.df) levene.test(score ~ factor(selfassess), computer.df)## data(computer.df) levene.test(score ~ factor(selfassess), computer.df)
Lists all case study R Markdown files shipped with the package and prints them as a formatted text table.
listCaseStudies() listCS() lcs()listCaseStudies() listCS() lcs()
Case studies are expected to live in inst/case_studies and to be named
using the pattern CS<chapter>_<number>.Rmd (e.g. CS9_2.Rmd).
The table has two columns: File (the case study identifier) and
Title (extracted from the YAML header). Case studies are listed in
numerical order, not alphabetical order.
The function invisibly returns a character vector of case study identifiers.
Invisibly returns a character vector of case study identifiers.
## Not run: listCaseStudies() ids = listCaseStudies() ## End(Not run)## Not run: listCaseStudies() ids = listCaseStudies() ## End(Not run)
Prices and ages of 124 Mazda cars collected from the Melbourne Age newspaper in 1991.
A data frame with 124 observations on 2 variables.
Numeric Price (Australian dollars)
Numeric Year of manufacture.
This data shows the monthly number of notifications meningococcal disease in New Zealand from January 1990 to December 2001.
A data frame with 144 observations on 3 variables.
Factor giving the month of notification.
Factor giving the year of notification.
Numeric number of notifications of meningococcal disease.
A random selection of 38 consummated mergers from the USA, 1982, giving the number of days between the date the merger was announced and the date the merger became effective.
A data frame with 38 observations on 1 variable.
Numeric number of days between the merger announcement and the effective date.
'modcheck()' is deprecated and is no longer exported. It plots four model checking plots: residuals versus fitted values, a normal Q-Q plot, a histogram of residuals with a normal distribution superimposed, and a Cook's distance plot.
modcheck(x, ...) ## S3 method for class 'lm' modcheck( x, plotOrder = 1:4, args = list(eovcheck = list(smoother = FALSE, twosd = FALSE, levene = FALSE, ...), normcheck = list(xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, ...), cooks20x = list(main = "Cook's Distance plot", xlab = "observation number", ylab = "Cook's distance", line = c(0.5, 0.1, 2), cex.labels = 1, axisOpts = list(xAxis = TRUE), ...)), parVals = list(mfrow = c(2, 2), xaxs = "r", yaxs = "r", pty = "s", mai = c(0.2, 0.2, 0.05, 0.05)), ... )modcheck(x, ...) ## S3 method for class 'lm' modcheck( x, plotOrder = 1:4, args = list(eovcheck = list(smoother = FALSE, twosd = FALSE, levene = FALSE, ...), normcheck = list(xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, ...), cooks20x = list(main = "Cook's Distance plot", xlab = "observation number", ylab = "Cook's distance", line = c(0.5, 0.1, 2), cex.labels = 1, axisOpts = list(xAxis = TRUE), ...)), parVals = list(mfrow = c(2, 2), xaxs = "r", yaxs = "r", pty = "s", mai = c(0.2, 0.2, 0.05, 0.05)), ... )
x |
a vector of observations, or the residuals from fitting a linear model. Alternatively, a fitted |
plotOrder |
the order of the plots: 1 for residuals versus fitted values, 2 for the normal Q-Q plot, 3 for the histogram, and 4 for Cook's distance plot. |
args |
a list containing three additional lists |
parVals |
the values that are set via |
... |
additional parameters. Included for future flexibility, but unsure how this might be used currently. |
Draws the selected model checking plots for teaching diagnostics. The function is called for its plotting side effects and does not provide a stable data return object.
modcheck(lm): Model checking plots
# Synthetic teaching example: an exponential growth curve set.seed(123) e = rnorm(100, 0, 0.1) x = rnorm(100) y = exp(5 + 3 * x + e) fit = lm(y ~ x, data = data.frame(x, y)) s20x:::modcheck(fit) # An exponential growth curve with the correct transformation fit = lm(log(y) ~ x, data = data.frame(x, y)) s20x:::modcheck(fit) # Peruvian Indians data data(peru.df) s20x:::modcheck(lm(BP ~ weight, data = peru.df))# Synthetic teaching example: an exponential growth curve set.seed(123) e = rnorm(100, 0, 0.1) x = rnorm(100) y = exp(5 + 3 * x + e) fit = lm(y ~ x, data = data.frame(x, y)) s20x:::modcheck(fit) # An exponential growth curve with the correct transformation fit = lm(log(y) ~ x, data = data.frame(x, y)) s20x:::modcheck(fit) # Peruvian Indians data data(peru.df) s20x:::modcheck(lm(BP ~ weight, data = peru.df))
Draw the teaching diagnostic plots used by older 's20x' workflows. 'modelcheck()' is retained as an exported compatibility helper for model checking, while newer teaching material may use focused diagnostic helpers such as [eovcheck()], [normcheck()], and [cooks20x()] directly.
modelcheck(x, ...) ## S3 method for class 'lm' modelcheck(x, which = 1:3, mar = c(3, 4, 1.5, 4), ...)modelcheck(x, ...) ## S3 method for class 'lm' modelcheck(x, which = 1:3, mar = c(3, 4, 1.5, 4), ...)
x |
The fitted model. |
which |
The plot(s) to be drawn. Residuals versus fitted values
( |
mar |
Margins applied to each selected plot. |
... |
any other arguments to pass to |
Draws diagnostic plots for teaching model checking. The function is called for its plotting side effects and does not provide a stable data return object.
modelcheck(lm): Model checking plots for linear models.
data(peru.df) lmFit = lm(BP ~ weight, data = peru.df) # Plot residuals versus fitted values only modelcheck(lmFit, 1) # Plot residuals versus fitted values, histogram, and Q-Q plot modelcheck(lmFit, 1:2) # Plot all diagnostics modelcheck(lmFit)data(peru.df) lmFit = lm(BP ~ weight, data = peru.df) # Plot residuals versus fitted values only modelcheck(lmFit, 1) # Plot residuals versus fitted values, histogram, and Q-Q plot modelcheck(lmFit, 1:2) # Plot all diagnostics modelcheck(lmFit)
Length of movements from 11 of Mozart's early symphonies and 11 of his late symphonies.
A data frame with 88 observations on 3 variables.
Numeric Time of each movement (seconds)
Factor Movement (M1, M2, M3, M4)
Factor Period that the symphony was written (early, late)
Calculates and prints the estimate, multiple 95% confidence intervals, unadjusted, Tukey and Bonferroni p-values for all possible differences in means in a one-way ANOVA.
multipleComp(fit, conf.level = 0.95, FUN = identity)multipleComp(fit, conf.level = 0.95, FUN = identity)
fit |
Output from the command [lm()]. |
conf.level |
Confidence level for the confidence interval, expressed as a percentage. |
FUN |
Optional function to be applied to estimates and confidence intervals. Typically used for back-transformation operations. |
Returns a list of estimates, confidence intervals and p-values.
## computer data data(computer.df) fit = lm(score ~ factor(selfassess), data = computer.df) multipleComp(fit) ## butterfat data data("butterfat.df") fit = lm(log(Butterfat) ~ Breed, data = butterfat.df) multipleComp(fit, FUN = exp)## computer data data(computer.df) fit = lm(score ~ factor(selfassess), data = computer.df) multipleComp(fit) ## butterfat data data("butterfat.df") fit = lm(log(Butterfat) ~ Breed, data = butterfat.df) multipleComp(fit, FUN = exp)
These data were collected to determine whether quick drying nail polish or regular nail polish dried faster. The time for each type of nail polish to dry was recorded.
A data frame with 60 observations on 2 variables.
Factor Type of polish (Regular, Quick)
Integer Time (in seconds) for the polish to dry.
Plots two plots side by side. First, it draws a normal Q-Q plot of the
residuals, along with a line with intercept equal to the mean of the
residuals and slope equal to the standard deviation of the residuals. If
shapiro.wilk = TRUE, the P-value from the Shapiro-Wilk test for
normality is shown in the top-left corner of the Q-Q plot. Second, it draws
a histogram of the residuals. A normal distribution is fitted and
superimposed over the histogram. Note: if you want to leave the
x-axis blank in the histogram then use xlab = c("Theoretical Quantiles", " ")
, i.e. leave a space between the quotes. If you do not leave a space,
information will be extracted from x.
normcheck(x, ...) ## Default S3 method: normcheck( x, xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, ... ) ## S3 method for class 'lm' normcheck( x, xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, ... ) ## S3 method for class 'tslm' normcheck( x, xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, residualType = "normalised", ... )normcheck(x, ...) ## Default S3 method: normcheck( x, xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, ... ) ## S3 method for class 'lm' normcheck( x, xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, ... ) ## S3 method for class 'tslm' normcheck( x, xlab = c("Theoretical Quantiles", ""), ylab = c("Sample Quantiles", ""), main = c("", ""), col = "light blue", bootstrap = FALSE, B = 5, bpch = 3, bcol = "lightgrey", shapiro.wilk = FALSE, whichPlot = 1:2, usePar = TRUE, residualType = "normalised", ... )
x |
the residuals from fitting a linear model. Alternatively, a fitted |
... |
additional arguments which are passed to both |
xlab |
a title for the x-axis of both the Q-Q plot and the histogram: see |
ylab |
a title for the y-axis of both the Q-Q plot and the histogram: see |
main |
a title for both the Q-Q plot and the histogram: see |
col |
a colour for the bars of the histogram. |
bootstrap |
if |
B |
the number of bootstrap samples to take. Five should usually be sufficient. |
bpch |
the plotting symbol used for the bootstrap samples. Legal values are the same as any legal
value for |
bcol |
the plotting colour used for the bootstrap samples. Legal values are the same as any legal
value for |
shapiro.wilk |
if |
whichPlot |
legal values are |
usePar |
if |
residualType |
for |
Draws the selected normality diagnostic plots. The function is called for its plotting side effects and does not provide a stable data return object.
normcheck(default): Testing for normality plot
normcheck(lm): Testing for normality plot
normcheck(tslm): Testing for normality plot
# Synthetic teaching example: an exponential growth curve set.seed(123) e = rnorm(100, 0, 0.1) x = rnorm(100) y = exp(5 + 3 * x + e) fit = lm(y ~ x) normcheck(fit) # An exponential growth curve with the correct transformation fit = lm(log(y) ~ x) normcheck(fit) # Same example as above except we use normcheck.default normcheck(residuals(fit)) # Peruvian Indians data data(peru.df) normcheck(lm(BP ~ weight, data = peru.df))# Synthetic teaching example: an exponential growth curve set.seed(123) e = rnorm(100, 0, 0.1) x = rnorm(100) y = exp(5 + 3 * x + e) fit = lm(y ~ x) normcheck(fit) # An exponential growth curve with the correct transformation fit = lm(log(y) ~ x) normcheck(fit) # Same example as above except we use normcheck.default normcheck(residuals(fit)) # Peruvian Indians data data(peru.df) normcheck(lm(BP ~ weight, data = peru.df))
Quarterly alcohol available for consumption in New Zealand from 1935 to 2021. The data give volumes of alcoholic beverages available for consumption, grouped into broad beverage categories.
A data frame with quarterly observations on 4 variables.
Integer Year.
Ordered factor giving the month at the end of the quarter.
Numeric volume available for consumption, in million litres.
Factor beverage category: 'Total beer', 'Total wine', or 'Total spirits'.
The 'month' variable gives the month ending the quarter. It should be treated in calendar order for plotting and summaries. For this quarterly data set the intended order is March, June, September, and December.
The 'category' variable has three levels:
Total beer available for consumption.
Total wine available for consumption.
Total spirits and spirit-based drinks available for consumption.
Stats NZ, Alcohol available for consumption: Year ended December 2021.
Monthly international passenger arrivals to New Zealand from January 1921 to February 2026. Missing monthly observations, if present in the source series, are retained as rows with missing 'arrivals.count' values.
A data frame with monthly observations on 3 variables.
Integer year.
Factor month abbreviation with levels given by 'month.abb'.
Integer number of international passenger arrivals.
Stats NZ Infoshare, table ITM049AA, Total passenger movements (monthly), Arrivals, Actual Counts. Last updated 14 April 2026.
Displays stripplot/boxplot of the reponse variable with intervals by factor levels. It is used as part of a one-way ANOVA analysis.
onewayPlot(x, ...) ## Default S3 method: onewayPlot( x, f, conf.level = 0.95, interval.type = "tukey", pooled = TRUE, strip = TRUE, vert = TRUE, verbose = FALSE, ylabel = deparse(terms(formula)[[2]]), flabel = deparse(terms(formula)[[3]]), ... ) ## S3 method for class 'formula' onewayPlot( formula, data = parent.frame(), conf.level = 0.95, interval.type = "tukey", pooled = TRUE, strip = TRUE, vert = TRUE, verbose = FALSE, ylabel = deparse(terms(formula)[[2]]), flabel = deparse(terms(formula)[[3]]), ... ) ## S3 method for class 'lm' onewayPlot(x, ..., ylabel = nms[1], flabel = nms[2])onewayPlot(x, ...) ## Default S3 method: onewayPlot( x, f, conf.level = 0.95, interval.type = "tukey", pooled = TRUE, strip = TRUE, vert = TRUE, verbose = FALSE, ylabel = deparse(terms(formula)[[2]]), flabel = deparse(terms(formula)[[3]]), ... ) ## S3 method for class 'formula' onewayPlot( formula, data = parent.frame(), conf.level = 0.95, interval.type = "tukey", pooled = TRUE, strip = TRUE, vert = TRUE, verbose = FALSE, ylabel = deparse(terms(formula)[[2]]), flabel = deparse(terms(formula)[[3]]), ... ) ## S3 method for class 'lm' onewayPlot(x, ..., ylabel = nms[1], flabel = nms[2])
x |
a vector of responses, a formula object or an lm object |
... |
optional arguments. |
f |
if x is a vector of responses then f contains the group labels for each observation in x. That is, the ith value in f says which group the ith observation of x belongs to. |
conf.level |
confidence level of the intervals. |
interval.type |
three options for intervals appearing on plot: 'hsd','lsd' or 'ci'. |
pooled |
two options: pooled or unpooled standard deviation used for plotted intervals. |
strip |
if strip=F, boxplots are displayed instead. |
vert |
if vert=F, horizontal stripplots are displayed instead (boxplots can only be displayed vertically). |
verbose |
if true, print intervals on console. |
ylabel |
can be used to replace variable name of y by another string. |
flabel |
can be used to replace variable name of f by another string. |
formula |
a symbolic description of the model to be fit. |
data |
an optional data frame in which to evaluate the formula. |
onewayPlot(default): One-way Analysis of Variance Plot
onewayPlot(formula): One-way Analysis of Variance Plot
onewayPlot(lm): One-way Analysis of Variance Plot
##see example in 'summary1way' ##computer data: data(computer.df) onewayPlot(score~selfassess, data = computer.df) ##apple data: data(apples.df) twosampPlot(Weight~Propagated, data = apples.df) ##oyster data: data(oysters.df) onewayPlot(log(Oysters)~Site, data = oysters.df) ##oyster data: data(oysters.df) oyster.fit = lm(log(Oysters)~Site, data = oysters.df) onewayPlot(oyster.fit)##see example in 'summary1way' ##computer data: data(computer.df) onewayPlot(score~selfassess, data = computer.df) ##apple data: data(apples.df) twosampPlot(Weight~Propagated, data = apples.df) ##oyster data: data(oysters.df) onewayPlot(log(Oysters)~Site, data = oysters.df) ##oyster data: data(oysters.df) oyster.fit = lm(log(Oysters)~Site, data = oysters.df) onewayPlot(oyster.fit)
Opens a case study .Rmd file for interactive use. The file shipped
inside the package is copied to dest_dir (so it is writable), then
opened in the RStudio editor when available (otherwise the system editor).
openCaseStudy(id, dest_dir = getwd(), overwrite = FALSE, ...) opencs(id, dest_dir = getwd(), overwrite = FALSE, ...) ocs(id, dest_dir = getwd(), overwrite = FALSE, ...)openCaseStudy(id, dest_dir = getwd(), overwrite = FALSE, ...) opencs(id, dest_dir = getwd(), overwrite = FALSE, ...) ocs(id, dest_dir = getwd(), overwrite = FALSE, ...)
id |
Case study identifier. Flexible formats are accepted, including
|
dest_dir |
Directory to copy the case study into. Defaults to the
current working directory. This legacy argument is retained for
compatibility; new code may use the camelCase |
overwrite |
Logical; overwrite an existing file in |
... |
Additional compatibility arguments. Currently supports
|
Invisibly returns the path to the copied file.
## Not run: openCaseStudy("2.1") openCaseStudy("2.1", destDir = tempdir()) ## End(Not run)## Not run: openCaseStudy("2.1") openCaseStudy("2.1", destDir = tempdir()) ## End(Not run)
Data from an experiment to determine the abundance of oysters recruiting from three sites in two different estuaries in New South Wales. One in Georges River and two in Port Stephens. The number of oysters was recorded for 10 cm by 10 cm panels over a two year period.
A data frame with 87 observations on 2 variables.
Numeric number of oysters on each experimental panel.
Factor giving the location of the experimental panels (GR = Georges River, PS1 = first Port Stephens site, PS2 = second Port Stephens site).
Plots pairwise scatter plots with histograms and correlations for the data frame.
pairs20x(x, na.rm = TRUE, ...)pairs20x(x, na.rm = TRUE, ...)
x |
a data frame. |
na.rm |
if TRUE then only complete cases will be displayed |
... |
optional argumments which are passed to the generic pairs function. |
Returns the plots.
'pairs', 'panel.smooth', 'panel.cor', 'panel.hist'
##peruvian indians data(peru.df) pairs20x(peru.df)##peruvian indians data(peru.df) pairs20x(peru.df)
A random sample of Peruvian Indians born in the Andes mountains, but who have since migrated to lower altitudes. The sample was collected to assess the long term effects of altitude on blood pressure.
A data frame with 39 observations on 5 variables.
Numeric Subject's age.
Numeric Number of years since migration.
Numeric Subject's weight (kg)
Numeric Subject's height (mm)
Numeric Subject's systolic blood pressure (mm Hg; standard clinical unit in New Zealand).
Teaching helper for linear-model predictions. It wraps
predict.lm and prints a compact table containing fitted values,
confidence intervals for the mean response, and prediction intervals for new
observations.
predict20x(object, newdata, cilevel = 0.95, digit = 3, print.out = TRUE, ...)predict20x(object, newdata, cilevel = 0.95, digit = 3, print.out = TRUE, ...)
object |
an |
newdata |
prediction data frame. |
cilevel |
confidence level for the intervals. |
digit |
number of decimal places to print. |
print.out |
if |
... |
optional arguments that are passed to |
This is not an S3 predict() method and is not intended to be a
drop-in replacement for base R prediction methods. It is a compatibility
helper for older teaching material that expects confidence and prediction
intervals to be printed together. The standard predict
interface is preferred for new work.
Note: newdata must be a data frame with the same column order and
data types as those used in fitting the model. This is stricter than the usual
predict.lm() interface and is kept for compatibility with the original
teaching wrapper.
Invisibly returns a list with components
printed data frame containing predictions, confidence intervals, and prediction intervals.
prediction values.
standard errors of predictions.
residual standard deviation.
residual degrees of freedom.
confidence level of the interval.
This function is deprecated because it is no longer used in class.
Prefer the standard predict method for new work.
predict, predict.lm, as.data.frame.
# Zoo data data(zoo.df) zoo.df = within(zoo.df, {day.type = factor(day.type)}) zoo.fit = lm(log(attendance) ~ time + sun.yesterday + nice.day + day.type + tv.ads, data = zoo.df) pred.zoo = data.frame(time = 8, sun.yesterday = 10.8, nice.day = 0, day.type = factor(3), tv.ads = 1.181) predict20x(zoo.fit, pred.zoo) # Peruvian Indians data data(peru.df) peru.fit = lm(BP ~ age + years + I(years^2) + weight + height, data = peru.df) pred.peru = data.frame(age = 21, years = 2, `I(years^2)` = 2, weight = 71, height = 1629) predict20x(peru.fit, pred.peru)# Zoo data data(zoo.df) zoo.df = within(zoo.df, {day.type = factor(day.type)}) zoo.fit = lm(log(attendance) ~ time + sun.yesterday + nice.day + day.type + tv.ads, data = zoo.df) pred.zoo = data.frame(time = 8, sun.yesterday = 10.8, nice.day = 0, day.type = factor(3), tv.ads = 1.181) predict20x(zoo.fit, pred.zoo) # Peruvian Indians data data(peru.df) peru.fit = lm(BP ~ age + years + I(years^2) + weight + height, data = peru.df) pred.peru = data.frame(age = 21, years = 2, `I(years^2)` = 2, weight = 71, height = 1629) predict20x(peru.fit, pred.peru)
Teaching helper for count predictions from a log-link generalised linear
model. It wraps predict.glm, constructs confidence intervals on
the link scale, exponentiates the fitted values and limits, rounds the result,
and optionally prints the returned table.
predictCount(object, newdata, cilevel = 0.95, digit = 3, print.out = TRUE, ...)predictCount(object, newdata, cilevel = 0.95, digit = 3, print.out = TRUE, ...)
object |
a |
newdata |
prediction data frame. |
cilevel |
confidence level for the intervals. |
digit |
number of decimal places to print. |
print.out |
if |
... |
optional arguments that are passed to |
This is not an S3 predict() method and is not intended to be a
drop-in replacement for base R prediction methods. It is a specialised
count-focused teaching wrapper. For a more general log-link or logit-link GLM
helper, see predictGLM.
Note: newdata must be a data frame with the same column order and
data types as those used in fitting the model. This stricter interface is kept
for compatibility with the original teaching wrapper.
Invisibly returns a data frame with three columns:
the predicted count on the response scale.
the lower confidence limit on the response scale.
the upper confidence limit on the response scale.
predict, predict.glm, predictGLM, as.data.frame.
Teaching helper for predictions from log-link and logit-link generalised
linear models. It wraps predict.glm with standard errors and
returns fitted values with confidence limits on either the link scale or the
response scale.
predictGLM(object, newdata, type = "link", cilevel = 0.95, quasit = FALSE, ...)predictGLM(object, newdata, type = "link", cilevel = 0.95, quasit = FALSE, ...)
object |
a |
newdata |
prediction data frame. |
type |
|
cilevel |
confidence level for the intervals. |
quasit |
if |
... |
optional arguments that are passed to |
This is not an S3 predict() method and is not intended to be a
drop-in replacement for base R prediction methods. It is the more general GLM
teaching helper in this package; predictCount remains a
specialised count-focused wrapper with rounded response-scale output.
Note: newdata must include all first-order terms used in the fitted
model. This simplified requirement reflects the teaching-wrapper interface
and is not a complete reproduction of predict.glm().
A data frame with columns fit, lwr, and upr containing
fitted values and confidence limits on the requested scale.
predict, predict.glm, predictCount.
This function is called by rowdistr.
propslsd.new(crosstablist, conf.level = 0.95, arrowlength = 0.1)propslsd.new(crosstablist, conf.level = 0.95, arrowlength = 0.1)
crosstablist |
A list produced by |
conf.level |
Confidence level of the intervals. |
arrowlength |
Length of the arrows. |
This is an internal legacy helper used by rowdistr(). It is
not exported and should not be called directly by users.
Data from an experiment to see if seeding clouds with Silver Nitrate effects the amount of rainfall.
A data frame with 50 observations on 3 variables.
Numeric amount of rain, measured in acre-feet (the volume of water required to cover one acre of land to a depth of one foot).
Numeric amount of rain expressed in megalitres (ML); rain_ML = rain * 1.23348184.
Factor indicating whether the clouds were seeded (seeded, unseeded).
Plots a scatter plot for the variables of the residuals and fitted values from the linear model, lmfit. A lowess smooth line for the underlying trend, as well as one standard deviation error bounds for the scatter about this trend, are added to this scatter plot. A test for a quadratic relationship between the residuals and the fitted values is also computed.
residPlot(lmfit, f = 0.5)residPlot(lmfit, f = 0.5)
lmfit |
an |
f |
the smoother span. This gives the proportion of points in the plot which influence the smooth at each value. Larger values give more smoothness. |
Returns the plot.
This is a legacy diagnostic plotting helper retained for compatibility
with older teaching material. New code should usually prefer the current
diagnostic workflow used by modelcheck().
# Peruvian Indians data data(peru.df) fit = lm(BP ~ age + years + weight + height, data = peru.df) residPlot(fit)# Peruvian Indians data data(peru.df) fit = lm(BP ~ age + years + weight + height, data = peru.df) residPlot(fit)
Produces summaries and plots from a cross-tabulation. The output produced depends on the parameter 'comp'. Columns relate to response categories and rows to different populations.
rowdistr( crosstablist, comp = c("basic", "within", "between"), conf.level = 0.95, plot = TRUE, suppressText = FALSE )rowdistr( crosstablist, comp = c("basic", "within", "between"), conf.level = 0.95, plot = TRUE, suppressText = FALSE )
crosstablist |
a list produced by 'crosstabs' or a matrix containing a 2-way table of counts (without marginal totals). |
comp |
three options: 'basic' (default), 'within', and 'between'. |
conf.level |
confidence level of the intervals. |
plot |
if |
suppressText |
if |
The 'basic' option (default) produces the response distribution for each row population together with comparative bar charts.
If comp = 'between' the resulting output displays how the probability of falling into a response class (column) differs between populations. Confidence intervals for differences in proportions are produced together with a set of barcharts with LSD intervals.
If comp = 'within' the resulting output shows the extent to which the component probabilities of the same row distribution differ. Separate Chi-square tests for uniformity are produced for each row distribution as are confidence intervals for differences in proportions within the same distribution.
Arguments plot and suppressText are really only used when
producing knitr or Sweave documents so that just the plot or just the text
can be displayed in the document.
Invisibly returns the matrix of row proportions printed by the
teaching summary when suppressText = FALSE. When
suppressText = TRUE, the function invisibly returns NULL
because no text summary is constructed. Plotting remains a side effect
controlled by plot.
data(body.df) z = crosstabs(~ ethnicity + married, data = body.df) rowdistr(z) rowdistr(z, comp = "between") rowdistr(z, comp = "within") ## from matrix of counts z = matrix(c(4, 3, 2, 6, 47, 20, 40, 62, 11, 8, 7, 22, 3, 0, 1, 10), 4, 4) rowdistr(z)data(body.df) z = crosstabs(~ ethnicity + married, data = body.df) rowdistr(z) rowdistr(z, comp = "between") rowdistr(z, comp = "within") ## from matrix of counts z = matrix(c(4, 3, 2, 6, 47, 20, 40, 62, 11, 8, 7, 22, 3, 0, 1, 10), 4, 4) rowdistr(z)
These data record the number of seeds (out of 100) that germinated when given different amounts of water. The seeds were either exposed to light or kept in the dark. Four identical boxes were used for each combination of water and light
A data frame with 48 observations on 3 variables.
Factor indicating whether the seeds were exposed to light (N = No, Y = Yes).
Integer amount of water, with higher levels corresponding to more water (1, 2, 3, 4, 5, 6).
Integer number of seeds that germinated, out of 100.
Weight measurements for sheep under combinations of copper and cobalt supplementation.
A data frame with 100 observations on 3 variables.
Integer Weight of sheep (kilograms, kg).
Factor indicating whether copper supplementation was given (No, Yes).
Factor indicating whether cobalt supplementation was given (No, Yes).
Calculates the skewness statistic of the data in 'x'. Values close to zero correspond to reasonably symmetric data, positive values of this measure indicate right-skewed data whereas negative values indicate left-skewness.
skewness(x, ...)skewness(x, ...)
x |
vector containing the data. |
... |
any other variables to be passed to |
Returns the value of the skewness.
## Merger data: data(mergers.df) skewness(mergers.df$mergerdays)## Merger data: data(mergers.df) skewness(mergers.df$mergerdays)
Male Egyptian skulls from five different epochs. Each skull has had four measurements taken of it, BH, Basibregmatic Height, BL, Basialveolar Length, MB, Maximum Breadth and NH, Nasal Height. It is of interest to investigate the change in shape over time. A gradual change, would indicate inbreeding of the populations. This data only includes the maximum breadth measurements.
A data frame with 150 observations on 2 variables.
Integer maximum breadth measurement of the skull.
Integer epoch year group for the skull.
A Handbook of Small Data Sets
Hand, D.J., Daly, F., Lunn, A.D., McConway, K.J. and Ostrowski, E. (1994). A Handbook of Small Data Sets. Boca Raton, Florida: Chapman and Hall/CRC.
Thomson, A. and Randall-Maciver, R. (1905). Ancient Races of the Thebaid. Oxford: Oxford University Press.
Weight and length measurements of 844 snapper (Pagrus auratus) caught in the Hauraki Gulf, near Auckland, New Zealand.
A data frame with 844 observations on 2 variables.
Numeric fork length in centimetres. Fork length is measured from the tip of the snout to the fork of the tail.
Numeric weight of the fish, in kilograms.
Russell Millar, University of Auckland.
Data from an experiment to examine the effects of different planting times on the yield of soya beans, given four different cultivars.
A data frame with 32 observations on 3 variables.
Numeric Yield of each plant.
Factor Cultivar used (cult1, cult2, cult3, cult4)
Factor Month of planting (Novemb, Decemb)
Littler, R. University of Waikato
'stripqq()' is deprecated and is no longer exported. It draws strip charts and normal quantile quantile plots of x for each value of the grouping variable g
stripqq(formula, ...) ## S3 method for class 'formula' stripqq(formula, data = NULL, ...)stripqq(formula, ...) ## S3 method for class 'formula' stripqq(formula, data = NULL, ...)
formula |
A symbolic specification of the form |
data |
An optional data frame in which to evaluate the formula |
... |
Optional arguments that are passed to the |
stripqq(formula): Strip charts and normal quantile-quantile plots
This is a legacy teaching helper retained for compatibility with older course material. New teaching material should prefer current diagnostic plotting workflows.
## Zoo data data(zoo.df) s20x:::stripqq(attendance~day.type, data = zoo.df)## Zoo data data(zoo.df) s20x:::stripqq(attendance~day.type, data = zoo.df)
Displays summary information for a one-way anova analysis. The lm object must come from a numerical response variable and a single factor. The output includes: (i) anova table; (ii) numeric summary; (iii) table of effects; (iv) plot of data with intervals.
summary1way( fit, digit = 5, conf.level = 0.95, inttype = "tukey", pooled = TRUE, print.out = TRUE, draw.plot = TRUE, ... )summary1way( fit, digit = 5, conf.level = 0.95, inttype = "tukey", pooled = TRUE, print.out = TRUE, draw.plot = TRUE, ... )
fit |
an lm object, i.e. the output from |
digit |
decimal numbers after the point. |
conf.level |
confidence level of the intervals. |
inttype |
three options for intervals appeared on plot: 'hsd','lsd' or 'ci'. |
pooled |
two options: pooled or unpooled standard deviation used for plotted intervals. |
print.out |
if |
draw.plot |
if |
... |
more options. |
Invisibly returns a list containing the one-way ANOVA summary components used in the printed teaching output. The list contains:
Df |
degrees of freedom for between groups, within groups, and total. |
Sum of Sq |
sum of squares for between groups, within groups, and total. |
Mean Sq |
mean squares for between groups and within groups. |
F value |
the one-way ANOVA F statistic. |
Pr(F) |
the P-value associated with the F test. |
Main Effect |
the grand mean of the response. |
Group Effects |
group deviations from the grand mean. |
The printed ANOVA table, numeric summary, effects table, and optional plot are the primary teaching interface. The returned list is invisible so classroom use can focus on the printed output while programmatic callers can still inspect the computed values.
summary2way, anova, aov, dummy.coef, onewayPlot
## Computer questionnaire data: data(computer.df) computer.df = within(computer.df, { selfassess = factor(selfassess) }) computer.fit = lm(score ~ selfassess, data = computer.df) result = summary1way(computer.fit) result## Computer questionnaire data: data(computer.df) computer.df = within(computer.df, { selfassess = factor(selfassess) }) computer.fit = lm(score ~ selfassess, data = computer.df) result = summary1way(computer.fit) result
Displays summary information for a two-way anova analysis. The lm object
must come from a numerical response variable and factors. The output depends
on the value of page:
summary2way( fit, page = c("table", "means", "effects", "interaction", "nointeraction"), digit = 5, conf.level = 0.95, print.out = TRUE, new = TRUE, all = FALSE, FUN = "identity", ... )summary2way( fit, page = c("table", "means", "effects", "interaction", "nointeraction"), digit = 5, conf.level = 0.95, print.out = TRUE, new = TRUE, all = FALSE, FUN = "identity", ... )
fit |
an lm object, i.e. the output from |
page |
options for output: |
digit |
the number of decimal places in the display. |
conf.level |
confidence level of the intervals. |
print.out |
if |
new |
if |
all |
Only applicable to |
FUN |
optional function to be applied to estimates and confidence intervals. Typically for backtransformation operations. |
... |
other arguments such as |
page = "table": ANOVA table.
page = "means": cell means matrix and numeric summary.
page = "effects": table of effects.
page = "interaction": interaction contrast tables.
page = "nointeraction": main-effect contrast tables.
'summary2way()' prints the requested teaching summary page and invisibly returns the current summary components. The returned list has the following components:
Df |
degrees of freedom for regression, residual and total. |
Sum of
Sq |
sum squares for regression, residual and total. |
Mean
Sq |
mean squares for regression and residual. |
F
value |
F-statistic value. |
Pr(F) |
The P-value associated with each F-test. |
Grand Mean |
The overall mean of the response variable. |
Row Effects |
The main effects for the first (row) factor. |
Col Effects |
The main effects for the second (column) factor. |
Interaction Effects |
The
interaction effects if an interaction model has been fitted,
otherwise |
results |
If |
.
summary1way, model.tables,
TukeyHSD
## Arousal data: data(arousal.df) arousal.fit = lm(arousal ~ gender * picture, data = arousal.df) summary2way(arousal.fit) ## Butterfat data: data("butterfat.df") fit = lm(log(Butterfat) ~ Breed + Age, data = butterfat.df) summary2way(fit, page = "nointeraction", FUN = exp)## Arousal data: data(arousal.df) arousal.fit = lm(arousal ~ gender * picture, data = arousal.df) summary2way(arousal.fit) ## Butterfat data: data("butterfat.df") fit = lm(log(Butterfat) ~ Breed + Age, data = butterfat.df) summary2way(fit, page = "nointeraction", FUN = exp)
Produces a table of summary statistics for the data. If the argument
group is missing, calculates a matrix of summary statistics for the
data in x. If group is present, the elements of group
are interpreted as group labels and the summary statistics are displayed for
each group separately.
summaryStats(x, ...) ## Default S3 method: summaryStats( x, group = rep("Data", length(x)), data.order = TRUE, digits = 2, ... ) ## S3 method for class 'formula' summaryStats(x, data = NULL, data.order = TRUE, digits = 2, ...) ## S3 method for class 'matrix' summaryStats(x, data.order = TRUE, digits = 2, ...)summaryStats(x, ...) ## Default S3 method: summaryStats( x, group = rep("Data", length(x)), data.order = TRUE, digits = 2, ... ) ## S3 method for class 'formula' summaryStats(x, data = NULL, data.order = TRUE, digits = 2, ...) ## S3 method for class 'matrix' summaryStats(x, data.order = TRUE, digits = 2, ...)
x |
either a single vector of values, a formula of the form
|
... |
Optional arguments that are passed to the summary statistic functions.
For example |
group |
a vector of group labels. |
data.order |
if |
digits |
the number of decimal places to display. |
data |
an optional data frame containing the variables in the model. |
A teaching summary is printed as a side effect. The returned value is invisible so that classroom use can focus on the printed summary while programmatic use can still save the result.
If x is a single variable and no grouping is supplied, an invisible
list is returned with the following named items:
min |
Minimum value. |
max |
Maximum value. |
mean |
Mean value. |
var |
Variance – the average of the squares of the deviations of the data values from the sample mean. |
sd |
Standard deviation – the square root of the variance. |
n |
Number of data values – size of the dataset. |
nMissing |
If there are missing values, and |
iqr |
Midspread (IQR) – the range spanned by the central half of the data; the interquartile range. |
skewness |
Skewness statistic – indicates how skewed the data set is. Positive values indicate right-skew data. Negative values indicate left-skew data. |
lq |
Lower quartile. |
median |
Median – the middle value when the batch is ordered. |
uq |
Upper quartile. |
If grouping is provided, either by using the group argument, by using
a formula, or by passing a matrix whose columns represent groups, the
function invisibly returns a data.frame with one row for each group
and columns containing the summary statistics.
summaryStats(default): Summary Statistics
summaryStats(formula): Summary Statistics
summaryStats(matrix): Summary Statistics
## STATS20x data: data(course.df) ## Single variable summary with(course.df, summaryStats(Exam)) ## Using a formula summaryStats(Exam ~ Stage1, course.df) ## Using a matrix courseMatrix = cbind(course.df$Exam, course.df$Assign, course.df$Test) summaryStats(courseMatrix) ## Saving and extracting the information sumStats = summaryStats(Exam ~ Degree, course.df) sumStats ## Just the BAs sumStats['BA', ] ## Just the means sumStats$mean## STATS20x data: data(course.df) ## Single variable summary with(course.df, summaryStats(Exam)) ## Using a formula summaryStats(Exam ~ Stage1, course.df) ## Using a matrix courseMatrix = cbind(course.df$Exam, course.df$Assign, course.df$Test) summaryStats(courseMatrix) ## Saving and extracting the information sumStats = summaryStats(Exam ~ Degree, course.df) sumStats ## Just the BAs sumStats['BA', ] ## Just the means sumStats$mean
Data from an experiment to assess the impact of three different teaching methods on language ability. 30 students were randomly allocated into three groups, one for each method. The students' IQ before instruction and a language test score after instruction were recorded.
A data frame with 30 observations on 3 variables.
Numeric Language test score after instruction.
Numeric Student's IQ.
Factor Teaching method (1, 2, 3)
Salary information for all salaried employees of the Technitron Company.
A data frame with 46 observations on 8 variables.
Numeric Annual Salary (dollars)
Numeric Number of years employed at Technitron.
Numeric Number of years prior experience.
Numeric Years of education after high school.
Numeric Company identification number.
Numeric Gender (0 = female, 1 = male)
Numeric Department employee works in (1 = Sales, 2 = Purchasing, 3 = Advertising, 4 = Engineering)
Numeric Number of employees supervised.
Data from an experiment to assess the effect of a new drug on the weight of the thyroid gland using 16 laboratory animals. The animals were randomly assigned into either a control group, or a treatment group, and each animal had its bodyweight recorded at the beginning of the experiment and its thyroid weight measured at the end of the experiment.
A data frame with 16 observations on 3 variables.
Numeric Weight of thyroid gland after 7 days (mg)
Numeric Animal body weight before experiment began (g)
Factor Animal's group (1 = control, 2 = drug)
Two random samples of households, one of households who purchase Crest toothpaste and one of households who do not. For each household the age is recorded of the person responsible for purchasing the toothpaste.
A data frame with 20 observations on 2 variables.
Numeric Age of the person in the household responsible for purchases of Crest.
Numeric Age of the person in the household responsible for purchases of other brands of toothpaste.
Plots a scatter plot for the variables x, y along with a lowess smooth for the underlying trend. One standard deviation error bounds for the scatter about this trend are also plotted.
trendscatter(x, ...) ## Default S3 method: trendscatter(x, y = NULL, f = 0.5, xlab = NULL, ylab = NULL, main = NULL, ...) ## S3 method for class 'formula' trendscatter( x, f = 0.5, data = NULL, xlab = NULL, ylab = NULL, main = NULL, ... )trendscatter(x, ...) ## Default S3 method: trendscatter(x, y = NULL, f = 0.5, xlab = NULL, ylab = NULL, main = NULL, ...) ## S3 method for class 'formula' trendscatter( x, f = 0.5, data = NULL, xlab = NULL, ylab = NULL, main = NULL, ... )
x |
the coordinates of the points in the scatter plot. Alternatively, a formula. |
... |
Optional arguments |
y |
the y coordinates of the points in the plot, ignored if |
f |
the smoother span. This gives the proportion of points in the plot which influence the smooth at each value. Larger values give more smoothness. |
xlab |
a title for the x axis: see |
ylab |
a title for the y axis: see |
main |
a title for the plot: see |
data |
an optional data frame containing the variables in the model. |
Returns the plot.
trendscatter(default): Trend and scatter plot
trendscatter(formula): Trend and scatter plot
# Synthetic teaching example: a simple polynomial set.seed(123) x = rnorm(100) e = rnorm(100) y = 2 + 3 * x - 2 * x^2 + 4 * x^3 + e trendscatter(y ~ x) # Synthetic teaching example: an exponential growth curve e = rnorm(100, 0, 0.1) y = exp(5 + 3 * x + e) trendscatter(log(y) ~ x) # Peruvian Indians data data(peru.df) trendscatter(BP ~ weight, data = peru.df) # Note: this usage is deprecated with(peru.df, trendscatter(weight, BP))# Synthetic teaching example: a simple polynomial set.seed(123) x = rnorm(100) e = rnorm(100) y = 2 + 3 * x - 2 * x^2 + 4 * x^3 + e trendscatter(y ~ x) # Synthetic teaching example: an exponential growth curve e = rnorm(100, 0, 0.1) y = exp(5 + 3 * x + e) trendscatter(log(y) ~ x) # Peruvian Indians data data(peru.df) trendscatter(BP ~ weight, data = peru.df) # Note: this usage is deprecated with(peru.df, trendscatter(weight, BP))
'tslm()' is a teaching-friendly wrapper for fitting linear models with optional AR(p) error structures. Students specify the mean model using an ordinary formula and add an 'ar(p)' term to request autoregressive errors.
tslm(formula, data, time, method = "REML", ...)tslm(formula, data, time, method = "REML", ...)
formula |
a model formula. Use 'ar(p)' in the right hand side to specify AR(p) errors, for example 'y ~ x + ar(1)'. |
data |
a data frame containing the variables in the model. |
time |
optional unquoted or quoted name of the time variable in 'data'. If omitted for an AR model, the row order of 'data' is used. |
method |
fitting method passed to [nlme::gls()] for AR models. Defaults to '"REML"'. |
... |
additional arguments passed to [stats::lm()] or [nlme::gls()]. |
When no 'ar(p)' term is present, 'tslm()' fits an ordinary [stats::lm()] model. When an 'ar(p)' term is present, 'tslm()' fits a [nlme::gls()] model with an AR(p) correlation structure using [nlme::corARMA()]. The 'ar(p)' term changes the error model, not the mean-model terms printed in the formula.
The formula describes the mean model, just as it does for [stats::lm()]. The special term 'ar(p)' is removed from the mean model before fitting and is used only to specify the correlation structure for the errors. For example, 'log(passengers) ~ t + month + ar(1)' fits a trend and seasonal mean model with AR(1) errors.
For AR-error models, 'time' should usually name the variable giving the time order of the observations. If 'time' is omitted, 'tslm()' fits the model using the row order of 'data' and gives a warning so that this assumption is visible.
Diagnostic methods for AR-error models use normalised residuals by default, because these residuals account for the fitted correlation structure. Use 'residualType = "response"' when the raw response residuals are required. '"normalised"' and '"normalized"' are both accepted for compatibility.
An object of class 'tslm', containing the original formula, the mean formula fitted internally, the AR order, the time variable if supplied, and the underlying fitted model.
[stats::lm()], [nlme::gls()], [nlme::corARMA()]
data(beer.df) fit = tslm(beer ~ t + ar(1), data = beer.df, time = t) coef(fit) data(airpass.df) fitAr = tslm(log(passengers) ~ t + month + ar(1), data = airpass.df, time = t ) summary(fitAr) anova(fitAr) plot(fitAr) plot(fitAr, residualType = "response")data(beer.df) fit = tslm(beer ~ t + ar(1), data = beer.df, time = t) coef(fit) data(airpass.df) fitAr = tslm(log(passengers) ~ t + month + ar(1), data = airpass.df, time = t ) summary(fitAr) anova(fitAr) plot(fitAr) plot(fitAr, residualType = "response")
Data for 455 days of attendance records for Auckland Zoo, from January 1, 1993. Note that only 440 values are given due to missing values. It was of interest to assess whether an advertising campaign was effective in increasing attendance.
A data frame with 440 observations on 6 variables.
Numeric Number of visitors.
Numeric Time in days since the start of the study.
Numeric Hours of sunshine the previous day.
Numeric Average spending on TV advertising in the previous week (1000s of dollars per day)
Factor Assessment based on number of hours of sunshine (0 = No, 1 = Yes)
Factor Type of day (1 = ordinary weekday, 2 = weekend day, 3 = school holiday weekday, 4 = public holiday)