やったこと

一行で書くと「細胞の理論生物学: ダイナミクスの視点から」の一章の前半30%くらいをまとめて行間を埋める等しました．これから少しだけ細かめに書いていきます．

イントロ

生物システムはその内部に「状態」という時間的に変化しうる状態を持っている．例えば，生態系における種の個体数や，多細胞生物を構成する各種の細胞の数や配置等である．これら「状態」の変化や多様化を記述するための手段として，力学系理論で用いられている状態空間という手法を使う．状態空間では，状態をいくつかの変数の組で表現する．例えば，状態が$k$個の変数の組で与えられるとすると，その状態を${x_1, x_2, \cdots, x_k}$で表す.こうすれば，状態を$k$次元空間の1点として表現することができる．

例1.1 たとえば(細胞) = (多くの化学成分からなる，その分子の数(濃度))で状態をあらわす．細胞内にある化学成分を$k$個用意し，成分1, 2, 3, $\cdots$, $k$の濃度の組$x_1, x_2, x_3, \cdots, x_k$で細胞の状態をあらわす．これらの成分はたとえばカルシウムの量，タンパク質の量，あるいはさまざまな代謝生産物などの量である．

適切な変数の数について一般的に答えることは困難なので，状態の何に着目したいかを考えることが重要である．

生物システムの変数をすべて組み込むことは現実的に不可能であるし，変数が多いモデルが良いモデルとなるとは限らない．では，どのように変数を選択すればよいのだろうか．教科書では熱力学の成功を例に挙げている．

熱力学は，分子の速度や位置といった情報を含む非常にミクロな変数ではなく，温度や体積といった(分子オーダーと比べて)マクロな変数を用いて，我々が必要なスケールの状態を記述することに成功している．もちろん，よりミクロな理解が必要であれば，例えば分子オーダーのモデルが必要である．つまり，モデル化する時は，知りたい現象のレイヤーについて考えて，使用する状態変数を考える必要がある．

モデルをたてる時，記述したモデルが無視した多くの変数に依存しないことを期待しなければならない．つまり，無視した変数によらない「普遍的な」ふるまいがあることを期待して，対象について理解する必要がある*1．

当然，変数の組によって表現がよいかどうかの違いはある．しかし，はじめからどの変数あベストかの問いをたててそれに悩みすぎてもおそらく正解は得られない．むしろ，いま考える現象に応じて適当に選び，足りなければ，さらに別な状態変数を考える，不要なものがあれば消す，といった態度で望むほうが生産的であろう．

力学系とは

一般に，状態は時間とともに変動する．そして，その変化は状態空間上の点の軌跡として表現できる，と考えるのが力学系の考え方である．

我々が考えたい現象を十分に記述できているモデルを考えると，状態の変化の仕方はその状態変数で与えられ，各点でどの方向にどの速さを示した矢印で表現できることが期待できる．これを数学的に表現するのが力学系である．

つまり，状態空間と，その状態空間の各点で時間発展の規則が与えられた時，初期条件が与えられれば状態変化(軌跡)を考えることができる．変化の矢印が各点$(x_1, x_2, \cdots, x_k)$の関数$f_i (i = 1, 2, \cdots, k)$で与えられるので，微分方程式を用いて記述すれば

$$\frac{dx_i}{dt} = f_i(x_1, x_2, \cdots, x_k)$$

と書ける．ただし，$f_i$は向きを持つベクトルである．このとき，時間発展は$k$次元空間で完全に決定される*2．

ここでは，この式を解くことよりも，それが$k$次元状態空間で，どういう矢印を構成し，それにそってどういう状態の流れがあるかを理解することが重要である．

例1.3 状態が化学成分の濃度であれば，その変化は化学反応による濃度変化で与えられるであろう．たとえば$X_1, X_2, X_3$の3成分を考えその濃度を$x_1, x_2, x_3$とする．ここで一定量の成分$A$から$X_1$がつくられ，一定量の成分$B$から$X_2$がつくられ，そして成分$X_1$と$X_2$が反応して$X_3$がつくられ(その合成レートをrとする)． $X_3$が$r_d$のレートで分解するとする．反応は一方向に進むとする．成分3の増加率は濃度$x_1$と$x_2$の積に比例し， $$\frac{dx_1}{dt} = a - rx_1 x_2$$ $$\frac{dx_2}{dt} = b - rx_1 x_2$$ $$\frac{dx_3}{dt} = rx_1 x_2 - r_d x_3$$ となる．もちろん，反応がどれだけ生じるかは，もとは分子の衝突なので，ぴったり上の値になるわけではなく，平均として，である．

数学の準備

力学系の式が与えられた時，その時間発展を調べるにはそれぞれの状態変数の増減を調べる必要がある．これは，$f_i(x_1, x_2, \cdots, x_k)$で与えられる．そこで，$f_i(x_1, x_2, \cdots, x_k) = 0$となる空間(面)を定めておくことが重要で，これはヌルクライン(nullcline)と呼ばれる．

