find_differentiated_clusters: Find Differentiated Clusters

View source: R/IDclust.R

find_differentiated_clustersR Documentation

Find Differentiated Clusters

Description

Find Differentiated Clusters

Usage

find_differentiated_clusters(object, ...)

## Default S3 method:
find_differentiated_clusters(
  object,
  differential_function = differential_ChromSCape,
  by = "IDcluster",
  logFC.th = log2(1.5),
  qval.th = 0.01,
  min_frac_cell_assigned = 0.1,
  limit = 5,
  FP_linear_model = NULL,
  cluster_of_origin = "Alpha",
  min_cluster_size = 30,
  verbose = TRUE,
  ...
)

## S3 method for class 'Seurat'
find_differentiated_clusters(
  object,
  differential_function = differential_edgeR_pseudobulk_LRT,
  by = "IDcluster",
  logFC.th = log2(1.5),
  qval.th = 0.01,
  limit = 5,
  cluster_of_origin = "Alpha",
  min_frac_cell_assigned = 0.1,
  min_cluster_size = 30,
  verbose = TRUE,
  ...
)

Arguments

object

A Seurat object containing scRNA dataset with 'IDcluster' column.

...

Additional parameters passed to the differential_function. See differential_ChromSCape() for more information on additional parameters for the default function.

differential_function

A function that take in entry a SingleCellExperiment object and parameters passed in ... and returns a data.frame containing the significantly differential features for each cluster. See differential_edgeR_pseudobulk_LRT for the default function.

by

A character specifying the name of the metadata column referencing the clusters.

logFC.th

A numeric specifying the log2 fold change of activation above/below which a feature is considered as significantly differential passed to the differential_function.

qval.th

A numeric specifying the adjusted p-value below which a feature is considered as significantly differential passed to the differential_function.

min_frac_cell_assigned

A numeric between 0 and 1 specifying the minimum percentage of the total cells in the SingleCellExperiment object that needs to be assigned. If a lower proportion is assigned, all cells are assigned to the cluster of origin.

limit

An integer specifying the minimum number of features required for a subcluster to be called 'true' subcluster.

FP_linear_model

Optional. A linear model (see stats::lm()) of the number of false positive expected for a given cluster size. The lm_list list of linear models present in this package gives default values accross multiple binsizes. (See calculate_FDR_scEpigenomics).

cluster_of_origin

A character specifying the name of the cluster of origin that will be concatenated before the name of true subclusters.

min_cluster_size

An integer specifying the minimum number of cells in a cluster to consider it as a 'true' subcluster.

verbose

A logical specifying wether to print.

Details

Find significantly differential features between the given set of clusters (within the 'IDcluster' column of the SingleCellExperiment). For each cluster, if enough differences are found, mark the cluster as a 'true' subcluster and gives it the alias 'cluster_of_origin_cluster'. The function will use by default ChromSCape::differential_activation() function to define differential features.

Value

A list containing :

  • "diffmat_n" - A data.frame containing the number of differential regions foun per cluster and the new assignations of the subclusters.

  • "res" - A data.frame containing the differential analysis.

  • "passing_min_pct_cell_assigned" - A boolean indicating if enough cells were assigned

See Also

differential_ChromSCape()

Examples

if(requireNamespace("Seurat", quietly=TRUE)){
# Find differentiated clusters in Seurat object using 
# edgeR_pseudobulk_LRT function
data("Seu")

DA = find_differentiated_clusters(
    Seu,
    differential_function = differential_edgeR_pseudobulk_LRT,
    logFC.th = log2(1.5),
    qval.th = 0.01,
    by = "seurat_clusters",
    limit = 5,
    cluster_of_origin = "Alpha",
    min_frac_cell_assigned = 0.1,
    verbose = TRUE
)

# Summary of differential genes per cluster
head(DA$diffmat_n)

# Differential analysis
head(DA$res)

# Did the clustering pass the minimum percent of cell assigned threshold ?
print(DA$passing_min_pct_cell_assigned)

# Find differentiated clusters in Seurat object using Seurat function
data("Seu")

DA = find_differentiated_clusters(
    Seu,
    differential_function = differential_Seurat,
    logFC.th = log2(1.5),
    qval.th = 0.01,
    by = "seurat_clusters",
    limit = 5,
    cluster_of_origin = "Alpha",
    min_frac_cell_assigned = 0.1,
    verbose = TRUE
)

# Find differentiated clusters in Seurat object using Seurat function, 
# passing additional arguments to differential_Seurat and thus 
# Seurat::FindAllMarkers funtion.
data("Seu")

DA = find_differentiated_clusters(
    Seu,
    differential_function = differential_Seurat,
    logFC.th = log2(1.5),
    qval.th = 0.01,
    by = "seurat_clusters",
    limit = 5,
    cluster_of_origin = "Alpha",
    min_frac_cell_assigned = 0.1,
    verbose = TRUE,
    test.use = "roc" # additional argument
)
}
# Find differentiated clusters in SingleCellExperiment object using 
# differential_ChromSCape function.
if(requireNamespace("ChromSCape", quietly=TRUE)){
data("scExp")
scExp = ChromSCape::find_clusters_louvain_scExp(scExp, 
resolution = 0.1)
DA = find_differentiated_clusters(
    scExp,
    differential_function = differential_ChromSCape,
    logFC.th = log2(5),
    qval.th = 0.01,
    by = "IDcluster",
    limit = 5,
    cluster_of_origin = "Alpha",
    min_frac_cell_assigned = 0.1,
    verbose = TRUE,
)
}

vallotlab/IDclust documentation built on July 5, 2024, 3:26 p.m.