意識が朦朧としてるから、拡散方程式から正規分布を出す (スラド日記 2018-04-25)

 ; 書類が……、というより、その前の調整でガリガリと……。

∂u/∂t = D ∂^2u/∂x^2, u(0, x)=0 when x!=0, ∫(-∞,∞)u(0,x)dx = 1 っていうのを解きたい。

要するにフーリエ変換するのだけど、 u≈Σ(k=-∞,∞)exp(i 2πk/T x)φ(t, 2πk/T) * 2π/T * 1/2π と書いてみる。

x で微分すると、 ∂u/∂x≈Σ(k=-∞,∞)(i 2πk/T)exp(i 2πk/T x)φ(t, 2πk/T) * 2π/T * 1/2π.

もう1度 x で微分すると、 ∂^2u/∂x^2≈Σ(k=-∞,∞)-(2πk/T)^2exp(i 2πk/T x)φ(t, 2πk/T) * 2π/T * 1/2π.

t で微分すると、 ∂u/∂t≈Σ(k=-∞,∞)(i 2πk/T)exp(i 2πk/T x)∂φ(t, 2πk/T)/∂t * 2π/T * 1/2π.

微分方程式から、∂φ(t, 2πk/T)/∂t ≈ - D(2πk/T)^2 * φ(t, 2πk/T).

2πk/T=ωとして、T→∞にすると∂φ(t, ω)/∂t = - Dω^2 * φ(t, ω).

ところで、のっけの式を t=0のところで、-T/2からT/2まで∫すると、 ∫(-T/2, T/2) u(0, x)dx ≈ φ(0, 0).

ということで、初期条件から、φ(0, 0)=1 なので、φ(t, ω)=exp(- Dω^2 t).

のっけの式に入れて、T→∞にすると、 u= 1/(2π) ∫(-∞,∞) exp(i ω x)exp(- Dω^2 t) dω.

あとは指数部分を平方完成して、地道に計算すると u = 1/(√(2π)√(2Dt)) * exp(- x^2 / (2 * 2Dt)).

ということで、分散 2Dt の正規分布.

とりあえず、正規分布の出し方でメジャーなやつは説明したはず。 次からストレスで何な状態になったら、何を書くのだろう……。

気分がはれないので、雑な中心極限定理 (スラド日記 2018-04-24)

 中心極限定理で正規分布を出すには、 キュムラント母関数を使うのが簡単。 だけど、2種類あるキュムラント母関数を出すには、 特性関数か積率母関数を出す必要がある。 という具合で、先が長い。

ということで、特性関数。

密度関数を波の合成と思えば、
f(x)≈Σ[k=-∞, ∞] 1/(2π) * exp(-ix* 2πk/T) * φ(2πk/T) * 2π/T.

逆向きにφ(2πk/T) ≈ ∫(-T/2, T/2) exp(ix * 2πk/T) f(x) dx.

そこでフーリエ変換とその逆変換で、φ(t) = ∫(-∞, ∞) exp(itx) f(x) dx, f(x) = 1/(2π) * ∫(-∞, ∞) exp(-itx) φ(t) dt.

前者が特性関数だけど、それを期待値の式とみると、φ(t) = E[exp(itX)].

マクローリン展開をすると、 φ(t)=E[1]/0!*(it)^0 + E[X]/1!*(it)^1+E[X^2]/2!*(it)^2+...

何度も微分して t=0 を代入すると、 φ^(k)(0) = i^k E[X^k] と原点まわりのモーメントが出せる。

独立な確率変数の和については、 φ_X+Y(t)= E[exp(it(X+Y))] = E[exp(itX) * exp(itY)] = E[exp(itX)] * E[exp(itY)] = φ_X(t) * φ_Y(t).

特性関数の対数が(第2)キュムラント母関数で、H_X(t)=logφ_X(t) = log E[exp(itX)].

独立な確率変数の和については、 H_X+Y(t)=H_X(t)+H_Y(t) と和なのがうれしい。

マクローリン展開して、 H_X(t) = 〈X〉_c/1! * (it) + 〈X^2〉_c / 2! * (it)^2 + 〈X^3〉_c / 3! * (it)^3 + ....

〈X^k〉_c のところが k 次のキュムラント。 何度も微分して 0 を代入すると出てくる。 しかも、〈X〉_c=μ, 〈X^2〉_c=σ^2.

