###############################
##                           ##
## CHAPTER 4: SINGLE SAMPLES ##
##                           ##
###############################

###########################################################################

##  AN EXAMPLE OF  ##
##  SINGLE-SAMPLE  ##
##       DATA      ##

mydata <- read.table(file.choose(), header = T, sep = ",") # load the data in from the Excel CSV. file
                                                           # 'extinction_data' and name it 'mydata'
View(mydata) # view the data

values <- mydata$percent # give the data to plot a convenient name

###########################################################################

##              ##
##  HISTOGRAMS  ##
##              ##

### Simple histogram

hist(values, # plot the data as a histogram
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include a title

## Playing with bin number

play <- read.table(file.choose(), header = T, sep = ",") # load the data in from the Excel csv. file
                                                         # 'bin_play' and name it 'play'
View(play) # view the data

numbers <- play$numbers # give the data to plot a convenient name

length(numbers) # check the sample size
                # n=36

1+3.322*(log10(36)) # calculate bins with Sturges' rule
                    # 6.17 bins

2*(36^(1/3)) # calculate bins with Rice rule
             # 6.60 bins

bins<-c(6.17,6.60) # create a list of the two suggestions
mean(bins) # calculate average 
           # 6.385 bins, rounds to 6

hist(numbers, # plot the data as a histogram
     breaks=6, # specify 6 bins
     main="6 bins", # give plot a title
     xlab="Numbers",ylab="Frequency") # set axis labels

par(mfrow=c(1,3)) # arrange 3 plots in 1 row with 3 columns
hist(numbers, breaks=2, main="2 bins", xlab="Numbers",ylab="Frequency") # histogram with 2 bins specified
hist(numbers, breaks=6, main="6 bins", xlab="Numbers",ylab="Frequency") # histogram with 6 bins specified
hist(numbers, breaks=20, main="20 bins", xlab="Numbers",ylab="Frequency") # histogram with 20 bins specified
par(mfrow=c(1,1)) # draw future plots in 1 row with 1 column

### Simple histogram with specified bins

length(values) # check the sample size
               # n=130

1+3.322*(log10(130)) # calculate bins with Sturges' rule
                     # 8.02 bins

2*(130^(1/3)) # calculate bins with Rice rule
              # 10.13 bins

bins<-c(8.02,10.13) # create a list of the two suggestions
mean(bins) # calculate average
           # 9.08 bins, rounds to 9

hist(values, breaks=6, # plot a histogram suggesting 6 bins 
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include a title

hist(values, breaks=c(0,10,20,30,40,50,60), # plot a histogram forcing 6 bins
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include a title

### Extreme values and histograms

nonzero <- subset(mydata, percent!="0") # create a new data set that excludes studies 
                                        # reporting 0% extinction, using the subset function 
                                        # (where ! means 'does not')
View(nonzero) # view the new data

newvalues <- nonzero$percent # give the data to plot a convenient name

hist(newvalues, # plot the data as a histogram
     breaks=12, # specify 12 bins
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include title

## BIGGER PICTURE 4.1: DEALING WITH A BROAD SPREAD OF VALUES ##

install.packages("MASS") # install package 'MASS'
library(MASS) # activate package
data<-mammals # use 'mammals' as our data
View(data) # view the data - brain and body weights of 62 spp. land mammals

par(mfrow=c(1,2)) # arrange 2 plots in 1 row with 2 columns
hist(data$body,
     xlab="Body weight (kg)",
     col="firebrick3",
     main=NULL)
mtext("a.", side=3,line=1,adj=-0.2,cex=1.5,font=2)

logbody<-log10(data$body)
hist(logbody,
     xlab=expression(paste("Log"[10]," body weight (kg)")),
     col="firebrick3",
     main=NULL)
mtext("b.", side=3,line=1,adj=-0.2,cex=1.5,font=2)

par(mfrow=c(1,1)) # draw future plots in 1 row with 1 column

### Refined histogram: customising axes, adding grid lines and colours

hist(values, # plot the data as a histogram
     breaks=c(-5,0,5,10,15,20,25,30,35,40,45,50,55,60), # force specific breaks to
                                                        # create an inital cell of 0 values
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include a title

hist(values, # plot the data as a histogram
     breaks=c(seq(-5,60,5)), # force specific breaks 
     xaxt="n", # do not draw on the x-axis
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include a title

xlabels <- c("0","0-5","5-10","10-15","15-20","20-25","25-30","30-35",
             "35-40","40-45","45-50","50-55","55-60") # create list of bin ranges
axis(side=1, # draw an x-axis at the bottom of the plot
     at=c(seq(-2.5,57.5,5)), # specify spacing of tick marks
     labels=xlabels, # label tick marks with list 'xlabels'
     tick=FALSE, # do not draw tick marks or axis line
     las=2, # draw labels perpendicular to x-axis
     line=-1) # position labels close to cells

hist(values, breaks=c(seq(-5,60,5)), # plot histogram with specified bins
     xaxt="n",yaxt="n", # do not draw axes
     ylab= "",xlab="", # do not include axis labels
     main=NULL) # do not include a title
abline(h=(seq(0,30,1)), col="lightgray") # add minor grid lines
abline(h=(seq(0,30,5)), col="darkgray") # add major grid lines
par(new=TRUE) # prepare to draw on top of existing plot
hist(values, # plot data as histogram
     col=c("mediumpurple4",rep("seagreen2",12)), # list colours of cells
     breaks=c(seq(-5,60,5)), # specify bins
     xaxt="n", # do not draw x-axis
     ylab= "Frequency",xlab="Predicted extinction rate %", # set axis labels
     main=NULL) # do not include a title
xlabels <- c("0","0-5","5-10","10-15","15-20","20-25","25-30","30-35",
             "35-40","40-45","45-50","50-55","55-60") # create list of bin ranges
axis(side=1, at=c(seq(-2.5,57.5,5)), labels=xlabels,tick=FALSE,las=2,line=-1) # draw on x-axis

###########################################################################

##            ##
##  BOXPLOTS  ##
##            ##

### Simple boxplot

boxplot(values, # plot data as boxplot
        ylab= "Predicted extinction rate %") # set y-axis label

### Refined boxplot: editing axes, adding colour and text

boxplot(values, # plot data as boxplot
        ylab= "Predicted extinction rate %", # set y-axis label
        ylim=c(-5,70), # define the bounds of the y-axis
        yaxs = "i", # fit neat tick labels based on ylim
        col="cadetblue1") # set box colour

summary(values) # check descriptive statistics

text(0.65,58.8, labels="outlier = 58.8") # add text to plot at set coordinates
text(0.65,6.7, labels="median = 6.7")
text(0.65,18, labels="IQR = 0.5 - 22.2")

abline(h=32,col="red",lty=2,lwd=2) # add a reference line

text(1.35,38,labels="test value",col="red") # label reference line