寻找基因间区域

Question

我想提取染色体的基因间坐标。我编写了大量代码，但由于我对这些包并不熟悉，所以我不确定是否遵循了这里的正确逻辑：

library(IRanges)
library(GenomicFeatures)
library(TxDb.Hsapiens.UCSC.hg19.knownGene)
txdb = transcriptsBy(TxDb.Hsapiens.UCSC.hg19.knownGene, by = "gene",use.names=TRUE)

#For example, only I am interested in intergenic coordinates in chromosome 1
seqlevels(txdb, force=TRUE) = c("chr1")

#Creates IRanges 

ir = IRanges(start=unlist(start(txdb)), end=unlist(end(txdb)), names=names(txdb)) 

# Reduce overlapping gene positions to continuous range
ir.red = reduce(ir) # Collapses ranges of overlapping genes to continuous range.

#Identify overlaps among the initial and reduced range data sets
overlap = findOverlaps(ir, ir.red) 

我必须照顾"+“和"-”吗？

Answer 1

最简单的方法是从gene()访问器开始，减少()，然后得到补语：

library(TxDb.Hsapiens.UCSC.hg19.knownGene)
genic <- genes(TxDb.Hsapiens.UCSC.hg19.knownGene)
genic <- reduce(genic, ignore.strand=T)
intergenic <- gaps(genic)
intergenic <- intergenic[strand(intergenic) == "*"] #This is important!!!

最后一步非常重要，否则每条染色体都会有另外两个条目(+和-各一个)。

Answer 2

虽然不能直接回答您的问题，但值得一提的是，您可以直接从UCSC浏览器http://genome.ucsc.edu/cgi-bin/hgTables获得这些数据。表浏览器允许您在两个轨道或它们的补充之间设置交叉，甚至可以将它们返回到表浏览器中的自定义轨道进行进一步操作，从而允许对多个轨道和表进行任意复杂的嵌套查询。为您的目的，步骤如下..。

1)选择你的基因组和兴趣的集合。

2)选择“基因和基因预测”作为分类，选择"UCSC基因“作为跟踪，选择"knownGene”作为表。

3)如果你想要单个染色体或染色体的一部分，在“区域”线上的方框中输入它，否则将其设置为“基因组”。

4)点击“交叉”线上的“创建”

5)在下一页中，选择group=“映射和测序”、track=“组装”、table=“组装(黄金)”，单击单选按钮“UCSC基因和汇编的基对交叉(和)”，选中“碱基对交叉/联合之前的互补UCSC基因”框，然后单击"submit“。

6)回到主Table屏幕上，输入输出选项所需的任何内容，然后按下"get output“以检索基因间区域。

狩猎愉快！