ネガティブ・フィードバックの 振動子をつくる

初期値

(lacI, cI, tetR, LacI, CI, TetR) = (0.2, 0.3, 0.1, 0.1, 0.5, 0.4) パラメータ

演習 9: ネガティブ・フィードバックの

function repressilator( input_args ) time = 0.001:1:200;

s0 = [0.2, 0.3, 0.1, 0.1, 0.5, 0.4]; % Initial values param = [20, 0, 0.2, 2]; % Constants

[t,time_course] = ode15s(@(t,s) ODE(t,s,param),time,s0);

plot_time_course(t,time_course);

plot_phase_plane(time_course);

end

MATLAB コード例 (1/2)

function dsdt = ODE(t,s,param)

lacI = s(1); cI = s(2); tetR = s(3);

LacI = s(4); CI = s(5); TetR = s(6);

alpha = param(1);

alpha_zero = param(2);

beta = param(3);

n = param(4);

dsdt(1,:) = dsdt(2,:) = dsdt(3,:) = dsdt(4,:) = dsdt(5,:) = dsdt(6,:) = end

MATLAB コード例 (2/2): 空欄を埋めなさい

正弦波状の振動を示す

function dsdt = ODE(t,s,param)

lacI = s(1); cI = s(2); tetR = s(3);

LacI = s(4); CI = s(5); TetR = s(6);

alpha = param(1);

alpha_zero = param(2);

beta = param(3);

n = param(4);

dsdt(1,:) = alpha_zero + alpha / ( 1 + CI^n ) - lacI;

dsdt(2,:) = alpha_zero + alpha / ( 1 + TetR^n ) - cI;

dsdt(3,:) = alpha_zero + alpha / ( 1 + LacI^n ) - tetR;

dsdt(4,:) = - beta * ( LacI - lacI );

dsdt(5,:) = - beta * ( CI - cI );

dsdt(6,:) = - beta * ( TetR - tetR );

end

解答例

• Plan A: ヌルクライン、ベクトル場 

• 計算がさほど難しくない

• ただし2次元限定

• Plan B: 固有値を調べる

• 計算が複雑になる

• しかし多次元でも使える

• Repressilator は 6 変数なので Plan B で行く

Hopf 分岐の存在を示すには

• 6変数のヤコビ行列

• 固有値を求める

• 固有値の実部が正（不安定）となる条件

固有値を調べて Hopf 分岐点の条件を導く

∂m˙ ₁

∂p₃ = −αnpⁿ₃⁻¹

(1 + p²₃)² = X

p

₃

= α

1 + p

ⁿ₃

+ α

₀

(β + 1) ² (2X + 4) − 3βX ² < 0

J =











−1 0 0 0 0 X

0 −1 0 X 0 0

0 0 −1 0 X 0

β 0 0 −β 0 0

0 β 0 0 −β 0

0 0 β 0 0 −β











ただし

p3は、方程式 の解

• パラメータβを振り、下の条件を満たす点を境にして、振動が発生もし くは消失することを確認しなさい。

演習 10: Hopf 分岐が起きることを確認しなさい

(β + 1) ² (2X + 4) − 3βX ² < 0

p

₃

= α

1 + p

ⁿ₃

+ α

₀

を満たします(∵固定点の座標であるため)。上で α=10, α0=0, n=2 と定義した に現れる p3は、方程式

ただし、 β以外のパラメータはα=10, α0=0, n=2とします。

X = − αnp

ⁿ₃⁻¹

(1 + p

²₃

)

²

また、

出力例

β=0.13 β=0.14

(β+1)

²

(2X+4)-3X

²

= 0.0231 > 0 (β+1)

²

(2X+4)-3X

²

= -0.0355 < 0

• Hill係数の増加に伴って不安定解の領域が広がる

• 漏れ(α0)がΚmに近づくと不安定領域は狭まる

• 結合のエネルギーでオペレータ占有率を予測

• 安定性について同様の結果

論文本文より

発展課題

紙と鉛筆で解く Toggle switch と Pitchfork 分岐

発展課題： Toggle switch と Pitchfork 分岐

• 双安定性を示す人工遺伝子回路のモデル

Gardner et al. (2000)

 



˙

x = a

1 + y ² − x

˙

y = a

1 + x ² − y

LacI (x) λCI (y)

x y

Pitchfork 分岐

a=2 x a>2 x

y y

分岐図が熊手のように 見えるので ”pitchfork”

Stable

Unstable

y

ẏ = 0

ẋ = 0

発展課題 : ヌルクラインとベクトル場を描く

1. Toggle switch モデルのヌルクラインを相平面上に描きなさい。

固定点の数はa の値によって1 個または3 個になりますが、ど ちらの場合でも構いません。

