佐藤 總夫:自然の数理と社会の数理 Ⅱ

作成日 : 2013-01-23
最終更新日 :

概要

微分方程式の使われる分野として、美術、経営、軍事、政治、生理学、疫学、生態学などから例を挙げて、 懇切丁寧に説明している。 第 I 巻に続き、この巻では疫学や生態学からの例を挙げている。 著者はまえがきで、I, II, III の三分冊にした、 と書いているが、第 III 巻は未だに出版されていない。

感想

名著である。 そして著者は既に亡くなっている。 残念なことだ。 (2009-01-12)。

軍拡モデル

著者は、第1話「軍縮はなぜ軍拡に変わるのか」で、リチャードソンのモデルを紹介している。 ルイス・フライ・リチャードソンはイギリスの数理物理学者である。リチャードソンは数値解析による天気予報、フラクタルの概念などで現在知られているが、 同書では戦争の数学的解析のモデルに関する業績で紹介されている。 さて、リチャードソンのモデルを途中まで著者に沿って紹介しよう。

まず、自国の軍備の増加率が相手国の軍備の大きさに比例する、と仮定する。ある時点での `X` 国の戦力を `x` で、 `Y` 国の戦力を `y` で表す。 微小な時間 `dt` における `x`, `y` の増加分をそれぞれ `dx`, `dy` とすると、適当な比例係数 `k` を用いて、次の連立微分方程式を得る。

`(dx)/(dt) = ky, quad (dy)/(dt) = kx`

のちの便宜のために、この連立微分方程式を次のように書き換える。

`{(x' = ky), (y' = kx):}`

この正の比例定数 `k` を以下防衛係数と呼ぶ。相手国の戦力が単位量だけ増えたとき、 それに応じて自国の戦力を為政者がどれぐらいの割合で増やすかを表す係数である。

この連立微分方程式を解くと、`t -> oo` で `x -> oo, y -> oo` となり、果てしない軍拡競争となる(その結果は戦争となる)。 しかし、実際には軍事といえども予算に限界があるのでこうはならない。モデルをどのように修正すべきか。

まず、自国と相手国の防衛係数が同じというのは現実とは異なるので、自国と相手国の防衛係数を異なるものとして、それぞれ `k, l (k != 0)` とする。 そして、自国の増大した防衛力の結果として、それ以上の防衛力を抑制すべきである。したがって、次のようなモデルが考えられる。

`{(x' = ky - alpha x), (y' = kx - beta y):}`

ここで `alpha, beta` は正の係数で、自国の戦力を維持するための経済的負担に対する疲労を表すものであるから、これを疲労係数と呼ぶ。

さらに、ナショナリストたちの野望の程度や自国民が抱く平和の安心の程度を表す因子を上の右辺の定数項`g, h`として導入する。これを不平因子と呼ぶ。 これを加えた式は次のとおりである。

`{(x' = ky - alpha x + g), (y' = kx - beta y + h):}`

不平因子が正のときは相対的に自国の軍備に不安を抱くとき、負のときは自国の軍備に安心するときと考える。

この最後のモデルがリチャードソンのモデルとして、同書の第2章以下で考察されている。 今、この文章を書いているのが 2016 年 11 月 13 日、 合衆国次期大統領にトランプが選ばれたばかりの時期であることに、 妙な符合を感じる。

糖尿病モデル

糖尿病とは、代謝異常の病気である。人間は食物を摂取して栄養を得て代謝によってエネルギー源に変える。 この代謝の過程で血液中のブドウ糖は一定になるように調節されているが、 この調節の過程に異常が生じるのが糖尿病である。

血糖を調節するシステムは複雑である。このシステムをモデル化する試みは多くあり、 その中でアッカーマン(Ackerman)は 4 個のパラメータで数理モデルを作るに至った。 以下はそのモデルをさらに簡便化して示す。 なお、このモデルは OGTT (後述)によって正常な人と、 軽症糖尿病または前糖尿病人とを識別する診断に限り役に立つ。ここで、 OGTT とは経口ぶどう糖負荷試験(Oral Glucose Tolerance Test)のことである。 もう一つ付言すると、以下はモデルの正当性を示すために血液中の血糖値を6回以上測定することを前提としているが、 臨床的には高々4回である。

