TCGAutils: Helper functions for working with TCGA datasets
September 14, 2020
Contents
- 1 Overview
- 2 Installation
- 3
curatedTCGADatautility functions- 3.1 obtaining TCGA as
MultiAssayExperimentobjects fromcuratedTCGAData - 3.2
sampleTables: what sample types are present in the data? - 3.3
splitAssays: separate the data from different tissue types - 3.4
getSubtypeMap: manually curated molecular subtypes - 3.5
getClinicalNames: key “level 4” clinical & pathological data
- 3.1 obtaining TCGA as
- 4 Converting Assays to SummarizedExperiment
- 5 Importing TCGA text data files to Bioconductor classes
- 6 Translating and interpreting TCGA identifiers
- 7 OncoPrint - oncoPrintTCGA
- 8 Reference data
- 9
sessionInfo
1 Overview
The TCGAutils package completes a suite of Bioconductor packages for
convenient access, integration, and analysis of The Cancer Genome Atlas.
It includes:
0. helpers for working with TCGA through the Bioconductor packages
MultiAssayExperiment (for coordinated representation and
manipulation of multi-omits experiments) and curatedTCGAData,
which provides unrestricted TCGA data as MultiAssayExperiment objects,
0. helpers for importing TCGA data as from flat data structures such as
data.frame or DataFrame read from delimited data structures provided by
the Broad Institute’s Firehose, GenomicDataCommons, and
0. functions for interpreting TCGA barcodes and for mapping between
barcodes and Universally Unique Identifiers (UUIDs).
2 Installation
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("TCGAutils")Required packages for this vignette:
library(TCGAutils)
library(curatedTCGAData)
library(MultiAssayExperiment)
library(RTCGAToolbox)
library(BiocFileCache)
library(rtracklayer)
library(R.utils)3 curatedTCGAData utility functions
Functions such as getSubtypeMap and getClinicalNames provide information
on data inside a MultiAssayExperiment object downloaded from
curatedTCGAData. sampleTables and splitAssays support useful
operations on these MultiAssayExperiment objects.
3.1 obtaining TCGA as MultiAssayExperiment objects from curatedTCGAData
For demonstration we download part of the Colon Adenocarcinoma (COAD) dataset
usingcuratedTCGAData via ExperimentHub. This command will download any
data type that starts with CN* such as CNASeq:
coad <- curatedTCGAData::curatedTCGAData(diseaseCode = "COAD",
assays = c("CNASeq", "Mutation", "miRNA*",
"RNASeq2*", "mRNAArray", "Methyl*"), dry.run = FALSE)For a list of all available data types, use dry.run = FALSE and an
asterisk * as the assay input value:
curatedTCGAData("COAD", "*")## Title DispatchClass
## 79 COAD_CNASeq-20160128 Rda
## 80 COAD_CNASNP-20160128 Rda
## 81 COAD_CNVSNP-20160128 Rda
## 83 COAD_GISTIC_AllByGene-20160128 Rda
## 84 COAD_GISTIC_Peaks-20160128 Rda
## 85 COAD_GISTIC_ThresholdedByGene-20160128 Rda
## 87 COAD_Methylation_methyl27-20160128_assays H5File
## 88 COAD_Methylation_methyl27-20160128_se Rds
## 89 COAD_Methylation_methyl450-20160128_assays H5File
## 90 COAD_Methylation_methyl450-20160128_se Rds
## 91 COAD_miRNASeqGene-20160128 Rda
## 92 COAD_mRNAArray-20160128 Rda
## 93 COAD_Mutation-20160128 Rda
## 94 COAD_RNASeq2GeneNorm-20160128 Rda
## 95 COAD_RNASeqGene-20160128 Rda
## 96 COAD_RPPAArray-20160128 Rda3.2 sampleTables: what sample types are present in the data?
The sampleTables function gives a tally of available
samples in the dataset based on the TCGA barcode information.
sampleTables(coad)## $`COAD_CNASeq-20160128`
##
## 01 10 11
## 68 55 13
##
## $`COAD_miRNASeqGene-20160128`
##
## 01 02
## 220 1
##
## $`COAD_mRNAArray-20160128`
##
## 01 11
## 153 19
##
## $`COAD_Mutation-20160128`
##
## 01
## 154
##
## $`COAD_RNASeq2GeneNorm-20160128`
##
## 01
## 191
##
## $`COAD_Methylation_methyl27-20160128`
##
## 01 11
## 165 37
##
## $`COAD_Methylation_methyl450-20160128`
##
## 01 02 06 11
## 293 1 1 38For reference in interpreting the sample type codes, see the sampleTypes
table:
data("sampleTypes")
head(sampleTypes)## Code Definition Short.Letter.Code
## 1 01 Primary Solid Tumor TP
## 2 02 Recurrent Solid Tumor TR
## 3 03 Primary Blood Derived Cancer - Peripheral Blood TB
## 4 04 Recurrent Blood Derived Cancer - Bone Marrow TRBM
## 5 05 Additional - New Primary TAP
## 6 06 Metastatic TM3.3 splitAssays: separate the data from different tissue types
TCGA datasets include multiple -omics for solid tumors, adjacent normal
tissues, blood-derived cancers and normals, and other tissue types, which may
be mixed together in a single dataset. The MultiAssayExperiment object
generated here has one patient per row of its colData, but each patient may
have two or more -omics profiles by any assay, whether due to assaying of
different types of tissues or to technical replication. splitAssays separates
profiles from different tissue types (such as tumor and adjacent normal) into
different assays of the MultiAssayExperiment by taking a vector of sample
codes, and partitioning the current assays into assays with an appended sample
code:
(tnmae <- splitAssays(coad, c("01", "11")))## Warning: Some 'sampleCodes' not found in assays## A MultiAssayExperiment object of 11 listed
## experiments with user-defined names and respective classes.
## Containing an ExperimentList class object of length 11:
## [1] 01_COAD_CNASeq-20160128: RaggedExperiment with 40530 rows and 68 columns
## [2] 11_COAD_CNASeq-20160128: RaggedExperiment with 40530 rows and 13 columns
## [3] 01_COAD_miRNASeqGene-20160128: SummarizedExperiment with 705 rows and 220 columns
## [4] 01_COAD_mRNAArray-20160128: SummarizedExperiment with 17814 rows and 153 columns
## [5] 11_COAD_mRNAArray-20160128: SummarizedExperiment with 17814 rows and 19 columns
## [6] 01_COAD_Mutation-20160128: RaggedExperiment with 62530 rows and 154 columns
## [7] 01_COAD_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 191 columns
## [8] 01_COAD_Methylation_methyl27-20160128: SummarizedExperiment with 27578 rows and 165 columns
## [9] 11_COAD_Methylation_methyl27-20160128: SummarizedExperiment with 27578 rows and 37 columns
## [10] 01_COAD_Methylation_methyl450-20160128: SummarizedExperiment with 485577 rows and 293 columns
## [11] 11_COAD_Methylation_methyl450-20160128: SummarizedExperiment with 485577 rows and 38 columns
## Features:
## experiments() - obtain the ExperimentList instance
## colData() - the primary/phenotype DataFrame
## sampleMap() - the sample availability DFrame
## `$`, `[`, `[[` - extract colData columns, subset, or experiment
## *Format() - convert into a long or wide DataFrame
## assays() - convert ExperimentList to a SimpleList of matricesThe MultiAssayExperiment package then provides functionality to
merge replicate profiles for a single patient (mergeReplicates()), which
would now be appropriate but would not have been appropriate before
splitting different tissue types into different assays, because that would
average measurements from tumors and normal tissues.
MultiAssayExperiment also defines the MatchedAssayExperiment class, which
eliminates any profiles not present across all assays and ensures identical
ordering of profiles (columns) in each assay. In this example, it will match
tumors to adjacent normals in subsequent assays:
(matchmae <- as(tnmae[, , c(4, 6, 7)], "MatchedAssayExperiment"))## harmonizing input:
## removing 853 sampleMap rows not in names(experiments)
## removing 260 colData rownames not in sampleMap 'primary'## A MatchedAssayExperiment object of 3 listed
## experiments with user-defined names and respective classes.
## Containing an ExperimentList class object of length 3:
## [1] 01_COAD_mRNAArray-20160128: SummarizedExperiment with 17814 rows and 138 columns
## [2] 01_COAD_Mutation-20160128: RaggedExperiment with 62530 rows and 138 columns
## [3] 01_COAD_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 138 columns
## Features:
## experiments() - obtain the ExperimentList instance
## colData() - the primary/phenotype DataFrame
## sampleMap() - the sample availability DFrame
## `$`, `[`, `[[` - extract colData columns, subset, or experiment
## *Format() - convert into a long or wide DataFrame
## assays() - convert ExperimentList to a SimpleList of matricesOnly about 12 participants have both a matched tumor and solid normal sample.
3.4 getSubtypeMap: manually curated molecular subtypes
Per-tumor subtypes are saved in the metadata of the colData
slot of MultiAssayExperiment objects downloaded from curatedTCGAData.
These subtypes were manually curated from the supplemental tables of all
primary TCGA publications:
getSubtypeMap(coad)## COAD_annotations COAD_subtype
## 1 Patient_ID patientID
## 2 msi MSI_status
## 3 methylation_subtypes methylation_subtype
## 4 mrna_subtypes expression_subtype
## 5 histological_subtypes histological_type3.5 getClinicalNames: key “level 4” clinical & pathological data
The curatedTCGAData colData contain hundreds of columns, obtained from
merging all unrestricted levels of clinical, pathological, and biospecimen data.
This function provides the names of “level 4” clinical/pathological variables,
which are the only ones provided by most other TCGA analysis tools.
Users may then use these variable names for subsetting or analysis, and may
even want to subset the colData to only these commonly used variables.
getClinicalNames("COAD")## [1] "years_to_birth"
## [2] "vital_status"
## [3] "days_to_death"
## [4] "days_to_last_followup"
## [5] "tumor_tissue_site"
## [6] "pathologic_stage"
## [7] "pathology_T_stage"
## [8] "pathology_N_stage"
## [9] "pathology_M_stage"
## [10] "gender"
## [11] "date_of_initial_pathologic_diagnosis"
## [12] "days_to_last_known_alive"
## [13] "radiation_therapy"
## [14] "histological_type"
## [15] "residual_tumor"
## [16] "number_of_lymph_nodes"
## [17] "race"
## [18] "ethnicity"Warning: some names may not exactly match the colData names in the object
due to differences in variable types. These variables are kept separate and
differentiated with x and y. For example, vital_status in this case
corresponds to two different variables obtained from the pipeline. One variable
is interger type and the other character:
class(colData(coad)[["vital_status.x"]])## [1] "integer"class(colData(coad)[["vital_status.y"]])## [1] "character"table(colData(coad)[["vital_status.x"]])##
## 0 1
## 355 102table(colData(coad)[["vital_status.y"]])##
## DECEASED LIVING
## 22 179Such conflicts should be inspected in this manner, and conflicts resolved by choosing the more complete variable, or by treating any conflicting values as unknown (“NA”).
4 Converting Assays to SummarizedExperiment
This section gives an overview of the operations that can be performed on
a given set of metadata obtained particularly from data-rich objects such
as those obtained from curatedTCGAData. There are several operations that
work with microRNA, methylation, mutation, and assays that have gene symbol
annotations.
4.1 CpGtoRanges
Using the methylation annotations in
IlluminaHumanMethylation450kanno.ilmn12.hg19 and the minfi package, we
look up CpG probes and convert to genomic coordinates with CpGtoRanges.
The function provides two assays, one with mapped probes and the other with
unmapped probes. Excluding unmapped probes can be done by setting the
unmapped argument to FALSE. This will run for both types of methylation
data (27k and 450k).
methcoad <- CpGtoRanges(coad)## Setting options('download.file.method.GEOquery'='auto')## Setting options('GEOquery.inmemory.gpl'=FALSE)## harmonizing input:
## removing 535 sampleMap rows not in names(experiments)4.2 mirToRanges
microRNA assays obtained from curatedTCGAData have annotated sequences
that can be converted to genomic ranges using the mirbase.db package.
The function looks up all sequences and converts them to (‘hg19’) ranges.
For those rows that cannot be found, an ‘unranged’ assay is introduced
in the resulting MultiAssayExperiment object.
mircoad <- mirToRanges(coad)## harmonizing input:
## removing 221 sampleMap rows not in names(experiments)4.3 qreduceTCGA
The qreduceTCGA function converts RaggedExperiment mutation data objects
to RangedSummarizedExperiment using org.Hs.eg.db and the qreduceTCGA
utility function from RaggedExperiment to summarize ‘silent’ and ‘non-silent’
mutations based on a ‘Variant_Classification’ metadata column in the original
object.
It uses ‘hg19’ transcript database (‘TxDb’) package internally to summarize
regions using qreduceAssay. The current genome build (‘hg18’) in the data
must be translated to ‘hg19’.
In this example, we first translate the genome strings in the RaggedExperiment
COAD_Mutation-20160128 dataset using translateBuild. We do this for the
entire vector by using vapply.
rag <- "COAD_Mutation-20160128"
genome(coad[[rag]]) <- translateBuild(genome(coad[[rag]]))
seqlevelsStyle(coad[[rag]]) <- "UCSC"
genome(coad[[rag]])## chr1 chr2 chr3 chr4 chr5 chr6 chr7 chr8 chr9 chr10 chr11
## "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18"
## chr12 chr13 chr14 chr15 chr16 chr17 chr18 chr19 chr20 chr21 chr22
## "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18" "hg18"
## chrX chrY
## "hg18" "hg18"Now we use liftOver from rtracklayer to translate ‘hg18’ builds
to ‘hg19’ using the downloaded chain file.
lifturl <-
"http://hgdownload.cse.ucsc.edu/goldenpath/hg18/liftOver/hg18ToHg19.over.chain.gz"
bfc <- BiocFileCache()
qfile <- bfcquery(bfc, "18to19chain", exact = TRUE)[["rpath"]]
cfile <-
if (length(qfile) && file.exists(qfile)) {
bfcquery(bfc, "18to19chain", exact = TRUE)[["rpath"]]
} else {
bfcadd(bfc, "18to19chain", lifturl)
}
chainfile <- file.path(tempdir(), gsub("\\.gz", "", basename(cfile)))
R.utils::gunzip(cfile, destname = chainfile, remove = FALSE)
chain <- suppressMessages(
rtracklayer::import.chain(chainfile)
)
ranges19 <- rtracklayer::liftOver(rowRanges(coad[[rag]]), chain)The same can be done to convert hg19 to hg38 (the same build that the
Genomic Data Commons uses) with the corresponding chain file:
liftchain <-
"http://hgdownload.cse.ucsc.edu/goldenpath/hg19/liftOver/hg19ToHg38.over.chain.gz"
bfc <- BiocFileCache()
q38file <- bfcquery(bfc, "19to38chain", exact = TRUE)[["rpath"]]
c38file <-
if (length(q38file) && file.exists(q38file)) {
bfcquery(bfc, "19to38chain", exact = TRUE)[["rpath"]]
} else {
bfcadd(bfc, "19to38chain", liftchain)
}
cloc38 <- file.path(tempdir(), gsub("\\.gz", "", basename(c38file)))
R.utils::gunzip(c38file, destname = cloc38, remove = FALSE)
chain38 <- suppressMessages(
rtracklayer::import.chain(cloc38)
)
## then use the liftOver function using the 'chain38' object
## as above
ranges38 <- rtracklayer::liftOver(unlist(ranges19), chain38)This will give us a list of ranges, each element corresponding to a single row
in the RaggedExperiment. We remove rows that had no matches in the liftOver
process and replace the ranges in the original RaggedExperiment with the
replacement method. Finally, we put the RaggedExperiment object back into the
MultiAssayExperiment.
re19 <- coad[[rag]][as.logical(lengths(ranges19))]
ranges19 <- unlist(ranges19)
genome(ranges19) <- "hg19"
rowRanges(re19) <- ranges19
# replacement
coad[["COAD_Mutation-20160128"]] <- re19
rowRanges(re19)## GRanges object with 62523 ranges and 0 metadata columns:
## seqnames ranges strand
## <Rle> <IRanges> <Rle>
## [1] chr20 1552407-1552408 +
## [2] chr1 161736152-161736153 +
## [3] chr7 100685895 +
## [4] chr7 103824453 +
## [5] chr7 104783644 +
## ... ... ... ...
## [62519] chr9 36369716 +
## [62520] chr9 37692640 +
## [62521] chr9 6007456 +
## [62522] chrX 123785782 +
## [62523] chrX 51487184 +
## -------
## seqinfo: 24 sequences from hg19 genome; no seqlengthsNow that we have matching builds, we can finally run the qreduceTCGA function.
coad <- qreduceTCGA(coad, keep.assay = TRUE)## ## 403 genes were dropped because they have exons located on both strands
## of the same reference sequence or on more than one reference sequence,
## so cannot be represented by a single genomic range.
## Use 'single.strand.genes.only=FALSE' to get all the genes in a
## GRangesList object, or use suppressMessages() to suppress this message.## 'select()' returned 1:1 mapping between keys and columns## Warning in `seqlevelsStyle<-`(`*tmp*`, value = value): more than one seqlevels
## style supplied, using the 1st one only4.4 symbolsToRanges
In the cases where row annotations indicate gene symbols, the symbolsToRanges
utility function converts genes to genomic ranges and replaces existing
assays with RangedSummarizedExperiment objects. Gene annotations are given
as ‘hg19’ genomic regions.
symbolsToRanges(coad)## 403 genes were dropped because they have exons located on both strands
## of the same reference sequence or on more than one reference sequence,
## so cannot be represented by a single genomic range.
## Use 'single.strand.genes.only=FALSE' to get all the genes in a
## GRangesList object, or use suppressMessages() to suppress this message.## 'select()' returned 1:1 mapping between keys and columns## 403 genes were dropped because they have exons located on both strands
## of the same reference sequence or on more than one reference sequence,
## so cannot be represented by a single genomic range.
## Use 'single.strand.genes.only=FALSE' to get all the genes in a
## GRangesList object, or use suppressMessages() to suppress this message.## 'select()' returned 1:1 mapping between keys and columns## harmonizing input:
## removing 363 sampleMap rows not in names(experiments)## A MultiAssayExperiment object of 11 listed
## experiments with user-defined names and respective classes.
## Containing an ExperimentList class object of length 11:
## [1] COAD_CNASeq-20160128: RaggedExperiment with 40530 rows and 136 columns
## [2] COAD_miRNASeqGene-20160128: SummarizedExperiment with 705 rows and 221 columns
## [3] COAD_Mutation-20160128: RaggedExperiment with 62523 rows and 154 columns
## [4] COAD_Methylation_methyl27-20160128: SummarizedExperiment with 27578 rows and 202 columns
## [5] COAD_Methylation_methyl450-20160128: SummarizedExperiment with 485577 rows and 333 columns
## [6] COAD_Mutation-20160128_simplified: RangedSummarizedExperiment with 22929 rows and 154 columns
## [7] COAD_CNASeq-20160128_simplified: RangedSummarizedExperiment with 22929 rows and 136 columns
## [8] COAD_mRNAArray-20160128_ranged: RangedSummarizedExperiment with 14310 rows and 172 columns
## [9] COAD_mRNAArray-20160128_unranged: SummarizedExperiment with 3504 rows and 172 columns
## [10] COAD_RNASeq2GeneNorm-20160128_ranged: RangedSummarizedExperiment with 17296 rows and 191 columns
## [11] COAD_RNASeq2GeneNorm-20160128_unranged: SummarizedExperiment with 3205 rows and 191 columns
## Features:
## experiments() - obtain the ExperimentList instance
## colData() - the primary/phenotype DataFrame
## sampleMap() - the sample availability DFrame
## `$`, `[`, `[[` - extract colData columns, subset, or experiment
## *Format() - convert into a long or wide DataFrame
## assays() - convert ExperimentList to a SimpleList of matrices5 Importing TCGA text data files to Bioconductor classes
A few functions in the package accept either files or classes such as
data.frame and FirehoseGISTIC as input and return standard Bioconductor
classes.
5.1 makeGRangesListFromExonFiles
The GenomicDataCommons package can be used to obtain ‘legacy’ exon
quantification files via:
Note. File downloads disabled on Windows due to long file names.
library(GenomicDataCommons)## Loading required package: magrittr##
## Attaching package: 'magrittr'## The following object is masked from 'package:R.utils':
##
## extract## The following object is masked from 'package:R.oo':
##
## equals##
## Attaching package: 'GenomicDataCommons'## The following object is masked from 'package:matrixStats':
##
## count## The following object is masked from 'package:S4Vectors':
##
## expand## The following object is masked from 'package:stats':
##
## filterqueso <- files(legacy = TRUE) %>%
filter( ~ cases.project.project_id == "TCGA-COAD" &
data_category == "Gene expression" &
data_type == "Exon quantification")
gdc_set_cache(directory = tempdir())## GDC Cache directory set to: /tmp/RtmpXtOsmWWe then use makeGRangesListFromExonFiles to create a GRangesList from
vectors of file paths. There are options to provide file names when file names
are too long to download (Windows OS). The nrows argument only keeps the
first 5 rows in each of the files read in due to invalid character exon ranges.
## FALSE until gdcdata works
qu <- manifest(queso)
qq <- gdcdata(qu$id[1:4])
makeGRangesListFromExonFiles(qq, nrows = 4)Note GRangesList objects must be converted to RaggedExperiment
class to incorporate them into a MultiAssayExperiment.
5.1.1 Work around for long file names on Windows
Due to file name length, Windows may not be able to read / display all files.
The workaround uses the fileNames argument from a character vector of file
names and will convert them to TCGA barcodes.
## Load example file found in package
pkgDir <- system.file("extdata", package = "TCGAutils", mustWork = TRUE)
exonFile <- list.files(pkgDir, pattern = "cation\\.txt$", full.names = TRUE)
exonFile## [1] "/tmp/Rtmpe8RRub/Rinsta68684f733f/TCGAutils/extdata/bt.exon_quantification.txt"## We add the original file prefix to query for the UUID and get the
## TCGAbarcode
filePrefix <- "unc.edu.32741f9a-9fec-441f-96b4-e504e62c5362.1755371."
## Add actual file name manually
makeGRangesListFromExonFiles(exonFile,
fileNames = paste0(filePrefix, basename(exonFile)))## Parsed with column specification:
## cols(
## exon = col_character(),
## raw_counts = col_double(),
## median_length_normalized = col_double(),
## RPKM = col_double()
## )## GRangesList object of length 1:
## $`TCGA-AA-3678-01A-01R-0905-07`
## GRanges object with 100 ranges and 3 metadata columns:
## seqnames ranges strand | raw_counts median_length_normalized
## <Rle> <IRanges> <Rle> | <numeric> <numeric>
## [1] chr1 11874-12227 + | 4 0.492918
## [2] chr1 12595-12721 + | 2 0.341270
## [3] chr1 12613-12721 + | 2 0.398148
## [4] chr1 12646-12697 + | 2 0.372549
## [5] chr1 13221-14409 + | 39 0.632997
## ... ... ... ... . ... ...
## [96] chr1 881782-881925 - | 179 1
## [97] chr1 883511-883612 - | 151 1
## [98] chr1 883870-883983 - | 155 1
## [99] chr1 886507-886618 - | 144 1
## [100] chr1 887380-887519 - | 158 1
## RPKM
## <numeric>
## [1] 0.322477
## [2] 0.449436
## [3] 0.523655
## [4] 1.097661
## [5] 0.936105
## ... ...
## [96] 35.4758
## [97] 42.2492
## [98] 38.8033
## [99] 36.6933
## [100] 32.2085
## -------
## seqinfo: 1 sequence from an unspecified genome; no seqlengths5.2 makeGRangesListFromCopyNumber
Other processed, genomic range-based data from TCGA data can be imported using
makeGRangesListFromCopyNumber. This tab-delimited data file of copy number
alterations from bladder urothelial carcinoma (BLCA) was obtained from the
Genomic Data Commons and is included in TCGAUtils as an example:
grlFile <- system.file("extdata", "blca_cnaseq.txt", package = "TCGAutils")
grl <- read.table(grlFile)
head(grl)## Sample Chromosome Start End Num_Probes
## 1 TCGA-BL-A0C8-01A-11D-A10R-02 14 70362113 73912204 NA
## 2 TCGA-BL-A0C8-01A-11D-A10R-02 9 115609546 131133898 NA
## 5 TCGA-BL-A13I-01A-11D-A13U-02 13 19020028 49129100 NA
## 6 TCGA-BL-A13I-01A-11D-A13U-02 1 10208 246409808 NA
## 9 TCGA-BL-A13J-01A-11D-A10R-02 23 3119586 5636448 NA
## 10 TCGA-BL-A13J-01A-11D-A10R-02 7 10127 35776912 NA
## Segment_Mean
## 1 -0.182879931
## 2 0.039675162
## 5 0.002085552
## 6 -0.014224752
## 9 0.877072555
## 10 0.113873871makeGRangesListFromCopyNumber(grl, split.field = "Sample")## GRangesList object of length 116:
## $`TCGA-BL-A0C8-01A-11D-A10R-02`
## GRanges object with 2 ranges and 0 metadata columns:
## seqnames ranges strand
## <Rle> <IRanges> <Rle>
## [1] 14 70362113-73912204 *
## [2] 9 115609546-131133898 *
## -------
## seqinfo: 24 sequences from an unspecified genome; no seqlengths
##
## $`TCGA-BL-A13I-01A-11D-A13U-02`
## GRanges object with 2 ranges and 0 metadata columns:
## seqnames ranges strand
## <Rle> <IRanges> <Rle>
## [1] 13 19020028-49129100 *
## [2] 1 10208-246409808 *
## -------
## seqinfo: 24 sequences from an unspecified genome; no seqlengths
##
## $`TCGA-BL-A13J-01A-11D-A10R-02`
## GRanges object with 2 ranges and 0 metadata columns:
## seqnames ranges strand
## <Rle> <IRanges> <Rle>
## [1] 23 3119586-5636448 *
## [2] 7 10127-35776912 *
## -------
## seqinfo: 24 sequences from an unspecified genome; no seqlengths
##
## ...
## <113 more elements>makeGRangesListFromCopyNumber(grl, split.field = "Sample",
keep.extra.columns = TRUE)## GRangesList object of length 116:
## $`TCGA-BL-A0C8-01A-11D-A10R-02`
## GRanges object with 2 ranges and 2 metadata columns:
## seqnames ranges strand | Num_Probes Segment_Mean
## <Rle> <IRanges> <Rle> | <logical> <numeric>
## [1] 14 70362113-73912204 * | <NA> -0.1828799
## [2] 9 115609546-131133898 * | <NA> 0.0396752
## -------
## seqinfo: 24 sequences from an unspecified genome; no seqlengths
##
## $`TCGA-BL-A13I-01A-11D-A13U-02`
## GRanges object with 2 ranges and 2 metadata columns:
## seqnames ranges strand | Num_Probes Segment_Mean
## <Rle> <IRanges> <Rle> | <logical> <numeric>
## [1] 13 19020028-49129100 * | <NA> 0.00208555
## [2] 1 10208-246409808 * | <NA> -0.01422475
## -------
## seqinfo: 24 sequences from an unspecified genome; no seqlengths
##
## $`TCGA-BL-A13J-01A-11D-A10R-02`
## GRanges object with 2 ranges and 2 metadata columns:
## seqnames ranges strand | Num_Probes Segment_Mean
## <Rle> <IRanges> <Rle> | <logical> <numeric>
## [1] 23 3119586-5636448 * | <NA> 0.877073
## [2] 7 10127-35776912 * | <NA> 0.113874
## -------
## seqinfo: 24 sequences from an unspecified genome; no seqlengths
##
## ...
## <113 more elements>5.3 makeSummarizedExperimentFromGISTIC
This function is only used for converting the FirehoseGISTIC class of the
RTCGAToolbox package. It allows the user to obtain thresholded
by gene data, probabilities and peak regions.
tempDIR <- tempdir()
co <- getFirehoseData("COAD", clinical = FALSE, GISTIC = TRUE,
destdir = tempDIR)## gdac.broadinstitute.org_COAD-TP.CopyNumber_Gistic2.Level_4.2016012800.0.0.tar.gzselectType(co, "GISTIC")## Dataset:COAD## FirehoseGISTIC object, dim: 24776 454class(selectType(co, "GISTIC"))## [1] "FirehoseGISTIC"
## attr(,"package")
## [1] "RTCGAToolbox"makeSummarizedExperimentFromGISTIC(co, "Peaks")## class: RangedSummarizedExperiment
## dim: 66 451
## metadata(0):
## assays(1): ''
## rownames(66): 23 24 ... 65 66
## rowData names(12): rowRanges Unique.Name ... V461 type
## colnames(451): TCGA-3L-AA1B-01A-11D-A36W-01
## TCGA-4N-A93T-01A-11D-A36W-01 ... TCGA-T9-A92H-01A-11D-A36W-01
## TCGA-WS-AB45-01A-11D-A40O-01
## colData names(0):5.4 mergeColData: expanding the colData of a MultiAssayExperiment
This function merges a data.frame or DataFrame into the
colData of an existing MultiAssayExperiment object. It will match
column names and row names to do a full merge of both data sets. This
convenience function can be used, for example, to add subtype information
available for a subset of patients to the colData. Here is a simplified
example of adding a column to the colData DataFrame:
race_df <- DataFrame(race_f = factor(colData(coad)[["race"]]),
row.names = rownames(colData(coad)))
mergeColData(coad, race_df)## A MultiAssayExperiment object of 9 listed
## experiments with user-defined names and respective classes.
## Containing an ExperimentList class object of length 9:
## [1] COAD_CNASeq-20160128: RaggedExperiment with 40530 rows and 136 columns
## [2] COAD_miRNASeqGene-20160128: SummarizedExperiment with 705 rows and 221 columns
## [3] COAD_mRNAArray-20160128: SummarizedExperiment with 17814 rows and 172 columns
## [4] COAD_Mutation-20160128: RaggedExperiment with 62523 rows and 154 columns
## [5] COAD_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 191 columns
## [6] COAD_Methylation_methyl27-20160128: SummarizedExperiment with 27578 rows and 202 columns
## [7] COAD_Methylation_methyl450-20160128: SummarizedExperiment with 485577 rows and 333 columns
## [8] COAD_Mutation-20160128_simplified: RangedSummarizedExperiment with 22929 rows and 154 columns
## [9] COAD_CNASeq-20160128_simplified: RangedSummarizedExperiment with 22929 rows and 136 columns
## Features:
## experiments() - obtain the ExperimentList instance
## colData() - the primary/phenotype DFrame
## sampleMap() - the sample availability DFrame
## `$`, `[`, `[[` - extract colData columns, subset, or experiment
## *Format() - convert into a long or wide DFrame
## assays() - convert ExperimentList to a SimpleList of matrices6 Translating and interpreting TCGA identifiers
6.1 Translation
The TCGA project has generated massive amounts of data. Some data can be
obtained with Universally Unique IDentifiers (UUID) and other
data with TCGA barcodes. The Genomic Data Commons provides a JSON API for
mapping between UUID and barcode, but it is difficult for many people to
understand. TCGAutils makes simple functions available for two-way
translation between vectors of these identifiers.
6.1.1 TCGA barcode to UUID
Here we translate the first two TCGA barcodes of the previous copy-number alterations dataset to UUID:
(xbarcode <- head(colnames(coad)[["COAD_CNASeq-20160128_simplified"]], 4L))## [1] "TCGA-A6-2671-01A-01D-1405-02" "TCGA-A6-2671-10A-01D-1405-02"
## [3] "TCGA-A6-2674-01A-02D-1167-02" "TCGA-A6-2674-10A-01D-1167-02"barcodeToUUID(xbarcode)## submitter_aliquot_ids aliquot_ids
## 56 TCGA-A6-2671-01A-01D-1405-02 82e23baf-da11-4175-bee0-81c0c0137d72
## 63 TCGA-A6-2671-10A-01D-1405-02 da65c9d3-62ac-4fb5-b452-1e9c551ba243
## 10 TCGA-A6-2674-01A-02D-1167-02 9fdfc199-b878-4049-994e-b5ca384678fb
## 32 TCGA-A6-2674-10A-01D-1167-02 dd75656d-6df6-4c53-972d-439791c908ac6.1.2 UUID to TCGA barcode
Here we have a known case UUID that we want to translate into a TCGA barcode.
UUIDtoBarcode("ae55b2d3-62a1-419e-9f9a-5ddfac356db4", from_type = "case_id")## case_id submitter_id
## 1 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 TCGA-B0-5117In cases where we want to translate a known file UUID to the associated TCGA
patient barcode, we can use UUIDtoBarcode.
UUIDtoBarcode("0001801b-54b0-4551-8d7a-d66fb59429bf", from_type = "file_id")## file_id associated_entities.entity_submitter_id
## 1 0001801b-54b0-4551-8d7a-d66fb59429bf TCGA-B0-5094-11A-01D-1421-08Translating aliquot UUIDs is also possible by providing a known aliquot UUID to
the function and giving a from_type, “aliquot_ids”:
UUIDtoBarcode("d85d8a17-8aea-49d3-8a03-8f13141c163b", from_type = "aliquot_ids")## portions.analytes.aliquots.aliquot_id portions.analytes.aliquots.submitter_id
## 1 d85d8a17-8aea-49d3-8a03-8f13141c163b TCGA-CV-5443-01A-01D-1510-01Additional UUIDs may be supported in future versions.
6.1.3 UUID to UUID
We can also translate from file UUIDs to case UUIDs and vice versa as long as
we know the input type. We can use the case UUID from the previous example to
get the associated file UUIDs using UUIDtoUUID. Note that this translation
is a one to many relationship, thus yielding a data.frame of file UUIDs for a
single case UUID.
head(UUIDtoUUID("ae55b2d3-62a1-419e-9f9a-5ddfac356db4", to_type = "file_id"))## case_id files.file_id
## 1 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 48c342b0-e7a2-4a7b-8556-55bcd8ad9ea0
## 2 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 db8ba5d3-76be-4a67-a575-803ba483b6f9
## 3 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 f580489b-55ea-43c5-9489-b54c13146992
## 4 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 bf72ffef-d8c4-423d-9c5a-7bb5c23b2f31
## 5 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 b36f4e88-89ca-40bf-b543-d0e3c08ad342
## 6 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 4c3a899f-be0f-454f-b5dc-e30e29314c49One possible way to verify that file IDs are matching case UUIDS is to
browse to the Genomic Data Commons webpage with the specific file UUID.
Here we look at the first file UUID entry in the output data.frame:
https://portal.gdc.cancer.gov/files/0ff55a5e-6058-4e0b-9641-e3cb375ff214
In the page we check that the case UUID matches the input.
6.2 Parsing TCGA barcodes
Several functions exist for working with TCGA barcodes, the main function being
TCGAbarcode. It takes a TCGA barcode and returns information about
participant, sample, and/or portion.
## Return participant barcodes
TCGAbarcode(xbarcode, participant = TRUE)## [1] "TCGA-A6-2671" "TCGA-A6-2671" "TCGA-A6-2674" "TCGA-A6-2674"## Just return samples
TCGAbarcode(xbarcode, participant = FALSE, sample = TRUE)## [1] "01A" "10A" "01A" "10A"## Include sample data as well
TCGAbarcode(xbarcode, participant = TRUE, sample = TRUE)## [1] "TCGA-A6-2671-01A" "TCGA-A6-2671-10A" "TCGA-A6-2674-01A" "TCGA-A6-2674-10A"## Include portion and analyte data
TCGAbarcode(xbarcode, participant = TRUE, sample = TRUE, portion = TRUE)## [1] "TCGA-A6-2671-01A-01D" "TCGA-A6-2671-10A-01D" "TCGA-A6-2674-01A-02D"
## [4] "TCGA-A6-2674-10A-01D"6.3 Sample select
Based on lookup table values, the user can select certain sample types from a vector of sample barcodes. Below we select “Primary Solid Tumors” from a vector of barcodes, returning a logical vector identifying the matching samples.
## Select primary solid tumors
TCGAsampleSelect(xbarcode, "01")## 01 10 01 10
## TRUE FALSE TRUE FALSE## Select blood derived normals
TCGAsampleSelect(xbarcode, "10")## 01 10 01 10
## FALSE TRUE FALSE TRUE6.4 data.frame representation of barcode
The straightforward TCGAbiospec function will take the information contained
in the TCGA barcode and display it in data.frame format with appropriate
column names.
TCGAbiospec(xbarcode)## submitter_id sample_definition sample vial portion analyte plate center
## 1 TCGA-A6-2671 Primary Solid Tumor 01 A 01 D 1405 02
## 2 TCGA-A6-2671 Blood Derived Normal 10 A 01 D 1405 02
## 3 TCGA-A6-2674 Primary Solid Tumor 01 A 02 D 1167 02
## 4 TCGA-A6-2674 Blood Derived Normal 10 A 01 D 1167 027 OncoPrint - oncoPrintTCGA
We provide a convenience function that investigates metadata within
curatedTCGAData objects to present a plot of molecular alterations
within a paricular cancer. MultiAssayExperiment objects are required to
have an identifiable ‘Mutation’ assay (using text search). The variantCol
argument identifies the mutation type column within the data.
Note. Functionality streamlined from the ComplexHeatmap package.
oncoPrintTCGA(coad, matchassay = rag)## 403 genes were dropped because they have exons located on both strands
## of the same reference sequence or on more than one reference sequence,
## so cannot be represented by a single genomic range.
## Use 'single.strand.genes.only=FALSE' to get all the genes in a
## GRangesList object, or use suppressMessages() to suppress this message.## 'select()' returned 1:1 mapping between keys and columns## All mutation types: Frame Shift Del, Frame Shift Ins, Intron, Missense
## Mutation, Nonsense Mutation
8 Reference data
The TCGAutils package provides several helper datasets for working with TCGA barcodes.
8.1 sampleTypes
As shown previously, the reference dataset sampleTypes defines sample codes
and their sample types (see ?sampleTypes for source url).
## Obtained previously
sampleCodes <- TCGAbarcode(xbarcode, participant = FALSE, sample = TRUE)
## Lookup table
head(sampleTypes)## Code Definition Short.Letter.Code
## 1 01 Primary Solid Tumor TP
## 2 02 Recurrent Solid Tumor TR
## 3 03 Primary Blood Derived Cancer - Peripheral Blood TB
## 4 04 Recurrent Blood Derived Cancer - Bone Marrow TRBM
## 5 05 Additional - New Primary TAP
## 6 06 Metastatic TM## Match codes found in the barcode to the lookup table
sampleTypes[match(unique(substr(sampleCodes, 1L, 2L)), sampleTypes[["Code"]]), ]## Code Definition Short.Letter.Code
## 1 01 Primary Solid Tumor TP
## 10 10 Blood Derived Normal NBSource: https://gdc.cancer.gov/resources-tcga-users/tcga-code-tables/sample-type-codes
8.2 clinicalNames - Firehose pipeline clinical variables
clinicalNames is a list of the level 4 variable names (the most commonly used
clinical and pathological variables, with follow-ups merged) from each
colData datasets in curatedTCGAData. Shipped curatedTCGAData
MultiAssayExperiment objects merge additional levels 1-3 clinical,
pathological, and biospecimen data and contain many more variables than the ones
listed here.
data("clinicalNames")
clinicalNames## CharacterList of length 33
## [["ACC"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["BLCA"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["BRCA"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["CESC"]] years_to_birth vital_status ... age_at_diagnosis clinical_stage
## [["CHOL"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["COAD"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["DLBC"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["ESCA"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["GBM"]] years_to_birth vital_status days_to_death ... race ethnicity
## [["HNSC"]] years_to_birth vital_status days_to_death ... race ethnicity
## ...
## <23 more elements>lengths(clinicalNames)## ACC BLCA BRCA CESC CHOL COAD DLBC ESCA GBM HNSC KICH KIRC KIRP LAML LGG LIHC
## 16 18 17 48 16 18 11 19 12 19 18 19 19 9 12 16
## LUAD LUSC MESO OV PAAD PCPG PRAD READ SARC SKCM STAD TGCT THCA THYM UCEC UCS
## 20 20 17 12 19 13 18 18 12 17 17 15 21 11 9 11
## UVM
## 149 sessionInfo
sessionInfo()## R version 4.0.2 (2020-06-22)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.4 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.11-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.11-bioc/R/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] GenomicDataCommons_1.12.0 magrittr_1.5
## [3] rhdf5_2.32.2 R.utils_2.10.1
## [5] R.oo_1.24.0 R.methodsS3_1.8.1
## [7] rtracklayer_1.48.0 BiocFileCache_1.12.1
## [9] dbplyr_1.4.4 RTCGAToolbox_2.18.0
## [11] curatedTCGAData_1.10.0 MultiAssayExperiment_1.14.0
## [13] SummarizedExperiment_1.18.2 DelayedArray_0.14.1
## [15] matrixStats_0.56.0 Biobase_2.48.0
## [17] GenomicRanges_1.40.0 GenomeInfoDb_1.24.2
## [19] IRanges_2.22.2 S4Vectors_0.26.1
## [21] BiocGenerics_0.34.0 TCGAutils_1.8.1
## [23] BiocStyle_2.16.0
##
## loaded via a namespace (and not attached):
## [1] circlize_0.4.10
## [2] AnnotationHub_2.20.2
## [3] RCircos_1.2.1
## [4] plyr_1.8.6
## [5] splines_4.0.2
## [6] BiocParallel_1.22.0
## [7] digest_0.6.25
## [8] mirbase.db_1.2.0
## [9] foreach_1.5.0
## [10] htmltools_0.5.0
## [11] magick_2.4.0
## [12] memoise_1.1.0
## [13] cluster_2.1.0
## [14] limma_3.44.3
## [15] ComplexHeatmap_2.4.3
## [16] Biostrings_2.56.0
## [17] readr_1.3.1
## [18] annotate_1.66.0
## [19] askpass_1.1
## [20] siggenes_1.62.0
## [21] prettyunits_1.1.1
## [22] colorspace_1.4-1
## [23] blob_1.2.1
## [24] rvest_0.3.6
## [25] rappdirs_0.3.1
## [26] xfun_0.17
## [27] dplyr_1.0.2
## [28] crayon_1.3.4
## [29] RCurl_1.98-1.2
## [30] jsonlite_1.7.1
## [31] TxDb.Hsapiens.UCSC.hg19.knownGene_3.2.2
## [32] genefilter_1.70.0
## [33] GEOquery_2.56.0
## [34] RaggedExperiment_1.12.0
## [35] survival_3.2-3
## [36] iterators_1.0.12
## [37] glue_1.4.2
## [38] zlibbioc_1.34.0
## [39] XVector_0.28.0
## [40] GetoptLong_1.0.2
## [41] Rhdf5lib_1.10.1
## [42] shape_1.4.5
## [43] HDF5Array_1.16.1
## [44] DBI_1.1.0
## [45] rngtools_1.5
## [46] Rcpp_1.0.5
## [47] xtable_1.8-4
## [48] progress_1.2.2
## [49] clue_0.3-57
## [50] bumphunter_1.30.0
## [51] bit_4.0.4
## [52] mclust_5.4.6
## [53] preprocessCore_1.50.0
## [54] httr_1.4.2
## [55] RColorBrewer_1.1-2
## [56] ellipsis_0.3.1
## [57] pkgconfig_2.0.3
## [58] reshape_0.8.8
## [59] XML_3.99-0.5
## [60] locfit_1.5-9.4
## [61] RJSONIO_1.3-1.4
## [62] tidyselect_1.1.0
## [63] rlang_0.4.7
## [64] later_1.1.0.1
## [65] AnnotationDbi_1.50.3
## [66] BiocVersion_3.11.1
## [67] tools_4.0.2
## [68] generics_0.0.2
## [69] RSQLite_2.2.0
## [70] ExperimentHub_1.14.2
## [71] evaluate_0.14
## [72] stringr_1.4.0
## [73] fastmap_1.0.1
## [74] yaml_2.2.1
## [75] org.Hs.eg.db_3.11.4
## [76] knitr_1.29
## [77] bit64_4.0.5
## [78] beanplot_1.2
## [79] scrime_1.3.5
## [80] purrr_0.3.4
## [81] nlme_3.1-149
## [82] doRNG_1.8.2
## [83] mime_0.9
## [84] nor1mix_1.3-0
## [85] xml2_1.3.2
## [86] biomaRt_2.44.1
## [87] compiler_4.0.2
## [88] png_0.1-7
## [89] curl_4.3
## [90] interactiveDisplayBase_1.26.3
## [91] tibble_3.0.3
## [92] stringi_1.5.3
## [93] GenomicFeatures_1.40.1
## [94] minfi_1.34.0
## [95] IlluminaHumanMethylation450kanno.ilmn12.hg19_0.6.0
## [96] lattice_0.20-41
## [97] Matrix_1.2-18
## [98] multtest_2.44.0
## [99] vctrs_0.3.4
## [100] pillar_1.4.6
## [101] lifecycle_0.2.0
## [102] BiocManager_1.30.10
## [103] GlobalOptions_0.1.2
## [104] data.table_1.13.0
## [105] bitops_1.0-6
## [106] httpuv_1.5.4
## [107] R6_2.4.1
## [108] bookdown_0.20
## [109] promises_1.1.1
## [110] codetools_0.2-16
## [111] MASS_7.3-53
## [112] assertthat_0.2.1
## [113] rjson_0.2.20
## [114] openssl_1.4.2
## [115] GenomicAlignments_1.24.0
## [116] Rsamtools_2.4.0
## [117] GenomeInfoDbData_1.2.3
## [118] hms_0.5.3
## [119] quadprog_1.5-8
## [120] grid_4.0.2
## [121] tidyr_1.1.2
## [122] base64_2.0
## [123] rmarkdown_2.3
## [124] DelayedMatrixStats_1.10.1
## [125] illuminaio_0.30.0
## [126] shiny_1.5.0