Perform composite test on results from mashr

The posterior probabilities for all genes and conditions is obtained as 1-lFSR. Let prob be an array storing results for one gene. The probability that _no_ conditions in the exclusion set are non-zero is prod(1 - prob[exclude]). The probability that _all_ conditions in the inclusion set are non-zero is prod(prob[include]). The probability that _at least one_ condition in the inclusion set is non-zero is 1 - prod(1 - prob[include]). The composite test is the product of the probabilties computed from the inclusion and exclusion sets.

Usage

compositePosteriorTest(
  x,
  include,
  exclude = NULL,
  test = c("at least 1", "all")
)

Arguments

x

"dreamlet_mash_result" from run_mash()

include

array of conditions in the inclusion set

exclude

array of conditions in the exclusion set. Defaults to NULL for no exclusion

test

evaluate the posterior probability of a non-zero effect in "at least 1" or "all" conditions

Details

Perform composite test evaluating the specificity of an effect. Evalute the posterior probability that an a non-zero effect present in _all_ or _at least one_ condition in the inclusion set, but _no conditions_ in the exclusion set.

See also

run_mash()

Examples

library(muscat)
library(mashr)
#> Loading required package: ashr
library(SingleCellExperiment)

data(example_sce)

# create pseudobulk for each sample and cell cluster
pb <- aggregateToPseudoBulk(example_sce[1:100, ],
  assay = "counts",
  cluster_id = "cluster_id",
  sample_id = "sample_id",
  verbose = FALSE
)

# voom-style normalization
res.proc <- processAssays(pb, ~group_id)
#>   B cells...
#> 0.048 secs
#>   CD14+ Monocytes...
#> 0.051 secs
#>   CD4 T cells...
#> 0.054 secs
#>   CD8 T cells...
#> 0.045 secs
#>   FCGR3A+ Monocytes...
#> 0.051 secs

# Differential expression analysis within each assay,
# evaluated on the voom normalized data
res.dl <- dreamlet(res.proc, ~group_id)
#>   B cells...
#> 0.031 secs
#>   CD14+ Monocytes...
#> 0.035 secs
#>   CD4 T cells...
#> 0.033 secs
#>   CD8 T cells...
#> 0.024 secs
#>   FCGR3A+ Monocytes...
#> 0.034 secs

# run MASH model
# This can take 10s of minutes on real data
# This small datasets should take ~30s
res_mash <- run_mash(res.dl, "group_idstim")

# Composite test based on posterior probabilities
# to identify effect present in *at least 1* monocyte type
# and *NO* T-cell type.
include <- c("CD14+ Monocytes", "FCGR3A+ Monocytes")
exclude <- c("CD4 T cells", "CD8 T cells")

# Perform composite test
prob <- compositePosteriorTest(res_mash, include, exclude)

# examine the lFSR for top gene
get_lfsr(res_mash$model)[which.max(prob), , drop = FALSE]
#>       B cells CD14+ Monocytes CD4 T cells CD8 T cells FCGR3A+ Monocytes
#> FBXO6      NA              NA          NA          NA      1.981564e-17

# Test if *all* cell types have non-zero effect
prob <- compositePosteriorTest(res_mash, assayNames(res.dl))