16 Relating Ranks (Ch. 16)

Soil <- LETTERS[1:17]
Prod <- c(2, 6, 3, 7, 4, 8, 8, 1, 3, 5, 1, 8, 7, 2, 4, 3, 6)
Density <- c(0.26, 1.35, 0.44, 1.26, 0.35, 2.3, 1.76, 0.31, 0.37, 0.78, 0.04, 
             1.62, 1.34, 0.47, 0.56, 0.48, 0.76)
KPNeolithic <- data.frame(Soil, Prod, Density)

16.1 Calculating Spearman’s Rank Correlation (\(\rho\))

We can use the rank() function to calculate the rankings of productivity and density, and begin to build Table 16.1.

KPNeolithic$ProdRank <- rank(KPNeolithic$Prod, ties.method = "average")
KPNeolithic$DensRank <- rank(KPNeolithic$Density, ties.method = "average")

KPNeolithic$d <- KPNeolithic$ProdRank - KPNeolithic$DensRank
KPNeolithic$d2 <- KPNeolithic$d ^ 2 

kable(KPNeolithic) %>% kableExtra::kable_classic(full_width = F)

Soil	Prod	Density	ProdRank	DensRank	d	d2
A	2	0.26	3.5	2	1.5	2.25
B	6	1.35	11.5	14	-2.5	6.25
C	3	0.44	6.0	6	0.0	0.00
D	7	1.26	13.5	12	1.5	2.25
E	4	0.35	8.5	4	4.5	20.25
F	8	2.30	16.0	17	-1.0	1.00
G	8	1.76	16.0	16	0.0	0.00
H	1	0.31	1.5	3	-1.5	2.25
I	3	0.37	6.0	5	1.0	1.00
J	5	0.78	10.0	11	-1.0	1.00
K	1	0.04	1.5	1	0.5	0.25
L	8	1.62	16.0	15	1.0	1.00
M	7	1.34	13.5	13	0.5	0.25
N	2	0.47	3.5	7	-3.5	12.25
O	4	0.56	8.5	9	-0.5	0.25
P	3	0.48	6.0	8	-2.0	4.00
Q	6	0.76	11.5	10	1.5	2.25

Calculating the values of t and T for soil productivity and settlement density is slightly trickier, since it involves relating multiple columns (e.g., “how many other rows have the same value of ‘Prod’ as this one?”). We can do this with a for() loop that, for each row, returns the total number of matching values in the relevant column.

i <- 1
for (i in 1:nrow(KPNeolithic)){
KPNeolithic$t_Prod[i] <- sum(KPNeolithic$Prod == KPNeolithic$Prod[i])
}

i <- 1 #need to reset i 
for (i in 1:nrow(KPNeolithic)){
KPNeolithic$t_Dens[i] <- sum(KPNeolithic$Density == KPNeolithic$Density[i])
}

Once we have values of t, we can calculate values of T: \(T = \dfrac{t^3-t}{12}\). ### why 12? why t^3###

KPNeolithic$T_Prod <- (KPNeolithic$t_Prod^3 - KPNeolithic$t_Prod) / 12
KPNeolithic$T_Dens <- (KPNeolithic$t_Dens^3 - KPNeolithic$t_Dens) / 12

We need the sums of \(d^2\), ‘T_Prod’, and ‘T_Dens’ in order to calculate Spearman’s rank correlation.

sum_d2 <- sum(KPNeolithic$d2)
sum_T_Prod <- sum(KPNeolithic$T_Prod)
sum_T_Dens <- sum(KPNeolithic$T_Dens)

These are 56.5, 17, and 0, respectively. With those we calculate a sum of squares for each of our variables, using the equation \(\sum{x}^2 = \dfrac{{n}^3 - n}{12} - \sum{T_x}\), where n is the number of samples, and \(\sum{T_x}\) is the sum of the T values for that variable, calculated above.

sum_x2_Prod <- ((nrow(KPNeolithic)^3 - nrow(KPNeolithic)) / 12) - sum_T_Prod

sum_x2_Dens <- ((nrow(KPNeolithic)^3 - nrow(KPNeolithic)) / 12) - sum_T_Dens

With those sums of squares we now have all the elements necessary to calculate Spearman’s rank correlation, using the equation \(r_s = \dfrac{\sum{{x}^2}+\sum{{y}^2}-\sum{{d}^2}}{2\sqrt{\sum{{x}^2}\sum{{y}^2}}}\).

rs <- (sum_x2_Prod + sum_x2_Dens - sum_d2) / (2 * sqrt(sum_x2_Prod * sum_x2_Dens)); rs

## [1] 0.929497

You’ll be pleased to know that R enables you to shortcut the process, calculating Spearman’s rank correlation using cor.test(). The value of \(\rho\) is a measure of correlation that can be interpreted as analogous to Pearson’s r (though the two should not be directly compared).

cor.test(KPNeolithic$Density, KPNeolithic$Prod, alternative="two.sided", method="spearman", exact = F)

## 
##  Spearman's rank correlation rho
## 
## data:  KPNeolithic$Density and KPNeolithic$Prod
## S = 56.929, p-value = 6.362e-08
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##       rho 
## 0.9302338

16.2 Significance

Note that cor.test() will return a p-value, equivalent to what Drennan describes on pp226-228.

16.3 Assumptions and Robust Methods

There are other tests of rank correlation, e.g. Kendall’s tau (\(\tau\)) rank correlation, also accessible using cor.test().

cor.test(KPNeolithic$Density, KPNeolithic$Prod, alternative="two.sided", method="kendall", exact = F)

## 
##  Kendall's rank correlation tau
## 
## data:  KPNeolithic$Density and KPNeolithic$Prod
## z = 4.3712, p-value = 1.236e-05
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.8053132

16.4 Practice

Data for the practice problems can be found in Teixeira.txt.