これら$k$個のヌルクラインの条件を同時に満たす点を固定点という．この点から出発すると，以降の時刻でその点から移動しない．

固定点では，それ以降，状態が変化しないことが分かった．では，固定点の近くではどういった振る舞いをするのだろうか．それを調べるには固定点付近の安定性について考える必要がある．これは，式を固定点の周りでテイラー展開して，(最低次で)線型化すれば判定できる．今，固定点を$(x^{*}_1, x^{*}_2, \cdots, x^{*}_k)$とし，固定点の定義より$f_j(x^*_1, x^*_2, \cdots, x^*_k) = 0 (j = 1, 2, \cdots, k)$を満たす．ここで，$x_i(t) = x^*_i + \delta x_i(t)$として，これを$f_j$に代入して，$x^*_i$の周りでテイラー展開して$\delta x_i(t)$の1次だけを残す．

$$ \frac{d\delta x_i(t)}{dt} = \frac{dx_i(t)}{dt}|_{(x_1, x_2, \cdots, x_k) = (x^*_1, x^*_2, \cdots, x^*_k)} $$ $$ = f_i((x^*_1, x^*_2, \cdots, x^*_k)) $$ $$ = (\delta x_1(t) \frac{\partial}{\partial x_1} + \delta x_2(t) \frac{\partial}{\partial x_2} + \cdots + \delta x_k(t) \frac{\partial}{\partial x_k}) \times f(x_1, x_2, \cdots, x_k) \quad(R(x^*_1, x^*_2, \cdots, x^*_k)は無視) $$

となる．ここで，$J_{ij} = \frac{\partial f_i}{\partial x_j}$は固定点周りのヤコビ行列である．ベクトルで書けば $$ \frac{d \delta \bm{x}(dt)}{dt} = \bm{J}\delta \bm{x}(t) $$ となる．

式は線型なので行列$\bm{J}$を対角化して解くことができる．対角化するための変換行列を$\bm{T}$として，$\bm{\Lambda} = \bm{TJT}^{-1}$が対角行列になるとし，その固有値を$\lambda_1, \lambda_2, \cdots, \lambda_k$とする．

ここで，$\bm{v} = \bm{T}\delta \bm{x}$とすれば，式に左から$\bm{T}$をかけ，$\bm{J}$と$\delta \bm{x}(t)$の間に$\bm{T}^{-1}\bm{T}$を挿入して，

$$ \frac{d \bm{v}(t)}{dt} = \Lambda \bm{v}(t) $$

を得る．ベクトルの成分ごとに書くと

$$ \frac{d v_i(t)}{dt} = \lambda_i v_i $$

となって，これを解くと$v_i(t) = v_i(0)\exp(\lambda_i t)$となる．ここで，

$$ \delta \bm{x} = \bm{T}^{-1} \cdot \bm{v} $$ $$ = \bm{T}^{-1} \cdot (v_1(t), v_2(t), \cdots, v_k(t))^\top $$ $$ = (\bm{T}^{-1}_1 \cdot \bm{v}, \bm{T}^{-1}_2 \cdot \bm{v}, \cdots, \bm{T}^{-1}_k \cdot \bm{v})^\top $$

であるから，$\delta \bm{x}$の各成分を見ると

$$ \delta x_i(t) = \bm{T}^{-1}_i \cdot \bm{v} $$ $$ = \sum_j T^{-1}_{ij}v_j(t) $$ $$ = \sum_j \exp(\lambda_j t) T^{-1}_{ij} v_j(0) $$ $$ = \sum_j \exp(\lambda_j t) T^{-1}_{ij} \sum_m T^{-1}_{jm} x_m(0) $$ $$ = \sum_{j, m} \exp(\lambda_j t) T^{-1}_{ij} T_{jm} \delta x_m(0) $$

となる．これから，$\delta \bm{x}(t)$の時間発展は$\exp(\lambda_j t)$の重ね合わせであわらせることがわかる．ここで，式は実数であるが，一般には固有値$\lambda_j$が複素数の場合もある．そのときは，複素共役の$\lambda_j^*$も固有値になり， $\lambda_j = a + i\omega (iは虚数単位，aは\lambda_jの実部，\omega は実数)$に対して

$$ \exp(\lambda_i t) = \exp((a \pm i\omega)t) $$ $$ = \exp(at)\exp(\pm i\omega t) $$

であるから，オイラーの公式より，$\delta x_i(t)$は$\exp(at)\cos(\omega t)，\exp(at)\sin(\omega t)$の重ね合わせで書けることが分かる．

