R/cutoff.R

In plbaldoni/mixNBHMM: Differential Peak Calling via Hidden Markov Model With Mixture of Negative Binomial Distributions

#' DEPRECATED, USE fdrControl INSTEAD!
#'
#' For a given FDR control threshold and a set of posterior probabilities, `cutoff()` classifies genomic windows pertaining to a particular HMM state.
#'
#' @param postprob numeric vector of posterior probabilities of genomic windows pertaining to a particular HMM state
#' @param alpha FDR threshold (default is 0.05)
#'
#' @return A logical vector indicating whether genomic windows pertain to a particular HMM state (\code{TRUE}) or not (\code{FALSE}), under a given FDR control threshold.
#'
#' @author Pedro L. Baldoni, \email{pedrobaldoni@gmail.com}
#' @references \url{https://biostats.bepress.com/jhubiostat/paper115/}
#' @references \url{https://doi.org/10.1093/biostatistics/5.2.155}
#'
#' @examples
#' data(ENCODE)
#' ENCODE <- createOffset(object = ENCODE,type = 'loess',span = 1)
#' \dontrun{ENCODE <- mixNBHMM(object = ENCODE,control = controlEM())}
#' \dontrun{cutoff(metadata(ENCODE)$prob$Differential,alpha = 0.05)}
#'
#' @import data.table
#'
#' @export
cutoff <- function(postprob,alpha = 0.05){
    if(any(postprob<0) | any(postprob>1)){stop('Posterior probabilities outside 0-1.')}
    notpp = FDR = Window = NULL

    DT <- data.table::data.table(notpp = 1-postprob,Window = 1:length(postprob),key = 'Window')
    return(DT[order(notpp),][,FDR := ((cumsum(notpp)/1:nrow(DT))<alpha)][order(Window),]$FDR)
}