2.1 Pre-requisites

• Basic knowledge of R syntax
• Familiarity with the GRanges and SummarizedExperiment classes
• Familiarity with ’omics data types including copy number and gene expression

2.2 Workshop Participation

Participants will have a chance to build a MultiAssayExperiment object from scratch, and will also work with more complex objects provided by the curatedTCGAData package.

2.3 R/Bioconductor packages used

To install the workshop dependencies, assuming you have already installed Bioconductor:

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

## Bioconductor 3.8 Stable Release
BiocManager::install(version = "3.8")

BiocManager::install(
    c("org.Hs.eg.db", "mirbase.db", "EnsDb.Hsapiens.v86")
)

BiocManager::install(
    "waldronlab/curatedTCGAWorkshop",
    dependencies = TRUE,
    build_vignettes = TRUE
)

Note, the first line installs experimental data packages not installed automatically by BiocManager.

Here is a list of Bioconductor packages required for this workshop:

• curatedTCGAData
• MultiAssayExperiment
• SummarizedExperiment
• GenomicRanges
• RaggedExperiment
• TCGAutils
• UpSetR
• mirbase.db
• EnsDb.Hsapiens.v86
• org.Hs.eg.db

library(curatedTCGAData)
library(MultiAssayExperiment)
library(SummarizedExperiment)
library(GenomicRanges)
library(RaggedExperiment)
library(GenomicDataCommons)
library(TCGAutils)
library(UpSetR)
library(mirbase.db)
library(EnsDb.Hsapiens.v86)
library(org.Hs.eg.db)
library(readxl)
library(dplyr)
library(kableExtra)

2.4 Time outline

2h total

3.1 Learning goals

• identify advantages of providing integrative data in a analysis-ready platform
• overview common data classes available in Bioconductor
• gain familiarity with available functionality in TCGAutils for the management and coordination of TGCA data

3.2 Learning objectives

• use curatedTCGAData to create on-the-fly TCGA MultiAssayExperiment datasets
• create a MultiAssayExperiment for TCGA or other multi’omics data
• explore functionality available in TCGAutils with curatedTCGAData objects

4.1 curatedTCGAData

Many tools exist for accessing and downloading The Cancer Genome Atlas (TCGA) data. These include but are not limited to RTCGAToolbox, GenomicDataCommons (package and website), TCGAbiolinks, cBioPortal, and Broad GDAC Firehose. These tools encompass a spectrum of strengths in ease-of-use, integration, and completeness of data. Few tools provide an integrative and user-friendly representation of TCGA data in a widely used analysis platform such as Bioconductor.

The curatedTCGAData experiment data package provides on-the-fly construction of TCGA datasets for 33 different cancer types from the Broad GDAC Firehose pipeline. Generally, it provides data using build hg19. curatedTCGAData facilitates access and integration of TCGA data by providing multi-’omics data objects using the MultiAssayExperiment data class. Where other platforms provide fragmented datasets, curatedTCGAData ensures that all data provided is matched and accounted for within the phenotypic metadata.

A list of available cancer types can be obtained from TCGAutils::diseaseCodes:

data("diseaseCodes")
knitr::kable(diseaseCodes, align = "l", escape = FALSE,
    caption = "List of available cancer types from curatedTCGAData") %>%
    kableExtra::kable_styling(
        bootstrap_options = c("hover", "striped", "responsive"),
        full_width = FALSE
    )

A descriptive table of available ’omics types in curatedTCGAData (thanks to Ludwig G. @lgeistlinger):

## from dataTypesHTML.R
knitr::kable(dataTypes, align = "l", escape = FALSE,
    caption = "Descriptions of data types available by Bioconductor data class") %>%
    kable_styling(bootstrap_options = c("hover", "striped", "responsive"),
        full_width = FALSE) %>% group_rows(index = setNames(iwidth, groups)) %>%
    footnote(symbol = footnote)

4.2 TCGAutils

TCGAutils is a companion package that enhances curatedTCGAData by allowing additional exploration and manipulation of samples and metadata in rows and columns.

Available operations in TCGAutils and MultiAssayExperiment enable user-friendly operations for subsetting, separating, converting, and reshaping of sample and feature TCGA data.

