忍者ブログ
バイオインフォマティックス技術者試験、情報処理試験など、IT系の試験を基礎から勉強します。また、Javaなどプログラミングを勉強します。

【バイオインフォ実習】第5回:NCBI DBから複数配列を自動検索・一括取得!マルチFASTAの活用とSeqIO.parseによるループ統計解析

前回(第4回)は、あらかじめ判明している1つのアクセッション番号(TP53)を指定して、NCBI APIから単一データを取得する方法を解説しました。

しかし、実際のデータ解析の現場では、「特定の検索条件にヒットする複数の配列を一括取得し、まとめて統計解析や比較解析を行う」というフローが頻出します。

今回は、NCBIデータベースへのキーワード検索(esearch)から、複数データのマルチFASTA一括ダウンロード(efetch)、そしてBiopythonの SeqIO.parse() を使ったループ統計処理までを完全自動化するパイプラインを構築します。

1. マルチFASTA(Multi-FASTA)とは?

マルチFASTAとは、1つのファイル内に複数のDNA・RNA・タンパク質配列を順番に格納したフォーマットです。

>NC_013993.1 Homo sp. Altai mitochondrion, complete genome
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT...
>NC_012920.1 Homo sapiens mitochondrion, complete genome
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT...
>NC_011137.1 Homo sapiens neanderthalensis mitochondrion, complete genome
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT...

単一FASTAとマルチFASTAの比較

  • 単一FASTA: 1ファイル = 1配列SeqIO.read() で読み込む)
  • マルチFASTA: 1ファイル = 複数の配列SeqIO.parse() でループ処理する)

複数種の比較解析(マルチプルアライメントや系統樹作成)や、特定の遺伝子群を一括処理する際に必須となる形式です。

2. 検索から一括取得・ループ解析までの完全コード

以下のスクリプトは、NCBIの nuccore データベースから「ヒト属のミトコンドリア完全ゲノム(RefSeq)」を検索し、ヒットした上位3件をマルチFASTA形式で自動取得して解析するコードです。

from Bio import Entrez, SeqIO
from Bio.SeqUtils import gc_fraction

# 1. 共通設定
Entrez.email = "your_email@example.com" # 自身のメールアドレスに変更してください

# 検索クエリ(ヒト属のRefSeqミトコンドリア完全ゲノム)
search_term = "Homo sapiens[Organism] AND mitochondrion[Title] AND srcdb_refseq[PROP]"
max_results = 3 # 取得件数上限

print(f"検索クエリ: {search_term}")
print("NCBI データベースを検索中...")

# 2. NCBI 検索(esearch): 条件に合う ID リストの取得
with Entrez.esearch(db="nuccore", term=search_term, retmax=max_results) as handle:
    search_results = Entrez.read(handle)

id_list = search_results["IdList"]
count = search_results["Count"]

print(f"検索ヒット件数: {count} 件")
print(f"取得対象 ID リスト (上位{len(id_list)}件): {id_list}\n")

if not id_list:
    print("該当するデータが見つかりませんでした。")
    exit()

# 3. NCBI 一括取得(efetch): 複数 ID を指定してマルチFASTAを受信
print("マルチFASTAデータを一括ダウンロード中...")

# id_list(GI番号等のリスト)をカンマ区切り文字列にして API に送る
with Entrez.efetch(db="nuccore", id=",".join(id_list), rettype="fasta", retmode="text") as handle:
    # 取得したマルチFASTAデータをローカルファイルに保存
    output_filename = "ncbi_multi_records.fasta"
    with open(output_filename, "w") as out_f:
        out_f.write(handle.read())

print(f"保存完了: {output_filename}\n")

# 4. SeqIO.parse() によるマルチFASTAのループ解析
print(f"=== {output_filename} の解析結果 ===")

# 複数配列が含まれるため SeqIO.parse を使用(イテレータ処理)
total_length = 0
record_count = 0

for record in SeqIO.parse(output_filename, "fasta"):
    record_count += 1
    seq_len = len(record.seq)
    total_length += seq_len
    
    # GC含有率の計算
    gc_val = gc_fraction(record.seq) * 100
    
    print(f"[{record_count}] ID: {record.id}")
    print(f" 概要: {record.description[:60]}...")
    print(f" 配列長: {seq_len:,} bp")
    print(f" GC率 : {gc_val:.2f}%")
    print("-" * 60)

# 全体の要約統計
avg_length = total_length / record_count if record_count > 0 else 0
print(f"解析完了: 合計 {record_count} 件 | 平均配列長: {avg_length:,.1f} bp")

3. 実行結果の確認

上記のプログラムを実行すると、アルタイ人(デニソワ人洞窟の個体)、現代人、ネアンデルタール人のミトコンドリア全ゲノムデータが自動で取得され、以下のようにパースされます。

検索クエリ: Homo sapiens[Organism] AND mitochondrion[Title] AND srcdb_refseq[PROP]
NCBI データベースを検索中...
検索ヒット件数: 3 件
取得対象 ID リスト (上位3件): ['292606408', '251831106', '196123578']

マルチFASTAデータを一括ダウンロード中...
保存完了: ncbi_multi_records.fasta

=== ncbi_multi_records.fasta の解析結果 ===
[1] ID: NC_013993.1
    概要: NC_013993.1 Homo sp. Altai mitochondrion, complete genome...
    配列長: 16,570 bp
    GC率 : 44.31%
------------------------------------------------------------
[2] ID: NC_012920.1
    概要: NC_012920.1 Homo sapiens mitochondrion, complete genome...
    配列長: 16,569 bp
    GC率 : 44.36%
------------------------------------------------------------
[3] ID: NC_011137.1
    概要: NC_011137.1 Homo sapiens neanderthalensis mitochondrion, complete...
    配列長: 16,565 bp
    GC率 : 44.39%
------------------------------------------------------------
解析完了: 合計 3 件 | 平均配列長: 16,568.0 bp

4. コードとデータ処理のポイント

一連のコードの中で、APIとBiopythonがどのように連動しているかをエンジニア目線で解説します。

  • esearch による内部ID(GI番号)の取得:
    Entrez.esearch を実行すると、search_results["IdList"]['292606408', '251831106', '196123578'] というGI番号(NCBI内部の識別ID)が動的に返されます。
  • カンマ区切りによる efetch の一括リクエスト:
    id=",".join(id_list) により、複数のGI番号を1つのリクエストとして送信しています。NCBIサーバーからまとめてマルチFASTA形式で返却されるため、通信回数を削減して効率的なダウンロードが可能です。
  • SeqIO.parse() によるメモリ効率の良いループ処理:
    SeqIO.read() は1ファイル1配列専用ですが、SeqIO.parse() は配列を1つずつ順番に読み込むイテレータ(Generator)として機能します。そのため、数百〜数千件の配列が含まれる大規模ファイルでもメモリを消費せずに高速処理できます。

5. ITエンジニア的まとめ

今回の実習のポイントは以下の通りです。

  • 自動検索と一括取得: esearch で得たIDリストをカンマ連結して efetch に送ることで、マルチFASTAの一括ダウンロードが可能。
  • イテレータ処理によるパイプライン化: 複数配列の解析には SeqIO.parse() を使用し、for ループで各レコードの塩基長やGC率を動的計算する。

次回は【バイオインフォ実習】第6回:GenBank形式ファイルからのアノテーション抽出とPandas DataFrame化に挑戦します!




PR