Skip to content

Commit

Permalink
Pre-pre-Berkeley meeting
Browse files Browse the repository at this point in the history
  • Loading branch information
gfalbery committed Sep 18, 2019
1 parent 7e97ba2 commit b028535
Show file tree
Hide file tree
Showing 12 changed files with 1,390 additions and 664 deletions.
54 changes: 14 additions & 40 deletions Greg ENM Code/00_Iceberg Species Dropping.R
Original file line number Diff line number Diff line change
Expand Up @@ -8,10 +8,11 @@ library(tidyverse); library(raster); library(parallel)
# Removing Marine Hosts ####

Panth1 <- read.delim("data/PanTHERIA_1-0_WR05_Aug2008.txt") %>%
dplyr::rename(Sp = MSW05_Binomial, hOrder = MSW05_Order)
dplyr::rename(Sp = MSW05_Binomial, hOrder = MSW05_Order, hFamily = MSW05_Family)
Panth1$Sp <- Panth1$Sp %>% str_replace(" ", "_")

Panth1 %>% filter(hOrder%in%c("Cetacea", "Sirenia")) %>% pull(Sp) ->
Panth1 %>% filter(hOrder%in%c("Cetacea", "Sirenia")|
hFamily%in%c("Phocidae", "Odobenidae", "Otariidae")) %>% pull(Sp) ->
MarineSp

# Land use ####
Expand Down Expand Up @@ -89,11 +90,21 @@ Root <- paste0("Iceberg Input Files/","MaxEnt","/01_Raw/Currents")
Files <- list.files(Root)
MaxEntSp <- Files %>% str_remove(".tif$")

# Dispersals ####

Dispersals <- read.csv("Iceberg Input Files/Data for dispersal_Corrected.csv", header = T)

Dispersals$Scientific_name <- Dispersals$Scientific_name %>% str_replace(" ", "_")
Dispersals <- Dispersals %>% filter(!is.na(Scientific_name), !is.na(disp50))

Dispersals %>% pull(Scientific_name) -> DispersalSp

# Combining them ####

PhylogeneticSp %>%
setdiff(MarineSp) %>%
setdiff(NonEutherianSp) %>%
intersect(DispersalSp) %>%
intersect(MaxEntSp) ->
FullMaxEntSp

Expand All @@ -102,6 +113,7 @@ FullMaxEntSp; length(FullMaxEntSp)
PhylogeneticSp %>%
setdiff(MarineSp) %>%
setdiff(NonEutherianSp)%>%
intersect(DispersalSp) %>%
setdiff(NullRangeBagRares) %>%
intersect(RangeBagSp) ->
FullRangeBagSp
Expand All @@ -112,41 +124,3 @@ OnlyRangeBags <- setdiff(FullRangeBagSp, FullMaxEntSp)

SpeciesList <- list(MaxEnt = FullMaxEntSp,
RangeBags = FullRangeBagSp)

stop()

# Dispersals ####

Dispersals <- read.csv("Iceberg Input Files/Data for dispersal.csv", header = T)

Dispersals <- Dispersals %>% filter(!is.na(Scientific_name), !is.na(disp50))
DispersalSp <- Dispersals$Scientific_name <- Dispersals$Scientific_name %>% str_replace(" ", "_")

Panth1 %>% mutate(DispersalKnown = as.numeric(Sp%in%DispersalSp)) %>%
filter(Sp%in%unlist(SpeciesList)) %>% ggregplot::ggMMplot("hOrder", "DispersalKnown")

lapply(NewEncountersList, function(a){

lapply(a, function(b){

AllSp <- c(b$Sp, b$Sp2)

table(AllSp%in%DispersalSp)/sum(table(AllSp%in%DispersalSp))

})

})


lapply(NewEncountersList, function(a){

lapply(a, function(b){

b %>% filter(!Sp%in%DispersalSp|!Sp2%in%DispersalSp) %>% nrow %>% magrittr::divide_by(nrow(b))

})

})