TCGAutils was developed specifically for TCGA data and for curatedTCGAData products. It provides convenience / helper functions in three major areas:

1. conversion / summarization of row annotations to genomic ranges
2. identification and separation of samples
3. translation and interpretation of TCGA identifiers

along with additional reference data sets to explore TCGA data.

## from TCGAutils_cheatHTML.R script
knitr::kable(tabInfo, align = "l", escape = FALSE,
    caption = "Summary of available functionality in TCGAutils") %>%
    kable_styling(bootstrap_options = c("hover", "striped", "responsive"),
        full_width = FALSE) %>% group_rows(index = setNames(iwidth, groups))

To better understand how it all fits together, this schematic shows the relationship among all as part of the curatedTCGAData pipeline.

Figure 1: Schematic of curatedTCGAData Pipeline

6.1 RaggedExperiment

RaggedExperiment is a flexible data representation for segmented copy number, somatic mutations such as represented in .vcf files, and other ragged array schema for genomic location data. Like the GRangesList class from GenomicRanges, RaggedExperiment can be used to represent differing genomic ranges on each of a set of samples. In fact, RaggedExperiment contains a GRangesList:

showClass("RaggedExperiment")
#> Class "RaggedExperiment" [package "RaggedExperiment"]
#> 
#> Slots:
#>                                                       
#> Name:       assays      rowidx      colidx    metadata
#> Class: GRangesList     integer     integer        list
#> 
#> Extends: "Annotated"

However, RaggedExperiment provides a flexible set of Assay methods to support transformation of such data to matrix format.

Figure 3: RaggedExperiment object schematic
Rows and columns represent genomic ranges and samples, respectively. Assay operations can be performed with (from left to right) compactAssay, qreduceAssay, and sparseAssay.

6.2 MultiAssayExperiment

MultiAssayExperiment is an integrative container for coordinating multi-omics experiment data on a set of biological specimens. As much as possible, its methods adopt the same vocabulary as SummarizedExperiment. A MultiAssayExperiment can contain any number of assays with different representations. Assays may be ID-based, where measurements are indexed identifiers of genes, microRNA, proteins, microbes, etc. Alternatively, assays may be range-based, where measurements correspond to genomic ranges that can be represented as GRanges objects, such as gene expression or copy number. For ID-based assays, there is no requirement that the same IDs be present for different experiments. For range-based assays, there is also no requirement that the same ranges be present for different experiments; furthermore, it is possible for different samples within an experiment to be represented by different ranges. The following data classes have been tested to work as elements of a MultiAssayExperiment:

1. matrix: the most basic class for ID-based datasets, could be used for example for gene expression summarized per-gene, microRNA, metabolomics, or microbiome data.
2. SummarizedExperiment and derived methods: described above, could be used for miRNA, gene expression, proteomics, or any matrix-like data where measurements are represented by IDs.
3. RangedSummarizedExperiment: described above, could be used for gene expression, methylation, or other data types referring to genomic positions.
4. ExpressionSet: Another rich representation for ID-based datasets, supported only for legacy reasons
5. RaggedExperiment: described above, for non-rectangular (ragged) ranged-based datasets such as segmented copy number, where segmentation of copy number alterations occurs and different genomic locations in each sample.
6. RangedVcfStack: For VCF archives broken up by chromosome (see VcfStack class defined in the GenomicFiles package)
7. DelayedMatrix: An on-disk representation of matrix-like objects for large datasets. It reduces memory usage and optimizes performance with delayed operations. This class is part of the DelayedArray package.

Note that any data class extending these classes, and in fact any data class supporting row and column names and subsetting can be used as an element of a MultiAssayExperiment.

Figure 4: MultiAssayExperiment object schematic
colData provides data about the patients, cell lines, or other biological units, with one row per unit and one column per variable. The experiments are a list of assay datasets of arbitrary class. The sampleMap relates each column (observation) in ExperimentList to exactly one row (biological unit) in colData; however, one row of colData may map to zero, one, or more columns per assay, allowing for missing and replicate assays. sampleMap allows for per-assay sample naming conventions. Metadata can be used to store information in arbitrary format about the MultiAssayExperiment. Green stripes indicate a mapping of one subject to multiple observations across experiments.

7.1 API cheat sheet

