CI type in diffCalc

Compare 95% CIs for \(D_{Jost}\), calculated by bootstraping across loci or individuals.

Calculate pairwise Fst and 95% CIs

library(diveRsity)
data(Test_data)
# bootstrapping individuals
ind <- diffCalc(Test_data, bs_pairwise = TRUE, boots = 999, para = TRUE)
# bootstrapping loci
loc <- diffCalc(Test_data, bs_pairwise = TRUE, boots = 999, para = TRUE, 
                ci_type = "loci")

Extract the relevant data to plot

plot_data <- data.frame(pop = c(as.character(ind$bs_pairwise$D$populations),
                                as.character(loc$bs_pairwise$D$populations)), 
                        actual = c(ind$bs_pairwise$D$actual,
                                   loc$bs_pairwise$D$actual),
                        low = c(ind$bs_pairwise$D$lower,
                                loc$bs_pairwise$D$lower),
                        hi = c(ind$bs_pairwise$D$upper,
                                loc$bs_pairwise$D$upper),
                        type = c(rep("ind", nrow(ind$bs_pairwise$D)),
                                 rep("loc", nrow(loc$bs_pairwise$D))))
pander::pandoc.table(head(plot_data))

## 
## ---------------------------------------------
##      pop        actual    low     hi    type 
## -------------- -------- ------- ------ ------
## pop1, vs pop2,  0.0027  -0.0122 0.0201  ind  
## 
## pop1, vs pop3,  0.1802  0.1492  0.2137  ind  
## 
## pop1, vs pop4,  0.1484  0.1225  0.1764  ind  
## 
## pop1, vs pop5,  0.2527  0.2223  0.2865  ind  
## 
## pop1, vs pop6,  0.1494  0.1202  0.1834  ind  
## 
## pop2, vs pop3,  0.1579   0.129  0.1901  ind  
## ---------------------------------------------

plot the results

library(ggplot2)
dodge <- position_dodge(width=0.7) 
p <- ggplot(data = plot_data, aes(x = pop, y = actual, ymax = max(actual),
                                  colour = type)) + 
  geom_point(position = dodge, size = 3) +
  geom_errorbar(aes(ymin = low, ymax = hi), width = 0, position = dodge,
                lwd = 1) +
  theme_bw() +
  theme(axis.text.x=element_text(angle = -90, size = 15, vjust = 0.3),
        axis.title=element_text(size = 20),
        axis.text.y=element_text(size = 15)) +
  scale_x_discrete(labels = levels(plot_data$pop)) +
  labs(x = "Pairwise Comparison", y = expression("D"["Jost"]),
       colour = "Bootstrap type")
p

Conclusions

All CIs contain the sample estimate. CIs calculated by bootstrapping individuals are tighter than loci CIs for D estimates > 0.05. The oposite is true below D = 0.05 (not sure why yet). One example (pop3 vs pop4) is seen where ind CIs fail to reject the null of panmixia, while loc CI does not.

What about other measures

\(G_{ST}\) (Nei & Chesser, 1983)

plot_data <- data.frame(pop = c(as.character(ind$bs_pairwise$gst$populations),
                                as.character(loc$bs_pairwise$gst$populations)),
                        actual = c(ind$bs_pairwise$gst$actual,
                                   loc$bs_pairwise$gst$actual),
                        low = c(ind$bs_pairwise$gst$lower,
                                loc$bs_pairwise$gst$lower),
                        hi = c(ind$bs_pairwise$gst$upper,
                                loc$bs_pairwise$gst$upper),
                        type = c(rep("ind", nrow(ind$bs_pairwise$gst)),
                                 rep("loc", nrow(loc$bs_pairwise$gst))))
pander::pandoc.table(head(plot_data))

## 
## --------------------------------------------
##      pop        actual   low     hi    type 
## -------------- -------- ------ ------ ------
## pop1, vs pop2,  0.0019  -0.001 0.0051  ind  
## 
## pop1, vs pop3,  0.0341  0.0283 0.0402  ind  
## 
## pop1, vs pop4,  0.0296  0.025  0.035   ind  
## 
## pop1, vs pop5,  0.0447  0.0386 0.0521  ind  
## 
## pop1, vs pop6,  0.0293  0.0244 0.0349  ind  
## 
## pop2, vs pop3,  0.0287  0.0237 0.0348  ind  
## --------------------------------------------

library(ggplot2)
dodge <- position_dodge(width=0.7) 
p <- ggplot(data = plot_data, aes(x = pop, y = actual, ymax = max(actual),
                                  colour = type)) + 
  geom_point(position = dodge, size = 3) +
  geom_errorbar(aes(ymin = low, ymax = hi), width = 0, position = dodge,
                lwd = 1) +
  theme_bw() +
  theme(axis.text.x=element_text(angle = -90, size = 15, vjust = 0.3),
        axis.title=element_text(size = 20),
        axis.text.y=element_text(size = 15)) +
  scale_x_discrete(labels = levels(plot_data$pop)) +
  labs(x = "Pairwise Comparison", y = expression("G"["ST"]),
       colour = "Bootstrap type")
p