The Student’s \(t\) Distribution

We use the \(t\) instead of the \(z\) distribution when either of the following conditions apply:

  1. \(\sigma\) is unknown
  2. \(\sigma\) is known but \(n < 30\)

\[t = \dfrac{\bar{Y} - \mu}{SE_{\bar{Y}}}\]

\[\text{where } SE_{\bar{Y}} = \dfrac{s}{\sqrt{n}}\]

\[\text{and } df=n-1\]

Working with the \(t\) Distribution

Just as we did with the \(z\) distribution, we can find the area above, below, and within the bounds of specified \(t-scores\) via the pt() command. Note that now we must specify the degrees of freedom.

pt(1.99, df=2, lower.tail=FALSE)
## [1] 0.09243554
pt(-1.99, df=2, lower.tail=TRUE)
## [1] 0.09243554
pt(-1.99, df=2, lower.tail=FALSE) - pt(1.99, df=2, lower.tail=FALSE)
## [1] 0.8151289

We can also calculate the \(t-score\) associated with a certain area under the curve and for a given degrees of freedom via the qt() command:

qt(0.995, df=79, lower.tail=TRUE)
## [1] 2.639505
qt(0.975, df=61, lower.tail=FALSE)
## [1] -1.999624

One-Sample \(t-tests\) and Confidence Intervals

The test-statistic is given by

\[t = \dfrac{\bar{Y} - \mu_{0}}{SE_{\bar{Y}}}\]

where \(\mu_{0}\) is the population mean specified in \(H_0\). As usual, you will have to specify \(H_0\) and \(H_A\), set \(\alpha\), set the decision rule, and then carry out the test. I’ve done this for one of the problems (see below).

Assumptions that must be met

The test is built upon some assumptions that must be met. In particular:

  1. Data represent a random sample
  2. The outcome is a numerical variable that is normally distributed

Testing Assumptions

Assumption (1) is assumed to be true since otherwise there would be no sense in running any kind of a statistical test. Assumption (2), however, must be tested via:

  1. Visual exploration of the data, and
  2. Formal tests of normality.

Visual “tests” of normality

These tests involve examining the box-plot for skew and outliers. With symmetric distributions, the t-test is quite accurate but with skew and/or outliers the test can be inaccurate in small samples. Which then begs the question: How small is “small”? A reasonable guide for a small sample is anything with \(n < 50\), though of course this rule will vary across some disciplines. If, for example, a discipline regularly deals with samples of 11, 15, 20 observations then they may be predisposed to accepting evidence based on these sample sizes.

The classic visual test is via a quantile-quantile plot, often referred to as a Q-Q plot. These plots essentially involve comparing the empirical distribution (i.e., what is visible in our sample) against the theoretical distribution we assume generates our sample data. The theoretical distribution will show up as 45-degree straight line with the empirical distribution superimposed on this 45-degree line. If too many observations deviate from the 45-degree line then we know our sample deviates from the distribution we assume the sample comes from.

set.seed(199166)
y = rnorm(100, 0, 1)
hist(y, col = "cornflowerblue")

boxplot(y, horizontal = TRUE, col = "cornflowerblue")

qqnorm(y); qqline(y, col = "darkred")

Notice how while most of the data are close to the line a few wander off. You also see the lone outlier on the high side (see the box-plot and the Q-Q plot). Of course it isn’t readily apparent when to decide that too many data are wandering off the line and hence that we have a problem. This is precisely why we square the visual evidence against formal tests.

Formal tests of normality

There are several formal tests of normality. Which one should be used is partly dependent upon our sample and partly upon our disciplines. Below I outline some considerations for when to use which test. Note that this only covers two of what are a large battery of tests, largely because these two are the most often used tests. Each of these tests has \(H_0\): The data are drawn from a normal distribution.

  1. The Shapiro-Wilk test – is used with sample sizes of 5,000 or less. If the p-value is \(\leq 0.05\) we end up rejecting the \(H_0\) of normality.
  2. The Anderson-Darling test – is used with sample sizes that exceed 5,000. If the p-value is \(\leq 0.05\) we end up rejecting the \(H_0\) of normality.
