ExcelでボックスのM検定の分布がどうなるか乱数で調べる
本日は、ExcelでボックスのM検定の分布がどうなるか乱数で調べたいと思います。ボックスのM検定とは、何グループかの分散共分散行列が等しいかどうかを検定する統計量です。判別分析などでは等分散の仮定が入ります。標準正規分布の乱数を生成させボックスのM検定を実行するプログラムを作成しました。下記のようになります。---------------------------Sub Macro1()'' Macro1 Macro'' Keyboard Shortcut: Ctrl+r' Dim x1(1 To 100, 1 To 2) As Double Dim x2(1 To 100, 1 To 2) As Double Dim i As Long Dim mean1(1 To 2) As Double Dim mean2(1 To 2) As Double Dim S1(1 To 2, 1 To 2) As Double Dim S2(1 To 2, 1 To 2) As Double Dim SW(1 To 2, 1 To 2) As Double Randomize Timer 'SEEDの変更 N = 30 N_ALL = 60 rho = 0.3 p = 2 g = 2 For k = 1 To 10000 '--- 一様乱数生成 --- Randomize For j = 1 To g For i = 1 To N x1(i, j) = Rnd() x1(i, j) = WorksheetFunction.Norm_S_Inv(x1(i, j)) x2(i, j) = Rnd() x2(i, j) = WorksheetFunction.Norm_S_Inv(x2(i, j)) x2(i, j) = rho * x1(i, j) + ((1 - (rho) ^ 2) ^ 0.5) * x2(i, j) Next i Next j '--- 平均 --- For j = 1 To g For i = 1 To N mean1(j) = mean1(j) + x1(i, j) mean2(j) = mean2(j) + x2(i, j) Next i Next j For i = 1 To g mean1(i) = mean1(i) / N mean2(i) = mean2(i) / N Next i '--- 不偏分散・共分散行列 S1 --- For i = 1 To N S1(1, 1) = S1(1, 1) + (x1(i, 1) - mean1(1)) * (x1(i, 1) - mean1(1)) S1(1, 2) = S1(1, 2) + (x1(i, 1) - mean1(1)) * (x2(i, 1) - mean2(1)) S1(2, 1) = S1(1, 2) S1(2, 2) = S1(2, 2) + (x2(i, 1) - mean2(1)) * (x2(i, 1) - mean2(1)) Next i ' 不偏化 For i = 1 To g S1(i, 1) = S1(i, 1) / (N - 1) S1(i, 2) = S1(i, 2) / (N - 1) Next i For i = 1 To N S2(1, 1) = S2(1, 1) + (x1(i, 2) - mean1(2)) * (x1(i, 2) - mean1(2)) S2(1, 2) = S2(1, 2) + (x1(i, 2) - mean1(2)) * (x2(i, 2) - mean2(2)) S2(2, 1) = S2(1, 2) S2(2, 2) = S2(2, 2) + (x2(i, 2) - mean2(2)) * (x2(i, 2) - mean2(2)) Next i ' 不偏化 For i = 1 To g S2(i, 1) = S2(i, 1) / (N - 1) S2(i, 2) = S2(i, 2) / (N - 1) Next i For i = 1 To p For j = 1 To p SW(i, j) = S1(i, j) * (N - 1) + S2(i, j) * (N - 1) Next j Next i 'プールした分散共分散行列 For i = 1 To p For j = 1 To p SW(i, j) = SW(i, j) / (N_ALL - g) Next j Next i tmp = SW(1, 1) * SW(2, 2) - SW(2, 1) * SW(1, 2) tmp = WorksheetFunction.Ln(tmp) tmp = tmp * (N_ALL - g) tmp2 = S1(1, 1) * S1(2, 2) - S1(2, 1) * S1(1, 2) tmp2 = WorksheetFunction.Ln(tmp2) tmp2 = tmp2 * (N - 1) tmp3 = S2(1, 1) * S2(2, 2) - S2(2, 1) * S2(1, 2) tmp3 = WorksheetFunction.Ln(tmp3) tmp3 = tmp3 * (N - 1) '統計量 M = tmp - tmp2 - tmp3 C = 1 - (2 * p * p + 3 * p - 1) / 6 / (p + 1) / (g - 1) * (1 / (N - 1) + 1 / (N - 1) - 1 / (N_ALL - g)) Box_m = M Cells(k, 1) = C * Box_mNext kEnd Sub-----------------------------------A列に10000個ボックスのM検定を実行した結果を出力します。A列を検証するためにExcelに式を組みました。D1からD4は定数です。セルD5=(D1-1)*D2*(D2+1)/2セルF3=G3-0.25でF39までコピーします。セルG3=G2+0.5でG39までコピーします。セルH2:H39を選択し=FREQUENCY(A:A,G2:G39)と入力しCtrl+Shift+Enterを押す。セルI3=CHISQ.DIST(F3,$D$5,FALSE)*10000*0.5でI39までコピーします。上の例題では自由度が3のカイ二乗分布に従います。表の意味は0.5から1の間に乱数が879個存在し理論値が880.16という意味です。グラフにするとデータを減らすとずれが大きくなるのは確認しました。もっとデータを大きくするともっとフィットするのかな。Excelでマクロを実行するには1.開発タブが無い『ExcelでVBAを使うには(開発タブを表示する)』本日は、ExcelでVBAを使うために開発タブを表示させる方法を解説したい。まず、開発タブが表示されているか確認します。表示されていなければ、「ファイル」…ameblo.jp2.マクロを有効にする『Excelでマクロを有効にするには』本日は、Excelでマクロを有効にする設定変更について解説します。「ファイル」を押します。「オプション」を押します。ポップアップ画面で「トラストセンター…ameblo.jp3.上記ソースコードの使い方『Excelでveq02020のVBAソースコードを使うには』本日は、Excelでveq02020のVBAソースコードを使い方を解説したい。「開発」タブに移動します。「マクロの記録」をクリックします。ポップアップ画…ameblo.jp新・金融政策入門 (岩波新書)Amazon(アマゾン)2000年頃からのFRBと日銀の歴史。ドキュメント 強権の経済政策 官僚たちのアベノミクス2 (岩波新書)Amazon(アマゾン)小説に近い。司馬遼太郎まではいかないが。よろしくお願いいたします。以上です!