print(paste(
"R version:",
getRversion()
))[1] "R version: 4.6.0"print(paste(
"Bioc version:",
BiocManager::version()
))[1] "Bioc version: 3.23"Working with Bioconductor objects
Why use Bioconductor? From a user perspective, the answer is clear: because many statisticians, bioinformaticians, and computer scientists have spent time writing methods and algorithms specifically for biological/genomic data. A reason for this, and why many people have contributed to this project, is that there is a shared infrastructure for common data types. This infrastructure is built up of object classes. An example of a class is GRanges (stands for “genomic ranges”), which is a way to specify a set of ranges in a particular genome, e.g. from basepair 101 to basepair 200 on chromosome 1 of the human genome (version 38). What’s an object? Well everything in R is an object, but usually when we talk about Bioconductor objects, we mean data structures containing many attributes, so more complex than a vector or matrix. And the objects have specific methods that help you either access the information in the object, run analyses on the object, plot the object, etc. Bioconductor also allows for inheritance, which means that you can define a class of object that inherits the structure and methods of a superclass on which it depends. This last point is mostly important for people who are developing new software for Bioconductor (maybe that’s you!)
Getting started with Bioconductor
Before we get started, you need to know how to install Bioconductor packages. The most important details are:
Bioconductor is a package repository, like CRAN
All Bioconductor packages should be installed following the instructions here: https://bioconductor.org/install (the only real exception is if you want to obtain Linux binaries, but for working on Windows or Mac, you should stick to the above instructions)
Bioconductor packages are linked in their versions, both to each other and to the version of R
Bioconductor’s installation function will look up your version of R and give you the appropriate versions of Bioconductor packages
If you want the latest version of Bioconductor, you need to use the latest version of R
Our version of R/Bioconductor
We will be using this version of R and Bioconductor:
How do you know if a package is a Bioconductor package? For one thing, you can just google the package name and you’ll see either CRAN or Bioconductor as a first result (packages must be in one or the other, they are not allowed to be on both repositories). But also, you can use Bioconductor’s installation function to install any packages, even ones on CRAN. By the way, you can install multiple packages at once by making a string vector: BiocManager::install(c("foo","bar"))
Why all this stress on versioning? This is because the packages in Bioconductor are highly interdependent, and also some are very dependent on R internals. So that the project can guarantee the code will run and not give errors on many systems (Linux, Mac and Windows have support for the majority of Bioconductor packages), new development is locked into cycles, such that a release of Bioconductor shouldn’t contain any two packages which conflict and could potentially cause errors.
Details: of course, Bioconductor is also a project, made up of people. There is a core team which is supported by an NIH grant, and developers who contribute to the open source Bioconductor packages. There are also yearly conferences (one in US, one in Europe, and one in Asia, etc.).
Working with Bioconductor objects
We will introduce the core Bioconductor objects this week. In this particular document, we will discuss one of the most important classes of object, which is the SummarizedExperiment, or SE.
SEs have the structure:
- a matrix of data, rows are genomic features, and columns are samples
- a table of data about the samples (columns)
- a table of data about the features (rows)
A diagram of this 3-part structure can be found here.
In SE, the 3 parts of the object are called 1) assay, 2) colData and 3) rowData or rowRanges.
Note: There was a class of object that came before the SE, called the ExpressionSet, which was used primarily to store microarray data. Here we will skip over the ExpressionSet, and just look at SEs.
It helps to start by making a small toy SE, to see how the pieces come together. (Often you won’t make an SE manually, but it will be downloaded from an external source, or generated by a function that you call, e.g. tximeta or some other data loading function.)
library(SummarizedExperiment)
col_data <- data.frame(sample=factor(1:6),
condition=factor(c("A","A","B","B","C","C")),
treated=factor(rep(0:1,3)))
col_data sample condition treated
1 1 A 0
2 2 A 1
3 3 B 0
4 4 B 1
5 5 C 0
6 6 C 1An important aspect of SEs is that the rows can optionally correspond to particular set of GRanges, e.g. a row of an SE could give the number of RNA-seq reads that can be assigned to a particular gene, and the row could also have metadata in the 3rd slot including, e.g. location of the gene in the genome. In this case, we use the rowRanges slot to specify the information.
If we don’t have ranges, we can just put a table on the “side” of the SE by specifying rowData.
I will show in the example though how to provide rowRanges. Let’s use the first 10 genes from the UCSC hg38 knownGene annotation. The following code loads a database, pulls out all the genes (as GRanges), removes extra “non-standard” chromosomes, and then subsets to the first 10 genes.
library(TxDb.Hsapiens.UCSC.hg38.knownGene)
txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene
suppressMessages({
g <- sort(genes(txdb))
})
g <- GenomeInfoDb::keepStandardChromosomes(g, pruning.mode="coarse")
row_ranges <- g[1:10]We will make up some simulated “expression” measurements, and then store these in the SE. I call list so I can name the matrix, otherwise it would not be named.
exprs <- matrix(rnorm(6 * 10), ncol=6, nrow=10)
se <- SummarizedExperiment(assay = list("exprs" = exprs),
colData = col_data,
rowRanges = row_ranges)
seclass: RangedSummarizedExperiment
dim: 10 6
metadata(0):
assays(1): exprs
rownames(10): 100287102 107985730 ... 400728 643837
rowData names(1): gene_id
colnames: NULL
colData names(3): sample condition treatedWe see this object has one named matrix. The object could have multiple matrices (as long as these are the same shape). In that case you could access the first with assay and in general by name, e.g. assay(se, "exprs") or equivalently assays(se)[["exprs"]] .
assayNames(se)[1] "exprs"Finally, if we wanted to add data onto the rows, for example, the score of a test on the matrix data, we use the metadata columns function, or mcols:
mcols(se)$score <- rnorm(10)
mcols(se)DataFrame with 10 rows and 2 columns
gene_id score
<character> <numeric>
100287102 100287102 1.1418449
107985730 107985730 1.2985157
100302278 100302278 1.6422630
79501 79501 -0.2848453
102725121 102725121 2.6633842
124903815 124903815 0.0929827
105378580 105378580 -0.4121722
124903817 124903817 0.0185722
400728 400728 0.0207893
643837 643837 -1.2392003Adding data to the column metadata is even easier, we can just use $:
se$librarySize <- runif(6,1e6,2e6)
colData(se)DataFrame with 6 rows and 4 columns
sample condition treated librarySize
<factor> <factor> <factor> <numeric>
1 1 A 0 1176159
2 2 A 1 1081057
3 3 B 0 1275764
4 4 B 1 1198264
5 5 C 0 1807873
6 6 C 1 1631031Using the ranges of a SE
How does this additional functionality of the rowRanges facilitate faster data analysis? Suppose we are working with another data set besides se and we find a region of interest on chromsome 1. If we want to pull out the expression data for that region, we just ask for the subset of se that overlaps. First we build the query region, and then use the GRanges function overlapsAny within single square brackets (like you would subset any matrix-like object:
query <- GRanges("chr1", IRanges(5000,35000))
se_sub <- se[overlapsAny(se, query), ]We could have equivalently used the shorthand code:
se_sub <- se[se %over% query,]We get just the overlapping ranges, and the corresponding rows of the SE:
rowRanges(se_sub)GRanges object with 3 ranges and 2 metadata columns:
seqnames ranges strand | gene_id score
<Rle> <IRanges> <Rle> | <character> <numeric>
100287102 chr1 11121-14421 + | 100287102 1.14184
107985730 chr1 28977-31109 + | 107985730 1.29852
100302278 chr1 30366-30503 + | 100302278 1.64226
-------
seqinfo: 25 sequences (1 circular) from hg38 genomeassay(se_sub) [,1] [,2] [,3] [,4] [,5] [,6]
100287102 1.2415408 -2.8962065 -0.5820229 0.2291434 1.1412485 2.853652
107985730 2.0337348 0.6713310 0.5457085 2.0416815 1.5675024 1.422555
100302278 -0.3548704 0.1838906 -1.0552120 -0.8388115 0.6943347 2.125127Another useful property is that we know metadata about the chromosomes, and the version of the genome. (If you were not yet aware, the basepair position of a given feature, say gene XYZ, will change between versions of the genome, as sequences are added or rearranged.)
seqinfo(se)Seqinfo object with 25 sequences (1 circular) from hg38 genome:
seqnames seqlengths isCircular genome
chr1 248956422 FALSE hg38
chr2 242193529 FALSE hg38
chr3 198295559 FALSE hg38
chr4 190214555 FALSE hg38
chr5 181538259 FALSE hg38
... ... ... ...
chr21 46709983 FALSE hg38
chr22 50818468 FALSE hg38
chrX 156040895 FALSE hg38
chrY 57227415 FALSE hg38
chrM 16569 TRUE hg38Downloading SE data
We previously introduced the computational project, called recount2, which performs a basic summarization of public data sets with gene expression data. We will use data from recount2 again.
This dataset contains RNA-seq samples from human airway epithelial cell cultures. The paper is here. The structure of the experiment was that, cell cultures from 6 asthmatic and 6 non-asthmatics donors were treated with viral infection or left untreated (controls). So we have 2 samples (control or treated) for each of the 12 donors.
library(here)here() starts at /Users/love/teach/compbio/compbio_srcurl <- "http://duffel.rail.bio/recount/SRP046226/rse_gene.Rdata"
file <- here("bioc","asthma.rda")
if (!file.exists(file)) download.file(url, file)
load(file)We use a custom function to produce a matrix which a count of RNA fragments for each gene (rows) and each sample (columns).
(Recount project calls these objects rse for RangedSummarizedExperiment, meaning it has rowRanges information.)
source(here("bioc","my_scale_counts.R"))
rse <- my_scale_counts(rse_gene)We can take a peek at the column data:
colData(rse)[,1:6]DataFrame with 24 rows and 6 columns
project sample experiment run read_count_as_reported_by_sra
<character> <character> <character> <character> <integer>
SRR1565926 SRP046226 SRS694613 SRX692912 SRR1565926 12866750
SRR1565927 SRP046226 SRS694614 SRX692913 SRR1565927 12797108
SRR1565928 SRP046226 SRS694615 SRX692914 SRR1565928 13319016
SRR1565929 SRP046226 SRS694616 SRX692915 SRR1565929 13725752
SRR1565930 SRP046226 SRS694617 SRX692916 SRR1565930 10882416
... ... ... ... ... ...
SRR1565945 SRP046226 SRS694632 SRX692931 SRR1565945 13791854
SRR1565946 SRP046226 SRS694633 SRX692932 SRR1565946 13480842
SRR1565947 SRP046226 SRS694634 SRX692933 SRR1565947 13166594
SRR1565948 SRP046226 SRS694635 SRX692934 SRR1565948 13320398
SRR1565949 SRP046226 SRS694636 SRX692935 SRR1565949 13002276
reads_downloaded
<integer>
SRR1565926 12866750
SRR1565927 12797108
SRR1565928 13319016
SRR1565929 13725752
SRR1565930 10882416
... ...
SRR1565945 13791854
SRR1565946 13480842
SRR1565947 13166594
SRR1565948 13320398
SRR1565949 13002276Parsing sample information
The information we are interested in is contained in the characteristics column (which is a character list).
class(rse$characteristics)[1] "CompressedCharacterList"
attr(,"package")
[1] "IRanges"rse$characteristics[1:3]CharacterList of length 3
[[1]] cell type: Isolated from human trachea-bronchial tissues passages: 2 disease state: asthmatic treatment: HRV16
[[2]] cell type: Isolated from human trachea-bronchial tissues passages: 2 disease state: asthmatic treatment: HRV16
[[3]] cell type: Isolated from human trachea-bronchial tissues passages: 2 disease state: asthmatic treatment: HRV16rse$characteristics[[1]][1] "cell type: Isolated from human trachea-bronchial tissues"
[2] "passages: 2"
[3] "disease state: asthmatic"
[4] "treatment: HRV16" We can pull out the 3 and 4 element using the sapply function and the square bracket function. I know this syntax looks a little funny, but it’s really just saying, use the single square bracket, pull out the third element (or fourth element).
rse$condition <- sapply(rse$characteristics, `[`, 3)
rse$treatment <- sapply(rse$characteristics, `[`, 4)
table(rse$condition, rse$treatment)
treatment: HRV16 treatment: Vehicle
disease state: asthmatic 6 6
disease state: non-asthmatic 6 6The following code I use to clean up the condition and treatment variables:
library(stringr)
rse$condition <- rse$condition |>
str_remove("disease state: ") |>
str_replace("-", ".") |>
factor()
rse$treatment <- rse$treatment |>
str_remove("treatment: ") |>
factor()Now we have:
table(rse$condition, rse$treatment)
HRV16 Vehicle
asthmatic 6 6
non.asthmatic 6 6Looking at the genes and their ranges
Let’s see what the rowRanges of this experiment look like:
rowRanges(rse)GRanges object with 58037 ranges and 3 metadata columns:
seqnames ranges strand | gene_id bp_length
<Rle> <IRanges> <Rle> | <character> <integer>
ENSG00000000003.14 chrX 100627109-100639991 - | ENSG00000000003.14 4535
ENSG00000000005.5 chrX 100584802-100599885 + | ENSG00000000005.5 1610
ENSG00000000419.12 chr20 50934867-50958555 - | ENSG00000000419.12 1207
ENSG00000000457.13 chr1 169849631-169894267 - | ENSG00000000457.13 6883
ENSG00000000460.16 chr1 169662007-169854080 + | ENSG00000000460.16 5967
... ... ... ... . ... ...
ENSG00000283695.1 chr19 52865369-52865429 - | ENSG00000283695.1 61
ENSG00000283696.1 chr1 161399409-161422424 + | ENSG00000283696.1 997
ENSG00000283697.1 chrX 149548210-149549852 - | ENSG00000283697.1 1184
ENSG00000283698.1 chr2 112439312-112469687 - | ENSG00000283698.1 940
ENSG00000283699.1 chr10 12653138-12653197 - | ENSG00000283699.1 60
symbol
<CharacterList>
ENSG00000000003.14 TSPAN6
ENSG00000000005.5 TNMD
ENSG00000000419.12 DPM1
ENSG00000000457.13 SCYL3
ENSG00000000460.16 C1orf112
... ...
ENSG00000283695.1 <NA>
ENSG00000283696.1 <NA>
ENSG00000283697.1 LOC101928917
ENSG00000283698.1 <NA>
ENSG00000283699.1 MIR4481
-------
seqinfo: 25 sequences (1 circular) from an unspecified genome; no seqlengthsseqinfo(rse)Seqinfo object with 25 sequences (1 circular) from an unspecified genome; no seqlengths:
seqnames seqlengths isCircular genome
chr1 <NA> <NA> <NA>
chr2 <NA> <NA> <NA>
chr3 <NA> <NA> <NA>
chr4 <NA> <NA> <NA>
chr5 <NA> <NA> <NA>
... ... ... ...
chr21 <NA> <NA> <NA>
chr22 <NA> <NA> <NA>
chrX <NA> <NA> <NA>
chrY <NA> <NA> <NA>
chrM <NA> TRUE <NA>The rowRanges here were determined by the quantification method that the recount2 authors used. We don’t know what the genome is from the seqinfo, but we can see from the recount2 website the following information:
The RangedSummarizedExperiment object for the counts summarized at the gene level using the Gencode v25 (GRCh38.p7, CHR) annotation.
library(Seqinfo)
si_obj <- Seqinfo::Seqinfo(genome="hg38")
seqinfo(rse) <- si_obj[ seqlevels(rse) ]
seqinfo(rse)Seqinfo object with 25 sequences (1 circular) from hg38 genome:
seqnames seqlengths isCircular genome
chr1 248956422 FALSE hg38
chr2 242193529 FALSE hg38
chr3 198295559 FALSE hg38
chr4 190214555 FALSE hg38
chr5 181538259 FALSE hg38
... ... ... ...
chr21 46709983 FALSE hg38
chr22 50818468 FALSE hg38
chrX 156040895 FALSE hg38
chrY 57227415 FALSE hg38
chrM 16569 TRUE hg38Now with seqinfo we get warnings (as we should) for shifting features past the limit of the chromosome.
rowRanges(rse)[1]GRanges object with 1 range and 3 metadata columns:
seqnames ranges strand | gene_id bp_length
<Rle> <IRanges> <Rle> | <character> <integer>
ENSG00000000003.14 chrX 100627109-100639991 - | ENSG00000000003.14 4535
symbol
<CharacterList>
ENSG00000000003.14 TSPAN6
-------
seqinfo: 25 sequences (1 circular) from hg38 genomeseqinfo(rse)["chrX"]Seqinfo object with 1 sequence from hg38 genome:
seqnames seqlengths isCircular genome
chrX 156040895 FALSE hg38shift(rowRanges(rse)[1], 55e6)GRanges object with 1 range and 3 metadata columns:
seqnames ranges strand | gene_id bp_length
<Rle> <IRanges> <Rle> | <character> <integer>
ENSG00000000003.14 chrX 155627109-155639991 - | ENSG00000000003.14 4535
symbol
<CharacterList>
ENSG00000000003.14 TSPAN6
-------
seqinfo: 25 sequences (1 circular) from hg38 genomeshift(rowRanges(rse)[1], 57e6)Warning in valid.GenomicRanges.seqinfo(x, suggest.trim = TRUE): GRanges object contains 1 out-of-bound range located on sequence chrX. Note that ranges located
on a sequence whose length is unknown (NA) or on a circular sequence are not considered
out-of-bound (use seqlengths() and isCircular() to get the lengths and circularity flags of the
underlying sequences). You can use trim() to trim these ranges. See ?`trim,GenomicRanges-method`
for more information.
Warning in valid.GenomicRanges.seqinfo(x, suggest.trim = TRUE): GRanges object contains 1 out-of-bound range located on sequence chrX. Note that ranges located
on a sequence whose length is unknown (NA) or on a circular sequence are not considered
out-of-bound (use seqlengths() and isCircular() to get the lengths and circularity flags of the
underlying sequences). You can use trim() to trim these ranges. See ?`trim,GenomicRanges-method`
for more information.GRanges object with 1 range and 3 metadata columns:
seqnames ranges strand | gene_id bp_length
<Rle> <IRanges> <Rle> | <character> <integer>
ENSG00000000003.14 chrX 157627109-157639991 - | ENSG00000000003.14 4535
symbol
<CharacterList>
ENSG00000000003.14 TSPAN6
-------
seqinfo: 25 sequences (1 circular) from hg38 genomeAdding gene metadata from a GTF file
One piece of metadata we might want to add to the rows of our SE is the gene biotype, e.g. whether a gene is protein-coding, a lncRNA, a pseudogene, etc. This information is available in the GENCODE annotation GTF file. Here we use rtracklayer to import the GTF directly, keeping only the gene-level entries. We will see other ways to work with GTF files and annotation databases later in the Bioconductor section of the course. The GTF is ~37 Mb.
library(rtracklayer)
# GENCODE v25
gtf_url <- "https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_25/gencode.v25.annotation.gtf.gz"
gtf_file <- here("bioc", "gencode.v25.gtf.gz")
if (!file.exists(gtf_file)) download.file(gtf_url, gtf_file)
genes_gtf <- rtracklayer::import(gtf_file, feature.type="gene")The feature.type argument limits the import to gene-level rows, so we avoid parsing the full multi-million line file. The result is ranges with metadata columns from the GTF, including gene_type. Both the GTF and our rse use versioned GENCODE gene IDs (e.g. ENSG00000000003.14), so we can match them directly.
m <- match(rownames(rse), genes_gtf$gene_id)
mcols(rse)$gene_type <- genes_gtf$gene_type[m]
tab <- table(mcols(rse)$gene_type)
head(sort(tab, decreasing = TRUE),10)
protein_coding processed_pseudogene lincRNA antisense
19950 10275 7539 5530
unprocessed_pseudogene misc_RNA snRNA miRNA
2663 2213 1900 1569
TEC snoRNA
1048 944 Visualizing count matrix data in a SE
We will discuss transformations and normalization in a following section, but here we will just use a transformation so that we can compute meaningful distances on count data. We build a DESeqDataSet and then specify the experimental design using a ~ and the variables that we expect to produce differences in the counts. (These variables are used to assess how much technical variability is in the data, but not used in the transformation function itself.)
library(DESeq2)
dds <- DESeqDataSet(rse, ~condition + treatment)converting counts to integer modeWe use this function, which implements a variance stabilizing transformation (more on this next time):
vsd <- vst(dds)We calculate the variance across all samples (on the transformed data):
library(matrixStats)
rv <- rowVars(assay(vsd))
o <- order(rv, decreasing=TRUE)[1:100]Finally, before plotting a heatmap, we extract the covariates that we want to annotated the top of the plot.
anno_col <- as.data.frame(colData(vsd)[,c("condition","treatment")])
anno_col condition treatment
SRR1565926 asthmatic HRV16
SRR1565927 asthmatic HRV16
SRR1565928 asthmatic HRV16
SRR1565929 asthmatic HRV16
SRR1565930 asthmatic HRV16
SRR1565931 asthmatic HRV16
SRR1565932 asthmatic Vehicle
SRR1565933 asthmatic Vehicle
SRR1565934 asthmatic Vehicle
SRR1565935 asthmatic Vehicle
SRR1565936 asthmatic Vehicle
SRR1565937 asthmatic Vehicle
SRR1565938 non.asthmatic HRV16
SRR1565939 non.asthmatic HRV16
SRR1565940 non.asthmatic HRV16
SRR1565941 non.asthmatic HRV16
SRR1565942 non.asthmatic HRV16
SRR1565943 non.asthmatic HRV16
SRR1565944 non.asthmatic Vehicle
SRR1565945 non.asthmatic Vehicle
SRR1565946 non.asthmatic Vehicle
SRR1565947 non.asthmatic Vehicle
SRR1565948 non.asthmatic Vehicle
SRR1565949 non.asthmatic VehicleThis code pull out the top of the transformed data by variance, and adds an annotation to the top of the plot. By default the rows and columns will be clustered by Euclidean distance. See ?pheatmap for more details on this function (it’s a very detailed manual page).
library(pheatmap)
pheatmap(assay(vsd)[o,],
annotation_col=anno_col,
show_rownames=FALSE,
show_colnames=FALSE)
We can also easily make a PCA plot with dedicated functions:
plotPCA(vsd, intgroup="treatment")using ntop=500 top features by variance
SingleCellExperiment
An example of a class that extends the SE is SingleCellExperiment. This is a special object type for looking at single cell data.
For more details, there is a free online book “Orchestrating Single Cell Analysis With Bioconductor” produced by a group within the Bioconductor Project, with lots of example analyses: OSCA.
Here we show a quick example of how this object extends the SE.
library(SingleCellExperiment)
sce <- as(rse, "SingleCellExperiment")
sceclass: SingleCellExperiment
dim: 58037 24
metadata(0):
assays(1): counts
rownames(58037): ENSG00000000003.14 ENSG00000000005.5 ... ENSG00000283698.1
ENSG00000283699.1
rowData names(4): gene_id bp_length symbol gene_type
colnames(24): SRR1565926 SRR1565927 ... SRR1565948 SRR1565949
colData names(23): project sample ... condition treatment
reducedDimNames(0):
mainExpName: NULL
altExpNames(0):There are special functions dedicated to scaling the samples (we will discuss this technical aspect soon):
library(scran)Loading required package: scuttlesce <- computeSumFactors(sce)
sizeFactors(sce) [1] 0.7672143 0.8205514 0.8686567 0.9479224 0.6484723 0.9815079 1.0797070 1.0569889 1.4377886
[10] 0.9465292 1.4759422 1.2630195 0.8889808 1.0524670 0.9677885 0.8086102 0.8806503 0.8999780
[19] 0.9505805 1.0430322 1.2527967 0.9908707 0.5208294 1.4491155Similarly, dedicated functions for transformations (the deprecation warning is ignorable for now, as of 2026 functions for single cell transformation are being migrated to new packages…)
sce <- logNormCounts(sce)Warning in .local(x, ...): 'normalizeCounts' is deprecated.
Use 'scrapper::normalizeCounts' instead.
See help("Deprecated")assayNames(sce)[1] "counts" "logcounts"And dedicated functions and new slots for reduced dimensions. (Again ignore the deprecation warnings…)
set.seed(1)
sce <- fixedPCA(sce, rank=5, subset.row=NULL)Warning in fixedPCA(sce, rank = 5, subset.row = NULL): 'fixedPCA' is deprecated.
Use 'scrapper::runPca.se' instead.
See help("Deprecated")reducedDimNames(sce)[1] "PCA"We can manually get at the PCs:
pca <- reducedDim(sce, "PCA")
plot(pca[,1:2])
But we can more easily use dedicated visualization functions:
library(scater)Warning: package 'scater' was built under R version 4.6.1plotReducedDim(sce, "PCA", color_by="treatment")
Other specialized objects in Bioconductor
Some other specialized objects that build on the SummarizedExperiment include:
- SpatialExperiment1 with best practices outlined here
- QFeatures2 for high-throughput mass spectrometry
- MultiAssayExperiment3 for multiple experimental assays performed on an overlapping set of specimens, e.g. RNA-seq, copy number, DNA methylation, etc.
- MethylationArray4 for DNA methylation arrays
- TreeSummarizedExperiment5 for data with hierarchical or tree-like structure on the rows or columns
sessionInfo()R version 4.6.0 (2026-04-24)
Platform: aarch64-apple-darwin23
Running under: macOS Tahoe 26.5.2
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/4.6/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.6/Resources/lib/libRlapack.dylib; LAPACK version 3.12.1
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
time zone: America/New_York
tzcode source: internal
attached base packages:
[1] stats4 stats graphics grDevices datasets utils methods base
other attached packages:
[1] scater_1.40.2 ggplot2_4.0.3
[3] scran_1.40.0 scuttle_1.22.0
[5] SingleCellExperiment_1.34.0 pheatmap_1.0.13
[7] DESeq2_1.52.0 rtracklayer_1.72.0
[9] stringr_1.6.0 here_1.0.2
[11] TxDb.Hsapiens.UCSC.hg38.knownGene_3.22.0 GenomicFeatures_1.64.0
[13] AnnotationDbi_1.74.0 SummarizedExperiment_1.42.0
[15] Biobase_2.72.0 GenomicRanges_1.64.0
[17] Seqinfo_1.2.0 IRanges_2.46.0
[19] S4Vectors_0.50.1 BiocGenerics_0.58.1
[21] generics_0.1.4 MatrixGenerics_1.24.0
[23] matrixStats_1.5.0
loaded via a namespace (and not attached):
[1] DBI_1.3.0 bitops_1.1-0 gridExtra_2.3.1
[4] rlang_1.3.0 magrittr_2.0.5 otel_0.2.0
[7] compiler_4.6.0 RSQLite_3.53.3 png_0.1-9
[10] vctrs_0.7.3 pkgconfig_2.0.3 crayon_1.5.3
[13] fastmap_1.2.0 XVector_0.52.0 labeling_0.4.3
[16] Rsamtools_2.28.0 rmarkdown_2.31 ggbeeswarm_0.7.3
[19] UCSC.utils_1.8.0 bit_4.6.0 xfun_0.60
[22] bluster_1.22.0 cachem_1.1.0 beachmat_2.28.0
[25] cigarillo_1.2.1 GenomeInfoDb_1.48.0 jsonlite_2.0.0
[28] blob_1.3.0 DelayedArray_0.38.2 BiocParallel_1.46.0
[31] irlba_2.3.7 parallel_4.6.0 cluster_2.1.8.3
[34] R6_2.6.1 stringi_1.8.9 RColorBrewer_1.1-3
[37] limma_3.68.4 Rcpp_1.1.2 knitr_1.51
[40] Matrix_1.7-6 igraph_2.3.3 tidyselect_1.2.1
[43] viridis_0.6.5 dichromat_2.0-1 abind_1.4-8
[46] yaml_2.3.12 codetools_0.2-20 curl_7.1.0
[49] lattice_0.22-9 tibble_3.3.1 withr_3.0.3
[52] KEGGREST_1.52.2 S7_0.2.2 evaluate_1.0.5
[55] Biostrings_2.80.1 pillar_1.11.1 BiocManager_1.30.27
[58] rprojroot_2.1.1 RCurl_1.98-1.19 scales_1.4.0
[61] glue_1.8.1 metapod_1.20.0 tools_4.6.0
[64] BiocIO_1.22.0 BiocNeighbors_2.6.0 ScaledMatrix_1.20.0
[67] locfit_1.5-9.12 GenomicAlignments_1.48.0 XML_3.99-0.23
[70] cowplot_1.2.0 grid_4.6.0 edgeR_4.10.1
[73] beeswarm_0.4.0 BiocSingular_1.28.0 vipor_0.4.7
[76] restfulr_0.0.17 cli_3.6.6 rsvd_1.0.5
[79] viridisLite_0.4.3 S4Arrays_1.12.0 dplyr_1.2.1
[82] gtable_0.3.6 digest_0.6.39 ggrepel_0.9.8
[85] SparseArray_1.12.2 dqrng_0.4.1 rjson_0.2.23
[88] htmlwidgets_1.6.4 farver_2.1.2 memoise_2.0.1
[91] htmltools_0.5.9 lifecycle_1.0.5 httr_1.4.8
[94] statmod_1.5.2 bit64_4.8.2