Convert transcripts to genes

convertTranscriptsToGenes(object, ...)

# S4 method for character
convertTranscriptsToGenes(object, tx2gene)

# S4 method for matrix
convertTranscriptsToGenes(object, tx2gene, aggregate = TRUE)

# S4 method for sparseMatrix
convertTranscriptsToGenes(object, tx2gene, aggregate = TRUE)

# S4 method for SummarizedExperiment
convertTranscriptsToGenes(object)

Arguments

object

Object.

tx2gene

Tx2Gene. Transcript-to-gene mappings.

aggregate

logical(1). For objects supporting dim(), aggregate counts to gene level and collapse the matrix.

...

Additional arguments.

Value

  • character: factor. Genes in the values, transcripts in the names.

  • matrix, sparseMatrix, SummarizedExperiment: Object containing counts collapsed to gene level by default (see aggregate argument).

Note

For objects containing a count matrix, the object rows will be collapsed to gene level using aggregateRows. This applies to our SummarizedExperiment method.

Updated 2019-08-11.

See also

Examples

data(SummarizedExperiment_transcripts, package = "acidtest") txse <- SummarizedExperiment_transcripts object <- txse t2g <- Tx2Gene(object) print(t2g)
#> Tx2Gene with 6 rows and 2 columns #> transcriptID geneID #> <character> <character> #> ENST00000494424 ENST00000494424 ENSG00000000003 #> ENST00000496771 ENST00000496771 ENSG00000000003 #> ENST00000612152 ENST00000612152 ENSG00000000003 #> ENST00000371584 ENST00000371584 ENSG00000000419 #> ENST00000371588 ENST00000371588 ENSG00000000419 #> ENST00000413082 ENST00000413082 ENSG00000000419
transcripts <- rownames(object) print(transcripts)
#> [1] "ENST00000494424" "ENST00000496771" "ENST00000612152" "ENST00000371584" #> [5] "ENST00000371588" "ENST00000413082"
## character ==== ## Returns as factor. x <- convertTranscriptsToGenes(transcripts, tx2gene = t2g) print(x)
#> ENST00000494424 ENST00000496771 ENST00000612152 ENST00000371584 ENST00000371588 #> ENSG00000000003 ENSG00000000003 ENSG00000000003 ENSG00000000419 ENSG00000000419 #> ENST00000413082 #> ENSG00000000419 #> Levels: ENSG00000000003 ENSG00000000419
str(x)
#> Factor w/ 2 levels "ENSG00000000003",..: 1 1 1 2 2 2 #> - attr(*, "names")= chr [1:6] "ENST00000494424" "ENST00000496771" "ENST00000612152" "ENST00000371584" ...
## matrix ==== ## Note that transcript IDs currently must be in the rows. counts <- counts(object) print(counts)
#> sample1 sample2 sample3 sample4 #> ENST00000494424 1 2 3 4 #> ENST00000496771 5 6 7 8 #> ENST00000612152 9 10 11 12 #> ENST00000371584 13 14 15 16 #> ENST00000371588 17 18 19 20 #> ENST00000413082 21 22 23 24
## Aggregate to gene level. x <- convertTranscriptsToGenes(counts, tx2gene = t2g, aggregate = TRUE)
#> Aggregating counts using 'sum()'. #> Groupings: #> ENST00000494424 ENST00000496771 ENST00000612152 ENST00000371584 ENST00000371588 #> ENSG00000000003 ENSG00000000003 ENSG00000000003 ENSG00000000419 ENSG00000000419 #> ENST00000413082 #> ENSG00000000419 #> Levels: ENSG00000000003 ENSG00000000419
#> sample1 sample2 sample3 sample4 #> ENSG00000000003 15 18 21 24 #> ENSG00000000419 51 54 57 60
#> sample1 sample2 sample3 sample4 #> 66 72 78 84
## Simply map to rownames. x <- convertTranscriptsToGenes(counts, tx2gene = t2g, aggregate = FALSE) print(x)
#> sample1 sample2 sample3 sample4 #> ENSG00000000003 1 2 3 4 #> ENSG00000000003 5 6 7 8 #> ENSG00000000003 9 10 11 12 #> ENSG00000000419 13 14 15 16 #> ENSG00000000419 17 18 19 20 #> ENSG00000000419 21 22 23 24
#> sample1 sample2 sample3 sample4 #> 66 72 78 84
## SummarizedExperiment ==== x <- convertTranscriptsToGenes(object)
#> Aggregating counts using 'sum()'. #> Groupings: #> ENST00000494424 ENST00000496771 ENST00000612152 ENST00000371584 ENST00000371588 #> ENSG00000000003 ENSG00000000003 ENSG00000000003 ENSG00000000419 ENSG00000000419 #> ENST00000413082 #> ENSG00000000419 #> Levels: ENSG00000000003 ENSG00000000419
#> class: SummarizedExperiment #> dim: 2 4 #> metadata(0): #> assays(1): counts #> rownames(2): ENSG00000000003 ENSG00000000419 #> rowData names(0): #> colnames(4): sample1 sample2 sample3 sample4 #> colData names(0):