Figure 5: The MultiAssayExperiment API for construction, access, subsetting, management, and reshaping to formats for application of R/Bioconductor graphics and analysis packages

8.1 Important clinical information

These objects contain most unrestricted TCGA assay and clinical / pathological data, as well as material curated from the supplements of published TCGA primary papers at the end of the colData columns:

dim(colData(acc))
#> [1]  92 822
tail(colnames(colData(acc)), 10)
#>  [1] "MethyLevel"       "miRNA.cluster"    "SCNA.cluster"    
#>  [4] "protein.cluster"  "COC"              "OncoSign"        
#>  [7] "purity"           "ploidy"           "genome_doublings"
#> [10] "ADS"

The TCGAutils::getClinicalNames function will display relevant clinical column names obtained from RTCGAToolbox that are also commonly found in other cancer types.

(acccol <- getClinicalNames("ACC"))
#>  [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] "gender"                              
#> [10] "date_of_initial_pathologic_diagnosis"
#> [11] "radiation_therapy"                   
#> [12] "histological_type"                   
#> [13] "residual_tumor"                      
#> [14] "number_of_lymph_nodes"               
#> [15] "race"                                
#> [16] "ethnicity"

all(acccol %in% names(colData(acc)))
#> [1] TRUE

head(colData(acc)[, acccol])
#> DataFrame with 6 rows and 16 columns
#>              years_to_birth vital_status days_to_death days_to_last_followup
#>                   <integer>    <integer>     <integer>             <integer>
#> TCGA-OR-A5J1             58            1          1355                    NA
#> TCGA-OR-A5J2             44            1          1677                    NA
#> TCGA-OR-A5J3             23            0            NA                  2091
#> TCGA-OR-A5J4             23            1           423                    NA
#> TCGA-OR-A5J5             30            1           365                    NA
#> TCGA-OR-A5J6             29            0            NA                  2703
#>              tumor_tissue_site pathologic_stage pathology_T_stage
#>                    <character>      <character>       <character>
#> TCGA-OR-A5J1           adrenal         stage ii                t2
#> TCGA-OR-A5J2           adrenal         stage iv                t3
#> TCGA-OR-A5J3           adrenal        stage iii                t3
#> TCGA-OR-A5J4           adrenal         stage iv                t3
#> TCGA-OR-A5J5           adrenal        stage iii                t4
#> TCGA-OR-A5J6           adrenal         stage ii                t2
#>              pathology_N_stage      gender
#>                    <character> <character>
#> TCGA-OR-A5J1                n0        male
#> TCGA-OR-A5J2                n0      female
#> TCGA-OR-A5J3                n0      female
#> TCGA-OR-A5J4                n1      female
#> TCGA-OR-A5J5                n0        male
#> TCGA-OR-A5J6                n0      female
#>              date_of_initial_pathologic_diagnosis radiation_therapy
#>                                         <integer>       <character>
#> TCGA-OR-A5J1                                 2000                no
#> TCGA-OR-A5J2                                 2004                no
#> TCGA-OR-A5J3                                 2008                no
#> TCGA-OR-A5J4                                 2000                no
#> TCGA-OR-A5J5                                 2000                no
#> TCGA-OR-A5J6                                 2006                no
#>                                 histological_type residual_tumor
#>                                       <character>    <character>
#> TCGA-OR-A5J1 adrenocortical carcinoma- usual type             r0
#> TCGA-OR-A5J2 adrenocortical carcinoma- usual type             r2
#> TCGA-OR-A5J3 adrenocortical carcinoma- usual type             r0
#> TCGA-OR-A5J4 adrenocortical carcinoma- usual type             r2
#> TCGA-OR-A5J5 adrenocortical carcinoma- usual type             r2
#> TCGA-OR-A5J6 adrenocortical carcinoma- usual type             r0
#>              number_of_lymph_nodes                      race
#>                          <integer>               <character>
#> TCGA-OR-A5J1                    NA                     white
#> TCGA-OR-A5J2                     0                     white
#> TCGA-OR-A5J3                     0                     white
#> TCGA-OR-A5J4                     2                     white
#> TCGA-OR-A5J5                    NA                     white
#> TCGA-OR-A5J6                    NA black or african american
#>                       ethnicity
#>                     <character>
#> TCGA-OR-A5J1                 NA
#> TCGA-OR-A5J2 hispanic or latino
#> TCGA-OR-A5J3 hispanic or latino
#> TCGA-OR-A5J4 hispanic or latino
#> TCGA-OR-A5J5 hispanic or latino
#> TCGA-OR-A5J6 hispanic or latino

