（その４）では

　　Φ(x)≒{1-exp(-2x^2/π)･f(x)}^(1/2)

　　Φ(x)≒{1-exp(-2x^2/π･g(x))}^(1/2)

が最良近似式となる関数f(x),g(x)を求めるために，変分法を用いて

　　f(x)=exp(2x^2/π)[1-{Φ(x)}^2]

　　　　=1+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

　　-2x^2/π･g(x)=log[1-{Φ(x)}^2]

　　　　　　　　 =-2/π･x^2+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

と計算されることを示した．

　

　しかし，変分法を用いなくても

　　exp(-2x^2/π)f(x)=[1-{Φ(x)}^2]

　　-2x^2/π･g(x)=log[1-{Φ(x)}^2]

は自明のことであるので，これだけでは変分法の有難味は少ない．ところが，この計算に際して，両者に式

　　+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

が共通して現れることがわかった．

　

　このことは，この種のデータを見慣れているひとにとって，数学的に異常あるいは不自然と感じられるのではなかろうか？　すなわち，多項式を指数関数中に引き込んでも引き込まなくても，同じ式が出現するというのである．今回のコラムでは，関数f(x)とg(x)の関係について取り上げ，これに引き続く項も等しいのかどうか，その真偽について確かめてみたい．

　

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

　

【１】宿題報告

　

　まず最初に，前回積み残しの宿題を片づけておきたい．統計数値表（日本規格協会）には，

　　f(x)=Σ(-1)^k･a2k･x^(2k)

　　a2k＝1/k!2^k･4/π∫(0,π/4)(sec^2θ-π/4)^kdθ

であって，実際に計算すると

　　a0=1,a2=0,a4=2(π-3)/3π^2,a6=(7π^2-60π+120)/45π^3,･･･

が得られるとある．この積分がどのようにして求められたものか，というのが前回の宿題である．

　

　周辺分布がともに平均０，分散１の正規分布となる２次元正規分布

　　p(x,y)dxdy=1/2π･exp(-(x^2+y^2)/2)dxdy

において，x=rcosθ，y=rsinθと極座標変換すると，ヤコビアンは

　　D(x,y)/D(r,θ)=r

であるから

　　p(x,y)dxdy=rexp(-r2/2)dr*1/2πdθ

　

　このことから，

　　{Φ(x)}^2=(2π)^(-1)∫(-x,x)∫(-x,x)exp{-(t1^2+t2^2)/2}dt1dt2

　　　　　　 <∫(0,r)rexp(-r^2/2)dr∫(0,2π)1/2πdθ

　　　　　　 =∫(0,r)rexp(-r^2/2)dr=1-exp(-r^2/2)

であることは，すでにみてきたとおりである．

　

　ここでｘ軸，(x,0)を通りｙ軸に平行な直線，直線ｙ＝ｘに囲まれた直角二等辺三角形領域で，２次元正規分布を積分することを考える．この積分を△で表すとx=rcosθより，

　　△=1/2π∫(0,π/4)∫(0,x/cosθ)rexp(-r^2/2)drdθ

　　　=1/2π∫(0,π/4){1-exp(-x^2/2cos^2θ)}dθ

と表される．このことさえわかれば，あとは簡単な式の変形だけである．

　

　この三角形領域での積分は正方形領域の積分の1/8であるから，

　　{Φ(x)}^2=8△

　　　　　　 =4/π∫(0,π/4){1-exp(-x^2/2･sec^2θ)}dθ

　　1-{Φ(x)}^2=4/π∫(0,π/4)exp(-x^2/2･sec^2θ)dθ

　

　したがって，

　　f(x)=exp(2x^2/π)[1-{Φ(x)}^2]

　　　　=exp(2x^2/π)4/π∫(0,π/4)exp(-x^2/2･sec^2θ)dθ

　　　　=4/π∫(0,π/4)exp(-x^2/2･sec^2θ+2x^2/π)dθ

　　　　=4/π∫(0,π/4)exp{-x^2(sec^2θ/2-2/π)}dθ

　　　　=4/π∫(0,π/4)Σ{-x^2(sec^2θ/2-2/π)}^k/k!dθ

　　　　=4/πΣ(-1)^k･x^(2k)/k!2^k∫(0,π/4)(sec^2θ-4/π)}^kdθ

となって，

　　f(x)=Σ(-1)^k･a2k･x^(2k)

　　a2k＝1/k!2^k･4/π∫(0,π/4)(sec^2θ-π/4)^kdθ

が得られるというわけである．

　

　具体的な係数値は，

　　a4=.0095642,a6=.0004240,a8=.0000425,･･･････････

と求められているが，これらの数値は後述する議論において必要となる．

　

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

　

【２】２つの級数は等しいか？

　

　２つの級数が等しいとは項ごとに等しいことを指すのであるが，前回のコラムにおいて，

　　f(x)=exp(2x^2/π)[1-{Φ(x)}^2]

　　　　=1+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

　　-2x^2/π･g(x)=log[1-{Φ(x)}^2]

　　　　　　　　 =-2/π･x^2+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

であること，両者に式

　　+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

