cellxgenedp 1.2.2
This package is available in Bioconductor version 3.15 and later. The following code installs cellxgenedp as well as other packages required for this vignette.
pkgs <- c("cellxgenedp", "zellkonverter", "SingleCellExperiment", "HDF5Array")
required_pkgs <- pkgs[!pkgs %in% rownames(installed.packages())]
BiocManager::install(required_pkgs)
Use the following pkgs
vector to install from GitHub (latest,
unchecked, development version) instead
pkgs <- c(
"mtmorgan/cellxgenedp", "zellkonverter", "SingleCellExperiment", "HDF5Array"
)
Load the package into your current R session. We make extensive use of the dplyr packages, and at the end of the vignette use SingleCellExperiment and zellkonverter, so load those as well.
suppressPackageStartupMessages({
library(zellkonverter)
library(SingleCellExperiment) # load early to avoid masking dplyr::count()
library(dplyr)
library(cellxgenedp)
})
cxg()
Provides a ‘shiny’ interfaceThe following sections outline how to use the cellxgenedp package
in an R script; most functionality is also available in the cxg()
shiny application, providing an easy way to identify, download, and
visualize one or several datasets. Start the app
cxg()
choose a project on the first tab, and a dataset for visualization, or one or more datasets for download!
Retrieve metadata about resources available at the cellxgene data
portal using db()
:
db <- db()
Printing the db
object provides a brief overview of the available
data, as well as hints, in the form of functions like collections()
,
for further exploration.
db
## cellxgene_db
## number of collections(): 117
## number of datasets(): 702
## number of files(): 2797
The portal organizes data hierarchically, with ‘collections’ (research studies, approximately), ‘datasets’, and ‘files’. Discover data using the corresponding functions.
collections(db)
## # A tibble: 117 × 15
## collec…¹ acces…² conso…³ conta…⁴ conta…⁵ curat…⁶ data_…⁷ descr…⁸ links name
## <chr> <chr> <lgl> <chr> <chr> <chr> <chr> <chr> <list> <chr>
## 1 43d4bb3… READ NA raymon… Raymon… Batuha… 1.0 Pertur… <list> Tran…
## 2 0434a9d… READ NA avilla… Alexan… Batuha… 1.0 The va… <list> Acut…
## 3 a925421… READ NA st9@sa… Sarah … Batuha… 1.0 During… <list> Sing…
## 4 3472f32… READ NA wongcb… Raymon… Batuha… 1.0 The re… <list> A si…
## 5 03cdc7f… READ NA tien.p… Tien P… Jason … 1.0 scRNA-… <lgl> Emph…
## 6 2902f08… READ NA lopes@… S. M. … Wei Kh… 1.0 The ov… <list> Sing…
## 7 83ed3be… READ NA tom.ta… Tom Ta… Jennif… 1.0 During… <list> Inte…
## 8 2b02dff… READ NA miriam… Miriam… Batuha… 1.0 Clinic… <list> Sing…
## 9 eb735cc… READ NA rv4@sa… Roser … Batuha… 1.0 Human … <list> Samp…
## 10 44531dd… READ NA tallul… Tallul… Jennif… 1.0 The cr… <list> Sing…
## # … with 107 more rows, 5 more variables: publisher_metadata <list>,
## # visibility <chr>, created_at <date>, published_at <date>,
## # updated_at <date>, and abbreviated variable names ¹collection_id,
## # ²access_type, ³consortia, ⁴contact_email, ⁵contact_name, ⁶curator_name,
## # ⁷data_submission_policy_version, ⁸description
datasets(db)
## # A tibble: 702 × 26
## dataset_id colle…¹ donor…² assay cell_…³ cell_…⁴ datas…⁵ devel…⁶ disease
## <chr> <chr> <list> <list> <int> <list> <chr> <list> <list>
## 1 d62144d1-6e98… 43d4bb… <list> <list> 68036 <list> https:… <list> <list>
## 2 07f8f239-2136… 0434a9… <list> <list> 59506 <list> https:… <list> <list>
## 3 605b0fa6-e441… a92542… <list> <list> 70325 <list> https:… <list> <list>
## 4 da2cbedd-f892… 3472f3… <list> <list> 19694 <list> https:… <list> <list>
## 5 9aca7e29-d19c… 03cdc7… <list> <list> 35699 <list> https:… <list> <list>
## 6 731aeab6-e1d6… 03cdc7… <list> <list> 3662 <list> https:… <list> <list>
## 7 bfdafcbc-d785… 03cdc7… <list> <list> 18386 <list> https:… <list> <list>
## 8 218af677-b80b… 2902f0… <list> <list> 20676 <list> https:… <list> <list>
## 9 9e11bf54-0ea6… 83ed3b… <list> <list> 71732 <list> https:… <list> <list>
## 10 7a637398-175f… 2b02df… <list> <list> 32458 <list> https:… <list> <list>
## # … with 692 more rows, 17 more variables: is_primary_data <chr>,
## # is_valid <lgl>, mean_genes_per_cell <dbl>, name <chr>, organism <list>,
## # processing_status <list>, published <lgl>, revision <int>,
## # schema_version <chr>, self_reported_ethnicity <list>, sex <list>,
## # suspension_type <list>, tissue <list>, tombstone <lgl>, created_at <date>,
## # published_at <date>, updated_at <date>, and abbreviated variable names
## # ¹collection_id, ²donor_id, ³cell_count, ⁴cell_type, ⁵dataset_deployments, …
files(db)
## # A tibble: 2,797 × 8
## file_id datas…¹ filen…² filet…³ s3_uri user_…⁴ created_at updated_at
## <chr> <chr> <chr> <chr> <chr> <lgl> <date> <date>
## 1 87ce2086-f2eb-4… d62144… "raw.h… RAW_H5… s3://… TRUE 1970-01-01 1970-01-01
## 2 8eabe485-be71-4… d62144… "local… H5AD s3://… TRUE 1970-01-01 1970-01-01
## 3 33dfa406-132a-4… d62144… "local… RDS s3://… TRUE 1970-01-01 1970-01-01
## 4 454e0b3c-e207-4… d62144… "" CXG s3://… TRUE 1970-01-01 1970-01-01
## 5 62deea2a-d722-4… 07f8f2… "raw.h… RAW_H5… s3://… TRUE 1970-01-01 1970-01-01
## 6 b50c48f1-ef5e-4… 07f8f2… "local… H5AD s3://… TRUE 1970-01-01 1970-01-01
## 7 59309f21-0db8-4… 07f8f2… "local… RDS s3://… TRUE 1970-01-01 1970-01-01
## 8 399c2655-f56a-4… 07f8f2… "" CXG s3://… TRUE 1970-01-01 1970-01-01
## 9 4b70a566-5409-4… 605b0f… "local… RDS s3://… TRUE 1970-01-01 1970-01-01
## 10 0a978dbb-6e45-4… 605b0f… "" CXG s3://… TRUE 1970-01-01 1970-01-01
## # … with 2,787 more rows, and abbreviated variable names ¹dataset_id,
## # ²filename, ³filetype, ⁴user_submitted
Each of these resources has a unique primary identifier (e.g.,
file_id
) as well as an identifier describing the relationship of the
resource to other components of the database (e.g.,
dataset_id
). These identifiers can be used to ‘join’ information
across tables.
facets()
provides information on ‘levels’ present in specific columnsNotice that some columns are ‘lists’ rather than atomic vectors like ‘character’ or ‘integer’.
datasets(db) |>
select(where(is.list))
## # A tibble: 702 × 11
## donor_id assay cell_…¹ devel…² disease organ…³ processing…⁴ self_…⁵ sex
## <list> <list> <list> <list> <list> <list> <list> <list> <list>
## 1 <list [9]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 2 <list> <list> <list> <list> <list> <list> <named list> <list> <list>
## 3 <list> <list> <list> <list> <list> <list> <named list> <list> <list>
## 4 <list [3]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 5 <list [6]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 6 <list [6]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 7 <list [6]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 8 <list [5]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 9 <list [5]> <list> <list> <list> <list> <list> <named list> <list> <list>
## 10 <list> <list> <list> <list> <list> <list> <named list> <list> <list>
## # … with 692 more rows, 2 more variables: suspension_type <list>,
## # tissue <list>, and abbreviated variable names ¹cell_type,
## # ²development_stage, ³organism, ⁴processing_status, ⁵self_reported_ethnicity
This indicates that at least some of the datasets had more than one
type of assay
, cell_type
, etc. The facets()
function provides a
convenient way of discovering possible levels of each column, e.g.,
assay
, organism
, self_reported_ethnicity
, or sex
, and the
number of datasets with each label.
facets(db, "assay")
## # A tibble: 32 × 4
## facet label ontology_term_id n
## <chr> <chr> <chr> <int>
## 1 assay 10x 3' v3 EFO:0009922 345
## 2 assay 10x 3' v2 EFO:0009899 169
## 3 assay Slide-seqV2 EFO:0030062 129
## 4 assay 10x 5' v1 EFO:0011025 46
## 5 assay Smart-seq2 EFO:0008931 37
## 6 assay Visium Spatial Gene Expression EFO:0010961 35
## 7 assay 10x multiome EFO:0030059 24
## 8 assay Drop-seq EFO:0008722 12
## 9 assay 10x 3' transcription profiling EFO:0030003 9
## 10 assay 10x 5' v2 EFO:0009900 9
## # … with 22 more rows
facets(db, "self_reported_ethnicity")
## # A tibble: 18 × 4
## facet label ontol…¹ n
## <chr> <chr> <chr> <int>
## 1 self_reported_ethnicity unknown unknown 354
## 2 self_reported_ethnicity European HANCES… 200
## 3 self_reported_ethnicity na na 184
## 4 self_reported_ethnicity Asian HANCES… 81
## 5 self_reported_ethnicity African American HANCES… 39
## 6 self_reported_ethnicity multiethnic multie… 25
## 7 self_reported_ethnicity Greater Middle Eastern (Middle Easter… HANCES… 21
## 8 self_reported_ethnicity Hispanic or Latin American HANCES… 16
## 9 self_reported_ethnicity African American or Afro-Caribbean HANCES… 10
## 10 self_reported_ethnicity East Asian HANCES… 4
## 11 self_reported_ethnicity African HANCES… 3
## 12 self_reported_ethnicity South Asian HANCES… 2
## 13 self_reported_ethnicity Chinese HANCES… 1
## 14 self_reported_ethnicity Eskimo HANCES… 1
## 15 self_reported_ethnicity Finnish HANCES… 1
## 16 self_reported_ethnicity Han Chinese HANCES… 1
## 17 self_reported_ethnicity Oceanian HANCES… 1
## 18 self_reported_ethnicity Pacific Islander HANCES… 1
## # … with abbreviated variable name ¹ontology_term_id
facets(db, "sex")
## # A tibble: 3 × 4
## facet label ontology_term_id n
## <chr> <chr> <chr> <int>
## 1 sex male PATO:0000384 610
## 2 sex female PATO:0000383 369
## 3 sex unknown unknown 54
Suppose we were interested in finding datasets from the 10x 3’ v3
assay (ontology_term_id
of EFO:0009922
) containing individuals of
African American ethnicity, and female sex. Use the facets_filter()
utility function to filter data sets as needed
african_american_female <-
datasets(db) |>
filter(
facets_filter(assay, "ontology_term_id", "EFO:0009922"),
facets_filter(self_reported_ethnicity, "label", "African American"),
facets_filter(sex, "label", "female")
)
Use nrow(african_american_female)
to find the number of datasets
satisfying our criteria. It looks like there are up to
african_american_female |>
summarise(total_cell_count = sum(cell_count))
## # A tibble: 1 × 1
## total_cell_count
## <int>
## 1 2657433
cells sequenced (each dataset may contain cells from several
ethnicities, as well as males or individuals of unknown gender, so we
do not know the actual number of cells available without downloading
files). Use left_join
to identify the corresponding collections:
## collections
left_join(
african_american_female |> select(collection_id) |> distinct(),
collections(db),
by = "collection_id"
)
## # A tibble: 8 × 15
## collect…¹ acces…² conso…³ conta…⁴ conta…⁵ curat…⁶ data_…⁷ descr…⁸ links name
## <chr> <chr> <lgl> <chr> <chr> <chr> <chr> <chr> <list> <chr>
## 1 c9706a92… READ NA hnaksh… Harikr… Jennif… 1.0 "Singl… <list> A si…
## 2 2f75d249… READ NA rsatij… Rahul … Jennif… 1.0 "This … <list> Azim…
## 3 b9fc3d70… READ NA bruce.… Bruce … Jennif… 1.0 "Numer… <list> A We…
## 4 62e8f058… READ NA chanj3… Joseph… Jennif… 1.0 "155,0… <list> HTAN…
## 5 a98b828a… READ NA markus… Markus… Jennif… 1.0 "Tumor… <lgl> HCA …
## 6 625f6bf4… READ NA a5wang… Allen … Jennif… 1.0 "Large… <list> Lung…
## 7 b953c942… READ NA icobos… Inma C… Jennif… 1.0 "Tau a… <list> Sing…
## 8 bcb61471… READ NA info@k… KPMP Jennif… 1.0 "Under… <list> An a…
## # … with 5 more variables: publisher_metadata <list>, visibility <chr>,
## # created_at <date>, published_at <date>, updated_at <date>, and abbreviated
## # variable names ¹collection_id, ²access_type, ³consortia, ⁴contact_email,
## # ⁵contact_name, ⁶curator_name, ⁷data_submission_policy_version, ⁸description
cellxgene
Discover files associated with our first selected dataset
selected_files <-
left_join(
african_american_female |> select(dataset_id),
files(db),
by = "dataset_id"
)
## Warning in left_join(select(african_american_female, dataset_id), files(db), : Each row in `x` is expected to match at most 1 row in `y`.
## ℹ Row 1 of `x` matches multiple rows.
## ℹ If multiple matches are expected, set `multiple = "all"` to silence this
## warning.
selected_files
## # A tibble: 88 × 8
## dataset_id file_id filen…¹ filet…² s3_uri user_…³ created_at updated_at
## <chr> <chr> <chr> <chr> <chr> <lgl> <date> <date>
## 1 24205601-0780-4… e1c842… "raw.h… RAW_H5… s3://… TRUE 1970-01-01 1970-01-01
## 2 24205601-0780-4… 15e9d9… "local… H5AD s3://… TRUE 1970-01-01 1970-01-01
## 3 24205601-0780-4… 0d3974… "" CXG s3://… TRUE 1970-01-01 1970-01-01
## 4 24205601-0780-4… e254f9… "local… RDS s3://… TRUE 1970-01-01 1970-01-01
## 5 d084f1e1-42f1-4… 59bd46… "local… RDS s3://… TRUE 1970-01-01 1970-01-01
## 6 d084f1e1-42f1-4… 3a2467… "" CXG s3://… TRUE 1970-01-01 1970-01-01
## 7 d084f1e1-42f1-4… d9f9d0… "local… H5AD s3://… TRUE 1970-01-01 1970-01-01
## 8 d084f1e1-42f1-4… bee803… "raw.h… RAW_H5… s3://… TRUE 1970-01-01 1970-01-01
## 9 0f6d4993-c092-4… c76643… "raw.h… RAW_H5… s3://… TRUE 1970-01-01 1970-01-01
## 10 0f6d4993-c092-4… f6d9f2… "local… H5AD s3://… TRUE 1970-01-01 1970-01-01
## # … with 78 more rows, and abbreviated variable names ¹filename, ²filetype,
## # ³user_submitted
The filetype
column lists the type of each file. The cellxgene service
can be used to visualize datasets that have CXG
files.
selected_files |>
filter(filetype == "CXG") |>
slice(1) |> # visualize a single dataset
datasets_visualize()
Visualization is an interactive process, so datasets_visualize()
will only open up to 5 browser tabs per call.
Datasets usually contain CXG
(cellxgene visualization), H5AD
(files produced by the python AnnData module), and Rds
(serialized
files produced by the R Seurat package). There are no public parsers
for CXG
, and the Rds
files may be unreadable if the version of
Seurat used to create the file is different from the version used to
read the file. We therefore focus on the H5AD
files. For
illustration, we download one of our selected files.
local_file <-
selected_files |>
filter(
dataset_id == "24205601-0780-4bf2-b1d9-0e3cacbc2cd6",
filetype == "H5AD"
) |>
files_download(dry.run = FALSE)
basename(local_file)
## [1] "15e9d9af-60dc-4897-88e4-89842655d6b8.H5AD"
These are downloaded to a local cache (use the internal function
cellxgenedp:::.cellxgenedb_cache_path()
for the location of the
cache), so the process is only time-consuming the first time.
H5AD
files can be converted to R / Bioconductor objects using
the zellkonverter package.
h5ad <- readH5AD(local_file, use_hdf5 = TRUE)
## Warning: 'X' matrix does not support transposition and has been skipped
h5ad
## class: SingleCellExperiment
## dim: 33234 31696
## metadata(3): default_embedding schema_version title
## assays(1): X
## rownames(33234): ENSG00000243485 ENSG00000237613 ... ENSG00000277475
## ENSG00000268674
## rowData names(4): feature_is_filtered feature_name feature_reference
## feature_biotype
## colnames(31696): CMGpool_AAACCCAAGGACAACC CMGpool_AAACCCACAATCTCTT ...
## K109064_TTTGTTGGTTGCATCA K109064_TTTGTTGGTTGGACCC
## colData names(34): donor_id self_reported_ethnicity_ontology_term_id
## ... self_reported_ethnicity development_stage
## reducedDimNames(3): X_pca X_tsne X_umap
## mainExpName: NULL
## altExpNames(0):
The SingleCellExperiment
object is a matrix-like object with rows
corresponding to genes and columns to cells. Thus we can easily
explore the cells present in the data.
h5ad |>
colData(h5ad) |>
as_tibble() |>
count(sex, donor_id)
## # A tibble: 7 × 3
## sex donor_id n
## <fct> <fct> <int>
## 1 female D1 2303
## 2 female D2 864
## 3 female D3 2517
## 4 female D4 1771
## 5 female D5 2244
## 6 female D11 7454
## 7 female pooled [D9,D7,D8,D10,D6] 14543
The Orchestrating Single-Cell Analysis with Bioconductor
online resource provides an excellent introduction to analysis and
visualization of single-cell data in R / Bioconductor. Extensive
opportunities for working with AnnData objects in R but using the
native python interface are briefly described in, e.g., ?AnnData2SCE
help page of zellkonverter.
The hca package provides programmatic access to the Human Cell Atlas data portal, allowing retrieval of primary as well as derived single-cell data files.
## R version 4.2.2 (2022-10-31)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.5 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.16-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.16-bioc/R/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB 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] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] cellxgenedp_1.2.2 dplyr_1.1.0
## [3] SingleCellExperiment_1.20.0 SummarizedExperiment_1.28.0
## [5] Biobase_2.58.0 GenomicRanges_1.50.2
## [7] GenomeInfoDb_1.34.8 IRanges_2.32.0
## [9] S4Vectors_0.36.1 BiocGenerics_0.44.0
## [11] MatrixGenerics_1.10.0 matrixStats_0.63.0
## [13] zellkonverter_1.8.0 BiocStyle_2.26.0
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.10 here_1.0.1 dir.expiry_1.6.0
## [4] lattice_0.20-45 png_0.1-8 rprojroot_2.0.3
## [7] digest_0.6.31 utf8_1.2.3 mime_0.12
## [10] R6_2.5.1 evaluate_0.20 httr_1.4.4
## [13] pillar_1.8.1 basilisk_1.10.2 zlibbioc_1.44.0
## [16] rlang_1.0.6 curl_5.0.0 jquerylib_0.1.4
## [19] Matrix_1.5-3 DT_0.27 reticulate_1.28
## [22] rmarkdown_2.20 htmlwidgets_1.6.1 RCurl_1.98-1.10
## [25] shiny_1.7.4 DelayedArray_0.24.0 compiler_4.2.2
## [28] httpuv_1.6.8 xfun_0.37 pkgconfig_2.0.3
## [31] htmltools_0.5.4 tidyselect_1.2.0 tibble_3.1.8
## [34] GenomeInfoDbData_1.2.9 bookdown_0.32 fansi_1.0.4
## [37] withr_2.5.0 later_1.3.0 bitops_1.0-7
## [40] rjsoncons_1.0.0 basilisk.utils_1.10.0 grid_4.2.2
## [43] xtable_1.8-4 jsonlite_1.8.4 lifecycle_1.0.3
## [46] magrittr_2.0.3 cli_3.6.0 cachem_1.0.6
## [49] XVector_0.38.0 promises_1.2.0.1 bslib_0.4.2
## [52] ellipsis_0.3.2 filelock_1.0.2 generics_0.1.3
## [55] vctrs_0.5.2 tools_4.2.2 glue_1.6.2
## [58] parallel_4.2.2 fastmap_1.1.0 yaml_2.3.7
## [61] BiocManager_1.30.19 knitr_1.42 sass_0.4.5