ESCO_simulation: Simulate Datasets by ESCO

View source: R/21-ESCO.R

ESCO_simulationR Documentation

Simulate Datasets by ESCO

Description

This function is used to simulate datasets from learned parameters by escoEstimate function in ESCO package.

Usage

ESCO_simulation(
  parameters,
  return_format,
  other_prior = NULL,
  verbose = FALSE,
  seed
)

Arguments

parameters

A object generated by ESCO::escoEstimate()

return_format

A character. Alternative choices: list, SingleCellExperiment, Seurat, h5ad. If you select h5ad, you will get a path where the .h5ad file saves to.

other_prior

A list with names of certain parameters. Some methods need extra parameters to execute the estimation step, so you must input them. In simulation step, the number of cells, genes, groups, batches, the percent of DEGs are usually customed, so before simulating a dataset you must point it out. See Details below for more information.

verbose

Logical. Whether to return messages or not.

seed

A random seed.

Details

In addtion to simulate datasets with default parameters, users want to simulate other kinds of datasets, e.g. a counts matrix with 2 or more cell groups. In ESCO, you can set extra parameters to simulate datasets.

The customed parameters you can set are below:

  1. nCells. You can directly set other_prior = list(nCells = 1000) to simulate 1000 cells.

  2. nGenes. You can directly set other_prior = list(nGenes = 5000) to simulate 5000 genes.

  3. nGroups. You can not directly set other_prior = list(nGroups = 3) to simulate 3 groups. Instead, you should set other_prior = list(prob.group = c(0.2, 0.3, 0.5)) where the sum of group probabilities must equal to 1.

  4. de.prob. You can directly set other_prior = list(de.prob = 0.2) to simulate DEGs that account for 20 percent of all genes.

  5. prob.group. You can directly set other_prior = list(prob.group = c(0.2, 0.3, 0.5)) to assign three proportions of cell groups. Note that the number of groups always equals to the length of the vector.

  6. If users want to simulate tree or trajectory structured data, set other_prior = list(type = "tree") or other_prior = list(type = "traj"). See Examples.

For more customed parameters in ESCO, please check ESCO::escoParams().

References

Tian J, Wang J, Roeder K. ESCO: single cell expression simulation incorporating gene co-expression. Bioinformatics, 2021, 37(16): 2374-2381. https://doi.org/10.1093/bioinformatics/btab116

Github URL: https://github.com/JINJINT/ESCO

Examples

## Not run: 
## Estimation
ref_data <- simmethods::data

estimate_result <- simmethods::ESCO_estimation(ref_data = ref_data,
                                               other_prior = NULL,
                                               verbose = TRUE,
                                               seed = 111)
# 1) Simulate with default parameters
simulate_result <- simmethods::ESCO_simulation(
  parameters = estimate_result[["estimate_result"]],
  other_prior = NULL,
  return_format = "list",
  verbose = TRUE,
  seed = 111
)

## counts
counts <- simulate_result[["simulate_result"]][["count_data"]]
dim(counts)
## cell information
col_data <- simulate_result[["simulate_result"]][["col_meta"]]
table(col_data$group)


# 2) Simulate two groups (20% proportion of DEGs)
simulate_result <- simmethods::ESCO_simulation(
  parameters = estimate_result[["estimate_result"]],
  other_prior = list(nCells = 1000,
                     nGenes = 2000,
                     de.prob = 0.2,
                     prob.group = c(0.3, 0.7)),
  return_format = "list",
  verbose = TRUE,
  seed = 111
)

## counts
counts <- simulate_result[["simulate_result"]][["count_data"]]
dim(counts)
## cell information
col_data <- simulate_result[["simulate_result"]][["col_meta"]]
table(col_data$group)/1000
## gene information
row_data <- simulate_result[["simulate_result"]][["row_meta"]]
table(row_data$de_gene)[2]/2000
### The result of ESCO contains the factors of different groups and uses can
### calculate the fold change by division. For example, the DEFactors of A gene
### in Group1 and Group2 are respectively 2 and 1, and the fold change of A gene
### is 2/1=2 or 1/2=0.5.
fc_group1_to_group2 <- row_data$DEFacGroup2/row_data$DEFacGroup1
table(fc_group1_to_group2 != 1)[2]/2000 ## de.prob = 0.2

# ----------------- Simulate tree or trajectory structured data -------------
# Load data
ref_data <- simmethods::data
# Estimate parameters
estimate_result <- simmethods::ESCO_estimation(ref_data = ref_data,
                                               other_prior = list(tree = TRUE),
                                               verbose = TRUE,
                                               seed = 10)

# (1) Simulate tree structured cell groups
simulate_result <- simmethods::ESCO_simulation(parameters = estimate_result[["estimate_result"]],
                                               other_prior = list(batchCells = 1000,
                                                                  nGenes = 1000,
                                                                  type = "tree"),
                                               return_format = "SingleCellExperiment",
                                               verbose = TRUE,
                                               seed = 111)
## plot
result <- scater::logNormCounts(simulate_result[["simulate_result"]])
result <- scater::runPCA(result)
plotPCA(result, colour_by = "group")


# (2) Simulate tree structured cell groups (specify de.prob and group.prob)
simulate_result <- simmethods::ESCO_simulation(parameters = estimate_result[["estimate_result"]],
                                               other_prior = list(batchCells = 1000,
                                                                  nGenes = 1000,
                                                                  type = "tree",
                                                                  group.prob = c(0.4, 0.3, 0.3),
                                                                  de.prob = 0.5,
                                                                  de.center = 2),
                                               return_format = "SingleCellExperiment",
                                               verbose = TRUE,
                                               seed = 111)
## plot
result <- scater::logNormCounts(simulate_result[["simulate_result"]])
result <- scater::runPCA(result)
plotPCA(result, colour_by = "group")


# (3) Simulate continous cell trajectory
simulate_result <- simmethods::ESCO_simulation(parameters = estimate_result[["estimate_result"]],
                                               other_prior = list(batchCells = 1000,
                                                                  nGenes = 1000,
                                                                  type = "traj",
                                                                  group.prob = c(0.4, 0.3, 0.3),
                                                                  de.prob = 0.5,
                                                                  de.center = 2),
                                               return_format = "SingleCellExperiment",
                                               verbose = TRUE,
                                               seed = 111)
## plot
result <- scater::logNormCounts(simulate_result[["simulate_result"]])
result <- scater::runPCA(result)
plotPCA(result, colour_by = "group")

## End(Not run)


duohongrui/simmethods documentation built on June 17, 2024, 10:49 a.m.