加载中…
个人资料
  • 博客等级:
  • 博客积分:
  • 博客访问:
  • 关注人气:
  • 获赠金笔:0支
  • 赠出金笔:0支
  • 荣誉徽章:
正文 字体大小:

BSgenome构建新的参考基因组

(2017-10-21 17:18:34)
标签:

bsgenome

fatotwobit

BSgeome自己有很多构建好的基因组但是有时候我们需要自己构建基因组的信息


参考文档: https://www.bioconductor.org/packages/devel/bioc/vignettes/BSgenome/inst/doc/BSgenomeForge.pdf


现在输入的文件是下载的dm6.fa,里面有果蝇的基因组的参考序列信息

 

####前期工作,在linux下工作

 

#1trans fa to twoBit formate. BSgeome不识别fa的格式. 需要下载ucscfaToTwoBit工具。

#下载地址:http://hgdownload.cse.ucsc.edu/admin/exe/linux.x86_64/faToTwoBit

 

faToTwoBit dm6.fa dm6.2bit

 

##2)获得每一条染色体的名字

less -S dm6.fa | grep ">" |awk '{print $1}' | sed 's/^>//g' >dm6.chromName.txt

 

#程序的关键是seed文件,参考已有的seed修个seed格式的文件。根据参考文档第四页的2.2 Prepare the BSgenome data packages seed file.

 

#修改,给出序列的文件名,比如 seqfile_name: dm6.2bit

##找出线粒体的名字,在seed文件中的 添加<</span>线粒体>的名字,并且手动删除dm6.chromName.txt中的<</span>线粒体>名字

 

 

#####R中的工作。

 

library(BSgenome)

 

#读取前面构建的染色体名字的文件,注意 seqnames应该和seed文件中的关键字一致

seqnames=read.table("dm6.chromName.txt")

seqnames=as.character(seqnames$V1)

 

#染色体名有<211000022280427>,这种纯数字的,这种chr是不能识别的

seqnames1=seqnames[grep("^211",seqnames,invert=T)]

 

 

#给出seed文件

dm6_seed="BSgenome.Dmelanogaster.Ensemble.dm6-seed"

 

#seqs_srcdir;destdir 序列文件所在以及输出的位置

forgeBSgenomeDataPkg(dm6_seed, seqs_srcdir=getwd(), destdir=getwd(), verbose=TRUE)

 

forgeBSgenomeDataPkg(dm6_seed, verbose=TRUE)


 

##参考 说明文档2.3 构建R包。

 

R CMD build BSgenome.Dmelanogaster.Ensemble.dm6

R CMD check BSgenome.Dmelanogaster.Ensemble.dm6_1.4.2.tar.gz

R CMD INSTALL BSgenome.Dmelanogaster.Ensemble.dm6_1.4.2.tar.gz

 

#大功告成。可以欢快的用 BSgenome.Dmelanogaster.Ensemble.dm6 这个自己构建的包了。



附上我自己写的seed文件,成功的关键就是seed文件的设置,需要理解里面的参数。


 

Package: BSgenome.Dmelanogaster.Ensemble.dm6
Title: Full genome sequences for Drosophila melanogaster (Ensemble version dm6)
Description: Full genome sequences for Drosophila melanogaster(furitfly) as provided by Ensemb
le (dm6, Nov. 2004) and stored in Biostrings objects.
Version: 1.4.2
organism: Drosophila melanogaster
common_name: Furitfly
provider: Ensemble
provider_version: dm6
release_date: Nov. 2004
release_name: Baylor College of Medicine HGSC v3.4
#source_url: http://hgdownload.cse.ucsc.edu/goldenPath/rn4/bigZips/
organism_biocview: Drosophila_melanogaster
BSgenomeObjname: Dmelanogaster
#seqnames: paste("chr", c(1:20, "X", "M", "Un", paste(c(1:20, "X", "Un"), "_random", sep="")),
 sep="")
#seqnames: set the sequence names in R first: seqnames_all
seqnames: seqnames1
circ_seqs: "dmel_mitochondrion_genome"
#SrcDataFiles: chromFa.tar.gz from http://hgdownload.cse.ucsc.edu/goldenPath/rn4/bigZips/
#seqs_srcdir: /fh/fast/morgan_m/BioC/BSgenomeForge/srcdata/BSgenome.Rnorvegicus.UCSC.rn4/seqs
#seqfile_name: the full of sequence name
seqfile_name: dm6.2bit




0

阅读 收藏 喜欢 打印举报/Report
  

新浪BLOG意见反馈留言板 欢迎批评指正

新浪简介 | About Sina | 广告服务 | 联系我们 | 招聘信息 | 网站律师 | SINA English | 产品答疑

新浪公司 版权所有