はじめに

等温滴定カロリメトリー（ITC）は、タンパク質とリガンドの結合を解析することができる手法のひとつです。濃度条件などを検討する過程で、結果のグラフをシミュレーションしてみたいと考えたので、Pythonを用いたITCデータシミュレーションの簡単なプログラムを作ってみました。
いろいろなパラメータをいじったときに、結果のプロットがどう変化するか、原理がイメージできていれば特に必要ないかもしれませんが、はじめて実験するとき私には難しかったので、プログラム上でなら自由に値をいじれるだろうから、遊んでみてイメージしよう、というモチベーションです。（ITCの解析用ソフトフェアにも同様の機能があるようです。）

プログラムの実装

基本的な考え方

ITCでは、滴定ごとにリガンドがタンパク質と結合する際の熱量変化を測定します。なお、ITCの測定の詳細は総説やメーカーの解説記事が詳しいのでそちらをご参照ください（参考文献）。この熱量は、以下の手順で考えます。

結合平衡式: タンパク質 (P) とリガンド (L) が結合して複合体 (PL) を形成する（今回は、結合比が1:1のシンプルな結合を考えます）

   P + L \rightleftharpoons PL

解離定数 ($K_d$): 結合の親和性を表す定数

   K_d = \frac{[P][L]}{[PL]}

ここで、$[P]_0$をタンパク質の初期濃度、$[L]$を滴下したリガンドの濃度とすると、

   K_d = \frac{([P]_0-[PL])([L]-[PL])}{[PL]}

となるので、結合状態の濃度$[PL]$は

[PL] = \frac{[P]_0 + [L] + K_d - \sqrt{([P]_0 + [L] + K_d)^2 - 4 [P]_0 [L]}}{2}

で求められます。
ITCで計測する熱量の変化は、

   dQ = d[PL]V_0\Delta H

の関係があります。実際の実験では、発生する熱の差分をとった$\Delta Q$が測定されます。

シミュレーション

1. 必要なライブラリのインポート

import numpy as np
import matplotlib.pyplot as plt

2. 入力パラメータの設定

以下のパラメータを設定します。

解離定数 $K_d$
結合に伴うエンタルピー変化 $\Delta H$
タンパク質濃度 $[P]_0$
反応セルの体積 $V_0$
結合部位（反応比） $n$
リガンド終濃度（最後まで滴下した場合）と滴定回数

# 結合モデルパラメータ
Kd_true = 1e-6           # 結合定数 (M)
DeltaH_true = -40e3      # 結合に伴うエンタルピー変化 (J/mol)
P0 = 20e-6               # 初期タンパク質濃度 (M)
V0 = 200e-6              # 反応セルの体積 (L)
n_true = 1.0             # 結合部位数 (1:1モデル)

# リガンド濃度の滴定ステップ
L_conc = np.linspace(0, 40e-6, 20)  # 滴定ごとのリガンド濃度

3. 熱量計算

累積熱量は、$\Delta H$と$PL$のモル数から、

Q = n \cdot \Delta H \cdot [PL] \cdot V_0

とし、差分熱量$\Delta Q$を滴定ごとに求めます。

# 結合モデル: 累積熱 Q の計算
def binding_heat(L, n, DeltaH, Kd, P0, V0):
    PL = (P0 + L + Kd - np.sqrt((P0 + L + Kd)**2 - 4 * P0 * L)) / 2  # 結合状態 (M)
    Q = n * DeltaH * PL * V0                                        # 熱量 (J)
    return Q

# 累積熱量 Q の計算
Q = binding_heat(L_conc, n_true, DeltaH_true, Kd_true, P0, V0)

# 差分熱量 ΔQ の計算
delta_Q = np.diff(Q)

4. 結果の可視化

シミュレーション結果をプロットします。

plt.figure(figsize=(8, 6))
plt.plot((L_conc[1:])/P0, delta_Q, 'o-', label='ΔQ', color='red')
plt.xlabel("molar ratio")
plt.ylabel("ΔQ (J)")
plt.legend()
plt.show()

結果の概要

まず、上記のパラメータで実行すると、次のようになりました。

では次に、$K_d$の値の違いによって、カーブがどのように変わるかを見てみます。

# Kdの値を3つ設定
Kd_values = [1e-7, 1e-6, 10e-6]

# グラフの描画
plt.figure(figsize=(8, 6))

# 各Kdの値に対してΔQを計算してプロット
for Kd in Kd_values:
    # 累積熱量 Q の計算
    Q = binding_heat(L_conc, n_true, DeltaH_true, Kd, P0, V0)

    # 差分熱量 ΔQ の計算
    delta_Q = np.diff(Q)

    # ΔQのプロット
    plt.plot((L_conc[1:])/P0, delta_Q, label=f'Kd = {Kd:.0e}', marker='o')

$K_d$が10 µMになってくると、シグモイドとしてはうまくプロットされづらくなってきそうだという気がします。この場合は、タンパク質・リガンドの濃度を上げる必要が出てきそうです。
では、この$K_d$において、適切なタンパク質・リガンド濃度を考えます。ここでは、Molar ratioが2になるまでをみることにして量比を固定し、次のように変えてみました。

# P0の値を3つ設定
P0_values = [20e-6, 60e-6, 100e-6]

# グラフの描画
plt.figure(figsize=(8, 6))

# 各Kdの値に対してΔQを計算してプロット
for P0 in P0_values:
    L_conc = np.linspace(0, 2 * P0, 20)  # 滴定ごとのリガンド濃度
    # 累積熱量 Q の計算
    Q = binding_heat(L_conc, n_true, DeltaH_true, Kd, P0, V0)

    # 差分熱量 ΔQ の計算
    delta_Q = np.diff(Q)

    # ΔQのプロット
    plt.plot((L_conc[1:])/P0, delta_Q, label=f'P0 = {P0:.0e}', marker='o')

タンパク質の初期濃度を100 µMにすると、$K_d$が10 µMであってもうまくプロットされそうな気がします。

参考文献

（最終閲覧はいずれも2024/12/19）

cytiva「原理：其は王道を歩む基礎体力（7）　等温滴定型熱量測定」https://www.cytivalifesciences.co.jp/technologies/biacore/road/road01_07.html
織田昌幸（2020）「等温滴定熱測定」蛋白質科学会アーカイブ, 1, e030 https://www.pssj.jp/archives/files/articles/030.pdf