Reference documentation in Firehose Broad GDAC:

browseURL("https://dx.doi.org/10.7908/C1RX9BD4")

8.2 Subtype Information

Using TCGAutils, getSubtypeMap shows the user the column names associated with published subtype information.

getSubtypeMap(acc)
#>          ACC_annotations     ACC_subtype
#> 1             Patient_ID          SAMPLE
#> 2  histological_subtypes       Histology
#> 3          mrna_subtypes         C1A/C1B
#> 4          mrna_subtypes         mRNA_K4
#> 5                   cimp      MethyLevel
#> 6      microrna_subtypes   miRNA cluster
#> 7          scna_subtypes    SCNA cluster
#> 8       protein_subtypes protein cluster
#> 9   integrative_subtypes             COC
#> 10     mutation_subtypes        OncoSign

(subtypeCols <- getSubtypeMap(acc)[["ACC_subtype"]])
#>  [1] "SAMPLE"          "Histology"       "C1A/C1B"         "mRNA_K4"        
#>  [5] "MethyLevel"      "miRNA cluster"   "SCNA cluster"    "protein cluster"
#>  [9] "COC"             "OncoSign"

## for older versions of TCGAutils
subtypeCols <- gsub("SAMPLE", "patientID", subtypeCols)

colData(acc)[, make.names(subtypeCols)]
#> DataFrame with 92 rows and 10 columns
#>                 patientID   Histology     C1A.C1B
#>               <character> <character> <character>
#> TCGA-OR-A5J1 TCGA-OR-A5J1  Usual Type         C1A
#> TCGA-OR-A5J2 TCGA-OR-A5J2  Usual Type         C1A
#> TCGA-OR-A5J3 TCGA-OR-A5J3  Usual Type         C1A
#> TCGA-OR-A5J4 TCGA-OR-A5J4  Usual Type          NA
#> TCGA-OR-A5J5 TCGA-OR-A5J5  Usual Type         C1A
#> ...                   ...         ...         ...
#> TCGA-PK-A5H9 TCGA-PK-A5H9  Usual Type         C1B
#> TCGA-PK-A5HA TCGA-PK-A5HA  Usual Type         C1B
#> TCGA-PK-A5HB TCGA-PK-A5HB  Usual Type         C1A
#> TCGA-PK-A5HC TCGA-PK-A5HC  Usual Type          NA
#> TCGA-P6-A5OG TCGA-P6-A5OG          NA          NA
#>                                           mRNA_K4        MethyLevel
#>                                       <character>       <character>
#> TCGA-OR-A5J1 steroid-phenotype-high+proliferation         CIMP-high
#> TCGA-OR-A5J2 steroid-phenotype-high+proliferation          CIMP-low
#> TCGA-OR-A5J3               steroid-phenotype-high CIMP-intermediate
#> TCGA-OR-A5J4                                   NA         CIMP-high
#> TCGA-OR-A5J5               steroid-phenotype-high CIMP-intermediate
#> ...                                           ...               ...
#> TCGA-PK-A5H9                steroid-phenotype-low          CIMP-low
#> TCGA-PK-A5HA                steroid-phenotype-low         CIMP-high
#> TCGA-PK-A5HB               steroid-phenotype-high         CIMP-high
#> TCGA-PK-A5HC                                   NA                NA
#> TCGA-P6-A5OG                                   NA                NA
#>              miRNA.cluster SCNA.cluster protein.cluster         COC
#>                <character>  <character>       <integer> <character>
#> TCGA-OR-A5J1       miRNA_1        Quiet              NA        COC3
#> TCGA-OR-A5J2       miRNA_1        Noisy               1        COC3
#> TCGA-OR-A5J3       miRNA_6  Chromosomal               3        COC2
#> TCGA-OR-A5J4       miRNA_6  Chromosomal              NA          NA
#> TCGA-OR-A5J5       miRNA_2  Chromosomal              NA        COC2
#> ...                    ...          ...             ...         ...
#> TCGA-PK-A5H9       miRNA_5        Quiet               3        COC1
#> TCGA-PK-A5HA       miRNA_5  Chromosomal               2        COC1
#> TCGA-PK-A5HB       miRNA_6        Noisy              NA        COC3
#> TCGA-PK-A5HC            NA  Chromosomal              NA          NA
#> TCGA-P6-A5OG            NA           NA              NA          NA
#>                 OncoSign
#>              <character>
#> TCGA-OR-A5J1         CN2
#> TCGA-OR-A5J2    TP53/NF1
#> TCGA-OR-A5J3         CN2
#> TCGA-OR-A5J4         CN1
#> TCGA-OR-A5J5    TP53/NF1
#> ...                  ...
#> TCGA-PK-A5H9    TP53/NF1
#> TCGA-PK-A5HA         CN2
#> TCGA-PK-A5HB    TP53/NF1
#> TCGA-PK-A5HC    TP53/NF1
#> TCGA-P6-A5OG          NA

