R/utilities.R

In tangchao7498/ICAS:

#' PSI of SJ
#'
#' @title PSI of SJ
#' @name PSIofSJ
#'
#' @param object an ICASDataSet
#' @param SJ a charactor or vector of SJ
#' @param cells which cells are output, NULL for all cells
#'
#' @export

PSIofSJ <- function(object, SJ, cells = NULL) {
  if(!is(object, "ICASDataSet"))
    stop("The object must be a ICASDataSet data")

  if(any(!is.character(SJ)))
    stop("SJ must be a character formated chr:start-end:strand")

  if(any(!SJ %in% row.names(counts(object))))
    warning("some SJ are not in row.names(counts(object))")

  if(any(!SJ %in% row.names(psi(object))))
    warning("some SJ are constitutive SJ")

  SJ <- SJ[SJ %in% row.names(counts(object))]

  if(length(SJ) == 0) {
    stop("all SJ are not in row.names(counts(object))")
  }

  if(!is.null(cells)) {
    cells <- cells[cells %in% colnames(counts(object))]
    if(length(cells) == 0) {
      stop("all cells are not in colnames(counts(object))")
    }
  } else {
    cells <- colnames(counts(object))
  }

  lapply(SJ, function(sj) {
    if(sj %in% row.names(psi(object))) {
      psi(object)[sj, cells]
    } else {
      MinMax(counts(object)[sj, cells], 0, 1)
    }
  }) -> res
  res <- do.call(rbind, res)
  rownames(res) <- SJ
  return(res)
}

#' Count of SJ
#'
#' @title Count of SJ
#' @name CountSJ
#'
#' @param object an ICASDataSet
#' @param SJ a charactor or vector of SJ
#' @param cells which cells are output, NULL for all cells
#'
#' @export

CountSJ <- function(object, SJ, cells = NULL) {
  if(!is(object, "ICASDataSet"))
    stop("The object must be a ICASDataSet data")

  if(any(!is.character(SJ)))
    stop("SJ must be a character formated chr:start-end:strand")

  if(any(!SJ %in% row.names(counts(object))))
    warning("some SJ are not in row.names(counts(object))")

  SJ <- SJ[SJ %in% row.names(counts(object))]

  if(length(SJ) == 0) {
    stop("all SJ are not in row.names(counts(object))")
  }

  if(!is.null(cells)) {
    cells <- cells[cells %in% colnames(counts(object))]
    if(length(cells) == 0) {
      stop("all cells are not in colnames(counts(object))")
    }
  } else {
    cells <- colnames(counts(object))
  }

  lapply(SJ, function(sj) {
    counts(object)[sj, cells]
  }) -> res
  res <- do.call(rbind, res)
  rownames(res) <- SJ
  return(res)
}

#' SJEvent
#'
#' @title AS type and event ID of SJ
#' @name SJEvent
#'
#' @param object an ICASDataSet
#' @param SJ a charactor or vector of SJ
#' @param HostGene Output host gene of SJ or not
#'
#' @importFrom S4Vectors mcols
#' @export

SJEvent <- function(object, SJ, HostGene = TRUE) {
  if(!is(object, "ICASDataSet"))
    stop("The object must be a ICASDataSet data")

  if(is.null(object@ASType))
    stop("The AS type of ASSJ of object have not yet been identified")

  if(HostGene) {
    if(is.null(object@HostGene)) {
      HostGene <- FALSE
      warning("The Host Gene of SJ of object have not yet been identified")
    }
  }

  if(any(!is.character(SJ)))
    stop("SJ must be a character formated chr:start-end:strand")

  if(any(SJ %in% setdiff(row.names(counts(object)), row.names(psi(object))))) {
    warning("some SJ are constitutive SJ")
  }

  SJ <- SJ[SJ %in% row.names(psi(object))]
  if(length(SJ) == 0) {
    stop("all SJ are not ASSJ")
  }

  restype <- S4Vectors::mcols(object@ASType)[match(SJ, as.character(object@ASType)), ]

  if(HostGene) {
    reshostg <- S4Vectors::mcols(object@HostGene)[match(SJ, as.character(object@HostGene)), ]
    restype <- cbind(restype, reshostg)
  }
  row.names(restype) <- SJ
  return(restype)
}


