Skip to content

Commit 9ddb373

Browse files
GLatombeGLatombe
committed
First commit
1 parent a431107 commit 9ddb373

24 files changed

+1387
-1
lines changed

‎.gitignore‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,4 @@
1+
.Rproj.user
2+
.Rhistory
3+
.RData
4+
.Ruserdata

‎Kernel-Integral-Method.Rproj‎

Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,13 @@
1+
Version: 1.0
2+
3+
RestoreWorkspace: Default
4+
SaveWorkspace: Default
5+
AlwaysSaveHistory: Default
6+
7+
EnableCodeIndexing: Yes
8+
UseSpacesForTab: Yes
9+
NumSpacesForTab: 2
10+
Encoding: UTF-8
11+
12+
RnwWeave: Sweave
13+
LaTeX: pdfLaTeX

‎README.md‎

Lines changed: 7 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1 +1,7 @@
1-
# Kernel-Integral-Method
1+
# Kernel-Integral-Method
2+
The main scripts calling all other functions are:
3+
- Run_FT_replicates.R
4+
- Test_abundance_all_figures.R
5+
- Test_abundance.R
6+
7+
The description of each function is included as a comment when they are called

‎Run_FT_replicates.R‎

Lines changed: 103 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,103 @@
1+
rm(list=ls(all.names = TRUE))
2+
gc()
3+
graphics.off()
4+
#quartz()
5+
6+
7+
##I admit I lost track of all the functions that are called, and not all packages may not be necessary, but better safe than sorry
8+
library("BAT") #version 2.0.1
9+
library("FD") #version 1.1-12
10+
library("ggplot2")
11+
library("gridExtra")
12+
library("labdsv")
13+
library("hypervolume") #version 2.0.11
14+
library("psych")
15+
library("StatMatch")
16+
library("TPD") #version 1.1.0
17+
library("geometry")
18+
library("mFD")
19+
library("MASS")
20+
library("rgl")
21+
library("reshape2")
22+
library("TreeTools")
23+
library("StatMatch")
24+
library("picante")
25+
library("alphahull")
26+
library("data.table")
27+
library("progress")
28+
library("abind")
29+
30+
##function to generate the theoretical data
31+
source("creation_data.R")
32+
33+
##Functions for the tree-based method
34+
source("treeFT.R") ##wrapper function calling all the others
35+
source("drop.tip.custom.R") ##adapt drop.tip so that the id of the nodes and tips is kept the same, to be able to compute the intersection of the trees. The resulting trees are not valid, but this will be fixed by reorder.phy.
36+
source("keep.tip.custom.R") ##like keep.tip(), but calls drop.tip.custom() instead
37+
source("reorder.phy.R") ##sort the id of the nodes and tips to make a valid phylo object
38+
source("merge.trees.R") ##compute the intersection between i trees
39+
40+
41+
##Functions for the functional-space methods
42+
source("all_FT_function.R") ##This is the function generating the theoretical point distributions and computing the different indices of functional turnover
43+
source("beta_BAT.R") ##wrapper function to compute convex hull indices
44+
source("hypervolume.bat.R") ##wrapper function to compute hypervolume-based indices
45+
source("for.multi.R") ##function to compute KIM V1 and V2 indices
46+
source("for.multi.nonunif.R") ##function to compute KIM V3 indices (no resampling of random points)
47+
source("hypervolume_gaussian.custom.R") ##compute point distributions using a modification of the hypervolume_gaussian() function from package hypervolume, to stop everything before resampling to uniform density
48+
source("sample_model_ellipsoid.custom.R") ##modification of sample_model_ellipsoid() to get the correct points
49+
source("predict_function_gaussian.custom.R") ##modification of predict_function_gaussian() needed for sample_model_ellipsoid.custom()
50+
source("estimate_bandwidth_silent.R") ##just to avoid some annoying message that is not an error
51+
52+
53+
54+
#creation of different square configurations
55+
configs <- c("overlap","disjoint","nested") ##relative position of the squares
56+
sizediffs <- c("same","different") ##size of the two squares
57+
npointss <- c(10,40,70,100) ##number of species-ponts
58+
59+
is.res4 <-0
60+
for(config in configs){
61+
is.res3 <-0
62+
for(sizediff in sizediffs){
63+
is.res2 <-0
64+
for(npoints in npointss){
65+
is.res1 <-0
66+
for(r in 1:50){ ##replicates
67+
if(is.res1==0){
68+
results1 <- all_FT_function(config=config,sizediff=sizediff,npoints=npoints)
69+
is.res1 <- 1
70+
}else{
71+
res <- all_FT_function(config=config,sizediff=sizediff,npoints=npoints)
72+
results1 <- abind(results1,res,along=4)
73+
}
74+
}
75+
if(is.res2==0){
76+
results2 <- results1
77+
is.res2 <- 1
78+
}else{
79+
results2 <- abind(results2,results1,along=5)
80+
}
81+
}
82+
if(is.res3==0){
83+
results3 <- results2
84+
is.res3 <- 1
85+
}else{
86+
results3 <- abind(results3,results2,along=6)
87+
}
88+
}
89+
if(is.res4==0){
90+
results4 <- results3
91+
is.res4 <- 1
92+
}else{
93+
results4 <- abind(results4,results3,along=7)
94+
}
95+
}
96+
97+
98+
# save(results4,file="wkspace.RData")
99+
100+
101+
102+
103+

