---
title: "General Helpers & Utilities"
date: 'Compiled: `r format(Sys.Date(), "%B %d, %Y")`'
output: rmarkdown::html_vignette
theme: united
df_print: kable
vignette: >
  %\VignetteIndexEntry{General Helpers & Utilities}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---
***

<style>
p.caption {
  font-size: 0.9em;
}
</style>
 
```{r setup, include=FALSE}
all_times <- list()  # store the time for each chunk
knitr::knit_hooks$set(time_it = local({
  now <- NULL
  function(before, options) {
    if (before) {
      now <<- Sys.time()
    } else {
      res <- difftime(Sys.time(), now, units = "secs")
      all_times[[options$label]] <<- res
    }
  }
}))
knitr::opts_chunk$set(
  tidy = TRUE,
  tidy.opts = list(width.cutoff = 95),
  message = FALSE,
  warning = FALSE,
  time_it = TRUE
)
```

# General Helper/Utility Functions  
scCustomize has several helper functions to simplify/streamline common tasks in scRNA-seq analysis.
Let's load packages and raw data object for this tutorial.
```{r init}
# Load Packages
library(ggplot2)
library(dplyr)
library(magrittr)
library(patchwork)
library(Seurat)
library(scCustomize)

pbmc <- pbmc3k.SeuratData::pbmc3k
```


```{r include=FALSE}
pbmc <- UpdateSeuratObject(pbmc)
```

We'll add some random meta data variables to pbmc data form use in this vignette
```{r}
pbmc$sample_id <- sample(c("sample1", "sample2", "sample3", "sample4", "sample5", "sample6"), size = ncol(pbmc), replace = TRUE)
pbmc$batch <- sample(c("Batch1", "Batch2"), size = ncol(pbmc), replace = TRUE)
```