確率変数に定数を足したとき、 H_X+b(t) = log E[exp(it(X+b))] = H_X(t) + itb
となり、〈X+b〉_c=〈X〉_c+b, 〈(X+b)^k〉_c=〈X^k〉_c (k>1).

確率変数を定数倍したとき H_aX(t)=E[exp(itaX)] = H_X(at) なので、 〈(aX)^k〉_c=a^k〈X^k〉_c.

なお、正規分布のばあい、 logφ(t)=log(exp(iμt-(σ^2)/2*t^2)=iμt-(σ^2)/2*t^2.

ようやく到着。

X₁, ..., X_nが独立で同一分布、 それぞれ平均μ、分散σ^2のとき、
X=((X₁+...+X_n)/n - μ)/(σ/√n)の分布が、 n→∞で、 平均 0、分散 1 に近付く。

Xを簡単にすると X = 1/(σ√n){X₁+...+X_n}-μ√n/σ.

平均こと第1次のキュムラントを計算すると、 1/(σ√n) * (n * μ) - μ√n/σ = 0.

分散こと第2次のキュラントを計算すると、 1/(σ^2 n) * (n * σ^2) = 1.

3次より上のキュムラントは、1/(σ√n)^k * (n * 〈X^k〉_c)→0.

ということで、キュムラント母関数が正規分布のそれに近付くので、 特性関数も正規分布のに近付き、 密度関数も近付く。

独立だったら、平均や分散が違っていても、 モーメントが高次でも定義できてれば、同じ感じ。

意識が茫洋とするので、ポアソン分布から雑に正規分布を出す(スラド日記 2018-05-02)

P(r)=λ^r * e^(-λ) / r! なので、 P(r-1) / P(r) = r / λ.

ポアソン分布は平均λ、分散λなので、t=(r-λ)/√λ と置き換える ことにする。 また、ΔP=P(r)-P(r-1) とする。このとき r = λ(1+t/√λ).

そうすると r を1だけ変えたとき、Δr = 1/√λ.

ということで、1- ΔP / P = 1+ tΔt. つまり、ΔP / P = - tΔt.

dP / P = -t dt という微分方程式を解いて、 ガウス積分で定数を決めてやれば、おしまい。

考える気力がないので、雑にカイ二乗分布から正規分布を出す (スラド日記 2018-05-22)

再生性があるので、自由度を大きくしていくと、 ゆっくり正規分布に近づく。

f(x; k) = x^(k/2 - 1) * e^(-x/2) / {2^(k/2) * Γ(k/2)} で対数をとって、

log f(x; k) = (k/2 - 1)log x - x/2 - log {2^(k/2) * Γ(k/2)}.

x で微分して、 f'(x) / f(x) = (k/2-1) * 1/x - 1/2.

t = (x - k) / √(2k) とおきかえる。 このとき、dx = √(2k) dt, x = k(1+√(2/k)).

f'(x)dx / f(x) = {(k/2 - 1) / k(1+√(2/k)) - 1/2} √(2k) dt.

k が大きいと、f'(x)dx / f(x) ≈ {1/2 * (1-√(2/k)) - 1/2}√(2k) dt = - tdt.

いつもの dp / p = - t dt が出てきたので、あとは同じ。

連続分布・分散有限でエントロピー最大が正規分布、雑なまとめ後編 (スラド日記 2018-04-22)

承前。 とりあえず最大の候補は出た。 変形の一部が思い付かなかったので、 式変形の元ネタ

まずは、正規分布のときのエントロピーを求めておく。

HN=-∫(1/√(2π)σ)*exp(-(x-μ)^2 / 2σ^2) log{(1/√(2π)σ)*exp(-(x-μ)^2 / 2σ^2)}dx
=log(√(2π)σ)∫(1/√(2π)σ) * exp(-(x-μ)^2 / 2σ^2)dx + 1/(2σ^2) ∫(x-μ)^2(1/√(2π)σ) * exp(-(x-μ)^2 / 2σ^2)dx
=log(√(2π)σ)+ 1/2.

差を計算する。

H - HN = -∫p log p dx - log(√(2π)σ)- 1/2
= -∫p log p dx - log(√(2π)σ)∫p dx - 1/(2σ^2)∫(x-μ)^2 p dx
= -∫p log p dx + ∫p log{1/√(2π)σ * exp(-(x-μ)^2 / 2σ^2)}dx
= ∫p log({1/√(2π)σ * exp(-(x-μ)^2 / 2σ^2)}/p) dx.

グラフの形から log x ≦ x - 1 なので、

H - HN ≦ ∫p ({1/√(2π)σ * exp(-(x-μ)^2 / 2σ^2)}/p - 1)dx
=∫{1/√(2π)σ * exp(-(x-μ)^2 / 2σ^2)}dx - ∫ p dx
= 1 - 1 = 0.

ということで、H ≦ H_N で、等号成立は p = 1/√(2π)σ * exp(-(x-μ)^2 / 2σ^2) のとき。

メジャーな出し方で、まだ書いてないのは、 拡散方程式と中心極限定理かな。

連続分布・分散有限でエントロピー最大が正規分布、雑なまとめ前編 (スラド日記 2018-04-21)

;; 編集しようとして、隣の削除押した……。

ちゃんとやると道具だけで大変なので、適当。

エントロピーは H=-∫p(x) log p(x) dx. 面倒なので、積分区間の-∞から∞は省略。 情報っぽく対数の底を2にするのも面倒なので、e。 あと (x) も略。

pをp+Δpにちょっと変えて、ΔH を考える。 Δp が小さいので、Δp^2 のあとは無視する。 log (1+Δp/p)≈Δp/p なので、

ΔH = -∫(p+Δp)log(p+Δp) dx + ∫p log p dx = - ∫(p log(1+Δp/p) +Δp(log p + log(1+Δp/p)))dx ≈ - ∫dx (1+log p)Δp.

制約は、確率の和が1の M0=∫p dx = 1, 平均がμの M1=∫x p dx = μ, 分散がσ^2 の M2=∫(x-μ)^2 p dx = σ^2. 最後のを原点まわりに書き直して、M'2=∫x^2 p dx = σ^2+μ^2.

これの差分もほとんど0で、 ΔM0=∫dx Δp ≈ 0, ΔM1=∫dx x Δp ≈ 0, ΔM'2=∫dx x^2 Δp ≈ 0.

ということは、 1 + log p = λ0+λ1x+λ2x^2 とすると、ΔH≈0 にできる (これが Lagrange の未定乗数法)。

p = exp(λ0-1+λ1x+λ2x^2) で指数部分を平方完成して、 制約を満たすようにすると、

p = 1/(√(2π)σ) * exp{-(x-μ)^2/(2σ^2)}.

正規分布が出てきたけど、これはまだ候補なので、後編に続く。

Adrainの正規分布の出し方、無茶だぜい (スラド日記 2018-04-15)

Robert Adrain https://en.wikipedia.org/wiki/Robert_Adrain 最小二乗法の人。 Gaussの導出より早いらしい。 元論文を誰かが打ち直したもの https://pdfs.semanticscholar.org/1900/ca31f56572c7fc588462c3fb5e10789e0ea5.pdf

雑なまとめ。棒の長さを測るつもり。 真の長さA、測定値a で誤差はα=a-A。 この棒のまっすぐ別の棒を繋いで、 2本目の棒でも真の長さB、測定値b で誤差はβ=b-B。 前提は:

• 誤差の密度関数は「相似」で 1/a*f(α/a)と1/b*f(β/b)。
• 誤差の合計がεのとき、ε=α+βで、このとき いちばんもっともらしいのは a:b=α:βのときのはず。
• こういう関数で一番簡単なのが欲しい。

P=1/a*f(α/a)と1/b*f(β/b)を最大にするけど、対数をとって、楽をする。 log P = log f(α/a) + log f(β/b) - log ab。 これを何かで微分して、 (log P)' = f'(α/a)/f(α/a) * α'/a + f'(β/b)/f(β/b) * β'/b.

最大になるときなので、ε=α+βから 0=α'+β' で (log P)' = 0。

そこから整理すると、f'(α/a)/f(α/a) /a = f'(β/b)/f(β/b) /b。

つまり、a:b = f'(α/a)/f(α/a) : f'(β/b)/f(β/b)。

a:b=α:βでもあるから、f'(α/a)/f(α/a) : f'(β/b)/f(β/b)=α:β。

そこで、何か定数 k があって、f'(α/a)/f(α/a)=kα, f'(β/b)/f(β/b)=kβ.

もう少し変形して、f'(α/a)/f(α/a)/a=kα/a, f'(β/b)/f(β/b)/b=kβ/b.

これらをみたす簡単な関係は (log f(t))'=kt で、これが微分方程式。