ローレンツアトラクタで「バタフライ効果」を測ったら、精度を10倍にしても予測できる時間は2.7しか伸びなかった
約9分で読めます

ローレンツアトラクタで「バタフライ効果」を測ったら、精度を10倍にしても予測できる時間は2.7しか伸びなかった

「蝶の羽ばたきが竜巻を起こす」というたとえは、気象学者エドワード・ローレンツが1963年に見つけた現象から来ている。彼が使っていた計算機で、途中経過を打ち直して計算を再開したら、まったく違う結果になった。打ち直した数字が小数点以下3桁で丸められていた、というのが原因だった。

そのとき彼が扱っていた方程式が、これだ。

dxdt=σ(y−x),dydt=x(ρ−z)−y,dzdt=xy−βz\frac{dx}{dt} = \sigma(y-x), \quad \frac{dy}{dt} = x(\rho - z) - y, \quad \frac{dz}{dt} = xy - \beta z

σ=10\sigma=10、ρ=28\rho=28、β=8/3\beta=8/3。3行しかない。ランダムな要素はどこにも入っていない。

「わずかな差が増幅される」ことは知っている。知りたいのは速さのほうだ。10億分の1のズレは、何秒後に目に見える差になるのか。そして、精度を上げれば予測できる時間はどれだけ伸びるのか。

実験に使ったコードの全文はGitHubに置いています。

まず、アトラクタを見る

4次のルンゲ・クッタ法で解いて、軌道を描く。

2つの渦を行き来しながら、同じ場所を二度は通らない。閉じた輪にならず、かといって発散もせず、有限の領域の中を永遠にさまよう。この形が蝶の羽に見えることから、バタフライ効果という名前がこちらの由来にもなっている。

10億分の1のズレが見えるまで

初期値の xx を 10−910^{-9} だけずらした2本目を用意して、並走させた。

左は2本の軌道の距離の時間変化を対数目盛で描いたグラフ。10のマイナス9乗から始まり、ほぼ直線的に(=指数的に)増えていき、t=23.8で距離1.0の破線を超える。右は2本の軌道を3次元で重ねた図で、蝶の形の同じアトラクタの上を青と赤の線が走っている

時間2本の距離
00.000000001
50.00000026
100.000035
150.0025
200.027
23.81.0

対数目盛でほぼ直線ということは、距離が指数関数的に増えているということだ。t=23.8t=23.8 で距離1.0——アトラクタ全体の大きさが40程度なので、この時点で「同じあたりを走っている2本」ではあるが、もう別の点だ。

動かして見るとこうなる。

ローレンツアトラクタ上を2本の軌道が走るアニメーション。序盤は青と赤の線が完全に重なっていて1本にしか見えないが、後半で赤い点が青い点から離れ始め、最後には反対側の羽を回っている。タイトルに現在の時刻と2点の距離が表示される

前半は完全に重なっていて1本の線にしか見えない。後半で急に分かれる。「じわじわ離れる」ではなく「しばらく何も起きず、あるとき急に」という見え方になるのは、指数関数がそう見えるからだ。

リアプノフ指数を、1回で測ってはいけなかった

この指数的な増大の速さが最大リアプノフ指数 λ\lambda で、ρ=28\rho=28 のローレンツ系では文献値0.906が知られている。上のグラフの傾きから求めると——0.760。合わない。

原因は測り方だった。1本の摂動ベクトルの伸び方は、アトラクタ上のどこにいるかで大きく変動する。1回の走行だけで傾きを取ると、その区間のたまたまの伸び方を拾ってしまう。

そこで、一定時間ごとに摂動の大きさを元に戻しながら伸び率を積み上げる方法(Benettin法)に変えた。

左は相関積分C(r)の両対数プロットで、きれいな直線になっており傾き2.050。右はリアプノフ指数の推定値が平均時間とともに収束していく様子で、最初は0.85から0.95の間で振れるが、時間が経つと文献値0.906の破線にぴたりと張り付く。下方に「1回だけ測ったときの0.760」の点線がある

0.9057。 文献値0.906とほぼ一致した。右のグラフのとおり、最初の100区間では0.889、500区間で0.917、2000区間で0.909と、行ったり来たりしながら収束する。

同じデータ、同じ現象でも、1回の観測から傾きを読むか、何度も測り直して平均するかで0.15違った。前者は0.906に対して16%低い。

精度を上げても、予測できる時間はほとんど伸びない

λ\lambda が分かると、実務的に一番効く量が出せる。初期誤差 δ\delta が許容誤差 Δ\Delta に育つまでの時間は

T=1λln⁡ΔδT = \frac{1}{\lambda}\ln\frac{\Delta}{\delta}

