快速从fasta序列文件中提取指定序列的方法

如何从fasta序列文件中提取指定序列？假如我有一个fasta序列文件，里面有>开头的行是ID信息，之后的内容是序列信息，如果有成千上万条序列，如何从中找到需要的序列？

>Sequence1
ATCGATCGATCGATCGATCGATCGATCG...
>Sequence2
ATCGATCGATCGATCGATCGATTAGCTA...
>Sequence3
ATCGATCGATCGATCGATCGATCGATCG...
>Sequence4
GCTAGCTAGCATCGATCGATCGTAGCTA...
......（此处省略10万字）

seqkit 是一个高效的命令行工具集，用于处理 FASTA 文件。你可以使用 seqkit grep 命令从 FASTA 文件中提取指定 ID 的序列。

安装方法

# 使用 conda 安装
conda install -c bioconda seqkit

# 使用 Homebrew 安装（适用于 macOS 和 Linux）
brew install seqkit

# 或者从 GitHub 下载二进制文件
wget https://github.com/shenwei356/seqkit/releases/download/v2.0.0/seqkit_linux_amd64.tar.gz
tar -zxvf seqkit_linux_amd64.tar.gz
mv seqkit /usr/local/bin/

提取指定 ID 的序列

假设你的输入 FASTA 文件名为 input.fasta，目标 ID 文件名为 ids.txt（每行一个 ID），输出文件名为 output.fasta。

准备包含目标 ID 的文件 ids.txt

Sequence3
Sequence4

注意不用添加>，只需要把像提取的序列名称输入即可。

提取指定序列

seqkit grep -f ids.txt input.fasta -o output.fasta

seqkit 会从 input.fasta 中提取 ids.txt 中列出的 ID 的序列，并将它们保存到 output.fasta 中。

解释：

seqkit grep: seqkit 的子命令，用于在序列文件中进行模式匹配。

-f ids.txt: 指定包含目标 ID 的文件，每行一个 ID。

input.fasta: 输入的 FASTA 文件。

-o output.fasta: 输出的 FASTA 文件

使用Python脚本提取指定序列

你可以使用 Python 脚本来处理 FASTA 文件，从中提取指定 ID 的序列并生成新文件。

def extract_fasta_sequences(input_file, output_file, target_ids):
    """
    从FASTA文件中提取指定ID的序列并保存到新文件中。

    参数：
    input_file -- 输入FASTA文件路径
    output_file -- 输出文件路径
    target_ids -- 目标ID列表
    """
    with open(input_file, 'r') as infile, open(output_file, 'w') as outfile:
        write_sequence = False
        for line in infile:
            if line.startswith('>'):
                # 检查ID是否在目标ID列表中
                sequence_id = line[1:].strip()
                if sequence_id in target_ids:
                    write_sequence = True
                    outfile.write(line)
                else:
                    write_sequence = False
            elif write_sequence:
                outfile.write(line)

if __name__ == "__main__":
    # 设置输入文件、输出文件和目标ID列表
    input_fasta_file = 'input.fasta'
    output_fasta_file = 'output.fasta'
    target_sequence_ids = ['id1', 'id2', 'id3']  # 替换为你的目标ID列表

    extract_fasta_sequences(input_fasta_file, output_fasta_file, target_sequence_ids)

脚本说明

extract_fasta_sequences 函数从输入的 FASTA 文件中提取指定 ID 的序列，并将它们写入到输出文件中。

input_file 是输入的 FASTA 文件路径。

output_file 是输出 FASTA 文件路径。

target_ids 是目标序列 ID 的列表。

使用方法

保存上述脚本到一个 Python 文件，例如 extract_fasta.py， 在同一目录下准备你的输入 FASTA 文件，例如 input.fasta， 编辑脚本中的 target_sequence_ids 列表，将其替换为你希望提取的序列 ID，最后运行如下代码：

python extract_fasta.py

运行脚本后，指定 ID 的序列将被提取并保存到 output.fasta 文件中，你可以根据需要修改文件路径和目标 ID 列表。

本文由 mdnice 多平台发布