R/prep.binary.lsn.mtx.R

In GRIN2: Genomic Random Interval (GRIN)

#' Prepare Binary Lesion Matrix
#'
#' @description
#' Constructs a binary lesion matrix representing the presence or absence of
#' specific lesion types affecting individual genes across patients. Each row
#' represents a gene-lesion type combination, and each column represents a
#' patient.
#'
#' @usage
#' prep.binary.lsn.mtx(ov.data,
#'                    min.ngrp = 0)
#'
#' @param ov.data GRIN gene-lesion overlap results, typically the output from
#' the `find.gene.lsn.overlaps` function.
#' @param min.ngrp Optional integer specifying the minimum number of patients
#' required in both the affected and unaffected groups for a gene-lesion type
#' combination to be retained. The default is `0`, which retains all observed
#' gene-lesion type combinations.
#'
#' @details
#' The function uses the gene-lesion overlaps in `ov.data$gene.lsn.hits` to
#' construct a binary matrix with gene-lesion type combinations as rows and
#' patients as columns.
#'
#' Each row is labeled using the format `gene_lesion.type` (for example,
#' `ENSG00000118513_gain`). For each gene-lesion type combination, a patient
#' receives a value of `1` if affected by that lesion type in the corresponding
#' gene and `0` otherwise.
#'
#' When `min.ngrp > 0`, a row is retained only when both the affected group
#' (`1`) and unaffected group (`0`) contain at least `min.ngrp` patients. This
#' can be useful when the resulting binary matrix is used for downstream
#' analyses that require a minimum number of patients in each comparison group.
#'
#' @return
#' A numeric binary matrix in which:
#' \itemize{
#'   \item Rows represent gene-lesion type combinations
#'   (`gene_lesion.type`).
#'   \item Columns represent patient IDs.
#'   \item Entries are `1` when the patient is affected by the specified lesion
#'   type in that gene and `0` otherwise.
#' }
#'
#' @export
#'
#' @references
#' Cao, X., Elsayed, A. H., & Pounds, S. B. (2023). Statistical Methods
#' Inspired by Challenges in Pediatric Cancer Multi-omics.
#'
#' @author
#' Abdelrahman Elsayed \email{abdelrahman.elsayed@stjude.org} and
#' Stanley Pounds \email{stanley.pounds@stjude.org}
#'
#' @seealso \code{\link{prep.gene.lsn.data}},
#' \code{\link{find.gene.lsn.overlaps}}
#'
#' @examples
#' data(lesion_data)
#' data(hg38_gene_annotation)
#'
#' # 1) Prepare gene-lesion input data
#' prep.gene.lsn <- prep.gene.lsn.data(lesion_data,
#'                                     hg38_gene_annotation)
#'
#' # 2) Identify gene-lesion overlaps
#' gene.lsn.overlap <- find.gene.lsn.overlaps(prep.gene.lsn)
#'
#' # 3) Create a binary lesion matrix requiring at least 5 patients
#' # in both the affected and unaffected groups
#' lsn.binary.mtx <- prep.binary.lsn.mtx(gene.lsn.overlap,
#'                                       min.ngrp = 5)
#'
prep.binary.lsn.mtx=function(ov.data,     # output of the find.gene.lsn.overlaps function
                             min.ngrp=0)  # minimum number of patients required in both affected and unaffected groups
{
  gene.lsn=ov.data$gene.lsn.hits

  # Order and index data by gene and lesion type

  ord=order(gene.lsn$gene,gene.lsn$lsn.type)
  gene.lsn=gene.lsn[ord,]
  m=nrow(gene.lsn)

  new.gene.lsn=which((gene.lsn$gene[-1]!=gene.lsn$gene[-m])|
                       (gene.lsn$lsn.type[-1]!=gene.lsn$lsn.type[-m]))

  row.start=c(1,new.gene.lsn+1)
  row.end=c(new.gene.lsn,m)

  gene.lsn.index=cbind.data.frame(gene=gene.lsn$gene[row.start],
                                  lsn.type=gene.lsn$lsn.type[row.start],
                                  row.start=row.start,
                                  row.end=row.end)

  k=nrow(gene.lsn.index)
  uniq.ID=unique(gene.lsn$ID)
  n=length(uniq.ID)

  gene.lsn.mtx=matrix(0,k,n)
  colnames(gene.lsn.mtx)=as.character(uniq.ID)
  rownames(gene.lsn.mtx)=paste0(gene.lsn.index$gene,"_",
                                gene.lsn.index$lsn.type)

  for (i in 1:k)
  {
    rows=(gene.lsn.index$row.start[i]:gene.lsn.index$row.end[i])
    ids=as.character(gene.lsn$ID[rows])
    gene.lsn.mtx[i,ids]=1
  }

  n.hit=rowSums(gene.lsn.mtx)
  min.n=pmin(n.hit,n-n.hit)

  keep.row=which(min.n>=min.ngrp)
  gene.lsn.mtx=matrix(gene.lsn.mtx[keep.row,],
                      length(keep.row),n)

  colnames(gene.lsn.mtx)=as.character(uniq.ID)
  rownames(gene.lsn.mtx)=paste0(gene.lsn.index$gene,"_",
                                gene.lsn.index$lsn.type)[keep.row]

  return(gene.lsn.mtx)
}