library(ggplot2)
theme_set(theme_bw())
setwd("/home/duncan/Dropbox/Public/Analyses/Mariana_bombus_iucn")
d<-read.csv("/home/duncan/Dropbox/Public/Analyses/Mariana_bombus_iucn/EOO_AOO_simplified.csv",sep=";")
str(d)
## 'data.frame':    16 obs. of  9 variables:
##  $ taxonid            : int  13152906 13340289 13340333 13340348 13340369 13340374 13340410 13340462 13342541 13356497 ...
##  $ assessmentid       : int  57047934 57345436 57346148 56956517 57347466 57347749 57349231 57349805 57366703 57368180 ...
##  $ region_name        : Factor w/ 1 level "Europe": 1 1 1 1 1 1 1 1 1 1 ...
##  $ friendly_name      : Factor w/ 16 levels "Bombus alpinus",..: 1 3 4 6 7 8 9 10 12 13 ...
##  $ taxonomic_authority: Factor w/ 15 levels "Curtis, 1835",..: 4 15 11 5 7 6 12 14 4 1 ...
##  $ category           : Factor w/ 3 levels "CR","EN","VU": 3 2 3 3 2 3 3 2 3 3 ...
##  $ criteria           : Factor w/ 10 levels "A2a","A2c","A2c+3c+4c",..: 8 7 3 2 3 6 10 6 2 9 ...
##  $ EOO                : int  3711722 14369 4910210 7921424 2128849 914838 962757 230166 12504130 1493904 ...
##  $ AOO                : int  1288 156 4336 8864 944 1432 420 368 17296 700 ...
d<-subset(d,category!="CR")

Plot EOO

g0<-ggplot(d,aes(x=category,y=EOO))
g1<-g0+geom_boxplot()+ stat_summary(fun.y = mean, geom = "point")
g1<-g1 + stat_summary(fun.data = mean_cl_boot, geom = "errorbar",col="red")
g1

plot of chunk unnamed-chunk-2

Plot AOO

g0<-ggplot(d,aes(x=category,y=AOO))
g1<-g0+geom_boxplot()+ stat_summary(fun.y = mean, geom = "point")
g1<-g1 + stat_summary(fun.data = mean_cl_boot, geom = "errorbar",col="red")
g1

plot of chunk unnamed-chunk-3

Logs for AOO

d$log10AOO<-log10(d$AOO)

g0<-ggplot(d,aes(x=category,y=log10AOO))
g1<-g0+geom_boxplot()+ stat_summary(fun.y = mean, geom = "point")
g1<-g1 + stat_summary(fun.data = mean_cl_boot, geom = "errorbar",col="red")
g1

plot of chunk unnamed-chunk-4

EOO boot

library(boot)
## 
## Attaching package: 'boot'
## 
## The following object is masked from 'package:survival':
## 
##     aml
## 
## The following object is masked from 'package:lattice':
## 
##     melanoma
## Calculate the differences in means

dif <- function(data, indices) {
  d <- data[indices,] #
  EOO_EN<-d$EOO[d$category=="EN"]
  EOO_VU<-d$EOO[d$category=="VU"]
  dif<-mean(EOO_EN)-mean(EOO_VU)
  return(dif)
}
# bootstrapping with 1000 replications
results <- boot(data=d, statistic=dif,
   R=1000)

## Confidence intervals. If they go from negative to positive (include zero) not significant

boot.ci(results)
## Warning: bootstrap variances needed for studentized intervals
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = results)
## 
## Intervals : 
## Level      Normal              Basic         
## 95%   (-6238213,   476387 )   (-6165392,   448324 )  
## 
## Level     Percentile            BCa          
## 95%   (-6124388,   489328 )   (-6649808,   202808 )  
## Calculations and Intervals on Original Scale

AOO boot

dif <- function(data, indices) {
  d <- data[indices,] 
  AOO_EN<-d$AOO[d$category=="EN"]
  AOO_VU<-d$AOO[d$category=="VU"]
  dif<-mean(AOO_EN)-mean(AOO_VU)
  return(dif)
}
# bootstrapping with 1000 replications
results <- boot(data=d, statistic=dif,
   R=1000)


boot.ci(results)
## Warning: bootstrap variances needed for studentized intervals
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = results)
## 
## Intervals : 
## Level      Normal              Basic         
## 95%   (-8232,  -637 )   (-7641,    89 )  
## 
## Level     Percentile            BCa          
## 95%   ( -8764,  -1033 )   (-10257,  -1649 )  
## Calculations and Intervals on Original Scale
## Some BCa intervals may be unstable