固定点周りでの安定性を調べていたことを思いだすと，ある$j$に対して固有値$\lambda_j$の実部$a$が正のとき，$\delta x_i(t)$は指数関数的に増加するから安定でない．つまり，(線型)安定性を持つには，$\bm{J}$の固有値すべての実部が負であることが必要である．

以上のことから，固定点周りでの振る舞いは，固有値の実部の符号と，虚部で決まることが分かった．ここでは，2変数の場合を例にとって議論する．

実部と虚部がともに負の場合 : 直線的にどの方向からも固定点に向かう．
複素共役なペアの固有値で実部が負の場合 : 角速度$\omega$で固定点へ向かう．
ともに実数で，符号が異なる場合 : 固有値が負の方の固有ベクトルの方向からは固定点に引き込まれ，正である固有ベクトルの方向には固定点から離れる方向に向かう．
複素共役なペアで固有値で実部が正の場合 : (ii)と反対，実部が正なので回転しながら固定点から離れる．
ともに実数で，正の場合 : (i)と反対に，固定点から直線的に離れる．
実部が0で虚部が非零の場合 : 固定点の周りで回転する．

例1.4 より具体的で，また後の実験の例でも何度か現れるものとして，タンパク質による発現の制御を考えてみよう．ここで，2つの遺伝子が互いに抑制しあう場合を考えよう．具体的には遺伝子1が発現してつくられたタンパク質1が遺伝子2の発現を抑え，遺伝子2が発現してつくられたタンパク質2が1の発現を抑える場合である．遺伝子1からmRNA1がつくられ，それからタンパク質1(その濃度を$p_1$とする)がつくられ，遺伝子2についても同様にそれからmRNA2を経てタンパク質2がつくられ，その濃度を$p_2$とする．そこで状態空間は$(p_1, p_2)$になる(なお，間にあったmRNAの濃度を考えないで，タンパク質だけにしたのは，mRNAは合成，分解の時間スケールが早く，つねにタンパク質の濃度に比例するとしたためである．
このような早い変数の消去については1.12節で述べる．後述するように，この系はタンパク質1が発現し(‘オン’)2が発現していない(‘オフ’)安定状態とその逆の安定状態がある．このシステムは，電気工学などで用いられている．つまみでスイッチのオンオフを切り換えるシステムにならって，トグル(toggle)・スイッチと呼ばれている．
一方で，合成されたタンパク質はある割合で分解していく．簡単のためにそれの速度定数はちらも同じとしよう．時間の単位を適当にとれば，その定数を1とすることができる．ここで各遺伝子が発現し，そのタンパク質が合成される割合は，相手側のタンパク質，つまりタンパク質1ならタンパク質2，タンパク質2ならタンパク質1により抑制されるとする．具体的にはそれぞれのタンパク質の合成率が$\frac{\alpha}{1 + p_2^2}，\frac{\alpha}{1 + p_1^2}$であらわせるとする．ここで2乗で減っていくとしたのは，たとえばタンパク質が2つくっついて2量体をつくり，それが遺伝子のプロモーターにくっついてそれがmRNAがつくられるのを阻害するというように，量の2次で抑制が起こる場合を考えたからである(これは例であって，他にもいろいろな場合が考えられる)．すると，この場合，$p_1，p_2$の時間発展の力学系は $$ \frac{dp_1}{dt} = \frac{\alpha}{1 + p_2^2} - p_1 $$ $$ \frac{dp_2}{dt} = \frac{\alpha}{1 + p_1^2} - p_2 $$ であらわされる．以上はタンパク質がその合成に互いに「抑制」する場合であるが，促進する場合であれば，$\frac{\alpha}{1 + p_2^2}$のような減少関数のかわりに$p_2$に対する増加関数が用いられる．

次回は，トグル・スイッチの固定点周りの安定性から始める．

*1 教科書曰く，物理では熱力学，流体力学での成功例がある．
*2 決定論方程式と呼ばれる
*3 ヌルクラインは各関数$f_i$ごとに定まるものであることに注意せよ
*4 固有値は縮退していないとして，変換行列$\bm{T}$を用いて対角化できるとした
*5 0になる場合は次回以降で検討する

参考にした書籍やページ

細胞の理論生物学: ダイナミクスの視点から. 金子邦彦, 澤井哲, 高木拓明, 古澤力.
EMANの物理学. 多変数関数のテイラー展開
工業大学生ももやまのうさぎ塾. うさぎでもわかる微分方程式　Part11 対角化を用いた連立微分方程式の解き方と指数行列
予備校のノリで学ぶ「大学の数学・物理」. 【大学数学】線形代数入門⑭(対角化：重解がある場合)

感想

前日から準備を始めて徹夜で資料を作ったせいで，とても疲れた．次回で一章を終わらせて，早めに生物パートに取り掛かりたいと思っているが，はたして…

tax_free/iGEM TokyoTechで勉強会してきた

やったこと

イントロ

力学系とは

数学の準備

参考にした書籍やページ

感想