179 changes: 179 additions & 0 deletions Greg ENM Code/00b_Unknown Dispersals.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,179 @@

# 00b_Assigning Dispersals to unknown generation times


# Iceberg pre-GAM spatial processing ####

# Rscript "01_Iceberg ENM Currents.R"

library(tidyverse); library(raster); library(parallel); library(sf); library(Matrix); library(magrittr); library(SpRanger); library(cowplot)

CORES <- 45

t1 <- Sys.time()

print("Dropping Species!")

Method = "MaxEnt"
# Method = "RangeBags"

source("00_Iceberg Species Dropping.R")

paste0("Iceberg Input Files/","MaxEnt","/01_Raw/Currents") %>%
list.files(full.names = T) %>%
append(paste0("Iceberg Input Files/","RangeBags","/01_Raw/Currents") %>% list.files(full.names = T)) ->
FullFiles

paste0("Iceberg Input Files/","MaxEnt","/01_Raw/Currents") %>%
list.files() %>% str_remove(".rds$") %>%
append(paste0("Iceberg Input Files/","RangeBags","/01_Raw/Currents") %>% list.files() %>% str_remove(".rds$")) ->
names(FullFiles)

Species <- SpeciesList %>% unlist %>% sort()

Files <- FullFiles[Species]

PredReps <- c("Currents", paste0("Futures", 1:4))

# Blanks
blank <- matrix(0,360*2,720*2) # proper resolution
blank <- raster(blank)
extent(blank) <- c(-180,180,-90,90)
projection(blank) <- CRS("+proj=longlat +datum=WGS84")

UniversalBlank <- raster("Iceberg Input Files/UniversalBlank.tif")
Land = which(raster::values(UniversalBlank)==0)
Sea = which(is.na(raster::values(UniversalBlank)))

# Grid areas
AreaRaster <- raster("Iceberg Input Files/LandArea.asc")
AreaValues <- raster::values(AreaRaster)

# Land Use Data
iucndat <- read.csv('Iceberg Input Files/IucnHabitatData.csv')
convcodes <- read.csv('Iceberg Input Files/IUCN_LUH_conversion_table.csv')

iucndat %>%
left_join(convcodes, by = c("code" = "IUCN_hab")) %>%
mutate(name = name %>% str_replace(" ", "_")) ->
Habitats

lapply(Species, function(a){

Habitats %>% filter(name == a) %>% pull(LUH)

}) -> HabitatList

names(HabitatList) <- Species

landuse2017 <- brick('Iceberg Input Files/landuse2017.grd')

# Continents ####
print("Continents!")

ContinentRaster <- raster("Iceberg Input Files/continents-greenland.tif") %>%
resample(blank, method = "ngb")

ContinentWhich <- lapply(1:6, function(a) which(values(ContinentRaster)==a))
names(ContinentWhich) <- c("Africa", "Eurasia", "Greenland", "NAm", "Oceania", "SAm")

# IUCN ranges for continent clipping ####
print("IUCN!")

load("~/LargeFiles/MammalStackFullMercator.Rdata")
load("Iceberg Input Files/IUCNBuffers.Rdata")

IUCNSp <- names(MammalStackFull) %>% intersect(Species)
MammalStackFull <- MammalStackFull[IUCNSp]

# Dispersals ####

Dispersals <- read.csv("Iceberg Input Files/Data for dispersal.csv", header = T)

Dispersals$Scientific_name <- Dispersals$Scientific_name %>% str_replace(" ", "_")
Dispersals <- Dispersals %>% filter(!is.na(Scientific_name))

ToBuffer <- intersect(Species, Dispersals %>% filter(!is.na(disp50)) %>% pull(Scientific_name))

ToFill <- setdiff(Species, ToBuffer)

write.csv(ToFill, file = "Iceberg Input Files/ImputedDispersals.csv")

# Assigning generation times ####

