linux下，如何從保存多條基因數據的.fa文件中提取特定的一條基因數據？

先描述個人項目內容：

1. 將50bp長的DNA序列進行單次比對（linux環境，算法gapmis已經寫好）；
2. 500萬個基因序列文件單次比對，會耗費大量I/O時間。爲此但願將1萬條基因數據保存在一個AT50_1_0.fasta文件中，每一條基因數據單獨保存爲一行，以下圖所示：
3. 依次提取各行數據，並調用比對算法gapmis，輸出每一行的比對結果。

---

「將包含一條基因數據的文件依次進行比對，轉化成二維數據進行比對」，直接目前存在的問題：

1. linux環境下，是否能夠編寫c循環程序：對.fa中的文件按行讀取？
2. 「.fa」文件格式說明：按照「>」標識來界定是否爲一條基因數據。若是隻有一個「>」，斷定只存在一條基因數據；
3. gapsmis程序在linux環境下的比對命令爲："./gapsmis -a a.fasta -b b.fasta" （將序列a與序列b進行比對）。換句話說：咱們須要修改gapsmis程序接口，將命令中的「文件名」輸入形式，轉換爲「基因字符串」輸入形式。

（具體解決思路，等我實現再補充，各位大佬見諒）

---

下面補充一個問題以及解決方案：如何經過基因ID文件從fasta文件中提取特定的基因？

概念介紹：

• .fa文件格式：

  「>基因ID
     基因序列」
  以下圖：AT50_0就是基因ID

---

解決思路：

• 訪問網址http://hgdownload.cse.ucsc.ed...
• 下載faSomeRecords腳本：faSomeRecords.txt
• 將txt後綴刪除（這裏爲了上傳方便），拷貝到指定的文件夾，運行命令行：
• $ chmod +x faSomeRecords #賦予文件可執行權限，爲Linux系統下執行

• 執行以下命令，進行數據處理

  ./faSomeRecords AT_50_1_a.fasta ID.txt query.fasta 
     #其中AT_50_1_a.fasta是原始的fasta文件，ID.txt列入了須要查找的基因ID，每行一個，query.fasta爲輸出文件，包含對應ID的序列信息。

---

linux操做圖：

AT_50_1_a.fasta是原始的fasta文件：

ID.txt列入了須要查找的基因ID：

query.fasta爲輸出文件：

最終從原始文件提取出基因ID爲：AT50_0,AT50_5的基因序列，並將結果保存在query.fasta文件中。

---

注：本文參考科學網陳振璽博客：http://blog.sciencenet.cn/blo... Tanvir Ahamed博士在biostarts中的回覆（https://www.biostars.org/p/12...），實際操做後整理所得。