美文网首页重点关注scATAC
BSgenome 构建自己的参考基因组

BSgenome 构建自己的参考基因组

作者: 尧小飞 | 来源:发表于2021-01-11 18:19 被阅读0次

    1. 目的

      最近在某个分析的时候,需要BSgenome object,然而我的物种不是常见的物种,是自己组装出来的多倍体物种,因此没法用已有的BSgenome参考基因组,需要自己构建BSgenome参考基因组,通过BSgenome相关介绍说明的时候,发现有一个How to forge a BSgenome data package说明,这说明了我们能够自己通过BSgenome包,构建自己的参考基因组,下面就是构建BSgenome参考基因组的过程

    2. 输入fasta的准备

      构建BSgenome参考基因组的方法有两种,一种是不需要masked sequences,另外一种是需要masked sequences文件,这里我选择第一种进行构建。其实两者最主要的差别在于是否提供masked sequences文件。

      构建BSgenome参考基因组具体需要什么文件呢?说明文档中也有详细的说明,其中最主要的文件就是一个fasta文件,这里的fasta文件可以是所有染色体合并一起的一个文件,也可以是单条染色体的fasta文件,文件格式可以是压缩文件格式,后续的处理几乎都是在fasta文件的基础上进行处理。

    输入fasta文件介绍
      前面提到了fasta文件需要进行处理,究竟如何处理呢?由于输入的fasta文件是TwoBit文件,因此需要下载faToTwoBit工具,对fasta文件进行转换,具体命令如下:
    wget  -c "http://hgdownload.cse.ucsc.edu/admin/exe/linux.x86_64/faToTwoBit"
    chmod 755 faToTwoBit
    faToTwoBit gene.fa gene.2bit
    

      这里faToTwoBit运行需要一定的时间,主要是与参考基因组的大小相关,生成的是二进制文件,生成twoBit后,就可以进行下一步的准备了。

    3. seed文件的准备

      构建BSgenome参考基因组另外一个需要准备的文件就是seed文件,其实这里的seed文件就是在开发R包的DESCRIPTION文件, 这里相当于构建了一个R包,R包需要添加一下DESCRIPTION信息,因此seed的文件,可以参考DESCRIPTION的格式说明,具体格式要求如下:

    seed文件的标准格式
      上图可以看出,seed文件,主要包括Package/Title/Description/Version/Author/Maintainer/License,其他的选填项可以不用填写。seed文件格式有很多种,可根据实际需要进行撰写;另外,也可以根据已经安装的BSgenome包来查看seed格式,具体的查看命令如下:
     library(BSgenome)
    seed_files <- system.file("extdata", "GentlemanLab", package="BSgenome")
    tail(list.files(seed_files, pattern="-seed$"))
    #[1] "BSgenome.Sscrofa.UCSC.susScr11-seed"        "BSgenome.Sscrofa.UCSC.susScr3-seed"         #"BSgenome.Sscrofa.UCSC.susScr3.masked-seed"  "BSgenome.Tguttata.UCSC.taeGut1-seed"    
    #这里可以看到有很多seed文件,可以随便cp一个文件,然后在此基础上就行即可。
    

      seed文件其他的一些可选项也有很多,比如这里的organism、common_name、genome、provider、release_date、source_url、organism_biocview等等,这些是为了更详细的记录一些参考基因组信息,一遍后续使用或者构建的时候了解相关信息,其他的一些说明,这里不再进行详细的介绍,可以直接通过How to forge a BSgenome data package进行相关信息的查询。

    可选项

    4. forgeBSgenomeDataPkg

      seed文件生成以后,就可以直接构建BSgenome包了,这里构建BSgenome包的具体命令如下:

    seed_file="gene.物种.seed"
    #seqs_srcdir;destdir 序列文件所在以及输出的位置
    forgeBSgenomeDataPkg(seed_file, seqs_srcdir=getwd(), destdir=getwd(), verbose=TRUE,unlink=TRUE)
    #forgeBSgenomeDataPkg(seed_file, verbose=TRUE)
    #unlink参数表示是否overwrite已有的目录,seqs_srcdir是twoBit的目录,destdir为生成包的目录,这里需要一定的时候。
    

    5. 包的构建以及安装

      前面forgeBSgenomeDataPkg运行成功以后,会在目的目录中生成R包的几个基础文件,如下图所示:

    forgeBSgenomeDataPkg结果目录
      这是一个基本的R包的结果目录,然后在此基础上,在linux命令进行包的build、check、install,具体的命令如下:
    R CMD build <pkgdir>
    R CMD check <tarball>
    R CMD INSTALL <tarball>
    

      一般来说,build没有问题的话,就不用再check了,我这边构建的包后进行check,结果报错,报错的是LaTeX有问题,然后尝试不管他,再直接安装,结果能够正常安装,并且正常导入包。因此这里的报错可以不用考虑,这应该是要生成pdf的说明文件报错,对后续其他的结果没有影响。

    LaTeX报错

    6. 后记

      这里构建BSgenome参考基因组是选择的最简单的方式,没有提供masked sequences文件,如果需要masked sequences文件,那么可以直接参考方说明文档。可以通过此方法,构建任何自己的参考基因组,进行其他相关的分析。

    参考文档

    BSgenome构建新的参考基因组

    2021年1月11日

    相关文章

      网友评论

        本文标题:BSgenome 构建自己的参考基因组

        本文链接:https://www.haomeiwen.com/subject/xvpdaktx.html