Mercurial > repos > iuc > virannot_rps2tsv
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 |
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 } |