が共通して現れることを示した．

　

　これに続く項も等しいのだろうか？　この疑問に応えるために，

　　-2x^2/π･g(x)=-2x^2/π+h(x)

　　h(x)=2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

とおいて，

　　f(x)=1+h(x)

が成り立つかどうかを調べてみることにした．

　

　　exp(-2x^2/π)f(x)=[1-{Φ(x)}^2]

の対数をとると

　　-2x^2/π+logf(x)=log[1-{Φ(x)}^2]=-2x^2/π･g(x)

-2x^2/πを移項すると，

　　logf(x)=-2x^2/π(g(x)-1)

となるが，右辺はh(x)より，

　　log{1+h(x)}=h(x)あるいはexp(h(x))=1+h(x)

が成り立てば，

　　f(x)=1+h(x)

が成り立つことになる．

　

　　exp(-2x^2/π)f(x)=[1-{Φ(x)}^2]

　　-2x^2/π･g(x)=log[1-{Φ(x)}^2]

から，級数h(x)がf(x)とg(x)とで項ごとに等しくなることが予想されたのだが，もし，それが正しいのならば，

　　exp(h(x))=1+h(x)+1/2{h(x)}^2+1/6{h(x)}^3+･･･=1+h(x)

とならなければならず，これは不合理である．

　

　このことから，まったく同一ではないことがわかったが，少なくともh(x)が小さい値をとるときには，近似的に

　　exp(h(x))≒1+h(x)

すなわち，

　　f(x)≒1+h(x)

となるということは直ちに理解されるであろう．

　　exp(x)=1+x

が正確に成り立つのはx=0のときだけであるが，ｘが小さければ，

　　exp(x)≒1+x

というわけである．

　

　したがって，Williams=山内の近似式

　　Φ(x)≒{1-exp(-2x^2/π)･(1+2(π-3)/3π^2･x^4)}^(1/2)

と，Cadwellの近似式

　　Φ(x)≒{1-exp(-2x^2/π+2(π-3)/3π^2･x^4)}^(1/2)

は（近似式的に）同値であるし，また，

　　Φ(x)≒[1-exp(-2x^2/π){1+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6}]^(1/2)

は

　　Φ(x)≒[1-exp(-2x^2/π+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6}]^(1/2)

と同値になる．

　

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

　

【３】Cadwellの近似式への疑惑

　

　それでは，

　　f(x)=1+h(x)

　　-2x^2/π･g(x)=-2x^2/π+h'(x)

とおいて，級数h(x)とh'(x)はどの項まで等しいのだろうか？　少なくともx^6の項まで等しくなることは手計算で確認済みである．

　

　　log{1+h(x)}=h'(x)

は成り立つので，左辺を展開することにすると，

　　h(x)-1/2･{h(x)}^2+1/3･{h(x)}^3-1/4{h(x)}^4+･･･=h'(x)

　

　ここで，h(x)は４次式なので，両者はx^8の項で食い違ってくることがわかるが，h'(x)のx^8の項の係数は，この式から

　　a8-1/2･a4^2

と計算される．

　

　前述したように

　　a4=.0095642,a6=.0004240,a8=.0000425

であるから，

　　a8-1/2･a4^2=-.0000032

となり，Cadwellの近似式は正確には，

　　Φ(x)≒{1-exp(-.63662*x^2+.0095642*x^4-.0004240*x^6-.0000032*x^8)}

となるはずである．

　

　しかし，このことは，Cadwellの近似式

　　Φ(x)≒{1-exp(-.63662*x^2+.0095642*x^4-.0005*x^6+.00002*x^8)}

に対して疑念を抱かざるを得ない結果でもあった．なぜなら，x^6の項ですでに係数の値が違っているし，x^8の項では符号の正負さえ異なっているからである．

　

　確認のため，

　　f(x)=exp(2x^2/π)[1-{Φ(x)}^2]

　　　　=1+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

および

　　-2x^2/π･g(x)=log[1-{Φ(x)}^2]

　　　　　　　　 =-2/π･x^2+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+･･･

のx^8の係数を求めてみることにした．これも手計算でやったので，保証はできないが，

　　f(x)=1+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+(9π^3-98π^2+420π-630)/315π^4･x^8-･･･

　　-2x^2/π･g(x)=-2/π･x^2+2(π-3)/3π^2･x^4-(7π^2-60π+120)/45π^3･x^6+(15π^3-280π^2-1400π+2100)/525π^4･x^8+･･･

となり，小生の計算のほうが数値的には合っているようである．

　

　実は，Cadwell（1951）の論旨には理解不能なところ，というよりCadwell自身の無理解と思われるフシが見受けられる．要するに，よくわかっていないのではないかと疑惑を抱かせる論文なのである．論文の再検討が必要となろうが，その点，統計数値表の計算は非常によくできていて，さすが米国においても参考文献として引用されるだけのことはあると思われた次第である．

　

　また，最初から数式処理ソフトを使っていれば，このような無駄な計算は避けて通れたのかもしれないが，たまにはソフトを用いずに計算して，先人の苦労を偲ぶことも意義のあることであろう．

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