shapiro.test(y)
## 
##  Shapiro-Wilk normality test
## 
## data:  y
## W = 0.99535, p-value = 0.9832
library(nortest)
ad.test(y)
## 
##  Anderson-Darling normality test
## 
## data:  y
## A = 0.17566, p-value = 0.9214

Of course, these tests are not foolproof. For example, consider the following example.

set.seed(244263)
x = rt(500000, 200)
ad.test(x)
## 
##  Anderson-Darling normality test
## 
## data:  x
## A = 1.5643, p-value = 0.0005048
qqnorm(x);  qqline(x, col = "darkred")

hist(x, col = "cornflowerblue")

boxplot(x, col = "cornflowerblue", horizontal = TRUE)

y = rnorm(1000, 0, 1)
shapiro.test(y)
## 
##  Shapiro-Wilk normality test
## 
## data:  y
## W = 0.99784, p-value = 0.2228
qqnorm(y);  qqline(y, col = "darkred")

hist(y, col = "cornflowerblue")

boxplot(y, col = "cornflowerblue", horizontal = TRUE)

In brief, even when visual tests show no major problems, formal tests in large samples will often end up rejecting the \(H_0\) even if only a few observations are deviant. In small samples the opposite is true: Even a large number of deviant observations will not lead to a rejection of \(H_0\).

So what should we do? Gain experience, and marry that to rules our discipline seems to follow. For humbler, class purposes, of course, we’ll stick to the results of the formal tests. If the formal test rejects Normality then we will conclude that the normality assumption is being violated, regardless of how trivial or large the problem seems to be in the Q-Q plots.

Normality is rejected, now what?

When the assumption of normality is rejected we have no choice but to try and transform the data. There are several transformations at our disposal but for now we’ll focus on one – the natural logarithm. If we convert the original variable via the natural logarithm and then the formal tests no longer allow us to reject the assumption of normality, we can then carry out the t-test with the log-transformed version of the outcome variable. If the log-transform doesn’t work to “ensure” normality then we may have to try other transformations. We’ll cover these in subsequent chapters. For now we’ll focus on testing for normality, both visually and with formal tests.

Practice Problems

Problem 9

Community ecologists draw “food webs” to describe the predator-prey relationships between all organisms living in an area. A theoretical model predicts that a measure of the structure of the food webs called “diet discontinuity” should be zero. Diet discontinuity is a measure of the relative number of predators whose prey are not ordered contiguously. Researchers have measured discontinuity scores for seven different food webs in nature. The values are given below (see R code). Assume that discontinuity has a normal distribution. Are these results consistent with the model’s prediction of a zero mean discontinuity score? Carry out an appropriate hypothesis test.

\(H_0:\) Mean discontinuity is zero \((\mu=0)\)
\(H_A:\) Mean discontinuity is not zero \((\mu \neq 0)\)
\(\alpha=0.05\)

y = c(0.35, 0.08, 0.004, 0.08, 0.32, 0.28, 0.17)

Testing for Normality

Don’t forget: The \(H_0\) is of Normality.

boxplot(y, col = "salmon", horizontal = TRUE, xlab = "Discontinuity Score")

qqnorm(y);  qqline(y, col = "darkred")

shapiro.test(y)
## 
##  Shapiro-Wilk normality test
## 
## data:  y
## W = 0.91513, p-value = 0.4325

The p-value is 0.4325 so we fail to reject the Null of Normality, and can proceed with the t-test.

The t-test

Note the structure of the t.test() command. The null hypothesis is specified via mu = ?, and then the alternative hypothesis plus the confidence interval are specified as usual.

t.test(y, mu = 0, alternative = "two.sided", conf.level = 0.95)
## 
##  One Sample t-test
## 
## data:  y
## t = 3.5925, df = 6, p-value = 0.01147
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  0.05849238 0.30836477
## sample estimates:
## mean of x 
## 0.1834286

The \(p-value < 0.05\) so we can reject \(H_0\). These data suggest that mean diet discontinuity of food webs is not \(0\). We can be about 95% confident that TRUE mean diet discontinuity falls in the interval given by 0.0584 and 0.3083.

Problem 13

