diff_dss_test: Calculates differential methylation statistics under general...

View source: R/diff_dss_test.R

diff_dss_testR Documentation

Calculates differential methylation statistics under general experimental design

Description

This function is a wrapper for DSS::DMLtest.multiFactor with the added feature of reporting methylation rates alongside the test results via the methylation_group_column and methylation_groups parameters. See documentation below.

Usage

diff_dss_test(
  bs,
  diff_fit,
  contrast,
  methylation_group_column = NA,
  methylation_groups = NA
)

Arguments

bs

a BSseq, the same used used to create diff_fit.

diff_fit

a list object output by diff_dss_fit().

contrast

a contrast matrix for hypothesis testing. The number of rows should match the number of columns design. Consult diff_fit$X to ensure the contrast correponds to the intended test.

methylation_group_column

Optionally, a column from diff_fit$design by which to group samples and capture methylation rates. This column can be a character, factor, or numeric. In the case of numeric the samples are grouped according to the top and bottom 25 percentiles of the covariate, and the mean methlyation for each group is calculated. If not a numeric, use the methylation_groups parameter to specify case and control.

methylation_groups

Optionally, a named character vector indicating the case and control factors of methylation_group_column by which to group samples and capture methylation rates. If specified, must also specify methylation_group_column.

Value

A GRanges object containing the following mcols:

stat:

The test statistic.

pvalue:

The p-value.

fdr:

The Benjamini-Hochberg adjusted p-values using p.adjust(method = 'BH').

If methylation_group_column is specified, also the following mcols:

meth_case:

Methylation estimate for case.

meth_control:

Methylation estimate for control.

meth_diff:

The difference meth_case - meth_control.

direction:

The group for which the locus is hyper-methylated. Note, this is not subject to significance thresholds.

Examples

data(BS.cancer.ex, package = 'bsseqData')

bs = filter_loci_by_group_coverage(
    bs = BS.cancer.ex,
    group_column = 'Type',
    c('cancer' = 2, 'normal' = 2))

small_test = bs[1:50]

diff_fit = diff_dss_fit(
    bs = small_test,
    design = bsseq::pData(bs),
    formula = '~ Type')

result = diff_dss_test(
    bs = small_test,
    diff_fit = diff_fit,
    contrast = matrix(c(0,1), ncol = 1)
)

result_with_meth = diff_dss_test(
    bs = small_test,
    diff_fit = diff_fit,
    contrast = matrix(c(0,1), ncol = 1),
    methylation_group_column = 'Type',
    methylation_groups = c('case' = 'cancer', 'control' = 'normal')
)


sartorlab/methylSig documentation built on Aug. 22, 2024, 5:58 p.m.