δ\delta の対数にしか依存しない。つまり、精度を10倍にしても予測時間は ln⁡10/λ\ln 10/\lambda しか伸びない。実測してみた。

横軸が初期のズレ(10のマイナス3乗からマイナス15乗まで、右へ行くほど精密)、縦軸が2本の距離が1.0に開くまでの時間の折れ線グラフ。3.6、8.7、18.5、23.8、28.5、33.9、34.0とほぼ等間隔に増えていき、直線的に見える

初期のズレ予測がもつ時間
10−310^{-3}3.6
10−710^{-7}18.5
10−1110^{-11}28.5
10−1510^{-15}34.0

精度を10倍にするごとに伸びるのは2.71(理論値 ln⁡10/λ=3.03\ln10/\lambda = 3.03)。10−310^{-3} から 10−1510^{-15} へ、1兆倍精密にして、伸びたのは10倍にも満たない。

倍精度浮動小数点の限界(およそ 10−1610^{-16})から始めても、34時間ぶんで別物になる。これ以上を望むなら、計算の精度ではなく方程式そのものが正しいことと初期状態の観測が必要になる。天気予報が2週間先で頭打ちになる理由が、この式の形の中にある。

アトラクタの次元は、素直に測れた

マンデルブロ集合の境界の次元は理論値2に対して1.09しか出せなかった。ローレンツアトラクタでも測ってみる。今回は箱を数えるのではなく、軌道上の点対のうち距離が rr 未満のものの割合 C(r)C(r) を使う相関次元という方法にした。

2.050。 文献値の約2.05とほぼ一致した(上の図の左)。両対数プロットが2桁半にわたってきれいな直線になっている。

同じ「フラクタル次元を測る」でも、マンデルブロ境界とは対照的にすんなり出た。アトラクタは軌道が実際に通る場所の集合なので、点を並べればそのまま標本になる。境界のように「1画素より細い枝」を追いかける必要がない。測りやすさは対象の複雑さではなく、標本の取り方で決まっていた。

羽の乗り換えは、一度も予測できなかった

最後に、素朴な疑問を試した。左右どちらの羽を回っているかを0と1で記録すると、記号の列ができる。この列の履歴から「次に乗り換えるか」を当てられるだろうか。

2万点ぶんで乗り換えは3.55%。つまり「乗り換えない」と言い続けるだけで96.45%当たる。

直前 kk 個の並びを文脈として、その文脈での乗り換え率が50%を超えたら「乗り換える」と予測する、という単純な予測器を組んだ。kk を1から16まで動かした結果は——

どの kk でも、乗り換えの再現率は0.0%だった。 50%を超える文脈がひとつも存在しない。予測器は最後まで「乗り換えない」としか言わなかった。

96.45%という正解率だけ見れば優秀に見える。中身は、何も予測していない。

まとめ

  • 10−910^{-9} のズレは t=23.8t=23.8 で距離1.0に育つ。増え方は指数的で、前半は完全に重なって見える
  • リアプノフ指数は1回の走行から傾きを取ると0.760。何度も測り直して平均すると0.9057(文献値0.906)
  • 精度を10倍にしても予測できる時間は2.71しか伸びない。10−1510^{-15} から始めても34時間ぶん
  • 相関次元は2.050(文献値約2.05)。マンデルブロ境界と違い素直に測れた
  • 羽の乗り換えは、履歴16個ぶんを使っても再現率0.0%。96.45%の正解率は「乗り換えない」と言い続けた結果

効いていたのは、式のどこかに ln⁡\ln が付いている、それだけだった。

予測できる時間は初期精度の対数でしか伸びない。だから精度への投資は、最初はよく効いて、すぐに効かなくなる。10−310^{-3} から 10−710^{-7} にすると15伸びたのに、10−1110^{-11} から 10−1510^{-15} にしても5.5しか伸びない。同じ「1万倍の改善」なのに。

仕事で見積もりの精度を上げようとするとき、私は暗黙にこれと逆のことを期待している気がする。調査を2倍丁寧にやれば、見通せる期間も2倍になるだろう、と。実際には、丁寧さは掛け算で効いても、見通せる期間には足し算でしか効かない。

そして羽の乗り換えのほうは、もっと直接的な戒めだった。96.45%当てる予測器を作って、乗り換えは一度も当てていない。知りたかったのは残り3.55%のほうだったのに、正解率という数字はそれを隠す。

関連記事

このブログでは、機械学習の仕組みを実際にコードで実装して確かめた実験をテーマ別に整理しています（拡散モデル・Transformer・強化学習・多様体学習・過学習など）。