Basic Statistics

The most common tools used in statistics are means, variances, correlations, and t-tests. These
are used in R with easy-to-use functions such as mean, var., cor, and t-test.
18.1 Summary Statistics
Summary statistics are a set of measures that describe the main features of a dataset, giving a
quick overview of its distribution, central tendency, and spread.
Common summary statistics include:
• Minimum (Min)
• First Quartile (1st Qu.)
• Median (Median)
• Mean (Mean)
• Third Quartile (3rd Qu.)
• Maximum (Max)
• Variance(var)
• Correlation(cor)
• T-test
Using summary() in R
The summary() function in R provides all these statistics at once for numeric vectors, factors,
or entire data frames.
Example 1: Numeric Vector
x<-sample(x=1:100)
> summary(x)
Min. 1st Qu. Median Mean 3rd Qu. Max.
1.00 25.75 50.50 50.50 75.25 100.00

Randomly set 20 elements using sample to set to NA

y[sample(x = 1:100, size = 20, replace = FALSE)] <- NA
sample(x = 1:100, size = 20, replace = FALSE)

• Randomly selects 20 unique numbers between 1 and 100.

• replace = FALSE ensures no repeats.
y[ ... ] <- NA
• Replaces the elements at the randomly chosen positions with NA.

> x[sample(x = 1:100, size = 20, replace = FALSE)] <- NA

>x

[1] 58 49 NA NA 17 55 96 32 45 62 88 NA 52 NA 78 85 100 NA 89 15 57 2 5
[24] 97 83 1 11 35 81 48 67 NA 14 NA NA 4 9 NA 91 66 99 7 92 75 NA NA
[47] 19 36 71 NA 29 6 72 21 50 90 68 33 24 25 NA NA 61 38 NA 87 94 69 77
[70] 8 12 18 98 NA 31 53 NA 37 82 84 NA 56 76 16 34 95 NA 3 74 80 13 46
[93] 93 51 22 40 73 20 39 NA

> mean(x)
[1] NA

> mean(x, [Link]=TRUE)

[1] 51.1125
To calculate Weighted mean
> grades=c(95,72,87,66)

1
 > weights=c(1/2,1/4,1/8,1/8)
> mean(grades)
[1] 80
> [Link](x=grades,w=weights)
[1] 84.625
Formula for Weighted Mean
The weighted mean is a type of average where different values contribute differently
according to their weights

y=sample(x=1:100)
>y
[1] 1 5 52 14 93 30 55 33 53 46 91 26 37 92 27 69 61 78 39 10 97 9 82
[24] 4 58 75 41 28 38 12 25 62 19 3 73 65 17 51 45 79 40 29 31 89 20 71
[47] 2 96 66 21 24 86 36 74 54 63 85 83 8 100 44 56 47 77 34 59 90 42 80
[70] 84 16 7 70 22 99 94 57 18 87 95 6 13 68 72 98 60 88 23 81 76 11 43
[93] 48 49 32 64 67 50 15 35

> var(y)
[1] 841.6667
Variance Formula

>y
[1] 1 5 52 14 93 30 55 33 53 46 91 26 37 92 27 69 61 78 39 10 97 9 82
[24] 4 58 75 41 28 38 12 25 62 19 3 73 65 17 51 45 79 40 29 31 89 20 71
[47] 2 96 66 21 24 86 36 74 54 63 85 83 8 100 44 56 47 77 34 59 90 42 80
[70] 84 16 7 70 22 99 94 57 18 87 95 6 13 68 72 98 60 88 23 81 76 11 43
[93] 48 49 32 64 67 50 15 35

> sd(y)
[1] 29.01149

2
>x
[1] 58 49 NA NA 17 55 96 32 45 62 88 NA 52 NA 78 85 100 NA 89 15 57 2 5
[24] 97 83 1 11 35 81 48 67 NA 14 NA NA 4 9 NA 91 66 99 7 92 75 NA NA
[47] 19 36 71 NA 29 6 72 21 50 90 68 33 24 25 NA NA 61 38 NA 87 94 69 77
[70] 8 12 18 98 NA 31 53 NA 37 82 84 NA 56 76 16 34 95 NA 3 74 80 13 46
[93] 93 51 22 40 73 20 39 NA

> sd(x)
[1] NA

> sd(x,[Link]=TRUE)
[1] 30.88771

> min(y)
[1] 1

> max(y)
[1] 100

> median(y)
[1] 50.5

> min(x)
[1] NA

> min(x,[Link]=TRUE)
[1] 1

> summary(x)
Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
1.00 21.75 51.50 51.11 78.50 100.00 20

> summary(y)
Min. 1st Qu. Median Mean 3rd Qu. Max.
1.00 25.75 50.50 50.50 75.25 100.00

> quantile(y)
0% 25% 50% 75% 100%
1.0 25.75 50.50 75.25 100.00

> quantile(y,probs=c(.25,.75))
25% 75%
25.75 75.25

> quantile(x)
Error in [Link](x) :
missing values and NaN's not allowed if '[Link]' is FALSE

> quantile(x,[Link]=TRUE)
0% 25% 50% 75% 100%

3
 1.00 21.75 51.50 78.50 100.00

> quantile(y,probs=c(.1,.25,.5,.75,.99))
10% 25% 50% 75% 99%
10.90 25.75 50.50 75.25 99.01

> quantile(x,probs=c(.1,.25,.5,.75,.99))
Error in [Link](x, probs = c(0.1, 0.25, 0.5, 0.75, 0.99)) :
missing values and NaN's not allowed if '[Link]' is FALSE

> quantile(x,probs=c(.1,.25,.5,.75,.99),[Link]=TRUE)
10% 25% 50% 75% 99%
8.90 21.75 51.50 78.50 99.21

18.2 Correlation and Covariance

