やりとりを繰り返すだけで貧富の差が生まれる？ ~「データの見えざる手」の分布~ part3
part2 の続き
さて、では、
$$Z(N,L) = \sum_{n_1 =0}^N \sum_{n_2=0}^N ... \sum_{n_L=0}^N \delta_{N,\sum_{l=1}^Ln_l}$$
を計算してみよう。
そのために、という式が、空間的に考えるとどんな図形になっているのかを考えよう。
すぐにはわからないので、L=2の場合を考える。
すると、これは
$$N=n_1 + n_2$$
を満たす場合の数なので、
下の図みたいに、の二次元で見た時に、(N,0)から(0,N)までを結ぶ直線上の格子点の数に一致する。
当然この場合はN+1通りである。
じゃあ、L=3の場合はどうなんだろうか？
こっちも図にしてみると、(N,0,0)、(0,N,0)、(0,0,N)の三点を結んで出来る正三角形上の格子点の数を数えることと同じだ。
さて、これで少しわかって来たぞ。
$$Z(N,L) = \sum_{n_1 =0}^N \sum_{n_2=0}^N ... \sum_{n_L=0}^N \delta_{N,\sum_{l=1}^Ln_l}$$
これは、L次元空間での、(N,0,0,...)、(0,N,0,....)、......、（...0,0,N）というような各軸上の原点からNの距離にある点を結んだ超平面に沿って、点を数えていくことに相当している。
・・・
・・・ここからどうしようか？？
一工夫が必要だ。
今回の計算では、Nが大きい場合を想定している。
Nが十分大きければ、図の格子点の数はこの超平面の面積とほとんど一致してしまうはずである。
そこで、ここは思い切って、和を積分で置き換えてしまおう。
さて、そうすると、
$$Z(N,L) \sim \int_0^Ndn_1\int_0^Ndn_2 \cdots \int_0^Ndn_L \delta(N - \sum_{l=1}^Ln_l)$$
と表せる。
ここで、 というのは、デルタ関数と呼ばれるもので、
の時のみ無限の値を持ち、それ意外が0、そして積分すると1になるという特殊な関数だ。（厳密には関数ではなく超関数だ。）
このデルタ関数を使えば、ジグザグに和を取っていく計算を、なめらかな関数の積分として表現しなおすことが出来る。
（もちろん、厳密には和の計算と積分の計算の答えは異なる。しかし、の場合には結果は変わらないだろう。ちゃんと証明してないけど。）
さて、では面積をどうやってもとめようか。
（実はここで、どうやってこれを求めれば良いかわからなくて、3日ほど費やした。きっと数学出来る人なら一瞬でやり方を思いつくのであろう・・・。）
超平面の面積は、超体積を微分すれば求まる
この超平面の面積を直接求めるのは大変だ。
しかし、実はこの超平面と軸が囲う超体積ならば簡単に求めることが出来る。
そして、面積はこの超体積を微分すれば求まるはずだ。
まずはL=2の場合でみてみよう。（この場合、超体積は面積、超平面の面積は線分の長さになる。）
この図の真ん中の直角二等辺三角形の面積は当然で、緑の線分の長さS(N)は
となる。
この緑の線分に並行な線分をT方向に積分しても、面積は求めることが出来る。
緑の線分の長さをTで表して、S(T)とおく。緑の線の中心の座標はなので、Tの長さは
となっているので、
$$\int^{\frac{N}{\sqrt{2}}}_0 S(T) dt = \frac{N^2}{2}$$
と書ける。
変数変換して、tをで表すと、
となるので、書き直すと、
$$\int^{N}_0 \frac{1}{\sqrt{2}} S(n_1) dn_1 = \frac{N^2}{2}$$
両辺をで微分すると、
$$S(N) = \sqrt{2} N$$
となって、ちゃんと答えが出てくる。
これを高次元に一般化すると、内部の超体積を求めて、原点から超平面までの距離で微分すれば、超平面の超面積を求める事ができる。
さて、L次元の場合で１辺がNの場合の体積は、
$$V_L(N) = \frac{N^L}{L!}$$
となることが分かっている。
実際、1辺をNとした場合に、2次元と3次元の場合の体積は、
$$V_2(N) = \frac{N^2}{2}$$
$$V_3(N) = \frac{N^3}{6}$$
となる。
の場合と同じように、高次元において、 超平面の中心の座標は (N/L, N/L, N/L, ... N/L) で表される。
なので、TをLとNで表すと、
$$T = \frac{N}{\sqrt{L}}$$
とかける。
そこで、これで変数変換を行うと
と書ける。
を代入してから、両辺Tで微分して左右を入れ替えると、
$$Z(N,L) = \frac{\sqrt{L}N^{L - 1}}{(L - 1)!}$$
となることがわかる。
（あ〜ようやく求まった。）
これで、元の問題に戻れる。
さて、本当に求めたかったのは、なんだったかというと、
「ボールがN個でマス目がLマスの場合に、ボールのやりとりを繰り返すと、
ある一つのマスのボールの確率分布はどうなるのか」ということだった。
NとLが大きい場合にシミュレーションして出てきた仮説は、平均的な玉の数をとして、
$$\rho(n_1) =\frac{1}{n} e^{-\frac{n_1}{n}}$$
という形になりそうだ、というものだった。
そして、厳密な計算の結果、この確率分布は
$$\rho(n_1) = \frac{Z(N - n_1,L - 1)}{Z(N,L)}$$
という式で表せるという所までは示せていた。
ここに、上で導いた、
$$Z(N,L) = \frac{\sqrt{L}N^{L - 1}}{(L - 1)!}$$
を代入すると
（数式を代入して、変形して、式にしたもの。3段かな）
と言うかたちになる。
ここで、マス目一つ辺りの平均的な玉の数をnとすると、なので、Lを消去すると、
と表される。
さて、ここでNを無限大に飛ばして、
$$\lim_{N \to \infty} (1+\frac{x}{N})^N = e^x$$
を用いると、
$$\rho(n_1) =\frac{1}{n} e^{-\frac{n_1}{n}}$$
と言うかたちになる。
おーー！！まじで導けた！！！
すごい！！！
ボルツマン分布との関係性
さて、ここまででわかった方もいるかもしれないが、実はこれは統計力学に出てくるボルツマン分布と全く同じものになっている。
マス目を分子、玉をエネルギーに見立てれば、この計算は全エネルギー一定の場合のボルツマン分布そのものなのだ。
矢野先生は本の中で、この分布は人間の活動の色々なところに現れる不思議な分布で、U分布と呼んでいた。
ただ、ここまで計算すれば、なんでそんな現象が現れるのかの理解はある程度出来る。
つまり、お金、土地、時間などどんなものでも、やりとりがとても多い場合には、ありうる場合が網羅的に生じる。
その結果として、全ての場合の数を考えた場合に理論的に出てくるU分布（ボルツマン分布）が現れるということだ。
まとめと直感的な説明
全体をまとめると、全体の量が決まったの何か（お金でもエネルギーでも）を、 たくさんのプレイヤー（人だったり分子だったり）がランダムにやりとりすると、たくさん持っているプレイヤーと何も持っていないプレイヤーが必然的に現れる。
そして、たくさん持っているプレイヤーの数は平均的な量（平均資産や平均エネルギー 温度）を係数にして、エクスポネンシャルに減っていくということだった。
なんでこんなことになるかって言うと、答えは簡単で、状況のパターンを網羅的に調べると、みんなが大体同じくらい持っている場合よりも、不平等にもっている場合のパターンの方が圧倒的に多いからだ。
単にそれだけなのだ。
eが現れる理由の直感的な説明はちょっとわからない。なんか自然界の法則がうまいこと出来ているからとしか言いようが無い気がする。
計算では導けるけど、直感的にはわからないものはたくさんある。これもその一つだと思う。
もしくは、わかりやすい説明出来る方がいたら是非教えてほしい。
おしまい
やりとりを繰り返すだけで貧富の差が生まれる？ ~「データの見えざる手」の分布~ part2
（part1からの続き）
ところで、なぜ自然対数の底:eが出てくるのだろう？
この分布は玉の数が0個になる場合はどれくらいの確率で起こるのか？ というようなことを表しているにすぎない。言ってしまえば、場合の数を計算しているだけのはずだ。
場合の数の計算で、eが現れるのは不思議だ。しかも玉の数Nやマスの数Lの値を変えても同じ様に現れる。
おそらくこれはNかLが大きい場合の近似式なのだろう。
ということは、場合の数を計算して、NかLを無限大に飛ばす操作をしたら、この式が導出できるはずだ。
やってみよう。
簡単な例での計算
詳しい計算に入る前に、まずは簡単な例で計算してみよう。
こうすると頭が整理されることが多い。
マス目が3マス、玉が3個という例を考える。
この場合、配置ごとの場合の数を考えると
3つのマス目に1個づつ玉が集まっている場合は1通り
1つのマス目に2個、もう1つのマス目に1個、残りのマス目に0個の場合は6通り
1つのマス目に3個集まって他2つが0個の場合は3通り
となっている。
なるほど、この例だと全部が平等に入っている時よりも、偏っている時の方が場合の数が多いということが分かった。
さらに、全通りを考えた時に、マスにいくつ入っている場合が最も多いのか数えると、
| 玉の数 | 場合の数 |
|---|---|
| 0 | |
| 1 | |
| 2 | |
| 3 | |
となっている。
なるほど。考えてみれば当たり前だが、0個のマス目の場合が一番多くて、玉の数が増えるに従って場合の数は減っていく。
なんとなくだけど、玉のやりとりを繰り返すと0個のマスが多くなる理由が見えてきたぞ。
ちなみに、この例を、前回と同じシミュレーションで回してみる。
マス目を3個、玉を3個用意して、ランダムに2つのマスを選んで、片方から片方に移すという操作を繰り返して、9操作毎にその時の分布を作る。
これを100000回繰り返して、その平均的な分布を求めてみると、
| 玉の数 | 確率（場合の数による） | 確率（シミュレーション） |
|---|---|---|
| 0 | 0.4 | 0.399937 |
| 1 | 0.3 | 0.300061 |
| 2 | 0.2 | 0.200069 |
| 3 | 0.1 | 0.099934 |
うむ。ちゃんと同じになっている。
つまり、この玉のやりとりのシミュレーションというのは、すべてのパターンをなるべく網羅するように、ランダムに動き続けている状況になっているということだ。
そしてある時点でのスナップショットを取ってきて、各マス目に入っているボールの数をヒストグラムにしているのが、図1で、スナップショットを平均してきれいな分布にしたのが、図2ということだ。
ところで、3個の場合は、なのだが、
分布を正規化して確率分布にしてみても、 の形になっていない。
やはり、この分布はNかLが大きい場合の近似式なのだろう。
さて、これでなんとなく、概観が分かった。
玉とマス目が多い場合
ここからが本題。マス目の数Lも玉の数Nもとても多い場合を考えてみよう。
全てのマス目について考えるのは大変なので、一つのマス目だけについて、玉がいくつの時に確率がいくつになるのかを計算する。
をマス目のインデックス、
をマス目
の玉の個数とする。
のマス目の玉の個数
の確率分布
を求めてみる。
まず、全てのありうる場合の数を式で表してみよう。
これは
$$\sum_{n_1 =0}^N \sum_{n_2=0}^N ... \sum_{n_L=0}^N \delta_{N,\sum_{l=1}^Ln_l}$$
と書くことが出来る。
ここではクロネッカーのデルタというもので、
$$
\delta_{i,j} =
\left\{
\begin{array}{}
1 & ( i = j) \\
0 & ( i \neq j)
\end{array}
\right.
$$
i,jを整数として、単にiとjが同じ値だったら1、そうでない場合は0になるという記号だ。
つまり上の式は、全てのマス目がバラバラに0個からN個の玉が入っている状況を考えてみて、その内全玉がN個の時のみカウントするという考え方だ。
さて、これを踏まえた上で、
確率分布を考えると、
$$\rho(n_1) = \frac{\sum_{n_2 =0}^N \sum_{n_3=0}^N ... \sum_{n_L=0}^N \delta_{N-n_1,\sum_{l=2}^Ln_l}}{\sum_{n_1 =0}^N \sum_{n_2=0}^N ... \sum_{n_L=0}^N \delta_{N,\sum_{l=1}^Ln_l}}$$
と書ける。
分母は全場合の数だが、分子ではの値が固定されていることに注意してほしい。
つまり、以外の全てのマス目がバラバラに0個からN個の玉が入っている状況を考えてみて、その内
以外の全ての玉が
個の時のみ許されるという考え方だ。（繰り返しになりますが。）
さて、の分母も分子も式としては同じ形をしている。
ということは、
$$Z(N,L) = \sum_{n_1 =0}^N \sum_{n_2=0}^N ... \sum_{n_L=0}^N \delta_{N,\sum_{l=1}^Ln_l}$$
という、玉がN個でマス目がLの時の全場合の数：が計算出来れば、
$$\rho(n_1) = \frac{Z(N -n_1,L -1)}{Z(N,L)}$$
なので、も計算出来るということだ。
後少しだ。が計算出来れば、多分すぐに終わりだ。
でも、この計算は意外と長い。なので次に 続く。
（続く）
やりとりを繰り返すだけで貧富の差が生まれる？ ~「データの見えざる手」の分布~
以前読んだ本「データの見えざる手」の第一章に面白いシミュレーションがあった。
何かというと、900個のマス目を用意して、その中に72000個の大量の玉をばらまく。
最初は、ランダムにばらまくので、どのマス目にも均等に大体80個づつ入っている。
ところが、「ランダムに2つのマス目を選んで片方から片方に玉を移動させる」という単純な動作を大量に繰り返すと、なんと玉の分布が偏り、1個も玉が入ってないマス目がたくさん現れる一方で、たくさんの玉が入っているマス目も現れるという。
つまり、平等な状態からスタートして、平等にやりとりを繰り返しても、貧富の差のようなものが必然的に生じてしまうということだった。
「ランダムな動作を繰り返すのだから、結果も乱雑になって、だいたいどれも同じような玉の数になる。」という直感を裏切られる、面白いシミュレーションだと思う。
では、U分布のような状態において、玉が無いマスとたくさんあるマスはどれくらいの確率で存在するのであろうか？
実はこれは指数関数的な分布をしているということが分かっているそうだ。
そして、このような何もない状態とたくさん持ってる状態が指数分布的に存在している状況は、日常の様々な現象の中に普遍的(Universal)に現れるらしく、本の中ではこれをU分布と名付けていた。
実際、この分布はリニアにプロットすると、下の右図のように、玉が少ない状態がとても多い。そして、この縦軸を対数プロットすると、直線になる。つまり指数分布するということがわかっているそうだ。
さて、せっかく面白いものに出会ったので、自分でも少しシミュレーションをして、色々考察をしてみようと思う。
まずは、 この玉のやり取りのシミュレーションを自分でもやってみて、実際にU分布（指数分布）になるのかを確かめてみた。
本の中と同じように、900個のマス目を用意して72000個の玉を最初にばらまいた後、玉のやり取りをランダムに何度も行う。
回のやり取りを1stepとして、100万stepやり取りを繰り返す。
各step毎に玉がいくつ入ったマスが何マスあるのかを計測して、最後に平均をとるだけだ。
（コードは一番下に貼っておく。）
ある瞬間の玉の数の分布と、最後に100万回分平均をとった分布の二つの結果を見てみた。
おお！ある瞬間での分布は少しいびつでも、100万回分の平均を取るとめちゃくちゃきれいな図になった！
縦軸をlogでとってみよう。ちゃんと直線になるのだろうか？
おお！ちゃんときれいな直線になっている！
まあ、玉の個数が大きすぎる所は、数回しか現れないので多少ガタついているが、十分直線といえるだろう。
ということは、やはりこの分布はきれいな指数分布になっているということだ。
さて、ここまでは本の内容を再現しただけだ。
もう少し深く分布の性質を見てみたい。
指数分布になっているということは、というような式に従うということだ。
このaとbはどんな形になっているのだろうか。
これを調べるために、玉の個数や、マス目の数を変えてプロットしてみよう。
そして、うまく当てはまりそうな、指数の係数aとbを見つけてみる。
まず、玉の個数とマス目の数を変えたものをプロットしてみると、
どれもきれいな指数分布になっている。
改めて上の確率分布の式を とする。
係数a,bに関係するのはマス目の数と玉の数しかないだろう。そして、そんなに複雑な式ではないはずだ（特に根拠はない）。
そこで、色々なパターンの式をプロットして、きれいにシミュレーションの直線に乗っかる式を探してみよう。
まず、なので、aの値は
の時の値と同じである。
玉の数とマス目の数を変えた時にどうなっているかというと、
玉の数=16000, マス目の数=400: a=0.025
玉の数=72000, マス目の数=900 : a=0.125
玉の数=500000, マス目の数=5000 : a=0.01
となっている。どうやら、となっていそうだ。
次に、に固定して、bの値を探索的に探してみる。
すると、マス目の数をL、玉の個数をNとして、横軸のマス目に入っている玉の個数をxとすると、
と言う形になりそうである。
実際プロットしてみると下のようにピッタリ合う。
ここで、 として、nを1マス辺りの平均的な玉の数とすると、
と、めちゃくちゃキレイな形になる！
これは、なんかある！そしてなんだか導ける気がする！
（続く）
以下コード
_______________________________
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "mt19937ar.h"
int shuffle(int *x , int *y) {
if(*x==0){
return 1;
} else {
(*x)--;
(*y)++;
}
return 0;
}
int MC_step(int *box , int N , int M){
int i;
int x,y;
for(i = 0;i < N*M;i++){
x = genrand_int31() % N;
y = genrand_int31() % N;
shuffle(&box[x] , &box[y]);
}
return 0;
}
void add_sum_box(int *box , double *sum_box, int N){
int i,j;
for(i = 0;i < N;i++){
sum_box[box[i]]++ ;
}
}
int main(){
int N = 900 ; // # of box
int M = 80; // # of ball per box
int MM = M*50;
int num_MC_steps = 1000000;
int box[N];
double sum_box[MM];
int i,j,k,l,m;
for(i = 0;i < N;i++){
box[i] = M;
}
for(i = 0;i < MM;i++){
sum_box[i] = 0;
}
for(i = 0;i < num_MC_steps;i++){
MC_step(box,N,M);
add_sum_box(box,sum_box,N);
}
// 1ループ辺り、sum_boxにはNが加えられていく。トータルで1にする。
for(i = 0;i < MM;i++){
sum_box[i] /= (num_MC_steps*N);
}
for(i = 0;i < MM;i++){
printf("%d,%f \n",i,sum_box[i]);
}
return 0;
}