8.3 Pan-Cancer Data

curatedTCGAData works to build datasets on-the-fly including Pan-Cancer datasets where assays from two or more particular cancer types need to be combined. Here we demonstrate a Pan-Cancer MultiAssayExperiment using ovarian and breast cancer data.

(
ovbrca <- curatedTCGAData(diseaseCode = c("OV", "BRCA"),
    assays = c("RNA*", "Mutation"), dry.run = FALSE)
)

#> A MultiAssayExperiment object of 6 listed
#>  experiments with user-defined names and respective classes. 
#>  Containing an ExperimentList class object of length 6: 
#>  [1] BRCA_Mutation-20160128: RaggedExperiment with 90490 rows and 993 columns 
#>  [2] BRCA_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 1212 columns 
#>  [3] BRCA_RNASeqGene-20160128: SummarizedExperiment with 20502 rows and 878 columns 
#>  [4] OV_Mutation-20160128: RaggedExperiment with 20219 rows and 316 columns 
#>  [5] OV_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 307 columns 
#>  [6] OV_RNASeqGene-20160128: SummarizedExperiment with 19990 rows and 299 columns 
#> Features: 
#>  experiments() - obtain the ExperimentList instance 
#>  colData() - the primary/phenotype DataFrame 
#>  sampleMap() - the sample availability DataFrame 
#>  `$`, `[`, `[[` - extract colData columns, subset, or experiment 
#>  *Format() - convert into a long or wide DataFrame 
#>  assays() - convert ExperimentList to a SimpleList of matrices

9.1 “Simplification” of curatedTCGAData objects

A number of helper functions are available for managing datasets from curatedTCGAData. These include:

mirToRanges(acc)
#> harmonizing input:
#>   removing 80 sampleMap rows not in names(experiments)
#> A MultiAssayExperiment object of 9 listed
#>  experiments with user-defined names and respective classes. 
#>  Containing an ExperimentList class object of length 9: 
#>  [1] ACC_CNVSNP-20160128: RaggedExperiment with 21052 rows and 180 columns 
#>  [2] ACC_GISTIC_AllByGene-20160128: SummarizedExperiment with 24776 rows and 90 columns 
#>  [3] ACC_GISTIC_Peaks-20160128: RangedSummarizedExperiment with 42 rows and 90 columns 
#>  [4] ACC_GISTIC_ThresholdedByGene-20160128: SummarizedExperiment with 24776 rows and 90 columns 
#>  [5] ACC_Mutation-20160128: RaggedExperiment with 20166 rows and 90 columns 
#>  [6] ACC_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 79 columns 
#>  [7] ACC_RPPAArray-20160128: SummarizedExperiment with 192 rows and 46 columns 
#>  [8] ACC_miRNASeqGene-20160128_ranged: RangedSummarizedExperiment with 1002 rows and 80 columns 
#>  [9] ACC_miRNASeqGene-20160128_unranged: SummarizedExperiment with 44 rows and 80 columns 
#> Features: 
#>  experiments() - obtain the ExperimentList instance 
#>  colData() - the primary/phenotype DataFrame 
#>  sampleMap() - the sample availability DataFrame 
#>  `$`, `[`, `[[` - extract colData columns, subset, or experiment 
#>  *Format() - convert into a long or wide DataFrame 
#>  assays() - convert ExperimentList to a SimpleList of matrices

