fit_models: Fits a model for each gene in a cell_data_set object.

View source: R/expr_models.R

fit_modelsR Documentation

Fits a model for each gene in a cell_data_set object.

Description

This function fits a generalized linear model for each gene in a cell_data_set. Formulae can be provided to account for additional covariates (e.g. day collected, genotype of cells, media conditions, etc).

Usage

fit_models(
  cds,
  model_formula_str,
  expression_family = "quasipoisson",
  reduction_method = "UMAP",
  cores = 1,
  clean_model = TRUE,
  verbose = FALSE,
  ...
)

Arguments

cds

The cell_data_set upon which to perform this operation.

model_formula_str

A formula string specifying the model to fit for the genes.

expression_family

Specifies the family function used for expression responses. Can be one of "quasipoisson", "negbinomial", "poisson", "binomial", "gaussian", "zipoisson", "zinegbinomial", or "mixed-negbinomial". Default is "quasipoisson".

reduction_method

Which method to use with clusters() and partitions(). Default is "UMAP".

cores

The number of processor cores to use during fitting.

clean_model

Logical indicating whether to clean the model. Default is TRUE.

verbose

Logical indicating whether to emit progress messages.

...

Additional arguments passed to model fitting functions.

Value

a tibble where the rows are genes and columns are

  • id character vector from rowData(cds)$id

  • gene_short_names character vector from rowData(cds)$gene_short_names

  • num_cells_expressed int vector from rowData(cds)$num_cells_expressed

  • gene_id character vector from row.names(rowData(cds))'

  • model GLM model list returned by speedglm

  • model_summary model summary list returned by summary(model)

  • status character vector of model fitting status: OK when model converged, otherwise FAIL

Examples

  
    cell_metadata <- readRDS(system.file('extdata',
                                         'worm_embryo/worm_embryo_coldata.rds',
                                         package='monocle3'))
    gene_metadata <- readRDS(system.file('extdata',
                                         'worm_embryo/worm_embryo_rowdata.rds',
                                         package='monocle3'))
    expression_matrix <- readRDS(system.file('extdata',
                                             'worm_embryo/worm_embryo_expression_matrix.rds',
                                             package='monocle3'))
   
    cds <- new_cell_data_set(expression_data=expression_matrix,
                             cell_metadata=cell_metadata,
                             gene_metadata=gene_metadata)

    cds <- preprocess_cds(cds, num_dim=50)
    cds <- align_cds(cds, alignment_group = "batch",
                     residual_model_formula_str = "~ bg.300.loading + bg.400.loading +
                     bg.500.1.loading + bg.500.2.loading + bg.r17.loading + bg.b01.loading +
                     bg.b02.loading")
    cds <- reduce_dimension(cds)
    ciliated_genes <- c("che-1", "hlh-17", "nhr-6", "dmd-6", "ceh-36", "ham-1")
    cds_subset <- cds[rowData(cds)$gene_short_name %in% ciliated_genes,]
    gene_fits <- fit_models(cds_subset, model_formula_str = "~embryo.time")
  


cole-trapnell-lab/monocle3 documentation built on May 24, 2022, 5:25 p.m.