Mercurial > repos > iuc > virannot_blast2tsv

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