PCRプライマーのTm値の計算機を作ってみよう③~pythonで実装編~

seq = 'ATGAGCCCTGAAGTG'
AT_count = seq.count('A') + seq.count('T')
print(AT_count)

GC_count = seq.count('G') + seq.count('C')
print(GC_count)
#7  →ATの総和
#8  →GCの総和

ATとGCのカウントができました。
あとはAT数に2をGC数に4を掛け、その和がWallace法でのTm値となります。

Wallace法によるTm値の算出

tm = AT_count*2 + GC_count*4
print(tm)
#46 →Wallace法で算出したTm値

最近接塩基対法による計算

では次に、最近接塩基対法によるTm値の算出を行います。
詳細は前回のブログをご参照ください。

パラメータ辞書の作成

まずは、最近接塩基対法のパラメータを作成します。
それぞれの近接ペアの[ΔH, ΔS]の順でリストを作成し、辞書化します。

#最近接塩基対法のパラメータを辞書型オブジェクトとして作成する
param = {
    'AA':[-9.1, -24],
    'TT':[-9.1, -24],
    'AT':[-8.6, -23.9],
    'TA':[-6, -16.9],
    'CA':[-5.8, -12.9],
    'TG':[-5.8, -12.9],
    'GT':[-6.5, -17.3],
    'AC':[-6.5, -17.3],
    'CT':[-7.8, -20.8],
    'AG':[-7.8, -20.8],
    'GA':[-5.6, -13.5],
    'TC':[-5.6, -13.5],
    'CG':[-11.9, -27.8],
    'GC':[-11.1, -26.7],
    'GG':[-11, -26.6],
    'CC':[-11, -26.6]
        }

隣接ペアの抽出

次に、プライマーの隣接する2塩基のペアをすべてリスト化します。
文字列もリスト同様スライスが使えますので、
前から順に2塩基ずつ取り出します。
ここはリスト内包表記を使うと良いでしょう！

seq = 'ATGAGCCCTGAAGTG'      
pair = [seq[i:i+2] for i in range(len(seq)-1)]
print(pair)
#['AT', 'TG', 'GA', 'AG', 'GC', 'CC', 'CC', 'CT', 'TG', 'GA', 'AA', 'AG', 'GT', 'TG']

ΔHとΔSの総和を求める

すべての隣接ペアのリストが作成できましたので、
先程作成した辞書型オブジェクトparamを参照し、
それぞれの隣接塩基ペアに対するΔHとΔSを取得してリスト化します。

#リスト内包表記を使ってすべてのペアに対するΔHをparamから抜き出します。
H_list = [param[i][0] for i in pair] 
#同様にΔSも抜き出し、リスト化します。
S_list = [param[i][1] for i in pair]
print(H_list)
print(S_list)
#[-8.6, -5.8, -5.6, -7.8, -11.1, -11, -11, -7.8, -5.8, -5.6, -9.1, -7.8, -6.5, -5.8] ΔHのリスト
#[-23.9, -12.9, -13.5, -20.8, -26.7, -26.6, -26.6, -20.8, -12.9, -13.5, -24, -20.8, -17.3, -12.9]　ΔSのリスト

それぞれのΔHおよびΔSの総和を求めます。
Pythonのリスト内の数字の総和は、
sum(リストオブジェクト)で求められます。

sumH = sum(H_list)
sumS = sum(S_list)
print(sumH)
print(sumS)
#-109.29999999999998 →sumHの値
#-273.2 →sumSの値

これで最近接塩基対法に必要な材料がそろいましたので、
最近接塩基対法の式に代入して計算したいと思います。

最近接塩基対法によるTm値の算出

Tm = (1000 x ΔH) /( -10.8 + ΔS + 1.987 x ln(Conc/4)) - 273.15 + 16.6log[Na+]でしたので、
これを計算していきます。
Na+の濃度は50mM、プライマーの濃度(Conc)は0.5 μMを使用します。

また、log10やlnの計算はmathライブラリを使います。

log10を使うとき：math.log10(真数)
lnを使うとき：math.log(真数)

import math
mol = 0.5*(10**(-6))
Na = 50*(10**(-3))

tm = (1000*H/(-10.8+S+1.987*math.log(mol/4)))-273.15+16.6*math.log10(Na)
 print(tm)

以上で最近接塩基法によるTm値の計算ができるようになりました。

LifeScience Hack

生物系創薬研究者がAI（誇大表示）を手に入れるまでの過程（Python、Deep Learning、ライフサイエンス）