R/add_cpg_islands.R

In riccardo-trozzo/BlissR: Package to analyse Bliss-seq data after BLISS-seq preprocessing pipeline

#' add cpg islands annotation
#'
#' This function allows you to add to the annotation table the cpg islands information of the gene
#' @param samples list of annotated peaks
#' @param genome hg38 or mm10, based on organism on sample table
#' @param cpg_annots the kind of annotation that you want, based also on the genome
#' @keywords cpg islands
#' @import data.table tidyverse annotatr GenomicRanges
#' @export
#' @examples
#' add_cpg_islands(samples, genome, cpg_annots)
add_cpg_islands <- function(samples, genome ="hg38", cpg_annots){
  #transform to GRanges object, which is needed by the package
  peaks_ranges <- lapply(samples, function(condition){
    t <- condition %>% as.data.frame() %>% distinct()
    t <- GRanges(t)
  })
  ###adding cpg islands information to annotation
  ###done using package annotatr
  cpg_annots <- build_annotations(genome = genome, annotations = cpg_annots)
  cpg_annotation <- list()
  for (condition in names(peaks_ranges)) {
    cpg_annotation[[condition]] <- annotate_regions(
      regions = peaks_ranges[[condition]],
      annotations = cpg_annots,
      ignore.strand = TRUE,
      quiet = FALSE
    )
    cpg_annotation[[condition]] <- as.data.frame(cpg_annotation[[condition]])
    names(cpg_annotation[[condition]])[names(cpg_annotation[[condition]]) == "seqnames"] <- "chrom"
    cpg_annotation[[condition]] <- cpg_annotation[[condition]] %>% select(chrom, start, end) %>% mutate(cpg_island = "yes")
    samples[[condition]] <- merge(samples[[condition]], cpg_annotation[[condition]], by = c("chrom", "start", "end"), all.x = TRUE)
  }
  rm(cpg_annots, cpg_annotation)
  return(samples)

}