To test the relationship between more than one variable, we use covariance, and
correlation is used.
What Is Correlation?
Correlation measures the strength and direction of the linear relationship between two
numeric variables.
• Range: -1 to +1
o +1 → perfect positive correlation
o -1 → perfect negative correlation
o 0 → no linear correlation
Positive correlation: as one variable increases, the other increases.
Negative correlation: as one variable increases, the other decreases.

Using cor() in R
Basic Syntax
cor(x, y, method = "pearson")
• x, y → numeric vectors
• method → correlation method:
o "pearson" → default, measures linear correlation
o "spearman" → rank-based correlation (non-parametric)
o "kendall" → rank correlation, less sensitive to ties

To understand the concept, we use the economics dataset from ggplot2.

[Link]("ggplot2")
library(ggplot2)
head(economics)
date pce pop psavert uempmed unemploy
1 1967-07-01 506.7 198712 12.6 4.5 2944
2 1967-08-01 509.8 198911 12.6 4.7 2945
3 1967-09-01 515.6 199113 11.9 4.6 2958
4 1967-10-01 512.2 199311 12.9 4.9 3143
5 1967-11-01 517.4 199498 12.8 4.7 3066
6 1967-12-01 525.1 199657 11.8 4.8 3018

> cor(economics$pce, economics$psavert)

4
[1] -0.7928546

To compute Pearson’s coefficient, we multiply deviations from the mean for X times those
for Y and divide by the product of the standard deviations. Here is the formula:

The correlation (r) always ranges between -1 and +1.

• When r = +1, it indicates a perfect positive correlation, meaning both variables
increase together in the same direction.
• When r = -1, it indicates a perfect negative correlation, meaning as one variable
increases, the other decreases.
• When r = 0, it indicates no linear relationship between the two variables.

> xpart<-economics$pce-mean(economics$pce)
> ypart<-economics$psavert-mean(economics$psavert)
> n<-(nrow(economics)-1)
> xsd<-sd(economics$pce)
> ysd<-sd(economics$psavert)
> sum(xpart*ypart)/(n*xsd*ysd)
[1] -0.7928546

> cor(economics[,c(2,4:6)])
pce psavert uempmed unemploy
pce 1.0000000 -0.7928546 0.7269616 0.6145176
psavert -0.7928546 1.0000000 -0.3251377 -0.3093769
uempmed 0.7269616 -0.3251377 1.0000000 0.8693097
unemploy 0.6145176 -0.3093769 0.8693097 1.0000000

To visualize the information using a plot. ggpairs from the GGally package is used to
visualize this.
[Link]("GGally")
library(GGally)

ggpairs(economics[, c(2,4:6)])

5
 Figure 18.1 Pairs plot of economics data showing the relationship between each pair of
variables as a scatterplot with the correlations printed as numbers.

This is similar to a small multiples plot except that each pane has different x and y axes. This
plot displays the original data, but it does not actually show the correlation. To show this, we
build a heatmap of the correlation numbers as shown in fig 18.2.
• High positive correlation shows a positive relationship between variables.
• High negative correlation indicates a relationship between the variables, and
• Non-zero correlation shows no strong relationship.
> library(reshape2)
> library(scales)
> econcor=cor(economics[,c(2,4:6)])
> econmelt=melt(econcor,varnames=c("x","y"),[Link]="Correlation")
> print(econmelt)
x y Correlation
1 pce pce 1.0000000
2 psavert pce -0.7928546
3 uempmed pce 0.7269616
4 unemploy pce 0.6145176
5 pce psavert -0.7928546
6 psavert psavert 1.0000000
7 uempmed psavert -0.3251377
8 unemploy psavert -0.3093769
9 pce uempmed 0.7269616
10 psavert uempmed -0.3251377
11 uempmed uempmed 1.0000000
12 unemploy uempmed 0.8693097
13 pce unemploy 0.6145176
14 psavert unemploy -0.3093769
15 uempmed unemploy 0.8693097
16 unemploy unemploy 1.0000000

> econmelt=econmelt[order(econmelt$Correlation),]
> print(econmelt)
x y Correlation
2 psavert pce -0.7928546

6
 5 pce psavert -0.7928546
7 uempmed psavert -0.3251377
10 psavert uempmed -0.3251377
8 unemploy psavert -0.3093769
14 psavert unemploy -0.3093769
4 unemploy pce 0.6145176
13 pce unemploy 0.6145176
3 uempmed pce 0.7269616
9 pce uempmed 0.7269616
12 unemploy uempmed 0.8693097
15 uempmed unemploy 0.8693097
1 pce pce 1.0000000
6 psavert psavert 1.0000000
11 uempmed uempmed 1.0000000
16 unemploy unemploy 1.0000000

1. Loading libraries
library(reshape2) # for melting the data
library(scales) # for muted() colors and scaling functions
• reshape2: provides the melt() function that reshapes a wide matrix/data frame into a
long format (needed for ggplot2 heatmaps).
• scales: provides helpers like muted("red"), color scales, and formatting functions.

2. Building the correlation matrix

econcor = cor(economics[, c(2, 4:6)])
• economics is a built-in dataset in ggplot2.
• economics[, c(2,4:6)] selects columns 2, 4, 5, and 6 (variables like population,
unemployment, etc.).
• cor() computes the correlation matrix between these variables.
Result: a symmetric square matrix with correlation values ranging from -1 to 1.

3. Melting into long format

econmelt = melt(econcor, varnames = c("x", "y"), [Link] = "Correlation")
print(econmelt)
• melt() converts the matrix into a long-format data frame with three columns:
o x → variable name from row
o y → variable name from column
o Correlation → correlation value
• This makes it compatible with ggplot2, where each (x, y) pair maps to a correlation
value.

4. Ordering by correlation
econmelt = econmelt[order(econmelt$Correlation), ]
print(econmelt)
• order() sorts rows by correlation values.
• This is helpful if you want to inspect strongest/weakest correlations in tabular form.
• For plotting, order doesn’t affect the heatmap (but it helps when printing the table).

5. Plotting the heatmap

