PythonでDNA配列を分析する
このトピックを終えたら
PythonでDNA配列のGC contentを計算し、FASTAファイルをパースし、コドンをアミノ酸に変換できるようになります。
変数:データを入れる器
実験ノートに「試料濃度:2.5 mg/mL」と記すように、Pythonでは変数に値を保存します。
gene_name = "BRCA1"sequence = "ATGCGATCGATCGATCG"gc_ratio = 0.529
print(f"遺伝子: {gene_name}")print(f"配列の長さ: {len(sequence)}bp")print(f"GC比率: {gc_ratio:.1%}")
assert gene_name == "BRCA1"assert len(sequence) == 17文字列:DNA配列を扱う
DNA配列はA、T、G、Cの4文字からなる文字列(string)です。Pythonの文字列メソッドで配列を分析できます。
sequence = "ATGCGATCGATCGATCG"
# 特定の塩基を数えるg_count = sequence.count("G")c_count = sequence.count("C")print(f"G: {g_count}個, C: {c_count}個")
# 相補配列を作る (A↔T, G↔C)complement_table = str.maketrans("ATGC", "TACG")complement = sequence.translate(complement_table)print(f"元の配列: {sequence}")print(f"相補配列: {complement}")print(f"逆相補: {complement[::-1]}")
assert complement == "TACGCTAGCTAGCTAGC"assert complement[::-1] == "CGATCGATCGATCGCAT"関数:GC Content計算機
関数(function)は実験プロトコルの1ステップを再利用可能にパッケージしたものです。一度作れば、どんな配列にも適用できます。
def calculate_gc_content(sequence: str) -> float: sequence = sequence.upper() gc_count = sequence.count("G") + sequence.count("C") return (gc_count / len(sequence)) * 100
# テストseq1 = "ATGCGATCGATCGATCG"seq2 = "AAAAAAAAAA"seq3 = "GGGGGGGGGG"
print(f"seq1 GC: {calculate_gc_content(seq1):.1f}%")print(f"seq2 GC: {calculate_gc_content(seq2):.1f}%")print(f"seq3 GC: {calculate_gc_content(seq3):.1f}%")
assert abs(calculate_gc_content(seq1) - 52.9) < 0.1assert calculate_gc_content(seq2) == 0.0assert calculate_gc_content(seq3) == 100.0リストとループ:複数の配列を一度に処理
96ウェルプレートのすべてのウェルに同じ処理をするように、forループは複数のデータに同じ作業を繰り返します。
def calculate_gc_content(sequence: str) -> float: sequence = sequence.upper() gc_count = sequence.count("G") + sequence.count("C") return (gc_count / len(sequence)) * 100
genes = { "BRCA1": "ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC", "TP53": "ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA", "EGFR": "ATGCGACCCTCCGGGACGGCCGGGGCAGCGCTCCTGGCGCTG",}
results = []for name, seq in genes.items(): gc = calculate_gc_content(seq) results.append(gc) print(f"{name}: GC={gc:.1f}%, 長さ={len(seq)}bp")
assert len(results) == 3assert all(0 <= gc <= 100 for gc in results)ディクショナリ:FASTAファイルをパースする
ディクショナリ(dictionary)は実験試料にラベルを貼るのと同じです。遺伝子名(キー)で配列(値)をすぐに見つけられます。
def parse_fasta(fasta_text: str) -> dict[str, str]: sequences: dict[str, str] = {} current_header = "" for line in fasta_text.strip().split("\n"): if line.startswith(">"): current_header = line[1:].strip() sequences[current_header] = "" else: sequences[current_header] += line.strip() return sequences
sample_fasta = """>BRCA1_humanATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC>TP53_humanATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA"""
result = parse_fasta(sample_fasta)print(f"パースされた配列数: {len(result)}")for name, seq in result.items(): print(f" {name}: {len(seq)}bp")
assert len(result) == 2assert result["BRCA1_human"] == "ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC"assert result["TP53_human"] == "ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA"コドン → アミノ酸変換
DNA配列の3文字(コドン)が1つのアミノ酸を指定します。この変換テーブルをディクショナリで作れます。
CODON_TABLE = { "ATG": "M", # メチオニン(開始コドン) "TTT": "F", "TTC": "F", # フェニルアラニン "TTA": "L", "TTG": "L", "CTT": "L", "CTC": "L", # ロイシン "GAT": "D", "GAC": "D", # アスパラギン酸 "GAA": "E", "GAG": "E", # グルタミン酸 "GCT": "A", "GCC": "A", # アラニン "TAA": "*", "TAG": "*", "TGA": "*", # 終止コドン}
def translate_sequence(dna: str) -> str: protein = [] for i in range(0, len(dna) - 2, 3): codon = dna[i:i+3] amino_acid = CODON_TABLE.get(codon, "?") if amino_acid == "*": break protein.append(amino_acid) return "".join(protein)
test_seq = "ATGGATTTTGAA"protein = translate_sequence(test_seq)print(f"DNA: {test_seq}")print(f"Protein: {protein}")
assert protein == "MDFE"やってみよう(Faded Example)
空欄を埋めて、GC content計算関数を完成させてください。
def gc_content(seq):seq = seq.upper()gc = seq.count('G') + seq.count('')return / len(seq) * 100
よくあるエラーと解決法
Q: IndentationError: expected an indented block
Pythonはインデント(字下げ)でコードブロックを区別します。def、for、ifの次の行は必ずスペース4つのインデントが必要です。
Q: KeyError: 'BRCA1'
ディクショナリに該当するキーがない時に発生します。dict.get("BRCA1", "なし")を使えば、キーがなくてもエラーなしでデフォルト値を返します。
Q: 配列に小文字が混ざっていてcountが合いません
.upper()でまず大文字に統一してから処理してください。実際のFASTAファイルでは、小文字はリピートマスクされた領域(repeat-masked region)を意味することがあります。
次の記事では、Gitでこの分析コードをバージョン管理する方法を学びます。