ガウス関数への回帰 update 20181115
単峰性関数のガウス関数による回帰分析
「ガウス関数に似た関数」では、関数を5つ取り上げ、ガウス関数（y=exp(-x^2)）に対してそれぞれの規格化相互相関係数（規格化相互相関関数の最大値を規格化相互相関係数と仮称した）が最大になる形を（係数を）定めました。
今回は逆に、得られた4つ（矩形は除く）の関数に回帰分析を使ってガウス関数を当てはめてみようというものです。
回帰分析に用いるガウス関数は、y=c*exp(-b*(x-a)^2)です。
ただし、この関数を回帰分析に直接使うのではなく、両辺の自然対数をとって得られる2次曲線を回帰分析に使います。
y=c*exp(-b*(x-a)^2)の両辺の自然対数をとります。
ln(ｙ)=ln(c)-b*(x-a)^2
=ln(c)-b*(x^2-2*a*x+a^2)
=(ln(c)-b*a^2)+2*a*b*x-b*x^2
ln(ｙ)＝Y、(ln(c)-b*a^2)＝α、2*a*b＝β、-b＝δ、とおくと
Y=α+β*x+δ*x^2 となります。
この2次曲線の2次多項式回帰を行うことができます。
α、β、δを使ってもともとの係数a、b、cが求まることを確認しておきましょう。
b＝-δ
a＝β/(2*b)=-β/(2*δ)
ln(c)＝α+b*a^2＝α+（-δ)*(-β/(2*δ))^2＝α-(β^2)/(4*δ)
∴c＝exp(α-(β^2)/(4*δ))
今、y=1/(1+2.2491689188*(x^2))、-3≦x≦3（0.1ステップ）の範囲で得られる61個のデータ(x1,y1)、(x1,y1)、(x2,y2)、・・・(x60,y60)、(x61,y61)のうち、y座標の値について自然対数をとり、（x1、Y1）、（x2,Y2）、・・・(x61,Y61)として、2次多項式回帰分析を行います。
【平均と分散から係数α、β、δを求める】
β＝(Ｓ(x･Y)*Ｓ(x^2･x^2)-Ｓ(x^2･Y)*Ｓ(x･x^2))/(Ｓ(x･x)*Ｓ(x^2･x^2)-(Ｓ(x･x^2))^2)
δ＝(Ｓ((x^2･Y)*Ｓ(x･x)-Ｓ(x･Y)*Ｓ(x･x^2))/(Ｓ(x･x)*Ｓ(x^2･x^2)-(Ｓ(x･x^2))^2)
α＝[Y]-β[x]-δ[x^2]
ここで、
【平均：[x]、[Y]、[x^2]】
[x]＝(∑xi)/n ,(i：1→ｎ)
[Y]＝(∑Yi)/n ,(i：1→ｎ)
[x^2]＝（∑(x^2)）/n ,(i：1→ｎ)
【分散：Ｓ(x･x)、Ｓ(x･Y)、Ｓ(x･x^2)、Ｓ(x^2･x^2)、Ｓ(x^2･Y)】
Ｓ(x･x)=(∑（xi-[x]）^2)/n =(∑xi^2)/n－[x]＾2
Ｓ(x･Y)=((∑(xi-[x])*(Yi-[Y]))/n =(∑(xi*Yi))/n-[x]*[Y]
Ｓ(x･x^2)=(∑(xi-[x])*((xi^2-[x^2]))/n =(∑xi^3)/n-[x]*[x^2]
Ｓ(x^2･x^2)=(∑(xi^2-[x^2])^2)/n =(∑xi^4)/n-[x^2]*[x^2]
Ｓ((x^2･Y)=(∑(xi^2-[x^2])*(Yi－[Y]))/n =(∑(xi^2*Yi))/n-[x^2]*[Y]
です。
今、-3≦x≦3（0.1ステップ）の範囲で、
y=1/(1+2.2491689188*(x^2))の、回帰曲線y=c*exp(-b*(x-a)^2)を計算します。
n＝61
[x]＝0.000000
[Y]＝-1.680151
[x^2]＝3.100000
(∑xi^2)/ｎ＝3.100000
(∑xi*Yi)/ｎ＝0.000000
(∑xi^3)/ｎ＝0.000000
(∑xi^4)/ｎ＝17.291800
(∑xi^2*Y)/ｎ＝-7.711895
S(x･x)＝(∑xi^2)/ｎ-[x]^2＝3.1000000000
S(x･Y)＝(∑xi*Yi)/ｎ-[x]*[Y]＝0.0000000000
S(x･x^2)＝(∑xi^3)/ｎ-[x]*[x^2]＝0.0000000000
S(x^2･x^2)＝(∑xi^4)/ｎ-[x^2]*[x^2]＝7.6818000000
S(x^2･Y)＝(∑xi^2*Y)/ｎ-[x^2]*[Y]＝-2.5034260881
β＝(S(x･Y)*S(x^2･x^2)-S(x^2･Y)*S(x･x^2))/(S(x･x)*S(x^2･x^2)-S(x･x^2)*S(x･x^2))＝0.0000000000
δ＝(S(x^2･Y)*S(x･x)-S(x･Y)*S(x･x^2))/(S(x･x)*S(x^2･x^2)-S(x･x^2)*S(x･x^2))＝-0.3258905580
α＝[Y]-β*[x]-δ*[x^2]＝-0.6698905997
b＝-δ＝0.3258905580
a＝-β/(2*δ)＝0.0000000000
c＝exp(α-(β^2)/(4*δ))＝0.5117645619
したがって、回帰曲線は
y=0.51173*exp(-1*(0.3259022*(x-0))^2)
また、【相関係数ｒ】は、
ｒ=sqrt(1-((∑(Yi-(α+β*xi+γ*x^2))^2)/(∑(Yi-[Y])^2))))
です。
2次多項式回帰分析で、どのような曲線をあてはめたのか確認します。
下のグラフにおいて、
曲線(赤)：y=1/(1+2.2491689188*(x^2))
曲線(青)：yの自然対数 Y=ln(y)
曲線(灰)：Y=ln(y)に対する2次多項式回帰曲線 Y＝α+β*x+δ*x^2
→Y＝-0.6698905997-0.3258905580*x^2
当初の目論見とだいぶ違うな～というのが正直な感想です。
ちなみに・・・
ちなみに、上に凸に2次曲線
y=1-1*0.542133027*(1*(x-0))^2、（-1.3≦x≦1.3、0.1ステップ）
については、
a=0.000000000
b=1.215002811
c=1.182414574
であり、
y=1.1824146*exp(-1*1.2150028*(x-0)^2)
が回帰曲線となります。
さらに、三角関数の余弦曲線
y=1.1824146*exp(-1*1.2150028*(x-0)^2)、（-1.7≦x≦1.7、0.1ステップ）
については、
a=0.000000000
b=1.494667803
c=1.346702851
であり、
y=1.3467029*exp(-1*(1.4946678*(x-0))^2)
が回帰曲線になります。
また、三角波
y=1*(-1*0.62538275*1*(abs(x)-0)+1)、（-1.5≦x≦1.5、0.1ステップ）
については
a=0.000000000
b=1.052484707
c=0.95818133
であり、
y=1*(-1*0.62538275*1*(abs(x)-0)+1)
が回帰曲線になります。
----- ここまで 20181115 -----