R/helper_functions.R

In cbg-ethz/timeseriesNEM: Map an Observational Time Series onto a Static NEM Signalling Network

#' Generate a random NEM.
#'
#' @param n The number of S-nodes
#' @param m The number of E-genes
#' @return A `nem` object.
#' @import nem
sample_random_NEM_and_data <- function(n = 10, m = 1000) {
  Phi <- sampleRndNetwork(paste0("S", 1:n), DAG = TRUE)
  Data <- sampleData(Phi, m, type = "density")
  rownames(Data$D) <- paste0("E", 1:m)
  hyper <- set.default.parameters(colnames(Data$D))
  hyper$type <- "CONTmLLMAP"
  nem.res <- nem(Data$D, "search", list(Phi), hyper)
  nem.res$pos <- nem.res$pos[[1]]
  nem.res$mappos <- nem.res$mappos[[1]]
  nem.res$mLL <- nem.res$mLL[[1]]
  nem.res$LLperGene <- nem.res$LLperGene[[1]]
  nem.res$para <- nem.res$para[[1]]
  nem.res$control <- hyper
  nem.res
}

#' Generate an effect profile from a static NEM.
#'
#' @param nemObject A `nem` object as generated by `sample_random_NEM_and_data`.
#' @param perturb The name of an S-node to perturb.
#' @import nem
generate_tsNEM_dataset <- function(nemObject, perturb = NULL) {
  # perturb <- which(rowSums(transitive.reduction(nemObject$graph)) != 0) %>%
  #   sample(size = 1) %>%
  #   names()
  sim.int <- sim.intervention(x = nemObject, int = perturb)
  affectedSgenes <- sim.int[["Sgenes.effected"]][["intervention 1"]]
  affectedEgenes <- sim.int[["Egenes.effected"]][["intervention 1"]]
  alpha <- sample(seq(0.1, 0.9, by = 0.1), size = 1)
  beta <- sample(5:50, size = 1)
  lambda <- sample(seq(0.01, 0.49, by = 0.01), size = 2, replace = TRUE)
  lambda <- c(1 - sum(lambda), lambda)
  m <- nrow(nemObject$pos)
  palt <- nem:::bum.ralt(m, c(alpha, beta), lambda)
  p0 <- runif(m)
  P <- ifelse(nemObject$selected %in% affectedEgenes, palt, p0)
  D <- matrix(log(nem:::bum.dalt(P, c(alpha, beta), lambda)), ncol = 1,
              dimnames = list(nemObject$selected, "V1"))
  list(nemObject = nemObject, logDensities = D, affectedSgenes = affectedSgenes)
}