- ベストアンサー
べき乗則を最小二乗法で求める
卒業研究で計測した分布をべき乗則 y = ax^k に類似させるべく最小二乗法を計算しようと思ってるのですが 自力で計算した結果 k = {loga×Σ(logy) - Σ(logy × logx)} / Σ(logx)^2 a = Σ(yx^k)/Σ{x^(2k)} で行き詰まって係数の計算までに至りません 最小二乗法で調べてみても簡単な一次式の例しか出ず参考になりません どなたか適切な計算方法について説明お願いします
- みんなの回答 (3)
- 専門家の回答
質問者が選んだベストアンサー
y = ax^kの両辺の対数をとって、lny=lna+klnx。 Y=lny,A=lna,X=lnxと置くとY=A+kXとなって、一次式の形になります。 各データから、Y(i)=lny(i),X(i)=lnx(i)を使って求めれば普通の最小二乗法でできます。
その他の回答 (2)
- arrysthmia
- ベストアンサー率38% (442/1154)
No.2 と比べれば判るように、貴方の式は、 k = { (log a)×Σ(log y) - Σ(log y × log x) } / Σ(log x)^2 が、 (∂/∂k) Σ{ (log y) - log(a x^k) }^2 = 0 を変形したもので、 log y の二乗誤差を最小にする a, k を求める式の一部、 a = Σ(y x^k) / Σ{ x^(2k) } が、 (∂/∂a) Σ( y - a x^k )^2 = 0 を変形したもので、 y の二乗誤差を最小にする a, k を求める式の一部、 になっています。どちらか一方に統一すれば、それなりの解が得られます。 どちらに統一すべきかは、応用上の都合で決まるので、 数学というより、その卒業研究の内容に従って考えるべきでしょう。 y の二乗誤差のほうを最小化することに決めた場合、 a, k の連立方程式は、非線形方程式になりますから、 反復法などを使って、近似解を求めるしかありません。 そちらが希望なら、ニュートン法などについて調べてみるとよいでしょう。 参考: http://akita-nct.jp/yamamoto/lecture/2007/5E_comp_app/nonlinear_equation/nonlinear_eq_html/node8.html#SECTION00082000000000000000
- sanori
- ベストアンサー率48% (5664/11798)
こんばんは。 No.1様のご回答のように、 y = ax^k から lny = klnx + lna (自然対数) あるいは、 logy = klogx + loga (常用対数) としてから、 logy と logx との関係を最小二乗法として考える方法が一つ。 この方法は、両対数グラフを作ったときの、直線と各点とを近づけたい場合に用います。 もう一つは、対数目盛りでないグラフ(曲線になります)において、 各点と曲線とを近づけるという目的で行う方法です。 誤差をεと置けば、 y + ε = ax^k ε = ax^k - y ε^2 = (ax^k - y)^2 = a^2・x^2k - 2ayx^k + y^2 データに番号をつければ、 εn^2 = a^2・xn^2k - 2aynxn^k + yn^2 よって、2乗誤差の合計Sは、 S = Σεn^2 = Σa^2・xn^(2k) - 2Σaynxn^k + Σyn^2 = a^2・Σxn^(2k) - 2a・Σynxn^k + Σyn^2 ここで、 xn、yn はデータなので、既知の定数です。 逆に、aとkは不明ですから、変数として考えることができます。 Sを最小にするには、Sがaやkに関して極値を取るようにすればよいのですから、 それは、Sをaやkで偏微分したときにゼロになるようにすればよいということです。 S = a^2・Σxn^(2k) - 2a・Σ(yn・xn^k) + Σyn^2 でしたから、 ∂S/∂a = 2a・Σxn^(2k) - 2・Σ(yn・xn^k) = 2(a・Σxn^(2k) - Σyn・xn^k) ∂S/∂k = a^2・Σ(2lnxn・xn^(2k)) - 2a・Σ(yn・lnxn・xn^k) = 2a(a・Σ(lnxn・xn^(2k)) - Σ(yn・lnxn・xn^k)) これらがゼロになればよいので、 a・Σxn^(2k) - Σ(yn・xn^k) = 0 a・Σ(lnxn・xn^(2k)) - Σ(yn・lnxn・xn^k) = 0 という連立方程式で、aとkが決まります。 この先は、計算していません。 以上、ご参考になりましたら幸いです。
お礼
なるほど、そう考えると単純ですね 私は対数は取ったのですが S = Σ(logy-loga-klogx)^2 のように変換せずそのまま最小二乗法をやってましたorz 参考にします、ありがとうございます。