R/varBP.R

In nghiavtr/BPSC: Beta-Poisson Model for Single-Cell RNA-seq Data Analyses

#' Get variance of a beta-Poisson model
#'
#' @param alp A non-negative value, the parameter of Beta function (alpha)
#' @param bet A non-negative value, the parameter of Beta function (beta)
#' @param lam1 A non-negative value, the parameter for scaling
#' @param lam2 A non-negative value, the parameter for smoothing, used for BP4 or BP5
#' @param prob0 A non-negative value in [0;1), the parameter for the portion of fixed density of zero in the distribution, used for BP5
#' @return The variance of a beta-Poisson model
#' @export
#' @examples
#' varBP(alp=0.6,bet=1.5,lam1=20,lam2=0.05)
varBP = function(alp,bet,lam1=1,lam2=1,prob0=0){
  if (missing(bet)){
  	par=alp;
  	alp=par[1]; bet=par[2]; lam1=par[3]; lam2=1; if (length(par)> 3) lam2=par[4]; prob0=0;if (length(par)> 4) prob0=par[5]
  }

  v4=lam2^2*lam1*alp/((alp+bet)) + (lam2*lam1)^2 *alp*bet/((alp+bet)^2 *(alp+bet+1))  
  v5=NULL;
  res=v4;
  if (prob0!=0) {
      res=v5;
      cat("\nWarning: We do not compute the variance for BP5 model!!!")
  }  
  return(res)
}