#' SJHostGene
#'
#' @title Host gene of SJ
#' @name SJHostGene
#'
#' @param object an ICASDataSet
#' @param SJ a charactor or vector of SJ
#'
#' @importFrom S4Vectors mcols
#' @export

SJHostGene <- function(object, SJ) {
  if(!is(object, "ICASDataSet"))
    stop("The object must be a ICASDataSet data")

  if(is.null(object@HostGene))
    stop("The Host Gene of SJ of object have not yet been identified")

  if(any(!is.character(SJ)))
    stop("SJ must be a character formated chr:start-end:strand")

  SJ <- SJ[SJ %in% as.character(object@HostGene)]

  if(length(SJ) == 0) {
    stop("all SJ are not in object@HostGene")
  }

  reshostg <- S4Vectors::mcols(object@HostGene)[match(SJ, as.character(object@HostGene)), ]
  row.names(reshostg) <- SJ
  return(reshostg)
}


#' plotSJ
#'
#' @title plotSJ
#' @name plotSJ
#'
#' @param object an ICASDataSet
#' @param SJ a charactor or vector of SJ
#' @param curvature The curvature for splice junction
#'
#' @import ggplot2
#' @importFrom data.table as.data.table
#' @export

plotSJ <- function(object, SJ, curvature = -0.1) {
  if(!is(object, "ICASDataSet"))
    stop("The object must be a ICASDataSet data")

  if(length(SJ) != 1) {
    stop("SJ must be one SJ")
  }

  if(!SJ %in% row.names(counts(object)))
    stop("SJ are not in object")

  if(!SJ %in% row.names(psi(object))) {
    warning("SJ are ASSJ")
  }

  event <- SJEvent(object = object, SJ = SJ, HostGene = TRUE)

  Mat5 <- data.table::as.data.table(as(strsplit(event$ID, "@")[[1]], "GRanges"))
  Mat6 <- data.table::as.data.table(as(SJ, "GRanges"))

  ggplot() +
    geom_curve(data = Mat5, aes(x = start, xend = end, y = 1, yend = 1), curvature = curvature) +
    geom_curve(data = Mat6, aes(x = start, xend = end, y = 1, yend = 1), curvature = curvature, color = 2) +
    scale_x_continuous(expand = c(0, 0)) +
    theme_classic() +
    labs(x = "Genomic Coordinates", title = event$gene_name, subtitle = event$ID, caption = event$AS) +
    theme(axis.line.y = element_blank(),
          axis.text.y = element_blank(),
          axis.ticks.y = element_blank(),
          axis.title.y = element_blank(),
          axis.text.x = element_text(size = 12),
          axis.title.x = element_text(size = 16))
}

#' FetchMeta
#'
#' @title FetchMeta
#' @name FetchMeta
#'
#' @param object an ICASDataSet
#' @param variable a charactor of colnames in colData of object
#' @param cells which cells are output, NULL for all cells
#'
#' @export

FetchMeta <- function(object, variable, cells = NULL) {
  if(!is(object, "ICASDataSet"))
    stop("The object must be a ICASDataSet data")

  if(any(!is.character(variable)))
    stop("variable must be a character")

  if(any(!variable %in% colnames(colData(object))))
    stop("some variable are not in colnames(colData(object))")

  variable <- variable[variable %in% colnames(colData(object))]

  if(length(variable) != 1) {
    stop("all variable are not in colnames(colData(object))")
  }

  if(!is.null(cells)) {
    cells <- cells[cells %in% colnames(counts(object))]
    if(length(cells) == 0) {
      stop("all cells are not in colnames(counts(object))")
    }
  } else {
    cells <- colnames(counts(object))
  }

  res <- colData(object)[cells, variable]
  names(res) <- row.names(colData(object)[cells, ])
  return(res)
}