Dispersals[Dispersals$Scientific_name%in%ToFill,"Carnivore"] # All there
Dispersals[Dispersals$Scientific_name%in%ToFill,"BodyMass.Value"] # All there
Dispersals[Dispersals$Scientific_name%in%ToFill,"GenerationLength_d"] # All NA's

# Importing nearest phylogenetic neighbour ####

if(!file.exists("Iceberg Input Files/FullSTMatrix.csv")){

library(geiger);library(ape);library(picante);library(dplyr)

STFull <- read.nexus("data/ele_1307_sm_sa1.tre")[[1]]
FullSTMatrix <- as.data.frame(cophenetic(STFull)) %>% as.matrix

write.csv(FullSTMatrix, file = "Iceberg Input Files/FullSTMatrix.csv", row.names = F)

} else{FullSTMatrix <- as.matrix(read.csv("Iceberg Input Files/FullSTMatrix.csv", header = T)) }

colnames(FullSTMatrix) %>% setdiff(ToFill, .)

rownames(FullSTMatrix) <- colnames(FullSTMatrix)

ThinSTMatrix <- FullSTMatrix[,intersect(colnames(FullSTMatrix),
Dispersals %>% filter(!is.na(disp50)) %>% pull(Scientific_name))]

for(a in 1:length(ToFill)){

Sp <- ToFill[a]

Vector <- ThinSTMatrix[Sp,]

MinVector <- which(Vector==min(Vector))

if(length(MinVector)>1){

MinVector %>% sample(1) %>% names ->
Sp2

}else{

MinVector %>% names ->
Sp2

}

FillGenerationLength <- Dispersals[Dispersals$Scientific_name == Sp2, "GenerationLength_d"] ->

Dispersals[Dispersals$Scientific_name == Sp, "GenerationLength_d"]

}

Dispersals[Dispersals$Scientific_name%in%ToFill,"GenerationLength_d"] # it worked!

# Colin's loop ####

for(i in 1:nrow(Dispersals)){

if(is.na(Dispersals$Carnivore[i])) {} else {

if(Dispersals$Carnivore[i]==1) {

Dispersals$disp[i] <- (40.7*(Dispersals$BodyMass.Value[i]/1000)^0.81)/(as.numeric(Dispersals$GenerationLength_d[i])/365.25)

} else {

Dispersals$disp[i] <- (3.31*(Dispersals$BodyMass.Value[i]/1000)^0.65)/(as.numeric(Dispersals$GenerationLength_d[i])/365.25)
}

}
}

Dispersals$disp50 <- 50*Dispersals$disp

write.csv(Dispersals,'Iceberg Input Files/Data for dispersal_Corrected.csv')

paste0("Iceberg Input Files/GretCDF/Currents/", ToFill,".rds") %>% file.remove()
paste0("Iceberg Input Files/GretCDF/Futures/", ToFill,".rds") %>% file.remove()
66 changes: 22 additions & 44 deletions Greg ENM Code/01_Iceberg ENM Currents.R
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,6 @@ t1 <- Sys.time()

print("Dropping Species!")

Method = "MaxEnt"
# Method = "RangeBags"

source("00_Iceberg Species Dropping.R")

paste0("Iceberg Input Files/","MaxEnt","/01_Raw/Currents") %>%
Expand Down Expand Up @@ -68,11 +65,11 @@ landuse2017 <- brick('Iceberg Input Files/landuse2017.grd')
# Continents ####
print("Continents!")

ContinentRaster <- raster("Iceberg Input Files/continents-greenland.tif") %>%
ContinentRaster <- raster("Iceberg Input Files/continents-madagascar.tif") %>%
resample(blank, method = "ngb")

