|
1
|
###################
|
|
2
|
##Time-lag analysis
|
|
3
|
#data format: species matrix
|
|
4
|
require(vegan)
|
|
5
|
|
|
6
|
#creating distance matrix
|
|
7
|
data.dm <- vegdist(spp_matrix, method="euclidean",
|
|
8
|
diag=FALSE, upper=FALSE, na.rm=TRUE)
|
|
9
|
as.matrix(data.dm)
|
|
10
|
|
|
11
|
#dividing each column from matrix into new data.frames and removing upper triangle values
|
|
12
|
c1 <- data.frame(data.dm[2:25,1])
|
|
13
|
c2 <- data.frame(data.dm[3:25, 2])...
|
|
14
|
|
|
15
|
#binding and renaming new columns
|
|
16
|
dm.final <- cbindX(c1,c2,c3...)
|
|
17
|
colnames(dm.final) <- c("x1", "x2", "x3"...)
|
|
18
|
|
|
19
|
#adding time-lag ID column
|
|
20
|
dm.final$time.lag <- c(1:yr#-1)
|
|
21
|
|
|
22
|
#melting matrix into long format so to do regression
|
|
23
|
require(reshape2)
|
|
24
|
|
|
25
|
dm.long <- melt(dm.final, idvars = time.lag,
|
|
26
|
measure.vars = 1:24,
|
|
27
|
value.name = "distance",
|
|
28
|
na.rm = TRUE)
|
|
29
|
#regression
|
|
30
|
require(car)
|
|
31
|
|
|
32
|
lm.dm <- lm(distance ~ time.lag, data=dm.long)
|
|
33
|
|
|
34
|
#obtaining p-value
|
|
35
|
#method= "pearson", "spearman"
|
|
36
|
cor.test(dm.long$time.lag, dm.long$distance, method="pearson")
|
|
37
|
|
|
38
|
#plotting regression
|
|
39
|
#see Plotting info below
|
|
40
|
par(mar=c(5, 5, 1, 1), xpd=TRUE)
|
|
41
|
plot(dm.long$time.lag, dm.long$distance, type="n",
|
|
42
|
xlab="distance",
|
|
43
|
ylab="time lag",
|
|
44
|
xlim=c(xMin,xMax), ylim=c(yMin,yMax))
|
|
45
|
|
|
46
|
points(dm.long$time.lag, dm.long$distance, cex = 1.0, pch=16, col="#99CC66")
|
|
47
|
|
|
48
|
#plotting regression line
|
|
49
|
require(plotrix)
|
|
50
|
ablineclip(lm.dm, x1=0, x2=25, type="l", lty= 1, lwd=3.5, col="#000000")
|
|
51
|
#a=yStart, b=yEnd, x1=xStart, x2= xEnd
|