ggplot(econmelt, aes(x = x, y = y)) +
geom_tile(aes(fill = Correlation)) +
scale_fill_gradient2(
low = muted("red"), # negative correlations

7
 mid = "white", # neutral
high = "steelblue", # positive correlations
guide = guide_colorbar(ticks = FALSE, barheight = 10),
limits = c(-1, 1) # fix the scale between -1 and 1
)+
theme_minimal() +
labs(x = NULL, y = NULL)
• geom_tile() draws squares for each (x,y) pair.
• scale_fill_gradient2() maps correlations to colors:
o red = strong negative correlation
o white = no correlation
o blue = strong positive correlation
• theme_minimal() gives a clean theme.
• labs(x=NULL, y=NULL) removes axis labels (since variable names are already
visible).

Figure 18.2 Heatmap of the correlation of the economics data. Diagonal -with correlation 1.
Red indicates a highly negative correlation, blue indicates a highly positive correlation, and
white indicates no correlation.

Handling Missing Data in Correlation Calculations in R

When calculating correlations using the cor() function in R, the dataset may contain missing
values (NA). By default, missing values can cause problems because correlation requires
complete numerical data.
To handle missing values, the use argument in the cor() function is used. It specifies how R
should treat missing observations.

Syntax of cor() with missing data options

cor(x, use = "method", method = "pearson")
• x → data frame or matrix of numeric values
• use → method for handling missing values
• method → type of correlation: "pearson", "spearman", or "kendall"

8
Options for use
1. "everything" (default)
• Assumes that the data has no missing values.
• If any NA is present, the result will be NA.
• Strictest option.
cor(df, use = "everything")

2. "[Link]"
• Assumes there are no missing values in the dataset.
• If missing values are found, R will throw an error.
• Safe only when data is already cleaned.
cor(df, use = "[Link]")

3. "[Link]"
• Uses only the rows that have no missing values in any variable.
• Produces a full correlation matrix without NAs.
• May reduce the dataset size, but ensures consistency.
cor(df, use = "[Link]")

4. "[Link]"
• Similar to "[Link]", but avoids unnecessary checks if there are no missing
values at all.
• Slightly faster than "[Link]".
cor(df, use = "[Link]")

5. "[Link]"
• Uses all available pairs of variables.
• If two variables have missing data in different places, R still computes the correlation
for the available pairs.
• Maximizes data usage but may result in a correlation matrix that is not mathematically
consistent (not positive definite).
cor(df, use = "[Link]")

Summary Table

Option Behavior with Missing Data

"everything" Returns NA if any missing values are present.

"[Link]" Assumes no missing data; gives error if NA is found.

"[Link]" Uses only rows with no missing values at all.

"[Link]" Same as [Link], but faster when no NAs exist.

"[Link]" Uses all available pairs; may cause inconsistencies.

Example
Step 1: Create the vectors
m <- c(9, 9, NA, 3, NA, 5, 8, 1, 10, 4)

9
n <- c(2, NA, 1, 6, 6, 4, 1, 1, 6, 7)
p <- c(8, 4, 3, 9, 10, NA, 3, NA, 9, 9)
q <- c(10, 10, 7, 8, 4, 2, 8, 5, 5, 2)
r <- c(1, 9, 7, 6, 5, 6, 2, 7, 9, 10)
• Each vector has 10 elements.
• Some of them contain NA values (m, n, p).

Step 2: Combine into a matrix

mat <- cbind(m, n, p, q, r)
print(mat)
The result:
m n p q r
[1,] 9 2 8 10 1
[2,] 9 NA 4 10 9
[3,] NA 1 3 7 7
[4,] 3 6 9 8 6
[5,] NA 6 10 4 5
[6,] 5 4 NA 2 6
[7,] 8 1 3 8 2
[8,] 1 1 NA 5 7
[9,] 10 6 9 5 9
[10,] 4 7 9 2 10

Step 3: Correlation with different use methods

1. "everything"
cor(mat, use = "everything")
• Default method: returns NA if there are missing values.
• Only q and r have full overlap (no shared NA), so their correlation is computed.
Result:
m n p q r
m 1 NA NA NA NA
n NA 1 NA NA NA
p NA NA 1 NA NA
q NA NA NA 1.0000000 -0.4242958
r NA NA NA -0.4242958 1.0000000

2. "[Link]"
cor(mat, use = "[Link]")
• Requires no missing values at all.
• Since m, n, and p contain NA, R throws an error:
Error in cor(mat, use = "[Link]") : missing observations in cov/cor

3. "[Link]"
cor(mat, use = "[Link]")
• Uses only rows with no missing values across all variables.
• Here, the complete rows are: 1, 4, 7, 9, 10 (5 rows).
• Correlation is computed on these.
Result:
m n p q r
m 1.0000000 -0.5228840 -0.2893527 0.2974398 -0.3459470

10
n -0.5228840 1.0000000 0.8090195 -0.7448453 0.9350718
p -0.2893527 0.8090195 1.0000000 -0.3613720 0.6221470
q 0.2974398 -0.7448453 -0.3613720 1.0000000 -0.9059384
r -0.3459470 0.9350718 0.6221470 -0.9059384 1.0000000

4. "[Link]"
cor(mat, use = "[Link]")
• Same as "[Link]".
• Gives identical result (just skips some checks internally).
Result:
m n p q r
m 1.0000000 -0.5228840 -0.2893527 0.2974398 -0.3459470
n -0.5228840 1.0000000 0.8090195 -0.7448453 0.9350718
p -0.2893527 0.8090195 1.0000000 -0.3613720 0.6221470
q 0.2974398 -0.7448453 -0.3613720 1.0000000 -0.9059384
r -0.3459470 0.9350718 0.6221470 -0.9059384 1.0000000

5. Manual selection of complete rows