‎Test_abundance.R‎

Lines changed: 110 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,110 @@
1+
rm(list=ls(all.names = TRUE))
2+
gc()
3+
graphics.off()
4+
quartz()
5+
6+
7+
##I admit I lost track of all the functions that are called, and not all packages may not be necessary, but better safe than sorry
8+
library("BAT") #version 2.0.1
9+
library("FD") #version 1.1-12
10+
library("ggplot2")
11+
library("gridExtra")
12+
library("labdsv")
13+
library("hypervolume") #version 2.0.11
14+
library("psych")
15+
library("StatMatch")
16+
library("TPD") #version 1.1.0
17+
library("geometry")
18+
library("mFD")
19+
library("MASS")
20+
library("rgl")
21+
library("reshape2")
22+
library("TreeTools")
23+
library("StatMatch")
24+
library("picante")
25+
library("alphahull")
26+
library("data.table")
27+
library("progress")
28+
library("abind")
29+
30+
##function to generate the theoretical data
31+
source("creation_data.R")
32+
33+
##Functions for the tree-based method
34+
source("treeFT.R") ##wrapper function calling all the others
35+
source("drop.tip.custom.R") ##adapt drop.tip so that the id of the nodes and tips is kept the same, to be able to compute the intersection of the trees. The resulting trees are not valid, but this will be fixed by reorder.phy.
36+
source("keep.tip.custom.R") ##like keep.tip(), but calls drop.tip.custom() instead
37+
source("reorder.phy.R") ##sort the id of the nodes and tips to make a valid phylo object
38+
source("merge.trees.R") ##compute the intersection between i trees
39+
40+
41+
##Functions for the functional-space methods
42+
source("all_FT_function.R") ##This is the function generating the theoretical point distributions and computing the different indices of functional turnover
43+
source("beta_BAT.R") ##wrapper function to compute convex hull indices
44+
source("hypervolume.bat.R") ##wrapper function to compute hypervolume-based indices
45+
source("for.multi.R") ##function to compute KIM V1 and V2 indices
46+
source("for.multi.nonunif.R") ##function to compute KIM V3 indices (no resampling of random points)
47+
source("hypervolume_gaussian.custom.R") ##compute point distributions using a modification of the hypervolume_gaussian() function from package hypervolume, to stop everything before resampling to uniform density
48+
source("sample_model_ellipsoid.custom.R") ##modification of sample_model_ellipsoid() to get the correct points
49+
source("predict_function_gaussian.custom.R") ##modification of predict_function_gaussian() needed for sample_model_ellipsoid.custom()
50+
source("estimate_bandwidth_silent.R") ##just to avoid some annoying message that is not an error
51+
52+
53+
54+
55+
56+
configs <- c("overlap") ##relative position of the squares
57+
sizediffs <- c("different") ##size of the two squares
58+
#npointss <- c(10) ##number of species-ponts
59+
npoints <- 2
60+
bandwith.fac <- 1
61+
62+
Site1AX <- c(1,1,5,5)
63+
Site1AY <- c(4,8,4,8)
64+
65+
Site3AX <- c(4,4,6,6)
66+
Site3AY <- c(3,5,3,5)
67+
68+
##different sizes
69+
par(mfrow=c(2,2))
70+
comms_similar <- creation_data(comm1 = data.frame(X = Site1AX, Y = Site1AY),comm2 = data.frame(X = Site3AX, Y= Site3AY),mode = "Similar",points_number = npoints,tol.abs = 1,tol.rel = NULL, plot=TRUE)
71+
comms_uniform <- creation_data(comm1 = data.frame(X = Site1AX, Y = Site1AY),comm2 = data.frame(X = Site3AX, Y= Site3AY),mode = "Uniform",points_number = npoints,tol.abs = 1,tol.rel = NULL, plot=TRUE)
72+
comms_diff <- creation_data(comm1 = data.frame(X = Site1AX, Y = Site1AY),comm2 = data.frame(X = Site3AX, Y= Site3AY),mode = "Different",points_number = npoints,tol.abs = 1,tol.rel = NULL, plot=TRUE)
73+
74+
75+
76+
##If you play with the comments, you can test different combinations: same abundance and trait variability for all species, or make on or both vary for each species. In practice, you will want to play with the values to get something that works for your data
77+
78+
#abundance <- 1000 ##same abundance for all species (1000 random points generated)
79+
abundance <- matrix(ceiling(1000*runif(npoints+4)),2,npoints+4) ## different abundance for all species (within [1,1000])
80+
# trait.var <- 1 ##same trait variability for all species and for the two traits (that will make circles with a radius of 1)
81+
trait.var <- list() ## different trait variability for all species and for each trait (within [0,1])
82+
trait.var[[1]] <- matrix(runif((npoints+4)*2),npoints+4,2)
83+
trait.var[[2]] <- matrix(runif((npoints+4)*2),npoints+4,2)
84+
85+
#Kernel-based function fully adjusted: same bandwidth for the 2 communities and no resampling of random points to uniform distribution (1000 points for each species, leading to changes in density in the functional space)
86+
kernel.res.nonunif <- data.frame(Config=c("Joint_different","Joint_uniform","Joint_similar"),Jac=numeric(3),Turn=numeric(3),Nest=numeric(3))
87+
kernel.res.diff.nonunif <- for.multi.nonunif(comms_diff[[1]], site1 = "Site1", site2 = "Site2", trait_data = comms_diff[[2]],res = 0.05, bandwith.fac = bandwith.fac,chunk.size=abundance,scales=trait.var)
88+
kernel.res.unif.nonunif <- for.multi.nonunif(comms_uniform[[1]], site1 = "Site1", site2 = "Site2", trait_data = comms_uniform[[2]],res = 0.05, bandwith.fac = bandwith.fac,chunk.size=abundance,scales=trait.var)
89+
kernel.res.simil.nonunif <- for.multi.nonunif(comms_similar[[1]], site1 = "Site1", site2 = "Site2", trait_data = comms_similar[[2]],res = 0.05, bandwith.fac = bandwith.fac,chunk.size=abundance,scales=trait.var)
90+
91+
kernel.res.nonunif[1,2:4] <- kernel.res.diff.nonunif$indices
92+
kernel.res.nonunif[2,2:4] <- kernel.res.unif.nonunif$indices
93+
kernel.res.nonunif[3,2:4] <- kernel.res.simil.nonunif$indices
94+
95+
plot.new()
96+
plot(comms_uniform[[2]],col=c(rep("black",npoints+4),rep("red",npoints+4)),xlim=c(0,7),ylim=c(2,9))
97+
plot(kernel.res.unif.nonunif$points1,pch=20,cex=0.1,xlim=c(0,7),ylim=c(2,9))
98+
points(kernel.res.unif.nonunif$points2,pch=20,cex=0.1,col="red")
99+
100+
101+
102+
103+
104+
105+
106+
107+
108+
109+
110+

0 commit comments

Comments
 (0)