** 本記事は、統計学(数学)に関する記載が多くなっています。ご了承ください。

先日、IGOR Pro のサポートページに、「Global Fit で係数の信頼区間を算出する方法」にという情報が掲載されました。コマンドによる操作とはいえ、わずか2行。操作はさほど面倒なものではありませんが、フィッティングを利用していると、時々疑問に感じるポイントが含まれてますので、この統計学的な背景について、ちょっとここで整理してみます。

▼ 信頼区間とは

信頼区間とは、母数がどのような範囲に存在するかを推定量と確率分布を使って表現したものです。例えば、直線(y=ax+b) の傾き a の 95%信頼区間が

P (1.5 ≤ a ≤ 2.5) = 0.95

と与えられるとき、傾き a は、同様な操作を何度も繰り返していくと、95%の確率で 1.5~2.5 の値になりますよ、というものです。逆に言うと 5% はこの範囲の外側の値になることがある、そう推定されるというものです。

あまり信頼区間が広いと、フィッティングで得られた傾きの値は信頼できませんよね?だって、この傾きを使った式で予測される値がばらついてしまいますもの。ましてや、傾きの信頼区間に 0 が含まれると、y=ax+b というモデル式を想定したことが、そもそも怪しいなぁ、となるわけです。

もうちょっと踏み込んでみます。

なぜ、傾き a の値がばらついてしまうんでしょう?y=ax+b という1次関数のパラメータ a と b は、最小二乗法によると、解析的に、下記のように唯一無二に定めることができます。

fig01.png

にもかかわらず、なぜ「ばらつき」を考慮する必要があるのか?ここで重要になるのは、母集団と標本という概念です。

▼ 母集団と標本

母集団と標本の話は、全数調査と標本調査を考えるとわかりやすいでしょう。たとえば、ある小学校の4年2組の平均身長を求めるとき、人数は高々40人ですから、全員の身長を測って、平均を算出することは可能です。これはまさに全数調査で、クラス全員(=母集団)の身長から導かれた平均値は、測定の誤差を無視すれば、真の平均値(母平均)です。

では、全国の小学校4年生の平均身長を求めるにはどうしたらよいでしょう?その気になれば、全国の4年生の身長データを集めることは可能かもしれませんが、容易ではありません。こういうときに利用する手法が、標本調査です。「標本」とは、「母集団」から抽出された部分集合ですので、全国の小学校4年生という「母集団」に対して、さきほどのある小学校の4年2組という集団は「標本」に相当しますね。隣の4年1組も1つの「標本」ですし、隣の小学校の4年生もまた別の「標本」です。これらの標本の平均身長は、いずれも全国の小学校4年生の平均身長という真の値に対する推定量であり、標本の取り方で、得られる推定量は当然異なる(=ばらつく)わけです。ですので、「標本」は、無作為に、偏りなく抽出されること(無作為抽出)が好ましいとされるわけです。

y=ax+b というモデル式のパラメータ a と b についても同じことが言えます。何か実験をしたとして、そこで得られるデータセット X, Y は、考えうるすべてのデータ、つまり母集団そのものではありませんので、何かしらの標本であると考えるべきでしょう。そのため、あるデータセットから得られるフィッティングパラメータ a と b は、真の a と b に対する推定量であり、「ばらつき」を考慮するわけですね。以下では、a と b の推定量を α、β とそれぞれ表現して区別することにします。


▼ 想定するばらつき

フィッティングでは、x のばらつきは通常考慮しません(注1)。ある xi に対して、yi だけがばらつきを持つもの(=確率変数)として考えます。ここで、誤差項 ei

fig02.png

とすると、ei も確率変数となります。この誤差項については、

  • 期待値 E (ei) = 0
  • 分散 V (ei) = σ2
  • 共分散 Cov (ei, ej) = 0 (i ≠ j)

を満たすものと仮定します。で、ここからはちょっと 力技  数学なんですが、yi を消去して、α を ei で表現すると、

fig03.png
fig04.png                                  (1)

となります。α の期待値E (α) を考えると、e の期待値 E (ei) が 0 であることから、第2項、3項はいずれもゼロとなるので、

fig05.png

が得られます。さらに、e の分散の仮定(V (ei) = σ2)から、α の分散 V (α) は、

fig06.png

となることがわかります。さて、誤差項 ei が正規分布 N (0, σ2) に従うとさらに仮定すると、α は ei の1次式(式 1)で与えられていますので、同様に、次の正規分布に従うことになります。

fig07.png