2. 「1」の結果にベクトル場の概略を描き加え、固定点が1 個の 場合はStable Node に、3 個の場合はStable Node 2 個とUnstable

Node 1 個になると予想できることを確かめなさい。

発展課題 : 固定点を求める下準備

3. y についてのヌルクラインの式をx についてのヌルクラインの 式に代入し、固定点の x 座標が x⁵ - ax⁴ + 2x³ - 2ax² + (1 + a²)x -

a = 0 を満たすことを示しなさい。

4. 「3」で示した固定点 x 座標に関する条件式 x⁵ - ax⁴ + 2x³ - 2ax² + (1 + a²) x - a = 0 を因数分解すると (x³ + x - a)(x² - ax + 1) = 0 となることを確認しなさい。

発展課題 : 固定点の数をしらべる (1/2)

5. 一般に、x³ + ax + b = 0 の形をした3 次方程式の解の性質は、

判別式 D = - 4a³ - 27b² で予測することができる。D > 0 ならば 3 つの実数解、D = 0ならば重解、D < 0 ならば1 つの実数解と 2 つの虚数解を持つ。この知見を利用して、因数分解によって 得られた項 x³ + x - a についての方程式 x³ + x - a = 0 には1 つの 実数解（ヌルクラインの交点）が常に存在することを確認し なさい。（ヒント：a は速度定数なので常に正）

発展課題 : 固定点の数をしらべる (2/2)

6. 因数分解によって得られたもう一つの項である x² - ax + 1 についての方

程式 x² - ax + 1 = 0 の解の判別式を求め、a > 2 ならば実数解2つ、a = 2

ならば実数解1つがそれぞれ存在し、a < 2ならば実数解が存在しないこ とを確認しなさい。

7. 「5」「6」の結果をまとめると、固定点 x 座標に関する条件式 (x³ + x -

a)(x² - ax + 1) = 0 より、 (x³ + x - a)の部分に由来する固定点が常に1つ存

在し、 (x² - ax + 1) の部分からはa = 2を境にして固定点が0〜2個生じ る。Toggle switchモデルにおける 0 < a < 2 、a = 2 、 a > 2 それぞれの場 合の固定点の個数を答えなさい。 a = 2 の場合の個数に注意すること。

発展課題 : 固定点の安定性 (1/2)

8. Toggle-switchモデルのヤコビ行列が次の形になることを示しな

さい。

9. a > 2 のとき、 x² - ax + 1 = 0 の解として得られる固定点につい て次の問いに答えなさい。

(i) 固定点の座標が であることを示 しなさい

(ii) ヤコビ行列の2つの固有値がともに負であることを確認しなさい

! − 1 −

