-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathSalmon_LCDist_Script.R
70 lines (53 loc) · 2.48 KB
/
Salmon_LCDist_Script.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
#Get least cost distances among sampling sites using marmap
#Set working directory
setwd("C:/Users/Nick/Desktop/Postdoc/Atlantic Salmon Project/Paper 1 - Clines/sPCA")
#libraries -------------
library(gdistance)
library(maps)
library(mapdata)
library(maptools)
library(marmap)
library(PBSmapping)
library(ggplot2)
library(dplyr)
library(data.table)
library(gtable)
library(vegan)
#source mapping functions
# source("c:/Users/rystanley/OneDrive/PostDoc/DFO/Cline Analysis/Cartesian Analysis/ScaleBar.R")
#Make colours
blues <- c("lightsteelblue4", "lightsteelblue3",
"lightsteelblue2", "lightsteelblue1")
greys <- c(grey(0.6), grey(0.93), grey(0.99))
#data ---------
#Load workspace
#Load coordinates as degrees Lat and Long
SpeciesCoordinates <- read.csv("Salmon_Coords.csv")
## Set map limits------
Lat.lim=c(42,62)
Long.lim=c(-72,-52)
#Download bathymetry layer for salmon sampling range - this gets the bathymetry for your area to then calculate lc dists
#Setting the resolution to a lower number (e.g. 1 is lowest) takes longer but grabs more detail
getNOAA.bathy(lon1 = Long.lim[1], lon2 = Long.lim[2], lat1 = Lat.lim[1], lat2 = Lat.lim[2],
resolution = 1,keep=TRUE) -> SalmonDepths
#plot the bathy map to see how it looks
plot(SalmonDepths,image = TRUE, land = TRUE, lwd = 0.03,
bpal = list(c(0, max(SalmonDepths), greys),
c(min(SalmonDepths), 0, blues)))
plot(SalmonDepths, lwd = 1, deep = 0, shallow = 0, step = 0, add = TRUE) # highlight coastline
points(SpeciesCoordinates$AdjustLong, SpeciesCoordinates$AdjustLat,pch=19,cex=1,col="black") #Add points
#Add text labels
text(x=SpeciesCoordinates$AdjustLong,y=SpeciesCoordinates$AdjustLat,labels=SpeciesCoordinates$Code, cex=0.6,pos=3,offset = 0.5,font = 2)
#### Calculate distance (km) between stations for each species -------------
#Transition object - this can take a few minutes and is used as the trans object in the lc.dist function below
trans1 <- trans.mat(SalmonDepths,min.depth = -10)
sites1<-SpeciesCoordinates[,5:6]
rownames(sites1)<-SpeciesCoordinates$Code
#calculate least cost distances - just use res="dist" to get a distance matrix
dist.salmon <- lc.dist(trans1,
sites1,
res="dist")
#res = "path" will give you a list of all distances between each pairwise site but can take a long time
salmon.lines<- lc.dist(trans1,
sites1,
res="path")