diff_methylsig {methylSig} | R Documentation |
The function calculates differential methylation statistics between two groups of samples using a beta-binomial approach to calculate differential methylation statistics, accounting for variation among samples within each group. The function can be applied to a BSseq
object subjected to filter_loci_by_coverage()
, filter_loci_by_snps()
, filter_loci_by_group_coverage()
or any combination thereof. Moreover, the function can be applied to a BSseq
object which has been tiled with tile_by_regions()
or tile_by_windows()
.
diff_methylsig( bs, group_column, comparison_groups, disp_groups, local_window_size = 0, local_weight_function, t_approx = TRUE, n_cores = 1 )
bs |
a |
group_column |
a |
comparison_groups |
a named |
disp_groups |
a named |
local_window_size |
an |
local_weight_function |
a weight kernel function. The default is the tri-weight kernel function defined as |
t_approx |
a |
n_cores |
an |
A GRanges
object containing the following mcols
:
Methylation estimate for case.
Methylation estimate for control.
The difference meth_case - meth_control
.
The group for which the locus is hyper-methylated. Note, this is not subject to significance thresholds.
The p-value from the t-test (t_approx = TRUE
) or the Chi-Square test (t_approx = FALSE
).
The Benjamini-Hochberg adjusted p-values using p.adjust(method = 'BH')
.
The dispersion estimate.
The log likelihood ratio.
Degrees of freedom used when t_approx = TRUE
.
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[seq(50)] diff_gr = diff_methylsig( bs = small_test, group_column = 'Type', comparison_groups = c('case' = 'cancer', 'control' = 'normal'), disp_groups = c('case' = TRUE, 'control' = TRUE), local_window_size = 0, t_approx = TRUE, n_cores = 1)