これがすごく重要で、α が正規分布に従うことから、α を 標準化した Zα は標準正規分布 N (0, 1) に従います。

fig08.png

ついに、「ばらつき」が姿を現しました。


▼ t 分布の登場、標準誤差、そして信頼区間

まだゴールではありません。a だけでなく σ2 は母数ですので、(多くの場合)知ることができません。このように、異なる未知数が残った状態では、a を推定することができませんので、母分散 σ2 を、標本から計算できる不偏分散 s2で置き換えることを考えます。

fig10.pngただし、fig09.png

di は残差(Residual)と呼ばれ、期待値 E (di) =0 という特徴があります。また、不偏分散の自由度は n-2 となります。


いよいよ佳境です。Zα の σ2 を s2 で置き換えた統計量を t 統計量といいます。

fig11.png        (2)

t 統計量は分母の不偏分散の自由度(n-2)の t 分布に従います。めでたく、a を含む統計量が確率分布に従うことが導かれました。

さて、t 統計量の分母は、係数 α の標準誤差(Standard Error)と呼ばれます。分散のルートですから、標準偏差と呼びたいところですが、推計量の分散の標準偏差は、標準誤差と呼ばれるのがルールのようです。

fig15.png


さらに進めます。 tα の 95%信頼区間は次のように推定できますね。 fig12.png

ここで、t0.025 (n-2) は、自由度 n-2 の t 分布で、上側/下側の累積確率が 0.025 となる値です(両側検定)。例えば、自由度 n-2 が 20 の場合、tα は次のように分布しますので、グラフの青い部分 -2.086~2.086 という範囲に収まることが分かります。

あと一歩。tα を代入して整理すると式 3 が得られます。

fig13.png

fig14.png (3)

式 3 のP()の中身に注目です。a の範囲になっています!!つまり、推定量α とその標準誤差 s.e.(α)、自由度 n-2 の t 分布の値が分かれば、その傾き a の範囲を導くことができるわけです。だから、α と標準誤差 s.e.(α) が必要なんですね。


▼ では、IGOR Pro で

信頼区間を得るために必要な情報はわかりました。IGOR Pro でフィッティングを行えば、α も s.e.(α) も得られます。それぞれ、既定では、

  • 係数:W_Coef
  • 標準誤差:W_Sigma

という名称のウェーブに、定義されたパラメータの順番に保持されます。あと必要なのは、t0.025(n-2) ですね。IGOR Pro では StudentT という関数がこの値を返してくれます。

StudentT(Prob, DegFree)

StudentT は -t ~ t の累積確率が Prob と等しくなる t を返しますので、0.025 ではなく 0.95 と指定して OK です。また、DegFree は不偏分散の自由度を指定しますので n-2 となります。そういうわけで、サポートページに記載されているように、係数の95% 信頼区間を得ることができるわけです。

論より証拠。次のコマンドを実行してみてください。

Make/O/N=20 Data=0.5*x+2+gnoise(1)	// line with Gaussian noise
Display Data
ModifyGraph mode(Data)=3,marker(Data)=19,msize(Data)=2,rgb(Data)=(0,0,65280)
CurveFit line Data /D
print "intercept = ", W_coef[0], "±", W_sigma[0]*StudentT(0.95, V_npnts-2)
print "slope = ", W_coef[1], "±", W_sigma[1]*StudentT(0.95, V_npnts-2)


▼ まとめ

今回は、統計学の教科書を丸写ししたような記事になってしまいましたが、標準誤差、信頼区間の関係、うまく伝わったでしょうか?Global Fitting パネルには、信頼区間を出力させるオプションはまだ実装されていませんが、きっと近い将来実装されるのではないかと思います。

--
井上@技術部でした

参考書籍:東京大学教養学部統計学教室 編(1997)「統計学入門」 東京大学出版

注1:別の機会にxデータのばらつきも考慮したフィッティングについて触れようと思います。

※この記事の内容は執筆者の個人的見解で、ヒューリンクスによる公式情報ではありません。[免責事項]

トラックバック

この記事へのトラックバックURL
http://blog.hulinks.co.jp/cgi/mt/mt-tb.cgi/476
内容に対しての関連性がみられないものは削除する場合があります

フィッティングでは、通常、x 方向にはばらつきなく、真の値が得られていることを仮... 続きを読む

コメントの投稿

Emailアドレスは表示されません。は必須項目です。
ヒューリンクス取り扱い製品の内容や購入に関するお問い合わせはヒューリンクスサイト連絡先へお願いいたします。投稿前にその他の注意事項もご覧ください。

HULINKS サイトの新着情報