CallPeaks.paramEsti: Parameter estimation in m6A peak calling with multiple...

In haowulab/TRESS: Toolbox for mRNA epigenetics sequencing analysis

Parameter estimation in m6A peak calling with multiple replicates.

Description

This function estimates all involved parameters in Bayesian hierarchical negative binomial model, which is built for read counts from candidate regions generated from multiple input\& IP replicates.

Usage

CallPeaks.paramEsti(mat, sf = NULL, cutoff = NULL,
                    update = "Joint",
                    trans = NULL,
                    optM = "L-BFGS-B",
                    myfscale = -1e+06)

Arguments

mat	A matrix containing read counts from all paired input \& input replicates. The order of samples are: input1, IP1, input2, IP2,...
sf	A vector of size factors for each sample. It can be provided by the users or estimated automatically from the data. Default is NULL.
cutoff	Background methylation level, which can be automatically estimated based on the background read counts in IP and input samples, or provided by users. Defauls is NULL.
update	A logical value indicating whether jointly estimating the nuisance parameter theta with dispersion parameter phi listed in the proposed model. Possible options are "OnlyPhi", "Iterative" and "Joint". "OnlyPhi" means only updating phi_i using R function optimize while fixing parameter theta as the plug-in moment estimator; "Iterative" means iteratively updating and phi using R function optimize; "Joint" means updating them together using R function optim. Default is "Joint".
optM	A charactor value to specify which optimization algorithm used in the R function optim. The options are: "Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN" and "Brent". Default is "L-BFGS-B". See more details in help pages of optim.
trans	Needed when **optM == "Nelder-Mead"**. It specifies which transformation function used in the estimation of dispersion and/or theta parameter(s) which are subjected to the nonnegative constraints. Possible options are "sin()" and "exp()". Default is NULL.
myfscale	A stop criteria in optim. Default is -1e+06.

Details

This function mainly involves three estimation procedures:

• Estimate methylation levels

• Estimate dispersion parameters and the variance of the estimated methylation levels

• Calculate test statistics and p-values. Also, it calculates a score used for peak ranking.

Value

mu	Estimation of methylation levels of all peaks.
mu.var	Estimated variance for estimated methylation level.
shrkPhi	Shrinkage estimator for dispersion parameter phi_i.
shrkTheta	Shrinkage estimator for parameter theta_i if update == "Joint" or "Iterative". Otherwise it would be a plug-in moment estimator.
stats	Wald-test statisitcs.
pvals	P-values derived from normal distribution based on the Wald-test statisitcs.
p.adj	Adjusted p-values using Benjamini-Hochberg procedure.
rSocre	A score used to ranke each region. The higher the score, the higher the rank would be.

Examples

### A toy example using basal samples from mouse cortex
data("Basal")
res = CallPeaks.paramEsti(
    mat = as.matrix(Basal$Candidates$Counts),
    sf = Basal$Bins$sf,
    cutoff = 0.5
    )

haowulab/TRESS documentation built on Aug. 27, 2022, 7:11 p.m.