symbolsToRanges(acc)
#> 'select()' returned 1:many mapping between keys and columns
#> 'select()' returned 1:1 mapping between keys and columns
#> 'select()' returned 1:many mapping between keys and columns
#> 'select()' returned 1:1 mapping between keys and columns
#> 'select()' returned 1:many mapping between keys and columns
#> 'select()' returned 1:1 mapping between keys and columns
#> harmonizing input:
#>   removing 259 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] ACC_CNVSNP-20160128: RaggedExperiment with 21052 rows and 180 columns 
#>  [2] ACC_GISTIC_Peaks-20160128: RangedSummarizedExperiment with 42 rows and 90 columns 
#>  [3] ACC_miRNASeqGene-20160128: SummarizedExperiment with 1046 rows and 80 columns 
#>  [4] ACC_Mutation-20160128: RaggedExperiment with 20166 rows and 90 columns 
#>  [5] ACC_RPPAArray-20160128: SummarizedExperiment with 192 rows and 46 columns 
#>  [6] ACC_GISTIC_AllByGene-20160128_ranged: RangedSummarizedExperiment with 19538 rows and 90 columns 
#>  [7] ACC_GISTIC_AllByGene-20160128_unranged: SummarizedExperiment with 5238 rows and 90 columns 
#>  [8] ACC_GISTIC_ThresholdedByGene-20160128_ranged: RangedSummarizedExperiment with 19538 rows and 90 columns 
#>  [9] ACC_GISTIC_ThresholdedByGene-20160128_unranged: SummarizedExperiment with 5238 rows and 90 columns 
#>  [10] ACC_RNASeq2GeneNorm-20160128_ranged: RangedSummarizedExperiment with 17527 rows and 79 columns 
#>  [11] ACC_RNASeq2GeneNorm-20160128_unranged: SummarizedExperiment with 2974 rows and 79 columns 
#> Features: 
#>  experiments() - obtain the ExperimentList instance 
#>  colData() - the primary/phenotype DataFrame 
#>  sampleMap() - the sample availability DataFrame 
#>  `$`, `[`, `[[` - extract colData columns, subset, or experiment 
#>  *Format() - convert into a long or wide DataFrame 
#>  assays() - convert ExperimentList to a SimpleList of matrices

genome(acc[["ACC_Mutation-20160128"]]) <-
    vapply(
        X = genome(acc[["ACC_Mutation-20160128"]]),
        FUN = TCGAutils::translateBuild,
        FUN.VALUE = character(1L)
    )

seqlevelsStyle(acc[["ACC_Mutation-20160128"]]) <- "UCSC"

seqlevelsStyle(acc[["ACC_Mutation-20160128"]])
#> [1] "UCSC"

rowRanges(acc[["ACC_Mutation-20160128"]])
#> GRanges object with 20166 ranges and 0 metadata columns:
#>           seqnames              ranges strand
#>              <Rle>           <IRanges>  <Rle>
#>       [1]     chr1            11561526      +
#>       [2]     chr1            12309384      +
#>       [3]     chr1            33820015      +
#>       [4]     chr1 152785074-152785097      +
#>       [5]     chr1           152800122      +
#>       ...      ...                 ...    ...
#>   [20162]     chr5 131007363-131007364      +
#>   [20163]     chr7   90894459-90894460      +
#>   [20164]     chr9 139581758-139581759      +
#>   [20165]    chr16   90095596-90095597      +
#>   [20166]    chr19   58385798-58385799      +
#>   -------
#>   seqinfo: 24 sequences from hg19 genome; no seqlengths

