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 (often genomic) data. A reason for this (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 class hierarchy, 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!)
Before we get started, you need to know how to install Bioconductor packages. The most important details are:
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.).
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 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 1
An 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 in the Ensembl
database for human. 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(EnsDb.Hsapiens.v86)
txdb <- EnsDb.Hsapiens.v86
g <- genes(txdb)
g <- 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)
se
## class: RangedSummarizedExperiment
## dim: 10 6
## metadata(0):
## assays(1): exprs
## rownames(10): ENSG00000223972 ENSG00000227232 ... ENSG00000238009 ENSG00000239945
## rowData names(6): gene_id gene_name ... symbol entrezid
## colnames: NULL
## colData names(3): sample condition treated
We 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 7 columns
## gene_id gene_name gene_biotype seq_coord_system symbol
## <character> <character> <character> <character> <character>
## ENSG00000223972 ENSG00000223972 DDX11L1 transcribed_unproces.. chromosome DDX11L1
## ENSG00000227232 ENSG00000227232 WASH7P unprocessed_pseudogene chromosome WASH7P
## ENSG00000278267 ENSG00000278267 MIR6859-1 miRNA chromosome MIR6859-1
## ENSG00000243485 ENSG00000243485 MIR1302-2 lincRNA chromosome MIR1302-2
## ENSG00000237613 ENSG00000237613 FAM138A lincRNA chromosome FAM138A
## ENSG00000268020 ENSG00000268020 OR4G4P unprocessed_pseudogene chromosome OR4G4P
## ENSG00000240361 ENSG00000240361 OR4G11P unprocessed_pseudogene chromosome OR4G11P
## ENSG00000186092 ENSG00000186092 OR4F5 protein_coding chromosome OR4F5
## ENSG00000238009 ENSG00000238009 RP11-34P13.7 lincRNA chromosome RP11-34P13.7
## ENSG00000239945 ENSG00000239945 RP11-34P13.8 lincRNA chromosome RP11-34P13.8
## entrezid score
## <list> <numeric>
## ENSG00000223972 100287596,100287102,727856,... -0.8416606
## ENSG00000227232 NA -0.4797526
## ENSG00000278267 102466751 0.2195530
## ENSG00000243485 105376912,100302278 1.1615049
## ENSG00000237613 654835,645520,641702 -0.5559217
## ENSG00000268020 NA -1.7221926
## ENSG00000240361 NA -2.2658853
## ENSG00000186092 79501 -0.4598345
## ENSG00000238009 NA -0.0350413
## ENSG00000239945 NA -0.3726303
Adding 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 1637581
## 2 2 A 1 1179905
## 3 3 B 0 1483914
## 4 4 B 1 1412396
## 5 5 C 0 1374390
## 6 6 C 1 1016842
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("1", IRanges(25000,40000))
se_sub <- se[overlapsAny(se, query), ]
We could have equivalently used the shorthand code:
se_sub <- se[se %over% query,]
We get just three ranges, and three rows of the SE:
rowRanges(se_sub)
## GRanges object with 3 ranges and 7 metadata columns:
## seqnames ranges strand | gene_id gene_name gene_biotype
## <Rle> <IRanges> <Rle> | <character> <character> <character>
## ENSG00000227232 1 14404-29570 - | ENSG00000227232 WASH7P unprocessed_pseudogene
## ENSG00000243485 1 29554-31109 + | ENSG00000243485 MIR1302-2 lincRNA
## ENSG00000237613 1 34554-36081 - | ENSG00000237613 FAM138A lincRNA
## seq_coord_system symbol entrezid score
## <character> <character> <list> <numeric>
## ENSG00000227232 chromosome WASH7P <NA> -0.479753
## ENSG00000243485 chromosome MIR1302-2 105376912,100302278 1.161505
## ENSG00000237613 chromosome FAM138A 654835,645520,641702 -0.555922
## -------
## seqinfo: 25 sequences (1 circular) from GRCh38 genome
assay(se_sub)
## [,1] [,2] [,3] [,4] [,5] [,6]
## ENSG00000227232 -0.7110513 -0.01052055 0.05751345 -0.66014829 -0.1810667 0.5245639
## ENSG00000243485 -0.1407036 -0.25654069 -2.05784819 0.09087783 2.1176562 -1.7655651
## ENSG00000237613 -0.6278446 1.00544937 0.32292541 -0.25451136 -1.5731432 -1.9969867
Another 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 GRCh38 genome:
## seqnames seqlengths isCircular genome
## 1 248956422 FALSE GRCh38
## 10 133797422 FALSE GRCh38
## 11 135086622 FALSE GRCh38
## 12 133275309 FALSE GRCh38
## 13 114364328 FALSE GRCh38
## ... ... ... ...
## 8 145138636 FALSE GRCh38
## 9 138394717 FALSE GRCh38
## MT 16569 TRUE GRCh38
## X 156040895 FALSE GRCh38
## Y 57227415 FALSE GRCh38
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.
url <- "http://duffel.rail.bio/recount/SRP046226/rse_gene.Rdata"
file <- "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.)
library(here)
## here() starts at /Users/milove/teach/compbio/compbio_src
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 13002276
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: HRV16
rse$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 6
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 seqlengths
seqinfo(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 could look this up from
the project website.
The following code I use to clean up the condition and treatment variables:
library(magrittr)
rse$condition %<>% (function(x) {
factor(sub("-",".", sub("disease state: (.*)","\\1",x) ))
})
rse$treatment %<>% (function(x) factor(sub("treatment: (.*)","\\1",x)))
Now we have:
table(rse$condition, rse$treatment)
##
## HRV16 Vehicle
## asthmatic 6 6
## non.asthmatic 6 6
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 mode
We 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 Vehicle
This 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")
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")
sce
## class: SingleCellExperiment
## dim: 58037 24
## metadata(0):
## assays(1): counts
## rownames(58037): ENSG00000000003.14 ENSG00000000005.5 ... ENSG00000283698.1
## ENSG00000283699.1
## rowData names(3): gene_id bp_length symbol
## 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: scuttle
sce <- 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.4491155
Similarly, dedicated functions for transformations:
sce <- logNormCounts(sce)
assayNames(sce)
## [1] "counts" "logcounts"
And dedicated functions and new slots for reduced dimensions:
set.seed(1)
sce <- fixedPCA(sce, rank=5, subset.row=NULL)
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)
plotReducedDim(sce, "PCA", color_by="treatment")
If you are interested in combining the “tidy” dplyr style of interacting with datasets with SE there is tidySummarizedExperiment, a Bioconductor package that integrates this class with these familiar verbs. See for example:
https://stemangiola.github.io/tidySummarizedExperiment/
library(tidySummarizedExperiment)
se <- SummarizedExperiment(assay = list("exprs" = exprs),
colData = col_data)
# filter to samples in condition A
se %>% filter(condition == "A")
# total over all genes, per condition
se %>%
group_by(condition) %>%
summarize(total=sum(exprs))
# the features where mean expression > .25
se %>%
group_by(.feature) %>%
mutate(mean_exprs = mean(exprs)) %>%
filter(mean_exprs > .25)
sessionInfo()
## R version 4.2.1 (2022-06-23)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur ... 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats4 stats graphics grDevices datasets utils methods base
##
## other attached packages:
## [1] scater_1.25.4 ggplot2_3.3.6 scran_1.25.0
## [4] scuttle_1.7.4 SingleCellExperiment_1.19.0 pheatmap_1.0.12
## [7] DESeq2_1.37.4 magrittr_2.0.3 here_1.0.1
## [10] EnsDb.Hsapiens.v86_2.99.0 ensembldb_2.21.2 AnnotationFilter_1.21.0
## [13] GenomicFeatures_1.49.5 AnnotationDbi_1.59.1 SummarizedExperiment_1.27.1
## [16] Biobase_2.57.1 GenomicRanges_1.49.0 GenomeInfoDb_1.33.3
## [19] IRanges_2.31.0 S4Vectors_0.35.1 BiocGenerics_0.43.0
## [22] MatrixGenerics_1.9.1 matrixStats_0.62.0 pkgdown_2.0.5
## [25] testthat_3.1.4 rmarkdown_2.14 devtools_2.4.3
## [28] usethis_2.1.6
##
## loaded via a namespace (and not attached):
## [1] BiocFileCache_2.5.0 igraph_1.3.4 lazyeval_0.2.2
## [4] splines_4.2.1 BiocParallel_1.31.10 digest_0.6.29
## [7] htmltools_0.5.2 viridis_0.6.2 fansi_1.0.3
## [10] memoise_2.0.1 ScaledMatrix_1.5.0 cluster_2.1.3
## [13] limma_3.53.4 remotes_2.4.2 Biostrings_2.65.1
## [16] annotate_1.75.0 prettyunits_1.1.1 colorspace_2.0-3
## [19] ggrepel_0.9.1 blob_1.2.3 rappdirs_0.3.3
## [22] xfun_0.31 dplyr_1.0.9 callr_3.7.0
## [25] crayon_1.5.1 RCurl_1.98-1.7 jsonlite_1.8.0
## [28] genefilter_1.79.0 survival_3.3-1 glue_1.6.2
## [31] gtable_0.3.0 zlibbioc_1.43.0 XVector_0.37.0
## [34] DelayedArray_0.23.0 pkgbuild_1.3.1 BiocSingular_1.13.1
## [37] scales_1.2.0 edgeR_3.39.6 DBI_1.1.3
## [40] Rcpp_1.0.8.3 viridisLite_0.4.0 xtable_1.8-4
## [43] progress_1.2.2 dqrng_0.3.0 bit_4.0.4
## [46] rsvd_1.0.5 metapod_1.5.0 httr_1.4.3
## [49] RColorBrewer_1.1-3 ellipsis_0.3.2 pkgconfig_2.0.3
## [52] XML_3.99-0.10 farver_2.1.1 sass_0.4.1
## [55] dbplyr_2.2.1 locfit_1.5-9.5 utf8_1.2.2
## [58] tidyselect_1.1.2 labeling_0.4.2 rlang_1.0.3
## [61] munsell_0.5.0 tools_4.2.1 cachem_1.0.6
## [64] cli_3.3.0 generics_0.1.3 RSQLite_2.2.14
## [67] evaluate_0.15 stringr_1.4.0 fastmap_1.1.0
## [70] yaml_2.3.5 processx_3.6.1 knitr_1.39
## [73] bit64_4.0.5 fs_1.5.2 purrr_0.3.4
## [76] KEGGREST_1.37.2 sparseMatrixStats_1.9.0 xml2_1.3.3
## [79] biomaRt_2.53.2 brio_1.1.3 compiler_4.2.1
## [82] rstudioapi_0.13 beeswarm_0.4.0 filelock_1.0.2
## [85] curl_4.3.2 png_0.1-7 statmod_1.4.36
## [88] tibble_3.1.7 geneplotter_1.75.0 bslib_0.3.1
## [91] stringi_1.7.6 highr_0.9 ps_1.7.1
## [94] lattice_0.20-45 bluster_1.7.0 ProtGenerics_1.29.0
## [97] Matrix_1.4-1 vctrs_0.4.1 pillar_1.7.0
## [100] lifecycle_1.0.1 jquerylib_0.1.4 BiocNeighbors_1.15.1
## [103] bitops_1.0-7 irlba_2.3.5 rtracklayer_1.57.0
## [106] R6_2.5.1 BiocIO_1.7.1 gridExtra_2.3
## [109] vipor_0.4.5 sessioninfo_1.2.2 codetools_0.2-18
## [112] assertthat_0.2.1 pkgload_1.3.0 rprojroot_2.0.3
## [115] rjson_0.2.21 withr_2.5.0 GenomicAlignments_1.33.0
## [118] Rsamtools_2.13.3 GenomeInfoDbData_1.2.8 parallel_4.2.1
## [121] hms_1.1.1 grid_4.2.1 beachmat_2.13.4
## [124] DelayedMatrixStats_1.19.0 ggbeeswarm_0.6.0 restfulr_0.0.15