Biopython を使用して fasta 形式で 複数の配列アラインメントを指定して、経験的なコドン置換マトリックスを構築しようとしています。
アラインされた配列が同じ長さである場合、 AlignInfoモジュールを使用した単一ヌクレオチド置換マトリックスは比較的単純なようです。これが私がpython2.7を使ってなんとかしたことです：
#!/usr/bin/env python
import os
import argparse
from Bio import AlignIO
from Bio.Align import AlignInfo
from Bio import SubsMat
import sys
version = "0.0.1 (23.04.20)"
name = "Aln2SubMatrix.py"
parser=argparse.ArgumentParser(description="Outputs a codon substitution matrix given a multi-alignment in FastaFormat. Will raise error if alignments contain dots (\".\"), so replace those with dashes (\"-\") beforehand (e.g. using sed)")
parser.add_argument('-i','--input', action = "store", dest = "input", required = True, help = "(aligned) input fasta")
parser.add_argument('-o','--output', action = "store", dest = "output", help = "Output filename (default = <Input-file>.codonSubmatrix")
args=parser.parse_args()
if not args.output:
args.output = args.input + ".codonSubmatrix" #if no outputname was specified set outputname based on inputname
def main():
infile = open(args.input, "r")
outfile = open(args.output, "w")
align = AlignIO.read(infile, "fasta")
summary_align = AlignInfo.SummaryInfo(align)
replace_info = summary_align.replacement_dictionary()
mat = SubsMat.SeqMat(replace_info)
print >> outfile, mat
infile.close()
outfile.close()
sys.stderr.write("\nfinished\n")
main()
同じ長さの配列 (aln.fa) を持つ fasta 形式の複数配列アラインメント ファイルを使用すると、出力は、アラインメントで観察されたヌクレオチド置換の数に対応するハーフ マトリックスです (ギャップ (-) が許容されることに注意してください)。
python Aln2SubMatrix.py -i aln.fa
- 0
a 860 232
c 596 75 129
g 571 186 75 173
t 892 58 146 59 141
- a c g t
私が目指しているのは、同様の経験的置換行列を計算することですが、複数の配列アラインメントに存在するすべてのヌクレオチド トリプレット (コドン) についてです。
以下を変更して、ヌクレオチド トリプレットを受け入れるために、AlignInfo モジュールの_pair_replacement関数を微調整しようとしました。
305～308行目
for residue_num in range(len(seq1)):
residue1 = seq1[residue_num]
try:
residue2 = seq2[residue_num]
に
for residue_num in range(0, len(seq1), 3):
residue1 = seq1[residue_num:residue_num+3]
try:
residue2 = seq2[residue_num:residue_num+3]
この段階で、アラインメントからコドンを取得できますが、アルファベットについて不満があります (モジュールは 1 文字のアルファベットしか受け入れませんか?)。
ご了承ください
(i) 3 つの可能な読み枠を説明する置換行列を取得したい
どんな助けでも大歓迎です。