ChIPseekerでゲノムの位置情報にアノテーションをつける

ゲノムの位置情報があったときに、さくっと近傍の遺伝子などを知りたいとする。
bedtoolsとアノテーションファイルを使えばわかるかもしれないが、それよりも全自動で注釈をさくってつけてくれた方が楽でしょう。

そんなときにChIPseekerが使えるぽい。

ChIPseekerは文字通り、Chip-seq用のツールですが、アノテーションをつけたい場合は汎用的に使えそうです。

## loading packages
library(ChIPseeker)

# Gene
library(TxDb.Hsapiens.UCSC.hg38.knownGene)
txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene

chr    <- c("chrX", "chrX", "chrX")
start  <- c(71349, 77124, 86364)
end    <- c(72835, 77836, 88378)
strand <- c("-", "+", "*")

# Genomic Range
peak <- GRanges(seqnames = chr, ranges = IRanges(start, end), strand = strand)

# Annotation
peakAnno <- annotatePeak(peak,
                         tssRegion=c(-3000,3000),
                         TxDb = txdb,
                         annoDb="org.Hs.eg.db")

result <- as.GRanges(peakAnno)

こんな結果が得られる。

seqnames	start	end	width	strand	annotation	geneChr	geneStart	geneEnd	geneLength	geneStrand	geneId	transcriptId	distanceToTSS	ENSEMBL	SYMBOL	GENENAME
1	chrX	71349	72835	1487	-	Distal Intergenic	23	276322	303353	27032	1	55344	ENST00000399012.6	-204973	ENSG00000182378	PLCXD1
2	chrX	77124	77836	713	+	Distal Intergenic	23	276322	303353	27032	1	55344	ENST00000399012.6	-198486	ENSG00000182378	PLCXD1
3	chrX	86364	88378	2015	*	Distal Intergenic	23	276322	303353	27032	1	55344	ENST00000399012.6	-187944	ENSG00000182378	PLCXD1

他にもいろいろ便利な関数がついているようです。

plotAnnoPie(peakAnno)

今回の場合だと、全部intergenicであまり意味のあるグラフになりませんが…

また同じ作者の clusterProfiler を使ってenrichment解析にもっていけるようになっているようです。

この記事は以上です。