library(raster); #For raster-based loading, calculations, and mapping
library(rgdal); #For reading the M polygon

setwd("~/Dropbox/ENMSeminar/Labs:Homeworks/Lab5/Lab5Data/WorldClimTerrestrial/");
  #Navigate to folder containing environmental data

#Loading a raster stack
envtList <- list.files(pattern = ".asc"); #Gets a list of .asc files
envtStack <- stack(envtList); #Reads in .asc files as a raster stack
crs(envtStack) <- "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0" #Defines projection of layers
plot(envtStack); #Plots all the layers of the raster stack object

#Get M shapefile
setwd("~/Dropbox/ENMSeminar/Labs:Homeworks/Lab5/Lab5Data/SyrmaticusSoemmerringii/");#Changing working directory to the folder where the shapefile of your M from created in Part One of the lab.
copperPheasant <- readOGR("./Syrmaticus soemmerringii.shp"); #Reads in shapefile; You will change this file to the shapefile you created in part 1 of the lab.
crs(copperPheasant) <- crs(envtStack); #Defines Copper Pheasant M projection as identical to envtStack

#Masking your environmental variables to M
copperPheasantTraining <- crop(envtStack, copperPheasant) #Cropping the raster stack
copperPheasantTraining <- mask(copperPheasantTraining, copperPheasant) #Masking the raster stack
writeRaster(copperPheasantTraining, filename = "SyrmaticusSoemmerringii",
  format = "ascii", bylayer = T, suffix=names(envtStack), NAFlag = "-9999", overwrite = T);
    #Saves all the layers in the stack, with "SyrmaticusSoemmerringii" as a prefix

#Removing correlated variables
pairs(copperPheasantTraining);
  #Shows correlations among variables. Inspect visually. Which layers are highly correlated?
reducedCopperPheasantTraining <- copperPheasantTraining[[-3:-6]];
  #Removes Bio 10 through Bio 13.
pairs(reducedCopperPheasantTraining);
reducedCopperPheasantTraining <- reducedCopperPheasantTraining[[-6:-8]];
  #Removes Bio 17 through Bio 19
pairs(reducedCopperPheasantTraining);
reducedCopperPheasantTraining <- reducedCopperPheasantTraining[[-9:-11]];
  #Removes Bio 5 through Bio 7
pairs(reducedCopperPheasantTraining);
  #Ok, now none of the variables have a correlation coefficient of 0.90. This is a somewhat arbitrary threshold. Use your best judgement to find a logical breakpoint for removing correlated variables from your environmental dataset.
writeRaster(reducedCopperPheasantTraining, filename = "reducedSyrmaticusSoemmerringii",
            format = "ascii", bylayer = T, suffix=names(reducedCopperPheasantTraining), NAFlag = "-9999", overwrite = T);
  #Saved reduced set if you don't want to hunt for the variables you want to get rid of.