BSgenome简介

一、BSgenome和BSgenome数据包

Bioconductor提供了某些物种的全基因组序列数据包，这些数据包是基于Biostrings构建的，称为BSgenome数据包。不同物种的BSgenome数据包都有类似的数据结构，可以用统一的方式进行处理。但是BSgenome数据包仅包含有数据，它们的处理的方法由另外一个软件包提供，即BSgenome包。

先安装BSgenome包（如果没有安装）：

if (!requireNamespace("BiocManager", quietly = TRUE))
  install.packages("BiocManager")

BiocManager::install("BSgenome")

载入BSgenome包，并查看当前版本提供的BSgenome数据包：

library(BSgenome)
(ag <- available.genomes())

#当前版本提供了18个物种的37个基因组包

unique(gsub("BSgenome\\.([^\\.]+).*", "\\1", ag))

获取拟南芥的BSgenome数据包（BSgenome.Athaliana.TAIR.TAIR9）：

biocLite("BSgenome.Athaliana.TAIR.TAIR9")

载入BSgenome.Athaliana.TAIR.TAIR9包，查看数据信息：

library(BSgenome.Athaliana.TAIR.TAIR9)
# BSgenome.Athaliana.TAIR.TAIR9只定义了一个对象
ls("package:BSgenome.Athaliana.TAIR.TAIR9")

#查看总体信息

Athaliana

#不载入序列就可以获取染色体名称和长度信息

seqnames(Athaliana)

seqlengths(Athaliana)

二、BSgenome方法

BSgenome数据包内包含的序列为DNAString，DNAStringSet或MaskedDNAString对象，完全可以使用Biostrings包的方法进行处理：

#可以使用多种方式获取某一序列

Athaliana$Chr1

Athaliana[[1]]

Athaliana[["Chr1"]]

下面使用sapply函数统计碱基组成和GC含量。
#统计碱基组成

sapply(seqnames(Athaliana), function(x) alphabetFrequency(Athaliana[[x]]))

#计算GC含量

sapply(seqnames(Athaliana), function(x) letterFrequency(Athaliana[[x]], letters = "GC", 
    as.prob = TRUE))

但是BSgenome包试图提供一些简便的方法来处理基因组水平的BSgenome类数据，例如类似于apply函数的bsapply函数。要使用bsapply函数得先构建BSParams类对象，用于设置以下参数：

X：将要处理的BSgenome对象
FUN：将要对X中每条染色体进行处理函数
exclude：字符向量，表示排除的染色体名称
simplify：逻辑向量，它的意义和与sapply的simplify参数一样。默认为FALSE，bsapply返回GenomeData类数据；如果设置为TRUE，函数的结果尽量使用表格（数据框）类型显示
maskList：逻辑向量，表示各染色体使用掩膜的应用状态
motifList：字符向量，对序列进行掩膜的motif
userMask：RangesList对象，每个元素将用于对响应染色体进行掩膜处理，为用户提供额外的掩膜。
invertUserMask：TRUE/FALSE，表示是否对userMask进行反转。 

FUN函数的其他参数可以在bsapply函数中以有名参数的方式进行设置。 下面代码的作用和前面sapply函数的应用是相同的：

op <- options()
options(digits = 4)
params <- new("BSParams", X = Athaliana, FUN = letterFrequency, simplify = TRUE)
bsapply(params, letters = "GC", as.prob = TRUE)

BSParams是S4类，对象数据可以修改，方便重复使用：

params@motifList <- "N"
bsapply(params, letters = "GC", as.prob = TRUE)

options(op)
params@FUN <- alphabetFrequency
bsapply(params, baseOnly = TRUE)

BSgenome包实现了BSgenome类对象的matchPWM，countPWM，vmatchPattern，vcountPattern，vmatchPDict和vcountPDict等操作，除继承了Biostrings中相应的方法外还增加了一些针对BSgenome类对象的参数设置（和bsapply函数参数类似）：

vmatchPattern("ACGTGKC", Athaliana, fix = "subject", exclude = c("ChrC", "ChrM"))

BSgenome包还提供了其他一些方法，比如用户自行构建BSgenome数据包和SNP的处理等。这些不在我的关心范围之内，没去了解。

转载自：https://blog.csdn.net/u014801157/article/details/24372467