R语言实现基因组信息的筛选

2019-12-04 12:46:57 浏览数 (1)

今天给大家介绍一个R语言中的数据对象TxDb,此对象可以完美支持sqlite数据库导入,并且减少了检索的耗时,主要用来存储大量的基因信息数据。目前在R中存在大量数据存储的包,具体的框架及数据包如图:

首先我们看下这种数据的类型的构建,其需要用到一个包GenomicFeatures。此包可以操作sqlite数据库,并将其转化为TxDb数据对象。

此包的安装需要用到Biconductor库,我们就不细讲了,直接上安装代码:

代码语言:javascript复制
BiocManager::install("GenomicFeatures")BiocManager::install("TxDb.Hsapiens.UCSC.hg19.knownGene")BiocManager::install("BSgenome.Hsapiens.UCSC.hg19")

首先我们看数据的载入,如果是sqlite数据集可以直接用loadDb函数,我们直接看下实例:

代码语言:javascript复制
samplefile <-system.file("extdata", "hg19_knownGene_sample.sqlite",package="GenomicFeatures")txdb <- loadDb(samplefile)

至此TxDb数据就通过sqlite数据库文件构建好了。同时包还带了对一些数据库的直接构建TxDb数据对象的函数:makeTxDbFromUCSC,makeTxDbFromBiomart, makeTxDbFromGFF。那么构建好的数据集我们怎么用呢,接下来,我们看一个实例:TxDb.Hsapiens.UCSC.hg19.knownGene。

接下来我们就直接通过实际操作,数据筛选给大家看下如何去通过操作这个数据包找到我们想要的东西。

#载入数据

代码语言:javascript复制
library(TxDb.Hsapiens.UCSC.hg19.knownGene)txdb <-TxDb.Hsapiens.UCSC.hg19.knownGene 
代码语言:javascript复制
#数据信息描述seqinfo(txdb)
代码语言:javascript复制
 seqlevels(txdb)
代码语言:javascript复制
columns(txdb)
代码语言:javascript复制
keytypes(txdb)
代码语言:javascript复制
transcripts(txdb)#转录本信息 
代码语言:javascript复制
strand( transcripts(txdb))
代码语言:javascript复制
PR <- promoters(txdb, upstream=2000,downstream=400)#  启动子区
代码语言:javascript复制
exons(txdb)#外显子区
代码语言:javascript复制
cds(txdb)#编码区
代码语言:javascript复制
#抽取一个序列信息seqlevels(txdb) <- "chr15"
代码语言:javascript复制
#基于keys检索数据keys <- c("100033416","100033417", "100033420")select(txdb, keys = keys,columns="TXNAME", keytype="GENEID")
代码语言:javascript复制
#检索结果展示多列的数据cols <- c("TXNAME","TXSTRAND", "TXCHROM")select(txdb, keys=keys, columns=cols,keytype="GENEID")
代码语言:javascript复制
#检索符合要求的转录本信息GR <- transcripts(txdb,filter=list(tx_chrom = "chr15", tx_strand = " "))
代码语言:javascript复制
#基于基因的转录本分组GRList <- transcriptsBy(txdb, by ="gene")#另外还有exonsBy, and cdsBy, tx代表转录本transcript 。
代码语言:javascript复制
#获取所有的分组names(GRList)

当然也可以进行操作序列数据,那就需要导入序列的数据集BSgenome.Hsapiens.UCSC.hg19:

代码语言:javascript复制
library(BSgenome.Hsapiens.UCSC.hg19) #获取所有转录本区域的DNA序列tx_seqs1 <-extractTranscriptSeqs(Hsapiens, TxDb.Hsapiens.UCSC.hg19.knownGene,use.names=TRUE)
代码语言:javascript复制
#获取所有转录本区域的DNA蛋白质序列translate(tx_seqs1)

那这样两个包基本就整合到了基因的基本信息,那如何联合使用呢,我们看下面这个实例:

代码语言:javascript复制
cds_seqs <-extractTranscriptSeqs(Hsapiens,cdsBy(txdb, by="tx", use.names=TRUE))translate(cds_seqs)

两个包的完美组合指定能获得你想要的信息。

0 人点赞