基于Biostrings的基因組數(shù)據(jù)包Bioconductor項(xiàng)目提供了包含給定生物體全基因組序列的數(shù)據(jù)包。這些軟件包稱為基于Biostrings的基因組數(shù)據(jù)包,因?yàn)樗鼈儼男蛄写鎯?chǔ)在Biostrings軟件包中定義的一些基本容器中,如DNAString,DNAStringSet或MaskedDNAString容器。不管它們包含的特定序列數(shù)據(jù)如何,所有的基于生物基因組的基因組數(shù)據(jù)包非常相似,可以以一致且簡(jiǎn)單的方式進(jìn)行操作。為了正常工作,他們都需要使用BSgenome軟件包。不同于基于Biostrings的基因組數(shù)據(jù)包,這個(gè)軟件包與是一個(gè)提供支持它們所需基礎(chǔ)的軟件包(這就是為什么基于Biostrings的基因組數(shù)據(jù)包也稱為BSgenome數(shù)據(jù)包)。 BSgenome軟件包本身需要Biostrings軟件包。 安裝BSgenome包(如果沒(méi)有安裝)source('https:///biocLite.R')
biocLite('BSgenome.Hsapiens.UCSC.hg19')
載入BSgenome包,并查看當(dāng)前版本提供的BSgenome數(shù)據(jù)包:library(BSgenome)
(ag <- available.genomes())
unique(gsub('BSgenome\\.([^\\.]+).*', '\\1', ag))
通過(guò)運(yùn)行結(jié)果可得知,當(dāng)前版本提供了24個(gè)物種基因組包 [1] 'Alyrata' 'Amellifera' 'Athaliana'
[4] 'Btaurus' 'Celegans' 'Cfamiliaris'
[7] 'Dmelanogaster' 'Drerio' 'Ecoli'
[10] 'Gaculeatus' 'Ggallus' 'Hsapiens'
[13] 'Mfascicularis' 'Mfuro' 'Mmulatta'
[16] 'Mmusculus' 'Osativa' 'Ptroglodytes'
[19] 'Rnorvegicus' 'Scerevisiae' 'Sscrofa'
[22] 'Tgondii' 'Tguttata' 'Vvinifera'
獲取Ecoli 的BSgenome數(shù)據(jù)包,bing載入BSgenome.Ecolide包,查看數(shù)據(jù)信息: biocLite('BSgenome.Ecoli.NCBI.20080805')
library(BSgenome.Ecoli.NCBI.20080805)
結(jié)果如下: [1] 'BSgenome.Ecoli.NCBI.20080805'
[2] 'Ecoli'
查看數(shù)據(jù)包的概況: E. coli genome: # organism: Escherichia coli (E. coli)
# provider: NCBI
# provider version: 2008/08/05 (發(fā)行日期)
# release date: NA
# release name: NA
# 13 sequences: (包含的序列)
# NC_008253 NC_008563 NC_010468
# NC_004431 NC_009801 NC_009800
# NC_002655 NC_002695 NC_010498
# NC_007946 NC_010473 NC_000913
# AC_000091
# (use 'seqnames()' to see all the
# sequence names, use the '$' or '[['
# operator to access a given sequence)
獲取每個(gè)染色體的名字還有其對(duì)應(yīng)的長(zhǎng)度: > seqnames(Ecoli)
[1] 'NC_008253' 'NC_008563' 'NC_010468'
[4] 'NC_004431' 'NC_009801' 'NC_009800'
[7] 'NC_002655' 'NC_002695' 'NC_010498'
[10] 'NC_007946' 'NC_010473' 'NC_000913'
[13] 'AC_000091'
> seqlengths(Ecoli)
NC_008253 NC_008563 NC_010468 NC_004431
4938920 5082025 4746218 5231428
NC_009801 NC_009800 NC_002655 NC_002695
4979619 4643538 5528445 5498450
NC_010498 NC_007946 NC_010473 NC_000913
5068389 5065741 4686137 4639675
AC_000091
4646332
查看其中一條染色體 Ecoli$NC_008253
4938920-letter 'DNAString' instance
seq: AGCTTTTCATTCTGAC...TTAGTAAGTGATTTTC
查看NC_008253序列中GC的數(shù)量 > letterFrequency(Ecoli$NC_008253,'GC')
G|C
2495020
查看NC_008253序列中GC的含量 > letterFrequency(Ecoli$NC_008253,'GC',as.prob=TRUE)
G|C
0.5051752
結(jié)合sapply函數(shù)統(tǒng)計(jì)堿基組成和GC含量 #統(tǒng)計(jì)堿基組成
> sapply(seqnames(Ecoli), function(x) alphabetFrequency(Ecoli[[x]]))
NC_008253 NC_008563 NC_010468
A 1222723 1256126 1164516
C 1251581 1285309 1204681
G 1243439 1283517 1209555
T 1221177 1256945 1167466
M 0 11 0
R 0 34 0
W 0 25 0
S 0 17 0
Y 0 24 0
K 0 16 0
V 0 0 0
H 0 0 0
D 0 0 0
B 0 0 0
N 0 1 0
- 0 0 0
+ 0 0 0
. 0 0 0
統(tǒng)計(jì)GC含量
> sapply(seqnames(Ecoli), function(x) letterFrequency(Ecoli[[x]], letters = 'GC',
+ as.prob = TRUE))
NC_008253.G|C NC_008563.G|C NC_010468.G|C NC_004431.G|C NC_009801.G|C
0.5051752 0.5054729 0.5086652 0.5047480 0.5062162
NC_009800.G|C NC_002655.G|C NC_002695.G|C NC_010498.G|C NC_007946.G|C
0.5081961 0.5038297 0.5053706 0.5049668 0.5060419
NC_010473.G|C NC_000913.G|C AC_000091.G|C
0.5078129 0.5078970 0.5079958
學(xué)習(xí)心得雖然我對(duì)R也不是特別熟悉,跟著教程一個(gè)一個(gè)代碼來(lái)跑,然后逐步理解也不難,希望大家也可以像我一樣做到。最后也給大家推薦Coursera上講解BSgenome專題的一個(gè)視頻,幫助大家進(jìn)一步理解這個(gè)包。文章鏈接:https://www./learn/bioconductor/lecture/E57YU/biostrings
|