cor(mat[c(1,4,7,9,10), ])
• Here you manually pick rows without NA.
• This produces the same result as "[Link]".
m n p q r
m 1.0000000 -0.5228840 -0.2893527 0.2974398 -0.3459470
n -0.5228840 1.0000000 0.8090195 -0.7448453 0.9350718
p -0.2893527 0.8090195 1.0000000 -0.3613720 0.6221470
q 0.2974398 -0.7448453 -0.3613720 1.0000000 -0.9059384
r -0.3459470 0.9350718 0.6221470 -0.9059384 1.0000000

Verification of [Link] Using Manual Row Selection in Correlation Calculations

identical(cor(mat, use="[Link]"), cor(mat[c(1,4,7,9,10), ]))

[1] TRUE
Explanation:
cor(mat, use="[Link]")

Automatically selects only the rows that have no NA values across all columns.

In your matrix, those rows are 1, 4, 7, 9, 10.

cor(mat[c(1,4,7,9,10), ])

Here, you manually select the same rows that are complete.

The correlation is then calculated on exactly the same data.

identical() - Checks if the two correlation matrices are exactly the same (same values, same
dimensions).

Returns TRUE → confirms that [Link] is equivalent to manually using only complete
rows.
6. [Link]

11
 • Instead of removing entire rows with any NA (like "[Link]"), each pair of
variables uses all rows where both values are available.
• This allows maximum usage of the available data and avoids discarding rows
unnecessarily.
• As a result, the correlation matrix can be more "complete" even if individual rows have
NAs.

cor(mat, use="[Link]")
Output:
m n p q r
m 1.00000000 -0.02511812 -0.3965859 0.4622943 -0.2001722
n -0.02511812 1.00000000 0.8717389 -0.5070416 0.5332259
p -0.39658588 0.87173889 1.0000000 -0.5197292 0.1312506
q 0.46229434 -0.50704163 -0.5197292 1.0000000 -0.4242958
r -0.20017222 0.53322585 0.1312506 -0.4242958 1.0000000

Correlation Between Specific Columns Using [Link] in R

Code Example
cor(mat[, c("m", "n")], use = "[Link]")
Explanation
1. Selecting Columns
o mat[, c("m", "n")] selects only the m and n columns from the matrix, producing
a smaller matrix with 2 columns and 10 rows.
2. Handling Missing Values with [Link]
o "[Link]" keeps only rows where both m and n have non-missing
values.
o In this example, rows 1, 4, 6, 7, 8, 9, 10 are complete. Rows with at least one
NA (2, 3, 5) are removed.
3. Correlation Calculation
o Correlation is computed using only the complete rows:
m n
m 1.00000000 -0.02511812
n -0.02511812 1.00000000
o cor(m, m) and cor(n, n) = 1
o cor(m, n) = -0.0251 → very weak negative correlation based on complete data.

Correlation Between Columns m and p Using [Link] in R

Code
cor(mat[, c("m", "p")], use = "[Link]")

Step 1: Selecting Columns

• mat[, c("m", "p")] selects only the m and p columns from the matrix.
• The resulting data has 2 columns and 10 rows, some containing NA.

Step 2: Handling Missing Values with [Link]

• "[Link]" keeps only rows where both m and p are not NA.
In the dataset:

12
 m = 9, 9, NA, 3, NA, 5, 8, 1, 10, 4
p = 8, 4, 3, 9, 10, NA, 3, NA, 9, 9
• Complete rows (both values non-NA) are: 1, 2, 4, 7, 9, 10
• Rows 3, 5, 6, 8 contain NA in at least one column → removed.

Step 3: Correlation Calculation

• Correlation is computed using only these 6 complete rows.
Result:
m p
m 1.0000000 -0.3965859
p -0.3965859 1.0000000
• cor(m, m) = 1 and cor(p, p) = 1
• cor(m, p) = -0.3966 → moderate negative correlation between m and p based on
available complete data.
Using the tips Dataset and Pairwise Plots
Load the dataset
data(tips, package = "reshape2")
• The tips dataset comes from the reshape2 package.
• It contains information about restaurant tips, including variables like:
o total_bill → total bill amount
o tip → tip amount
o sex → sex of the customer
o smoker → smoker status
o day → day of the week
o time → lunch/dinner
o size → number of people in the party

View the first few rows

head(tips)
Example output:
total_bill tip sex smoker day time size
16.99 1.01 Female No Sun Dinner 2
10.34 1.66 Male No Sun Dinner 3
21.01 3.50 Male No Sun Dinner 3
23.68 3.31 Male No Sun Dinner 2
24.59 3.61 Female No Sun Dinner 4
25.29 4.71 Male No Sun Dinner 4
• head() shows the first 6 rows.
• Useful to quickly inspect the structure of the data.

Pairwise plots using GGally::ggpairs

GGally::ggpairs(tips)
• GGally::ggpairs() creates a pairwise plot matrix:
o Each variable is plotted against every other variable.
o Diagonal: usually histograms or density plots of each variable.
o Lower triangle: scatterplots showing relationships between numeric variables.
o Upper triangle: often shows correlation coefficients or smoothed plots.

13
 o Categorical variables are handled automatically (e.g., colored points,
boxplots).
• Helps visually explore correlations and distributions for all variables in one plot.

Using RXKCD to Fetch XKCD Comics

Load the Package
library(RXKCD)
• RXKCD is an R package for working with XKCD comics.
• XKCD is a popular webcomic focused on science, math, and technology.
• The package allows you to:
o Fetch comic data
o Download images
o Access comic metadata (title, alt text, date, etc.)
Fetch a Specific Comic
getXKCD(WHICH="552")
• getXKCD() is the main function to retrieve an XKCD comic.
• WHICH specifies the comic number:
o "552" → fetches comic number 552.
• Output of getXKCD() is usually a list containing:
o num → comic number
o title → comic title
o img → URL of the comic image
o alt → the alt-text that appears when you hover over the image
o day, month, year → date of publication

Example output:
$img
[1] "[Link]

$title

14
[1] "Pressure"
$alt
[1] "Not to be confused with 'force per unit area.'"
$num
[1] 552
$day
[1] 5
$month
[1] 9
$year
[1] 2009

