Mercurial > repos > iuc > ruvseq

annotate get_deseq_dataset.R @ 4:3a083c78896e 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/ruvseq commit 30117fce22f3771c9c0c13e70c3ad14b694de6e2
author iuc
date Fri, 21 Apr 2023 14:09:17 +0000
parents fed9d0350d72
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
1 get_deseq_dataset <- function(sample_table, header, design_formula, tximport, txtype, tx2gene) {
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
2
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
3 dir <- ""
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
4
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
5 has_header <- !is.null(header)
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
6 use_txi <- !is.null(tximport)
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
7 if (use_txi) {
1
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
8 if (is.null(tx2gene)) stop("A transcript-to-gene map or a GTF/GFF3 file is required for tximport")
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
9 if (tolower(file_ext(tx2gene)) == "gff") {
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
10 gff_file <- tx2gene
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
11 } else {
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
12 gff_file <- NULL
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
13 tx2gene <- read.table(tx2gene, header = has_header)
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
14 }
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
15 }
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
16
4
3a083c78896e planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 30117fce22f3771c9c0c13e70c3ad14b694de6e2
iuc
parents: 2
diff changeset
17 if (!use_txi && has_header) {
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
18 countfiles <- lapply(as.character(sample_table$filename), read.delim, row.names = 1)
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
19 tbl <- do.call("cbind", countfiles)
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
20 colnames(tbl) <- rownames(sample_table) # take sample ids from header
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
21
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
22 # check for htseq report lines (from DESeqDataSetFromHTSeqCount function)
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
23 old_special_names <- c(
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
24 "no_feature",
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
25 "ambiguous",
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
26 "too_low_aQual",
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
27 "not_aligned",
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
28 "alignment_not_unique"
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
29 )
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
30 special_rows <- (substr(rownames(tbl), 1, 1) == "_") | rownames(tbl) %in% old_special_names
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
31 tbl <- tbl[!special_rows, , drop = FALSE]
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
32
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
33 dds <- DESeqDataSetFromMatrix(
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
34 countData = tbl,
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
35 colData = subset(sample_table, select = -filename),
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
36 design = design_formula
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
37 )
4
3a083c78896e planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 30117fce22f3771c9c0c13e70c3ad14b694de6e2
iuc
parents: 2
diff changeset
38 } else if (!use_txi && !has_header) {
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
39
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
40 # construct the object from HTSeq files
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
41 dds <- DESeqDataSetFromHTSeqCount(
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
42 sampleTable = sample_table,
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
43 directory = dir,
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
44 design = design_formula
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
45 )
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
46 colnames(dds) <- row.names(sample_table)
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
47
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
48 } else {
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
49 # construct the object using tximport
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
50 library("tximport")
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
51 txi_files <- as.character(sample_table$filename)
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
52 labs <- row.names(sample_table)
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
53 names(txi_files) <- labs
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
54 if (!is.null(gff_file)) {
1
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
55 # first need to make the tx2gene table
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
56 # this takes ~2-3 minutes using Bioconductor functions
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
57 suppressPackageStartupMessages({
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
58 library("GenomicFeatures")
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
59 })
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
60 txdb <- makeTxDbFromGFF(gff_file)
1
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
61 k <- keys(txdb, keytype = "TXNAME")
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
62 tx2gene <- select(txdb, keys = k, columns = "GENEID", keytype = "TXNAME")
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
63 # Remove 'transcript:' from transcript IDs (when gff_file is a GFF3 from Ensembl and the transcript does not have a Name)
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
64 tx2gene$TXNAME <- sub("^transcript:", "", tx2gene$TXNAME) # nolint
1
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
65 }
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
66 try(txi <- tximport(txi_files, type = txtype, tx2gene = tx2gene))
1
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
67 if (!exists("txi")) {
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
68 # Remove version from transcript IDs in tx2gene...
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
69 tx2gene$TXNAME <- sub("\\.[0-9]+$", "", tx2gene$TXNAME) # nolint
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
70 # ...and in txi_files
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
71 txi <- tximport(txi_files, type = txtype, tx2gene = tx2gene, ignoreTxVersion = TRUE)
1
c24765926774 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 9ed3d83cc447ee897af867362bf1dd67af8a11c2
iuc
parents: 0
diff changeset
72 }
2
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
73 dds <- DESeqDataSetFromTximport(
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
74 txi,
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
75 subset(sample_table, select = -c(filename)),
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
76 design_formula
fed9d0350d72 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 4daa375d022673d2437d609b1865b78c64b04415"
iuc
parents: 1
diff changeset
77 )
0
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
78 }
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
79 return(dds)
61dffb20b6f9 planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit b95582cea8320d5488056a9576474f79cec53be8
iuc
parents:
diff changeset
80 }