R/plotCost.R

In Shicheng-Guo/edmr: Empirical differentially methylated regions calculation

#' plot cost function and the determined distance cutoff
#' @export
#' @param mixmdl object generated by \code{myDiff.to.mixmdl}
#' @param main title for the plot
#' @examples
#' library(edmr)
#' library(mixtools)
#' data(edmr)
#' 
#' myMixmdl=myDiff.to.mixmdl(chr22.myDiff)
#' plotCost(myMixmdl, main="cost function")
plotCost=function(mixmdl, ...){
  xlim=c(0,ceiling(max(mixmdl$x)))
  f = function(mixmdl,x){
    if(mixmdl$lambda[1]> mixmdl$lambda[2]) {i=1;k=2}
    else {i=2; k=1}
    p1=pnorm(x, mixmdl$mu[i],mixmdl$sigma[i], lower.tail = F)
    p2=pnorm(x, mixmdl$mu[k],mixmdl$sigma[k], lower.tail = T)
    cost=sum(mixmdl$lambda[i]*p1,mixmdl$lambda[k]*p2)
    return(cost)
  }
  x=seq(min(xlim), max(xlim), by=0.1)
  y=unlist(lapply(x, function(i){f(mixmdl,i)}))
  plot(x,y, type="n",xlab="x = log 2 distance", ylab="C(x)", ...)
  lines(x, y, col="blue", lwd=1.5)
  abline(v=get.break_point(mixmdl), col="red")
}