Note:
To get the details, visit the page
• Go to [Link]
• Each comic URL contains the comic number:
Example: [Link]

18.3 T-Tests in Statistics

Definition
A t-test is a statistical test used to compare means of one or two groups to determine if they
are significantly different from each other. It is based on the Student’s t-distribution, which
is used when the sample size is small and population standard deviation is unknown.
Types of T-Tests
1. One-Sample T-Test
o Compares the mean of a single sample to a known value (e.g., population
mean).
Formula:

Example in R:
x <- c(12, 14, 15, 13, 16)
[Link](x, mu=14) # test if mean is 14

2. Two-Sample T-Test
o Compares the means of two independent groups.
o Formula:

15
 o

Example in R:
group1 <- c(12, 14, 15)
group2 <- c(10, 11, 13)
[Link](group1, group2) # independent samples

3. Paired T-Test
o Compares means of two related groups (e.g., before and after treatment).
o Formula:

Example in R:
before <- c(85, 90, 88)
after <- c(88, 92, 89)
[Link](before, after, paired=TRUE)

Steps to Perform a T-Test

1. State the Hypotheses:
o Null Hypothesis (H0H_0H0): No difference in means
o Alternative Hypothesis (H1H_1H1): Means are different
2. Choose the Significance Level (α\alphaα)
o Commonly α=0.05\alpha = 0.05α=0.05
3. Calculate the Test Statistic
o Use the formula for the appropriate t-test
4. Determine the Critical Value / P-value
o Compare t-value with critical t from tables, or use p-value
5. Make Decision
o If p<αp < \alphap<α → reject H0H_0H0
o If p≥αp \ge \alphap≥α → fail to reject H0H_0H0

16
Assumptions of T-Tests
1. Data are continuous (interval/ratio).
2. Data are approximately normally distributed.
3. Observations are independent.
4. For two-sample t-tests, population variances are equal (can be relaxed with [Link]
= FALSE in R).

Advantages
• Simple and widely used.
• Can handle small sample sizes.
• Helps compare sample mean with population or another sample.

Disadvantages
• Sensitive to non-normal data in small samples.
• Assumes independent observations (except for paired t-test).
• Not suitable for categorical data.

Examples:
[Link](tips$tip, alternative = "[Link]", mu = 2.50)
This is performing a one-sample t-test. Here's what each part means:
[Link]()
This is the built-in R function for performing a t-test, which is used to compare means.
Depending on the input, it can do:
• One-sample t-test
• Two-sample t-test
• Paired t-test

tips$tip
This is the data vector we are testing. In this case:
tips is a dataset (probably the famous “tips” dataset from restaurants).
tips$tip refers to the column containing the tip amounts.

alternative = "[Link]"
This specifies the alternative hypothesis:
"[Link]" → Tests whether the mean of the sample is not equal to the hypothesized mean.
"less" → Tests if the mean is less than the hypothesized mean.
"greater" → Tests if the mean is greater than the hypothesized mean.

Since you chose "[Link]", you are checking:

mu = 2.50

What the test returns

17
head(tips)
total_bill tip sex smoker day time size
1 16.99 1.01 Female No Sun Dinner 2
2 10.34 1.66 Male No Sun Dinner 3
3 21.01 3.50 Male No Sun Dinner 3
4 23.68 3.31 Male No Sun Dinner 2
5 24.59 3.61 Female No Sun Dinner 4
6 25.29 4.71 Male No Sun Dinner 4

> unique(tips$sex)
[1] Female Male
Levels: Female Male
> unique(tips$day)
[1] Sun Sat Thur Fri
Levels: Fri Sat Sun Thur
> [Link](tips$tip,alternative="[Link]",mu=2.50)

One Sample t-test

data: tips$tip
t = 5.6253, df = 243, p-value = 5.08e-08
alternative hypothesis: true mean is not equal to 2.5
95 percent confidence interval:
2.823799 3.172758
sample estimates:
mean of x
2.998279

randT=rt(30000,df=NROW(tips)-1)

tipttest=[Link](tips$tip,alternative="[Link]",mu=2.50)

18
ggplot([Link](x=randT))+geom_density(aes(x=x),fill="grey",color="grey")+geom_vline(
xintercept=tipttest$statistic)+geom_vline(xintercept=mean(randT)+c(-
2,2)*sd(randT),linetype=2)

t-distribution and t-statistic for tip data. The dashed lines are two standard deviations from
mean in either directions. The thick black line, the t-statistic, is so far outside the distribution
that we must reject the null hypothesis and conclude that the true mean is not $2.50.

Explanation
1. Generate a random t-distribution
randT = rt(30000, df = NROW(tips) - 1)
• rt() generates random numbers from a t-distribution.
• 30000 → we are simulating 30,000 random t-values (a large sample to approximate
the distribution).
• df = NROW(tips)-1 → degrees of freedom is the number of rows in tips dataset minus
1.
o If tips has 244 rows, then df = 243.
This gives us a "null distribution" of t-values under the assumption that the null hypothesis is
true.

2. Perform a one-sample t-test

tipttest = [Link](tips$tip, alternative = "[Link]", mu = 2.50)
• Performs a two-sided one-sample t-test on the tip column in the tips dataset.
• Null hypothesis (H₀): The true mean tip = 2.50
• Alternative hypothesis (H₁): The true mean tip ≠ 2.50
• The result (tipttest) contains:
o tipttest$statistic → observed t-value from your data.
o tipttest$[Link], confidence interval, etc.

3. Plot the t-distribution + observed t-value

