R/ChIP.promoter.R

#' ChIP.promoter
#'
#' Get a GRanges with enhancer-positions
#'
#' @param H3K4me3 A BED-df with H3k27ac-peak positions.
#' @param TSS A data-frame from biomart, containing a bed-like structure of TSS's
#' @param ignorePeak A BED-df with peak positions of an mark not associated with promoters (e.g. H3KK27ac). Can be a list.
#' @param TSSpad Search-space for TSS-H3K4me3
#' @param padding Extend peaks.
#' @examples promoters <- ChIP.promoter(H3K4me3 = H3K4me3,TSS = TSS , ignorePeak = list(H3K4me1, H3K27ac), TSSpad = 1e3)
#'enhancers <- ChIP.enhancer(ignorePeak = H3K4me3,H3K27ac = H3K27ac, activePeak = H3K4me1)
#' @return A genomicRanges-object
#' @export
ChIP.promoter <- function(H3K4me3, TSS = NULL, ignorePeak = NULL, TSSpad = 2e3, padding = 100, verbose = F){
  require(RHWlib)
  require(GenomicRanges)
  require(data.table)
  ME3 <- df2gr(H3K4me3, keepMeta = F)
  tss <- NULL
  ignore <- NULL

  promsFiltered <- NULL # me3 + tss + ignore
  proms <- NULL # me3 + tss
  me3Filtered <- NULL # me3 + ignore

  if(is.list(ignorePeak)){
    for(i in 1:length(ignorePeak)){
      ignorePeak[[i]] <- ignorePeak[[i]][,1:3]
    }
    ignorePeak <- rbindlist(ignorePeak)
    }



  # If no TSS and ignore, warn that the fuction will just output H3K4me3!
  if(is.null(TSS)){
    if(is.null(ignorePeak)){
      warning("No peaks to ignore and no TSS. Outputting the H3K4me3 peaks!")
    } else {
      ignore <- df2gr(ignorePeak, keepMeta = F)
      warning("No TSS given to denote promoters.")
    }
  } else { # active is given
    tss <- df2gr(TSS, keepMeta = F)
    if(is.null(ignorePeak)){
      warning("No peaks to ignore. Promoters could also potentially overlap enhancers!")
    } else {
      ignore <- df2gr(ignorePeak, keepMeta = F)
      if(verbose){message("Both TSS-sites and non-promoter peaks given.\nLet's get to work!")}
    }
  }


  # Get Promoters
  if(!is.null(tss)){ # Cool, intersect TSS with ME3

    PThits <- suppressWarnings(findOverlaps(ME3,tss, maxgap = TSSpad))
    proms <- ME3[PThits@from]


    if(!is.null(ignore)){ # ignore is also there!

      IPThits <-  suppressWarnings(findOverlaps(proms,ignore, maxgap = padding))
      promsFiltered <- proms[ - IPThits@from]

    }

  } else { # no tss....


    if(!is.null(ignore)){ # ignore is  there!

      MIhits <- suppressWarnings(findOverlaps(ME3,ignore, maxgap = padding))
      me3Filtered <- ME3[ - MIhits@from]

    }
  }

  # eerst proberen beide
  # dan degene met ignore
  # dan degene met tss
  # als laatste me3

  if(!is.null(promsFiltered)){
    return(unique(promsFiltered))
  } else if(!is.null(me3Filtered)){
    return(unique(me3Filtered))
  } else if(!is.null(proms)){
    return(unique(proms))
  } else {
    return(unique(ME3))
  }
  # I have to test this output... IGV!
}
robinweide/RHWlib documentation built on May 7, 2019, 8:03 a.m.