qreduceTCGA(acc[ , , "ACC_Mutation-20160128"])
#> 'select()' returned 1:1 mapping between keys and columns
#> harmonizing input:
#>   removing 655 sampleMap rows not in names(experiments)
#>   removing 2 colData rownames not in sampleMap 'primary'
#> Warning in .Seqinfo.mergexy(x, y): The 2 combined objects have no sequence levels in common. (Use
#>   suppressWarnings() to suppress this warning.)
#> harmonizing input:
#>   removing 90 sampleMap rows not in names(experiments)
#> A MultiAssayExperiment object of 1 listed
#>  experiment with a user-defined name and respective class. 
#>  Containing an ExperimentList class object of length 1: 
#>  [1] ACC_Mutation-20160128_simplified: RangedSummarizedExperiment with 22942 rows and 90 columns 
#> Features: 
#>  experiments() - obtain the ExperimentList instance 
#>  colData() - the primary/phenotype DataFrame 
#>  sampleMap() - the sample availability DataFrame 
#>  `$`, `[`, `[[` - extract colData columns, subset, or experiment 
#>  *Format() - convert into a long or wide DataFrame 
#>  assays() - convert ExperimentList to a SimpleList of matrices

(simpa <- TCGAutils::simplifyTCGA(acc))

sampleTables(acc)
#> $`ACC_CNVSNP-20160128`
#> 
#> 01 10 11 
#> 90 85  5 
#> 
#> $`ACC_GISTIC_AllByGene-20160128`
#> 
#> 01 
#> 90 
#> 
#> $`ACC_GISTIC_Peaks-20160128`
#> 
#> 01 
#> 90 
#> 
#> $`ACC_GISTIC_ThresholdedByGene-20160128`
#> 
#> 01 
#> 90 
#> 
#> $`ACC_miRNASeqGene-20160128`
#> 
#> 01 
#> 80 
#> 
#> $`ACC_Mutation-20160128`
#> 
#> 01 
#> 90 
#> 
#> $`ACC_RNASeq2GeneNorm-20160128`
#> 
#> 01 
#> 79 
#> 
#> $`ACC_RPPAArray-20160128`
#> 
#> 01 
#> 46

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

splitAssays(acc)
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Blood Derived Normal' samples
#> Selecting 'Solid Tissue Normal' samples
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Primary Solid Tumor' samples
#> Selecting 'Primary Solid Tumor' samples
#> A MultiAssayExperiment object of 10 listed
#>  experiments with user-defined names and respective classes. 
#>  Containing an ExperimentList class object of length 10: 
#>  [1] 01_ACC_CNVSNP-20160128: RaggedExperiment with 21052 rows and 90 columns 
#>  [2] 10_ACC_CNVSNP-20160128: RaggedExperiment with 21052 rows and 85 columns 
#>  [3] 11_ACC_CNVSNP-20160128: RaggedExperiment with 21052 rows and 5 columns 
#>  [4] 01_ACC_GISTIC_AllByGene-20160128: SummarizedExperiment with 24776 rows and 90 columns 
#>  [5] 01_ACC_GISTIC_Peaks-20160128: RangedSummarizedExperiment with 42 rows and 90 columns 
#>  [6] 01_ACC_GISTIC_ThresholdedByGene-20160128: SummarizedExperiment with 24776 rows and 90 columns 
#>  [7] 01_ACC_miRNASeqGene-20160128: SummarizedExperiment with 1046 rows and 80 columns 
#>  [8] 01_ACC_Mutation-20160128: RaggedExperiment with 20166 rows and 90 columns 
#>  [9] 01_ACC_RNASeq2GeneNorm-20160128: SummarizedExperiment with 20501 rows and 79 columns 
#>  [10] 01_ACC_RPPAArray-20160128: SummarizedExperiment with 192 rows and 46 columns 
#> Features: 
#>  experiments() - obtain the ExperimentList instance 
#>  colData() - the primary/phenotype DataFrame 
#>  sampleMap() - the sample availability DataFrame 
#>  `$`, `[`, `[[` - extract colData columns, subset, or experiment 
#>  *Format() - convert into a long or wide DataFrame 
#>  assays() - convert ExperimentList to a SimpleList of matrices

UUIDtoBarcode("ae55b2d3-62a1-419e-9f9a-5ddfac356db4", from_type = "case_id")
#>                                case_id submitter_id
#> 1 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 TCGA-B0-5117