_(1+y^2ay2₎²

−

_(1+x^2ax2₎²

− 1

"

!a ±√

a² −4

2 , a∓√

a² −4 2

"

発展課題 : 固定点の安定性 (2/2)

10. x³ + x - a = 0 をの解として得られる固定点 x 座標は

となる。この x について次の問いに答えなさい。

(i) a = 2 のとき、 x = 1 である。また、x は a に対して単調に増

加する。これらの事実に基づき、0 < a < 2 のとき 0 < x < 1

および a > 2 のとき x > 1 であることを示しなさい。

(ii) 0 < x < 1 のとき（すなわち 0 < a < 2 のとき） Δ > 0 であるこ とおよび、 x > 1 のとき（すなわち a > 2 のとき） Δ < 0 であ

x = ³

!"

"

#a 2 +

$%a 2

&2

+ '1

3 (3

− 1

3 ³ )

a

2 + *+_a

2

,2

+ +₁

3

,3

発展課題 : Pitchfork 分岐図を描く

11. 「7.」から「10.」までの結果により、固定点の数と安定性は下の 表の通りであるとわかった。この表に基づいて横軸にa、縦軸に 固定点のx座標をとり、分岐図の概略を描きなさい。

0<a<2

^x3 + x - a = 0

1 1 Node x 1

0<a<2

x² - ax + 1 = 0

0 1 Node x 1

a=2

^x3 + x - a = 0

1 1

-a=2

x² - ax + 1 = 0

1 1

-a>2

^x3 + x - a = 0

1 3 Node x 2 Saddle x 1 a>2

x² - ax + 1 = 0

2 3 Node x 2

Saddle x 1

実数解の数 総計

1 の解答 (1/2)

y = a 1 + x

²

ẏ = 0 のヌルクラインの式は  である。これを次のように微分し、

√1

x

3

y’

---y’’ - 0 +

y^! = 0 ·(1 + x²) − a ·2x (1 +x²)²

= − 2ax

(1 + x²)²

y^!! = −2a(1 +x²)² −(−2ax) ·2(1 +x²)·2x (1 +x²)⁴

= −2a(1 +x²){(1 +x²) −2x ·2x} (1 +x²)⁴

= −2a(1 +x²)(1 −√

3x)(1 +√ 3x) (1 +x²)⁴

増減表をつくる。

y

a

ẏ = 0

1 の解答 (2/2)

ẋ = 0 については、 y = x について対称な曲線が得られる。よって

2つのヌルクラインは下の図のようになる。

x y

a

ẏ = 0 ẋ = 0

a

y = x

2 の解答

x

y

^{ẋ = 0}

a x

y

ẏ = 0

a

y

a

ẋ = 0

3 の解答

y = a 1 + x

²

x = a

1 + y

²

に を代入して、以下のように式変形。

x = a

1 + !

a 1+x²

"

2

= a(1 + x

²

)

²

(1 + x

²

)

²

+ a

²

(1 + x

²

)

²

x + a

²

x = a(1 + x

²

)

²

x + 2x

³

+ x

⁵

+ a

²

x = a + 2ax

²

+ ax

⁴

x

⁵

− ax

⁴

+ 2x

³

− 2ax

²

+ (1 + a

²

)x − a = 0

4 の解答

(x

³

+ x − a)(x

²

− ax + 1) = 0

x

⁵

− ax

⁴

+ x

³

+ x

³

− ax

²

+ x − ax

²

+ a

²

x − a = 0

x

⁵

− ax

⁴

+ 2x

³

− 2ax

²

+ (1 + a

²

)x − a = 0

5-7 の解答

5. x³ + x - a = 0 の判別式を計算すると、D = - 4 - 27a² < 0 となる。よっ て、この部分は常に1つの実数解と2つの虚数解を持つ。

6. x² - ax + 1 = 0 の解の判別式を計算すると、D = a² - 4 となる。a > 2 のとき、 D > 0 となるので実数解2つ。a = 2 ならば D = 0 より実数解 1つ。 a < 2 ならば D < 0 より実数解なし。

7.

0<a<2

^x3 + x - a = 0

1 1

0<a<2

x² - ax + 1 = 0

0 1

a=2

^x3 + x - a = 0

1 1

a=2

x² - ax + 1 = 0

1 1

x³ + x - a = 0

1

実数解の数 総計

3次の部分も2次 の部分もともに x=1が解

8 の解答

 



˙

x = a

1 + y

²

− x

˙

y = a

1 + x

²

− y より、

∂

∂x x ˙ = − 1 ^∂

∂y x ˙ = − 2ay (1 + y

²

)

²

∂

∂y y ˙ = − 1

∂

∂y x ˙ = − 2ax (1 + x

²

)

²

! − 1 −

_(1+y^2ay2₎²

1

"

よって

9-(ii) の解答

∆ =1 − 4a²xy

(1 + x²)²(1 + y²)²

= 1 − 4a²xy

(1 + x² +y² + x²y²)²

xy = (a ± √

a² −4)(a ∓ √

a² −4)

4 = 1

より

∆ =1 − 4a²

(2 + x² + y²)²

= 1 − 4a²

(2 + ^4a2₄⁻⁸)²

= 1 − 4a² a⁴

= 1 − 4

> 0 (∵ a > 2)

τ<0, Δ>0, τ

²

-4Δ>0 より、

2 つの固有値はともに

負の実数

10 の解答

(i) まず、0 < a < 2 の区間の端点における x の値を求める。x の式に

a = 0 を代入すると、x = 0 、 また、問題文より a = 2 のとき、 x

= 1 である。x は a に対して単調に増加するので、0 < a < 2 の区間 ではx は 0 より小さい値にならず、 1 より大きな値にもならな い。よって 0 < x < 1 。

(ii) この固定点は y=x 上にあるから、

∆ =1 − 4a²xy

(1 + x²)²(1 + y²)²

= 1 − 4a²x² (1 + x²)⁴ d∆

dx = 1 − 8a²x(1 + x²)³(x −1)² (1 + x²)⁴

これを微分すると となるので、Δはxの

11 の解答

x

2 a 1

実線は Node 、点線は Saddle

x³ + x - a = 0 由来の固定点

x² - ax + 1 = 0 由来の固定点 x² - ax + 1 = 0 由来の固定点

まとめ

• 連続値の微分方程式で離散的な分岐現象が起きるのはなぜか?

• 答え

• パラメータの値によってヌルクラインの交点の数が変わるから (Saddle-Node, Pitchfork)

• パラメータの値によってヤコビ行列の固有値の符号が変わるから (Hopf)

ドキュメント内 Underlying mechanisms of biochemical oscillations (ページ 102-133)