已知转录本test.ID,根据biomart包批量获得ensembl数据库中对应序列信息,并以fa格式输出

该方法需要魔法,我已进行了调试,可以成功运行

转录本test.ID样式(ENSMUST00000040608.4),以xlsx形式储存

输入数据格式

library(Biostrings)
library(biomaRt)
library(openxlsx)
library(seqinr)

rm(list = ls())
genelist <- read.xlsx("待测基因testID汇总.xlsx")
genelist <- genelist$test_id
#移除genelist的小数点
genelist<-unlist(lapply(genelist, FUN = function(x) {return(strsplit(x, split = ".",fixed = T)[[1]][1])}))

mart <- useMart("ensembl", dataset="mmusculus_gene_ensembl") # 小鼠的数据库
Ensemble_to_seq <- function(x) {cat("Getting sequence for gene", x, "\n")  # 添加打印语句
                                seq = biomaRt::getSequence(id = x, 
                                type = "ensembl_transcript_id", 
                                seqType = "cdna", 
                                mart = mart)
  seq<-as.data.frame(seq)
  seq$"Length"<-lapply(seq[,1],function(y){return(nchar(y))})
  return(seq)
}

outTab <- data.frame()
for (i in genelist){outTab<-rbind(outTab,Ensemble_to_seq(i))} # 此步需要魔法

# 保存为fa格式文件
for (i in 1:nrow(outTab)) {
  # 定义序列
  seq <- DNAStringSet(outTab$cdna[i])
  # 定义描述信息
  descr <- paste(outTab$ensembl_transcript_id[i], "length=", outTab1$Length[i], sep = "|")
  # 输出序列
  write.fasta(seq, names = outTab$ensembl_transcript_id[i], file.out = "output.fa", open = "a")
}

 最后是输出结果展示

  • 2
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值