モデルの導出過程は省略し、結果だけを述べると、t 時間後の血糖値は次の通りである:
`G(t) = G_0 + R exp (-alpha t) sin omega t` 。
`omega_0 = sqrt(omega^2 + alpha^2), T_0 = (2pi) / omega`

最初にお断りするが、同書では、指数関数と三角関数の積の係数は `R // omega` であり、 本来のモデルはこれが正しいが、モデルの当てはめに限って話を進めるため、`R // omega` を改めて `R` とする。

さて本書で取り上げられている3例のパラメータの推定値を次の表に示す。

表1 パラメータの推定値
G0 (mg/100ml)R (mg/100ml) α (min-1×103) ω (min-1×103) ω0 (min-1×103) T0 (min)
A67.810013.941.543.8 144
B78.71,60649.05.8349.3127
C86.42054.8216.717.4361

同書によれば、T0が 4時間超であれば軽症糖尿病患者、 4時間未満であれば正常人、という判断が示されている。

;

;

;

私の関心事は、グラフの外形からこれらのパラメータの概算値が求められないか、 ということである。というのは、最小二乗法の計算を行なうときの精度にかかわるからである。 たいていは初期パラメータから反復計算によって最適なパラメータを計算すれば、 初期パラメータの概算値が適切であれば少ない時間で妥当なパラメータが得られるが、 概算値が不適切であれば収束に時間がかかったり、収束したパラメータが不適切な値にいくことも考えられるからだ。 A 、B、C で共通なのは、初期値 G(0) G_0 として考えてよいだろうが。 他のパラメータは A、B、C それぞれのタイプで概算値の求め方は異なるだろう。

まず、A の形を考える。周期も明確だし、 包絡線によって指数関数の係数を決められそうだ。 さて、周期は 180 分のように見える。だから、`2 pi = 180 omega ` を満たす `omega` とすればよい。このとき、`omega = 0.0349 = 34.9 * 10^-3` だから、 表1の `omega` と10%の誤差の範囲内である。だから、これを使えばいいだろう。 残りの `R` と `alpha` はこの形のままだと難しい。そこで、`G(t)` を `t` に関して微分した式
`G'(t) = R exp (-alpha t) (-alpha sin omega t + omega cos omega t)` 。
を考える。この式が 0 になるのは右のカッコ内がゼロになるときで、このときグラフの傾きは 0 である。 従って、グラフが極値をとるような t を読み取って omega といっしょに計算すればよい。
グラフをみると、そのような t は 30 分のようだ。
`-alpha sin omega t + omega cos omega t = 0`
`alpha = (omega cos omega t) / (sin omega t) = omega / (tan omega t)`
`omega = 0.0349, t = 30 ` を代入すると、
`alpha = 0.0349 / (tan 0.0349 * 30) = 0.020`
表1の `alpha` (=0.0139) の2倍弱あるが、仕方ない。
最後に R は、グラフから、t = 30 のとき、G(t) = 130 であることを利用して求める。
なお、G(0) もグラフから 70 と読み取る。
`130 = 70 + R exp(-0.02 * 30) sin ( 0.0349 * 30)`
`R = (130 - 70)/ (exp(-0.02 * 30) sin ( 0.0349 * 30)) = 126`
表1と比べても、30 % の誤差ならばまあまあだろうか。

B はあとで考えることにして、C を考えることにする。周期は 180 * 2 = 360 (分) のように見えるので、 これをもとに `omega` が求められるはずだ(以下省略)。

問題は B である。この形は周期 `omega` がよくわからない。 そこで、まずは `omega` を求めるのを一度諦める。ただし、 `omega` が小さいことを考えると、`sin omega t ~~ omega t` という近似が成り立つと考え、 このもとで、`G'(t)` を計算すると次のようになる:
`G(t) ~~ G_0 + R omega exp (-alpha t) t` 。
`G'(t) ~~ R omega (1 - alpha t ) exp (-alpha t) omega t` これから、極値をとる `t = t_0` を求めれば、`alpha` が求められることになる。 グラフから、`t_0` は 20 分あたりであろう。すると、`alpha = 50 ` ぐらいだろうか。 `alpha` がこれぐらい精度よく当たればよい。 またこの alpha から G(t)の式に戻って、積 `R omega` を求めることができる。 しかし、このままでは `R` と `omega` を分離して求めることができない。 ここで、試みに `R omega = const` の条件でグラフを描いてみた。

