getRegionNearGenes {ELMER}R Documentation

Identifies nearest genes to a region

Description

Auxiliary function for GetNearGenes This will get the closest genes (n=numFlankingGenes) for a target region (TRange) based on a genome of refenrece gene annotation (geneAnnot). If the transcript level annotation (tssAnnot) is provided the Distance will be updated to the distance to the nearest TSS.

Usage

getRegionNearGenes(
  TRange = NULL,
  numFlankingGenes = 20,
  geneAnnot = NULL,
  tssAnnot = NULL
)

Arguments

TRange

A GRange object contains coordinate of targets.

numFlankingGenes

A number determine how many gene will be collected from each

geneAnnot

A GRange object contains gene coordinates of for human genome.

tssAnnot

A GRange object contains tss coordinates of for human genome.

Value

A data frame of nearby genes and information: genes' IDs, genes' symbols,

Author(s)

Tiago C Silva (maintainer: tiagochst@usp.br)

Examples

geneAnnot <-  ELMER:::get.GRCh("hg38",as.granges = TRUE)
tssAnnot <-  getTSS(genome = "hg38")
probe <- GenomicRanges::GRanges(seqnames = c("chr1","chr2"), 
range=IRanges::IRanges(start = c(16058489,236417627), end= c(16058489,236417627)), 
name= c("cg18108049","cg17125141"))
names(probe) <- c("cg18108049","cg17125141")
NearbyGenes <- getRegionNearGenes(numFlankingGenes = 20,
                             geneAnnot = geneAnnot,
                             TRange = probe,
                             tssAnnot = tssAnnot)

[Package ELMER version 2.18.0 Index]