把fasta序列读入到R里面去
原贴来自于生信技能树论坛： http://www.biotrainee.com/thread-806-1-1.html把fasta序列读入到R里面去~
fasta是什么，我就不多说了！
你一定会遇到这个需求，把fasta序列读入到R里面，至于读进去变成一个字符串还是一个list还是一个对象，是后话！
本地读取
我这里先给几个本地读取的方法：
library(Biostrings) ;
s = readDNAStringSet("nm.fasta")
readBStringSet(filepath, format="fasta", nrec=-1L, skip=0L, seek.first.rec=FALSE, use.names=TRUE)
library(seqinr);
fastafile<-read.fasta(file = "proteins.fasta", seqtype = "AA",as.string = TRUE)
library(ape);
read.dna();
read.FASTA在线获取
然后有几个包的在线获取方法也很赞！
用rentrez包也可以完成同样的任务.
首先用 entrez_link()根据基因的entrez ID来得到转录本的ID，然后根据转录本的ID用entrez_fetch()获取序列~
gene_ids <- c(351, 11647)
linked_seq_ids <- entrez_link(dbfrom="gene", id=gene_ids, db="nuccore")
linked_transripts <- linked_seq_ids$links$gene_nuccore_refseqrna
head(linked_transripts)
all_recs <- entrez_fetch(db="nuccore", id=linked_transripts, rettype="fasta")读入到R里面是一个什么东西，就需要你仔细探究了。
class(all_recs) ## 只是一个纯粹的 character 而已
cat(strwrap(substr(all_recs, 1, 500)), sep="\n")
write(all_recs, file="my_transcripts.fasta")
temp <- tempfile() ## 编程技巧，把变量写入临时文件~
write(all_recs, temp)
parsed_recs <- ape::read.dna(all_recs, temp)还是我前面说的，
http://www.biotrainee.com/thread-778-1-1.html
R语言初学笔记之吾日三省吾身！
其实你看到这一个简单的帖子，是因为我懒，但是对初学者来说，这里面的内容非常之多，值得学习的地方更是数不胜数，希望大家勇敢的学习这些技巧，然后在下面回帖写出自己的感悟！
其它推荐：http://www.biotrainee.com/thread-805-1-1.html ，R语言的protr包计算多条蛋白序列相似度 里面说到了用这个包的readFASTA函数也可以直接读取url的序列。
本文参与 腾讯云自媒体分享计划，分享自微信公众号。
原始发表：2018-12-04，如有侵权请联系 cloudcommunity@tencent.com 删除
评论
登录后参与评论
推荐阅读
目录