19
ggplot([Link](x=randT)) +
geom_density(aes(x=x), fill="grey", color="grey") +
geom_vline(xintercept=tipttest$statistic) +
geom_vline(xintercept=mean(randT)+c(-2,2)*sd(randT), linetype=2)
Explanation of layers:
• geom_density(...)
→ Plots the simulated t-distribution (the null distribution of t-values).
• geom_vline(xintercept = tipttest$statistic)
→ Adds a vertical line at the observed t-value from the real data.
This shows where your observed test statistic lies relative to the null distribution.
• geom_vline(xintercept = mean(randT) + c(-2,2)*sd(randT), linetype=2)
→ Adds two dashed vertical lines at:
mean of null distribution±2×standard deviation\text{mean of null distribution} \pm 2 \times
\text{standard deviation}mean of null distribution±2×standard deviation
This roughly corresponds to the ±2 standard deviation region (like a 95% interval under
normal approximation).

To conclude:
This code simulates a t-distribution under the null hypothesis, runs a one-sample t-test on
tips, and plots:
• The null distribution of t-values.
• The observed t-statistic from your data.
• Approximate 95% cutoff lines.
That way, you can visually check whether your observed t-value is far into the tails (→
small p-value) or close to the null distribution center.

We conduct a one-sided t-test to see if the mean is greater than 2.50.

[Link](tips$tip,alternative="greater",mu=2.50)

One Sample t-test

data: tips$tip
t = 5.6253, df = 243, p-value = 2.54e-08
alternative hypothesis: true mean is greater than 2.5
95 percent confidence interval:
2.852023 Inf
sample estimates:
mean of x
2.998279

20
Here’s a visual explanation of your one-sample t-test:
• The blue dashed line represents your calculated t-value (5.6253).
• The green dashed line is the critical t-value for a 95% confidence, one-tailed test.
• The red shaded area is the rejection region where you would reject the null hypothesis.
Since the blue line (t-value) is far to the right of the critical value, it clearly falls in the rejection
region, confirming that the mean tip is significantly greater than 2.50.

18.3.2 Two-Sample T-Test

The Shapiro-Wilk test is a statistical test used to check whether a dataset is normally
distributed. In other words, it tests if your data roughly follows a bell-shaped curve.
Here’s a detailed breakdown:

1. Purpose
• Many statistical tests (like t-tests, ANOVA) assume normality of the data.
• The Shapiro-Wilk test helps you check if this assumption is valid.

2. Hypotheses
Hypothesis Meaning
H₀ (null) The data is normally distributed
H₁ (alternative) The data is not normally distributed

3. Test Statistic
• Denoted as W.
• Ranges between 0 and 1.
o W close to 1: Data is close to normal.
o W far from 1: Data deviates from normality.

21
4. p-value
• If p-value > 0.05 → Fail to reject H₀ → Data is approximately normal.
• If p-value ≤ 0.05 → Reject H₀ → Data is not normal.

5. R Example
[Link](tips$tip)
• Suppose it returns W = 0.89781, p-value = 8.2e-12
• Interpretation:
o W = 0.89781 → Some deviation from normality
o p-value = 8.2e-12 → Extremely small → Reject H₀
o Conclusion: Tip data is not normally distributed

aggregate(tip~sex,data=tips,var)
sex tip
1 Female 1.344428
2 Male 2.217424

This calculates the variance of tips for each sex.

Interpretation:
• Variance of tips for Females ≈ 1.34
• Variance of tips for Males ≈ 2.22
The male tips are more spread out (higher variability) than female tips.

> [Link](tips$tip[tips$sex=="Female"])

Shapiro-Wilk normality test

data: tips$tip[tips$sex == "Female"]

W = 0.95678, p-value = 0.005448
• Null hypothesis (H₀): The data is normally distributed.
• Alternative hypothesis (H₁): The data is not normally distributed.
• p-value, which is extremely small (< 0.05).
Interpretation:
• We reject H₀. The tip data is not normally distributed.

> [Link](tips$tip[tips$sex=="Male"])

Shapiro-Wilk normality test

data: tips$tip[tips$sex == "Male"]

W = 0.87587, p-value = 3.708e-10

Meaning:
• Null hypothesis (H₀): Male customers’ tips are normally distributed.
• Alternative hypothesis (H₁): Male customers’ tips are not normally distributed.
• Test statistic (W) = 0.87587 → Far from 1, indicates departure from normality.
• p-value = 3.708e-10 (< 0.05) → Reject H₀.
Conclusion: The distribution of tips given by male customers is not normal.

22
> ggplot(tips,aes(x=tip,fill=sex))+geom_histogram(binwidth=.5,alpha=1/2)

What this does:

• aes(x=tip, fill=sex) → plots tips on the x-axis, colors bars by sex (Male vs Female).
• geom_histogram(binwidth=0.5) → groups tips into bins of width 0.5.
• alpha=1/2 → makes the bars semi-transparent (so overlapping areas can be seen).
Interpretation:
• You’ll see two overlapping histograms (Male and Female).
• The shape of the Male histogram should confirm the Shapiro result — not perfectly
bell-shaped, possibly skewed.
• Comparing Male vs Female:
o Male tips tend to be more spread out (higher variance, as you saw earlier).
o Female tips are usually more concentrated.

Since the data do not appear to be normally distributed, standard F-test nor bartlett test will
suffice. Hence Ansari-Bradley test is used to examine the equality of variances.

[Link](tip~sex,tips)

Ansari-Bradley test

data: tip by sex

AB = 5582.5, p-value = 0.376
alternative hypothesis: true ratio of scales is not equal to 1

This test indicates that the variances are equal, that is we can use the standard two sample t-
test.
[Link](tip~sex,data=tips,[Link]=TRUE)

23
 Two Sample t-test

data: tip by sex

t = -1.3879, df = 242, p-value = 0.1665
alternative hypothesis: true difference in means between group Female and group Male is not
equal to 0
95 percent confidence interval:
-0.6197558 0.1074167
sample estimates:
mean in group Female mean in group Male
2.833448 3.089618

Based on the test, the results are nor significant and we conclude that male and female
dinners tip roughly equal.
To check if the two means are within the two standard deviations of each other.

[Link]("plyr")
library(plyr)

