accCost2 <- function(x, fromCoords) {
fromCells <- cellFromXY(x, fromCoords)
tr <- transitionMatrix(x)
tr <- rBind(tr, rep(0, nrow(tr)))
tr <- cBind(tr, rep(0, nrow(tr)))
startNode <- nrow(tr)
adjP <- cbind(rep(startNode, times = length(fromCells)), fromCells)
tr[adjP] <- Inf
adjacencyGraph <- graph.adjacency(tr, mode = "directed", weighted = TRUE)
E(adjacencyGraph)$weight <- 1/E(adjacencyGraph)$weight
return(shortest.paths(adjacencyGraph, v = startNode, mode = "out")[-startNode])
}
Suitabilities <- list()
for(j in 1:length(BinSpp)){
Suitability <- list()
for (i in 1:nlayers(BinSpp[[j]])){
temp <- data.frame(Suitability = values(BinSpp[[j]][[i]]), ID = 1:length(values(BinSpp[[j]][[i]])), Time = i-1)
Suitability[[i]] <- temp[complete.cases(temp),]
}
Suitabilities[[j]]<- do.call("rbind", Suitability)
Suitabilities[[j]]$Spp <- names(BinSpp)[j]
}
Suitability <- do.call("rbind", Suitabilities)
s <- Suitability %>% group_by(ID) %>% summarise(SUMA = sum(Suitability)) %>% filter(SUMA > 0)
Suitability <- Suitability[Suitability$ID %in% s$ID,]
Spps <- unique(Suitability$Spp)
Suitability <- Suitability[,c(4,1,2,3)]
Suitabilities <- list()
for (i in Spps){
Suitabilities[[i]] <- dplyr::filter(Suitability, Spp == i)
temp <- split(Suitabilities[[i]], Suitabilities[[i]]$Time)
Suitabilities[[i]] <- do.call(cbind, lapply(1:length(temp), function(i){
if (i == 1){
setNames(data.frame(paste("[",temp[[i]][["Spp"]],",",temp[[i]][["ID"]], ",*","]", sep = ""), temp[[i]]["Time"], temp[[i]]["Suitability"]),
c("V", paste("T", i, sep = ""), i-1))
} else if (i == length(temp)){
setNames(data.frame(temp[[i]]["Time"], temp[[i]]["Suitability"], rep("\n", NROW(temp[[i]]))),
c(paste("T", i, sep = ""), i-1, "line"))
} else {
setNames(data.frame(temp[[i]]["Time"], temp[[i]]["Suitability"]),
c(paste("T", i, sep = ""), i-1))
}
}))
}
Suitability <-do.call("rbind", Suitabilities)
conns <- list()
for(j in 1:length(BinSpp)){
Raster <- sum(BinSpp[[j]])
Raster[values(Raster) > 0] = 1
Raster[values(Raster) == 0] = NA
h16 <- transition(Raster, transitionFunction=function(x){1},16,symm=FALSE)
h16 <- geoCorrection(h16, scl=FALSE)
ID <-c(1:ncell(Raster))[!is.na(values(Raster))]
B <- xyFromCell(Raster, cell = ID)
connections <- list()
#For each pair of cells in B
for (i in 1:nrow(B)){
#Create a temporal raster for each row with the distance from cell xy to all other cells
temp <- accCost2(h16,B[i,])
index <- which(temp < Dist)
connections[[i]] <- cbind(ID[i], index, temp[index])
}
#Get everything together as a large data frame
connections <- do.call("rbind", connections)
connections <- as.data.frame(connections)
colnames(connections) <- c("from", "to", "dist")
connections$Sp <- names(BinSpp)[j]
conns[[j]] <- connections
}
connections <- conns
connections <- do.call("rbind", connections)
Nchains <- data.frame(Spp = Spps, Nchains = nchains, Space = "\n")
Cost <- values(costlayer)[unique(unique(connections$to), unique(connections$to))]
sink(paste0(name, ".dat"))
cat(c("set V :=", unique(unique(connections$to), unique(connections$to)), ";"))
cat("\n")
cat("\n")
cat(c("set SP :=", names(BinSpp), ";"))
cat("\n")
cat("\n")
cat(c("set c :=", Cost, ";"))
cat("\n")
cat("\n")
cat(c("set E :=", paste0("(",unique(unite_(connections, col = "V", sep = ",", from = c("from", "to"))$V), ")"), ";"))
cat("\n")
cat("\n")
cat(paste0("param T:= ", (nlayers(BinSpp[[1]])-1),";"))
cat("\n")
cat("\n")
cat("param u :=")
cat("\n")
cat(do.call(paste, Suitability))
cat(";")
cat("\n")
cat("param nchains := ")
cat("\n")
cat(do.call(paste, Nchains))
cat(";")
cat("\n")
sink()
return(list(connections = connections, Suitability = Suitability))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.