;

`omega = 6` と `omega = 12` ではそれほど外見が違いがない。 仕方がないので、このような外形であれば、`omega = 6` と決め打ち、 グラフから求めた`R omega = const` の定数から `R` を逆算するしかないだろう。

これまで述べてきたのはあくまで反復法の初期値をどうやってそれらしく決めるかであるから、 多少最小二乗の値からずれても、収束さえうまくいけばいい。 このような、非線形モデルの最小二乗法によるデータのあてはめは、 ガウス・ニュートン法で述べる。

伝染病の伝播

第7話「伝染病の伝播」については、 感染症と微分方程式という記事で、 この内容を自分の勉強用に書き写したので、そこを参照されたい。

ミクシンスキーの演算子法

あとがきで、著者はこのように言っている。

第6話について一言.定数係数の線型微分方程式を解くだけならば, ラプラス変換を使うよりもミクシンスキーの演算子法がよいと思う. それが当節の常識であろう. (中略)本文でミクシンスキーの演算子法を説かなかったのは, 伝統に従ったまでのこと. それにアッカーマンたちの論文・著作に合わせることなどで他意はない. (中略)ところが,初心者にこの演算子法を教える場合, もっともめんどうなところはティッチマーシュの定理の証明であった.証明なんかあとまわし, (中略)わりきって突進しても, この演算子法を十分に使うことはできる.ただ,数学畑に育った人はそうもゆくまい. しかし,幸いにも,最近になって特別な場合ではあるが, その定理の証明はいちじるしく簡素化された.(後略)

このティッチマーシュの定理について同書に書かれていないのは残念である。 ティッチマーシュの定理(あるいはティッチマーシュのたたみ込み定理)とは、 関数 `f(t)` および `g(t)` が `int_0^t f(t-s)g(s)ds -= 0` を満たすならば `f(t) -= 0` または `g(t) -= 0`.すなわち
`f * g = 0` ならば `f = 0` または `g = 0` という定理である。

吉田耕作と加藤敏夫の 「応用数学」 には `f(t), g(t)` ともに `[0, oo)` で積分可能なとき、すなわち、
`int_0^oo "|"f(t)"|" dt < oo, quad int_0^oo "|"g(t)"|" dt < oo`
に限定した証明がのっている。 この限定が、著者の言う「特別な場合」で、 このときの証明が「簡素化された」証明にあたるのかどうかは不明である。

なお、ティッチマーシュの定理の本格的な証明は、 京都大学数理解析研究所の数学入門公開講座 バックナンバー(講義ノート) (www.kurims.kyoto-u.ac.jp) の 1988 年 8 月 2 日 - 8 月 11 日(第11回)の松浦 重武先生による「演算子法の話」の資料にある。

誤植

場所備考
p.70 下から 2 行目(i) から (vi) まで(i) から (iv) まで
p.202 図5 のグラフの正の解 `zeta_1``zeta_2`
p.202 下から 3 行目 `z_-oo = zeta_1``z_-oo = - zeta_1`
p.206 (31)式`y_0/N = 1 - (N/rho)^-1(1-log{:N/rho:})` `y_0/N = 1 - (N/rho)^-1(1+log{:N/rho:})`

数式記述

このページの数式は ASCIIMathML で記述している。

書誌情報

書 名自然の数理と社会の数理 II
著 者佐藤 總夫
発行日1987 年 6 月 10 日(第 1 版第 1 刷)
発行元日本評論社
定 価3400 円(本体)
サイズA5 版 ページ
ISBN
NDC

まりんきょ学問所数学の本 > 佐藤 總夫:自然の数理と社会の数理 II


MARUYAMA Satosi