- ベストアンサー
データ列の曲線によるフィッティング
データ列(xi,yi)(i=1,...,n)を関数 y=α(x+γ)^2 exp(-β(x+γ)) (α>0,β>0,γ>0) でフィッティングしたいです。対数をとって普通に最小２乗法で解こうとしたら得られる連立方程式が線形でなくて解けませんでした。 どうしたらいいのでしょうか？
- continuous
- お礼率15% (3/20)
- 数学・算数
- 回答数3
- ありがとう数0
- みんなの回答 （3）
- 専門家の回答
質問者が選んだベストアンサー
完全な意味での最小二乗法ではないけれど、無理矢理線形の問題に帰着して簡単に計算するstomachman我流のやりかたがあります。（超平面法と勝手に呼んでいます。）この方法は、「データがモデルに良く合う」という時に限り「最小二乗解とほぼ一致する」という性質があります。 モデルをデータと区別して、 f(t) = α((t+γ)^2) exp(-β(t+γ)) と書きましょう。 両辺をtで微分すると f'(t) = 2α(t+γ) exp(-β(t+γ))-βα((t+γ)^2) exp(-β(t+γ)) (t+γ)f'(t) = 2α(t+γ)^2) exp(-β(t+γ))-(t+γ)βα((t+γ)^2) exp(-β(t+γ)) (t+γ)f'(t) = 2f(t)-(t+γ)βf(t) (t+γ)f'(t) = (2-βγ)f(t)-βtf(t) という微分方程式が得られます。 この両辺をtで積分すると、(t=x1～xiまで） ∫(t+γ)f'(t)dt = (2-βγ)∫f(t)dt-β∫tf(t)dt 左辺を部分積分して (xi+γ)f(xi)-(x1+γ)f(x1)-∫f(t)dt = (2-βγ)∫f(t)dt-β∫tf(t)dt 整理すると xif(xi)-x1f(x1) = (3-βγ)∫f(t)dt-β∫tf(t)dt-γ(f(xi)-f(x1)) です。 ここで、f(xi)の代わりにyiを使っても、「データとモデルが良く合う」のなら大差は出ない筈。だから ∫f(t)dtの代わりにyをx1～xiまで数値積分したものJi、 ∫tf(t)dtの代わりにxyをx1～xiまで数値積分したものKi、を用いることにして、残差をεiと書くことにし、 xiyi-x1y1 = (3-βγ)Ji-βKi-γ(yi-y1)+εi (i=1,2,....,n) 改めて Yi = xiyi-x1y1 Li = yi - y1 A = (3-βγ) B = -β C = -γ とおくと、 Yi = A Ji + B Ki + C Li + εi (i=1,2,....,n) という線形問題になったわけで、これは簡単に解けます。 ところで、A,B,Cはβとγだけで表される一方、αが出てきませんね。 つまり ・A,B,Cからβ,γを求めると１通りに決まりません。これは気にしないで、たとえばAを無視してβ,γを計算します。ついでに、こうして得たβ,γを使って(3-βγ)を計算してみます。もしこれがAとまるで合わないようなら、「データがモデルに良く合う」という条件が満たされていない。yiに含まれるノイズが大きすぎるんです。 ・αを求めるには、 f(t) = α((t+γ)^2) exp(-β(t+γ)) を再び使います。γとβは既に分かっているから、簡単な線形最小二乗法の問題ですね。 より高精度の最小二乗解を求めるには、超平面法で得た解を初期値にして、非線形最小二乗法（ガウス・ニュートン法など）を使って改良すると良いです。（教科書としてはいつも、中川・小柳「最小二乗法による実験データ解析」東京大学出版会をお薦めしています。）
その他の回答 (2)
最小二乗法と非線型２情報の話しは既に出ていますので．ＯＲの話しにしましょう。 最適化問題で非線型問題がどうしても融けない場合に使います。 まず．でたらめにアルファ・べーた・ガンマーの値を決めます。 簡単にするのでしたらば対数で１ｅ－８から１え＋８の間が均等になるような点（１．１０．１００とか標準数（ＪＩＳＺ番号忘れ）とかあります）を選び．残差を求めます。 そのなかで．もっとも残差が少ない．数値の組み合わせを選びます。 次に．どれか一つをちょっと変えて計算してみます。「ちょっと」の程度ですが．切りの良い（２進数で考えてください）１／１２８ぐらいではじめます。もし．最初の計算よりも残差が少なくなったらばその解をえらびます。もし．残差がおおきくなってしまったらは．「ちょっと」の程度を更に１／２にします。計算していて．計算機Ｅの問題が出てきたら計算を止めます。 この考え方を基にして．３つの変数の解の付近（目安としては前回の計算した範囲程度）を１２８分割して演算をしてその中のもっとも残差の少ない点を解とします。 最低でも倍精度演算が必須です。誤差関数を求めるために８倍精度とか無限倍精度（文字演算で行うために計算がやたら遅いのでやりたくない）とかを使う方法もあります。 ＯＲ問題の亜種です。私は場当たり法と呼んでいます。メッシュで計算する法ですから．解の巣問題や鞍部問題（１の方の参考文献を読んでください）に対して（計算を間違えなければ）強力な計算方法です。しかし．計算時間がやたらかかるので．評判はよろしくないです。８－１０変数で．１００点．８０４８６で８時間なんてざらです（真の最適解とは限らないので．初期値をちょっとずらして再計算するとすぐに１－２ヶ月使ってしまいます）から。
- siegmund
- ベストアンサー率64% (701/1090)
これは非線形最小二乗法の問題になりますが， 過去いくつかのスレッドがあり，処方もいろいろ回答されています． たとえば， http://oshiete1.goo.ne.jp/kotaeru.php3?q=179919 http://oshiete1.goo.ne.jp/kotaeru.php3?q=190632 http://oshiete1.goo.ne.jp/kotaeru.php3?q=98446 http://oshiete1.goo.ne.jp/kotaeru.php3?q=97271 http://oshiete1.goo.ne.jp/kotaeru.php3?q=33318 など． 他にも質問検索で「最小二乗法」「最小自乗法」とやるとかなりヒットします．
関連するQ&A
- ガウスフィッティングについて
時間 t に対するデータ列 yi (i=1, 2, …) に 最小二乗法を用いて ガウス関数 y(t) = a + b exp( -(t-c)^2 / 2σ^2 ) を 当てはめたいのです。 これまでは、エクセルなどのソフトを用いて 未知数 a, b, c, σを求めていたのですが、 これをソフトなどを用いずに理論式で解くにはどうすれば良いでしょうか？ ひとつのやり方として、 ガウス関数の右辺の a を左辺に移行して、両辺の対数を取り 対数データ列 ln(yi - a) との誤差を求めて最小二乗法… という手を考えたのですが、 これだと、誤差の総和を求める際、Σの中にパラメータの a が ln(yi - a) という形で残ってしまい、 誤差の偏微分で得られるパラメータに関する方程式が、 非線形な連立方程式となって解くことができません。 (説明が分かり難くて申し訳ありません) その他にもいろいろ考えたのですが、 どうにもうまくいきません。 とにかく、パラメータの a が邪魔で、 これを b, c, σを用いずに先に求める、 という方向で考えているのですが、 何か良い方法がありましたら、よろしくお願いします。
- ベストアンサー
- 数学・算数
- 実験データのフィッティングについて
i 個の測定データ(x[i],y[i]) を，最小二乗法などを用いて下記の式にフィッティングさせ、AとBを求めたいと思うのですが、私の勉強不足で線形最小二乗法（グラフにプロットして、切片と傾きから求める方法）で解く方法が分かりません。 Y = Alog{x/(x-B)} (x>B) 考え方だけでも構いませんので，どうかご教授下さい。よろしくお願いいたします。 また、最小二乗法に関する大学学部生程度のレベルの教科書的な本がありましたら、教えて下さい。よろしくお願いします。
- ベストアンサー
- その他（学問・教育）
- 最小二乗平面
ある複数の空間座標(x1,y1,z1)～(xn,yn,zn)(nは3以上)から、平面近似式である最小二乗平面の方程式を求める関数を作ろうと考えています。 平面方程式はz=ax+by+c(a,b,cが定数)であらわされ、引数を座標と座標個数n、戻り値をa,b,cにします。 ここ(http://oshiete1.goo.ne.jp/qa2802443.html)を参考に 最小二乗平面の連立方程式を解くコードを書いたのですが、 どうも答えが合いません。どなたかご教授願えないでしょうか？ 開発環境はC++Builder2007です。 ↓の数式をコードにしましたが、コードが間違っているのか、 数式自体がダメなのかさっぱりわかりません。 //与えられるn個の3次元座標(xi,yi,zi)から平面方程式を求める //平面方程式:z = ax + by + c //最小二乗平面を求める連立方程式は下記のようになる。 // aΣxi^2 + bΣxiyi + cΣxi = Σxizi // aΣxiyi + bΣyi^2 + cΣyi = Σyizi // aΣxi + bΣyi + cn = Σzi //これを行列で解く // |Σxi^2 Σxiyi Σxi | |a| = |Σxizi| // |Σxiyi Σyi^2 Σyi | |b| = |Σyizi| // |Σxi Σyi n | |c| = |Σzi | //ここで // |Σxi^2 Σxiyi Σxi | // A = |Σxiyi Σyi^2 Σyi | // |Σxi Σyi n | // // |Σxizi| // B = |Σyizi| // |Σzi | // // |a| // C = |b| // |c| // //とすると // // C = B・A^-1 // //で求めることができる
- 締切済み
- C・C++・C#
- 指数関数のカーブフィッティング
時間ｔに対する１ｃｈデータ列ｙがありまして、それを ｙ＝ａ ｅｘｐ（ｂ ｔ） ＋ ｃ に対して客観的に、できれば自動的にフィッティングして、ａ，ｂ，ｃを求めたいです。 これがただの1次関数の最小二乗法ならわかりますし、ｃが既知なら１次関数の応用で、というところまでもわかります。恥ずかしながら渡井には非線形最小二乗法を一般論で理解して解けるような気がしません。 Excelを使った最小二乗法手順説明サイト http://szksrv.isc.chubu.ac.jp/lms/lms2.html のような方法か、 C/C++のプログラム http://www.sist.ac.jp/~suganuma/kougi/other_lecture/SE/predict/predict.htm#2 のようなアルゴリズムの説明をいただけると大変ありがたいです。 よろしくお願いします。
- ベストアンサー
- 物理学
- データのフィッティングについて
データを誤差関数erf(x)=(2/√π)×∫_0^x exp(-u^2)duでフィッティングしてパラメータxを求める方法が分かりません。webで調べたところ、scilabでフィッティングできるようなのでインストールしたのですが、詳しい解説書がなくわかりません。scilabに限らず、このようなフィッティングの手順の分かる方、詳しく教えていただければと思います。
- 締切済み
- その他([技術者向] コンピューター)
- Ｒを使った非線形最小二乗について
はじめまして。 Ｒを使った非線形最小二乗のフィッティングについて、ご教授お願いします。 （カテゴリーが違うようでしたらご指摘お願いします） 二つのデータの組Ｙ＝（y1,y2,y3,・・・）、Ｚ＝(z1,z2,z3,・・・)について、これらをf(x) = α / { 1 + β exp(-γx) }の関数でフィッティングしたいと考えています。 ここで、f(x)の未知パラメータα、βはＹとＺで同じ値を持ち、γはＹとＺに固有のパラメータとして、フィッティングしようとしています。 以上のことをＲのnls関数を用いて行いたいが、どのようにすればよいか。 以下の流れでフィッティングしてみましたが、上手くいきませんでした。 (1) Ｙをf(x)でフィッティングして、α、β、γを求める。 (2) (1)のα、βを用いてＺをf(x)でフィッティングする。 そもそも、αとβを共通の値と考えること自体が間違いというご指摘があるかもしれませんが、 あくまでもαとβが共通の値と考えたとき、ＹとＺをf(x)で最も上手くフィッティング（Ｙの残差とＺの残差の和が最小になる）したいと考えています。 よろしくお願いいたします。
- 締切済み
- 数学・算数
- 最小二乗法について
最小二乗法では二乗和の誤差 Σ[i=1～n]{Yi-(α+βXi )}^2 （iは添え字です） を最小化するα,βを推定することを考えますが、 これは単純にα,βで偏微分してそれを0とおいて 連立方程式を解くだけでよいのですか？ といいますのも、2変数関数の極値を求める場合、 Hessianを計算して判別しますよね？ ただ一階偏導関数が0になるからといって、 そこで極値をとるとは限らない気がしたので… それとも最小二乗法の場合は必ずとるようになっているのでしょうか？ 手元の本には、 「この二乗和は非負値なので、αとβで偏微分したものを0とするα,βが上式を最小にする値である」 とあるのですが、一般に非負値だとこの ようなことが言えるのでしょうか？
- ベストアンサー
- 数学・算数
- 最小二乗法のn次曲線について
最小二乗法のn次曲線について Pn(x)=anx^n+an-1x^(n-1)+……+a1x+a0 の時、最小二乗誤差E2は、 E2=Σ(i=1,m)(yi-Pn(xi))^2 =Σ(i=1,m)yi^2-2Σ(i=1,m)Pn(xi)yi+Σ(i=1,m)(Pn(xi))^2 =Σ(i=1,m)yi^2-2Σ(i=1,m)(Σ(j=0,n)ajxi^j)yi+Σ(i=1,m)(Σ(j=0,n)ajxi^j)^2 ここまではわかるんですが、次の式になる理由が分りません。 E2=Σ(i=1,m)yi^2-2Σ(j=0,n)aj(Σ(i=1,m)yixi^j)+Σ(j=0,n)ajΣ(k=0,n)ajak(Σ(i=1,m)xi^(j+k)) 一番後ろの項Σ(j=0,n)ajΣ(k=0,n)ajak(Σ(i=1,m)xi^(j+k))はどうやったらでてくるんでしょうか？ なんでいきなりkがでてくるんでしょうか？ jとiの組を２乗してるんだからkというのがでてくるのは変だとおもうんですが、どういう考え方なんでしょうか？
- ベストアンサー
- 数学・算数
- 最小二乗法の分散の求め方
http://oshiete1.goo.ne.jp/qa3077638.htmlに関連しての質問です。 例えば、ｙ＝Ｘβ+εに関して最小二乗解を求めると b = [ nΣ(xi yi) - (Σxi) (Σyi)]/[ nΣ(xi^2) - (Σxi)^2 ] となります。ここから分散を求めるためにはどうすればよいのでしょうか？教科書を引っ張ってみると求め方の行列の式しか書いていなくいまいちピンときません（確かに計算すれば正しい結果を得られるようですが）。具体的にこの式だけを使って分散を求めるということはできないのですか？
- 締切済み
- 数学・算数
