BioHaskell：阅读 FASTA 文件

Question

使用BioHaskell，如何读取包含氨基酸序列的FASTA文件？

我希望能够：

• 获取String序列列表
• 从 FASTA 注释（假定是唯一的）到序列字符串中获取一个 Map String String（来自 Data.Map）
• 在 BioHaskell 中实现的算法中使用序列。

注意：这个问题故意不显示研究工作，因为它立即以问答方式回答。

Answer 1

提取原始序列字符串

从现在开始，我们假设文件aa.fa 包含一些氨基酸FASTA 序列。让我们从一个提取序列列表的简单示例开始。

import Bio.Sequence.Fasta (readFasta)
import Bio.Sequence.SeqData (seqdata)
import qualified Data.ByteString.Lazy.Char8 as LB

main = do
    sequences <- readFasta "aa.fa"
    let listOfSequences = map (LB.unpack . seqdata) sequences :: [String]
    -- Just for show, we will print one sequence per line here
    -- This will basically execute putStrLn for each sequence
    mapM_ putStrLn listOfSequences

readFasta 返回IO [Sequence Unknown]。基本上这意味着没有关于序列是否包含氨基酸或核苷酸的信息。

请注意，我们在这里使用LB.unpack 而不是show，因为show 在生成的String 的开头和结尾添加了双引号（"）。使用LB.unpack 有效，因为在当前的BioHaskell 0.5.3 版本中，SeqData 只是定义为惰性ByteString。

我们可以使用castToAmino 或castToNuc 来解决这个问题：

转换为 AA/核苷酸序列

let aaSequences = map castToAmino sequences :: [Sequence Amino]

请注意，这些函数当前（BioHaskell 版本 0.5.3）不执行任何有效性检查。您可以在 BioHaskell 算法中使用[Sequence Amino] 或[Sequence Nuc]。

通过 FASTA 标头查找序列

我们现在假设我们的aa.fa 包含一个序列

>abc123
MGLIFARATNA...

现在，我们将从 FASTA 文件构建一个Map String String（在本例中我们将使用Data.Map.Strict）。我们可以使用这个映射来查找序列。

查找将产生Maybe String。此示例中的预期行为是在找到序列时打印该序列，如果在 Map 中找不到任何内容，则不打印任何内容。

由于Data.Maybe 是一个Monad，我们可以使用Data.Foldable.mapM_ 来完成这个任务。

import Bio.Sequence.Fasta (readFasta)
import Bio.Sequence.SeqData (Sequence, seqdata, seqheader)
import qualified Data.ByteString.Lazy.Char8 as LB
import Data.Foldable (mapM_)
import qualified Data.Map.Strict as Map

-- | Convert a Sequence to a String tuple (sequence label, sequence)
sequenceToMapTuple :: Sequence a -> (String, String)
sequenceToMapTuple s = (LB.unpack $ seqheader s, LB.unpack $ seqdata s)

main = do
    sequences <- readFasta "aa.fa"
    -- Build the sequence map (by header)
    let sequenceMap = Map.fromList $ map sequenceToMapTuple sequences
    -- Lookup the sequence for the "abc123" header
    mapM_ print $ Map.lookup "abc123" sequenceMap

编辑：感谢@GabrielGonzalez 的建议，最后一个示例现在使用Data.Foldable.mapM_ 而不是Data.Maybe.fromJust