"""Tools to analyze single-cell RNA"""
from pipen.utils import mark
from ..core.proc import Proc
from ..core.config import config
class SeuratLoading(Proc):DOCS
"""Seurat - Loading data
Deprecated, should be superseded by SeuratPreparing
Input:
metafile: The metadata of the samples
A tab-delimited file
Two columns are required:
- `Sample` to specify the sample names.
- `RNAData` to assign the path of the data to the samples
The path will be read by `Read10X()` from `Seurat`
Output:
rdsfile: The RDS file with a list of Seurat object
Envs:
qc: The QC filter for each sample.
This will be passed to `subset(obj, subset=<qc>)`.
For example
`nFeature_RNA > 200 & nFeature_RNA < 2500 & percent.mt < 5`
"""
input = "metafile:file"
output = "rdsfile:file:{{in.metafile | stem}}.seurat.RDS"
envs = {"qc": ""}
lang = config.lang.rscript
script = "file://../scripts/scrna/SeuratLoading.R"
class SeuratPreparing(Proc):DOCS
"""Load, prepare and apply QC to data, using `Seurat`
This process will -
- Prepare the seurat object
- Apply QC to the data (gene, cell, and contamination)
- Integrate the data from different samples
See also
- <https://satijalab.org/seurat/articles/pbmc3k_tutorial.html#standard-pre-processing-workflow-1)>
- <https://satijalab.org/seurat/articles/integration_introduction>
This process will read the scRNA-seq data, based on the information provided by
`SampleInfo`, specifically, the paths specified by the `RNAData` column.
The RNAData column should contain the path to the 10X, ParseBio data or HIVE data, either a directory or
a file If the path is a directory, the function will look for barcodes.tsv.gz, features.tsv.gz
and matrix.mtx.gz. The directory should be loaded by Seurat::Read10X, Seurat::ReadParseBio or
the HIVE data. Sometimes, there may be prefix in the file names, e.g. "'prefix'.barcodes.tsv.gz",
which is also supported. If the path is a file ending with ".loom", it will be loaded by
SeuratDisk::Connect() and converted to a Seurat object. Otherwise, if the path is a file,
it should be a h5 file that can be loaded by `Seurat::Read10X_h5()`.
Each sample will be loaded individually and then merged into one `Seurat` object, and then perform QC.
In order to perform QC, some additional columns are added to the meta data of the `Seurat` object. They are:
- `precent.mt`: The percentage of mitochondrial genes.
- `percent.ribo`: The percentage of ribosomal genes.
- `precent.hb`: The percentage of hemoglobin genes.
- `percent.plat`: The percentage of platelet genes.
For integration, two routes are available:
- [Performing integration on datasets normalized with `SCTransform`](https://satijalab.org/seurat/articles/seurat5_integration#perform-streamlined-one-line-integrative-analysis)
- [Using `NormalizeData` and `FindIntegrationAnchors`](https://satijalab.org/seurat/articles/seurat5_integration#layers-in-the-seurat-v5-object)
/// Note
When using `SCTransform`, the default Assay will be set to `SCT` in output, rather than `RNA`.
If you are using `cca` or `rpca` interation, the default assay will be `integrated`.
///
/// Note
From `biopipen` v0.23.0, this requires `Seurat` v5.0.0 or higher.
///
Input:
metafile: The metadata of the samples
A tab-delimited file
Two columns are required:
`Sample` to specify the sample names.
`RNAData` to assign the path of the data to the samples
The path will be read by `Read10X()` from `Seurat`, or the path
to the h5 file that can be read by `Read10X_h5()` from `Seurat`.
It can also be an RDS or qs2 file containing a `Seurat` object.
Note that it must has a column named `Sample` in the meta.data to specify the sample names.
Output:
outfile: The qs2 file with the Seurat object with all samples integrated.
Note that the cell ids are prefixied with sample names.
Envs:
ncores (type=int): Number of cores to use.
Used in `future::plan(strategy = "multicore", workers = <ncores>)`
to parallelize some Seurat procedures.
mutaters (type=json): The mutaters to mutate the metadata to the cells.
These new columns will be added to the metadata of the Seurat object and
will be saved in the output file.
See <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
min_cells (type=int): The minimum number of cells that a gene must be
expressed in to be kept. This is used in `Seurat::CreateSeuratObject()`.
Futher QC (`envs.cell_qc`, `envs.gene_qc`) will be performed after this.
It doesn't work when data is loaded from loom files or RDS/qs2 files.
min_features (type=int): The minimum number of features that a cell must
express to be kept. This is used in `Seurat::CreateSeuratObject()`.
Futher QC (`envs.cell_qc`, `envs.gene_qc`) will be performed after this.
It doesn't work when data is loaded from loom files or RDS/qs2 files.
features (type=auto): Features to rename. It can be a dict with the keys as
the original feature names and the values as the new feature names, or
file path to a two-column tab-delimited file (no header) with the first column as
the original feature names and the second column as the new feature names.
cell_qc: Filter expression to filter cells, using
`tidyrseurat::filter()`.
It can also be a dictionary of expressions, where the names of the list are
sample names.
You can have a default expression in the list with the name "DEFAULT" for
the samples that are not listed.
Available QC keys include `nFeature_RNA`, `nCount_RNA`,
`percent.mt`, `percent.ribo`, `percent.hb`, and `percent.plat`.
/// Tip | Example
Including the columns added above, all available QC keys include
`nFeature_RNA`, `nCount_RNA`, `percent.mt`, `percent.ribo`, `percent.hb`,
and `percent.plat`. For example:
```toml
[SeuratPreparing.envs]
cell_qc = "nFeature_RNA > 200 & percent.mt < 5"
```
will keep cells with more than 200 genes and less than 5%% mitochondrial
genes.
///
gene_qc (ns): Filter genes.
`gene_qc` is applied after `cell_qc`.
- min_cells: The minimum number of cells that a gene must be
expressed in to be kept.
- excludes: The genes to exclude. Multiple genes can be specified by
comma separated values, or as a list.
/// Tip | Example
```toml
[SeuratPreparing.envs]
gene_qc = { min_cells = 3 }
```
will keep genes that are expressed in at least 3 cells.
///
ccs_args (ns): Arguments for `RunSeuratCellCycleScoring()`.
When "S.Score" and/or "G2M.Score" are specified in `envs.SCTransform.vars-to-regress` when `envs.use_sct = TRUE`
or in `envs.ScaleData.vars-to-regress` when `envs.use_sct = FALSE`, this will be automatically enabled
to calculate the cell cycle scores and regress them out in the following transformation step.
Otherwise if this is empty, cell cycle scoring and the normalization before cell cycle scoring will be skipped.
- trans_args (ns): Arguments for `RunSeuratTransformation()` to normalize the data before cell cycle scoring.
`use_sct` will be `FALSE` by default (not using `envs.use_sct`).
See <https://github.com/satijalab/seurat/issues/7694>.
- <more>: See <https://pwwang.github.io/biopipen.utils.R/reference/RunSeuratTransformation.html>.
- <more>: See <https://pwwang.github.io/biopipen.utils.R/reference/RunSeuratCellCycleScoring.html>.
qc_plots (type=json): The plots for QC metrics.
It should be a json (or python dict) with the keys as the names of the plots and
the values also as dicts with the following keys:
* kind: The kind of QC. Either `gene`, `cell` (default), or `contam`/`contamination`.
* devpars: The device parameters for the plot. A dict with `res`, `height`, and `width`.
* more_formats: The formats to save the plots other than `png`.
* save_code: Whether to save the code to reproduce the plot.
* other arguments passed to
[`biopipen.utils::VizSeuratCellQC`](https://pwwang.github.io/biopipen.utils.R/reference/VizSeuratCellQC.html)
when `kind` is `cell`,
[`biopipen.utils::VizSeuratGeneQC`](https://pwwang.github.io/biopipen.utils.R/reference/VizSeuratGeneQC.html)
when `kind` is `gene`, or
[`biopipen.utils::RunSeuratContamination`](https://pwwang.github.io/biopipen.utils.R/reference/RunSeuratContamination.html)
when `kind` is `contam`/`contamination`.
use_sct (flag): Whether use SCTransform routine to integrate samples or not.
Before the following procedures, the `RNA` layer will be split by samples.
If `False`, following procedures will be performed in the order:
* [`NormalizeData`](https://satijalab.org/seurat/reference/normalizedata).
* [`FindVariableFeatures`](https://satijalab.org/seurat/reference/findvariablefeatures).
* [`ScaleData`](https://satijalab.org/seurat/reference/scaledata).
See <https://satijalab.org/seurat/articles/seurat5_integration#layers-in-the-seurat-v5-object>
and <https://satijalab.org/seurat/articles/pbmc3k_tutorial.html>
If `True`, following procedures will be performed in the order:
* [`SCTransform`](https://satijalab.org/seurat/reference/sctransform).
See <https://satijalab.org/seurat/articles/seurat5_integration#perform-streamlined-one-line-integrative-analysis>
no_integration (flag): Whether to skip integration or not.
By default, if there are multiple samples, integration will be performed. If `no_integration` is `True`, the samples will be merged without integration.
If there is only one sample, integration will be skipped regardless of the value of `no_integration`.
NormalizeData (ns): Arguments for [`NormalizeData()`](https://satijalab.org/seurat/reference/normalizedata).
`object` is specified internally, and `-` in the key will be replaced with `.`.
- <more>: See <https://satijalab.org/seurat/reference/normalizedata>
FindVariableFeatures (ns): Arguments for [`FindVariableFeatures()`](https://satijalab.org/seurat/reference/findvariablefeatures).
`object` is specified internally, and `-` in the key will be replaced with `.`.
- <more>: See <https://satijalab.org/seurat/reference/findvariablefeatures>
ScaleData (ns): Arguments for [`ScaleData()`](https://satijalab.org/seurat/reference/scaledata).
`object` and `features` is specified internally, and `-` in the key will be replaced with `.`.
You can specify `features` to scale specific features, or set it to `"__all__"` to scale all features.
- <more>: See <https://satijalab.org/seurat/reference/scaledata>
RunPCA (ns): Arguments for [`RunPCA()`](https://satijalab.org/seurat/reference/runpca).
`object` and `features` is specified internally, and `-` in the key will be replaced with `.`.
- npcs (type=int): The number of PCs to compute.
For each sample, `npcs` will be no larger than the number of columns - 1.
- <more>: See <https://satijalab.org/seurat/reference/runpca>
SCTransform (ns): Arguments for [`SCTransform()`](https://satijalab.org/seurat/reference/sctransform).
`object` is specified internally, and `-` in the key will be replaced with `.`.
- return-only-var-genes: Whether to return only variable genes.
- min_cells: The minimum number of cells that a gene must be expressed in to be kept.
A hidden argument of `SCTransform` to filter genes.
If you try to keep all genes in the `RNA` assay, you can set `min_cells` to `0` and
`return-only-var-genes` to `False`.
See <https://github.com/satijalab/seurat/issues/3598#issuecomment-715505537>
- <more>: See <https://satijalab.org/seurat/reference/sctransform>
IntegrateLayers (ns): Arguments for [`IntegrateLayers()`](https://satijalab.org/seurat/reference/integratelayers).
`object` is specified internally, and `-` in the key will be replaced with `.`.
When `use_sct` is `True`, `normalization-method` defaults to `SCT`.
- method (choice): The method to use for integration.
- CCAIntegration: Use `Seurat::CCAIntegration`.
- CCA: Same as `CCAIntegration`.
- cca: Same as `CCAIntegration`.
- RPCAIntegration: Use `Seurat::RPCAIntegration`.
- RPCA: Same as `RPCAIntegration`.
- rpca: Same as `RPCAIntegration`.
- HarmonyIntegration: Use `Seurat::HarmonyIntegration`.
- Harmony: Same as `HarmonyIntegration`.
- harmony: Same as `HarmonyIntegration`.
- FastMNNIntegration: Use `Seurat::FastMNNIntegration`.
- FastMNN: Same as `FastMNNIntegration`.
- fastmnn: Same as `FastMNNIntegration`.
- scVIIntegration: Use `Seurat::scVIIntegration`.
- scVI: Same as `scVIIntegration`.
- scvi: Same as `scVIIntegration`.
- <more>: See <https://satijalab.org/seurat/reference/integratelayers>
contam_correction (choice): The tool used to perform contamination correction.
If None or not specified, no contamination correction will be performed.
- decontX: Use `decontX` to perform contamination correction.
See: https://www.camplab.net/decontx/
- scCDC: Use `scCDC` to perform contamination correction.
See: https://github.com/ZJU-UoE-CCW-LAB/scCDC
decontX (ns): Arguments for `decontX()`.
- <more>: See <https://rdrr.io/bioc/celda/man/decontX.html>
scCDC (ns): Arguments for `scCDC` functions:
- Detection (ns): arguments for `scCDC::ContaminationDetection()`
- <more>: See https://github.com/ZJU-UoE-CCW-LAB/scCDC/blob/main/R/Contamination_Detection.R#L273
- Correction (ns): arguments for `scCDC::ContaminationCorrection()`
- <more>: See https://github.com/ZJU-UoE-CCW-LAB/scCDC/blob/main/R/Contamination_Correction.R#L147
- Quantification (ns): arguments for `scCDC::ContaminationQuantification()`
- <more>: See https://github.com/ZJU-UoE-CCW-LAB/scCDC/blob/main/R/Contamination_Quantification.R#L60
keep_contam_assay (flag): Whether to keep the "Contaminated" (original) assay after QC is finished.
If kept, we can use it to visualize some marker expressions for comparisons.
Note that when `False` (default), the `Contaminated` assay is dropped per-sample right after
contamination correction (before samples are merged) to reduce the memory usage.
Contamination-expression QC plots (with `metric` of `expr`/`expression`) in `qc_plots`
require this to be `True`.
doublet_detector (choice): The doublet detector to use.
- none: Do not use any doublet detector.
- DoubletFinder: Use `DoubletFinder` to detect doublets.
- doubletfinder: Same as `DoubletFinder`.
- scDblFinder: Use `scDblFinder` to detect doublets.
- scdblfinder: Same as `scDblFinder`.
DoubletFinder (ns): Arguments to run [`DoubletFinder`](https://github.com/chris-mcginnis-ucsf/DoubletFinder).
See also <https://demultiplexing-doublet-detecting-docs.readthedocs.io/en/latest/DoubletFinder.html>.
- PCs (type=int): Number of PCs to use for 'doubletFinder' function.
- doublets (type=float): Number of expected doublets as a proportion of the pool size.
- pN (type=float): Number of doublets to simulate as a proportion of the pool size.
- ncores (type=int): Number of cores to use for `DoubletFinder::paramSweep`.
Set to `None` to use `envs.ncores`.
Since parallelization of the function usually exhausts memory, if big `envs.ncores` does not work
for `DoubletFinder`, set this to a smaller number.
scDblFinder (ns): Arguments to run [`scDblFinder`](https://rdrr.io/bioc/scDblFinder/man/scDblFinder.html).
- dbr (type=float): The expected doublet rate.
- ncores (type=int): Number of cores to use for `scDblFinder`.
Set to `None` to use `envs.ncores`.
- <more>: See <https://rdrr.io/bioc/scDblFinder/man/scDblFinder.html>.
cache (type=auto): Whether to cache the information at different steps.
If `True`, the seurat object will be cached in the job output directory, which will be not cleaned up when job is rerunning.
The cached seurat object will be saved as `<signature>.<kind>.RDS` file, where `<signature>` is the signature determined by
the input and envs of the process.
See <https://github.com/satijalab/seurat/issues/7849>, <https://github.com/satijalab/seurat/issues/5358> and
<https://github.com/satijalab/seurat/issues/6748> for more details also about reproducibility issues.
To not use the cached seurat object, you can either set `cache` to `False` or delete the cached file at
`<signature>.RDS` in the cache directory.
Note that caching saves full snapshots of the seurat object at step boundaries, which transiently
increases the peak memory (on both save and load). Set `cache` to `False` on memory-constrained
runs with big data to avoid the extra peaks.
Requires:
r-seurat:
- check: {{proc.lang}} <(echo "library(Seurat)")
r-future:
- check: {{proc.lang}} <(echo "library(future)")
r-bracer:
- check: {{proc.lang}} <(echo "library(bracer)")
""" # noqa: E501
input = "metafile:file"
output = "outfile:file:{{in.metafile | stem}}.seurat.qs"
lang = config.lang.rscript
envs_depth = 4
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"min_cells": 0,
"min_features": 0,
"features": None,
"cell_qc": None, # "nFeature_RNA > 200 & percent.mt < 5",
"gene_qc": {"min_cells": 0, "excludes": []},
"ccs_args": {},
"qc_plots": {
"Violin Plots": {
"kind": "cell",
"plot_type": "violin",
"devpars": {"res": 100, "height": 600, "width": 1200},
},
"Scatter Plots": {
"kind": "cell",
"plot_type": "scatter",
"devpars": {"res": 100, "height": 800, "width": 1200},
},
"Ridge Plots": {
"kind": "cell",
"plot_type": "ridge",
"devpars": {"res": 100, "height": 800, "width": 1200},
},
"Distribution of number of cells a gene is expressed in": {
"kind": "gene",
"plot_type": "histogram",
"devpars": {"res": 100, "height": 1200, "width": 1200},
},
},
"use_sct": False,
"no_integration": None,
"NormalizeData": {},
"FindVariableFeatures": {},
"ScaleData": {},
"RunPCA": {},
"SCTransform": {
"verbose": True,
},
"IntegrateLayers": {"method": "harmony"},
"contam_correction": None,
"decontX": {},
"scCDC": {"Detection": {}, "Quantification": {}, "Correction": {}},
"keep_contam_assay": False,
"doublet_detector": "none",
"DoubletFinder": {"PCs": 10, "pN": 0.25, "doublets": 0.075, "ncores": 1},
"scDblFinder": {"dbr": 0.075, "ncores": 1},
"cache": config.path.tmpdir,
}
script = "file://../scripts/scrna/SeuratPreparing.R"
plugin_opts = {
"report": "file://../reports/common.svelte",
}
class SeuratClustering(Proc):DOCS
"""Determine the clusters of cells without reference using Seurat FindClusters
procedure.
Input:
srtobj: The seurat object loaded by SeuratPreparing
Output:
outfile: The seurat object with cluster information at `seurat_clusters` or
the name specified by `envs.ident`
Envs:
ncores (type=int;order=-100): Number of cores to use.
Used in `future::plan(strategy = "multicore", workers = <ncores>)`
to parallelize some Seurat procedures.
See also: <https://satijalab.org/seurat/articles/future_vignette.html>
ident: The name in the metadata to save the cluster labels.
A shortcut for `envs["FindClusters"]["cluster.name"]`.
python: The path to the python executable where `umap-learn` is installed.
Needed for `RunUMAP` when `umap.method` is set to `umap-learn`.
RunUMAP (ns): Arguments for [`RunUMAP()`](https://satijalab.org/seurat/reference/runumap).
`object` is specified internally, and `-` in the key will be replaced with `.`.
`dims=N` will be expanded to `dims=1:N`; The maximal value of `N` will be the minimum of `N` and the number of columns - 1 for each sample.
You can also specify `features` instead of `dims` to use specific features for UMAP. It can be
a list with the following fields: `order` (the order of the markers to use for UMAP, e.g. "desc(abs(avg_log2FC))", and
`n` (the number of total features to use for UMAP, e.g. 30).
If `features` is a list, it will run `biopipen.utils::RunSeuratDEAnalysis` to get the markers
for each group, and then select the top `n`/`ngroups` features for each group based on the `order` field.
If `features` is a numeric value, it will be treated as the `n` field
in the list above, with the default `order` being "desc(abs(avg_log2FC))".
- dims (type=int): The number of PCs to use
- reduction: The reduction to use for UMAP.
If not provided, `sobj@misc$integrated_new_reduction` will be used.
- <more>: See <https://satijalab.org/seurat/reference/runumap>
RunPCA (ns): Arguments for [`RunPCA()`](https://satijalab.org/seurat/reference/runpca).
FindNeighbors (ns): Arguments for [`FindNeighbors()`](https://satijalab.org/seurat/reference/findneighbors).
`object` is specified internally, and `-` in the key will be replaced with `.`.
- reduction: The reduction to use.
If not provided, `sobj@misc$integrated_new_reduction` will be used.
- <more>: See <https://satijalab.org/seurat/reference/findneighbors>
FindClusters (ns): Arguments for [`FindClusters()`](https://satijalab.org/seurat/reference/findclusters).
`object` is specified internally, and `-` in the key will be replaced with `.`.
The cluster labels will be saved in cluster names and prefixed with "c".
The first cluster will be "c1", instead of "c0".
- resolution (type=auto): The resolution of the clustering. You can have multiple resolutions as a list or as a string separated by comma.
Ranges are also supported, for example: `0.1:0.5:0.1` will generate `0.1, 0.2, 0.3, 0.4, 0.5`. The step can be omitted, defaulting to 0.1.
The results will be saved in `<ident>_<resolution>`.
The final resolution will be used to define the clusters at `<ident>`.
- <more>: See <https://satijalab.org/seurat/reference/findclusters>
cache (type=auto): Where to cache the information at different steps.
If `True`, the seurat object will be cached in the job output directory, which will be not cleaned up when job is rerunning.
Set to `False` to not cache the results.
Requires:
r-seurat:
- check: {{proc.lang}} <(echo "library(Seurat)")
r-tidyr:
- check: {{proc.lang}} <(echo "library(tidyr)")
r-dplyr:
- check: {{proc.lang}} <(echo "library(dplyr)")
""" # noqa: E501
input = "srtobj:file"
output = "outfile:file:{{in.srtobj | stem}}.qs"
lang = config.lang.rscript
envs = {
"ncores": config.misc.ncores,
"ident": "seurat_clusters",
"python": config.lang.python,
"RunPCA": {},
"RunUMAP": {},
"FindNeighbors": {},
"FindClusters": {"resolution": 0.8},
"cache": config.path.tmpdir,
}
script = "file://../scripts/scrna/SeuratClustering.R"
class SeuratSubClustering(Proc):DOCS
"""Find clusters of a subset of cells.
It's unlike [`Seurat::FindSubCluster`], which only finds subclusters of a single
cluster. Instead, it will perform the whole clustering procedure on the subset of
cells. One can use metadata to specify the subset of cells to perform clustering on.
For the subset of cells, the reductions will be re-performed on the subset of cells,
and then the clustering will be performed on the subset of cells. The reduction
will be saved in `object@reduction$<casename>.<reduction>` of the original object and the
clustering will be saved in the metadata of the original object using the casename
as the column name.
Input:
srtobj: The seurat object in RDS or qs/qs2 format.
Output:
outfile: The seurat object with the subclustering information in qs/qs2 format.
Envs:
ncores (type=int;order=-100): Number of cores to use.
Used in `future::plan(strategy = "multicore", workers = <ncores>)`
to parallelize some Seurat procedures.
python: The path to the python executable where `umap-learn` is installed.
Needed for `RunUMAP` when `umap.method` is set to `umap-learn`.
mutaters (type=json): The mutaters to mutate the metadata to subset the cells.
The mutaters will be applied in the order specified.
See <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
subset: An expression to subset the cells, will be passed to
[`tidyseurat::filter()`](https://stemangiola.github.io/tidyseurat/reference/filter.html).
RunPCA (ns): Arguments for [`RunPCA()`](https://satijalab.org/seurat/reference/runpca).
`object` is specified internally as the subset object, and `-` in the key will be replaced with `.`.
- <more>: See <https://satijalab.org/seurat/reference/runpca>
RunUMAP (ns): Arguments for [`RunUMAP()`](https://satijalab.org/seurat/reference/runumap).
`object` is specified internally as the subset object, and `-` in the key will be replaced with `.`.
`dims=N` will be expanded to `dims=1:N`; The maximal value of `N` will be the minimum of `N` and the number of columns - 1 for each sample.
You can also specify `features` instead of `dims` to use specific features for UMAP. It can be
a list with the following fields: `order` (the order of the markers to use for UMAP, e.g. "desc(abs(avg_log2FC))", and
`n` (the number of total features to use for UMAP, e.g. 30).
If `features` is a list, it will run `biopipen.utils::RunSeuratDEAnalysis` to get the markers
for each group, and then select the top `n`/`ngroups` features for each group based on the `order` field.
If `features` is a numeric value, it will be treated as the `n` field
in the list above, with the default `order` being "desc(abs(avg_log2FC))".
- dims (type=int): The number of PCs to use
- reduction: The reduction to use for UMAP.
If not provided, `sobj@misc$integrated_new_reduction` will be used.
- <more>: See <https://satijalab.org/seurat/reference/runumap>
FindNeighbors (ns): Arguments for [`FindNeighbors()`](https://satijalab.org/seurat/reference/findneighbors).
`object` is specified internally, and `-` in the key will be replaced with `.`.
- reduction: The reduction to use.
If not provided, `object@misc$integrated_new_reduction` will be used.
- <more>: See <https://satijalab.org/seurat/reference/findneighbors>
FindClusters (ns): Arguments for [`FindClusters()`](https://satijalab.org/seurat/reference/findclusters).
`object` is specified internally, and `-` in the key will be replaced with `.`.
The cluster labels will be prefixed with "s". The first cluster will be "s1", instead of "s0".
- resolution (type=auto): The resolution of the clustering. You can have multiple resolutions as a list or as a string separated by comma.
Ranges are also supported, for example: `0.1:0.5:0.1` will generate `0.1, 0.2, 0.3, 0.4, 0.5`. The step can be omitted, defaulting to 0.1.
The results will be saved in `<casename>_<resolution>`.
The final resolution will be used to define the clusters at `<casename>`.
- <more>: See <https://satijalab.org/seurat/reference/findclusters>
cache (type=auto): Whether to cache the results.
If `True`, the seurat object will be cached in the job output directory, which will be not cleaned up when job is rerunning.
Set to `False` to not cache the results.
cases (type=json): The cases to perform subclustering.
Keys are the names of the cases and values are the dicts inherited from `envs` except `mutaters` and `cache`.
If empty, a case with name `subcluster` will be created with default parameters.
The case name will be passed to `biopipen.utils::SeuratSubCluster()` as `name`.
It will be used as the prefix for the reduction name, keys and cluster names.
For reduction keys, it will be `toupper(<name>)` + "PC_" and `toupper(<name>)` + "UMAP_".
For cluster names, it will be `<name>` + "." + resolution.
And the final cluster name will be `<name>`.
Note that the `name` should be alphanumeric and anything other than alphanumeric will be removed.
""" # noqa: E501
input = "srtobj:file"
output = "outfile:file:{{in.srtobj | stem}}.qs"
lang = config.lang.rscript
envs_depth = 1
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"python": config.lang.python,
"subset": None,
"RunPCA": {},
"RunUMAP": {},
"FindNeighbors": {},
"FindClusters": {"resolution": 0.8},
"cache": config.path.tmpdir,
"cases": {},
}
script = "file://../scripts/scrna/SeuratSubClustering.R"
class SeuratClusterStats(Proc):DOCS
"""Statistics of the clustering.
Including the number/fraction of cells in each cluster, the gene expression values
and dimension reduction plots. It's also possible to perform stats on
TCR clones/clusters or other metadata for each T-cell cluster.
Examples:
### Clustree Plot
```toml
[SeuratClusterStats.envs.clustrees."Clustree Plot"]
prefix = "seurat_clusters"
devpars = {height = 500}
```
{: width="80%" }
### Number of cells in each cluster (Bar Chart)
```toml
[SeuratClusterStats.envs.stats."Number of cells in each cluster (Bar Chart)"]
plot_type = "bar"
x_text_angle = 90
```
{: width="80%" }
### Number of cells in each cluster by Sample (Bar Chart)
```toml
[SeuratClusterStats.envs.stats."Number of cells in each cluster by Sample (Bar Chart)"]
plot_type = "bar"
group_by = "Sample"
x_text_angle = 90
```
{: width="80%" }
### Number of cells in each cluster by Diagnosis
```toml
[SeuratClusterStats.envs.stats."Number of cells in each cluster by Diagnosis"]
plot_type = "bar"
group_by = "Diagnosis"
frac = "group"
x_text_angle = 90
swap = true
position = "stack"
```
{: width="80%" }
### Number of cells in each cluster by Diagnosis (Circos Plot)
```toml
[SeuratClusterStats.envs.stats."Number of cells in each cluster by Diagnosis (Circos Plot)"]
plot_type = "circos"
group_by = "Diagnosis"
```
{: width="80%" }
### Number of cells in each cluster by Diagnosis (Sankey Plot)
```toml
[SeuratClusterStats.envs.stats."Number of cells in each cluster by Diagnosis (Sankey Plot)"]
plot_type = "sankey"
group_by = ["seurat_clusters", "Diagnosis"]
links_alpha = 0.6
devpars = {width = 800}
```
{: width="80%" }
### Number of cells in each cluster by Sample (Spider Plot)
```toml
[SeuratClusterStats.envs.stats."Number of cells in each cluster by Sample (Spider Plot)"]
plot_type = "spider"
group_by = "Diagnosis"
palette = "Set1"
```
{: width="80%" }
### Number of genes detected in each cluster
```toml
[SeuratClusterStats.envs.ngenes."Number of genes detected in each cluster"]
plot_type = "violin"
add_box = true
add_point = true
```
{: width="80%" }
### Feature Expression in Clusters (Violin Plots)
```toml
[SeuratClusterStats.envs.features_defaults]
features = ["CD3D", "CD4", "CD8A", "MS4A1", "CD14", "LYZ", "FCGR3A", "NCAM1", "KLRD1"]
[SeuratClusterStats.envs.features."Feature Expression in Clusters (Violin Plots)"]
plot_type = "violin"
ident = "seurat_clusters"
```
{: width="80%" }
### Feature Expression in Clusters (Ridge Plots)
```toml
# Using the same features as above
[SeuratClusterStats.envs.features."Feature Expression in Clusters (Ridge Plots)"]
plot_type = "ridge"
ident = "seurat_clusters"
flip = true
```
{: width="80%" }
### Feature Expression in Clusters by Diagnosis
```toml
# Using the same features as above
[SeuratClusterStats.envs.features."Feature Expression in Clusters by Diagnosis"]
plot_type = "violin"
group_by = "Diagnosis"
ident = "seurat_clusters"
comparisons = true
sig_label = "p.signif"
```
{: width="80%" }
### Feature Expression in Clusters (stacked)
```toml
# Using the same features as above
[SeuratClusterStats.envs.features."Feature Expression in Clusters (stacked)"]
plot_type = "violin"
ident = "seurat_clusters"
add_bg = true
stack = true
add_box = true
```
{: width="80%" }
### CD4 Expression on UMAP
```toml
[SeuratClusterStats.envs.features."CD4 Expression on UMAP"]
plot_type = "dim"
feature = "CD4"
highlight = "seurat_clusters == 'c1'"
```
{: width="80%" }
### Feature Expression in Clusters by Diagnosis (Heatmap)
```toml
[SeuratClusterStats.envs.features."Feature Expression in Clusters by Diagnosis (Heatmap)"]
# Grouped features
features = {"T cell markers" = ["CD3D", "CD4", "CD8A"], "B cell markers" = ["MS4A1"], "Monocyte markers" = ["CD14", "LYZ", "FCGR3A"], "NK cell markers" = ["NCAM1", "KLRD1"]}
plot_type = "heatmap"
ident = "Diagnosis"
columns_split_by = "seurat_clusters"
name = "Expression"
devpars = {height = 560}
```
{: width="80%" }
### Feature Expression in Clusters by Diagnosis (Heatmap with annotations)
```toml
# Using the default features
[SeuratClusterStats.envs.features."Feature Expression in Clusters by Diagnosis (Heatmap with annotations)"]
ident = "seurat_clusters"
cell_type = "dot"
plot_type = "heatmap"
name = "Expression Level"
dot_size = "nanmean"
dot_size_name = "Percent Expressed"
add_bg = true
rows_split_by = "Diagnosis"
cluster_rows = false
flip = true
palette = "YlOrRd"
column_annotation = ["percent.mt", "VDJ_Presence"]
column_annotation_type = {"percent.mt" = "violin", VDJ_Presence = "pie"}
column_annotation_params = {"percent.mt" = {show_legend = false}}
devpars = {width = 1400, height = 900}
```
{: width="80%" }
### Dimensional reduction plot
```toml
[SeuratClusterStats.envs.features."Dimensional reduction plot"]
label = true
```
{: width="80%" }
### Dimensional reduction plot (with marks)
```toml
[SeuratClusterStats.envs.dimplots."Dimensional reduction plot (with marks)"]
add_mark = true
mark_linetype = 2
```
{: width="80%" }
### Dimensional reduction plot (with hex bins)
```toml
[SeuratClusterStats.envs.dimplots."Dimensional reduction plot (with hex bins)"]
hex = true
hex_bins = 50
```
{: width="80%" }
### Dimensional reduction plot (with Diagnosis stats)
```toml
[SeuratClusterStats.envs.dimplots."Dimensional reduction plot (with Diagnosis stats)"]
stat_by = "Diagnosis"
stat_plot_type = "ring"
stat_plot_size = 0.15
```
{: width="80%" }
### Dimensional reduction plot by Diagnosis
```toml
[SeuratClusterStats.envs.dimplots."Dimensional reduction plot by Diagnosis"]
facet_by = "Diagnosis"
highlight = true
theme = "theme_blank"
```
{: width="80%" }
Input:
srtobj: The seurat object loaded by `SeuratClustering`
Output:
outdir: The output directory.
Different types of plots will be saved in different subdirectories.
For example, `clustree` plots will be saved in `clustrees` subdirectory.
For each case in `envs.clustrees`, both the png and pdf files will be saved.
Envs:
ncores (type=int): Number of cores to use for reading and writing the data.
mutaters (type=json): The mutaters to mutate the metadata to subset the cells.
The mutaters will be applied in the order specified.
You can also use the clone selectors to select the TCR clones/clusters.
See <https://pwwang.github.io/scplotter/reference/clone_selectors.html>.
You can also use key `<newcol>:ident` to set the `<newcol>` as the default
ident for the stats.
See also <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
cache (type=auto): Whether to cache the plots.
Currently only plots for features are supported, since creating the those
plots can be time consuming.
If `True`, the plots will be cached in the job output directory, which will
be not cleaned up when job is rerunning.
clustrees_defaults (ns): The parameters for the clustree plots.
- devpars (ns): The device parameters for the clustree plot.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- more_formats (type=list): The formats to save the plots other than `png`.
- save_code (flag): Whether to save the code to reproduce the plot.
- prefix (type=auto): string indicating columns containing clustering information.
The trailing dot is not necessary and will be added automatically.
When `TRUE`, clustrees will be plotted when there is `FindClusters` or
`FindClusters.*` in the `obj@commands`.
The latter is generated by `SeuratSubClustering`.
This will be ignored when `envs.clustrees` is specified
(the prefix of each case must be specified separately).
- <more>: Other arguments passed to `scplotter::ClustreePlot`.
See <https://pwwang.github.io/scplotter/reference/ClustreePlot.html>
clustrees (type=json): The cases for clustree plots.
Keys are the names of the plots and values are the dicts inherited from `env.clustrees_defaults` except `prefix`.
There is no default case for `clustrees`.
stats_defaults (ns): The default parameters for `stats`.
This is to do some basic statistics on the clusters/cells. For more comprehensive analysis,
see <https://pwwang.github.io/scplotter/reference/CellStatPlot.html>.
The parameters from the cases can overwrite the default parameters.
- subset: An expression to subset the cells, will be passed to `tidyrseurat::filter()`.
- devpars (ns): The device parameters for the clustree plot.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- descr: The description of the plot, showing in the report.
- more_formats (type=list): The formats to save the plots other than `png`.
- save_code (flag): Whether to save the code to reproduce the plot.
- save_data (flag): Whether to save the data used to generate the plot.
- <more>: Other arguments passed to `scplotter::CellStatPlot`.
See <https://pwwang.github.io/scplotter/reference/CellStatPlot.html>.
stats (type=json): The number/fraction of cells to plot.
Keys are the names of the plots and values are the dicts inherited from `env.stats_defaults`.
ngenes_defaults (ns): The default parameters for `ngenes`.
The default parameters to plot the number of genes expressed in each cell.
- more_formats (type=list): The formats to save the plots other than `png`.
- subset: An expression to subset the cells, will be passed to `tidyrseurat::filter()`.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
ngenes (type=json): The number of genes expressed in each cell.
Keys are the names of the plots and values are the dicts inherited from `env.ngenes_defaults`.
features_defaults (ns): The default parameters for `features`.
- features (type=auto): The features to plot.
It can be either a string with comma separated features, a list of features, a file path with `file://` prefix with features
(one per line), or an integer to use the top N features from `VariantFeatures(srtobj)`.
It can also be a dict with the keys as the feature group names and the values as the features, which
is used for heatmap to group the features.
- order_by (type=auto): The order of the clusters to show on the plot.
An expression passed to `dplyr::arrange()` on the grouped meta data frame (by `ident`).
For example, you can order the clusters by the activation score of
the cluster: `desc(mean(ActivationScore, na.rm = TRUE))`, suppose you have a column
`ActivationScore` in the metadata.
You may also specify the literal order of the clusters by a list of strings (at least two).
- subset: An expression to subset the cells, will be passed to `tidyrseurat::filter()`.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- descr: The description of the plot, showing in the report.
- more_formats (type=list): The formats to save the plots other than `png`.
- save_code (flag): Whether to save the code to reproduce the plot.
- save_data (flag): Whether to save the data used to generate the plot.
- <more>: Other arguments passed to `scplotter::FeatureStatPlot`.
See <https://pwwang.github.io/scplotter/reference/FeatureStatPlot.html>
features (type=json): The plots for features, include gene expressions, and columns from metadata.
Keys are the titles of the cases and values are the dicts inherited from `env.features_defaults`.
dimplots_defaults (ns): The default parameters for `dimplots`.
You can also set `dims` to `[1, 2, 3]` to make a 3D dimensional reduction plot.
Note that the default `n.components` for `RunUMAP` is 2, so if you want to make a 3D UMAP plot,
you need to set `n.components` to 3 in `RunUMAP` parameters in `SeuratSubClustering` or `SeuratClustering`.
- group_by: The identity to use.
- split_by: The column name in metadata to split the cells into different plots.
- subset: An expression to subset the cells, will be passed to `tidyrseurat::filter()`.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- reduction (choice): Which dimensionality reduction to use.
- dim: If `ident` is from subclustering, `<group_by>.<reduction>` will be used.
Otherwise, `scplotter:::default_dimreduc(obj)` will be used to determine
the default dimensionality reduction to use based on the reductions available in the object.
- auto: Same as `dim`
- umap: Use `Seurat::UMAPPlot`.
- tsne: Use `Seurat::TSNEPlot`.
- pca: Use `Seurat::PCAPlot`.
- <more>: See <https://pwwang.github.io/scplotter/reference/CellDimPlot.html>
dimplots (type=json): The dimensional reduction plots.
Keys are the titles of the plots and values are the dicts inherited from `env.dimplots_defaults`. It can also have other parameters from
[`scplotter::CellDimPlot`](https://pwwang.github.io/scplotter/reference/CellDimPlot.html).
Requires:
r-seurat:
- check: {{proc.lang}} -e "library(Seurat)"
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}.cluster_stats"
lang = config.lang.rscript
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"cache": config.path.tmpdir,
"clustrees_defaults": {
"devpars": {"res": 100},
"more_formats": [],
"save_code": False,
"prefix": True,
},
"clustrees": {},
"stats_defaults": {
"subset": None,
"descr": None,
"devpars": {"res": 100},
"more_formats": [],
"save_code": False,
"save_data": False,
},
"stats": {
"Number of cells in each cluster (Bar Chart)": {
"plot_type": "bar",
"x_text_angle": 90,
},
"Number of cells in each cluster by Sample (Bar Chart)": {
"plot_type": "bar",
"group_by": "Sample",
"x_text_angle": 90,
},
},
"ngenes_defaults": {
"subset": None,
"more_formats": [],
"add_box": True,
"devpars": {"res": 100, "height": 800, "width": 1000},
},
"ngenes": {
"Number of genes expressed in each cluster": {},
},
"features_defaults": {
"features": None,
"order_by": None,
"subset": None,
"devpars": {"res": 100},
"descr": None,
"more_formats": [],
"save_code": False,
"save_data": False,
},
"features": {},
"dimplots_defaults": {
"group_by": None, # use default ident
"split_by": None,
"subset": None,
"reduction": "dim",
"devpars": {"res": 100},
},
"dimplots": {
"Dimensional reduction plot": {
"label": True,
},
},
}
script = "file://../scripts/scrna/SeuratClusterStats.R"
plugin_opts = {
"report": "file://../reports/common.svelte",
"report_paging": 8,
}
class HdWGCNA(Proc):DOCS
"""Weighted gene co-expression network analysis (WGCNA) for single-cell RNA-seq data.
Wrapping [`hdWGCNA`](https://smorabit.github.io/hdWGCNA/index.html) (v0.4.12).
For a comprehensive introduction and usage guide, please refer to
the [hdWGCNA tutorials](https://smorabit.github.io/hdWGCNA/articles/index.html).
The analysis includes the following steps:
1. **Setup**: `SetupForWGCNA` to select the genes for WGCNA.
2. **Data aggregation** (choose one):
- **Metacells** (recommended): `MetacellsByGroups` + `NormalizeMetacells`, with optional
`ScaleMetacells`, `RunPCAMetacells`, `RunHarmonyMetacells` and `RunUMAPMetacells`.
- **Pseudobulk**: `AggregatePseudobulk` + `NormalizeCounts` (DESeq2 VST).
3. **Network construction**: `SetDatExpr` + `TestSoftPowers` + `ConstructNetwork`, or
`SetMultiExpr` + `TestSoftPowersConsensus` + `ConstructNetwork` for a **consensus network**
(`use_consensus`).
4. **Module analysis**: `ModuleEigengenes` + `ModuleConnectivity` (kME), optionally
`ResetModuleNames`/`ResetModuleColors`, and `GetHubGenes` for hub genes.
5. **Downstream analyses**: module visualization (`plots`), differential module expression
(`dmes`), module-trait correlation (`module_trait_corr`), functional enrichment (`enrich`),
module preservation/projection (`module_preservations`), TF regulatory networks
(`tf_network`) and motif overlap (`motifs`).
The envs that configure a single hdWGCNA function are named after the function itself
(e.g. `SetupForWGCNA`, `MetacellsByGroups`), and the keys within them are the argument
names of the corresponding R function, with `-` used in place of `.` in the argument
names (e.g. `group.by` -> `group-by`). The envs that bundle several functions
(`plots`, `dmes`, `module_trait_corr`, `enrich`, `module_preservations`, `tf_network`,
`motifs`) are case-based, where the keys are the titles of the cases and the values
inherit from the corresponding `*_defaults` env.
Packages:
**Core** (required for any analysis, checked by `Requires`):
[Seurat](https://satijalab.org/seurat/), [hdWGCNA](https://github.com/smorabit/hdWGCNA),
[WGCNA](https://cran.r-project.org/package=WGCNA), [igraph](https://igraph.org/),
[patchwork](https://patchwork.data-imaginist.com/),
[tidyseurat](https://tidyseurat.github.io/), [glue](https://glue.tidyverse.org/) and
[biopipen.utils](https://github.com/pwwang/biopipen.utils) (provides `RunEnrichment`/
`VizEnrichment`).
**Advanced** (only when the corresponding analyses are enabled):
- **Enrichment** (`enrich`): [enrichit](https://github.com/pwwang/enrichit)
(via `biopipen.utils::RunEnrichment`, requires internet access for the
Enrichr API).
- **GSEA** (`gsea`): [fgsea](https://bioconductor.org/packages/fgsea/);
[msigdbr](https://cran.r-project.org/package=msigdbr) only when
`genesets` is given as a dict of arguments for `msigdbr::msigdbr()`
(e.g. `{species: "human", collection: "C5", subcollection: "BP"}`) —
a GMT file needs no extra package.
- **Motif scan** (`tf_network.motif_scan` and `motifs`):
[TFBSTools](https://bioconductor.org/packages/TFBSTools/),
[motifmatchr](https://bioconductor.org/packages/motifmatchr/),
a JASPAR package (e.g. [JASPAR2020](https://bioconductor.org/packages/JASPAR2020/)),
a BSgenome package (e.g. `BSgenome.Hsapiens.UCSC.hg38`, set via
`species_genome`), and an EnsDb package (e.g. `EnsDb.Hsapiens.v86`,
set via `ensdb_package`).
- **TF network** (`tf_network.construct`): [xgboost](https://cran.r-project.org/package=xgboost).
- **Regulon scores** (`tf_network.regulon_scores`, and `ModuleExprScore`
with `method: "UCell"`): [UCell](https://github.com/carmonalab/UCell).
- **Regulon enrichment** (`tf_network.enrich_regulons`): [enrichit](https://github.com/pwwang/enrichit)
(via `biopipen.utils::RunEnrichment`, requires internet access for the
Enrichr API).
- **Regulatory heatmap** (`tf_network.regulatory_heatmap`): [forcats](https://forcats.tidyverse.org/)
(hdWGCNA uses `fct_rev()` in the heatmap's aes without importing forcats).
- **Module preservation** (`module_preservations.netrep`): [NetRep](https://cran.r-project.org/package=NetRep).
- **Metacells** (`RunHarmonyMetacells`): [harmony](https://cran.r-project.org/package=harmony).
- **Pseudobulk** (`NormalizeCounts` with `method: "VST"`): [DESeq2](https://bioconductor.org/packages/DESeq2).
Input:
srtobj: The Seurat object to analyze
Output:
outdir: The output directory
Envs:
ncores (type=int): The number of cores to use.
mutaters (type=json): The columns to add to the metadata of the object.
See `SeuratPreparing` for more details.
subset: An expression to subset the cells, will be passed to `tidyrseurat::filter()`.
seed (type=int): The random seed.
cache: The directory to cache the TOM matrices (they are large). `true` to use the job
output directory. By default, the TOM matrices are cached in the output directory.
use_pseudobulk (flag): Whether to aggregate the cells into pseudobulk
(via `AggregatePseudobulk` + `NormalizeCounts`) instead of metacells.
use_consensus (flag): Whether to construct a consensus network
(requires `SetMultiExpr` to be set).
ref_srtobj: The reference Seurat object (RDS/qs file) for module preservation/projection.
Relative paths are resolved against the directory of the input object.
Use `"self"` to test the preservation of the modules against the query
object itself (self-preservation test).
SetupForWGCNA (type=json): The arguments for [`hdWGCNA::SetupForWGCNA`](https://smorabit.github.io/hdWGCNA/reference/SetupForWGCNA.html).
- wgcna_name: The name of the WGCNA network, will be used throughout the analysis.
- gene_select: The gene selection method.
- fraction: The fraction of cells that a gene needs to be expressed in to be
included when `gene_select = "fraction"`.
MetacellsByGroups (type=json): The arguments for [`hdWGCNA::MetacellsByGroups`](https://smorabit.github.io/hdWGCNA/reference/MetacellsByGroups.html).
- group-by: The metadata columns to group the cells by.
- ident-group: The column to use as the identity of the metacells.
NormalizeMetacells (type=json): The arguments for [`hdWGCNA::NormalizeMetacells`](https://smorabit.github.io/hdWGCNA/reference/NormalizeMetacells.html).
ScaleMetacells: The arguments for [`hdWGCNA::ScaleMetacells`](https://smorabit.github.io/hdWGCNA/reference/ScaleMetacells.html).
`None` to skip.
RunPCAMetacells: The arguments for [`hdWGCNA::RunPCAMetacells`](https://smorabit.github.io/hdWGCNA/reference/RunPCAMetacells.html).
`None` to skip.
RunHarmonyMetacells: The arguments for [`hdWGCNA::RunHarmonyMetacells`](https://smorabit.github.io/hdWGCNA/reference/RunHarmonyMetacells.html).
`None` to skip. Requires `group-by-vars` to be set.
RunUMAPMetacells: The arguments for [`hdWGCNA::RunUMAPMetacells`](https://smorabit.github.io/hdWGCNA/reference/RunUMAPMetacells.html).
`None` to skip.
AggregatePseudobulk (type=json): The arguments for [`hdWGCNA::AggregatePseudobulk`](https://smorabit.github.io/hdWGCNA/reference/AggregatePseudobulk.html).
- replicate_col: The column in metadata for the replicate (sample) of each cell.
- group_col: The column in metadata for the group of each cell.
NormalizeCounts (type=json): The arguments for [`hdWGCNA::NormalizeCounts`](https://smorabit.github.io/hdWGCNA/reference/NormalizeCounts.html).
SetDatExpr (type=json): The arguments for [`hdWGCNA::SetDatExpr`](https://smorabit.github.io/hdWGCNA/reference/SetDatExpr.html).
Note that in the pseudobulk mode, `mat` and `layer` are set automatically.
TestSoftPowers (type=json): The arguments for [`hdWGCNA::TestSoftPowers`](https://smorabit.github.io/hdWGCNA/reference/TestSoftPowers.html).
SetMultiExpr: The arguments for [`hdWGCNA::SetMultiExpr`](https://smorabit.github.io/hdWGCNA/reference/SetMultiExpr.html).
Required when `use_consensus` is true. When `use_pseudobulk` is true,
`mat` and `layer` are set automatically, and `multi-group-by` should be
set to the column in the pseudobulk metadata that defines the groups
(e.g. `cell_type`).
TestSoftPowersConsensus (type=json): The arguments for [`hdWGCNA::TestSoftPowersConsensus`](https://smorabit.github.io/hdWGCNA/reference/TestSoftPowersConsensus.html).
ConstructNetwork (type=json): The arguments for [`hdWGCNA::ConstructNetwork`](https://smorabit.github.io/hdWGCNA/reference/ConstructNetwork.html).
- soft_power: The soft power to use. `None` to use the one determined by `TestSoftPowers`.
- tom_name: The name of the TOM matrix. Defaults to the `wgcna_name`.
- minModuleSize: The minimum number of genes in a module.
- mergeCutHeight: The cut height for merging similar modules.
ModuleEigengenes (type=json): The arguments for [`hdWGCNA::ModuleEigengenes`](https://smorabit.github.io/hdWGCNA/reference/ModuleEigengenes.html).
ModuleConnectivity (type=json): The arguments for [`hdWGCNA::ModuleConnectivity`](https://smorabit.github.io/hdWGCNA/reference/ModuleConnectivity.html).
ResetModuleNames: The arguments for [`hdWGCNA::ResetModuleNames`](https://smorabit.github.io/hdWGCNA/reference/ResetModuleNames.html).
`None` to skip. e.g. `{new_name: "M"}` to rename the modules to `M1`, `M2`, ...
ResetModuleColors: The arguments for [`hdWGCNA::ResetModuleColors`](https://smorabit.github.io/hdWGCNA/reference/ResetModuleColors.html).
`None` to skip.
GetHubGenes (type=json): The arguments for [`hdWGCNA::GetHubGenes`](https://smorabit.github.io/hdWGCNA/reference/GetHubGenes.html).
plots_defaults (ns): The default parameters for the cases in `plots`.
- kind: The kind of the plot. See `plots` for the available kinds.
- devpars (ns): The device parameters for the plots.
- more_formats (type=list): The formats to save the plots other than `png`.
plots (type=json): The plots to make. Keys are the titles of the plots and values
are the dicts inherited from `plots_defaults`.
- kind: The kind of the plot, one of `soft_powers` (the soft power selection plots,
SFT and fit, saved as `<title>.sft.png` and `<title>.fit.png`), `dendrogram`
(the gene dendrogram), `kmes` (the kME module membership plots, arguments for
`PlotKMEs`), `module_umap` (the module UMAP, arguments for `RunModuleUMAP` at the
case level except `umap_plot_args`, and the arguments for `ModuleUMAPPlot` in the
sub-dict `umap_plot_args`), `module_network` (the module network plots,
`ModuleNetworkPlot` writes the plots to `<outdir>/plots/ModuleNetworks.<title>`),
`hub_gene_network`, `module_feature` (arguments for `ModuleFeaturePlot`),
`module_radar` (arguments for `ModuleRadarPlot`), `module_correlogram`
(arguments for `ModuleCorrelogram`), `consensus_compare` (only when
`use_consensus` is true: re-runs the standard non-consensus workflow on the
same data under `standard_name`, then plots both the consensus and standard
module colors on the consensus dendrogram via `WGCNA::plotDendroAndColors`;
accepts `standard_name` (required), `soft_power` for the standard network
(given, the soft power test is skipped), and other `SetDatExpr`/`ConstructNetwork`
arguments, plus `main` for the plot title).
- descr: The description of the plot, showing in the report.
ModuleExprScore: The arguments for [`hdWGCNA::ModuleExprScore`](https://smorabit.github.io/hdWGCNA/reference/ModuleExprScore.html).
Runs before `dmes` so that the module scores can be used as the features for DMEs.
`None` to skip.
dmes_defaults (ns): The default parameters for the cases in `dmes`.
- mode (choice): `find_all` to use `FindAllDMEs`, `find` to use `FindDMEs`.
- group-by: The column in metadata to group the cells by (for `find_all`).
- barcodes1: The barcodes of the first group to compare (for `find`).
- barcodes2: The barcodes of the second group to compare (for `find`).
- lollipop (flag): Whether to make the lollipop plot.
- volcano (flag): Whether to make the volcano plot.
- devpars (ns): The device parameters for the plots.
dmes (type=json): The differential module expression analyses. Keys are the titles of the
cases and values are the dicts inherited from `dmes_defaults`. Other arguments of
`FindAllDMEs`/`FindDMEs` can also be set, plus `pvalue` for the p-value column
used by the lollipop plot (default: `p_val_adj`).
module_trait_corr_defaults (ns): The default parameters for the cases in `module_trait_corr`.
- traits: The columns in metadata to use as the traits (a character vector).
- group-by: The column in metadata to group the cells by.
- devpars (ns): The device parameters for the plots.
module_trait_corr (type=json): The module-trait correlation analyses. Keys are the titles
of the cases and values are the dicts inherited from `module_trait_corr_defaults`.
enrich_defaults (ns): The default parameters for the cases in `enrich`.
- dbs: The databases to use for the enrichment analysis (via [enrichit](https://github.com/pwwang/enrichit)).
- enrich_plots (type=json): The plots to make for each module. Keys are the titles of
the plots and values are the arguments for
[`plotthis::VizEnrichment`](https://pwwang.github.io/plotthis/reference/VizEnrichment.html),
including `db`, `plot_type`, `top_term`, `metric`, etc.
- devpars (ns): The device parameters for the plots.
enrich (type=json): The functional enrichment analyses for the module genes. Keys are the
titles of the cases and values are the dicts inherited from `enrich_defaults`.
The enrichment is done by `enrichit::RunEnrichment` and the plots by `plotthis::VizEnrichment`.
Requires internet access.
gsea_defaults (ns): The default parameters for the cases in `gsea`.
- genesets: The gene sets to test against. A dict of arguments for
[`msigdbr::msigdbr()`](https://igordot.github.io/msigdbr/reference/msigdbr.html)
(e.g. `{species: "human", collection: "C5", subcollection: "BP"}`), or a GMT file.
- metric: The significance metric used to rank the terms, `p.adjust` or `pvalue`.
- cutoff: The significance cutoff for the terms (on the `metric` column).
- top_term (type=int): The number of top terms to show in the summary plot, and
the default number of terms for the running-score plots.
- gsea_plots (type=json): The running-score plots. Keys are the titles and values are
the arguments for
[`plotthis::VizGSEA`](https://pwwang.github.io/plotthis/reference/GSEAPlot.html)
with `plot_type: "gsea"` (the summary plot is always made),
including `gs` (the gene set IDs to plot; default: the top significant terms),
`top_term`, etc.
- devpars (ns): The device parameters for the plots.
gsea (type=json): The gene set enrichment analyses of the module gene ranks (kME) for
each module. Keys are the titles of the cases and values are the dicts inherited
from `gsea_defaults`. Other arguments of
[`biopipen.utils::RunGSEA`](https://pwwang.github.io/biopipen.utils/reference/RunGSEA.html)
(e.g. `minSize`, `maxSize`) can also be set. The analysis is done by `fgsea` and
the plots by `plotthis::VizGSEA` (`GSEASummaryPlot`/`GSEAPlot`).
module_preservations_defaults (ns): The default parameters for the cases in `module_preservations`.
- project_modules (flag): Whether to project the modules to the reference object
(via `ProjectModules`).
- preserve (flag): Whether to test the module preservation
(via `ModulePreservation`).
- project_args (ns): The arguments for `ProjectModules`.
- preserve_args (ns): The arguments for `ModulePreservation`.
- plot (flag): Whether to make the preservation plot.
- lollipop (flag): Whether to make the lollipop plot.
- netrep (ns): The network reproduction test (via `ModulePreservationNetRep`). `{}` to skip.
- args: The arguments for `ModulePreservationNetRep`.
- topology_heatmap (type=json): The topology heatmaps. Keys are the titles and values
are the arguments for `ModuleTopologyHeatmap`, with `mod` required.
- topology_barplot (type=json): The topology barplots. Keys are the titles and values
are the arguments for `ModuleTopologyBarplot`, with `mod` required.
- devpars (ns): The device parameters for the plots.
module_preservations (type=json): The module preservation/projection analyses. Keys are the
titles of the cases and values are the dicts inherited from `module_preservations_defaults`.
Requires `ref_srtobj` to be set.
tf_network (ns): The TF regulatory network analysis. `{}` to skip.
- motif_scan (ns): The arguments for `MotifScan`. The `pfm` and `EnsDb` are constructed
from the JASPAR and Ensembl databases automatically.
- species_genome: The genome build name, e.g. `hg38` or `mm39`.
- ensdb_package: The name of the EnsDb package, e.g. `EnsDb.Hsapiens.v86`.
- jaspar_package: The name of the JASPAR package, default `JASPAR2020`.
- collection: The JASPAR collection, e.g. `CORE`.
- tax_group: The JASPAR taxonomic group, e.g. `vertebrates`.
- all_versions: Whether to include all versions of the JASPAR matrices.
- construct (ns): The arguments for `ConstructTFNetwork`, with `model_params` required.
`nthread` in `model_params` is set to `ncores` automatically if not set.
- assign (ns): The arguments for `AssignTFRegulons`.
- regulon_scores (ns): The arguments for `RegulonScores`. `{}` to skip.
- plots (type=json): The TF network plots. Keys are the titles and values are the
arguments for `TFNetworkPlot`, with `selected_tfs` required.
- regulon_bar_plots (type=json): The regulon bar plots. Keys are the titles and values
are the arguments for `RegulonBarPlot`, with `selected_tf` required.
- differential_regulons (type=json): The differential regulon analyses. Keys are the
titles and values are the arguments for `FindDifferentialRegulons`, with
`barcodes1` and `barcodes2` required.
- enrich_regulons (ns): The regulon enrichment of the positive/negative TF targets
against the Enrichr databases, using enrichit (via `biopipen.utils::RunEnrichment`,
requires internet access). `None` to skip.
- dbs: The Enrichr databases to test against.
- depth: The depth of the TF regulatory network to include.
- use_regulons: Whether to use only the regulon genes (instead of all targets).
- min_genes: The minimum number of target genes required to run the enrichment.
- regulatory_heatmap (ns): The arguments for `ModuleRegulatoryHeatmap`. `None` to skip.
- regulatory_network (ns): The arguments for `ModuleRegulatoryNetworkPlot`. `None` to skip.
motifs (ns): The motif overlap analysis. `{}` to skip. Requires `tf_network.motif_scan` to be run.
- overlap_bar_plot (ns): The arguments for `MotifOverlapBarPlot`. The plots are written
to `<outdir>/plots/MotifOverlap`.
Requires:
r-seurat:
- check: {{proc.lang}} <(echo "library(Seurat)")
r-hdwgcn:
- check: {{proc.lang}} <(echo "library(hdWGCNA)")
r-wgcna:
- check: {{proc.lang}} <(echo "library(WGCNA)")
r-igraph:
- check: {{proc.lang}} <(echo "library(igraph)")
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}.hdwgcn"
lang = config.lang.rscript
envs_depth = 4
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"subset": None,
"seed": 8525,
"cache": None,
"use_pseudobulk": False,
"use_consensus": False,
"ref_srtobj": None,
"SetupForWGCNA": {
"wgcna_name": "biopipen",
"gene_select": "fraction",
"fraction": 0.05,
},
"MetacellsByGroups": {
"group-by": "seurat_clusters",
"ident-group": "seurat_clusters",
},
"NormalizeMetacells": {},
"ScaleMetacells": None,
"RunPCAMetacells": None,
"RunHarmonyMetacells": None,
"RunUMAPMetacells": None,
"AggregatePseudobulk": {
"replicate_col": "Sample",
"group_col": "seurat_clusters",
},
"NormalizeCounts": {"method": "VST", "assay_name": "counts"},
"SetDatExpr": {},
"TestSoftPowers": {},
"SetMultiExpr": None,
"TestSoftPowersConsensus": {},
"ConstructNetwork": {},
"ModuleEigengenes": {},
"ModuleConnectivity": {},
"ResetModuleNames": None,
"ResetModuleColors": None,
"GetHubGenes": {"n_hubs": 10},
"plots_defaults": {
"kind": None,
"devpars": {"res": 100},
"more_formats": [],
},
"plots": {
"Soft Powers": {"kind": "soft_powers"},
"Dendrogram": {"kind": "dendrogram"},
},
"ModuleExprScore": None,
"dmes_defaults": {
"mode": "find_all",
"group-by": "seurat_clusters",
"barcodes1": None,
"barcodes2": None,
"lollipop": True,
"volcano": False,
"devpars": {"res": 100},
},
"dmes": {},
"module_trait_corr_defaults": {
"traits": None,
"group-by": "seurat_clusters",
"devpars": {"res": 100},
},
"module_trait_corr": {},
"enrich_defaults": {
"dbs": None,
"enrich_plots": {},
"devpars": {"res": 100},
},
"enrich": {},
"gsea_defaults": {
"genesets": None,
"metric": "p.adjust",
"cutoff": 0.05,
"top_term": 10,
"gsea_plots": {},
"devpars": {"res": 100},
},
"gsea": {},
"module_preservations_defaults": {
"project_modules": True,
"preserve": True,
"project_args": {},
"preserve_args": {},
"plot": True,
"lollipop": False,
"netrep": {},
"devpars": {"res": 100},
},
"module_preservations": {},
"tf_network": {},
"motifs": {},
}
script = "file://../scripts/scrna/HdWGCNA.R"
plugin_opts = {
"report": "file://../reports/common.svelte",
"report_paging": 2,
}
class ModuleScoreCalculator(Proc):DOCS
"""Calculate the module scores for each cell
The module scores are calculated by
[`biopipen.utils::RunModuleScoring()`](https://pwwang.github.io/biopipen.utils.R/reference/RunModuleScoring.html)
with the scoring method specified by `env.defaults.method` (or per module
by `method` in the module dict):
- `seurat`: [`Seurat::AddModuleScore()`](https://satijalab.org/seurat/reference/addmodulescore),
the default. The module scores are calculated as the average expression
levels of each program on single cell level, subtracted by the
aggregated expression of control feature sets. All analyzed features
are binned based on averaged expression, and the control features are
randomly selected from each bin. (Tirosh I, et al. 2016. Dissecting the
multicellular ecosystem of metastatic melanoma by single-cell RNA-seq.
*Science* 352(6282):189-196.
<https://www.science.org/doi/10.1126/science.aad0501>)
- `ucell`: [`UCell::AddModuleScore_UCell()`](https://bioconductor.org/packages/release/bioc/html/UCell.html)
(Andreatta M, Carmona SJ. 2021. UCell: Robust and scalable single-cell
gene signature scoring. *Comput Struct Biotechnol J* 19:3796-3798.
<https://doi.org/10.1016/j.csbj.2021.06.043>). Missing genes are
imputed with expression 0 (with a warning).
- `aucell`: [`AUCell::AUCell_calcAUC()`](https://bioconductor.org/packages/release/bioc/html/AUCell.html)
(Aibar S, et al. 2017. SCENIC: single-cell regulatory network inference
and clustering. *Nat Methods* 14:1083-1086.
<https://doi.org/10.1038/nmeth.4463>)
- `ssgsea`: [`GSVA::gsva()`](https://bioconductor.org/packages/release/bioc/html/GSVA.html)
with `method = "ssgsea"` (Barbie DA, et al. 2009. Systematic RNA
interference reveals that oncogenic KRAS-driven cancers require TBK1.
*Nature* 462:108-112. <https://doi.org/10.1038/nature08460>)
- `jasmine`: (Noureen N, et al. 2022. Integrated analysis of telomerase
enzymatic activity unravels an association with cancer stemness and
proliferation. *eLife* 11:e71994. <https://doi.org/10.7554/eLife.71994>)
- `scse`: (Pont F, et al. 2019. Single-cell signature explorer for
personalized transcriptomics studies and drug discovery. *Nucleic Acids
Res* 47(19):e90. <https://doi.org/10.1093/nar/gkz601>)
- `scps`: (the scPS benchmarking study:
<https://academic.oup.com/nargab/article/6/3/lqae124/7770961>)
Scores from different methods are not comparable with each other — only
the column names are consistent.
Input:
srtobj: The seurat object loaded by `SeuratClustering`
Output:
rdsfile: The seurat object with module scores added to the metadata.
Envs:
defaults (ns): The default parameters for `modules`.
- method (choice): The scoring method to use, one of `seurat`,
`ucell`, `aucell`, `ssgsea`, `jasmine`, `scse` or `scps`.
Can be overridden per module.
- features: The features (genes) to calculate the scores.
A comma-separated string of genes, e.g.
`"HAVCR2,ENTPD1,LAYN,LAG3"`, yields one score column named by
the module key. A list of gene vectors, e.g.
`["HAVCR2","ENTPD1"]`, yields one column per element, named
`{key}1`, `{key}2`, ... if unnamed, or `{key}_{name}` if
named. You can also specify `cc.genes`,
`cc.genes.updated.2019` or `cc.genes.mouse` (or use
`kind: "cc"`, with `features` defaulting to `cc.genes`) to
calculate cell cycle scores. Three columns will be added to
the metadata: `{key}_S.Score`, `{key}_G2M.Score` and
`{key}_Phase`. This works for all methods. Use one of the
reserved no-prefix keys (`"_"`, `"-"`, `"*"` or `"#"`) as the
module key to keep the plain `S.Score`, `G2M.Score` and
`Phase` names. For diffusion map modules (`kind: "dm"`), this
is the number of components to keep (default 2).
- nbin (type=int): Number of bins of aggregate expression levels
for all analyzed features. Only for the `seurat` method.
- ctrl (type=int): Number of control features selected from
the same bin per analyzed feature. Only for the `seurat`
method.
- k (flag): Use feature clusters returned from `DoKMeans`.
Only for the `seurat` method.
- assay: The assay to use (for tools that accept it).
- seed (type=int): Set a random seed. Only for the `seurat` method.
- search (flag): Search for symbol synonyms for features that
don't match features in object? Only for the `seurat` method.
- <more>: Other parameters, passed to the underlying tool of the
`method`. For `seurat`, they go to `Seurat::AddModuleScore()`
or `Seurat::CellCycleScoring()` (see
<https://satijalab.org/seurat/reference/addmodulescore> and
<https://satijalab.org/seurat/reference/cellcyclescoring>).
For `ucell`: `maxRank`, `w_neg` and `slot` (default
`"counts"`, not `layer`). For `aucell`: `aucMaxRank`,
`plotStats`. For `ssgsea`: `kcdf`, `verbose`, `min.sz`,
`max.sz`, `tau`, etc. For `scse`/`scps`: `layer` (default
`"data"`). `scps` uses the bundled scPS implementation, which
runs PCA on the `scale.data` of the object, so the signature
genes must be scaled first (`ScaleData(features = ...)` or
`SCTransform`); it also requires the
[`GSEABase`](https://bioconductor.org/packages/release/bioc/html/GSEABase.html)
package. For diffmap modules (`kind: "dm"`): `n_pcs` (use PCA
embeddings instead of assay data) and other arguments passed
to
[`destiny::DiffusionMap()`](https://www.rdocumentation.org/packages/destiny/versions/2.0.4/topics/DiffusionMap%20class).
The `agg` and `keep` parameters from old versions are removed
and ignored.
ncores (type=int): The number of cores to use for reading and writing the seurat object.
modules (type=json): The modules to calculate the scores.
Keys are the names of the expression programs and values are the
dicts inherited from `env.defaults`.
Here are some examples -
>>> {
>>> "CellCycleMouse": {"features": "cc.genes.mouse"},
>>> "CellCycle": {"kind": "cc", "features": "cc.genes.updated.2019"},
>>> "TcellState": {
>>> "features": {
>>> "Exhaustion": ["HAVCR2", "ENTPD1", "LAYN", "LAG3"],
>>> "Activation": ["IFNG"]
>>> },
>>> "method": "ucell",
>>> "maxRank": 500
>>> },
>>> "Proliferation": {"features": "STMN1,TUBB"},
>>> "DC": {"kind": "dm"}
>>> }
For `CellCycle`, the columns `CellCycle_S.Score`,
`CellCycle_G2M.Score` and `CellCycle_Phase` will be added to the
metadata.
For `TcellState`, the columns `TcellState_Exhaustion` and
`TcellState_Activation` will be added to the metadata, one for each
program in the named `features` list (the list values are gene
vectors, not comma-separated strings).
For `DC`, a diffusion map will be calculated with
[`destiny`](https://bioconductor.org/packages/release/bioc/html/destiny.html)
(regardless of `method`), and the first 2 components will be added
as the `DC` reduction as well as the `DC_1` and `DC_2` columns to
the metadata. `dm` is a shortcut for `diffmap`/`diffusion_map`.
You can later plot the diffusion map by using
`reduction = "DC"` in `env.dimplots` in `SeuratClusterStats`.
This requires
[`SingleCellExperiment`](https://bioconductor.org/packages/release/bioc/html/SingleCellExperiment.html)
and [`destiny`](https://bioconductor.org/packages/release/bioc/html/destiny.html) R packages.
post_mutaters (type=json): The mutaters to mutate the metadata after
calculating the module scores.
The mutaters will be applied in the order specified.
This is useful when you want to create new scores based on the
calculated module scores.
""" # noqa: E501
input = "srtobj:file"
output = "rdsfile:file:{{in.srtobj | stem}}.qs"
lang = config.lang.rscript
envs = {
"defaults": {
"method": "seurat",
"features": None,
"nbin": 24,
"ctrl": 100,
"k": False,
"assay": None,
"seed": 8525,
"search": False,
},
"ncores": config.misc.ncores,
"modules": {
# "CellCycle": {"features": "cc.genes.updated.2019"},
# "Exhaustion": {"features": "HAVCR2,ENTPD1,LAYN,LAG3"},
# "Activation": {"features": "IFNG"},
# "Proliferation": {"features": "STMN1,TUBB"},
},
"post_mutaters": {},
}
script = "file://../scripts/scrna/ModuleScoreCalculator.R"
@mark(DOCS
deprecated=(
"[{proc.name}] is deprecated, "
"use [SeuratClusterStats] or [ClonalStats] instead."
)
)
class CellsDistribution(Proc):
"""Distribution of cells (i.e. in a TCR clone) from different groups
for each cluster
This generates a set of pie charts with proportion of cells in each cluster
Rows are the cells identities (i.e. TCR clones or TCR clusters), columns
are groups (i.e. clinic groups).
Examples:
```toml
[CellsDistribution.envs.mutaters]
# Add Patient1_Tumor_Expanded column with CDR3.aa that
# expands in Tumor of patient 1
Patient1_Tumor_Expanded = '''
expanded(., region, "Tumor", subset = patient == "Lung1", uniq = FALSE)
'''
[CellsDistribution.envs.cases.Patient1_Tumor_Expanded]
cells_by = "Patient1_Tumor_Expanded"
cells_orderby = "desc(CloneSize)"
group_by = "region"
group_order = [ "Tumor", "Normal" ]
```
Input:
srtobj: The seurat object in RDS format
Output:
outdir: The output directory.
The results for each case will be saved in a subdirectory.
Envs:
mutaters (type=json): The mutaters to mutate the metadata
Keys are the names of the mutaters and values are the R expressions
passed by `dplyr::mutate()` to mutate the metadata.
cluster_orderby: The order of the clusters to show on the plot.
An expression passed to `dplyr::summarise()` on the grouped data frame (by `seurat_clusters`).
The summary stat will be passed to `dplyr::arrange()` to order the clusters. It's applied on the whole meta.data before grouping and subsetting.
For example, you can order the clusters by the activation score of
the cluster: `desc(mean(ActivationScore, na.rm = TRUE))`, suppose you have a column
`ActivationScore` in the metadata.
group_by: The column name in metadata to group the cells for the columns of the plot.
group_order (list): The order of the groups (columns) to show on the plot
cells_by: The column name in metadata to group the cells for the rows of the plot.
If your cell groups have overlapping cells, you can also use multiple columns, separated by comma (`,`).
These columns will be concatenated to form the cell groups. For the overlapping cells, they will be
counted multiple times for different groups. So make sure the cell group names in different columns
are unique.
cells_order (list): The order of the cells (rows) to show on the plot
cells_orderby: An expression passed to `dplyr::arrange()` to order the cells (rows) of the plot.
Only works when `cells-order` is not specified.
The data frame passed to `dplyr::arrange()` is grouped by `cells_by` before ordering.
You can have multiple expressions separated by semicolon (`;`). The expessions will be parsed by `rlang::parse_exprs()`.
4 extra columns were added to the metadata for ordering the rows in the plot:
* `CloneSize`: The size (number of cells) of clones (identified by `cells_by`)
* `CloneGroupSize`: The clone size in each group (identified by `group_by`)
* `CloneClusterSize`: The clone size in each cluster (identified by `seurat_clusters`)
* `CloneGroupClusterSize`: The clone size in each group and cluster (identified by `group_by` and `seurat_clusters`)
cells_n (type=int): The max number of groups to show for each cell group identity (row).
Ignored if `cells_order` is specified.
subset: An expression to subset the cells, will be passed to `dplyr::filter()` on metadata.
This will be applied prior to `each`.
descr: The description of the case, will be shown in the report.
hm_devpars (ns): The device parameters for the heatmaps.
- res (type=int): The resolution of the heatmaps.
- height (type=int): The height of the heatmaps.
- width (type=int): The width of the heatmaps.
devpars (ns): The device parameters for the plots of pie charts.
- res (type=int): The resolution of the plots
- height (type=int): The height of the plots
- width (type=int): The width of the plots
each: The column name in metadata to separate the cells into different plots.
prefix_each (flag): Whether to prefix the `each` column name to the
value as the case/section name.
section: The section to show in the report. This allows different cases to be put in the same section in report.
Only works when `each` is not specified.
overlap (list): Plot the overlap of cell groups (values of `cells_by`) in different cases
under the same section.
The section must have at least 2 cases, each case should have a single `cells_by` column.
cases (type=json;order=99): If you have multiple cases, you can specify them here.
Keys are the names of the cases and values are the options above except `mutaters`.
If some options are not specified, the options in `envs` will be used.
If no cases are specified, a default case will be used with case name `DEFAULT`.
Requires:
r-seurat:
- check: {{proc.lang}} -e "library(Seurat)"
r-dplyr:
- check: {{proc.lang}} -e "library(dplyr)"
r-tidyr:
- check: {{proc.lang}} -e "library(tidyr)"
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}.cells_distribution"
lang = config.lang.rscript
envs = {
"mutaters": {},
"cluster_orderby": None,
"group_by": None,
"group_order": [],
"cells_by": None,
"cells_order": [],
"cells_orderby": None,
"cells_n": 10,
"subset": None,
"descr": None,
"devpars": {},
"hm_devpars": {},
"each": None,
"prefix_each": True,
"section": "DEFAULT",
"overlap": [],
"cases": {},
}
script = "file://../scripts/scrna/CellsDistribution.R"
plugin_opts = {
"report": "file://../reports/scrna/CellsDistribution.svelte",
"report_paging": 8,
}
class SeuratMetadataMutater(Proc):DOCS
"""Mutate the metadata of the seurat object
Input:
srtobj: The seurat object loaded by SeuratPreparing
metafile: Additional metadata
A tab-delimited file with columns as meta columns and rows as
cells.
Output:
outfile: The seurat object with the additional metadata
Envs:
ncores (type=int): The number of cores to use for reading and writing the seurat object.
mutaters (type=json): The mutaters to mutate the metadata.
The key-value pairs will be passed the `dplyr::mutate()` to mutate the metadata.
See <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
subset: An expression to subset the cells, will be passed to `dplyr::filter()`.
This will be applied after mutating the metadata.
Requires:
r-seurat:
- check: {{proc.lang}} <(echo "library(Seurat)")
r-tibble:
- check: {{proc.lang}} <(echo "library(tibble)")
r-dplyr:
- check: {{proc.lang}} <(echo "library(dplyr)")
""" # noqa: E501
input = "srtobj:file, metafile:file"
output = "outfile:file:{{in.srtobj | stem}}.qs"
lang = config.lang.rscript
envs = {"mutaters": {}, "subset": None, "ncores": config.misc.ncores}
script = "file://../scripts/scrna/SeuratMetadataMutater.R"
@mark(deprecated="[{proc.name}] is deprecated, use [SeuratClusterStats] instead.")DOCS
class DimPlots(Proc):
"""Seurat - Dimensional reduction plots
Input:
srtobj: The seruat object in RDS format
configfile: A toml configuration file with "cases"
If this is given, `envs.cases` will be overriden
name: The name of the job, used in report
Output:
outdir: The output directory
Envs:
cases: The cases for the dim plots
Keys are the names and values are the arguments to
`Seurat::Dimplots`
"""
input = "srtobj:file, configfile:file, name:var"
output = "outdir:dir:{{in.srtobj | stem}}.dimplots"
lang = config.lang.rscript
script = "file://../scripts/scrna/DimPlots.R"
envs = {"cases": {"Ident": {"group.by": "ident"}}}
plugin_opts = {
"report": "file://../reports/scrna/DimPlots.svelte",
"report_toc": False,
}
class MarkersFinder(Proc):DOCS
"""Find markers between different groups of cells
When only `group_by` is specified as identity column in
`envs.cases`, the markers will be found for all the clusters.
You can also find the differentially expressed genes between
any two groups of cells by setting `group_by` to a different
column name in metadata. Follow `envs.cases` for more details.
Input:
srtobj: The seurat object loaded by `SeuratPreparing`
If you have your `Seurat` object prepared by yourself, you can also
use it here, but you should make sure that the object has been processed
by `PrepSCTFindMarkers` if data is not normalized using `SCTransform`.
Output:
outdir: The output directory for the markers and plots
Envs:
ncores (type=int): Number of cores to use for parallel computing for some `Seurat` procedures.
* Used in `future::plan(strategy = "multicore", workers = <ncores>)` to parallelize some Seurat procedures.
* See also: <https://satijalab.org/seurat/articles/future_vignette.html>
mutaters (type=json): The mutaters to mutate the metadata.
You can also use the clone selectors to select the TCR clones/clusters.
See <https://pwwang.github.io/scplotter/reference/clone_selectors.html>.
You can also use key `<newcol>:ident` to set the `<newcol>` as the default
ident for the stats.
See also <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
group_by: The column name in metadata to group the cells.
If only `group_by` is specified, and `ident_1` and `ident_2` are
not specified, markers will be found for all groups in this column
in the manner of "group vs rest" comparison.
`NA` group will be ignored.
If `None`, `Seurat::Idents(srtobj)` will be used, which is usually
`"seurat_clusters"` after unsupervised clustering.
ident_1: The first group of cells to compare
When this is empty, the comparisons will be expanded to each group v.s. the rest of the cells in `group_by`.
ident_2: The second group of cells to compare
If not provided, the rest of the cells are used for `ident_2`.
each: The column name in metadata to separate the cells into different
cases.
When this is specified, the case will be expanded for each value of
the column in metadata. For example, when you have `envs.cases."Cluster Markers".each = "Sample"`,
then the case will be expanded as `envs.cases."Cluster Markers - Sample1"`, `envs.cases."Cluster Markers - Sample2"`, etc.
You can specify `allmarker_plots` and `overlaps` to plot the markers for all cases in the same plot and plot the overlaps of the markers
between different cases by values in this column.
`"ident"` can be used as an alias for the default identity column (e.g. `"seurat_clusters"`).
dbs (list): The dbs to do enrichment analysis for significant markers.
You can use built-in dbs in `enrichit`, or provide your own gmt files.
See also <https://pwwang.github.io/enrichit/reference/FetchGMT.html>.
The built-in dbs include:
* "BioCarta" or "BioCarta_2016"
* "GO_Biological_Process" or "GO_Biological_Process_2025"
* "GO_Cellular_Component" or "GO_Cellular_Component_2025"
* "GO_Molecular_Function" or "GO_Molecular_Function_2025"
* "KEGG", "KEGG_Human", "KEGG_2021", or "KEGG_2021_Human"
* "Hallmark", "MSigDB_Hallmark", or "MSigDB_Hallmark_2020"
* "Reactome", "Reactome_Pathways", or "Reactome_Pathways_2024"
* "WikiPathways", "WikiPathways_2024", "WikiPathways_Human", or "WikiPathways_2024_Human"
You can also fetch more dbs from <https://maayanlab.cloud/Enrichr/#libraries>.
sigmarkers: An expression passed to `dplyr::filter()` to filter the
significant markers for enrichment analysis.
Available variables are `p_val`, `avg_log2FC`, `pct.1`, `pct.2` and
`p_val_adj`. For example, `"p_val_adj < 0.05 & abs(avg_log2FC) > 1"`
to select markers with adjusted p-value < 0.05 and absolute log2
fold change > 1.
enrich_style (choice): The style of the enrichment analysis.
The enrichment analysis will be done by `EnrichIt()` from [`enrichit`](https://pwwang.github.io/enrichit/).
Two styles are available:
- enrichr: `enrichr` style enrichment analysis (fisher's exact test will be used).
- clusterprofiler: `clusterProfiler` style enrichment analysis (hypergeometric test will be used).
- clusterProfiler: alias for `clusterprofiler`
assay: The assay to use.
subset: An expression to subset the cells for each case.
error (flag): Stop the job if errors happen.
Helpful when no/not enough markers are found or no pathways are enriched.
If `False`, empty results will be returned.
cache (type=auto): Where to cache the results.
If `True`, cache to `outdir` of the job. If `False`, don't cache.
Otherwise, specify the directory to cache to.
rest (ns): Rest arguments for `Seurat::FindMarkers()`.
Use `-` to replace `.` in the argument name. For example,
use `min-pct` instead of `min.pct`.
- <more>: See <https://satijalab.org/seurat/reference/findmarkers>
allmarker_plots_defaults (ns): Default options for the plots for all markers when `ident_1` is not specified.
To reproduce what [`Seurat::DoHeatmap()`](https://satijalab.org/seurat/reference/doheatmap) does, you can
`each = ':seurat_cluster'` to select the top N (default: 20) markers all together (instead of selecting top N markers for each cluster),
and set `plot_type = "heatmap"` and `cell_type = "bars"` to plot the heatmap of the top N markers.
You may also want to use `order_by` to order the markers and `select = N` to select the top N markers.
- plot_type: The type of the plot.
See <https://pwwang.github.io/biopipen.utils.R/reference/VizDEGs.html>.
Available types are `violin`, `box`, `bar`, `ridge`, `dim`, `heatmap` and `dot`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: Other arguments passed to [`biopipen.utils::VizDEGs()`](https://pwwang.github.io/biopipen.utils.R/reference/VizDEGs.html).
allmarker_plots (type=json): All marker plot cases.
The keys are the names of the cases and the values are the dicts inherited from `allmarker_plots_defaults`.
allenrich_plots_defaults (ns): Default options for the plots to generate for the enrichment analysis.
- plot_type: The type of the plot.
- descr: A description of the plot to be shown above the plot image.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- db (type=str): The database(s) to apply for the plot. If not specified, all databases will be applied.
`dbs` can also be used as an alias for this parameter.
- <more>: See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.html>.
allenrich_plots (type=json): Cases of the plots to generate for the enrichment analysis.
The keys are the names of the cases and the values are the dicts inherited from `allenrich_plots_defaults`.
The cases under `envs.cases` can inherit this options.
marker_plots_defaults (ns): Default options for the plots to generate for the markers.
- plot_type: The type of the plot.
See <https://pwwang.github.io/biopipen.utils.R/reference/VizDEGs.html>.
Available types are `violin`, `box`, `bar`, `ridge`, `dim`, `heatmap` and `dot`.
There are two additional types available - `volcano_pct` and `volcano_log2fc`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: Other arguments passed to [`biopipen.utils::VizDEGs()`](https://pwwang.github.io/biopipen.utils.R/reference/VizDEGs.html).
If `plot_type` is `volcano_pct` or `volcano_log2fc`, they will be passed to
[`scplotter::VolcanoPlot()`](https://pwwang.github.io/plotthis/reference/VolcanoPlot.html).
marker_plots (type=json): Cases of the plots to generate for the markers.
Plot cases. The keys are the names of the cases and the values are the dicts inherited from `marker_plots_defaults`.
The cases under `envs.cases` can inherit this options.
enrich_plots_defaults (ns): Default options for the plots to generate for the enrichment analysis.
- plot_type: The type of the plot.
See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.html>.
Available types are `bar`, `dot`, `lollipop`, `network`, `enrichmap` and `wordcloud`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.html>.
enrich_plots (type=json): Cases of the plots to generate for the enrichment analysis.
The keys are the names of the cases and the values are the dicts inherited from `enrich_plots_defaults`.
The cases under `envs.cases` can inherit this options.
overlaps_defaults (ns): Default options for investigating the overlapping of significant markers between different cases or comparisons.
This means either `ident_1` should be empty, so that they can be expanded to multiple comparisons.
- sigmarkers: The expression to filter the significant markers for each case.
If not provided, `envs.sigmarkers` will be used.
- descr: A description of the plot to be shown above the plot image.
- plot_type (choice): The type of the plot to generate for the overlaps.
- venn: Use `plotthis::VennDiagram()`.
- upset: Use `plotthis::UpsetPlot()`.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: More arguments pased to `plotthis::VennDiagram()`
(<https://pwwang.github.io/plotthis/reference/venndiagram1.html>)
or `plotthis::UpsetPlot()`
(<https://pwwang.github.io/plotthis/reference/upsetplot1.html>)
overlaps (type=json): Cases for investigating the overlapping of significant markers between different cases or comparisons.
The keys are the names of the cases and the values are the dicts inherited from `overlaps_defaults`.
There are two situations that we can perform overlaps:
1. If `ident_1` is not specified, the overlaps can be performed between different comparisons.
2. If `each` is specified, the overlaps can be performed between different cases, where in each case, `ident_1` must be specified.
cases (type=json): If you have multiple cases for marker discovery, you can specify them
here. The keys are the names of the cases and the values are the above options. If some options are
not specified, the default values specified above (under `envs`) will be used.
If no cases are specified, the default case will be added with the default values under `envs` with the name `Marker Discovery`.
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem0}}.markers"
lang = config.lang.rscript
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"group_by": None,
"ident_1": None,
"ident_2": None,
"each": None,
"dbs": ["KEGG_2021_Human", "MSigDB_Hallmark_2020"],
"sigmarkers": "p_val_adj < 0.05",
"enrich_style": "enrichr",
"assay": None,
"error": True,
"subset": None,
"cache": config.path.tmpdir,
"rest": {},
"allmarker_plots_defaults": {
"plot_type": None,
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"allmarker_plots": {},
"allenrich_plots_defaults": {
"plot_type": "heatmap",
"devpars": {"res": 100},
},
"allenrich_plots": {},
"marker_plots_defaults": {
"plot_type": None,
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"marker_plots": {
"Volcano Plot (diff_pct)": {"plot_type": "volcano_pct"},
"Volcano Plot (log2FC)": {"plot_type": "volcano_log2fc"},
"Dot Plot": {"plot_type": "dot"},
},
"enrich_plots_defaults": {
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"enrich_plots": {
"Bar Plot": {"plot_type": "bar", "ncol": 1, "top_term": 10},
},
"overlaps_defaults": {
"sigmarkers": None,
"plot_type": "venn",
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"overlaps": {},
"cases": {},
}
order = 5
script = "file://../scripts/scrna/MarkersFinder.R"
plugin_opts = {
"report": "file://../reports/scrna/MarkersFinder.svelte",
"report_paging": 8,
}
class TopExpressingGenes(Proc):DOCS
"""Find the top expressing genes in each cluster
Input:
srtobj: The seurat object in RDS or qs/qs2 format
Output:
outdir: The output directory for the tables and plots
Envs:
ncores (type=int): Number of cores to use for reading and writing the seurat object.
mutaters (type=json): The mutaters to mutate the metadata.
You can also use the clone selectors to select the TCR clones/clusters.
See <https://pwwang.github.io/scplotter/reference/clone_selectors.html>.
You can also use key `<newcol>:ident` to set the `<newcol>` as the default
ident for the stats.
See also <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
ident: The group of cells to find the top expressing genes.
The cells will be selected by the `group_by` column with this
`ident` value in metadata.
If not provided, the top expressing genes will be found for all
groups of cells in the `group_by` column.
group_by: The column name in metadata to group the cells.
each: The column name in metadata to separate the cells into different
cases.
`"ident"` can be used as an alias for the default identity column (e.g. `"seurat_clusters"`).
dbs (list): The dbs to do enrichment analysis for significant markers.
You can use built-in dbs in `enrichit`, or provide your own gmt files.
See also <https://pwwang.github.io/enrichit/reference/FetchGMT.html>.
The built-in dbs include:
* "BioCarta" or "BioCarta_2016"
* "GO_Biological_Process" or "GO_Biological_Process_2025"
* "GO_Cellular_Component" or "GO_Cellular_Component_2025"
* "GO_Molecular_Function" or "GO_Molecular_Function_2025"
* "KEGG", "KEGG_Human", "KEGG_2021", or "KEGG_2021_Human"
* "Hallmark", "MSigDB_Hallmark", or "MSigDB_Hallmark_2020"
* "Reactome", "Reactome_Pathways", or "Reactome_Pathways_2024"
* "WikiPathways", "WikiPathways_2024", "WikiPathways_Human", or "WikiPathways_2024_Human"
You can also fetch more dbs from <https://maayanlab.cloud/Enrichr/#libraries>.
n (type=int): The number of top expressing genes to find.
enrich_style (choice): The style of the enrichment analysis.
The enrichment analysis will be done by `EnrichIt()` from [`enrichit`](https://pwwang.github.io/enrichit/).
Two styles are available:
- enrichr: `enrichr` style enrichment analysis (fisher's exact test will be used).
- clusterprofiler: `clusterProfiler` style enrichment analysis (hypergeometric test will be used).
- clusterProfiler: alias for `clusterprofiler`
enrich_plots_defaults (ns): Default options for the plots to generate for the enrichment analysis.
- plot_type: The type of the plot.
See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.html>.
Available types are `bar`, `dot`, `lollipop`, `network`, `enrichmap` and `wordcloud`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.htmll>.
enrich_plots (type=json): Cases of the plots to generate for the enrichment analysis.
The keys are the names of the cases and the values are the dicts inherited from `enrich_plots_defaults`.
The cases under `envs.cases` can inherit this options.
subset: An expression to subset the cells for each case.
error (flag): Stop the job if errors happen.
Helpful when no/not enough markers are found or no pathways are enriched.
If `False`, empty results will be returned.
cases (type=json): If you have multiple cases, you can specify them
here. The keys are the names of the cases and the values are the
above options except `mutaters`. If some options are
not specified, the default values specified above will be used.
If no cases are specified, the default case will be added with
the default values under `envs` with the name `Top Expressing Genes`.
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}.top_expressing_genes"
lang = config.lang.rscript
script = "file://../scripts/scrna/TopExpressingGenes.R"
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"ident": None,
"group_by": None,
"each": None,
"dbs": ["KEGG_2021_Human", "MSigDB_Hallmark_2020"],
"n": 250,
"subset": None,
"error": False,
"enrich_style": "enrichr",
"enrich_plots_defaults": {
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"enrich_plots": {
"Bar Plot": {"plot_type": "bar", "ncol": 1, "top_term": 10},
},
"cases": {},
}
plugin_opts = {
"report": "file://../reports/common.svelte",
"report_paging": 8,
}
class ExprImputation(Proc):DOCS
"""This process imputes the dropout values in scRNA-seq data.
It takes the Seurat object as input and outputs the Seurat object with
imputed expression data.
Reference:
- [Linderman, George C., Jun Zhao, and Yuval Kluger. "Zero-preserving imputation of scRNA-seq data using low-rank approximation." BioRxiv (2018): 397588.](https://www.nature.com/articles/s41467-021-27729-z)
- [Li, Wei Vivian, and Jingyi Jessica Li. "An accurate and robust imputation method scImpute for single-cell RNA-seq data." Nature communications 9.1 (2018): 997.](https://www.nature.com/articles/s41467-018-03405-7)
- [Dijk, David van, et al. "MAGIC: A diffusion-based imputation method reveals gene-gene interactions in single-cell RNA-sequencing data." BioRxiv (2017): 111591.](https://www.cell.com/cell/abstract/S0092-8674(18)30724-4)
Input:
infile: The input file in RDS/qs format of Seurat object
Output:
outfile: The output file in RDS format of Seurat object
Note that with rmagic and alra, the original default assay will be
renamed to `RAW` and the imputed RNA assay will be
renamed to `RNA` and set as default assay.
Envs:
tool (choice): Either alra, scimpute or rmagic
- alra: Use RunALRA() from Seurat
- scimpute: Use scImpute() from scimpute
- rmagic: Use magic() from Rmagic
ncores (type=int): Number of cores to use for reading and writing,
the Seurat object, and for scimpute and rmagic if applicable.
scimpute_args (ns): The arguments for scimpute
- drop_thre (type=float): The dropout threshold
- kcluster (type=int): Number of clusters to use
- refgene: The reference gene file
rmagic_args (ns): The arguments for rmagic
- python: The python path where magic-impute is installed.
- threshold (type=float): The threshold for magic imputation.
Only the genes with dropout rates greater than this threshold (No. of
cells with non-zero expression / total number of cells) will be imputed.
alra_args (type=json): The arguments for `RunALRA()`
Requires:
r-scimpute:
- if: {{proc.envs.tool == "scimpute"}}
- check: {{proc.lang}} <(echo "library(scImpute)")
r-rmagic:
- if: {{proc.envs.tool == "rmagic"}}
- check: |
{{proc.lang}} <(\
echo "\
tryCatch(\
{ setwd(dirname(Sys.getenv('CONDA_PREFIX'))) }, \
error = function(e) NULL \
); \
library(Rmagic)\
"\
)
magic-impute:
- if: {{proc.envs.tool == "rmagic"}}
- check: {{proc.envs.rmagic_args.python}} -c "import magic")
r-dplyr:
- if: {{proc.envs.tool == "scimpute"}}
- check: {{proc.lang}} <(echo "library(dplyr)")
r-seurat:
- check: {{proc.lang}} <(echo "library(Seurat)")
r-seuratwrappers:
- if: {{proc.envs.tool == "alra"}}
- check: {{proc.lang}} <(echo "library(SeuratWrappers)")
""" # noqa: E501
input = "infile:file"
output = "outfile:file:{{in.infile | stem}}.imputed.qs"
lang = config.lang.rscript
envs = {
"tool": "alra",
"ncores": config.misc.ncores,
"rmagic_args": {"python": config.exe.magic_python, "threshold": 0.5},
"scimpute_args": {
"drop_thre": 0.5,
"kcluster": None,
"refgene": config.ref.refgene,
},
"alra_args": {},
}
script = "file://../scripts/scrna/ExprImputation.R"
class SCImpute(Proc):DOCS
"""Impute the dropout values in scRNA-seq data.
Deprecated. Use `ExprImputation` instead.
Input:
infile: The input file for imputation
Either a SeuratObject or a matrix of count/TPM
groupfile: The file to subset the matrix or label the cells
Could be an output from ImmunarchFilter
Output:
outfile: The output matrix
Envs:
infmt: The input format.
Either `seurat` or `matrix`
"""
input = "infile:file, groupfile:file"
output = [
"outfile:file:{{in.infile | stem | replace: '.seurat', ''}}." "{{envs.outfmt}}"
]
lang = config.lang.rscript
envs = {
"infmt": "seurat", # or matrix
"outfmt": "txt", # or csv, rds
"type": "count", # or TPM
"drop_thre": 0.5,
"kcluster": None, # or Kcluster
"ncores": config.misc.ncores,
"refgene": config.ref.refgene,
}
script = "file://../scripts/scrna/SCImpute.R"
class SeuratFilter(Proc):DOCS
"""Filtering cells from a seurat object
Input:
srtobj: The seurat object in RDS
filters: The filters to apply. Could be a file or string in TOML, or
a python dictionary, with following keys:
- mutaters: Create new columns in the metadata
- filter: A R expression that will pass to
`subset(sobj, subset = ...)` to filter the cells
Output:
outfile: The filtered seurat object in RDS
Envs:
invert: Invert the selection?
Requires:
r-seurat:
- check: {{proc.lang}} <(echo "library('Seurat')")
r-dplyr:
- check: {{proc.lang}} <(echo "library('dplyr')")
"""
input = "srtobj:file, filters:var"
output = "outfile:file:{{in.srtobj | stem}}.filtered.RDS"
lang = config.lang.rscript
envs = {"invert": False}
script = "file://../scripts/scrna/SeuratFilter.R"
class SeuratSubset(Proc):DOCS
"""Subset a seurat object into multiple seruat objects
Input:
srtobj: The seurat object in RDS
subsets: The subsettings to apply. Could be a file or string in TOML, or
a python dictionary, with following keys:
- <name>: Name of the case
mutaters: Create new columns in the metadata
subset: A R expression that will pass to
`subset(sobj, subset = ...)`
groupby: The column to group by, each value will be a case
If groupby is given, subset will be ignored, each value
of the groupby column will be a case
Output:
outdir: The output directory with the subset seurat objects
Envs:
ignore_nas: Ignore NA values?
Requires:
r-seurat:
- check: {{proc.lang}} <(echo "library('Seurat')")
r-dplyr:
- check: {{proc.lang}} <(echo "library('dplyr')")
"""
input = "srtobj:file, subsets:var"
output = "outdir:dir:{{in.srtobj | stem}}.subsets"
envs = {"ignore_nas": True}
lang = config.lang.rscript
script = "file://../scripts/scrna/SeuratSubset.R"
class SeuratSplit(Proc):DOCS
"""Split a seurat object into multiple seruat objects
Input:
srtobj: The seurat object in RDS
by: The metadata column to split by
Output:
outdir: The output directory with the subset seurat objects
Envs:
by: The metadata column to split by
Ignored if `by` is given in the input
recell: Rename the cell ids using the `by` column
A string of R function taking the original cell ids and `by`
"""
input = "srtobj:file, by:var"
output = "outdir:dir:{{in.srtobj | stem}}.subsets"
envs = {
"by": None,
"recell": None, # "function(cellid, by) {}",
}
lang = config.lang.rscript
script = "file://../scripts/scrna/SeuratSplit.R"
class Subset10X(Proc):DOCS
"""Subset 10X data, mostly used for testing
Requires r-matrix to load matrix.mtx.gz
Input:
indir: The input directory
Output:
outdir: The output directory
Envs:
seed: The seed for random number generator
nfeats: The number of features to keep.
If <=1 then it will be the percentage of features to keep
ncells: The number of cells to keep.
If <=1 then it will be the percentage of cells to keep
feats_to_keep: The features/genes to keep.
The final features list will be `feats_to_keep` + `nfeats`
"""
input = "indir:dir"
output = "outdir:dir:{{in.indir | stem}}"
envs = {
"seed": 8525,
"nfeats": 0.1,
"ncells": 0.1,
"feats_to_keep": [],
}
lang = config.lang.rscript
script = "file://../scripts/scrna/Subset10X.R"
class SeuratTo10X(Proc):DOCS
"""Write a Seurat object to 10X format
using `write10xCounts` from `DropletUtils`
Input:
srtobj: The seurat object in RDS
Output:
outdir: The output directory.
When `envs.split_by` is specified, the subdirectories will be
created for each distinct value of the column.
Otherwise, the matrices will be written to the output directory.
Envs:
version: The version of 10X format
"""
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}"
envs = {"version": "3", "split_by": None}
lang = config.lang.rscript
script = "file://../scripts/scrna/SeuratTo10X.R"
class ScFGSEA(Proc):DOCS
"""Gene set enrichment analysis for cells in different groups using `fgsea`
This process allows us to do Gene Set Enrichment Analysis (GSEA) on the expression data,
but based on variaties of grouping, including the from the meta data and the
scTCR-seq data as well.
The GSEA is done using the
[fgsea](https://bioconductor.org/packages/release/bioc/html/fgsea.html) package,
which allows to quickly and accurately calculate arbitrarily low GSEA P-values
for a collection of gene sets.
The fgsea package is based on the fast algorithm for preranked GSEA described in
[Subramanian et al. 2005](https://www.pnas.org/content/102/43/15545).
For each case, the process will generate a table with the enrichment scores for
each gene set, and GSEA plots for the top gene sets.
Examples:
### The summary and GSEA plots
{: width="80%"}
{: width="80%"}
### Summary plot for all subsets or idents
If you use `each` to separate the cells into different subsets, this is useful to
make a summary plot for all subsets. Or if you don't specify `ident_1`, the summary plot for all idents in `group_by` will be generated.
```toml
[ScFGSEA.envs]
group_by = "Diagnosis"
ident_1 = "Colitis"
ident_2 = "Control"
each = "seurat_clusters"
[ScFGSEA.envs.alleach_plots.Heatmap]
plot_type = "heatmap"
group_by = "Diagnosis"
```
{: width="80%"}
Input:
srtobj: The seurat object in RDS format
Output:
outdir: The output directory for the results and plots
Envs:
ncores (type=int): Number of cores for parallelization
Passed to `nproc` of `fgseaMultilevel()`.
cache (type=auto): Where to cache the results.
If `True`, cache to `outdir` of the job. If `False`, don't cache.
Otherwise, specify the directory to cache to.
mutaters (type=json): The mutaters to mutate the metadata.
The key-value pairs will be passed the `dplyr::mutate()` to mutate the metadata.
You can also use the clone selectors to select the TCR clones/clusters.
See <https://pwwang.github.io/scplotter/reference/clone_selectors.html>.
You can also use key `<newcol>:ident` to set the `<newcol>` as the default
ident for the stats.
See also <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
group_by: The column name in metadata to group the cells.
ident_1: The first group of cells to compare
ident_2: The second group of cells to compare, if not provided, the rest of the cells that are not `NA`s in `group_by` column are used for `ident_2`.
assay: The assay to use. If not provided, the default assay will be used.
each: The column name in metadata to separate the cells into different subsets to do the analysis.
`"ident"` can be used as an alias for the default identity column (e.g. `"seurat_clusters"`).
subset: An expression to subset the cells.
error (flag): Stop the job if errors happen.
Helpful when no/not enough markers are found or no pathways are enriched.
If `False`, empty results will be returned.
gmtfile: The pathways in GMT format, with the gene names/ids in the same format as the seurat object.
You can use built-in dbs in `enrichit`, or provide your own gmt files.
See also <https://pwwang.github.io/enrichit/reference/FetchGMT.html>.
The built-in dbs include:
* "BioCarta" or "BioCarta_2016"
* "GO_Biological_Process" or "GO_Biological_Process_2025"
* "GO_Cellular_Component" or "GO_Cellular_Component_2025"
* "GO_Molecular_Function" or "GO_Molecular_Function_2025"
* "KEGG", "KEGG_Human", "KEGG_2021", or "KEGG_2021_Human"
* "Hallmark", "MSigDB_Hallmark", or "MSigDB_Hallmark_2020"
* "Reactome", "Reactome_Pathways", or "Reactome_Pathways_2024"
* "WikiPathways", "WikiPathways_2024", "WikiPathways_Human", or "WikiPathways_2024_Human"
You can also fetch more dbs from <https://maayanlab.cloud/Enrichr/#libraries>.
method (choice): The method to do the preranking.
- signal_to_noise: Signal to noise.
The larger the differences of the means (scaled by the standard deviations);
that is, the more distinct the gene expression is in each phenotype and the more the gene
acts as a "class marker".
- s2n: Alias of signal_to_noise.
- abs_signal_to_noise: The absolute value of signal_to_noise.
- abs_s2n: Alias of abs_signal_to_noise.
- t_test: T test.
Uses the difference of means scaled by the standard deviation and number of samples.
- ratio_of_classes: Also referred to as fold change.
Uses the ratio of class means to calculate fold change for natural scale data.
- diff_of_classes: Difference of class means.
Uses the difference of class means to calculate fold change for nature scale data
- log2_ratio_of_classes: Log2 ratio of class means.
Uses the log2 ratio of class means to calculate fold change for natural scale data.
This is the recommended statistic for calculating fold change for log scale data.
top (type=auto): Do gsea table and enrich plot for top N pathways.
If it is < 1, will apply it to `padj`, selecting pathways with `padj` < `top`.
eps (type=float): This parameter sets the boundary for calculating the p value.
See <https://rdrr.io/bioc/fgsea/man/fgseaMultilevel.html>
alleach_plots_defaults (ns): Default options for the plots to generate for all pathways.
- plot_type: The type of the plot, currently either dot or heatmap (default)
- descr: A description of the plot to be shown above the plot image.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: See <https://pwwang.github.io/biopipen.utils.R/reference/VizGSEA.html>.
alleach_plots (type=json): Cases of the plots to generate for all pathways.
The keys are the names of the cases and the values are the dicts inherited from `alleach_plots_defaults`.
minsize (type=int): Minimal size of a gene set to test. All pathways below the threshold are excluded.
maxsize (type=int): Maximal size of a gene set to test. All pathways above the threshold are excluded.
rest (type=json;order=98): Rest arguments for [`fgsea()`](https://rdrr.io/bioc/fgsea/man/fgsea.html)
See also <https://rdrr.io/bioc/fgsea/man/fgseaMultilevel.html>
cases (type=json;order=99): If you have multiple cases, you can specify them here.
The keys are the names of the cases and the values are the above options except `mutaters`.
If some options are not specified, the default values specified above will be used.
If no cases are specified, the default case will be added with the name `GSEA`.
Requires:
bioconductor-fgsea:
- check: {{proc.lang}} -e "library(fgsea)"
r-seurat:
- check: {{proc.lang}} -e "library(seurat)"
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{(in.casefile or in.srtobj) | stem0}}.fgsea"
lang = config.lang.rscript
envs = {
"mutaters": {},
"ncores": config.misc.ncores,
"cache": config.path.tmpdir,
"assay": None,
"group_by": None,
"ident_1": None,
"ident_2": None,
"each": None,
"subset": None,
"error": False,
"gmtfile": "KEGG_2021_Human",
"method": "s2n",
"top": 20,
"minsize": 10,
"maxsize": 100,
"eps": 0,
"alleach_plots_defaults": {
"plot_type": "heatmap",
"devpars": {"res": 100},
},
"alleach_plots": {},
"rest": {},
"cases": {},
}
script = "file://../scripts/scrna/ScFGSEA.R"
plugin_opts = {
"report": "file://../reports/common.svelte",
"report_paging": 8,
}
class CellTypeAnnotation(Proc):DOCS
"""Annotate the cell clusters. Currently, the following ways are supported:
1. Pass the cell type annotation directly (at cluster-level or cell-level)
2. Use [`ScType`](https://github.com/IanevskiAleksandr/sc-type) (cluster-level, marker-based)
3. Use [`scCATCH`](https://github.com/ZJUFanLab/scCATCH) (cluster-level, marker-based)
4. Use [`hitype`](https://github.com/pwwang/hitype) (cell-level or cluster-level, marker-based)
5. Use [`celltypist`](https://github.com/Teichlab/celltypist) (cell-level, model-based)
6. Use [`scSorter`](https://pmc.ncbi.nlm.nih.gov/articles/PMC7898451/) (cell-level, marker-based)
7. Use [`SCINA`](https://github.com/jcao89757/SCINA) (cell-level, marker-based)
8. Use [`SingleR`](https://github.com/dviraran/SingleR) (cluster-level, model-based)
9. Use [`scHDeepInsight`](https://github.com/shangruJia/scHDeepInsight) (cell-level, model-based)
10. Use [`LLMCelltype`](https://github.com/pwwang/LLMCelltype) (cluster-level, LLM-based, a fork of [`gptcelltype`](https://github.com/Winnie09/GPTCelltype))
11. Use [`cellassign`](https://github.com/Irrationone/cellassign) (cell-level, marker-based)
12. Use [`scBERT`](https://github.com/TencentAILabHealthcare/scBERT) (cell-level, model-based)
13. Use [`CelliD`](https://github.com/RausellLab/CelliD) (cell-level, marker-based)
14. Use [`scAgentType`](https://github.com/sathyasjali/scAgentType) (cluster-level, LLM-based, agentic)
15. Use [`garnett`](https://cole-trapnell-lab.github.io/garnett/) (cell-level, model-based)
16. Use [`UCell`](https://github.com/carmonalab/UCell) (cell-level, marker-based scoring)
17. Use [`AUCell`](https://github.com/aertslab/AUCell) (cell-level, marker-based scoring)
18. Use [`GSVA`](https://bioconductor.org/packages/release/bioc/html/GSVA.html) (cell-level, marker-based scoring)
19. Use [`singscore`](https://bioconductor.org/packages/release/bioc/html/singscore.html) (cell-level, marker-based scoring)
20. Use [`scmap`](https://bioconductor.org/packages/release/bioc/html/scmap.html) (cell-level or cluster-level, reference-based)
21. Use [`CHETAH`](https://bioconductor.org/packages/release/bioc/html/CHETAH.html) (cell-level, reference-based)
22. Use [`scClassify`](https://bioconductor.org/packages/release/bioc/html/scClassify.html) (cell-level, reference-based)
23. Use [`scPred`](https://github.com/powellgenomicslab/scPred) (cell-level, reference-based)
24. Use [`Azimuth`](https://github.com/satijalab/azimuth) (cluster-level, reference-based)
25. Use [`SCSA`](https://github.com/bioinfo-ibms-pumc/SCSA) (cluster-level, marker-based, python)
26. Use [`MACA`](https://github.com/ImXman/MACA) (cell-level, marker-based, python, a modernized fork should be used, https://github.com/pwwang/MACA/tree/modernize)
27. Use [`scMapNet`](https://github.com/Yuz7/scMapNet) (cell-level, marker-based, python)
28. Use [`mLLMCelltype`](https://github.com/cafferychen777/mLLMCelltype) (cluster-level, LLM-based)
29. Use [`LICT`](https://github.com/Glowworm-cell/LICT) (cluster-level, LLM-based, a modernized fork should be used, https://github.com/pwwang/LICT/tree/modernize)
30. Use [`MapQuery`](https://satijalab.org/seurat/reference/mapquery) (cell-level, reference-based)
31. Use [`ScInfeR`](https://github.com/swainasish/ScInfeR) (cluster-level, marker-based; graph-based,
in two rounds: it calls the clusters from the marker table first, then the cells within a cluster
when two cell types score within 0.1 of each other; the marker table is either a universal one
in `envs.scinfer.db`, or one of ScInfeR's own ScInfeRDB tissues in `envs.scinfer.tissue`)
(Swain, Shekhawat & Yadav, 2025, doi:[10.1093/bib/bbaf253](https://doi.org/10.1093/bib/bbaf253);
there is a published correction, doi:[10.1093/bib/bbaf337](https://doi.org/10.1093/bib/bbaf337))
The tools can be divided into two categories:
- Cluster-level tools: annotate the clusters, each cluster being assigned one cell type.
These include `sctype`, `sccatch`, `scsorter`, `singler`, `azimuth`, `llmcelltype`,
`mllmcelltype`, `lict`, `scagenttype`, `scsa`, `scinfer`, and `direct`.
- Cell-level tools: annotate the cells, each cell being assigned one cell type.
These include `scina`, `hitype`, `cellassign`, `cellid`, `scbert`, `schdeepinsight`,
`llmcelltype`, `cell`, `garnett`, `ucell`, `aucell`, `gsva`, `singscore`, `scmap`,
`cheetah`, `scclassify`, `scpred`, `mapquery`, `maca`, and `scmapnet`.
For cell-level tools, a cluster-level annotation can also be generated by specifying
`envs.ident` (or `envs.cases.X.ident` for a specific case), where each cluster is
assigned the cell type by majority vote of the cells in the cluster.
- `hitype` can run at either level: cell-level when `envs.ident` is not set,
or cluster-level when `envs.ident` is set (hitype then scores each cluster
as a whole and assigns its cell type directly, not by majority vote).
- `celltypist` is a special cell-level tool, which can also generate cluster-level
annotations via its over-clustering mechanism (`envs.celltypist.over_clustering`
or `envs.ident`).
`scsa` annotates the clusters like `sctype` and `sccatch`, but its wrapper maps
the cluster labels back onto the cells, so a per-cell column is produced as well
(the same shape as `azimuth`: a cluster-level mapping plus per-cell labels).
The tools can also be divided by their input types:
- Marker-based tools: take marker genes for the cell types, including
`sctype`, `hitype`, `sccatch`, `scsorter`, `scina`, `cellassign`, and `cellid`.
They all accept the universal marker format (see the note below).
`scinfer` is marker-based as well: it takes the universal marker table through
`envs.scinfer.db`, or, like `sctype` and `sccatch` with their bundled databases,
a database of its own through `envs.scinfer.tissue` (ScInfeRDB: ScInfeR's own
per-tissue cell type → marker → weight set, fetched at call time).
The scorer family (`ucell`, `aucell`, `gsva`, and `singscore`) is also
marker-based: it scores the *same* marker table for every cell with equal
weights (no learned weights), and assigns each cell the cell type with the
highest score, so a `weight` column in the marker table is ignored.
- Model/reference-based tools: take a trained model or a reference object,
including `celltypist`, `scbert`, `singler`, `schdeepinsight`, `garnett`,
`scmap`, `cheetah`, `scclassify`, `scpred`, `azimuth`, and `mapquery`.
`mapquery` consumes a labeled reference and annotates the *cells* of the
query; like the other cell-level tools, the cluster-level view comes from
`envs.ident` (or `envs.cases.X.ident`), not from the tool itself.
A `garnett` classifier is trained from marker genes and expression data
(it is model-based, not a plain marker-table input like the marker-based tools).
- Direct-annotation tools: take the cell types directly via `envs.cell_types`,
including `direct` and `cell`.
`scsa`, `maca`, and `scmapnet` are python-based tools, driven through wrapper
scripts under `biopipen/scripts/scrna/` like the other python-based tools
(`celltypist`, `scbert`, `schdeepinsight`, `scagenttype`); their python
dependencies are not installed by biopipen, see the `envs.<tool>` below.
The annotated cell types will be saved to a new column (`envs.anno_col`, default: `CellType`)
in the metadata, so that the downstream processes will use the annotated cell types
(the identity will be set to the annotation column unless `envs.set_ident` is `False`).
/// Note
The original identity column (e.g. `seurat_clusters`) is never modified.
For `tool` set to `direct`, if `envs.cell_types` is not specified or is an empty list,
the original cell types will be kept and nothing will be changed.
///
If you are using a cluster-level tool (or a cell-level tool with `envs.ident`), a text file
containing the mapping from the original identity to the new cell types will be generated
and saved to `cluster2celltype.tsv` under the job output directory.
The per-cell annotations from cell-level tools will be saved to `cell2celltype.tsv`
under the job output directory.
/// Note
The following envs are deprecated and will be removed in future versions:
`sctype_tissue`, `sctype_db`, `hitype_tissue`, `hitype_db`, `scsorter_db`, `scsorter_args`,
`scina_db`, `scina_args`, `singler_db`, `singler_args`, `schdeepinsight_ref`,
`schdeepinsight_args`, `llmcelltype_args`, `cellassign_db`, `cellassign_args`,
`scbert_ref`, `scbert_model`, `scbert_label_dict`, `scbert_args`, `cellid_db`,
`cellid_args`, `sccatch_args`, `celltypist_args`, `newcol`, and `backup_col`.
Use the corresponding `envs.<tool>` namespace instead (e.g. `envs.sctype_db` →
`envs.sctype.db`). The deprecated envs still work, with a warning, and take precedence
over the new-style envs when both are provided.
`envs.newcol` is replaced by `envs.anno_col`, and `envs.backup_col` is no longer
needed (the original identity column is never modified).
///
/// Note
### Universal marker format
The marker-based tools (`sctype`, `hitype`, `sccatch`, `scsorter`, `scina`,
`cellassign`, and `cellid`) accept a universal marker table in addition to
their native formats. The table can be a TSV, CSV, or an RDS/qs/qs2 file
containing a data.frame, in long format with one row per gene per cell type:
- `cell_type` (required): the cell type.
- `gene` (required): the marker gene.
- `direction`: `positive`/`negative` (aliases: `pos`/`neg`/`+`/`-`).
For `sctype`/`hitype`, negative markers are used as down-regulated markers
(`geneSymbolmore2`). For `scsorter`, negative markers become negative
`Weight`s. For the other tools (`scina`, `cellassign`, `cellid`,
`sccatch`), negative markers cannot be represented and are ignored
(only positive markers are used).
- `weight`: a numeric weight. Used by `scsorter` (as the `Weight`
column) and, when present, by `hitype` (hitype >= 0.0.6), where the
weights are used as-is in the scoring. Ignored by the other tools.
- `species`, `cancer`, `tissue`: optional. When the matching env
(`envs.<tool>.species`/`cancer`/`tissue`) is set, only the rows with the
given value are kept (an error is raised if the table has no such column
or no rows match). For `sctype`/`hitype`, `tissue` also becomes the
`tissueType` column, and it can be used to filter a native ScType xlsx/TSV
database as well (the other two columns only exist in universal tables).
- `level`: an integer, only used by `sctype`/`hitype`.
Column aliases are auto-detected: `celltype`/`cellType`/`Type` → `cell_type`,
`marker`/`Marker`/`gene_symbol` → `gene`, `sign` → `direction`, and
`tissueType` → `tissue`.
You can specify the columns after `#` in the file path, for example,
`file:///path/to/markers.tsv#cell_type,gene,direction`, using column names
(aliases allowed) or 1-based indices. The `file://` prefix is optional.
A file without `cell_type`/`gene` columns is treated as the tool's native
format (e.g. a ScType xlsx for `sctype`/`hitype`, a per-cell-type-column CSV
for `scina`, a named-list RDS for `scina`/`cellassign`/`cellid`, an RDS
data.frame for `sccatch`, etc.).
///
Examples:
```toml
[CellTypeAnnotation.envs]
tool = "direct"
cell_types = ["CellType1", "CellType2", "-", "CellType4"]
```
The cell types will be assigned as:
```
0 -> CellType1
1 -> CellType2
2 -> 2
3 -> CellType4
```
Input:
sobjfile: The single-cell object in RDS/qs/qs2/h5ad format.
Output:
outfile: The rds/qs/qs2/h5ad file of seurat object with cell type annotated.
A text file containing the mapping from the old identity to the new cell types
will be generated and saved to `cluster2celltype.tsv` under the job output directory.
Another text file containing the per-cell annotations will be generated and saved
to `cell2celltype.tsv` under the job output directory, with the cell barcodes in
the first column (`Cell`) and one column per case that produces cell-level annotations.
Note that the identity of the output Seurat object will be set to the annotation
column (`envs.anno_col`) when `envs.set_ident` is `True` (see `envs.set_ident`).
Envs:
tool (choice): The tool to use for cell type annotation.
- sctype (cluster-level): Use `scType` to annotate cell types.
See <https://github.com/IanevskiAleksandr/sc-type>
- hitype (cell-level): Use `hitype` to annotate cell types.
Runs at cell-level by default; set `envs.ident` for cluster-level
annotations (hitype then scores each cluster as a whole and
assigns it one cell type). See <https://github.com/pwwang/hitype>
- sccatch (cluster-level): Use `scCATCH` to annotate cell types.
See <https://github.com/ZJUFanLab/scCATCH>
- celltypist (cell-level): Use `celltypist` to annotate cell types.
It can also generate cluster-level annotations via its over-clustering
mechanism (`envs.celltypist.over_clustering` or `envs.ident`).
See <https://github.com/Teichlab/celltypist>
- scsorter (cluster-level): Use `scSorter` to annotate cell types.
See <https://github.com/pwwang/scSorter>, an optimized version of
<https://pmc.ncbi.nlm.nih.gov/articles/PMC7898451/>, which increases the
speed of the original scSorter.
- scina (cell-level): Use `SCINA` to annotate cell types.
See <https://github.com/jcao89757/SCINA>
- singler (cluster-level): Use `SingleR` to annotate cell types.
See <https://github.com/dviraran/SingleR>
- garnett (cell-level): Use `Garnett` to annotate cell types with a
pre-trained classifier (trained from marker genes and expression data).
See <https://cole-trapnell-lab.github.io/garnett/>
- schdeepinsight (cell-level): Use `scHDeepInsight` to annotate cell types.
See <https://github.com/shangruJia/scHDeepInsight>
- llmcelltype (cluster-level): Use `LLMCelltype` to annotate cell types
with LLMs. See <https://github.com/pwwang/LLMCelltype>
It is the model provider agnostic version of
[`gptcelltype`](https://github.com/Winnie09/GPTCelltype).
- cellassign (cell-level): Use `cellassign` to annotate cell types
with a probabilistic model.
See <https://github.com/Irrationone/cellassign>
- scbert (cell-level): Use `scBERT` to annotate cell types with a
BERT-based transformer model.
See <https://github.com/TencentAILabHealthcare/scBERT>
- cellid (cell-level): Use `CelliD` to annotate cell types with MCA-based
per-cell gene signature enrichment.
See <https://github.com/RausellLab/CelliD>
- scagenttype (cluster-level): Use `scAgentType` to annotate cell types
with an agentic LLM workflow.
See <https://github.com/sathyasjali/scAgentType>
- ucell (cell-level): Use `UCell` to annotate cell types by scoring the
marker table with equal weights.
See <https://github.com/carmonalab/UCell>
- aucell (cell-level): Use `AUCell` to annotate cell types by scoring the
marker table with equal weights.
See <https://github.com/aertslab/AUCell>
- gsva (cell-level): Use `GSVA` to annotate cell types by scoring the
marker table with equal weights.
See <https://bioconductor.org/packages/release/bioc/html/GSVA.html>
- singscore (cell-level): Use `singscore` to annotate cell types by scoring
the marker table with equal weights.
See <https://bioconductor.org/packages/release/bioc/html/singscore.html>
- scmap (cell-level): Use `scmap` to transfer the cell types of a
reference object.
See <https://bioconductor.org/packages/release/bioc/html/scmap.html>
- cheetah (cell-level): Use `CHETAH` to transfer the cell types of a
reference object.
See <https://bioconductor.org/packages/release/bioc/html/CHETAH.html>
- scclassify (cell-level): Use `scClassify` to transfer the cell types of
a reference object.
See <https://bioconductor.org/packages/release/bioc/html/scClassify.html>
- scpred (cell-level): Use `scPred` to transfer the cell types of a
reference object.
See <https://github.com/powellgenomicslab/scPred>
- azimuth (cell-level): Use `Azimuth` to transfer the annotation levels of
a published reference.
See <https://github.com/satijalab/azimuth>
- scsa (cluster-level): Use `SCSA` to annotate cell types.
See <https://github.com/bioinfo-ibms-pumc/SCSA>
- maca (cell-level): Use `MACA` to annotate cell types.
See <https://github.com/ImXman/MACA>
- scmapnet (cell-level): Use `scMapNet` to annotate cell types.
See <https://github.com/Yuz7/scMapNet>
- mllmcelltype (cluster-level): Use `mLLMCelltype` to annotate cell types
with LLMs. See <https://github.com/cafferychen777/mLLMCelltype>
- lict (cluster-level): Use `LICT` to annotate cell types with LLMs.
See <https://github.com/Glowworm-cell/LICT>
- mapquery (cell-level): Use `Seurat::MapQuery()` to transfer the cell
types of a reference object.
See <https://satijalab.org/seurat/reference/mapquery>
- direct (cluster-level): Directly assign cell types
- cell (cell-level): Directly assign cell types, but at cell-level instead of cluster-level.
assay: The assay to use for the analysis. If not specified, the default assay will be used.
Will not be inherited by the cases under `envs.cases`.
To use a different assay for a case, specify it in the case args.
This will be also used to convert Seurat object to h5ad if the input is Seurat object
and the output is h5ad.
layer: The layer in the assay to use for the analysis. If not specified, the default layer will be used.
ident: The column name in metadata to use as the clusters.
For cluster-level tools, this is required, and if not specified,
the identity column will be used when input is rds/qs/qs2 (supposing we have a Seurat object).
If input data is h5ad, this is required to run cluster-based annotation tools.
For cell-level tools, if specified, a cluster-level annotation will also be generated
by majority vote of the cells in each cluster, and saved to `cluster2celltype.tsv`.
For `hitype`, this switches to cluster-level scoring (cell-level is the default
when it is not set).
For `celltypist`, this is a shortcut to set `over_clustering` in `envs.celltypist`
(see `envs.celltypist.over_clustering`).
`"ident"` can be used as an alias for the identity column.
To set it for a specific case, use `envs.cases.X.ident`.
anno_col (type=str): The name of the column to store the annotated cell types (default: `CellType`).
For cluster-level tools (or cell-level tools with `envs.ident`), the annotation
column stores the cluster-level cell types.
For cell-level tools, the per-cell annotations are also saved to a tool-specific
column (e.g. `scina_celltype`), and if `envs.ident` is specified, `anno_col` stores
the majority-vote results of the clusters.
For the default case (`DEFAULT`), the column is named as `anno_col`; for other cases,
the case name is prefixed to the column name unless `envs.add_prefix` is `False`.
set_ident (flag): Whether to set the identity of the output Seurat object to the annotation column.
If all cases have `set_ident` set to `False`, the original identity is kept.
If multiple cases have `set_ident` set to `True`, a warning is given and the last case wins.
Can be set per case via `envs.cases.X.set_ident` (default: True).
subset (type=str): An expression to subset the cells for each case.
It is passed to `tidyrseurat::filter()` on the metadata of the Seurat object.
When set, the tool only receives the cells in the subset, and the annotation
columns of the case are `NA` for the cells outside it. This makes
complementary cases cover all cells (e.g., two `mapquery` runs on the
bone-marrow and blood cells), which can be merged back with
`envs.post_mutaters`.
Can be set per case via `envs.cases.X.subset`.
sctype (ns): The arguments for `sctype` if `tool` is `sctype`.
- tissue: The tissue to use for `sctype`.
Available tissues should be the first column (`tissueType`) of `db`.
If not specified, all rows in `db` will be used.
- cancer: Filter the markers by the `cancer` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `cancer` column.
- species: Filter the markers by the `species` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `species` column.
- db: The database to use for sctype.
Check examples at <https://github.com/IanevskiAleksandr/sc-type/blob/master/ScTypeDB_full.xlsx>
Can also be a universal marker table (see the note above).
hitype (ns): The arguments for `hitype` if `tool` is `hitype`.
- tissue: The tissue to use for `hitype`.
Available tissues should be the first column (`tissueType`) of `db`.
If not specified, all rows in `db` will be used.
- cancer: Filter the markers by the `cancer` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `cancer` column.
- species: Filter the markers by the `species` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `species` column.
- db: The database to use for hitype.
Compatible with `sctype.db`.
See also <https://pwwang.github.io/hitype/articles/prepare-gene-sets.html>
You can also use built-in databases, including `hitypedb_short`, `hitypedb_full`, and `hitypedb_pbmc3k`.
Can also be a universal marker table (see the note above).
When the table has a `weight` column (e.g. trained by
`HitypeWeightTrainer`), the weights are used as-is for
scoring (hitype >= 0.0.6).
- norm: The normalization method for `hitype::hitype_score()`.
One of "sqrt", "weight", "none" (default: "sqrt").
"weight" is recommended when scoring with learned weights.
- use_sensitivity: Whether to weight markers by their
sensitivity (default: `True`). `False` is recommended when
scoring with learned weights.
- threshold: The assignment threshold passed to `RunHitype`
(default: `0.0`).
scsorter (ns): The arguments for `scSorter::RunScSorter()` if `tool` is `scsorter`.
- db: The database to use for scSorter. It will be loaded and passed to the `anno`
argument of `RunScSorter()`. It could be either:
* A TSV file with cell type annotations, with columns `Type`, `Marker`, and `Weight`.
* A RDS/qs2 file of the annotation data frame with the same columns as above.
Can also be a universal marker table (see the note above).
You can also use `#` followed by the column names (aliases allowed)
or 1-based indices to specify the columns, for example,
`file:///path/to/scsorter_db.tsv#celltype,marker,weight`.
A third column is used as `Weight` only if it is named `weight` (or an alias of it).
- assay: The assay to use for `RunScSorter()`.
If not specified, `envs.assay` will be used.
- tissue: Filter the markers by the `tissue` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `tissue` column.
- cancer: Filter the markers by the `cancer` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `cancer` column.
- species: Filter the markers by the `species` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `species` column.
- <more>: Other arguments for [`scSorter::RunScSorter()`](https://github.com/pwwang/scSorter/blob/9baae9f0e0904ddbf3f9bb5dacb9227503a8ce3e/R/scSorter.R#L73).
scina (ns): The arguments for `SCINA::SCINA()` if `tool` is `scina`.
- db (type=str): The path to the SCINA signature file.
It can be an RDS file containing a named list of
signature genes (the names are the cell types and the
values are the marker gene symbols), or a CSV file
with the markers for each cell type in a column.
Can also be a universal marker table (see the note above).
- tissue: Filter the markers by the `tissue` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `tissue` column.
- cancer: Filter the markers by the `cancer` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `cancer` column.
- species: Filter the markers by the `species` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `species` column.
- max_iter (type=int): Maximum number of EM iterations (default: 100).
- convergence_n (type=int): Stop if assignment stays stable for N consecutive rounds (default: 10).
- convergence_rate (type=float): Fraction of cells with stable assignment for convergence (default: 0.99).
- sensitivity_cutoff (type=float): Cutoff (0-1) for removing signatures of absent cell types (default: 1).
- rm_overlap (flag): Whether to remove genes shared between multiple signatures (default: TRUE).
- allow_unknown (flag): Whether to allow unknown cells (default: TRUE).
- <more>: Other arguments for [`SCINA::SCINA()`](https://rdrr.io/cran/SCINA/man/SCINA.html).
singler (ns): The arguments for `SingleR::SingleR()` if `tool` is `singler`.
Both the [Bioconductor](https://bioconductor.org/packages/SingleR)
and [CRAN](https://github.com/dviraran/SingleR) versions are
auto-detected.
- db (type=str): The path to the SingleR reference file.
It can be an RDS, qs, or qs2 file containing a reference object,
supporting:
* `SummarizedExperiment` (e.g., from the `celldex` package).
References can be obtained via `celldex::HumanPrimaryCellAtlasData()`,
`celldex::BlueprintEncodeData()`, `celldex::MonacoImmuneData()`,
`celldex::DatabaseImmuneCellExpressionData()`,
`celldex::NovershternHematopoieticData()`, `celldex::ImmGenData()`,
`celldex::MouseRNAseqData()`.
* `Seurat` object. Labels are auto-detected from metadata.
Save with `saveRDS()` or `biopipen.utils::write_obj()`.
Both the Bioconductor and CRAN versions of SingleR are supported
and auto-detected at runtime.
- label (type=str): The metadata/colData column name for
reference labels. Auto-detected from `label.main`,
`label.fine`, `label.ont`, `label` in order.
- <more>: See the SingleR documentation for your version:
[`Bioconductor`](https://rdrr.io/bioc/SingleR/man/SingleR.html)
or [`CRAN`](https://github.com/dviraran/SingleR).
garnett (ns): The arguments for `garnett::classify_cells()`
if `tool` is `garnett`. The process takes a **trained** classifier
and predicts the cell types for each cell — no training step is
performed. A classifier can be trained with the
`GarnettClassifierTrainer` process from a Seurat object and marker
genes, or you can use a pre-trained one from
<https://cole-trapnell-lab.github.io/garnett/classifiers/>.
The Seurat object is converted to a monocle3 `cell_data_set` for
classification, and cells that cannot be confidently classified
are labeled `Unknown`.
- classifier (type=str): The path to the trained classifier
file (an RDS file containing a `garnett_classifier` object,
saved by `garnett::save_classifier()`).
- db (type=str): The name of the installed annotation package
used to convert the gene IDs of the expression data to the
gene IDs the classifier was trained on (default: `none`, i.e.
no conversion). The official pre-trained classifiers are
trained on ENSEMBL gene IDs, so for expression data with gene
symbols you need e.g. `db: org.Hs.eg.db` together with
`cds_gene_id_type: SYMBOL`.
- cds_gene_id_type (choice): The gene ID type of the expression
data, used when `db` is not `none` (default: `custom`).
Supported types include `custom`, `SYMBOL`, `ENSEMBL`, `ENTREZID`.
- assay (type=str): The assay to use for classification.
If not specified, `envs.assay` will be used.
The assay must contain raw counts.
- cluster_extend (flag): Use Garnett's own cluster mode: the
clusters are passed as the `garnett_cluster` column of the
`cell_data_set`, and `classify_cells()` extends the labels over
each cluster (`cluster_ext_type`) instead of the engine
aggregating the per-cell labels by majority vote
(default: FALSE).
- <more>: Other arguments for
[`garnett::classify_cells()`](https://github.com/cole-trapnell-lab/garnett/blob/monocle3/R/classify.R),
e.g. `rank_prob_ratio`,
`cluster_extend_max_frac_unknown`, `cluster_extend_max_frac_incorrect`,
`return_type_levels`, and `verbose`.
schdeepinsight (ns): The arguments for scHDeepInsight
if `tool` is `schdeepinsight`. Gated by default: the
`SCHdeepinsight` package (0.3.5) is not installed, and the tool
additionally needs a reference file and a pretrained checkpoint.
To un-gate it, run `pip install SCHdeepinsight` (plus
`pip install git+https://github.com/alok-ai-lab/pyDeepInsight.git`)
in `python` and download the checkpoint from
<https://github.com/shangruJia/scHDeepInsight>.
- ref (type=str): The path to the scHDeepInsight
reference RDS file. The bundled `reference.rds` from the
scHDeepInsight repo provides immune cell reference.
See <https://github.com/shangruJia/scHDeepInsight>.
- batch_size (type=int): Batch size for CNN prediction
(default: 128).
- python (type=str): Path to Python executable with
`SCHdeepinsight` installed.
- assay (type=str): Assay to use for h5ad conversion.
llmcelltype (ns): The arguments for
`LLMCelltype::llmcelltype()` if `tool` is `llmcelltype`.
- api_key (type=str): OpenAI API key
- model (type=str): GPT model (required, e.g.
'gpt-4', 'gpt-4o').
- base_url (type=str): Custom base URL for
OpenAI-compatible providers.
- tissuename (type=str): Tissue name for context.
- sigmarkers (type=str): A expression to filter the result from `RunSeuratDEAnalysis()`, e.g. `avg_log2FC > 0.25 & p_val_adj < 0.05`,
to be used for generating the marker gene list for LLM.
- assay (type=str): Assay to use for
`FindAllMarkers()`.
- <more>: Additional args passed to
[`biopipen.utils::RunSeuratDEAnalysis()`](https://pwwang.github.io/biopipen.utils.R/reference/RunSeuratDEAnalysis.html).
cellassign (ns): The arguments for
`cellassign::cellassign()` if `tool` is `cellassign`.
- db (type=str): The path to the marker gene info
file for `cellassign`. Supports:
* RDS/qs2 file: a binary gene×celltype matrix or a
named list (cell type → vector of marker genes)
* CSV/TSV file: with columns `gene` and `cell_type`
Can also be a universal marker table (see the note above).
- tissue: Filter the markers by the `tissue` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `tissue` column.
- cancer: Filter the markers by the `cancer` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `cancer` column.
- species: Filter the markers by the `species` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `species` column.
- python (type=str): Path to Python with `tensorflow` installed.
- assay (type=str): Assay to extract raw counts from.
- min_delta (type=int): Min log-fold change for marker
overexpression (default: 2).
- B (type=int): Number of RBF dispersion bases
(default: 10).
- shrinkage (flag): Hierarchical shrinkage on delta
(default: TRUE).
- n_batches (type=int): Data subsample batches
(default: 1).
- learning_rate (type=float): ADAM learning rate
(default: 0.1).
- max_iter_em (type=int): Max EM iterations
(default: 20).
- verbose (flag): Print progress (default: TRUE).
- <more>: Additional args to
`cellassign::cellassign()`.
scbert (ns): The arguments for scBERT inference
if `tool` is `scbert`.
- ref (type=str): The path to the scBERT repo
directory (containing `performer_pytorch/`).
- model (type=str): The path to the fine-tuned
model checkpoint (.pth file).
- label_dict (type=str): The path to the label
dictionary pickle file (maps class indices to cell
type names).
- python (type=str): Path to Python with scBERT
dependencies (torch, scanpy, etc.).
- bin_num (type=int): Number of bins for
expression embedding (default: 5).
- gene_num (type=int): Number of genes expected
by the model (default: 16906).
- seed (type=int): Random seed (default: 2021).
- pos_embed (flag): Use Gene2vec positional
encoding (default: TRUE).
- novel_type (flag): Enable novel cell type
detection (default: FALSE).
- unassign_thres (type=float): Confidence
threshold for unassigned cells (default: 0.5).
- <more>: Additional args to the wrapper script.
cellid (ns): The arguments for CelliD
if `tool` is `cellid`.
- db (type=str): The path to the marker gene set
file for `cellid`. Supports:
* RDS/qs2 file: a named list (cell type → vector
of marker genes)
* CSV/TSV file: with columns `gene` and `cell_type`
Can also be a universal marker table (see the note above).
- tissue: Filter the markers by the `tissue` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `tissue` column.
- cancer: Filter the markers by the `cancer` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `cancer` column.
- species: Filter the markers by the `species` column of a universal
marker table (see the note above). Only works with a universal
marker table that has a `species` column.
- nmcs (type=int): Number of MCA components
(default: 50).
- n_features (type=int): Top n features per
cell for hypergeometric test (default: 200).
- dims (type=auto): MCA dimensions to use
(default: seq(nmcs)).
- min_size (type=int): Min overlapping genes
(default: 10).
- log_trans (flag): -log10 transform p-values
(default: TRUE).
- p_adjust (flag): Benjamini-Hochberg correction
(default: TRUE).
- group_gsea (flag): Use CelliD's own cluster mode
(`CelliD::RunGroupGSEA()`): the clusters are annotated from the
gene-set enrichment scores of the whole group (the best-scoring
pathway by NES) instead of the engine aggregating the per-cell
hypergeometric test by majority vote (default: FALSE).
cell_types (type=auto): The cell types to use for direct or cell-level annotation.
For `direct`, the cell types will be assigned to the clusters in the order of the original identities.
If given as a list (array), you can use `"-"` or `""` as the placeholder for the clusters that
you want to keep the original cell types. If the length of `cell_types` is shorter than the number of
clusters, the remaining clusters will be kept as the original cell types.
You can also use `NA` to remove the clusters from downstream analysis
(the cells in these clusters will be removed from the Seurat object).
If given as a dict (map), the keys are the original cluster names and the values are the new cell types.
For `cell`, it must be a TSV file with cell-level annotations.
You can specify the column names after the `#`. For example,
`file:///path/to/cell_types.tsv#cell_id,cell_type` will use `cell_id` as the cell id column
to match the cell ids in the Seurat object, and `cell_type` as the cell type column to assign the cell types.
Multiple cell type columns can be specified, and the first one will be used as the annotation column
(the others will be added as additional annotation columns).
You can also use 1-based column index to specify the columns, for example, `file:///path/to/cell_types.tsv#1,3`
will use the first column as the cell id column and the third column as the cell type column.
If cells in the Seurat object are not found in the cell type file, `NA`s will be assigned to those cells.
If no columns are specified, the first two columns will be used as the cell id and cell type columns.
Prefix `file://` is optional.
For `cell`, it must be a TSV file with cell-level annotations.
You can specify the column names after the `#`. For example,
`file:///path/to/cell_types.tsv#cell_id,cell_type` will use `cell_id` as the cell id column
to match the cell ids in the Seurat object, and `cell_type` as the cell type column to assign the cell types.
Multiple cell type columns can be specified, and the first one will be used as the new identity column.
You can also use 1-based column index to specify the columns, for example, `file:///path/to/cell_types.tsv#1,3`
will use the first column as the cell id column and the third column as the cell type column.
If cells in the Seurat object are not found in the cell type file, `NA`s will be assigned to those cells.
If not columns are specified, the first two columns will be used as the cell id and cell type columns.
Prefix `file://` is optional.
more_cell_types (type=json): The additional cell type annotations to add to the metadata.
The keys are the new column names and the values are the cell types lists.
The cell type lists work the same as `cell_types` above.
This is useful when you want to keep multiple annotations of cell types.
sccatch (ns): The arguments for `scCATCH::findmarkergene()` if `tool` is `sccatch`.
- species: The specie of cells.
When `marker` is a custom (universal) marker table, only the rows
matching the value are kept; otherwise it is used to filter the
built-in scCATCH database.
- cancer: If the sample is from cancer tissue, then the cancer type may be defined.
Defaults to "Normal" if not `if_use_custom_marker`.
When `marker` is a custom (universal) marker table, only the rows
matching the value are kept; otherwise it is used to filter the
built-in scCATCH database.
- tissue: Tissue origin of cells must be defined.
When `marker` is a custom (universal) marker table, only the rows
matching the value are kept; otherwise it is used to filter the
built-in scCATCH database.
- marker: The marker genes for cell type identification.
Can also be a universal marker table (see the note above).
An error is raised if `species`, `cancer`, or `tissue` is set but
the table has no such column or no rows match.
- if_use_custom_marker (flag): Whether to use custom marker genes.
When `marker` is provided, this is set to `True` automatically, and
`species`, `cancer`, and `tissue` filter the custom markers if set.
- <more>: Other arguments for [`scCATCH::findmarkergene()`](https://rdrr.io/cran/scCATCH/man/findmarkergene.html).
You can pass an RDS file to `sccatch.marker` to work as custom marker. If so,
`if_use_custom_marker` will be set to `TRUE` automatically.
celltypist (ns): The arguments for `celltypist::celltypist()` if `tool` is `celltypist`.
- model: The path to model file.
- python: The python path where celltypist is installed.
- majority_voting: When true, it refines cell identities within local subclusters after an over-clustering approach
at the cost of increased runtime.
- over_clustering (type=auto): The column name in metadata to use as clusters for majority voting.
Set to `False` to disable over-clustering.
When `in.sobjfile` is rds/qs/qs2 (supposing we have a Seurat object), the default ident is used by default.
Otherwise, it is False by default.
- assay: When converting a Seurat object to AnnData, the assay to use.
If input is h5seurat, this defaults to RNA.
If input is Seurat object in RDS, this defaults to the default assay.
scagenttype (ns): The arguments for the scAgentType annotation agent
if `tool` is `scagenttype`. It annotates each cluster via an agentic
LLM workflow (ReAct) and needs Python >=3.10 with the `scagenttype`
package installed, e.g.
`pip install "scagenttype[llm] @ git+https://github.com/sathyasjali/scAgentType.git"`.
- python: The python path where `scagenttype` is installed.
- api: The LLM provider: `openai`, `anthropic`, or `google` (default: `openai`).
- api_key: The API key for the provider.
When not set, the key is read from the environment variable of
the provider (`OPENAI_API_KEY`, `ANTHROPIC_API_KEY`, or `GOOGLE_API_KEY`).
- model: The model to use. Defaults are chosen per provider when not set.
- base_url: Custom API base URL (e.g. for proxies).
Passed to the subprocess via `OPENAI_BASE_URL` or `ANTHROPIC_BASE_URL`.
- tissue: The tissue of the data, e.g. `Human peripheral blood`.
Folded into `tissue_context` when `tissue_context` is not set.
- species: The species of the data.
Folded into `tissue_context` when `tissue_context` is not set.
- assay: When converting a Seurat object to AnnData, the assay to use.
- <more>: Other arguments for [`AnnotationAgent()`](https://github.com/sathyasjali/scAgentType),
e.g. `tissue_context`, `n_markers`, `max_react_steps`,
`confidence_threshold`, `self_consistency_n`, `cache_dir`,
`enable_cellxgene`.
ucell (ns): The arguments for [`UCell::AddModuleScore_UCell()`](https://github.com/carmonalab/UCell)
if `tool` is `ucell`. UCell scores the universal marker table with
equal (unit) weights, and each cell gets the cell type of its highest
score (see the marker-based note above).
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- assay (type=str): The assay to score. If not specified,
the default assay will be used.
- maxRank (type=int): The number of top-ranked genes used to score
a signature (default: 1000).
- w_neg (type=float): The weight of the negative markers
(default: 1).
- name (type=str): The suffix of the score columns added to the
metadata (default: `_UCell`).
- tissue: Filter the markers by the `tissue` column of the universal
marker table (see the note above).
- cancer: Filter the markers by the `cancer` column of the universal
marker table (see the note above).
- species: Filter the markers by the `species` column of the universal
marker table (see the note above).
aucell (ns): The arguments for [`AUCell::AUCell_calcAUC()`](https://github.com/aertslab/AUCell)
if `tool` is `aucell`. The negative-direction markers of the marker
table are dropped (a ranked list has no direction).
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- assay (type=str): The assay to score. If not specified,
the default assay will be used.
- aucMaxRank (type=int): The number of top-ranked genes to
calculate the AUC on (default: 5% of the genes).
- normAUC (flag): Whether to normalize the AUC to the maximum
possible AUC (default: TRUE).
- tissue (type=str): Filter the marker table by tissue, e.g.
`Immune system`.
- cancer (type=str): Filter the marker table by cancer, e.g.
`Breast cancer`.
- species (type=str): Filter the marker table by species, e.g.
`Human`.
gsva (ns): The arguments for [`GSVA::gsva()`](https://bioconductor.org/packages/release/bioc/html/GSVA.html)
if `tool` is `gsva`. The negative-direction markers of the marker
table are dropped (a gene set has no direction).
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- assay (type=str): The assay to score. If not specified,
the default assay will be used.
- kcdf (choice): The kernel to use for the enrichment scores.
`Gaussian` (default), `Poisson`, or `none`.
- minSize (type=int): The minimum number of genes in a gene set
(default: 1).
- maxSize (type=int): The maximum number of genes in a gene set
(default: `Inf`, i.e. no limit).
- tissue (type=str): Filter the marker table by tissue, e.g.
`Immune system`.
- cancer (type=str): Filter the marker table by cancer, e.g.
`Breast cancer`.
- species (type=str): Filter the marker table by species, e.g.
`Human`.
singscore (ns): The arguments for [`singscore::simpleScore()`](https://bioconductor.org/packages/release/bioc/html/singscore.html)
if `tool` is `singscore`. Unlike the other scorers, the direction of
the markers is used as-is (`upSet`/`downSet`).
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- assay (type=str): The assay to score. If not specified,
the default assay will be used.
- subSamples (type=int): The number of random subsets of the ranked
genes to score. Cells left out are assigned `NA`.
- centerScore (flag): Whether to center the scores to a [0, 1]
range.
- tissue (type=str): Filter the marker table by tissue, e.g.
`Immune system`.
- cancer (type=str): Filter the marker table by cancer, e.g.
`Breast cancer`.
- species (type=str): Filter the marker table by species, e.g.
`Human`.
scmap (ns): The arguments for `scmap::scmapCluster()` (`scmap::scmapCell()`
when `use_cell_index` is `True`) if `tool` is `scmap`. The reference
is a `SingleCellExperiment` (or a Seurat object, converted on the fly)
with a `logcounts` assay.
- db (type=str): The path to the reference (required).
- assay (type=str): The assay of the reference and of the object
(default: `RNA`).
- cluster_col (type=str): The column of the reference's `colData`
holding the cell types (default: `cell_type1`).
- features (type=auto): The features used for the scmap index.
Defaults to all the features of the reference.
- threshold (type=float): The similarity threshold below which a
cell is labeled `unassigned` (default: 0.5).
- use_cell_index (flag): Use the scmap-cell (k-means) index instead
of the scmap-cluster index (default: FALSE).
cheetah (ns): The arguments for [`CHETAH::CHETAHclassifier()`](https://bioconductor.org/packages/release/bioc/html/CHETAH.html)
if `tool` is `cheetah`. The reference is a `SingleCellExperiment`
(or a Seurat object, converted on the fly) holding counts.
- db (type=str): The path to the reference (required).
- assay (type=str): The assay of the reference and of the object
(default: `RNA`).
- input_c (type=str): The name of the assay of the input to use.
- thresh (type=float): The confidence threshold below which a cell
is labeled `Unassigned` (default: 0.1).
- n_genes (type=int): The number of genes used for the classification.
- pc_thresh (type=float): The threshold for the principal components.
- only_pos (flag): Whether to use only the positively correlated
genes. Defaults to the CHETAH default.
- ref_ct (type=str): The reference's cell-type column name
(auto-detected when not given).
- label (type=str): The label vector/column used when the reference
carries several.
scclassify (ns): The arguments for `scClassify::scClassify()` /
`scClassify::predict_scClassify()` if `tool` is `scclassify`.
The reference is either a labelled matrix (a list with `exprsMat` and
`cellTypes`), or a pre-trained model from
`scClassify::train_scClassify()`.
- db (type=str): The path to the reference (required).
- assay (type=str): The assay to classify (default: `RNA`).
- algorithm (choice): The algorithm to use. `WKNN` (default), `KNN`,
or `DWKNN` (the `lr` of the scClassify docs is not supported by
the installed scClassify).
- topN (type=int): The number of features used by the model
(default: 50).
- prob_threshold (type=float): The probability threshold below which
a cell is labeled `unassigned`.
- parallel (flag): Whether to run scClassify in parallel.
scpred (ns): The arguments for [`scPred::scPredict()`](https://github.com/powellgenomicslab/scPred)
if `tool` is `scpred`. The reference is a labelled Seurat object, and
the feature space is extracted and the model trained on it at runtime.
- db (type=str): The path to the reference (required).
- assay (type=str): The assay to predict.
- pvar (type=str): The metadata column of the reference holding the
cell types (default: `cell_type`).
- model (type=str): The model to train (default: `svmRadial`).
- reduction (type=str): The reduction used for the feature space
(default: `pca`).
- threshold (type=float): The prediction-probability threshold
passed to `scPredict()` (default: 0.55).
azimuth (ns): The arguments for [`Azimuth::RunAzimuth()`](https://github.com/satijalab/azimuth)
if `tool` is `azimuth`. The reference's own annotation levels are
transferred to the query, and each cluster gets the majority call of
its cells.
- ref (type=str): The reference name (looked up through `SeuratData`,
e.g. `pbmcref`, downloaded when not installed) or the path to a
directory holding `ref.Rds` + `idx.annoy` (required).
- db (type=str): An alias of `ref`.
- assay (type=str): The assay to transfer.
- anno_col_in (type=str): The `predicted.*` column of the Azimuth
result to use as the annotation. Defaults to
`predicted.celltype.l1` when present, otherwise the first
`predicted.*` column.
- dims: Ignored with a warning: the installed Azimuth reads the
dimensionality off the reference's own annoy index.
- k.anchor: Ignored with a warning, see `dims`.
scsa (ns): The arguments for [`SCSA`](https://github.com/bioinfo-ibms-pumc/SCSA)
if `tool` is `scsa`. SCSA is not on CRAN/Bioconductor/PyPI and its
`SCSA.py` does not run on current numpy/pandas, so the wrapper
lives in `biopipen/scripts/scrna/scsa-wrapper.py` and ports its
scoring; the marker table is all it needs (SCSA's own reference
database is not used).
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- foldchange (type=float): The minimum fold change of a marker for
it to be used (default: 2.0, SCSA's own `-f`).
- pvalue (type=float): The maximum adjusted p-value of a marker for
it to be used (default: 0.05, SCSA's own `-p`). A cluster with
no marker left is not annotated.
- python (type=str): The python path with the SCSA dependencies
(`pandas`, `numpy`, `scanpy`) installed
(default: the same python as the pipeline).
maca (ns): The arguments for [`MACA`](https://github.com/ImXman/MACA)
if `tool` is `maca`. MACA pins `scanpy==1.6.0` and `anndata==0.7.5`,
so it needs its own environment. The MACA script lives in
`biopipen/scripts/scrna/maca-wrapper.py`.
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- python (type=str): The python path with MACA installed
(default: the same python as the pipeline).
- n_pcs (type=int): The number of principal components.
- res (type=auto): The Louvain resolutions, e.g. `[1, 2, 3]`
(default: MACA's own `[1, 2, 3]`).
- n_neis (type=auto): The numbers of neighbors, e.g. `[5, 10]`
(default: MACA's own `[5, 10]`).
- freq (type=float): The frequency threshold of the cluster mapping
(default: 0.5).
- use_weight (flag): Whether to weight the markers by their order in
the marker table (default: FALSE).
scmapnet (ns): The arguments for [`scMapNet`](https://github.com/Yuz7/scMapNet)
if `tool` is `scmapnet`. scMapNet turns each cell into a treemap image
and classifies it with a vision transformer. Its pipeline (treemap
image generation + `main_finetune.py`) is driven by
`biopipen/scripts/scrna/scmapnet-wrapper.py`.
The pre-trained weights are a manual download and are licensed
CC BY-NC 4.0 (non-commercial).
- db (type=str): The path to the marker table (required).
Must be a universal marker table (see the note above).
- python (type=str): The python path with the scMapNet dependencies
(`torch`, `timm`, `torchvision`) installed
(default: the same python as the pipeline).
- scmapnet_dir (type=str): The path to the cloned scMapNet repo
(required), i.e. the directory holding `main_finetune.py` and
`generate_image_script.sh`.
- weights (type=str): The path to the checkpoint used for the
prediction (required). The pre-trained weights are not part of
the repo, see <https://github.com/Yuz7/scMapNet>. The checkpoint
must be fine-tuned on the cell types of the marker table, and
the cell types are indexed in their sorted order.
- organ (type=str): The organ of the cells. scMapNet builds the
treemap images from an organ → cell type → gene hierarchy, so
this is the `organ` of every marker in the table.
mllmcelltype (ns): The arguments for
`mLLMCelltype::annotate_cell_types()` if `tool` is `mllmcelltype`.
One call is made per run, with the top markers of every cluster.
- tissue (type=str): The tissue of the cells (required),
e.g. `human PBMC`.
- model (type=str): The LLM to use (default: `gpt-5.5`).
The model decides the provider, and so which API key is needed.
- api_key (type=str): The API key of the provider. When not set, the
key is read from the environment variable of the provider
(`OPENAI_API_KEY` or `ANTHROPIC_API_KEY`); without one, only the
prompt is built and an error is raised instead of annotating.
- top_gene_count (type=int): The number of top markers per cluster
put into the prompt (default: 10).
- base_urls (type=str): The base URLs of the providers, for
OpenAI-compatible endpoints. A single URL is the request URL
mLLMCelltype posts to, so a host-only one (or a `.../v1`) is
completed with the `/chat/completions` path; a named list of
per-provider URLs is used as it is.
- return_reasoning (flag): Whether to return the reasoning of the
model along with the cell types (default: FALSE).
lict (ns): The arguments for `LICT::LLMCellType()` if `tool` is `lict`.
LICT queries every provider at once and combines the answers; the
providers without a key are skipped.
- species (type=str): The species of the cells (default: `Human`).
- tissue (type=str): The tissue of the cells.
- topgenenumber (type=int): The number of top markers per cluster
put into the prompt (default: 30).
- validate (flag): Whether to validate the markers
(default: TRUE).
- percent (type=float): The percentage of the markers to use
(default: 0.5).
- keys (type=json): The API keys of the providers, as a map from the
environment variable name to its value, e.g.
`{ "openai_api_key": "<key>" }`. They are set in the R session
(and inherited by the child processes) before the call.
- provider (type=str): The provider whose answer to use. Defaults to
the first provider that answered.
mapquery (ns): The arguments for [`Seurat::MapQuery()`](https://satijalab.org/seurat/reference/mapquery)
if `tool` is `mapquery`. The reference's cell types are transferred to
the cells of the query; as for the other cell-level tools, each cluster
also gets the majority call of its cells when `envs.ident` is set.
- db (type=str): The reference Seurat object file (required).
- use (type=str): The metadata column of the reference holding the
cell types to transfer (required).
- ident_name (type=str): The metadata column of the query that receives
the transferred labels (kept distinct from `ident` so the clustering
column is not overwritten) (default: `predicted.id`).
- refnorm (choice): The normalization method the reference used; the
same method is used for the query (default: `auto`).
- LogNormalize: Normalize the query with `NormalizeData()`.
- SCT: Normalize the query with `SCTransform()`.
- SCTransform: An alias of `SCT`.
- auto: Automatically detect the normalization method.
- skip_if_normalized (flag): Whether to skip the normalization if the
query is already normalized (default: TRUE).
- ncores (type=int): Number of cores to use for the mapping.
- map_query (ns): The arguments for `Seurat::MapQuery()`.
- find_transfer_anchors (ns): The arguments for `Seurat::FindTransferAnchors()`.
- sctransform (ns): The arguments for `Seurat::SCTransform()`.
- normalize_data (ns): The arguments for `Seurat::NormalizeData()`.
scinfer (ns): The arguments for `ScInfeR::predict_celltype_scRNA_seurat()`
if `tool` is `scinfer`. ScInfeR is a graph-based annotator (GPL-3,
GitHub-only, not on CRAN): it scores the cell types of the marker table
against the top cluster markers, and resolves a cluster cell by cell
(over the reduction's neighborhood graph) when a second cell type scores
within 0.1 of the best, so the per-cell labels are saved to
`cell2celltype.tsv` along with the cluster mapping.
One of `db` and `tissue` is required, and the `tissue` route needs
network access.
- db (type=str): The path to the marker table (required, unless
`tissue` is set). Must be a universal marker table (see the note
above); a `weight` column is used as the marker weight
(default: 1), and the negative-direction markers are dropped
(ScInfeR has no direction, so a negative marker would be counted
as a marker of the cell type).
- tissue (type=str): The tissue of the cells. With `db`, it filters
the marker table by its `tissue` column (see the note above).
Without `db`, it selects one of ScInfeR's own tissues: `Bladder`,
`Bone marrow`, `Brain`, `Breast`, `Eye`, `Heart`, `Intestine`,
`Kidney`, `Liver`, `Lungs`, `Pancreas`, `PBMC`, or `Skin`.
ScInfeRDB is not part of the package, so the sheet of the tissue
is downloaded from GitHub at call time -- this route needs
network access, and an unknown tissue is an error.
- subtype_present (flag): Whether the marker table carries subtypes,
i.e. whether the cell types of `envs.scinfer.subtype_info` are
rows of the marker table as well (default: FALSE). They are then
held out of the round-1 scoring of their parent cell type, and
each cell gets its subtype in the second round.
- subtype_info (type=json): The subtypes of each cell type, as a named
list mapping a cell type of the marker table to its subtypes,
e.g. `{ "T-cell": ["CD8+ T-cell", "CD4+ T-cell"], "B cell": ["Naive B cell", "Memory B cell"] }`.
Only read when `subtype_present` is `True`; a cell type without an
entry (or without the list at all) keeps its round-1 label.
- assay (type=str): The assay of the object to annotate
(default: `RNA`).
- slot (type=str): The slot of the assay (default: `counts`).
Ignored by ScInfeR in this version: it reads the default layer of
the assay whichever slot is asked for.
- reduction (type=str): The reduction ScInfeR walks to resolve the
cells of a mixed cluster (default: `umap`). The object must carry
it; the run stops if it does not.
- own_weightage (type=float): The weight of a cell's own expression
(against that of its neighbors) in the per-cell prediction, from
0 to 1 (default: 0.5).
- n_neighbor (type=int): The number of neighbors whose expression is
averaged in the per-cell prediction, from 5 to 30 (default: 10).
merge (flag): Whether to merge the clusters with the same cell types.
Otherwise, a suffix will be added to the cell types (ie. `.1`, `.2`, etc).
post_mutaters (type=json): Mutaters to apply to the metadata after all cases
are embedded and after the identity is set (see `envs.set_ident`).
They are applied via `biopipen.utils::MutateSeuratMeta()`.
The key with the `:ident` suffix creates the column and sets it as the
identity, e.g. `CellType:ident = "coalesce(BM_CellType, PB_CellType)"`
merges two per-case columns into `CellType` and makes it the identity.
Not inherited by the cases.
add_prefix (flag): Whether to add a prefix to the new column names in metadata.
Only used when in non-default cases. The prefix will be the case name followed by `_`.
cases (type=json): Run multiple cases of cell type annotation on the same Seurat object.
The keys are the prefix of column names added to the metadata (unless `add_prefix` is `False`),
and the values will inherit the above options (including `subset`).
The default case is `DEFAULT`, meaning no prefix will be added to the column names.
The identity of the output Seurat object will be set according to `envs.set_ident` of each case
(see `envs.set_ident`).
If any case requires h5ad conversion (e.g., `celltypist`), the object is pre-converted
once and shared across those tools.
ncores (type=int): Number of cores to use for parallel execution of multiple cases.
When > 1, cases are run in parallel via `mclapply`. This is not inherited by individual cases.
outtype (choice): The output file type. Currently only works for `celltypist`.
An RDS file will be generated for other tools.
- input: Use the same file type as the input.
- rds: Use RDS file.
- qs: Use qs2 file.
- qs2: Use qs2 file.
- h5ad: Use AnnData file.
Requires:
r-HGNChelper:
- if: {{proc.envs.tool == 'sctype'}}
- check: {{proc.lang}} -e "library(HGNChelper)"
r-seurat:
- if: {{proc.envs.tool == 'sctype'}}
- check: {{proc.lang}} -e "library(Seurat)"
r-dplyr:
- if: {{proc.envs.tool == 'sctype'}}
- check: {{proc.lang}} -e "library(dplyr)"
r-openxlsx:
- if: {{proc.envs.tool == 'sctype'}}
- check: {{proc.lang}} -e "library(openxlsx)"
r-monocle3:
- if: {{proc.envs.tool == 'garnett'}}
- check: {{proc.lang}} -e "library(monocle3)"
r-garnett:
- if: {{proc.envs.tool == 'garnett'}}
- check: {{proc.lang}} -e "library(garnett)"
r-seuratwrappers:
- if: {{proc.envs.tool == 'garnett'}}
- check: {{proc.lang}} -e "library(SeuratWrappers)"
""" # noqa: E501
input = "sobjfile:file"
output = (
"outfile:file:"
"{{in.sobjfile | stem}}.annotated."
"{{- ext0(in.sobjfile) if envs.outtype == 'input' else envs.outtype -}}"
)
lang = config.lang.rscript
envs = {
"tool": "hitype",
"assay": None,
"layer": None,
"ident": None,
"anno_col": "CellType",
"set_ident": True,
"subset": None,
"cell_types": [],
"more_cell_types": None,
"sctype": { # tested
"tissue": None,
"cancer": None,
"species": None,
"db": config.ref.sctype_db,
},
"hitype": { # tested
"tissue": None,
"cancer": None,
"species": None,
"db": None,
"norm": "sqrt",
"use_sensitivity": True,
"threshold": 0.0,
},
"scsorter": { # tested
"db": None,
"assay": None,
"tissue": None,
"cancer": None,
"species": None,
},
"scina": { # tested
"db": None,
"tissue": None,
"cancer": None,
"species": None,
},
"singler": { # tested
"db": None,
},
"garnett": { # tested locally
"classifier": None,
"db": "none",
"cds_gene_id_type": "custom",
"assay": None,
"cluster_extend": False,
},
"schdeepinsight": { # test needed
"ref": None,
},
"llmcelltype": { # tested locally
"cache": config.path.tmpdir,
"sigmarkers": "p_val_adj < 0.05",
},
"cellassign": { # test needed
"db": None,
"python": config.lang.python,
"tissue": None,
"cancer": None,
"species": None,
},
"scbert": { # models are not available, test needed
"ref": None,
"model": None,
"label_dict": None,
},
"cellid": { # tested
"db": None,
"tissue": None,
"cancer": None,
"species": None,
"group_gsea": False,
},
"sccatch": { # tested
"species": None,
"cancer": None,
"tissue": None,
"marker": None,
"if_use_custom_marker": False,
},
"celltypist": { # tested
"model": None,
"python": config.lang.python,
"majority_voting": True,
"over_clustering": None,
"assay": None,
},
"scagenttype": { # tested locally
"python": config.lang.python,
"api": "openai",
"api_key": None,
"model": None,
"base_url": None,
"tissue": None,
"species": None,
"assay": None,
},
"ucell": { # tested
"db": None, # required
"assay": None,
"maxRank": 1000,
"w_neg": 1,
"name": "_UCell",
"tissue": None,
"cancer": None,
"species": None,
},
"aucell": { # tested
"db": None, # required
"assay": None,
"aucMaxRank": None,
"normAUC": True,
"tissue": None,
"cancer": None,
"species": None,
},
"gsva": { # tested
"db": None, # required
"assay": None,
"kcdf": "Gaussian",
"minSize": 1,
"maxSize": None,
"tissue": None,
"cancer": None,
"species": None,
},
"singscore": { # tested
"db": None, # required
"assay": None,
"subSamples": None,
"centerScore": None,
"tissue": None,
"cancer": None,
"species": None,
},
"scmap": { # tested
"db": None, # required
"assay": None,
"cluster_col": "cell_type1",
"features": None,
"threshold": 0.5,
"use_cell_index": False,
},
"cheetah": { # tested
"db": None, # required
"assay": None,
"input_c": None,
"thresh": 0.1,
"n_genes": None,
"pc_thresh": None,
"only_pos": False,
"ref_ct": None,
"label": None,
},
"scclassify": { # test needed
"db": None, # required
"assay": None,
"algorithm": None,
"topN": 50,
"prob_threshold": None,
"parallel": None,
},
"scpred": { # test needed
"db": None, # required
"assay": None,
"pvar": "cell_type",
"model": "svmRadial",
"reduction": "pca",
"threshold": 0.55,
},
"azimuth": { # test needed
"ref": None, # required (or `db`)
"db": None,
"assay": None,
"anno_col_in": None,
"dims": None, # ignored with a warning
"k.anchor": None, # ignored with a warning
},
"scsa": { # tested
"db": None, # required
"foldchange": 2.0,
"pvalue": 0.05,
"python": config.lang.python,
},
"maca": { # tested
"db": None, # required
"python": config.lang.python,
"n_pcs": None,
"res": None,
"n_neis": None,
"freq": 0.5,
"use_weight": False,
},
"scmapnet": { # test needed
"db": None, # required
"python": config.lang.python,
"scmapnet_dir": None, # required
"weights": None, # required
"organ": None,
},
"mllmcelltype": { # need to pass local tests
"tissue": None, # required
"model": "gpt-5.5",
"api_key": None,
"top_gene_count": 10,
"base_urls": None,
"return_reasoning": False,
},
"lict": {
"species": "Human",
"tissue": None,
"topgenenumber": 30,
"validate": True,
"percent": 0.5,
"keys": None,
"provider": None,
},
"mapquery": { # tested
"db": None, # required
"use": None, # required
"ident_name": "predicted.id",
"refnorm": "auto",
"skip_if_normalized": True,
"ncores": 1,
"map_query": {},
"find_transfer_anchors": {},
"sctransform": {},
"normalize_data": {},
},
"scinfer": { # test needed
"db": None, # required (or `tissue`)
"tissue": None, # required (or `db`); fetched from GitHub
"subtype_present": False,
"subtype_info": False,
"assay": "RNA",
"slot": "counts",
"reduction": "umap",
"own_weightage": 0.5,
"n_neighbor": 10,
},
"add_prefix": None,
"merge": False,
"post_mutaters": {},
"cases": {},
"ncores": 1,
"outtype": "input",
}
script = "file://../scripts/scrna/CellTypeAnnotation.R"
class GarnettClassifierTrainer(Proc):DOCS
"""Train a [`Garnett`](https://cole-trapnell-lab.github.io/garnett/) classifier
from a Seurat object and marker genes.
The classifier is trained with `garnett::train_cell_classifier()` (garnett
>= 0.2.22, monocle3 branch) and can be fed to `CellTypeAnnotation` with
`tool = "garnett"` / `envs.garnett.classifier` for classification. Cell
types whose marker genes are absent from the expression data are dropped
with a warning, and cell types with too few cells
(`envs.min_observations`) are dropped with a message.
Input:
srtobj: The Seurat object in RDS/qs/qs2 format with raw counts.
markerfile: The marker file. Either a garnett-native marker file or a
universal marker table (see below). The detection is by content: a
file whose first non-comment line starts with `>` is garnett-native.
A garnett-native marker file is a text file with one block per
cell type, headed by `> <cell type>`, followed by rules such as
`expressed: <marker genes>` and `not expressed: <genes>`. See the
[garnett docs](https://cole-trapnell-lab.github.io/garnett/docs/#1b-train-your-own-classifier)
for the full grammar. It is used as is, so `subtype of:`
hierarchies and `expressed above:`/`below:`/`between:` rules are
supported. `envs.species`/`envs.cancer`/`envs.tissue` cannot be
applied to this format (it has no such columns) and will raise an
error if set.
A universal marker table is a TSV/CSV (or an RDS/qs/qs2 file
containing a data.frame) in long format with one row per gene per
cell type, with `cell_type` (required) and `gene` (required)
columns. It can also have an optional `direction` column
(`positive`/`negative`, aliases: `pos`/`neg`/`+`/`-`), where
negative markers become `not expressed:` rules of the converted
garnett marker file, and optional `species`/`cancer`/`tissue`
columns, used only when the matching env is set, which keeps the
rows with the given value (an error is raised if the table has no
such column or no rows match). Column aliases are auto-detected:
`celltype`/`cellType`/`Type` → `cell_type`,
`marker`/`Marker`/`gene_symbol` → `gene`, `sign` → `direction`,
and `tissueType` → `tissue`. A table without `cell_type`/`gene`
columns is an error. The converted garnett marker file is saved as
`<stem>.markers.txt` in the job output directory for inspection.
Output:
outfile: The classifier file (RDS with a `garnett_classifier` object).
Use it with `CellTypeAnnotation` (`tool = "garnett"`) as
`envs.garnett.classifier`.
Envs:
species (type=str): Filter the markers by the `species` column of a
universal marker table (see `markerfile` above). An error if the
marker file is garnett-native.
cancer (type=str): Filter the markers by the `cancer` column of a
universal marker table (see `markerfile` above). An error if the
marker file is garnett-native.
tissue (type=str): Filter the markers by the `tissue` column of a
universal marker table (see `markerfile` above). An error if the
marker file is garnett-native.
db (type=str): The name of the installed annotation package
(e.g. `org.Hs.eg.db`) to convert the gene IDs of the markers and
the expression data to a common ID space (default: `none`, i.e.
no conversion; then the gene IDs of the markers must match those
of the expression data, and the classifier is said to have
`custom` gene IDs — classify with the garnett tool of
`CellTypeAnnotation` using `db: none` and the default
`cds_gene_id_type: custom`).
cds_gene_id_type (choice): The gene ID type of the expression data
(default: `custom`). Supported types include `SYMBOL`, `ENSEMBL`,
`ENTREZID`.
marker_file_gene_id_type (choice): The gene ID type of the marker genes
(default: `custom`).
classifier_gene_id_type (choice): The gene ID type of the gene models
(default: `custom`).
assay (type=str): The assay of the Seurat object to use for training.
It must contain raw counts. If not set, the default assay is used.
min_observations (type=int): The minimum number of cells a cell type
must have in the training data for its model to be trained
(default: 8). Cell types with fewer cells are dropped with a
message.
max_training_samples (type=int): The maximum number of cells used to
train the model of each cell type (default: 500).
num_unknown (type=int): The maximum number of cells that do not match
any cell type to include in the training data as the `Unknown`
class (default: 500). Cells matching a cell type are all used for
training (up to `envs.max_training_samples`).
propogate_markers (flag): Whether to add the markers of parent cell
types to their subtypes (default: `True`). Note the spelling
`propogate` is garnett's.
lambdas (type=json): The regularization parameters for the
elastic-net models to consider (default: `None`, i.e. let
`cv.glmnet` choose). Usually no tuning is needed, but garnett's
default grid starts at `1e+05` and the cross-validation may select
that point, which gives an all-zero (degenerate) model
(see `envs.fail_on_degenerate`); e.g.
`10^seq(1, -3, length.out = 30)` avoids it.
cores (type=int): Number of CPU cores to use for training (default: 1).
seed (type=int): The random seed to make the training reproducible
(default: `None`, i.e. not set). Note that the `Unknown` class
cells and the training samples are drawn at random. A warning is
emitted when it is not set.
fail_on_degenerate (flag): Whether to fail, instead of merely warning,
when the fitted node models are all degenerate, i.e. have all-zero
coefficients (default: `False`). This happens when the
cross-validation selects the first point of garnett's hard-coded
lambda grid (`unique(c(1e+05, 50000, seq(10000, 100, by = -200),
...))`, fitted with `standardize = FALSE`), which yields a null
model predicting the class priors, so every cell ends up as
`Unknown` (see `envs.lambdas` for the remedy). The per-node
diagnostics are always reported and written to
`<outfile>.diagnostics.tsv`.
Requires:
r-seurat:
- check: {{proc.lang}} -e "library(Seurat)"
r-monocle3:
- check: {{proc.lang}} -e "library(monocle3)"
r-garnett:
- check: {{proc.lang}} -e "library(garnett)"
r-seuratwrappers:
- check: {{proc.lang}} -e "library(SeuratWrappers)"
""" # noqa: E501
input = "srtobj:file, markerfile:file"
output = "outfile:file:{{in.srtobj | stem0}}.classifier.RDS"
lang = config.lang.rscript
envs = {
"species": None,
"cancer": None,
"tissue": None,
"db": "none",
"cds_gene_id_type": "custom",
"marker_file_gene_id_type": "custom",
"classifier_gene_id_type": "custom",
"assay": None,
"min_observations": 8,
"max_training_samples": 500,
"num_unknown": 500,
"propogate_markers": True,
"lambdas": None,
"cores": 1,
"seed": 8525,
"fail_on_degenerate": False,
}
script = "file://../scripts/scrna/GarnettClassifierTrainer.R"
class HitypeWeightTrainer(Proc):DOCS
"""Train [`hitype`](https://pwwang.github.io/hitype/) marker weights from a
Seurat object and marker genes.
The weights are learned with `hitype::train_weights()` (hitype >= 0.0.6;
the Seurat object is used as is — `data` layer, no scaling) and the output
is a **weighted universal marker table** (`cell_type`/`gene`/`direction`/
`weight`/`level`, one row per marker gene) that can be fed to
`CellTypeAnnotation` with `tool = "hitype"` / `envs.hitype.db`, where the
weights are used as-is for scoring.
The markers to train are either given by `in.markerfile` or, when not
given, first discovered from the data with `hitype::find_markers()` (over
the cell types of `envs.ident` or the current Idents) and then trained —
both branches converge to the same weighted output. The cell types of the
markers must be a subset of the cell types in `envs.ident`/the Idents (an
error otherwise); marker genes absent from the expression data are
dropped by hitype with a warning.
Input:
srtobj: The Seurat object in RDS/qs/qs2 format.
markerfile: The marker file to train on, optional. When not
given, the markers are found from the data first (see
`envs.find_markers`). It is an optional input; note the trade-off
that marker content changes do not invalidate cached jobs (same
as `CellTypeAnnotation`'s `envs.<tool>.db`). The formats below
are auto-detected by `hitype::gs_prepare()` (hitype >= 0.0.6):
A universal marker table — a TSV/CSV or an RDS/qs/qs2 file of
a data.frame, long format, one row per gene per cell type, with
`cell_type` and `gene` columns (required) and optional
`direction`/`species`/`cancer`/`tissue`/`level` columns
(column aliases auto-detected: `celltype`/`cellType`/`Type` →
`cell_type`, `marker`/`Marker`/`gene_symbol` → `gene`,
`sign` → `direction`, `tissueType` → `tissue`). A `weight`
column, if any, is ignored — the weights are learned from the
data.
A native hitype/ScType db-format file — a TSV with
`cellName`/`geneSymbolmore1`/... columns, or a native ScType
xlsx — is used as is.
The cell types of the markers must be a subset of the cell
types in `envs.ident`/the Idents.
Output:
outfile: The trained marker weights (TSV, universal marker format
with a numeric `weight` column). Use it with `CellTypeAnnotation`
(`tool = "hitype"`) as `envs.hitype.db`.
Envs:
ident (type=str): The metadata column of the Seurat object whose
values are the cell types (the labels to learn from; default:
the current Idents).
assay (type=str): The assay to use. If not set, the default assay
is used.
level (type=int): The marker level to train/find (default: 1). For
markers with multiple `level`s, only `envs.level` is trained;
the output has a single `level`.
species (type=str): Filter the markers by the `species` column of a
universal marker table given by `in.markerfile`. An error if the
markers are a native db-format file (which has no such columns)
or discovered from the data.
cancer (type=str): Filter the markers by the `cancer` column of a
universal marker table given by `in.markerfile`. Same constraints
as `envs.species`.
universal marker table given by `in.markerfile`. Same constraints
as `envs.species`. Note the trained output carries no
`tissue`/`species`/`cancer` columns, so a downstream
`envs.hitype.tissue` on it (in `CellTypeAnnotation`) will error.
`envs.hitype.tissue` on it (in `CellTypeAnnotation`) will error.
find_markers (ns): Arguments for `hitype::find_markers()`, used when
`in.markerfile` is not given. Defaults mirror the function's.
`level` comes from `envs.level`.
- method (choice): `fc`/`seurat`/`presto` (default: `fc`). The
default `fc` needs no extra packages; `seurat` requires a
Seurat object (as given), `presto` requires the presto
package.
- top (type=int): Number of markers to return per cell type
(default: 20).
- min_log2fc (type=float): Minimum log2 fold change for a gene
to be kept as a marker (default: 0.25).
- min_pct (type=float): Minimum fraction of cells in the cell
type expressing the gene (default: 0.1).
- pos_only (flag): Only keep genes higher in the cell type than
in the rest of the cells (default: `True`).
- top (type=auto): Number of markers per direction; a scalar or
[pos, neg] (default: [10, 10]).
- against (type=json): Reference group for one-vs-group contrast:
null (all other cells), a list of cell types, or "nearest"
(most-correlated type; default: null).
- max_pct_out (type=float): Drop positive candidates expressed
in more than this fraction of all other cells (pan-lineage
guard; default: 0.75).
train_weights (ns): Arguments for `hitype::train_weights()`.
Defaults mirror the function's (including `seed`). `level` comes
from `envs.level`; `scaled` is forced `False` (Seurat input);
`clusters` comes from `envs.ident`/the Idents.
- method (choice): `glmnet`/`lr`/`rf`/`xgb`/`lrp`/`correlation`/
`uniform` (default: `glmnet`).
- data_split (type=json): Fractions for training/validation/
testing (default: `[0.7, 0.2, 0.1]`; with two elements, no
testing set is used).
- epochs (type=int): Number of epochs to train, `lrp` method
only (default: 20).
- batch_size (type=int): The batch size, `lrp` method only
(default: 32).
- run_weights_on_test (flag): Whether to run the weights on the
test set (default: `True`). Requires three `data_split`
elements.
- cv_folds (type=int): Number of cross-validation folds for
weight estimation; > 1 averages weights across folds for
stability (default: 1).
- pos_only (flag): Only keep positive weights (default: `False`).
- seed (type=int): The random seed for reproducibility
(default: 8525).
Requires:
r-seurat:
- check: {{proc.lang}} -e "library(Seurat)"
r-hitype:
- check: {{proc.lang}} -e "library(hitype)"
""" # noqa: E501
input = "srtobj:file, markerfile:file"
output = "outfile:file:{{in.srtobj | stem0}}.hitype.tsv"
lang = config.lang.rscript
envs = {
"ident": None,
"assay": None,
"level": 1,
"species": None,
"cancer": None,
"tissue": None,
"find_markers": {
"method": "fc",
"top": [10, 10],
"min_log2fc": 0.25,
"min_pct": 0.1,
"pos_only": True,
"against": None,
"max_pct_out": 0.75,
},
"train_weights": {
"method": "glmnet",
"data_split": [0.7, 0.2, 0.1],
"epochs": 20,
"batch_size": 32,
"run_weights_on_test": True,
"cv_folds": 1,
"pos_only": False,
"seed": 8525,
},
}
script = "file://../scripts/scrna/HitypeWeightTrainer.R"
@mark(DOCS
deprecated=(
"[{proc.name}] is deprecated, "
"use `CellTypeAnnotation` with `envs.tool = 'mapquery' instead."
)
)
class SeuratMap2Ref(Proc):
"""Map the seurat object to reference
See: <https://satijalab.org/seurat/articles/integration_mapping.html>
and <https://satijalab.org/seurat/articles/multimodal_reference_mapping.html>
Details:
### Preparing a Seurat reference for mapping
#### Step 0: Create the Seurat reference object
Start from raw counts.
```r
reference <- CreateSeuratObject(counts = reference_counts)
```
At this stage the object should contain at least:
* RNA counts
* cell barcodes
#### Step 1: Choose normalization strategy
Two main normalization strategies are supported for reference mapping.
##### Option A — LogNormalize workflow
Recommended for:
* standard scRNA-seq
* single modality datasets
* smaller references
```r
reference <- NormalizeData(reference, normalization.method = "LogNormalize", scale.factor = 10000)
reference <- FindVariableFeatures(reference, selection.method = "vst", nfeatures = 2000)
reference <- ScaleData(reference)
```
This produces the normalized expression matrix used for PCA.
##### Option B — SCTransform workflow
Recommended for:
* large datasets
* heterogeneous samples
* multimodal references (e.g. CITE-seq)
```r
reference <- SCTransform(reference, verbose = FALSE)
```
Important notes:
* Creates an **SCT assay**
* Variable features and scaling are performed automatically
* Mapping later requires `normalization.method = "SCT"`
#### Step 2: Choose dimensional reduction
The reference must contain a dimensional reduction used for mapping.
##### Option A — PCA (standard references)
Used with LogNormalize.
```r
reference <- RunPCA(reference, verbose = FALSE)
```
Typical usage:
* 30-50 PCs
##### Option B — SPCA (supervised PCA)
Used with SCTransform references and often in multimodal workflows.
```r
reference <- RunSPCA(reference, assay = "SCT")
```
SPCA learns a projection supervised by a cell‑cell similarity graph and is commonly used in reference atlases.
#### Step 3: Compute neighbors and clustering (optional but recommended)
Precomputing neighbors allows faster anchor finding.
##### PCA reference
```r
reference <- FindNeighbors(reference, reduction = "pca", dims = 1:30)
reference <- FindClusters(reference, resolution = 0.5)
```
##### SPCA / multimodal reference
```r
reference <- FindMultiModalNeighbors(
reference,
reduction.list = list("spca"),
dims.list = list(1:30)
)
```
This creates a **weighted nearest neighbor (WNN)** graph.
#### Step 4: Compute UMAP and store the model
To allow `MapQuery` to project new cells into the same UMAP space, the model must be saved.
```r
reference <- RunUMAP(
reference,
reduction = "pca",
dims = 1:30,
return.model = TRUE
)
```
For WNN references:
```r
reference <- RunUMAP(
reference,
nn.name = "weighted.nn",
reduction.name = "wnn.umap",
return.model = TRUE
)
```
Storing the model enables **ProjectUMAP / MapQuery** to reuse the trained embedding.
#### Step 5: Annotate the reference
Reference mapping transfers metadata labels.
Add cell type annotations to metadata:
```r
reference$celltype <- annotated_celltypes
reference$celltype_l1 <- broad_labels
reference$celltype_l2 <- fine_labels
```
Any metadata field can later be transferred.
#### Step 6: Save the reference
```r
saveRDS(reference, "reference.rds")
```
Or save it in [qs2](https://github.com/qsbase/qs2) format for faster loading:
```r
biopipen.utils::save_obj(reference, "reference.qs")
```
This allows the reference to be reused across multiple mapping runs.
### Prepare the query dataset (can be done with `SeuratPreparing` process)
The **query must use the same normalization method as the reference**.
If the query is not normalized in the same way as the reference, you can specify arguments
in `envs.NormalizeData` or `envs.SCTransform` for this process to reproduce the same normalization.
You can also specify `envs.skip_if_normalized = false` to force re‑normalization of the query dataset.
### Find transfer anchors (arguments specified in `envs.FindTransferAnchors`)
Anchors link cells between the query and reference.
#### LogNormalize reference
```r
anchors <- FindTransferAnchors(
reference = reference,
query = query,
reference.reduction = "pca",
dims = 1:30,
normalization.method = "LogNormalize"
)
```
#### SCTransform reference
```r
anchors <- FindTransferAnchors(
reference = reference,
query = query,
reference.reduction = "spca",
dims = 1:30,
normalization.method = "SCT"
)
```
Optional useful parameters:
* `reference.assay = "SCT"`
* `recompute.residuals = TRUE`
* `reference.neighbors = "pca.nn"` (reuse neighbor index if precomputed)
### Map the query (arguments specified in `envs.MapQuery`)
```r
query <- MapQuery(
anchorset = anchors,
query = query,
reference = reference,
refdata = list(celltype = "celltype"),
reference.reduction = "pca",
reduction.model = "umap"
)
```
`MapQuery` performs:
1. `TransferData`
2. `IntegrateEmbeddings`
3. `ProjectUMAP`
The query cells will:
* receive predicted cell type labels
* be projected into the reference UMAP
The reference UMAP reduction is reused and saved in the query object as `ref.umap`.
### Decision summary
| Reference type | Normalization | Reduction | Typical use |
| ----------------------- | ------------- | ---------- | --------------------------------- |
| Standard scRNA‑seq | LogNormalize | PCA | Single modality datasets |
| Large / atlas reference | SCTransform | SPCA | Large heterogeneous datasets |
| Multimodal reference | SCTransform | SPCA + WNN | CITE‑seq / multimodal integration |
---
### Important rules
1. Query and reference **must use the same normalization strategy**.
2. The **reference dimensional reduction must already exist**.
3. Metadata labels in the reference are required for label transfer.
4. Store a **UMAP model (`return.model = TRUE`)** to enable projection.
5. Precomputing neighbors improves performance for repeated mapping.
6. Use `ref.umap` for consistent visualization of mapped query cells.
### Practical advice
For high-quality reference atlases:
* integrate multiple datasets first
* curate annotations carefully
* use ~30-50 PCs
* keep metadata hierarchy (broad → fine labels)
The reference effectively acts as a **pretrained cell atlas** that new datasets can be projected onto.
Input:
sobjfile: The seurat object
Output:
outfile: The rds file of seurat object with cell type annotated.
Note that the reduction name will be `ref.umap` for the mapping.
To visualize the mapping, you should use `ref.umap` as the reduction name.
Envs:
ncores (type=int;order=-100): Number of cores to use.
When `split_by` is used, this will be the number of cores for each object to map to the reference.
When `split_by` is not used, this is used in `future::plan(strategy = "multicore", workers = <ncores>)`
to parallelize some Seurat procedures.
See also: <https://satijalab.org/seurat/archive/v3.0/future_vignette.html>
mutaters (type=json): The mutaters to mutate the metadata.
This is helpful when we want to create new columns for `split_by`.
See <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>.
use: A column name of metadata from the reference
(e.g. `celltype.l1`, `celltype.l2`) to transfer to the query as the
cell types (ident) for downstream analysis. This field is required.
If you want to transfer multiple columns, you can use
`envs.MapQuery.refdata`.
ident: The name of the ident for query transferred from `envs.use` of the reference.
ref: The reference seurat object file.
Either an RDS file or a h5seurat file that can be loaded by
`Seurat::LoadH5Seurat()`.
The file type is determined by the extension. `.rds` or `.RDS` for
RDS file, `.h5seurat` or `.h5` for h5seurat file.
refnorm (choice): Normalization method the reference used. The same method will be used for the query.
- LogNormalize: Using [`NormalizeData`](https://satijalab.org/seurat/reference/normalizedata).
- SCTransform: Using [`SCTransform`](https://satijalab.org/seurat/reference/sctransform).
- SCT: Alias of SCTransform.
- auto: Automatically detect the normalization method.
If the default assay of reference is `SCT`, then `SCTransform` will be used.
split_by: The column name in metadata to split the query into multiple objects.
This helps when the original query is too large to process.
skip_if_normalized: Skip normalization if the query is already normalized.
Since the object is supposed to be generated by `SeuratPreparing`, it is already normalized.
However, a different normalization method may be used.
If the reference is normalized by the same method as the query, the normalization can be skipped.
Otherwise, the normalization cannot be skipped.
The normalization method used for the query set is determined by the default assay.
If `SCT`, then `SCTransform` is used; otherwise, `NormalizeData` is used.
You can set this to `False` to force re-normalization (with or without the arguments previously used).
SCTransform (ns): Arguments for [`SCTransform()`](https://satijalab.org/seurat/reference/sctransform)
- do-correct-umi (flag): Place corrected UMI matrix in assay counts layer?
- do-scale (flag): Whether to scale residuals to have unit variance?
- do-center (flag): Whether to center residuals to have mean zero?
- <more>: See <https://satijalab.org/seurat/reference/sctransform>.
Note that the hyphen (`-`) will be transformed into `.` for the keys.
NormalizeData (ns): Arguments for [`NormalizeData()`](https://satijalab.org/seurat/reference/normalizedata)
- normalization-method: Normalization method.
- <more>: See <https://satijalab.org/seurat/reference/normalizedata>.
Note that the hyphen (`-`) will be transformed into `.` for the keys.
FindTransferAnchors (ns): Arguments for [`FindTransferAnchors()`](https://satijalab.org/seurat/reference/findtransferanchors)
- normalization-method (choice): Name of normalization method used.
- LogNormalize: Log-normalize the data matrix
- SCT: Scale data using the SCTransform method
- auto: Automatically detect the normalization method.
See `envs.refnorm`.
- reference-reduction: Name of dimensional reduction to use from the reference if running the pcaproject workflow.
Optionally enables reuse of precomputed reference dimensional reduction.
- <more>: See <https://satijalab.org/seurat/reference/findtransferanchors>.
Note that the hyphen (`-`) will be transformed into `.` for the keys.
MapQuery (ns): Arguments for [`MapQuery()`](https://satijalab.org/seurat/reference/mapquery)
- reference-reduction: Name of reduction to use from the reference for neighbor finding
- reduction-model: `DimReduc` object that contains the umap model.
- refdata (type=json): Extra data to transfer from the reference to the query.
- <more>: See <https://satijalab.org/seurat/reference/mapquery>.
Note that the hyphen (`-`) will be transformed into `.` for the keys.
cache (type=auto): Whether to cache the information at different steps.
If `True`, the seurat object will be cached in the job output directory, which will be not cleaned up when job is rerunning.
The cached seurat object will be saved as `<signature>.<kind>.RDS` file, where `<signature>` is the signature determined by
the input and envs of the process.
See <https://github.com/satijalab/seurat/issues/7849>, <https://github.com/satijalab/seurat/issues/5358> and
<https://github.com/satijalab/seurat/issues/6748> for more details also about reproducibility issues.
To not use the cached seurat object, you can either set `cache` to `False` or delete the cached file at
`<signature>.RDS` in the cache directory.
plots (type=json): The plots to generate.
The keys are the names of the plots and the values are the arguments for the plot.
The arguments will be passed to `biopipen.utils::VizSeuratMap2Ref()` to generate the plots.
The plots will be saved to the output directory.
See <https://pwwang.github.io/biopipen.utils.R/reference/VizSeuratMap2Ref.html>.
Requires:
r-seurat:
- check: {{proc.lang}} -e "library(Seurat)"
""" # noqa: E501
input = "sobjfile:file"
output = "outfile:file:{{in.sobjfile | stem}}.qs"
lang = config.lang.rscript
envs_depth = 3
envs = {
"ncores": config.misc.ncores,
"use": None,
"ident": "seurat_clusters",
"mutaters": {},
"ref": None,
"refnorm": "auto",
"split_by": None,
"skip_if_normalized": True,
"SCTransform": {
"do-correct-umi": False,
"do-scale": False,
"do-center": True,
},
"NormalizeData": {
"normalization-method": "LogNormalize",
},
"FindTransferAnchors": {
# "reference-reduction": "spca",
},
"MapQuery": {
# "reference-reduction": "spca",
# "reduction-model": "wnn.umap",
"refdata": {
# "celltype-l1": "celltype.l1",
# "celltype-l2": "celltype.l2",
# "predicted_ADT": "ADT",
},
},
"cache": config.path.tmpdir,
"plots": {
"Mapped Identity": {
"features": "{ident}:{use}",
},
"Mapping Score": {
"features": "{ident}.score",
},
},
}
script = "file://../scripts/scrna/SeuratMap2Ref.R"
plugin_opts = {"report": "file://../reports/common.svelte"}
class RadarPlots(Proc):DOCS
"""Radar plots for cell proportion in different clusters.
This process generates the radar plots for the clusters of T cells.
It explores the proportion of cells in different groups (e.g. Tumor vs Blood)
in different T-cell clusters.
Examples:
Let's say we have a metadata like this:
| Cell | Source | Timepoint | seurat_clusters |
| ---- | ------ | --------- | --------------- |
| A | Blood | Pre | 0 |
| B | Blood | Pre | 0 |
| C | Blood | Post | 1 |
| D | Blood | Post | 1 |
| E | Tumor | Pre | 2 |
| F | Tumor | Pre | 2 |
| G | Tumor | Post | 3 |
| H | Tumor | Post | 3 |
With configurations:
```toml
[RadarPlots.envs]
by = "Source"
```
Then we will have a radar plots like this:
We can use `each` to separate the cells into different cases:
```toml
[RadarPlots.envs]
by = "Source"
each = "Timepoint"
```
Then we will have two radar plots, one for `Pre` and one for `Post`:
Using `cluster_order` to change the order of the clusters and show only the first 3 clusters:
```toml
[RadarPlots.envs]
by = "Source"
cluster_order = ["2", "0", "1"]
breaks = [0, 50, 100] # also change the breaks
```
/// Attention
All the plots used in the examples are just for demonstration purpose. The real plots will have different appearance.
///
Input:
srtobj: The seurat object in RDS or qs/qs2 format
Output:
outdir: The output directory for the plots
Envs:
mutaters (type=json): Mutaters to mutate the metadata of the
seurat object. Keys are the column names and values are the
expressions to mutate the columns. These new columns will be
used to define your cases.
by: Which column to use to separate the cells in different groups.
`NA`s will be ignored. For example, If you have a column named `Source`
that marks the source of the cells, and you want to separate the cells
into `Tumor` and `Blood` groups, you can set `by` to `Source`.
The there will be two curves in the radar plot, one for `Tumor` and
one for `Blood`.
each: A column with values to separate all cells in different cases
When specified, the case will be expanded to multiple cases for
each value in the column.
If specified, `section` will be ignored, and the case name will
be used as the section name.
prefix_each (flag): Whether to prefix the `each` column name to the values as the
case/section name.
breakdown: An additional column with groups to break down the cells
distribution in each cluster. For example, if you want to see the
distribution of the cells in each cluster in different samples. In
this case, you should have multiple values in each `by`. These values
won't be plotted in the radar plot, but a barplot will be generated
with the mean value of each group and the error bar.
test (choice): The test to use to calculate the p values.
If there are more than 2 groups in `by`, the p values will be calculated
pairwise group by group. Only works when `breakdown` is specified and
`by` has 2 groups or more.
- wilcox: Wilcoxon rank sum test
- t: T test
- none: No test will be performed
order (list): The order of the values in `by`. You can also limit
(filter) the values we have in `by`. For example, if column `Source`
has values `Tumor`, `Blood`, `Spleen`, and you only want to plot
`Tumor` and `Blood`, you can set `order` to `["Tumor", "Blood"]`.
This will also have `Tumor` as the first item in the legend and `Blood`
as the second item.
colors: The colors for the groups in `by`. If not specified,
the default colors will be used.
Multiple colors can be separated by comma (`,`).
You can specify `biopipen` to use the `biopipen` palette.
ident: The column name of the cluster information.
cluster_order (list): The order of the clusters.
You may also use it to filter the clusters. If not given,
all clusters will be used.
If the cluster names are integers, use them directly for the order,
even though a prefix `Cluster` is added on the plot.
breaks (list;itype=int): breaks of the radar plots, from 0 to 100.
If not given, the breaks will be calculated automatically.
direction (choice): Direction to calculate the percentages.
- inter-cluster: the percentage of the cells in all groups
in each cluster (percentage adds up to 1 for each cluster).
- intra-cluster: the percentage of the cells in all clusters.
(percentage adds up to 1 for each group).
section: If you want to put multiple cases into a same section
in the report, you can set this option to the name of the section.
Only used in the report.
subset: The subset of the cells to do the analysis.
bar_devpars (ns): The parameters for `png()` for the barplot
- res (type=int): The resolution of the plot
- height (type=int): The height of the plot
- width (type=int): The width of the plot
devpars (ns): The parameters for `png()`
- res (type=int): The resolution of the plot
- height (type=int): The height of the plot
- width (type=int): The width of the plot
cases (type=json): The cases for the multiple radar plots.
Keys are the names of the cases and values are the arguments for
the plots (`each`, `by`, `order`, `breaks`, `direction`,
`ident`, `cluster_order` and `devpars`).
If not cases are given, a default case will be used, with the
key `DEFAULT`.
The keys must be valid string as part of the file name.
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}.radar_plots"
lang = config.lang.rscript
script = "file://../scripts/scrna/RadarPlots.R"
envs = {
"mutaters": {},
"by": None,
"each": None,
"prefix_each": True,
"order": None,
"colors": "biopipen",
"ident": "seurat_clusters",
"cluster_order": [],
"breakdown": None,
"test": "wilcox",
"breaks": [],
"direction": "intra-cluster",
"section": "DEFAULT",
"subset": None,
"bar_devpars": {
"res": 100,
"width": 1200,
"height": 800,
},
"devpars": {
"res": 100,
"width": 1200,
"height": 1000,
},
"cases": {},
}
plugin_opts = {
"report": "file://../reports/scrna/RadarPlots.svelte",
}
@mark(deprecated=True)DOCS
class MetaMarkers(Proc):
"""Find markers between three or more groups of cells, using one-way ANOVA
or Kruskal-Wallis test.
Sometimes, you may want to find the markers for cells from more than 2 groups.
In this case, you can use this process to find the markers for the groups and
do enrichment analysis for the markers. Each marker is examined using either
one-way ANOVA or Kruskal-Wallis test.
The p values are adjusted using the specified method. The significant markers
are then used for enrichment analysis using
[enrichr](https://maayanlab.cloud/Enrichr/) api.
Other than the markers and the enrichment analysis as outputs, this process also
generates violin plots for the top 10 markers.
Input:
srtobj: The seurat object loaded by `SeuratPreparing`
Output:
outdir: The output directory for the markers
Envs:
ncores (type=int): Number of cores to use to parallelize for genes
mutaters (type=json): The mutaters to mutate the metadata
The key-value pairs will be passed the `dplyr::mutate()` to mutate the metadata.
group-by: The column name in metadata to group the cells.
If only `group-by` is specified, and `idents` are
not specified, markers will be found for all groups in this column.
`NA` group will be ignored.
idents: The groups of cells to compare, values should be in the `group-by` column.
each: The column name in metadata to separate the cells into different cases.
prefix_each (flag): Whether to add the `each` value as prefix to the case name.
dbs (list): The dbs to do enrichment analysis for significant markers.
You can use built-in dbs in `enrichit`, or provide your own gmt files.
See also <https://pwwang.github.io/enrichit/reference/FetchGMT.html>.
The built-in dbs include:
* "BioCarta" or "BioCarta_2016"
* "GO_Biological_Process" or "GO_Biological_Process_2025"
* "GO_Cellular_Component" or "GO_Cellular_Component_2025"
* "GO_Molecular_Function" or "GO_Molecular_Function_2025"
* "KEGG", "KEGG_Human", "KEGG_2021", or "KEGG_2021_Human"
* "Hallmark", "MSigDB_Hallmark", or "MSigDB_Hallmark_2020"
* "Reactome", "Reactome_Pathways", or "Reactome_Pathways_2024"
* "WikiPathways", "WikiPathways_2024", "WikiPathways_Human", or "WikiPathways_2024_Human"
You can also fetch more dbs from <https://maayanlab.cloud/Enrichr/#libraries>.
subset: The subset of the cells to do the analysis.
An expression passed to `dplyr::filter()`.
p_adjust (choice): The method to adjust the p values, which can be used to filter the significant markers.
See also <https://rdrr.io/r/stats/p.adjust.html>
- holm: Holm-Bonferroni method
- hochberg: Hochberg method
- hommel: Hommel method
- bonferroni: Bonferroni method
- BH: Benjamini-Hochberg method
- BY: Benjamini-Yekutieli method
- fdr: FDR method of Benjamini-Hochberg
- none: No adjustment
sigmarkers: An expression passed to `dplyr::filter()` to filter the
significant markers for enrichment analysis. The default is `p.value < 0.05`.
If `method = 'anova'`, the variables that can be used for filtering are:
`sumsq`, `meansq`, `statistic`, `p.value` and `p_adjust`.
If `method = 'kruskal'`, the variables that can be used for filtering are:
`statistic`, `p.value` and `p_adjust`.
section: The section name for the report.
Worked only when `each` is not specified.
Otherwise, the section name will be constructed from `each` and `group-by`.
If `DEFAULT`, and it's the only section, it not included in the case/section names.
method (choice): The method for the test.
- anova: One-way ANOVA
- kruskal: Kruskal-Wallis test
cases (type=json): If you have multiple cases, you can specify them
here. The keys are the names of the cases and the values are the
above options except `ncores` and `mutaters`. If some options are
not specified, the default values specified above will be used.
If no cases are specified, the default case will be added with
the default values under `envs` with the name `DEFAULT`.
""" # noqa: E501
input = "srtobj:file"
output = "outdir:dir:{{in.srtobj | stem}}.meta_markers"
lang = config.lang.rscript
script = "file://../scripts/scrna/MetaMarkers.R"
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"group-by": None,
"idents": None,
"each": None,
"subset": None,
"prefix_each": True,
"p_adjust": "BH",
"dbs": ["KEGG_2021_Human", "MSigDB_Hallmark_2020"],
"sigmarkers": "p_adjust < 0.05",
"section": "DEFAULT",
"method": "anova",
"cases": {},
}
plugin_opts = {
"report": "file://../reports/scrna/MetaMarkers.svelte",
"report_paging": 8,
}
class Seurat2AnnData(Proc):DOCS
"""Convert seurat object to AnnData
Input:
sobjfile: The seurat object file, in RDS or qs/qs2 format
Output:
outfile: The AnnData file
Envs:
assay: The assay to use for AnnData.
If not specified, the default assay will be used.
"""
input = "sobjfile:file"
output = "outfile:file:{{in.sobjfile | stem}}.h5ad"
lang = config.lang.rscript
script = "file://../scripts/scrna/Seurat2AnnData.R"
envs = {"assay": None}
class AnnData2Seurat(Proc):DOCS
"""Convert AnnData to seurat object
Input:
adfile: The AnnData .h5ad file
Output:
outfile: The seurat object file in RDS or qs/qs2 format
Envs:
assay: The assay to use to convert to seurat object.
ncores (type=int): Number of cores to use for data loading and saving.
ident: The column name in `adfile.obs` to use as the identity
for the seurat object.
If not specified, no identity will be set.
dotplot_check (type=auto): Whether to do a check with a dot plot.
(`scplotter::FeatureStatPlot(plot_type = "dot", ..)` will be used)
to see if the conversion is successful.
Set to `False` to disable the check.
If `True`, top 10 variable genes will be used for the check.
You can give a list of genes or a string of genes with comma (`,`) separated
to use for the check.
"""
input = "adfile:file"
output = "outfile:file:{{in.adfile | stem}}.qs"
lang = config.lang.rscript
envs = {
"assay": "RNA",
"ident": None,
"dotplot_check": True,
"ncores": config.misc.ncores,
}
script = "file://../scripts/scrna/AnnData2Seurat.R"
class ScSimulation(Proc):DOCS
"""Simulate single-cell data using splatter.
See <https://www.bioconductor.org/packages/devel/bioc/vignettes/splatter/inst/doc/splatter.html#2_Quickstart>
Input:
seed: The seed for the simulation
You could also use string as the seed, and the seed will be
generated by `digest::digest2int()`.
So this could also work as a unique identifier for the simulation (ie. Sample ID).
Output:
outfile: The output Seurat object/SingleCellExperiment in qs/qs2 format
Envs:
ngenes (type=int): The number of genes to simulate
ncells (type=int): The number of cells to simulate
nspikes (type=int): The number of spike-ins to simulate
When `ngenes`, `ncells`, and `nspikes` are not specified, the default
params from `mockSCE()` will be used. By default, `ngenes = 2000`,
`ncells = 200`, and `nspikes = 100`.
outtype (choice): The output file type.
- seurat: Seurat object
- singlecellexperiment: SingleCellExperiment object
- sce: alias for `singlecellexperiment`
method (choice): which simulation method to use. Options are:
- single: produces a single population
- groups: produces distinct groups (eg. cell types), or
- paths: selects cells from continuous trajectories (eg. differentiation processes)
params (ns): Other parameters for simulation.
The parameters are initialized `splitEstimate(mockSCE())` and then
updated with the given parameters.
See <https://rdrr.io/bioc/splatter/man/SplatParams.html>.
Hyphens (`-`) will be transformed into dots (`.`) for the keys.
""" # noqa: E501
input = "seed:var"
output = "outfile:file:simulatied_{{in.seed}}.RDS"
lang = config.lang.rscript
envs = {
"ngenes": None,
"ncells": None,
"nspikes": None,
"outtype": "seurat",
"method": "single",
"params": {},
}
script = "file://../scripts/scrna/ScSimulation.R"
class CellCellCommunication(Proc):DOCS
"""Cell-cell communication inference
This is implemented based on [LIANA](https://liana-py.readthedocs.io/en/latest/index.html),
which is a Python package for cell-cell communication inference and provides a list of existing
methods including [CellPhoneDB](https://github.com/ventolab/CellphoneDB),
[Connectome](https://github.com/msraredon/Connectome/), log2FC,
[NATMI](https://github.com/forrest-lab/NATMI),
[SingleCellSignalR](https://github.com/SCA-IRCM/SingleCellSignalR), Rank_Aggregate, Geometric Mean,
[scSeqComm](https://gitlab.com/sysbiobig/scseqcomm), and [CellChat](https://github.com/jinworks/CellChat).
You can also try `python -c 'import liana; liana.mt.show_methods()'` to see the methods available.
Note that this process does not do any visualization. You can use `CellCellCommunicationPlots`
to visualize the results.
Reference:
- [Review](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9184522/).
- [LIANA](https://www.biorxiv.org/content/10.1101/2023.08.19.553863v1).
Input:
sobjfile: The seurat object file in RDS or AnnData file.
Output:
outfile: The output file with the 'liana_res' data frame.
Stats are provided for both ligand and receptor entities, more specifically: ligand and receptor are
the two entities that potentially interact. As a reminder, CCC events are not limited to secreted signalling,
but we refer to them as ligand and receptor for simplicity.
Also, in the case of heteromeric complexes, the ligand and receptor columns represent the subunit with minimum
expression, while *_complex corresponds to the actual complex, with subunits being separated by _.
source and target columns represent the source/sender and target/receiver cell identity for each interaction, respectively
* `*_props`: represents the proportion of cells that express the entity.
By default, any interactions in which either entity is not expressed in above 10%% of cells per cell type
is considered as a false positive, under the assumption that since CCC occurs between cell types, a sufficient
proportion of cells within should express the genes.
* `*_means`: entity expression mean per cell type.
* `lr_means`: mean ligand-receptor expression, as a measure of ligand-receptor interaction magnitude.
* `cellphone_pvals`: permutation-based p-values, as a measure of interaction specificity.
A typical output will look like this:
| ligand | ligand_complex | ligand_props | ligand_trimean | mat_max | receptor | receptor_complex | receptor_props | receptor_trimean | source | target | lr_probs | cellchat_pvals | mag_score | spec_score |
|--------|---------------|--------------|----------------|---------|----------|------------------|----------------|------------------|--------|--------|----------|----------------|-----------|------------|
| VIM | VIM | 1.00 | 0.36 | 8.73 | CD44 | CD44 | 0.77 | 0.16 | c7 | c3 | 0.10 | 0.00 | 0.10 | 0.00 |
| MIF | MIF | 0.97 | 0.22 | 8.73 | CXCR4 | CD74_CXCR4 | 0.87 | 0.26 | c5 | c6 | 0.10 | 0.00 | 0.10 | 0.00 |
| HLA-B | HLA-B | 1.00 | 0.44 | 8.73 | KLRD1 | KLRD1 | 0.73 | 0.13 | c9 | c2 | 0.10 | 0.00 | 0.10 | 0.00 |
| HMGB1 | HMGB1 | 0.99 | 0.26 | 8.73 | CXCR4 | CXCR4 | 0.81 | 0.21 | c2 | c7 | 0.10 | 0.00 | 0.10 | 0.00 |
| CD48 | CD48 | 0.94 | 0.20 | 8.73 | CD2 | CD2 | 0.99 | 0.28 | c7 | c8 | 0.10 | 0.00 | 0.10 | 0.00 |
| HLA-C | HLA-C | 1.00 | 0.38 | 8.73 | CD8B | CD8B | 0.73 | 0.15 | c1 | c9 | 0.10 | 0.00 | 0.10 | 0.00 |
| LGALS1 | LGALS1 | 0.95 | 0.17 | 8.73 | CD69 | CD69 | 0.99 | 0.34 | c10 | c5 | 0.10 | 0.00 | 0.10 | 0.00 |
Envs:
method (choice): The method to use for cell-cell communication inference.
- CellPhoneDB: Use CellPhoneDB method.
Magnitude Score: lr_means; Specificity Score: cellphone_pvals.
- Connectome: Use Connectome method.
- log2FC: Use log2FC method.
- NATMI: Use NATMI method.
- SingleCellSignalR: Use SingleCellSignalR method.
- Rank_Aggregate: Use Rank_Aggregate method.
- Geometric_Mean: Use Geometric Mean method.
- scSeqComm: Use scSeqComm method.
- CellChat: Use CellChat method.
- cellphonedb: alias for `CellPhoneDB`
- connectome: alias for `Connectome`
- log2fc: alias for `log2FC`
- natmi: alias for `NATMI`
- singlesignaler: alias for `SingleCellSignalR`
- rank_aggregate: alias for `Rank_Aggregate`
- geometric_mean: alias for `Geometric_Mean`
- scseqcomm: alias for `scSeqComm`
- cellchat: alias for `CellChat`
subset: An expression in string to subset the cells.
When a `.rds` file is provided for `in.sobjfile`, you can provide an expression in `R`,
which will be passed to `base::subset()` in `R` to subset the cells.
But you can always pass an expression in `python` to subset the cells.
See <https://anndata.readthedocs.io/en/latest/tutorials/notebooks/getting-started.html#subsetting-using-metadata>.
You should use `adata` to refer to the AnnData object. For example, `adata.obs.groups == "g1"` will subset the cells
with `groups` equal to `g1`.
subset_using (choice): The method to subset the cells.
- python: Use python to subset the cells.
The expression will be evaluated by `adata.obs.query()` to subset the cells, e.g.: `groups == "g1"`.
- r: Use R to subset the cells.
A subset expression will be passed to `dplyr::filter()` applying to the metadata of the Seurat object.
Note that this is only available when the input is a Seurat object in RDS format and
only available for `envs.subset_using`, but not for `cases.<case>.subset_using`.
- R: alias for `r`
split_by (auto): The column names in metadata to split the cells to run the method separately.
The results will be combined together with this column in the final output.
Multiple columns can be provided as a list and the data will be split by the combination of the columns
(the column values are concatenated only for the splitting).
In the final output, each column will be recovered with its original values.
assay: The assay to use for the analysis.
Only works for Seurat object.
seed (type=int): The seed for the random number generator.
ncores (type=int): The number of cores to use.
groupby: The column name in metadata to group the cells.
Typically, this column should be the cluster id.
If provided input is a Seurat object, the default identity will be used by default.
Otherwise, it is recommended to provide this parameter.
"seurat_clusters" will be used with a warning if the input is in AnnData format and
this parameter is not provided.
group_by: alias for `groupby`
species (choice): The species of the cells.
- human: Human cells, the 'consensus' resource will be used.
- mouse: Mouse cells, the 'mouseconsensus' resource will be used.
expr_prop (type=float): Minimum expression proportion for the ligands and
receptors (+ their subunits) in the corresponding cell identities. Set to 0
to return unfiltered results.
min_cells (type=int): Minimum cells (per cell identity if grouped by `groupby`)
to be considered for downstream analysis.
n_perms (type=int): Number of permutations for the permutation test.
Relevant only for permutation-based methods (e.g., `CellPhoneDB`).
If `0` is passed, no permutation testing is performed.
rscript: The path to the Rscript executable used to convert RDS file to AnnData.
if `in.sobjfile` is an RDS file, it will be converted to AnnData file (h5ad).
You need `Seurat`, `SeuratDisk` and `digest` installed.
cases (type=json): Different cases for the analysis.
The keys are the names of the cases, which will be saved in the `Case` column in the output file.
The values are the arguments for from `envs` to override the default values, except
`method`, `species` and `rscript`, to keep the output from different cases with the same columns.
Note that `subset` and `subset_using` will not be inherited from `envs`, and only `subset_using = "python"`
is supported for the cases. When using subset in the cases, the data is subsetted after the default subset in `envs` is applied.
If no cases are given, a default case will be used, with the arguments from `envs`, and
no `Case` column will be added in the output file.
<more>: Other arguments for the method.
The arguments are passed to the method directly.
See the method documentation for more details and also
`help(liana.mt.<method>.__call__)` in Python.
""" # noqa: E501
input = "sobjfile:file"
output = "outfile:file:{{in.sobjfile | stem}}-ccc.txt"
lang = config.lang.python
envs = {
"method": "cellchat",
"assay": None,
"seed": 1337,
"subset": None,
"subset_using": "auto",
"split_by": None,
"ncores": config.misc.ncores,
"groupby": None,
"group_by": None,
"species": "human",
"expr_prop": 0.1,
"min_cells": 5,
"n_perms": 1000,
"rscript": config.lang.rscript,
"cases": {},
}
script = "file://../scripts/scrna/CellCellCommunication.py"
class CellCellCommunicationPlots(Proc):DOCS
"""Visualization for cell-cell communication inference.
Examples:
### Network Plot
```toml
[CellCellCommunicationPlots.envs.cases."Cell-Cell Communication Network"]
plot_type = "network"
legend-position = "none"
theme = "theme_blank"
theme_args = {add_coord = false}
```
{: width="80%"}
### Circos Plot
{: width="80%"}
### Heatmap Plot
{: width="80%"}
### Cell-Cell Communication Interaction (Box Plot)
```toml
[CellCellCommunicationPlots.envs.cases."Cell-Cell Communication Interaction (Box Plot)"]
plot_type = "box"
x_text_angle = 90
method = "interaction"
```
{: width="80%"}
Input:
cccfile: The output file from `CellCellCommunication`
Output:
outdir: The output directory for the plots.
Envs:
subset: An expression to pass to `dplyr::filter()` to subset the ccc data.
magnitude: The column name in the data to use as the magnitude of the
communication. By default, the second last column will be used.
See `li.mt.show_methods()` for the available methods in LIANA. or
<https://liana-py.readthedocs.io/en/latest/notebooks/basic_usage.html#Tileplot>
specificity: The column name in the data to use as the specificity of the communication.
By default, the last column will be used. If the method doesn't have a specificity, set it to None.
ligand_expr: The name of the column containing the mean (or otherwise summarized) expression of the ligand.
Only used when plot_type = "linkedheatmap".
If not specified, `ligand_means` will be tried first, then `ligand_trimean`.
receptor_expr: The name of the column containing the mean (or otherwise summarized) expression of the receptor.
Only used when plot_type = "linkedheatmap".
If not specified, `receptor_means` will be tried first, then `receptor_trimean`.
devpars (ns): The parameters for the plot.
- res (type=int): The resolution of the plot
- height (type=int): The height of the plot
- width (type=int): The width of the plot
more_formats (type=list): The additional formats to save the plots.
descr: The description of the plot.
cases (type=json): The cases for the plots.
The keys are the names of the cases and the values are the arguments for
the plots. The arguments include the ones inherited from `envs`.
You can have a special `plot_type` `"table"` to generate a table for the
ccc data to save as a text file and show in the report.
If no cases are given, a default case will be used, with the
key `Cell-Cell Communication`.
<more>: Other arguments passed to
[scplotter::CCCPlot](https://pwwang.github.io/scplotter/reference/CCCPlot.html)
""" # noqa: E501
input = "cccfile:file"
output = "outdir:dir:{{in.cccfile | stem}}_plots"
lang = config.lang.rscript
envs = {
"subset": None,
"magnitude": None,
"specificity": None,
"ligand_expr": None,
"receptor_expr": None,
"devpars": {"res": 100},
"more_formats": [],
"descr": "Cell-cell communication plot",
"cases": {},
}
script = "file://../scripts/scrna/CellCellCommunicationPlots.R"
plugin_opts = {
"report": "file://../reports/common.svelte",
}
class ScVelo(Proc):DOCS
"""Velocity analysis for single-cell RNA-seq data
This process is implemented based on the Python package `scvelo` (v0.3.3).
Note that it doesn't work with `numpy>=2`.
Input:
sobjfile: The seurat object file in RDS or h5seurat format or AnnData file.
Output:
outfile: The output object with the velocity embeddings and information.
In either RDS, h5seurat or h5ad format, depending on the `envs.outtype`.
There will be also plots generated in the output directory
(parent directory of `outfile`).
Note that these plots will not be used in the report, but can be used as
supplementary information for the velocity analysis.
To visualize the velocity embeddings, you can use the `SeuratClusterStats`
process with `v_reduction` provided to one of the `envs.dimplots`.
Envs:
ncores (type=int): Number of cores to use.
group_by: The column name in metadata to group the cells.
Typically, this column should be the cluster id.
If provided input is a Seurat object, the default identity will be used by
default. Otherwise, it is recommended to provide this parameter.
"seurat_clusters" will be used with a warning if the input is in AnnData
format and this parameter is not provided.
mode (type=list): The mode to use for the velocity analysis.
It should be a subset of `['deterministic', 'stochastic', 'dynamical']`,
meaning that we can perform the velocity analysis in multiple modes.
fitting_by (choice): The mode to use for fitting the velocities.
- stochastic: Stochastic mode
- deterministic: Deterministic mode
min_shared_counts (type=int): Minimum number of counts
(both unspliced and spliced) required for a gene.
n_neighbors (type=int): The number of neighbors to use for the velocity graph.
n_pcs (type=int): The number of PCs to use for the velocity graph.
denoise (flag): Whether to denoise the data.
denoise_topn (type=int): Number of genes with highest likelihood selected to
infer velocity directions.
kinetics (flag): Whether to compute the RNA velocity kinetics.
kinetics_topn (type=int): Number of genes with highest likelihood selected to
infer velocity directions.
calculate_velocity_genes (flag): Whether to calculate the velocity genes.
top_n (type=int): The number of top features to plot.
rscript: The path to the Rscript executable used to convert RDS file to AnnData.
if `in.sobjfile` is an RDS file, it will be converted to AnnData file
(h5ad). You need `Seurat`, `SeuratDisk` and `digest` installed.
outtype (choice): The output file type.
- <input>: The same as the input file type.
- h5seurat: h5seurat file
- h5ad: h5ad file
- qs: qs/qs2 file
- qs2: qs2 file
- rds: RDS file
"""
input = "sobjfile:file"
output = (
"outfile:file:{{in.sobjfile | stem}}-scvelo."
"{{ext0(in.sobjfile) if envs.outtype == '<input>' else envs.outtype}}"
)
lang = config.lang.python
envs = {
"ncores": config.misc.ncores,
"group_by": None,
"mode": ["deterministic", "stochastic", "dynamical"],
"fitting_by": "stochastic",
"min_shared_counts": 30,
"n_neighbors": 30,
"n_pcs": 30,
"denoise": False,
"denoise_topn": 3,
"kinetics": False,
"kinetics_topn": 100,
"calculate_velocity_genes": False,
"top_n": 6,
"rscript": config.lang.rscript,
"outtype": "<input>",
}
script = "file://../scripts/scrna/ScVelo.py"
class Slingshot(Proc):DOCS
"""Trajectory inference using Slingshot
This process is implemented based on the R package `slingshot`.
Input:
sobjfile: The seurat object file in RDS or qs format.
Output:
outfile: The output object with the trajectory information.
The lineages are stored in the metadata of the seurat object at
columns `LineageX`, where X is the lineage number. The `BranchID`
column contains the branch id for each cell.
One can use
`scplotter::CellDimPlot(object, lineages = c("Lineage1", "Lineage2", ...))`
to visualize the trajectories.
Envs:
group_by: The column name in metadata to group the cells.
Typically, this column should be the cluster id.
Default is the default identity of the seurat object.
reduction: The nonlinear reduction to use for the trajectory analysis.
dims (type=auto): The dimensions to use for the analysis.
A list or a string with comma separated values.
Consecutive numbers can be specified with a colon (`:`) or a dash (`-`).
Or a single number greater than 1, which will be expanded to `1:number`.
If `None`, all dimensions will be used.
start: The starting group for the Slingshot analysis.
end: The ending group for the Slingshot analysis.
reverse (flag): Logical value indicating whether to reverse the pseudotime variable.
align_start (flag): Whether to align the starting pseudotime values at the maximum pseudotime.
seed (type=int): The seed for the random number generator.
ncores (type=int): The number of cores to use for reading and writing the seurat object.
subset: An expression in string to subset the cells.
split_by: The column name in metadata to split the cells to run the method separately.
After run, the results will be combined together with this column in the final output.
The lineages will be all stored in the columns `<prefix>_LineageX` and `<prefix>_BranchID`.
If a lineage is not found in a split, the corresponding columns will be filled with `NA`.
cases: A dictionary of cases to run the analysis.
The keys are the names of the cases, which will be served as the
prefix to add to the column names of the resulting pseudotime variable.
For example, if the case name is `case1`, the resulting pseudotime variable
will be stored in the column `case1_LineageX` and `case1_BranchID`.
The values are the arguments and will be inherited from the `envs` above, except for `cases`.
The default case will be added with the default values under `envs` with an empty prefix.
""" # noqa: E501
input = "sobjfile:file"
output = "outfile:file:{{in.sobjfile | stem}}.qs"
lang = config.lang.rscript
envs = {
"ncores": config.misc.ncores,
"group_by": None,
"reduction": None,
"dims": None,
"start": None,
"end": None,
"subset": None,
"split_by": None,
"cases": {},
"reverse": False,
"align_start": False,
"seed": 8525,
}
script = "file://../scripts/scrna/Slingshot.R"
class LoomTo10X(Proc):DOCS
"""Convert Loom file to 10X format
Input:
loomfile: The Loom file
Output:
outdir: The output directory for the 10X format files,
including the `matrix.mtx.gz`, `barcodes.tsv.gz` and `features.tsv.gz`
files.
"""
input = "loomfile:file"
output = "outdir:dir:{{in.loomfile | stem}}.10X"
lang = config.lang.rscript
script = "file://../scripts/scrna/LoomTo10X.R"
class PseudoBulkDEG(Proc):DOCS
"""Pseduo-bulk differential gene expression analysis
This process performs differential gene expression analysis, instead of
on single-cell level, on the pseudo-bulk data, aggregated from the single-cell data.
Input:
sobjfile: The seurat object file in RDS or qs/qs2 format.
Output:
outdir: The output containing the results of the differential gene expression
analysis.
Envs:
ncores (type=int): Number of cores to use for parallelization.
mutaters (type=json): Mutaters to mutate the metadata of the
seurat object. Keys are the new column names and values are the
expressions to mutate the columns. These new columns can be
used to define your cases.
You can also use the clone selectors to select the TCR clones/clusters.
See <https://pwwang.github.io/scplotter/reference/clone_selectors.html>.
You can also use key `<newcol>:ident` to set the `<newcol>` as the default
ident for the stats.
See <https://pwwang.github.io/biopipen.utils.R/reference/MutateSeuratMeta.html>
each: The column name in metadata to separate the cells into different cases.
When specified, the case will be expanded to multiple cases for
each value in the column.
`"ident"` can be used as an alias for the default identity column (e.g. `"seurat_clusters"`).
cache (type=auto): Where to cache the results.
If `True`, cache to `outdir` of the job. If `False`, don't cache.
Otherwise, specify the directory to cache to.
subset: An expression in string to subset the cells.
aggregate_by: The column names in metadata to aggregate the cells.
layer: The layer to pull and aggregate the data.
assay: The assay to pull and aggregate the data.
error (flag): Stop the job if errors happen.
Helpful when no/not enough markers are found or no pathways are enriched.
If `False`, empty results will be returned.
group_by: The column name in metadata to group the cells.
ident_1: The first identity to compare.
ident_2: The second identity to compare.
If not specified, the rest of the identities will be compared with `ident_1`.
paired_by: The column name in metadata to mark the paired samples.
For example, subject. If specified, the paired test will be performed.
dbs (list): The databases to use for enrichment analysis.
You can use built-in dbs in `enrichit`, or provide your own gmt files.
See also <https://pwwang.github.io/enrichit/reference/FetchGMT.html>.
The built-in dbs include:
* "BioCarta" or "BioCarta_2016"
* "GO_Biological_Process" or "GO_Biological_Process_2025"
* "GO_Cellular_Component" or "GO_Cellular_Component_2025"
* "GO_Molecular_Function" or "GO_Molecular_Function_2025"
* "KEGG", "KEGG_Human", "KEGG_2021", or "KEGG_2021_Human"
* "Hallmark", "MSigDB_Hallmark", or "MSigDB_Hallmark_2020"
* "Reactome", "Reactome_Pathways", or "Reactome_Pathways_2024"
* "WikiPathways", "WikiPathways_2024", "WikiPathways_Human", or "WikiPathways_2024_Human"
You can also fetch more dbs from <https://maayanlab.cloud/Enrichr/#libraries>.
sigmarkers: An expression passed to `dplyr::filter()` to filter the
significant markers for enrichment analysis.
The default is `p_val_adj < 0.05`.
If `tool = 'DESeq2'`, the variables that can be used for filtering
are: `baseMean`, `log2FC`, `lfcSE`, `stat`, `p_val`, `p_val_adj`.
If `tool = 'edgeR'`, the variables that can be used for filtering
are: `logCPM`, `log2FC`, `LR`, `p_val`, `p_val_adj`.
enrich_style (choice): The style of the enrichment analysis.
- enrichr: Use `enrichr`-style for the enrichment analysis.
- clusterProfiler: Use `clusterProfiler`-style for the enrichment analysis.
allmarker_plots_defaults (ns): Default options for the plots for all markers when `ident_1` is not specified.
- plot_type: The type of the plot.
See <https://pwwang.github.io/scplotter/reference/FeatureStatPlot.html>.
Available types are `violin`, `box`, `bar`, `ridge`, `dim`, `heatmap` and `dot`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- order_by: an expression to order the markers, passed by `dplyr::arrange()`.
- genes: The number of top genes to show or an expression passed to `dplyr::filter()` to filter the genes.
- <more>: Other arguments passed to [`scplotter::FeatureStatPlot()`](https://pwwang.github.io/scplotter/reference/FeatureStatPlot.html).
allmarker_plots (type=json): All marker plot cases.
The keys are the names of the cases and the values are the dicts inherited from `allmarker_plots_defaults`.
allenrich_plots_defaults (ns): Default options for the plots to generate for the enrichment analysis.
- plot_type: The type of the plot.
- descr: A description of the plot to be shown above the plot image.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- db (type=str): The database(s) to apply for the plot. If not specified, all databases will be applied.
`dbs` can also be used as an alias for this parameter.
- <more>: See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.html>.
allenrich_plots (type=json): Cases of the plots to generate for the enrichment analysis.
The keys are the names of the cases and the values are the dicts inherited from `allenrich_plots_defaults`.
The cases under `envs.cases` can inherit this options.
marker_plots_defaults (ns): Default options for the plots to generate for the markers.
- plot_type: The type of the plot.
See <https://pwwang.github.io/scplotter/reference/FeatureStatPlot.html>.
Available types are `violin`, `box`, `bar`, `ridge`, `dim`, `heatmap` and `dot`.
There are two additional types available - `volcano_pct` and `volcano_log2fc`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- order_by: an expression to order the markers, passed by `dplyr::arrange()`.
- genes: The number of top genes to show or an expression passed to `dplyr::filter()` to filter the genes.
- <more>: Other arguments passed to [`scplotter::FeatureStatPlot()`](https://pwwang.github.io/scplotter/reference/FeatureStatPlot.html).
If `plot_type` is `volcano_pct` or `volcano_log2fc`, they will be passed to
[`scplotter::VolcanoPlot()`](https://pwwang.github.io/plotthis/reference/VolcanoPlot.html).
marker_plots (type=json): Cases of the plots to generate for the markers.
Plot cases. The keys are the names of the cases and the values are the dicts inherited from `marker_plots_defaults`.
The cases under `envs.cases` can inherit this options.
enrich_plots_defaults (ns): Default options for the plots to generate for the enrichment analysis.
- plot_type: The type of the plot.
See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.html>.
Available types are `bar`, `dot`, `lollipop`, `network`, `enrichmap` and `wordcloud`.
- descr: A description of the plot to be shown above the plot image.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: See <https://pwwang.github.io/scplotter/reference/EnrichmentPlot.htmll>.
enrich_plots (type=json): Cases of the plots to generate for the enrichment analysis.
The keys are the names of the cases and the values are the dicts inherited from `enrich_plots_defaults`.
The cases under `envs.cases` can inherit this options.
overlaps_defaults (ns): Default options for investigating the overlapping of significant markers between different cases or comparisons.
This means either `ident_1` should be empty, so that they can be expanded to multiple comparisons.
- sigmarkers: The expression to filter the significant markers for each case.
If not provided, `envs.sigmarkers` will be used.
- descr: A description of the plot to be shown above the plot image.
- plot_type (choice): The type of the plot to generate for the overlaps.
- venn: Use `plotthis::VennDiagram()`.
- upset: Use `plotthis::UpsetPlot()`.
- more_formats (type=list): The extra formats to save the plot in.
- save_code (flag): Whether to save the code to generate the plot.
- devpars (ns): The device parameters for the plots.
- res (type=int): The resolution of the plots.
- height (type=int): The height of the plots.
- width (type=int): The width of the plots.
- <more>: More arguments pased to `plotthis::VennDiagram()`
(<https://pwwang.github.io/plotthis/reference/venndiagram1.html>)
or `plotthis::UpsetPlot()`
(<https://pwwang.github.io/plotthis/reference/upsetplot1.html>)
overlaps (type=json): Cases for investigating the overlapping of significant markers between different cases or comparisons.
The keys are the names of the cases and the values are the dicts inherited from `overlaps_defaults`.
There are two situations that we can perform overlaps:
1. If `ident_1` is not specified, the overlaps can be performed between different comparisons.
2. If `each` is specified, the overlaps can be performed between different cases, where in each case, `ident_1` must be specified.
tool (choice): The method to use for the differential expression analysis.
- DESeq2: Use DESeq2 for the analysis.
- edgeR: Use edgeR for the analysis.
plots_defaults (ns): The default parameters for the plots.
- <more>: Parameters passed to `biopipen.utils::VizBulkDEGs()`.
See: <https://pwwang.github.io/biopipen.utils.R/reference/VizBulkDEGs.html>
plots (type=json): The parameters for the plots.
The keys are the names of the plots and the values are the parameters
for the plots. The parameters will override the defaults in `plots_defaults`.
If not specified, no plots will be generated.
cases (type=json): The cases for the analysis.
The keys are the names of the cases and the values are the arguments for
the analysis. The arguments include the ones inherited from `envs`.
If no cases are specified, a default case will be added with
the name `DEG Analysis` and the default values specified above.
""" # noqa: E501
input = "sobjfile:file"
output = "outdir:dir:{{in.sobjfile | stem}}.pseudobulk_deg"
lang = config.lang.rscript
script = "file://../scripts/scrna/PseudoBulkDEG.R"
envs = {
"ncores": config.misc.ncores,
"mutaters": {},
"cache": config.path.tmpdir,
"each": None,
"subset": None,
"aggregate_by": None,
"layer": "counts",
"assay": None,
"error": True,
"group_by": None,
"ident_1": None,
"ident_2": None,
"paired_by": None,
"tool": "DESeq2",
"dbs": ["KEGG_2021_Human", "MSigDB_Hallmark_2020"],
"sigmarkers": "p_val_adj < 0.05",
"enrich_style": "enrichr",
"allmarker_plots_defaults": {
"plot_type": None,
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
"order_by": "desc(abs(log2FC))",
"genes": 10,
},
"allmarker_plots": {},
"allenrich_plots_defaults": {
"plot_type": "heatmap",
"devpars": {"res": 100},
},
"allenrich_plots": {},
"marker_plots_defaults": {
"plot_type": None,
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
"order_by": "desc(abs(log2FC))",
"genes": 10,
},
"marker_plots": {
"Volcano Plot": {"plot_type": "volcano"},
},
"enrich_plots_defaults": {
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"enrich_plots": {
"Bar Plot": {"plot_type": "bar", "ncol": 1, "top_term": 10},
},
"overlaps_defaults": {
"sigmarkers": None,
"plot_type": "venn",
"more_formats": [],
"save_code": False,
"devpars": {"res": 100},
},
"overlaps": {},
"cases": {},
}
plugin_opts = {
"report": "file://../reports/common.svelte",
"report_paging": 8,
}
class CellSNPLite(Proc):DOCS
"""Genotyping bi-allelic SNPs on single cells using cellsnp-lite.
The output from cellsnp-lite can be directly used for downstream analysis such as -
* Donor deconvolution in multiplexed single-cell RNA-seq data (e.g., with vireo).
* Allele-specific CNV analysis in single-cell or spatial transcriptomics data (e.g., with Numbat, XClone, or CalicoST).
* Clonal substructure discovery using single cell mitochondrial variants (e.g., with MQuad).
Here we only support model `1a`/`2a` in cellsnp-lite, which is designed for a single bam file as input.
For model `1b`/`2b`, which is designed for multiple bam files as input (e.g., one per cell), you can still
run with this process, but only one bam file is allowed.
See <https://github.com/single-cell-genetics/cellsnp-lite> for more details about cellsnp-lite.
Input:
crdir: The cellranger output directory or the directory containing
the bam file and barcode file.
It should contain the `outs/possorted_genome_bam.bam` file and
the `outs/filtered_feature_bc_matrix/barcodes.tsv.gz` file.
Output:
outdir: The output directory for cellsnp-lite results.
Envs:
ncores (type=int): The number of cores to use.
Will pass to `-p` option in cellsnp-lite.
regionsVCF: A vcf file listing all candidate SNPs, for fetch each variants.
genotype (flag): Whether to perform genotyping.
If `False`, only the allele counts will be computed.
gzip (flag): Whether to gzip the output files.
<more>: Other arguments passed to cellsnp-lite.
See <https://cellsnp-lite.readthedocs.io/en/latest/main/manual.html#full-parameters> for more details.
""" # noqa: E501
input = "crdir:dir"
output = """
outdir:dir:
{%- if basename(in.crdir) == 'outs' -%}
{{in.crdir | dirname | basename}}
{%- else -%}
{{in.crdir | basename}}
{%- endif -%}
.cellsnp
""" # noqa: E501
lang = config.lang.python
envs = {
"cellsnp_lite": config.exe.cellsnp_lite,
"ncores": config.misc.ncores,
"regionsVCF": None,
"genotype": False,
"gzip": True,
}
script = "file://../scripts/scrna/CellSNPLite.py"
class MQuad(Proc):DOCS
"""Clonal substructure discovery using single cell mitochondrial variants with MQuad.
MQuad uses a Mixture Model for Mitochondrial Mutation detection in single-cell omics data.
MQuad is a tool that detects mitochondrial mutations that are informative for clonal substructure inference. It uses a binomial mixture model to assess the heteroplasmy of mtDNA variants among background noise.
Input:
cellsnpout: The output directory from `CellSNPLite` process, which should contain
AD and DP sparse matrices (.mtx) or the vcf file.
Output:
outdir: The output directory for MQuad results.
Envs:
ncores (type=int): The number of cores to use.
It will be passed to `--nproc` option in MQuad.
seed (type=int): The seed for the random number generator.
It will be passed to `--randSeed` option in MQuad.
<more>: Other arguments passed to MQuad.
See <https://github.com/single-cell-genetics/MQuad/blob/main/mquad/mquad_CLI.py> for more details.
""" # noqa: E501
input = "cellsnpout:dir"
output = "outdir:dir:{{in.cellsnpout | stem}}.mquad"
lang = config.lang.python
envs = {
"mquad": config.exe.mquad,
"ncores": config.misc.ncores,
"seed": 8525,
}
script = "file://../scripts/scrna/MQuad.py"
plugin_opts = {
"report": "file://../reports/scrna/MQuad.svelte",
}
class MQuadVcf(Proc):DOCS
"""Clonal substructure discovery using single cell mitochondrial variants with MQuad.
MQuad uses a Mixture Model for Mitochondrial Mutation detection in single-cell omics data.
MQuad is a tool that detects mitochondrial mutations that are informative for clonal substructure inference. It uses a binomial mixture model to assess the heteroplasmy of mtDNA variants among background noise.
Input:
cellvcf: The output cell vcf from `CellSNPLite` process.
Output:
outdir: The output directory for MQuad results.
Envs:
ncores (type=int): The number of cores to use.
It will be passed to `--nproc` option in MQuad.
seed (type=int): The seed for the random number generator.
It will be passed to `--randSeed` option in MQuad.
<more>: Other arguments passed to MQuad.
See <https://github.com/single-cell-genetics/MQuad/blob/main/mquad/mquad_CLI.py> for more details.
""" # noqa: E501
input = "cellvcf:file"
output = "outdir:dir:{{in.cellvcf | dirname | stem}}.mquad"
lang = config.lang.python
envs = {
"mquad": config.exe.mquad,
"ncores": config.misc.ncores,
"seed": 8525,
}
script = "file://../scripts/scrna/MQuadVcf.py"
plugin_opts = {
"report": "file://../reports/scrna/MQuad.svelte",
}
class MQuadMerge(Proc):DOCS
"""Merge multiple MQuad results for multiple samples.
We will merge the passed_ad.mtx, passed_dp.mtx and passed_variant_names.txt files
from multiple samples
Input:
mquaddirs: The output directories from `MQuad` process for multiple samples.
Output:
outdir: The output directory for merged MQuad results.
This can be later used as input to `VireoSNP` process.
"""
input = "mquaddirs:dirs"
output = "outdir:dir:{{in.mquaddirs | first | stem | append: '-etc'}}.mquadmerged"
lang = config.lang.python
script = "file://../scripts/scrna/MQuadMerge.py"
class VireoSNP(Proc):DOCS
"""Demultiplexing of single-cell RNA-seq data using vireoSNP.
VireoSNP is a Bayesian method for demultiplexing pooled single-cell RNA-seq data
using natural genetic variations (SNPs) without requiring genotype reference.
Refers to <https://github.com/single-cell-genetics/vireo/blob/master/examples/vireoSNP_clones.ipynb> for more details.
Input:
cellsnpout: The output directory from `CellSNPLite` process, which should
contain AD and DP sparse matrices (.mtx) or the vcf file.
To investigate the clonal substructure using mitochondrial variants, run
`MQuad` first to select the informative variants, and then use the
filtered vcf file from MQuad as input to this process.
Output:
outdir: The output directory for vireoSNP results.
Envs:
seed (type=int): The seed for the random number generator.
ncores (type=int): The number of cores to use for model fitting for
different number of clones.
n_init (type=int): The number of random initializations to perform.
n_clones (type=auto): The number of clones in the pooled single-cell RNA-seq data.
(Refered as `n_donor` in vireoSNP documentation.)
If a 2-element list of provided, we will try to estimate the best number of clones
between the two values using the elbow plot method.
min_iter (type=int): The minimum number of iterations to perform.
max_iter (type=int): The maximum number of iterations to perform.
""" # noqa: E501
input = "cellsnpout:dir"
output = "outdir:dir:{{in.cellsnpout | stem}}.vireoSNP"
lang = config.lang.python
envs = {
"seed": 8525,
"ncores": config.misc.ncores,
"n_init": 50,
"n_clones": [2, 10],
"min_iter": 30,
"max_iter": 100,
}
script = "file://../scripts/scrna/VireoSNP.py"
plugin_opts = {
"report": "file://../reports/scrna/VireoSNP.svelte",
}