制作fasta文件(標準化)

該教材是如何制作固定每行70堿基的fasta文件,主要設計如何將特別長的堿基序列均勻切割成每行70個堿基的文件。代碼如下:
如何構建fasta文件及其index索引文件
fasta文件更改header

rm(list = ls())
suppressMessages(library(Biostrings))
suppressMessages(library(tidyverse))

trim_fa<-function(fa_dir,prefix,out_dir){
  trimSeq <- function(Seq, binwidth = 70){
    N_line = nchar(Seq) / binwidth 
    if(N_line > as.integer(N_line)){
      N =  as.integer(N_line) + 1
    } else{
      N =  as.integer(N_line)
    }
    
    out_seq<-vector()
    bases = strsplit(Seq, "")[[1]]
    for(i in 0:(N-1)){
      xseq = paste0(bases[(1 + binwidth*i):(binwidth + i*binwidth)], collapse = "") 
      if(i != (N-1)){
        xseq = paste0(bases[(1 + binwidth*i):(binwidth + i*binwidth)], collapse = "") 
      } else{
        xseq = paste0(bases[(1 + binwidth*i):length(bases)], collapse = "") 
      }
      out_seq<-paste0(out_seq,xseq, "\n")
    }
    return(out_seq)
  }
  fa<-readDNAStringSet(fa_dir)
  seq_name = names(fa)
  sequence = paste(fa)
  df <- data.frame(seq_name, sequence)
  df$sequence<-trimSeq(Seq = as.character(df$sequence))
  df$seq_name<-paste0(">",df$seq_name)
  D <- do.call(rbind, lapply(seq(nrow(df)), function(i) t(df[i, ])))
  write.table(D,paste0(out_dir,prefix,".fa"), row.names = FALSE, col.names = FALSE, quote = FALSE)
}
整整齊齊
最后編輯于
?著作權歸作者所有,轉載或內容合作請聯(lián)系作者
【社區(qū)內容提示】社區(qū)部分內容疑似由AI輔助生成,瀏覽時請結合常識與多方信息審慎甄別。
平臺聲明:文章內容(如有圖片或視頻亦包括在內)由作者上傳并發(fā)布,文章內容僅代表作者本人觀點,簡書系信息發(fā)布平臺,僅提供信息存儲服務。

相關閱讀更多精彩內容

友情鏈接更多精彩內容