集団遺伝学に現れる退化放物弓偏微分方程式の近似一般解の構成法
城西大学理学部天野
–男
(KAZUO AMANO)
アブストラクト
. 著者の目的は、集団遺伝学に現れる退化放物型偏微分方程
式の近似的な–般解の構成法を与え、かつその方法の計算機上へのインプリ
メンテーションの仕方を具体的に説明することである。著者の結果は、退化
する偏微分方程式へ数式処理的なアプローチを行う際の、理論的な基盤を与
えると同時に、偏微分方程式に対する数式処理の新しい可能性を示唆してい
る。
1.
はじめに
良く知られているように、
与えられた偏微分方程式が簡単な場合には、
われわれはその厳
密な
–
般解を求積法で求めることができる。
しかしながら、
このような事は極めて稀である。
現在では、
ほとんどの研究者は、偏微分方程式の厳密な
–
般解を求めることにはあまり関心が
なく、主に関数解析的手法による理論解の研究、 および数値解法の研究に興味が移ってしまっ
ているように見受けられる。
われわれの目的は、近似一般解の可能性を示すことである。
じっさい、
ある種の初期値問
題の近似解が、極めて単純な代数的表現をもつことを、 われわれは証明する。具体的には、
近
似解が初期条件のシンボリックな和差積によって表されることを証明する。 われわれの近似
般解の構成法は、古典的な解の構成法とはまったく異なっている。具体的には、数値
$-$
数式ハ
イブリッド法を用いる。
われわれの方法は、
LISP
と
$\mathrm{C}$言語に依存せざるを得ない。
というの
は、
求める近似一般解の公式のサイズが、
5.4
メガバイトを越えてしまうからである。
このよ
うな公式は、古典的な紙と鉛筆による計算には、
あまりに大きすぎる。
このノートの目的は、集団遺伝学に現れる退化放物型偏微分
(1.1)
$\frac{\partial u}{\partial i}=\frac{1}{2}\frac{\partial^{2}}{\partial x^{2}}(V(X)u)-\frac{\partial}{\partial x}(M(x)u)$$(t>0,0<x<1)$
,
に対する初期値問題の近似一般解を構成することである。
ここで、
$V(x)= \frac{x(1-x)}{2N},$
$M(x)$
は
$x$の多項式を表し、
かつ
$N\in\{1,2,3, \cdots\}$
は
effective population number
を表す
あるクラスの
$M(x)$
に対して、
Crow-Kimura
は
(1.1)
の幾つかの厳密解を得た
(cf. [2]).
しかしながら、彼らの方法は、
職人芸的なテクニックを用いており、一般化することは不可能
である。
このノートにおいてわれわれは、一般の
$M(x)$
に対して Y
(1.1)
の近似一般解を構成
2.
補助定理
補題 2.1. 自然数
$n$と
$C^{n+1}$
級の関数
$u(t, x)$
に対して、
$u(t+h, x+k)= \sum_{\nu=0}^{n}\frac{1}{\nu!}(h\frac{\partial}{\partial t}+k\frac{\partial}{\partial x})^{\nu}u(t, x)$
(2.1)
$+ \int_{0}^{1}\frac{(1-\theta)^{n}}{n!}(h\frac{\partial}{\partial t}+k\frac{\partial}{\partial x})^{n+}1u(t+\theta h, x+\theta k)d\theta$
,
補題
2.1
は、
テイラーの定理より従う。
微分作用素
$A$
を
(2.4)
$Au= \frac{1}{2}\frac{\partial^{2}}{\partial x^{2}}(V(X)u)-\frac{\partial}{\partial x}(M(X)u)$,
で定義する。良く知られているように、
(2.4)
は集団遺伝学において重要な役割を演じる
$(cf.$
,
for example, [2]
$)$.
じっさい、
Mendelian population
の中で、 互いに対立する
2
つの遺伝子
$\alpha_{1}$
と
$\alpha_{2}$が与えられたとき、
$\alpha_{1}$の遺伝子頻度の確率密度関数
$u(t, x)$
は、偏微分方程式
$\partial u$$\overline{\partial t}$
(2.5)
$-=Au$
in
$(0, NT)\cross(0,1)$
.
をみたす。
ここで、
$T$
は正の定数である。
補題 22.
(2.6)
$a(x)= \frac{x(1-x)}{2}$
,
$b(x)= \frac{1-2x}{2}-NM(x)$
,
$c(x)=-NM’(x)$
,
とおき、
かつ微分作用素
$B$
を
(2.7)
$B= \frac{1}{2}a(x)\frac{\partial^{2}}{\partial x^{2}}+b(x)\frac{\partial}{\partial x}+c(x)$.
で定義する。 このとき、
$v(t, x)$
が方程式
$\partial v$(2.8)
$-=Bv$
in
$(0, T)\cross(0,1)$
$\overline{\partial l}$の解であれば、
(2.9)
$u(t, x)= \exp(_{-}\frac{t}{2N})v(\frac{t}{N}, X)$
は
(2.5)
をみたす。
補題
22
の証明は、単純な直接計算なので、省略する。
補題 22 により、偏微分方程式
(2.5)
が
(2.8)
に帰着された。
(2.8)
は
(2.5)
よりも、比較
的扱い易い。
じっさい
Y
(2.9)
において
fixation term
$\exp(-t/2N)$
が既に分離されており
$\text{、}$補題
2.3.
$u(t, x)\in C^{4}([0, \infty)\cross[0,1])$
は、偏微分方程式
$\partial u$
$\overline{\partial t}$
(2.10)
–$=Bu$
$(t>0,0<x<1)$ .
の解とする。 このとき、
$0<h\ll 1,0\leq x\pm\sqrt{a(x)}h\leq 1,0\leq x+b(x)h^{2}\leq 1$
かっ
$t-h^{2}\geq 0$
であれば、
(2.11)
$u(t, x)=Lu(t, X)-h^{4}Ru(t, X)$
がなりたつ。
ここで、
$Lu(t, X)$
(2.12)
$= \frac{1}{6}u(t, x+\sqrt{a(x)}h)+\frac{1}{6}u(t, x-\sqrt{a(x)}h)$
$+ \frac{1}{3}u(t, x+b(x)h^{2})+\frac{1}{3}u(t-h^{2}, X)+\frac{h^{2}}{3}c(x)u(t, X)$
かつ
$Ru(t, X)$
(2.13)
$= \frac{1}{3}\int_{0}^{1}\{\frac{a^{2}(x)(1-\theta)^{3}}{12}(\frac{\partial^{4}u}{\partial x^{4}}(t, x+\sqrt{a(x)}\theta h)+\frac{\partial^{4}u}{\partial x^{4}}(t, x-\sqrt{a(x)}\theta h))$$+b^{2}(X)(1- \theta)\frac{\partial^{2}u}{\partial x^{2}}(t, x+b(X)\theta h^{2})+(1-\theta)\frac{\partial^{2}u}{\partial t^{2}}(t-\theta h2, X)\}d\theta$
.
証明.
(2.1)
より、
$u(t, x\pm\sqrt{a(x)}h)$
$=u(t, x) \pm ha^{1/2}(X)\frac{\partial u}{\partial x}(t, x)$
(2.14)
$+ \frac{h^{2}}{2}a(_{X})\frac{\partial^{2}u}{\partial x^{2}}(t, X)\pm\frac{h^{3}}{6}a3/2(X)\frac{\partial^{3}u}{\partial x^{3}}(t, X)$
$+ \frac{h^{4}}{6}a^{2}(x)\int_{0}^{1}(1-\theta)^{3_{\frac{\partial^{4}u}{\partial x^{4}}}}(t, x\pm\theta\sqrt{a(x)}h)d\theta$
.
従って、
$\frac{1}{2}u(t, x+\sqrt{a(x)}h)+\frac{1}{2}u(t, x-\sqrt{a(x)}h)$
(2.15)
$=u(t, x)+ \frac{h^{2}}{2}a(x)\frac{\partial^{2}u}{\partial x^{2}}(i, X)$(2.1)
より、
$u(t, x+b(X)h^{2})$
$\partial u$
(2.16)
$=u(t, X)+h2b(X)\overline{\partial_{X}}(t, X)$
$+h^{4}b^{2}(x) \int_{0}^{1}(1-\theta)\frac{\partial^{2}u}{\partial x^{2}}(t, x+\theta b(x)h^{2})d\theta$
かつ
$u(t-h^{2}, x)$
(2.17)
$=u(t, x)-h^{2} \frac{\partial u}{\partial t}(t, X)$$+h^{4} \int_{0}^{1}(1-\theta)\frac{\partial^{2}u}{\partial t^{2}}(t-\theta h^{2}, X)d\theta$
.
方、
(2.18)
$\frac{1}{2}a(x)\frac{\partial^{2}u}{\partial x^{2}}+b(x)\frac{\partial u}{\partial x}-\frac{\partial u}{\partial t}=-c(x)u$が
(2.10)
から従う。
(2.18)
を
(2.15)
$)$(2.16), (2.17) と合成することにより、
$\frac{1}{2}u(t, x+\sqrt{a(x)}h)+\frac{1}{2}u(t, x-\sqrt{a(x)}h)+u(t, x+b(x)h^{2})+u(t-h^{2}, X)$
$=3u(t, x)-h2C(X)u(t, X)$
(2.19)
$+h^{4} \{_{\frac{a^{2}(x)}{12}\int^{1}}0(1-\theta)^{3}(\frac{\partial^{4}u}{\partial x^{4}}(t, x+\theta\sqrt{a(x)}h)+\frac{\partial^{4}u}{\partial x^{4}}(t, x-\theta\sqrt{a(x)}h))d\theta$$+b^{2}(x) \int_{0}^{1}(1-\theta)\frac{\partial^{2}u}{\partial x^{2}}(t, x+\theta b(x)h^{2})d\theta+\int_{0}^{1}(1-\theta)\frac{\partial^{2}u}{\partial l^{2}}(t-\theta h^{2}, x)d\theta\}$
を得る。
I
補題
23
により、
$u(t, x)$
が
$C^{4}$級の有界な解であれば、
$u(t, x)$
は次のように書き替え
られる
:
$u(t, x)= \frac{1}{6}u(t, x+\sqrt{a(x)}h)+\frac{1}{6}u(t, x-\sqrt{a(x)}h)$
(2.20)
$+ \frac{1}{3}u(t, x+b(x)h^{2})+\frac{1}{3}u(t-h^{2}, X)+\frac{h^{2}}{3}c(x)u(t, x)+O(h^{4})$
この式を再帰的に用いることにより、
われわれは、退化放物型方程式
(2.10)
の近似一般解を得
る。
3.
近似一般解の構成
この節で、
われわれは
(2.10)
の近似一般解を構成し、 その誤差評価を与える
(
定理
3.1,
3.2).
さらに、
この結果が簡単な
symbolic list operation
に翻訳され得ることを、
実証する
(Reduce Program).
応用上の
–
般性を失うことなく、
われわれは
(3.1)
$b(x)\{$
$\geq 0$
in a neighborhood
of
$0$in
$[0,1]$
$\leq 0$
in a neighborhood
of
1 in
$[0,1]$
と仮定する。 この節を通して、
$u(t, x)\in C^{4}([0, \infty)\cross[0,1])$
は
(2.10)
の古典解をあらわし Y
かつ
$h$は十分小さな正の定数とする。
(2.20)
において、
$(h^{2}/3)c(X)u(t, X)$
を
$\frac{h^{2}}{3}c(x)(\frac{1}{6}u(t, x+\sqrt{a(x)}h)+\frac{1}{6}u(t, x-\sqrt{a(x)}h)$
$+ \frac{1}{3}u(t, X+b(X)h^{2})+\frac{1}{3}u(t-h2, X)+\frac{h^{z}}{3}c(X)u(t, x)+O(h^{4}))$
で置き換えれば、
$u(t, x)=(1+ \frac{n^{B}}{3}c(X))(_{\overline{6}}\perp u(t, X+\sqrt{a(x)}h)+_{6}\perp-u(t, x-\sqrt{a(x)}h)$
(3.2)
$+. \frac{1}{3}u(t, x+b(x)h^{2})+\frac{1}{3}u(t-h^{2}, x))+\frac{h^{4}}{9}c^{2}(x)u(t, X)+O(h4)$
.
が得られる。
ここで、
$|c^{2}(X)|=N^{2}|M’(X)|^{2}$
は集団遺伝学においては比較的大きな値になる
ので、
$(h^{4}/9)C^{2}(X)u(t, X)$
を
$O(h^{4})$
で置き換えるべきではない。
そこで
Y
(3.2)
に対しても
う
–
度同様の変形を行うと、
$u(t, X)=(1+ \frac{h^{p}}{3}c(X)+\frac{h^{\mathrm{I}}}{9}C2(X))(\frac{\perp}{6}u(t, x+\sqrt{a(x)}h)+u(t, x\overline{6}\perp-\sqrt{a(x)}h)$
(3.3)
$+ \frac{1}{\mathrm{q}}u(t, X+b(X)h^{2})+.\frac{1}{\mathrm{q}}u(t-h2, x))+O(h^{4})$
;
が従う。
この等式は、本節において重要な役割を演じる。
関数
$\epsilon(x)$を
[
$h$if
$h^{2}\leq x\leq 1-h^{2}$
or
$x=0,1$
(3.4)
$\epsilon(x)=\{$
$\frac{x}{\sqrt{a(x)}}$if
$0<x<h^{2}$
$\frac{1-x}{\sqrt{a(x)}}$if
$1-h^{2}<x<1$
.
で定義する。
ここで、
(2.6), (3.1),
(3.4)
より、
$0\leq x\leq 1$
に対して、
(3.5)
$x\pm\sqrt{a(x)}\epsilon(x),$
$x+b(x)\epsilon 2(X)\in[0,1]$
がなりたつ。
補題
3.1.
$0<h\ll 1$
に対して、
(3.6)
$h^{2}\leq x\leq 1-h^{2}$
$\Rightarrow$ $\frac{h^{2}}{4}\leq x\pm\sqrt{a(x)}h\leq 1-\frac{h^{2}}{4}$.
補題 3.2.
$0<h\ll 1$
に対して、
(3.7)
$\frac{h^{2}}{4}\leq y\leq 1-\frac{h^{2}}{4}$ $\Rightarrow$ $\epsilon(y)\geq\frac{h}{\sqrt{2}}$差分作用素
$M$
を
(3.8)
$Mf(t, X)=\{$
$(1+ \frac{\epsilon^{2}(_{X)}}{3}c(_{X)}+\frac{\epsilon^{4}(x)}{9}c^{2}(_{X}))$
$\cross(\frac{1}{6,\rceil}f(t, x+\sqrt{a(x)}\epsilon(x))+1\frac{1}{6}f(t, x-\sqrt{a(x)}\epsilon(x\backslash ))$
$+ \frac{\perp}{3}f(t, x+b(x)\epsilon^{2}(x))+\frac{\perp}{3}f(t-\mathcal{E}^{2}(X), X))$
if
$t\geq h^{2}$
$f(t, x)$
otherwise.
で定義する。
さらに、関数
$Pk(t, X;s, y),$
$k=0,1,2,$
$\cdots$を次のように定義する :
(3.9)
$p_{0}(t, x;S, y)=\{$
1
if
$(s, y)=(t, X)$
$0$otherwise,
(3.10)
$p_{1}(t, x;s, y)=\{$
$\frac{1}{6}(1+\frac{\epsilon^{2}(_{X)}}{3}\mathrm{c}(x)+\frac{\epsilon^{4}(_{X)}}{9}c^{2}(x))$
if
$(s.y)=(t, x\pm\sqrt{a(x)}\epsilon(x))$
$\frac{1}{3}(1+\frac{\epsilon^{2}(_{X)}}{3}c(x)+\frac{\epsilon^{4}(_{X)}}{9}c^{2}(x)\mathrm{I}$
if
$(s.y)=(t, x+b(x)\epsilon 2(x))$
$\frac{1}{3}(1+\frac{\epsilon^{2}(_{X)}}{3}c(x)+\frac{\epsilon^{4}(_{X)}}{9}c^{2}(x))$
if
$(s.y)=(t-\epsilon^{2}(x), x)$
$0$
otherwise
かつ
$pk+1(t, x;s, y)(k\geq 1)$ は、
$D=[0, \infty]\mathrm{x}[0,1]$
上の任意の関数
$f(t, x)$
に対して、
(3.11)
$\int_{D}f(s, y)pk+1(t, x;ds, dy)=\int_{D}Mf(s, y)pk(t, x;ds, dy)$
$($
$i.e.$
,
$\sum_{(s,y)\in D}f(s, y)p_{k+1}(t, x;s, y)=\sum_{)(s,y\in D}Mf(s, y)p_{k}(t, x;s, y)$
$)$
(3.9)
より、
(3.12)
$u(t, x)= \int_{D}p_{0}(t, x;dS, dy)u(s, y)$
.
$\epsilon^{4}(x)=O(h^{4})$
なので、
(3.3)
と
(3.10)
より、
(3.13)
$u(t, x)= \int_{D}p_{1}(t, x;ds, dy)u(s, y)+O(h^{4})$
.
さらに、
(3.3), (3.11), (3.13)
より、
$u(t, X)= \int_{D}p_{1}(t, x;dS, dy)u(S, y)+o(h^{4})$
(3.14)
$= \int_{D}p_{1}(t, X;ds, dy)(Mu(s, y)+O(h^{4}))+O(h^{4})$
$= \int_{D}p_{2}(t, x;ds, dy)u(S, y)+O(h4)2$
.
われわれは、
このような計算を帰納的に繰り返すことができる。
以上より、つぎの定理を得る。
定理
31.
$u(t, x)\in C^{4}(.[0, \infty)\mathrm{X}[0,1])$
は退化放物型方程式
$(\mathit{2}, \mathit{1}\mathit{0})$の古典解とする。
このと
き、任意の
$k=0,1,2,$
$\cdots$に対して、
(3.15)
$u(t, x)= \int_{0\leq s\leq}\iota,$
$0\leq y\leq 1p_{k}(i, x;ds, dy)u(S, y)+O(h^{4})k$
が成り立つ。
ここで、
$O(h^{4})$
は
$k$には依存しない o
(3.15)
において、
$u(s, y),$
$0\leq s<h^{2}$
を初期条件
$\phi(y)$で置き換えることにより、次の
定理が得られる。
定理
32.
$u(t, x)\in C^{4}([0, \infty)\mathrm{x}[0,1])$
は退化放物型方程式
(2.10)
の古典解であると仮定
し、 かつ
$u(t, x)$
は初期条件
$u(0, x)=\phi(x)$
$(0<x<1)$
を満足するものと仮定する。
このとき、任意の
$k\gg 1$
に対して、
$u(t, X)= \int_{0\leq s<h^{2}},0\leq y\leq 1p_{k}(t, x;ds, dy)\emptyset(y)+o(h2)$
(3.16)
$+O(1) \sum k_{0}(\frac{1}{3})^{l}(\frac{2}{3})k-l+O(h^{4})k$
が成り立つ。
ただし、
$k_{0}$は
$2t/h^{2}$
以上の最小整数とする。
ここで、
$O(1),$
$O(h^{2})$
と
$O(h^{4})$
は
$k$に依存しない。
したがって、
目的の近似一般解
(3.17)
$u(t, x) \sim\int_{0\leq s<h}2,0\leq y\leq 1p_{k}(t, x;d_{S}, dy)\emptyset(y)$
$(k\gg 1, h^{4}k\ll 1)$
定理
32
の証明
.
定理
3.1
より、
$u(t, X)= \int_{0\leq s<h^{2}},0\leq y\leq 1p_{k}(t, x;dS, dy)\emptyset(y)+o(h2)$
(3.18)
$+ \int_{h^{2}\leq s\leq 1},0\leq y\leq 1pk(t, x;ds, dy)u(_{S}, y)+O(h^{4})k$
.
$p_{k}(t, x;s, y)$
の定義を見れば、補題
3.1,
32
から
$| \int_{h^{2}\leq s\leq 1,0}\leq y\leq 1|pk(i,x;ds,dy)u(s,y)$
$\leq C_{1}\int_{h^{2}\leq s\leq}1,0\leq y\leq 1p_{k(y}t,$
$X;ds,$
$d)$
$\leq C_{2}\sum_{=l0}^{0}k(\frac{1}{3})^{f}(\frac{2}{3})^{k-f}$
が従う。
ここで、
$C_{i}(i=1,2)$
は
$k$に依存しない正の定数をあらわす。
I
線形結合
$\int_{0\leq S}\leq b,$
$0\leq y\leq 1p_{k(S}t,$
$x;d,$
$dy)u(_{S}, y)=q_{1}u(t_{1,1}X)+q2u(t2, x2)+\cdots+q_{n}u(t_{n}, X_{n})$
をリスト
$((q_{1}t_{1}x_{1})(q_{2}t_{2}x_{2}).
.
. (q_{n}t_{n}x_{n}))$
と同
–
視すれば、定理
32
の主要部分は次のようなリスト演算に翻訳される。
Reduce Program.
procedure
hybrid-method
(
$\mathrm{t},$ $\mathrm{x}$,
n);
begin;
list-in
$:=\{\{1, \mathrm{t}, \mathrm{x}\}\}$;
list-tmp
$:=$
$\{\}$;
while
$\mathrm{n}>0$do
$<<$
while
length(list-in)
$>0$
do
$<<$
tmp
$:=\mathrm{f}$irst
(list-in);
$\mathrm{p}$ $:=\mathrm{f}$
irst
(
$\mathrm{t}\mathrm{m}\mathrm{p}^{)}$;
$\mathrm{s}$ $:=\mathrm{f}$
irst
(rest
$(\mathrm{t}\mathrm{m}\mathrm{p}^{)})$;
$\mathrm{y}$ $:=\mathrm{f}$
irst
(rest (rest
$(\mathrm{t}\mathrm{m}\mathrm{p}^{)})$);
if domain-p(s)
neq
$0$then
$<<$
list-tmp
$:=$
cons
(
$\{\mathrm{q}/6,$ $\mathrm{s},$ $\mathrm{y}+\mathrm{s}\mathrm{q}\mathrm{r}\mathrm{t}(\mathrm{a}(\mathrm{y}))*\mathrm{h}\}$,
list-tmp);
list-tmp
$:=$
cons
(
$\{\mathrm{q}/6,$ $\mathrm{s},$ $\mathrm{y}-\mathrm{s}\mathrm{q}\mathrm{r}\mathrm{t}(\mathrm{a}(\mathrm{y}))*\mathrm{h}\}$, list-tmp);
list-tmp
$:=$
cons
(
$\{\mathrm{q}/3,$ $\mathrm{s},$$\mathrm{y}+\mathrm{b}(\mathrm{y})*\mathrm{h}**2\}$
, list-tmp);
list-tmp
$:=$
cons
(
$\{\mathrm{q}/3,$ $\mathrm{s}-\mathrm{h}**2,$ $\mathrm{y}\}$, list-tmp)
$>>$
else
list-tmp
$:=$cons
(
$\{\mathrm{p},$$\mathrm{s},$ $\mathrm{y}\}$
,
list-tmp);
list-in
:
rest(list-in)
$>>$
;
list-in
$:=$
list-tmp;
list-tmp
$:=$
$\{\}$;
$\mathrm{n}$ $:=\mathrm{n}-\perp$$>>$
;
return list-in
end;
4.
数値実験
(3.17)
を用いて、
われわれは次の
2 っの問題を解く :
(4.1)
$\{$$\frac{\partial u}{\partial t}=\frac{1}{4N}\frac{\partial^{2}}{\partial x^{2}}(x(1-x)u)$
$(0<i\leq 4N, 0<x<1)$
$u(0, x)\sim\delta(x-0.5)$
,
(4.2)
$\{$$\frac{\partial u}{\partial t}=\frac{1}{4N}\frac{\partial^{2}}{\partial x^{2}}(x(1-x)u)$
$(0<t\leq 4N, 0<x<1)$
$u(\mathrm{O}, x)\sim\delta(x-0.1)$
.
(4.1)
と
(4.2)
の近似解のグラフは Y
Figure 1
と
Figure
2 のようになる。
われわれの近
似解は木村の厳密解と酷似しており、木村の結果と同様の締釘にわれわれも到達する
$(cf,$
$[2$
,
Chapter 8]).
’.
$\cdot$.
.
$\cdot$ .Figure
1
$u$