Can a human swim faster in water or in syrup? An experiment was done in which researchers filled one pool with water mixed with syrup and another pool with normal water. They had 18 swimmers swim laps in both pools in random order. The data are presented as the relative speed of each swimmer – the speed in syrup divided by the speed in water. If the syrup has no effect then the relative speed should have a mean \(=1\).

\(H_0:\) Relative speed is the same in syrup as in water \((\mu=1)\)
\(H_A:\) Relative speed is not the same in syrup as in water \((\mu \neq 1)\)
\(\alpha = 0.05\)

syrup = read.csv(url("http://whitlockschluter.zoology.ubc.ca/wp-content/data/chapter11/chap11q13SyrupSwimming.csv"))

head(syrup); names(syrup)
##   relativeSpeedInSyrup
## 1                 1.08
## 2                 0.94
## 3                 1.07
## 4                 1.04
## 5                 1.02
## 6                 1.01
## [1] "relativeSpeedInSyrup"
hist(syrup$relativeSpeedInSyrup, xlab = "Relative Speed", col = "cornflowerblue", main = "Distribution of Relative Speed")

boxplot(syrup$relativeSpeedInSyrup, xlab = "Relative Speed", col = "cornflowerblue", main = "Distribution of Relative Speed", horizontal = TRUE)

qqnorm(syrup$relativeSpeedInSyrup);  qqline(syrup$relativeSpeedInSyrup, col = "darkred")

shapiro.test(syrup$relativeSpeedInSyrup)
## 
##  Shapiro-Wilk normality test
## 
## data:  syrup$relativeSpeedInSyrup
## W = 0.94257, p-value = 0.3205

The histogram and the box-plot show skewed distributions with average relative speed around 1.02. The Shapiro-Wilk test does not allow us to reject the Null of Normality so we may proceed with the t-test.

t.test(syrup$relativeSpeedInSyrup, mu = 1, conf.level = 0.95, alternative = "two.sided")
## 
##  One Sample t-test
## 
## data:  syrup$relativeSpeedInSyrup
## t = 1.1695, df = 17, p-value = 0.2583
## alternative hypothesis: true mean is not equal to 1
## 95 percent confidence interval:
##  0.9906203 1.0327130
## sample estimates:
## mean of x 
##  1.011667

With the \(p-value > 0.05\) we are unable to reject \(H_0\); the data do not allow us to reject the notion that average relative speed is 1.

They want the 99% confidence interval so we’ll use the Rmisc package rather than running the t.test() with conf.level = 0.99; Rmisc would be more efficient.

library(Rmisc)
## Loading required package: lattice
## Loading required package: plyr
CI(syrup$relativeSpeedInSyrup, ci = 0.99)
##     upper      mean     lower 
## 1.0405778 1.0116667 0.9827555

The 99% confidence interval is \([0.9827, 1.0405]\).

Problem 16

In the Northern Hemisphere dolphins swim predominantly in a counterclockwise direction while sleeping. A group of researchers wanted to know whether the same was true of dolphins in the Southern Hemisphere. They watched eight sleeping dolphins and recorded the percentage of time they swam clockwise.

  1. What is the mean percentage of clockwise swimming for the Southern Hemisphere dolphins?
swim = read.csv(url("http://whitlockschluter.zoology.ubc.ca/wp-content/data/chapter11/chap11q16DolphinsClockwise.csv"))
head(swim); names(swim)
##   percentClockwise
## 1             77.7
## 2             84.8
## 3             79.4
## 4             84.0
## 5             99.6
## 6             93.6
## [1] "percentClockwise"
mean(swim$percentClockwise)
## [1] 88.2125

The mean percentage of time spent swimming clockwise is 88.2125

  1. and c. What are the 95% and 99% confidence intervals for the mean?
library(Rmisc)
CI(swim$percentClockwise, ci = 0.95)
##    upper     mean    lower 
## 94.97821 88.21250 81.44679
CI(swim$percentClockwise, ci = 0.99)
##    upper     mean    lower 
## 98.22529 88.21250 78.19971

These are \([81.4467, 94.9782]\) and \([78.1997, 98.2252]\), respectively.

  1. and e. What is the best estimate of the standard deviation and the median?
sd(swim$percentClockwise)
## [1] 8.092755
median(swim$percentClockwise)
## [1] 87.1

The standard deviation is 8.0927 and the median is 87.1