ContinentWhich <- lapply(1:6, function(a) which(values(ContinentRaster)==a))
names(ContinentWhich) <- c("Africa", "Eurasia", "Greenland", "NAm", "Oceania", "SAm")
ContinentWhich <- lapply(1:max(values(ContinentRaster), na.rm = T), function(a) which(values(ContinentRaster)==a))
names(ContinentWhich) <- c("Africa", "Eurasia", "Greenland", "Madagascar", "NAm", "Oceania", "SAm")

# IUCN ranges for continent clipping ####
print("IUCN!")
Expand All @@ -85,16 +82,21 @@ MammalStackFull <- MammalStackFull[IUCNSp]

# Dispersals ####

Dispersals <- read.csv("Iceberg Input Files/Data for dispersal.csv", header = T)
Dispersals <- read.csv("Iceberg Input Files/Data for dispersal_Corrected.csv", header = T)

Dispersals <- Dispersals %>% filter(!is.na(Scientific_name), !is.na(disp50))
Dispersals$Scientific_name <- Dispersals$Scientific_name %>% str_replace(" ", "_")

ToBuffer <- intersect(Species, Dispersals$Scientific_name)

# 01_Processing Currents ####
# Adding in exception for bats ####

Panth1 %>% filter(hOrder == "Chiroptera") %>% pull(Sp) ->
BatSpecies

ToBuffer <- setdiff(ToBuffer, BatSpecies)

Root <- paste0("Iceberg Input Files/", Method,"/01_Raw/","Currents")
# 01_Processing Currents ####

# Setting raster standards ####

Expand Down Expand Up @@ -134,7 +136,9 @@ mclapply(1:length(ToProcess), function(i){
GretCDF <- data.frame(
X = seq(from = XMin, to = XMax, length.out = NCol) %>% rep(NRow),
Y = seq(from = YMax, to = YMin, length.out = NRow) %>% rep(each = NCol),

Climate = as.numeric(!is.na(values(RasterLista)))

)

# IUCN Buffer clipping ####
Expand Down Expand Up @@ -182,6 +186,13 @@ mclapply(1:length(ToProcess), function(i){
SpWhich <- which(!is.na(values(r1)))
ContinentsInhabited <- unique(values(ContinentRaster)[SpWhich])

if(Sp == "Canis_lupus"){

ContinentsInhabited <-
ContinentsInhabited %>% setdiff(1)

}

if(length(ContinentsInhabited)>0){

GretCDF$Continent <- 0
Expand Down Expand Up @@ -278,41 +289,8 @@ mclapply(1:length(ToProcess), function(i){

}, mc.preschedule = F, mc.cores = CORES)

stop()

t2 <- Sys.time()

print(t2 - t1)


# Troubleshooting ####

iucndat <- read.csv('Iceberg Input Files/IucnHabitatData.csv')
convcodes <- read.csv('Iceberg Input Files/IUCN_LUH_conversion_table.csv')
iucndat %>%
left_join(convcodes, by = c("code" = "IUCN_hab")) %>%
mutate(name = name %>% str_replace(" ", "_")) ->
Habitats

list.files(paste0("~/Albersnet/Iceberg Input Files/GretCDF/Currents")) %>% str_remove(".rds") ->
Species

lapply(Species, function(a){
Habitats %>% filter(name == a) %>% pull(LUH)
}) -> HabitatList
names(HabitatList) <- Species
HabitatList[sapply(HabitatList, length)>0]

Sp <- "Microtus_daghestanicus"
Sp <- "Acinonyx_jubatus"
Sp <- "Lycaon_pictus"

CDF <- readRDS(paste0("~/Albersnet/Iceberg Input Files/GretCDF/Currents/",Sp,".rds")) %>%
as.matrix %>% as.data.frame()

CDF %>% colSums()

CDF %>% gather(key = "Key", value = "Value", -X, -Y) %>%
mutate(Key = factor(Key, levels = unique(Key))) %>%
ggplot(aes(X, Y, fill = Value)) +
coord_fixed() + theme_cowplot() +
geom_tile() + facet_wrap(~Key)

Loading

0 comments on commit b028535

Please sign in to comment.