カオス時系列データのべき乗則の確認(python)



スポンサードリンク

これまでに、カオス時系列データに関する勉強として、レスラー方程式やローレンツ方程式を扱ってきましたが、その中で得られた時系列データのデータの振れ幅(x(t)-x(t-1))について何となく気になりました。

これまでのカオス時系列データに関する記事は以下です。

>>>ローレンツ方程式をRunge-Kutta法(python)で計算

>>>Euler法とRunge-Kutta法のカオス時系列データの比較

今回は得られた時系列データが理論にあっているか(べき乗則)の確認として、べき分布を表示してみました。

念のため、今回計算に用いた時系列データは以下のグラフです。

ついでに時系列データはローレンツ方程式をRunge-Kutta法(python)で計算の記事にあるコードを用いれば、簡単に計算できます。

べき乗則とは?

統計モデルの一つ。もっとも一般的な定義は
\begin{eqnarray}
f(x)=ax^{k}+o(x^{k})
\end{eqnarray}
oはランダウの記号、kはスケーリング指数です。

両対数グラフにすると、線形になるのが特徴。

ついでに片対数グラフで線形になるのは、指数関数です。

そしてこのべき乗則は自然現象や社会現象でもみられます。有名なところで言えば、株式市場の株価や地震の発生回数とその規模などもこれに従います。

なので、例えば株式市場で言えば、ちっちゃな下落や上昇はたくさんあるが歴史的な株価の暴落は数える程しか起きないと言うのがそうであります。

地震で言えば、震度1や2の地震は頻繁に発生するが、震度6や7の地震はか数えるほどしか起きないということです。

要は、事象の発生回数とその規模の相関関係を示す法則と言えます。

計算方法

結果を表示する前に今回の計算の流れ(方法)を述べます。計算に用いたコードはこの記事の最後に載せます。

1、Δt間での変化量ごとにグルーピングする(x(t)-x(t-1)の大きさごとにいくつかに分類)。地震の例で言えば、これが震度1、震度2、震度3のような大きさの定義になります。

2、Δtの分け方については、基本は1ステップとし、比較のためにもう少し変化の勾配がきつくなるようにΔtを設定した。

3、全ての時系列に対して行い、規模ごとに分類し配列にぶち込む。

4、Δx/Δt=1,1.1の時を最小の規模(地震で言えば震度1)とし、その1.5倍を次の規模(震度2)、2倍(震度3)、3倍と増やしていき、最小規模の6倍の規模までカウントするようにします。

5、なお今回の計算では、Δx/Δt、事象規模の間隔(1.5倍、2倍…等の分け方)については直感的に決めました。もしかしたら最適な決め方があるのかもしれません。

上の事象規模の分類の仕方ですと、最小規模の1.1〜1.4倍の事象も1.5倍の事象として分類していることになります。

なので精度を挙げるなら、ここの間隔(地震で言えば震度の分類)をより細かくするべきでしょう。


スポンサードリンク

べき分布の結果

以下に計算結果を示します。


Δt=1と1.1の時のデータです。もう少し線形なグラフが得られるかなと予測したのですが、全然違いました。ガウス分布の両対数グラフに近いかもしれません。

上でも述べたように、発生規模の分類の仕方に問題があるかもしれません。

あとローレンツ力学系はカオス系なので、初期値鋭敏性を持っています。

なので違う初期条件や計算時間をもっと長くしたりするなどしたら、結果はまた変わってくるかもしれません。

ついでに片対数(縦軸対数)グラフです。

こっちの方が線形ですね笑

ついでに、両軸とも対数にしていないグラフは以下です。




スポンサードリンク

計算に用いたコード

先に謝ります。ひどいです。いや本当にひどいコード。

配列の作成、Δtの場合分けをもっとうまくやればこんなに長くなりません。

奮闘しましたが、書けなかったのでなんとも原始的なやり方になってしまいました。

事象の規模の分類のところもどんだけelif使うんだよって感じですが笑

もう少しpython 勉強します笑。

べき分布の計算をする際に役立ててもらえれば幸いです。

記事が役に立ったらシェア!

投稿者:

中村 俊

中村 俊

1993/09/04生まれ。機械系大学院を休学し、ベンチャーでインターンしている最中。直近では、デカルトの「方法序説」に感銘を受けた。 趣味:読書、web開発の勉強、異分野の論文読んだり、記事書いたり。 最終的には経営者か研究者になりたい。

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です

CAPTCHA