tipSummary=ddply(tips,"sex",summarize,[Link]=mean(tip),[Link]=sd(tip),Lower=[Link]
-2*[Link]/sqrt(NROW(tip)),Upper=[Link]+2*[Link]/sqrt(NROW(tip)))
tipSummary

This line is doing group-wise summary statistics with ddply() from the plyr package.
Step-by-Step Explanation
1. ddply(tips, "sex", summarize, ...)
• ddply() splits the data frame tips by the factor sex ("Male" and "Female").
• Then applies the functions inside summarize to each group.
• Finally returns a new data frame with results for each group.

2. [Link] = mean(tip)
• Calculates the average tip separately for males and females.

3. [Link] = sd(tip)
• Computes the standard deviation of tips for each sex.

4. Lower and Upper

These calculate an approximate 95% confidence interval for the mean tip.
So:
• Lower = lower bound of CI
• Upper = upper bound of CI

sex [Link] [Link] Lower Upper

1 Female 2.833448 1.159495 2.584827 3.082070
2 Male 3.089618 1.489102 2.851931 3.327304
Interpretation
• Female average tip ≈ 2.8 with a tighter spread.
• Male average tip ≈ 3.1 with larger variability.

24
 • The confidence intervals show the range where the true mean tip is likely to fall for each group.
• If the intervals overlap → difference may not be statistically significant.
• If they don’t overlap → stronger evidence of a real difference.

>ggplot(tipSummary,aes(x=[Link],y=sex))+geom_point()+geom_errorbarh(aes(xmin=Low
er,xmax=Upper),height=.2)

18.3.3 Paired Two-Sample T-Test

• A Paired t-test compares the means of two related groups (not independent).
• Instead of comparing the groups directly, it looks at the differences within each pair.
We use a Paired t-test when the two samples are dependent / matched:
Examples:
• Before vs After measurements on the same subjects (e.g., blood pressure before and
after treatment).
• Left vs Right measurements from the same person (e.g., left eye vs right eye).
• Father vs Son heights (each father is naturally paired with his son).
• Two methods tested on the same participants (e.g., exam scores of students under two
teaching methods).
In all these cases, observations are linked one-to-one.

[Link]("UsingR")
library("UsingR")
data([Link],package="UsingR")
head([Link])

fheight sheight
1 65.04851 59.77827
2 63.25094 63.21404
3 64.95532 63.34242
4 65.75250 62.79238
5 61.13723 64.28113
6 63.02254 64.24221

25
[Link]([Link]$fheight,[Link]$sheight,paired=TRUE)
• Performs a paired t-test between:
o [Link]$fheight → fathers’ heights
o [Link]$sheight → sons’ heights
• paired=TRUE means we are comparing matched pairs (each father with his own
son).

Logic of a Paired t-test

• A paired t-test checks whether the mean difference between paired observations is
significantly different from 0.
Hypotheses:
• H₀ (null): Mean difference = 0
(fathers and sons have the same average height)
• H₁ (alternative): Mean difference ≠ 0
(there is a difference in average height)

Paired t-test

data: [Link]$fheight and [Link]$sheight

t = -11.789, df = 1077, p-value < 2.2e-16
alternative hypothesis: true mean difference is not equal to 0
95 percent confidence interval:
-1.1629160 -0.8310296
sample estimates:
mean difference
-0.9969728
Interpretation
1. Mean Difference
• mean difference = -0.997
• This is dˉ=average(fheight−sheight)\bar{d} = \text{average}(fheight -
sheight)dˉ=average(fheight−sheight).
• Negative → On average, sons are about 1 inch taller than their fathers.

2. t-value
• t = -11.789
• This is very large in magnitude → the difference is highly significant.
• It means the observed difference is almost 12 standard errors away from 0.

3. Degrees of Freedom
• df = 1077 → 1,078 pairs (father-son).
• df=n−1df = n-1df=n−1, where nnn = number of pairs.

4. p-value
• p-value < 2.2e-16 (essentially 0).
• Much smaller than 0.05 → Reject H₀.
• Strong evidence of a real difference in heights between fathers and sons.

5. Confidence Interval
• 95% CI: [−1.163,−0.831][-1.163, -0.831][−1.163,−0.831]

26
 • Since the interval does not include 0, the difference is statistically significant.
• Interpretation: With 95% confidence, sons are between 0.83 and 1.16 inches taller
than fathers.

The paired t-test shows that sons are significantly taller than their fathers, by about 1 inch
on average.

heightdiff=[Link]$fheight- [Link]$sheight
ggplot([Link],aes(x=fheight-sheight)) +
geom_density() +
geom_vline(xintercept=mean(heightdiff)) +
geom_vline(xintercept=mean(heightdiff) + 2*c(-1,1)*sd(heightdiff)/sqrt(nrow([Link])),
linetype=2)

Interpretation
heightdiff=[Link]$fheight- [Link]$sheight
• This takes the difference between father’s height (fheight) and son’s height
(sheight) for each pair.
• So heightdiff = vector of differences.
• Positive → father taller than son.
• Negative → son taller than father.
ggplot creates a plot using [Link] dataset.
• aes(x = fheight - sheight) → the x-axis is the difference in heights.
• geom_density() → draws a smooth density curve (like a smoothed
histogram) showing the distribution of these differences.
geom_vline(xintercept = mean(heightdiff))
• mean(heightdiff) → average difference in height.
• A solid vertical line is drawn at this value.
• Shows whether, on average, fathers are taller or shorter than sons.
geom_vline( xintercept = mean(heightdiff) + 2 * c(-1, 1) * sd(heightdiff) /
sqrt(nrow([Link])), linetype = 2 )
• sd(heightdiff) → standard deviation of differences.
• sqrt(nrow([Link])) → square root of sample size (for standard error).
• 2 * ... → multiplying by 2 gives approximately a 95% confidence interval for
the mean (by normal approximation).
• c(-1,1) → gives both sides (mean - margin, mean + margin).
• linetype = 2 → dashed vertical lines.

