Mercurial > repos > ecology > ecoregion_brt_analysis

annotate nb_clust_G.R @ 3:d63d74194100 draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit ced658540f05bb07e1e687af30a3fa4ea8e4803c
author ecology
date Wed, 28 May 2025 10:10:52 +0000
parents 7d77be8fab92
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
1 # Script to determine the optimal number of clusters thanks to the optimization of the SIH index and to produce the files needed in the next step of clustering
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
2
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
3 #load packages
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
4 library(cluster)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
5 library(dplyr)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
6 library(tidyverse)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
7
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
8 #load arguments
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
9 args = commandArgs(trailingOnly=TRUE)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
10 if (length(args)==0)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
11 {
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
12 stop("This tool needs at least one argument")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
13 }else{
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
14 enviro <- args[1]
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
15 taxa_list <- args[2]
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
16 preds <- args[3]
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
17 max_k <- as.numeric(args[4])
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
18 metric <- args[5]
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
19 sample <- as.numeric(args[6])
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
20 }
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
21
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
22 #load data
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
23
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
24 env.data <- read.table(enviro, sep="\t", header = TRUE, dec = ".", na.strings = "-9999")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
25
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
26 ##List of modelled taxa used for clustering
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
27 tv <- read.table(taxa_list, dec=".", sep=" ", header=F, na.strings = "NA")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
28 names(tv) <- c("a")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
29
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
30 ################Grouping of taxa if multiple prediction files entered ################
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
31
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
32 data_split = str_split(preds,",")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
33 data.bio = NULL
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
34
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
35 for (i in 1:length(data_split[[1]])) {
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
36 data.bio1 <- read.table(data_split[[1]][i], dec=".", sep="\t", header=T, na.strings = "NA")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
37 data.bio <- rbind(data.bio,data.bio1)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
38 remove(data.bio1)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
39 }
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
40
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
41 names(data.bio) <- c("lat", "long", "pred", "taxon")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
42
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
43 #keep selected taxa
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
44 data.bio <- data.bio[which(data.bio$taxon %in% tv$a),]
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
45
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
46 write.table(data.bio,file="data_bio.tsv",sep="\t",quote=F,row.names=F)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
47
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
48 #format data
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
49
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
50 test3 <- matrix(data.bio$pred , nrow = nrow(env.data), ncol = nrow(data.bio)/nrow(env.data))
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
51 test3 <- data.frame(test3)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
52 names(test3) <- unique(data.bio$taxon)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
53
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
54 write.table(test3, file="data_to_clus.tsv", sep="\t",quote=F,row.names=F)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
55
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
56 #Max number of clusters to test
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
57 max_k <- max_k
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
58
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
59 # Initialization of vectors to store SIH indices
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
60 sih_values <- rep(0, max_k)
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
61
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
62 # Calculation of the SIH index for each number of clusters
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
63 for (k in 2:max_k) {
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
64 # Clara execution
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
65 clara_res <- clara(test3, k, metric =metric, samples = sample, sampsize = min(nrow(test3), (nrow(data.bio)/nrow(test3))+2*k))
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
66 # Calculation of the SIH index
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
67 sih_values[k] <- clara_res$silinfo$avg.width
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
68 }
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
69
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
70 # Plot SIH Index Chart by Number of Clusters
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
71 png("Indices_SIH.png")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
72 plot(2:max_k, sih_values[2:max_k], type = "b", xlab = "Number of clusters", ylab = "SIH index")
7d77be8fab92 planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/Ecoregionalization_workflow commit e03df85746a3b61a382a5ee7e3357a8bf42a5097
ecology
parents:
diff changeset
73 dev.off()