## Add Mitochondrial & Ribosomal Percentages
As discussed in [QC Plotting vignette](https://samuel-marsh.github.io/scCustomize/articles/QC_Plots.html) one the first steps after creating object if often to calculate and add mitochondrial and ribosomal count percentages per cell/nucleus.

### `Add_Mito_Ribo()`
scCustomize contains easy wrapper function to automatically add both Mitochondrial and Ribosomal percentages to meta.data slot.  If you are using mouse, human, marmoset, zebrafish, rat, drosophila, or rhesus macaque data all you need to do is specify the `species` parameter.
```{r eval=FALSE}
# These defaults can be run just by providing accepted species name
pbmc <- Add_Mito_Ribo(object = pbmc, species = "human")
```

To view list of accepted values for default species names simply set `list_species_names = TRUE`.

```{r eval=FALSE}
Add_Mito_Ribo(list_species_names = TRUE)
```


```{r echo=FALSE}
accepted_names <- Add_Mito_Ribo(object = pbmc, list_species_names = TRUE)

accepted_names %>%
  kableExtra::kbl(row.names = TRUE) %>%
  kableExtra::kable_styling(bootstrap_options = c("bordered", "condensed", "responsive", "striped")) 
```

#### Non-default species
However custom prefixes can be used for non-human/mouse/marmoset/rat/zebrafish/drosophila/macaque species with different annotations. Simply specify `species = other` and supply feature lists or regex patterns for your species of interest. 
*NOTE: If desired please submit issue on GitHub for additional default species.  Please include regex pattern or list of genes for both mitochondrial and ribosomal genes and I will add additional built-in defaults to the function.*
```{r eval=FALSE}
# Using gene name patterns
pbmc <- Add_Mito_Ribo(object = pbmc, species = "other", mito_pattern = "regexp_pattern", ribo_pattern = "regexp_pattern")

# Using feature name lists
mito_gene_list <- c("gene1", "gene2", "etc")
ribo_gene_list <- c("gene1", "gene2", "etc")
pbmc <- Add_Mito_Ribo(object = pbmc, species = "other", mito_features = mito_gene_list, ribo_features = ribo_gene_list)

# Using combination of gene lists and gene name patterns
pbmc <- Add_Mito_Ribo(object = pbmc, species = "Human", mito_features = mito_gene_list, ribo_pattern = "regexp_pattern")
```

### Warning Messages
The added benefit of `Add_Mito_Ribo` is that it will return informative warnings if no Mitochondrial or Ribosomal features are found using the current species, features, or pattern specification.
```{r message=TRUE, warning=TRUE, error=TRUE}
# For demonstration purposes we can set `species = mouse` for this object of human cells
pbmc <- Add_Mito_Ribo(object = pbmc, species = "mouse")
```

```{r include=FALSE}
pbmc <- pbmc3k.SeuratData::pbmc3k
pbmc <- UpdateSeuratObject(object = pbmc)
```


```{r message=TRUE, warning=TRUE, error=TRUE}
# Or if providing custom patterns/lists and features not found
pbmc <- Add_Mito_Ribo(object = pbmc, species = "other", mito_pattern = "^MT-", ribo_pattern = "BAD_PATTERN")
```

`Add_Mito_Ribo` will also return warnings if columns are already present in `@meta.data` slot and prompt you to provide override if you want to run the function.
```{r include=FALSE}
pbmc <- pbmc3k.SeuratData::pbmc3k
pbmc <- UpdateSeuratObject(object = pbmc)
pbmc <- Add_Mito_Ribo(object = pbmc, species = "human")
```

```{r message=TRUE, warning=TRUE, error=TRUE}
pbmc <- Add_Mito_Ribo(object = pbmc, species = "human")
```



### LIGER Objects: `Add_Mito_Ribo()`
scCustomize `Add_Mito_Ribo` also works seemlessly with LIGER objects.  
```{r eval=FALSE}
liger_obj <- Add_Mito_Ribo(object = liger_obj, species = "human")
```


## Add Cell Complexity/Novelty QC Metrics  
In addition to metrics like number of features and UMIs it can often be helpful to analyze the complexity of expression within a single cell.  scCustomize provides functions to add two of these metrics to meta data.  

### Cell Complexity (log10(nFeature) / log10(nCount))
scCustomize contains easy shortcut function to add a measure of cell complexity/novelty that can sometimes be useful to filter low quality cells.  The metric is calculated by calculating the result of log10(nFeature) / log10(nCount).  
```{r eval = FALSE}
# These defaults can be run just by providing accepted species name
pbmc <- Add_Cell_Complexity(object = pbmc)
```
*NOTE: The function also works seemlessly with LIGER objects.*


### Add Top Percent Expression QC Metric  
Additionally, (or alternatively), scCustomize contains another metric of complexity which is the top percent expression.  The user supplies an integer value for `num_top_genes` (default is 50) which species the number of genes and the function returns percentage of counts occupied by top XX genes in each cell.
```{r eval = FALSE}
# These defaults can be run just by providing accepted species name
pbmc <- Add_Top_Gene_Pct_Seurat(seurat_object = pbmc, num_top_genes = 50)
```


## Add All Cell QC Metrics with Single Function
To simplify the process of adding cell QC metrics scCustomize contains a wrapper function which can be customized to add all or some of the available QC metrics.  The default parameters of the function `Add_Cell_QC_Metrics` will add:  

* Mitochondrial and Ribosomal Percentages (default and custom species).
* Cell Complexity (log10(nFeature) / log10(nCount).
* Top XX Gene Percentage.
* Percentage of counts for IEG (human and mouse only).
* OXPHOS, APOP, and DNA Repair pathways (supported species only).
* Cell Cycle Scoring (Human only).

```{r eval=FALSE}
pbmc <- Add_Cell_QC_Metrics(seurat_object = pbmc, species = "human")
```



## Meta Data  
scCustomize contains a set of functions to aid in use of meta data both within and outside of objects.  

```{r, include = FALSE}
pbmc <- pbmc3k.SeuratData::pbmc3k
pbmc <- UpdateSeuratObject(object = pbmc)

pbmc$sample_id <- sample(c("sample1", "sample2", "sample3", "sample4", "sample5", "sample6"), size = ncol(pbmc), replace = TRUE)
pbmc$batch <- sample(c("Batch1", "Batch2"), size = ncol(pbmc), replace = TRUE)
```



### Extract all meta data (cell-level) from object  
`Fetch_Meta()` functions as simple getter function to obtain meta data from object and return data.frame.  
```{r eval = FALSE}
meta_data <- Fetch_Meta(object = pbmc)

head(meta_data, 10)
```

```{r echo=FALSE}
meta_data <- Fetch_Meta(object = pbmc)

head(meta_data, 10) %>%
  kableExtra::kbl(row.names = TRUE) %>%
  kableExtra::kable_styling(bootstrap_options = c("bordered", "condensed", "responsive", "striped")) 
```

### Extract sample-level meta data from object  
While cell-level meta data is helpful in some situations often all that is required is sample-level meta data.  This can easily be extracted and filtered using `Extract_Sample_Meta()`.  

```{r}
sample_meta <- Extract_Sample_Meta(object = pbmc, sample_name = "sample_id")
```

```{r echo=FALSE}
sample_meta %>%
  kableExtra::kbl(row.names = TRUE) %>%
  kableExtra::kable_styling(bootstrap_options = c("bordered", "condensed", "responsive", "striped")) 
```

#### Remove columns that are not relevant on sample level  
As you can see by default `Extract_Sample_Meta` removes a default set of columns (see documentation) which do not provide meaningful sample-level information (e.g., nFeature_RNA). However, you may want to remove other columns too.  This can be achieved using either positive or negative selection using `variables_include` or `variables_exclude` parameters.

```{r}
sample_meta <- Extract_Sample_Meta(object = pbmc, sample_name = "sample_id", variables_exclude = c("nFeature_RNA", "nCount_RNA", "seurat_annotations", "orig.ident"))
```


```{r echo=FALSE}
sample_meta %>%
  kableExtra::kbl(row.names = TRUE) %>%
  kableExtra::kable_styling(bootstrap_options = c("bordered", "condensed", "responsive", "striped")) 
```


#### Merge with sample-level summary statistic data  
While some original number columns are not valid at sample-level it can be valuable to get summary information for those variables.  This can be achieved by merging outputs with `Median_Stats` function.

```{r}
sample_meta <- Extract_Sample_Meta(object = pbmc, sample_name = "sample_id", variables_exclude = c("nFeature_RNA", "nCount_RNA", "seurat_annotations", "orig.ident"))

sample_median <- Median_Stats(seurat_object = pbmc, group_by_var = "sample_id")

sample_merged <- right_join(x = sample_meta, y = sample_median)
```

```{r echo=FALSE}
sample_merged %>%
  kableExtra::kbl(row.names = TRUE) %>%
  kableExtra::kable_styling(bootstrap_options = c("bordered", "condensed", "responsive", "striped")) 
```


### Add sample-level meta data to object  
scCustomize provides easy function to add sample-level meta data to object without the need to first convert it to cell-level meta data.  This makes adding meta data from summary or supplemental tables to cell level object data very easy.

In order to add meta data you will need to specify:  

* Seurat Object with `@meta.data` column that matches sample-level meta data (often "orig.ident").  
* Sample-level meta data/file.  Can be either a data.frame or tibble in environment or a path to file to read.  
* Joining parameters.  Define the column names the specify matching information in Seurat object and sample-level meta data.  

This is example command:
```{r, eval = FALSE}
obj <- Add_Sample_Meta(seurat_object = obj, meta_data = sample_meta, join_by_seurat = "orig.ident", join_by_meta = "sample_id")
```


### Add feature meta data  
Starting in Seurat V5 each assay now possess it's own meta.data slot which is feature-level meta data.  During course of normal analysis this is where information on variable features is stored.  However, we can also use it to store alternate feature names, in most cases this is Ensembl IDs matching the symbols used in object creation/analysis.

scCustomize provides the function `Add_Alt_Feature_ID()` to automatically match and add these features using the same files used in object creation.  Users only need to supply either path to the features.tsv.gz file or the hdf5 file produced from Cell Ranger output.

```{r eval=FALSE}
# Using features.tsv.gz file
obj <- Add_Alt_Feature_ID(seurat_object = obj,
features_tsv = "sample01/outs/filtered_feature_bc_matrix/features.tsv.gz", assay = "RNA")

# Using hdf5 file
obj <- Add_Alt_Feature_ID(seurat_object = obj,
hdf5_file = "sample01/outs/outs/filtered_feature_bc_matrix.h5"", assay = "RNA")
```

*NOTE:* If using features.tsv.gz file the file from either filtered or raw outputs can be used as they are identical.  

*NOTE:* If using hdf5 file the file from either filtered_feature_bc or raw_feature_bc can be used as the features slot is identical. Though it is faster to load filtered_feature_bc file due to droplet filtering.  


## Check for Features/Genes
scCustomize also makes forward-facing a number of utilities that are used internally in functions but may also have utility on their own.

### `Gene_Present()` to check for features.
`Gene_Present` is fairly basic function to check if feature exists in data.  It can be used with Seurat or LIGER objects as well as generic data formats (Matrix, data.frame, tibble).

In addition to some warning messages `Gene_Present` returns a list with 3 entries when run:  

* found_features: features found in the data.  
* bad_features: features not found in the data.  
* wrong_case_found_features: features found but in different case than present in input gene list.  
  - *If `bad_features` > 0 then `Gene_Present` will convert `the gene list `bad_features` to all upper case and to sentence case and check against all possible features to see if wrong case was provided.*  

```{r message=TRUE, warning=TRUE}
# Example gene list with all examples (found genes, wrong case (lower) and misspelled (CD8A forgetting to un-shift when typing 8))
gene_input_list <- c("CD14", "CD3E", "Cd4", "CD*A")

genes_present <- Gene_Present(data = pbmc, gene_list = gene_input_list)
```

Now let's look at the output:
```{r message=TRUE, warning=TRUE}
genes_present
```

### Turn warnings/messages off.
By default `Gene_Present` has 3 sets of warnings/messages it prints to console when it finds issues.  If using the function yourself on its own or wrapped inside your own function and you prefer no messages each of these can be toggled using optional parameters.  

* `case_check_msg` prints and list of features if alternate case features are found in data.  
* `omit_warn` prints warning and list of all features not found in data.  
* `print_msg` prints message if all features in `gene_list` are found in data.  

### Check for updated gene symbols
In order to keep run times down and support offer greater support for offline use `Gene_Present` does not include a check for updated gene symbols.  If you're dataset is from human cells/donors you can simply supply the not found features from `Gene_Present` to Seurat's `UpdateSymbolList` function.

```{r message=TRUE, warning=TRUE}
gene_input_list <- c("CD14", "CD3E", "Cd4", "CD*A", "SEPT1")

genes_present <- Gene_Present(data = pbmc, gene_list = gene_input_list)

check_symbols <- UpdateSymbolList(symbols = genes_present[[2]], verbose = TRUE)
```


## Merging & Splitting Data/Objects

### Merging raw data
It can often be advantageous to merge raw data before creating analysis objects vs creating lots of objects and merging them all later.  scCustomize features a modified version of the internal LIGER function `MergeSparseDataAll()`.

`Merge_Sparse_Data_All()` will combine a list of sparse matrices and return single sparse matrix.  Additionally, by specifying the `add_cell_ids` parameter you can specify a prefix to be added to the barcodes from each entry in the list (using "_" as delimiter).

This function can be especially useful when combined with any of the scCustomize's `Read_` data functions which automatically return named lists of matrices and the ability to specify sample `orig.ident` when creating Seurat objects.  See [Read & Write Vignette](https://samuel-marsh.github.io/scCustomize/articles/Read_and_Write_Functions.html) for more info on the data import functions.
```{r eval=FALSE}
# Read in data
GEO_10X <- Read10X_GEO(data_dir = "assets/GSE152183_RAW_Marsh/")

# Merge data and add sample prefix
GEO_10X_merged <- Merge_Sparse_Data_All(matrix_list = GEO_10X, add_cell_ids = names(GEO_10X))

# Create Seurat Object and specify orig.ident location
GEO_10X_Seurat <- Seurat::CreateSeuratObject(counts = GEO_10X_merged, names.field = 1, names.delim = "_", min.features = 200, min.cells = 5)
```

### Merging List of Seurat Objects
Sometimes it can be advantageous to create a list of multiple Seurat Objects in order to run similar pipeline on all objects in loop.  

To facilitate ease in merging such lists into single object scCustomize contains simple wrapper `Merge_Seurat_List` that uses `purrr::reduce()` to merge all objects in list into single combined object
```{r eval=FALSE}
list_of_objects <- list(obj1, obj2, obj2, ..., obj10)

merged_seurat <- Merge_Seurat_List(list_seurat = list_of_objects)

# Can also add sample specific ids to each object during the merge
cell_ids <- c("sample1", "sample2", "sample3", ..., "sample10")

merged_seurat <- Merge_Seurat_List(list_seurat = list_of_objects, add.cell.ids = cell_ids)
```


### Splitting Seurat V5 Layers  
Seurat V5 objects now have the ability to split within the object into layers.  However, I find that the syntax to do this is not the most intuitive and can be simplified with a new simple wrapper function: `Split_Layers()`

```{r message=TRUE, warning=TRUE}
pbmc <- Split_Layers(seurat_object = pbmc, split.by = "sample_id")
```

`Split_Layers()` defaults to "RNA" assay but can be used for any assay present in object (users should check whether splitting assay other than "RNA" is valid before proceeding).


## Storing Misc Information in Seurat Objects
Seurat objects contain an extra empty slot that can be used to store any extra information desired.  
scCustomize contains two functions `Store_Misc_Info_Seurat` and a wrapper around that function `Store_Palette_Seurat` to make this process easy.

```{r eval=FALSE}
# Data can be vectors or data.frames
misc_info <- "misc_vector_dataframe_list_etc"

# Add data to the @misc slot in Seurat Object
pbmc <- Store_Misc_Info_Seurat(seurat_object = pbmc, data_to_store = misc_info, data_name = "misc_info_name")
```

### Storing Lists
If you are storing a list in the `@misc` slot there is additional parameter that dictates whether to store the information as a list or whether to store each entry in the list separately.
```{r eval=FALSE}
# Create list
misc_info <- list("misc_item1", "misc_item2", etc)

# Store the list directly
pbmc <- Store_Misc_Info_Seurat(seurat_object = pbmc, data_to_store = misc_info, data_name = "misc_info_name", list_as_list = TRUE)

# Store each entry in list as separate entity in `@misc` slot
pbmc <- Store_Misc_Info_Seurat(seurat_object = pbmc, data_to_store = misc_info, data_name = "misc_info_name", list_as_list = FALSE)
```

### Storing Color Palettes
One of the most common times I use this function is to store color palettes associated with clustering or subclustering.  
To make it easier to remember function call in this situation scCustomize contains a wrapper function `Store_Palette_Seurat`.
```{r eval=FALSE}
# Data can be vectors or data.frames
annotated_color_palette <- c("color1", "color2", "color3", "etc")

# Add data to the @misc slot in Seurat Object
pbmc <- Store_Palette_Seurat(seurat_object = pbmc, palette = annotated_color_palette, palette_name = "Round01_Color_Pal")

# Then you can easily call that palette (with tab completion) when plotting without ever needing to reload the palette in current environment
DimPlot(object = pbmc, cols = pbmc@misc$Round01_Color_Pal)
```

## Modify Cell Barcodes in Raw Data
Sometimes, especially with public data, you may want to modify the cell barcode names before creating analysis object.

scCustomize contains a selection of functions to simplify this process:  

### Replace Suffixes  
`Replace_Suffix` can be used on single matrix/data.frame or list of matrices/data.frames to modify to remove suffixes
```{r eval=FALSE}
# For single object
data_mod <- Replace_Suffix(data = raw_data, current_suffix = "-1", new_suffix = "-2")

# For list of objects containing same suffix
raw_data_list <- list(raw_data1, raw_data2, raw_data3, ..., raw_data10)

new_suffixes <- c("-1", "-2", "-3", ..., "-10")

data_mod <- Replace_Suffix(data = raw_data_list, current_suffix = "-1", new_suffix = new_suffixes)


# For list of objects containing different suffixes
raw_data_list <- list(raw_data1, raw_data2, raw_data3, ..., raw_data10)

old_suffixes <- c("-A", "-B", "-C", ..., "-J")
new_suffixes <- c("-1", "-2", "-3", ..., "-10")

data_mod <- Replace_Suffix(data = raw_data_list, current_suffix = old_suffixes, new_suffix = new_suffixes)

```

### Strip Suffixes
`Replace_Suffix` can also be used to strip suffixes from data
```{r eval=FALSE}
# For single object
data_mod <- Replace_Suffix(data = raw_data, current_suffix = "-1", new_suffix = "")
```


### Change Prefix/Suffix Delimiters
scCustomize has 3 functions to facilitate changing the type of delimiters present in cell barcodes.  

* `Change_Delim_Prefix()` Change just the suffix delimiter.
* `Change_Delim_Suffix()` Change just the prefix delimiter.
* `Change_Delim_All()` Change all delimiters.

These functions all take identical inputs and can be applied to either single matrix/data.frames or lists of matrices/data.frames.

```{r eval=FALSE}
data_mod <- Change_Delim_Prefix(data = raw_data, current_delim = ".", new_delim = "_")

data_mod <- Change_Delim_Suffix(data = raw_data, current_delim = ".", new_delim = "_")

data_mod <- Change_Delim_All(data = raw_data, current_delim = ".", new_delim = "_")
```