So these two dashed lines mark the 95% confidence interval for the mean difference in
height.

27
Interpretation of the Plot
• The density curve shows the distribution of (father’s height − son’s height).
• The solid line marks the average difference.
• The two dashed lines show the 95% confidence interval of the mean difference.
• If the CI does not include 0, it suggests fathers and sons have a statistically significant
difference in average height.

18.4 ANOVA
ANOVA (Analysis of Variance) is a statistical test used to compare the means of
three or more groups to see if at least one group mean is significantly different.
It is an extension of the t-test:
• t-test → compares 2 groups.
• ANOVA → compares 3 or more groups.
Hypotheses in ANOVA
• H₀ (null): All group means are equal.
• H₁ (alternative): At least one group mean is different.

Types of ANOVA
1. One-way ANOVA → One independent variable (factor) with ≥ 3 groups.
2. Two-way ANOVA → Two independent variables (factors), may include interaction
effects.
3. Repeated Measures ANOVA → Same subjects measured under different conditions.

ANOVA in R
1. One-way ANOVA Example
Suppose we check if tip amount differs by day in the tips dataset:
model <- aov(tip ~ day-1,tips)
summary(model)

Df Sum Sq Mean Sq F value Pr(>F)

day 4 2203.0 550.8 290.1 <2e-16 ***
Residuals 240 455.7 1.9

28
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

tipintercept=aov(tip~day,tips)
model$coefficients

Here:
• tip ~ day means you are modeling tip as a function of the categorical variable day.
• In R, categorical variables (factors) are converted into dummy variables in the model
matrix.

dayFri daySat daySun dayThur

2.734737 2.993103 3.255132 2.771452
• This looks like the mean tip for each day directly.
• If you calculate mean(tip[tips$day == "Fri"]) you’ll get 2.734737, etc.
• This is just the group means, not regression coefficients.

tipintercept$coefficients
(Intercept) daySat daySun dayThur
2.73473684 0.25836661 0.52039474 0.03671477

Here’s what’s happening:

1. (Intercept) is the mean of the baseline category. In R, the first level of a factor is the
baseline. By default, day is alphabetically ordered, so Fri is baseline: (Intercept) =
mean(tip on Fri) = 2.734737.
2. daySat = 0.25836661 means the difference between Sat and Fri:
mean(tip on Sat)−mean(tip on Fri)=2.993103−2.734737≈0.258367

3. Similarly:
o daySun = 0.52039474 → difference Sun vs Fri
o dayThur = 0.03671477 → difference Thur vs Fri
So aov() stores the regression coefficients relative to the baseline, not the raw means.

tipsbyday=ddply(tips,"day",plyr::summarize,[Link]=mean(tip),[Link]=sd(tip),Length=NR
OW(tip),tfrac=qt(p=.90,df=Length-1),Lower=[Link] -
tfrac*[Link]/sqrt(Length),Upper=[Link]+tfrac*[Link]/sqrt(Length))
Step-by-step explanation
1. ddply(tips, "day", ...)
o ddply splits the tips data frame by the variable day.
o For each day (Fri, Sat, Sun, Thur), it will compute the summaries you define.
2. plyr::summarize
o This tells ddply to create a summary data frame with new columns.
3. [Link] = mean(tip)
o Computes the average tip for that day.
4. [Link] = sd(tip)
o Computes the standard deviation of tips for that day.

29
 5. Length = NROW(tip)
o Length stores the number of observations (tips) for that day.
6. tfrac = qt(p = .90, df = Length - 1)
o qt() gives the t-value for the 90th percentile (one-sided) of the t-distribution
with df = Length - 1 degrees of freedom.
o This is used to calculate a 90% confidence interval.
7. Lower = [Link] - tfrac * [Link] / sqrt(Length)
o The lower bound of the 90% confidence interval.
8. Upper = [Link] + tfrac * [Link] / sqrt(Length)
o The upper bound of the 90% confidence interval.

summary(tipsbyday)

day [Link] [Link] Length tfrac Lower

Fri :1 Min. :2.735 Min. :1.020 Min. :19.00 Min. :1.291 Min. :2.424
Sat :1 1st Qu.:2.762 1st Qu.:1.181 1st Qu.:51.25 1st Qu.:1.293 1st Qu.:2.531
Sun :1 Median :2.882 Median :1.238 Median :69.00 Median :1.294 Median :2.667
Thur:1 Mean :2.939 Mean :1.281 Mean :61.00 Mean :1.303 Mean :2.708
3rd Qu.:3.059 3rd Qu.:1.338 3rd Qu.:78.75 3rd Qu.:1.304 3rd Qu.:2.843
Max. :3.255 Max. :1.631 Max. :87.00 Max. :1.330 Max. :3.072
Upper
Min. :2.976
1st Qu.:3.028
Median :3.132
Mean :3.170
3rd Qu.:3.274
Max. :3.438

ggplot(tipsbyday,aes(x=[Link],y=day))+geom_point() +
geom_errorbarh(aes(xmin=Lower,xmax=Upper),height=.3)

30
1. nrow(tips)
• tips is a data frame (like a table).
• nrow() returns the number of rows in a data frame or matrix.
• Output:
[1] 244
The tips dataset has 244 rows.

2. nrow(tips$tip)
• tips$tip extracts just the tip column, which is a vector of length 244.
• nrow() works only on 2D objects (matrices, data frames).
• Since a vector has no rows/columns, nrow() returns:
NULL

3. NROW(tips$tip)
• NROW() is more general than nrow().
• It works on vectors, matrices, or data frames.
• For a vector, NROW(x) = length(x).
• Output:
[1] 244
Because the tip vector has 244 elements.

Key Difference
• nrow() → Only for 2D objects (data frames, matrices). Returns NULL for vectors.
• NROW() → More flexible. Works for vectors too, returns their length.

Conclusion
• ANOVA is a powerful tool for comparing means of multiple groups.
• It generalizes the t-test and uses the F-statistic.

31