R/SER.R

In statwangz/SuSiE: A simple new approach to variable selection in regression, with application to genetic fine-mapping

#' The single effect regression model
#'
#' @param X SNPs, a n*p matrix.
#' @param y trait, a n*1 vector.
#' @param pi0 hyperparameter of the prior distribution of gamma.
#' @param sig2 hyperparameter, the variance of the trait y.
#' @param sig02 hyperparameter, the variance of the coefficient.
#'
#' @return a tibble consists of the result of the function BSLR on the SNPs matrix X.
#' @export
SER <- function(X, y, pi0 = NULL, sig2 = 1, sig02 = 0.1) {
  n <- length(y)
  p <- ncol(X)

  if (is.null(pi0)) {
    pi0 <- rep(1 / p, p)
  }

  bslr <- apply(X, MARGIN = 2, FUN = BSLR, y = y, sig2 = sig2, sig02 = sig02) %>%
    t() %>%
    as_tibble()

  log_bf_mdf <- bslr$log_bf - max(bslr$log_bf)
  bf_mdf <- exp(log_bf_mdf)
  alpha <- (pi0 * bf_mdf) / sum(pi0 * bf_mdf)

  ser <- tibble(alpha, mu1 = bslr$mu1, sig12 = bslr$sig12)

  return(ser)
}