annotate seek_otu.R @ 0:bbaa89f070f4 draft

planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
author iuc
date Mon, 04 Mar 2024 19:56:16 +0000
parents
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
1 #!/usr/bin/env Rscript
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
2
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
3 ## Redirect R error handling to stderr.
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
4 options(show.error.messages = FALSE, error = function() {
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
5 cat(geterrmessage(), file = stderr())
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
6 q("no", 1, FALSE)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
7 })
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
8
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
9 ## Avoid crashing Galaxy with a UTF8 error on German LC settings
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
10 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
11
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
12 args <- commandArgs(trailingOnly = TRUE)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
13 if (length(args) == 0) {
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
14 stop("Arguments missing for Rscrpit", call. = FALSE)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
15 } else {
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
16 # percentage of identity
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
17 id_threshold <- as.numeric(args[3])
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
18 # get input data (matrix)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
19 data <- read.csv(args[1], header = FALSE, sep = ",", row.names = 1)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
20 # remove last 2 columns
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
21 data_length <- length(data)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
22 # create matrix
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
23 mat <- as.matrix(data[, 1:data_length], fill = TRUE)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
24 # create coordinate matrix
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
25 d <- as.dist(1 - mat)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
26 # create tree
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
27 hc <- hclust(d, method = "single")
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
28 # assign otu based on identity value
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
29 otu <- cutree(hc, h = -id_threshold)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
30 # group contigs by otu
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
31 # Print results to output file
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
32 output <- args[2]
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
33 # unique is used to know the number of different otu
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
34 for (i in unique(otu)) {
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
35 # retrieve contigs belonging to the same otu
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
36 clust <- which(otu == i)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
37 # write otu number and number of contigs in this otu
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
38 cat(
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
39 paste("OTU_", i, ",", length(clust), ",", sep = ""),
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
40 file = output, append = TRUE
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
41 )
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
42 for (n in names(clust)) {
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
43 # write contigs name
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
44 cat(paste(gsub(" ", "", n), ",", sep = ""), file = output, append = TRUE)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
45 }
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
46 cat("\n", sep = "", file = output, append = TRUE)
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
47 }
bbaa89f070f4 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/virAnnot commit 3a3b40c15ae5e82334f016e88b1f3c5bbbb3b2cd
iuc
parents:
diff changeset
48 }