---
title: "Summarizing cell and gene metadata"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Summarizing cell and gene metadata}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
options(width = 88, max.print = 256)
```

<!--`
THIS VIGNETTE IS BASED ON:
https://github.com/chanzuckerberg/cellxgene-census/blob/main/api/python/notebooks/analysis_demo/comp_bio_summarize_axis_query.ipynb
-->

This vignette provides examples for basic axis metadata handling. The CZ CELLxGENE Census stores `obs` (cell) and `var` (gene) metadata in a SOMA DataFrame objects via the [TileDB-SOMA](https://github.com/single-cell-data/TileDB-SOMA) API ([documentation](https://tiledbsoma.readthedocs.io/en/latest/)), which can be queried and read as an R data frame. The Census also has a convenience package which simplifies opening the census.

R data frames are in-memory objects. Take care that queries are small enough for results to fit in memory.

**Contents**

1. Opening the Census
2. Summarizing cell metadata
    - Example: Summarize all cell types
    - Example: Summarize a subset of cell types, selected with a `value_filter`
3. Full Census metadata stats

## Opening the Census

The `cellxgene.census` R package contains a convenient API to open any version of the Census (by default, the newest stable version).

```{r,message=FALSE,warning=FALSE}
library("cellxgene.census")
census <- open_soma()
```

If you open the Census, you should close it with `census$close()`. This can be automated using `on.exit(census$close(), add = TRUE)` immediately after `census <- open_soma()`.

You can learn more about the cellxgene.census methods by accessing their corresponding documentation. For example `?cellxgene.census::open_soma`.

## Summarizing cell metadata

Once the Census is open you can use its TileDB-SOMA methods as it is itself a `SOMACollection`. You can thus access the metadata `SOMADataFram`e objects encoding cell and gene metadata.

Tips:

* You can read an *entire* `SOMADataFrame` into R using `as.data.frame(soma_df$read()$concat())`.
* Queries will be much faster if you request only the DataFrame columns required for your analysis (e.g. `column_names = c("soma_joinid", "cell_type_ontology_term_id")`).
* You can also further refine query results by using a `value_filter`, which will filter the census for matching records.

### Example: Summarize all cell types

This example reads the cell metadata (obs) into an R data frame to summarize in a variety of ways.

```{r}
human <- census$get("census_data")$get("homo_sapiens")

# Read obs into an R data frame (tibble).
obs_df <- human$obs$read(column_names = c("cell_type"))
obs_df <- as.data.frame(obs_df$concat())

# Find all unique values in the cell_type column.
unique_cell_type <- unique(obs_df$cell_type)

cat(
  "There are",
  length(unique_cell_type),
  "cell types in the Census! The first few are: ",
  paste(head(unique_cell_type), collapse = ", ")
)
```

### Example: Summarize a subset of cell types, selected with a `value_filter`

This example utilizes a SOMA "value filter" to read the subset of cells with `tissue_ontology_term_id` equal to `UBERON:0002048` (lung tissue), and summarizes the query result.

```{r}
# Read cell_type terms for cells which have a specific tissue term
LUNG_TISSUE <- "UBERON:0002048"

obs_df <- human$obs$read(column_names = c("cell_type"), value_filter = paste0("tissue_ontology_term_id == '", LUNG_TISSUE, "'"))
obs_df <- as.data.frame(obs_df$concat())

# Find all unique values in the cell_type column as an R data frame.
unique_cell_type <- unique(obs_df$cell_type)
cat(
  "There are ",
  length(unique_cell_type),
  " cell types in the Census where tissue_ontology_term_id == ",
  LUNG_TISSUE,
  "!\nThe first few are:",
  paste(head(unique_cell_type), collapse = ", "),
  "\n"
)

# Report the 10 most common
top_10 <- sort(table(obs_df$cell_type), decreasing = TRUE)[1:10]
cat(
  "The top 10 cell types where tissue_ontology_term_id ==",
  LUNG_TISSUE,
  "are: ",
  paste(names(top_10), collapse = ", ")
)
```

You can also define much more complex value filters. For example:

* combine terms with `&` and `|`
* use the `%in%` operator to query on multiple values

```{r}
# You can also do more complex queries, such as testing for inclusion in a list of values
obs_df <- human$obs$read(
  column_names = c("cell_type_ontology_term_id"),
  value_filter = "tissue_ontology_term_id %in% c('UBERON:0002082', 'UBERON:OOO2084', 'UBERON:0002080')"
)

obs_df <- as.data.frame(obs_df$concat())

# Summarize
top_10 <- sort(table(obs_df$cell_type_ontology_term_id), decreasing = TRUE)[1:10]
print(top_10)
```

## Full Census metadata stats

This example queries all organisms in the Census, and summarizes the diversity of various metadata labels.

```{r}
cols_to_query <- c(
  "cell_type_ontology_term_id",
  "assay_ontology_term_id",
  "tissue_ontology_term_id"
)

total_cells <- 0
for (organism in census$get("census_data")$names()) {
  print(organism)

  obs_df <- census$get("census_data")$get(organism)$obs$read(column_names = cols_to_query)
  obs_df <- as.data.frame(obs_df$concat())

  total_cells <- total_cells + nrow(obs_df)
  for (col in cols_to_query) {
    cat("  Unique ", col, " values: ", length(unique(obs_df[[col]])), "\n")
  }
}
cat("Complete Census contains ", total_cells, " cells.")
```

### Close the census

After use, the census object should be closed to release memory and other resources.

```{r}
census$close()
```

This also closes all SOMA objects accessed via the top-level `census`. Closing can be automated using `on.exit(census$close(), add = TRUE)` immediately after `census <- open_soma()`.