(xcode <- head(colnames(acc)[["ACC_Mutation-20160128"]], 4))
#> [1] "TCGA-OR-A5J1-01A-11D-A29I-10" "TCGA-OR-A5J2-01A-11D-A29I-10"
#> [3] "TCGA-OR-A5J3-01A-11D-A29I-10" "TCGA-OR-A5J4-01A-11D-A29I-10"
barcodeToUUID(xcode)
#>           submitter_aliquot_ids                          aliquot_ids
#> 32 TCGA-OR-A5J1-01A-11D-A29I-10 352062e7-9b06-41cd-880c-38fe268c9bf3
#> 5  TCGA-OR-A5J2-01A-11D-A29I-10 d97c6076-7e4f-4dbe-85de-12bc0d84d8e8
#> 16 TCGA-OR-A5J3-01A-11D-A29I-10 7691e5bf-7a03-4e6a-9873-74f2c8390a41
#> 39 TCGA-OR-A5J4-01A-11D-A29I-10 4d76daa9-a336-435c-aeef-47344d97ac5b

head(UUIDtoUUID("ae55b2d3-62a1-419e-9f9a-5ddfac356db4", to_type = "file_id"))
#>                                case_id                        files.file_id
#> 1 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 d6625424-1503-4735-a7cc-2a3de606278f
#> 2 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 bc62b7e8-1b54-438b-944b-3e0655cdf6ac
#> 3 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 fad3f577-6e5e-4039-b35b-69269f42d488
#> 4 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 48c342b0-e7a2-4a7b-8556-55bcd8ad9ea0
#> 5 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 e8eb689b-c5e1-4a21-acec-5f2ac25dbf97
#> 6 ae55b2d3-62a1-419e-9f9a-5ddfac356db4 8f345db1-686b-4e72-9b90-704978cde526

library(GenomicDataCommons)
fquery <- files(legacy = TRUE) %>%
    GenomicDataCommons::filter(~ cases.project.project_id == "TCGA-ACC" &
        data_category == "Gene expression" &
        data_type == "Exon quantification")

fnames <- results(fquery)$file_name[1:6]

filenameToBarcode(fnames)

9.2 Other TCGA data types

Helper functions to add TCGA exon files (from legacy archive), copy number and GISTIC copy number calls to MultiAssayExperiment objects are also available in TCGAutils.

10.1 How many samples have data for each combination of assays?

Solution

The built-in upsetSamples creates an “upset” Venn diagram to answer this question:

upsetSamples(acc)

In this dataset only 44 samples have all 5 assays, 33 are missing RNA gene expression data, 2 are missing reverse-phase protein array (RPPA), 12 have only mutations and RPPA, etc.

10.2 Kaplan-meier plot stratified by pathology_N_stage

Create a Kaplan-meier plot, using pathology_N_stage as a stratifying variable.

Solution

The colData provides clinical data for things like a Kaplan-Meier plot for overall survival stratified by nodal stage.

library(survival)
Surv(acc$days_to_death, acc$vital_status)

And remove any patients missing overall survival information:

accsurv <- acc[, complete.cases(acc$days_to_death, acc$vital_status), ]

library(survminer)

fit <- survfit(Surv(days_to_death, vital_status) ~ pathology_N_stage, data = colData(accsurv))
ggsurvplot(fit, data = colData(accsurv), risk.table = TRUE)

10.3 Multivariate Cox regression including RNA-seq, copy number, and pathology

Choose the EZH2 gene for demonstration. This subsetting will drop assays with no row named EZH2:

wideacc <- wideFormat(acc["EZH2", , ],
    colDataCols=c("vital_status", "days_to_death", "pathology_N_stage"))

wideacc$y <- Surv(wideacc$days_to_death, wideacc$vital_status)

head(wideacc)

Perform a multivariate Cox regression with EZH2 copy number (gistict), log2-transformed EZH2 expression (RNASeq2GeneNorm), and nodal status (pathology_N_stage) as predictors:

coxph(
    Surv(days_to_death, vital_status) ~
        gistict_EZH2 + log2(RNASeq2GeneNorm_EZH2) + pathology_N_stage,
    data=wideacc
)

We see that EZH2 expression is significantly associated with overal survival (p < 0.001), but EZH2 copy number and nodal status are not. This analysis could easily be extended to the whole genome for discovery of prognostic features by repeated univariate regressions over columns, penalized multivariate regression, etc.

For further detail, see the main MultiAssayExperiment vignette.

back to top