Mathematica で多変量線形回帰を実行する方法を教えてください.

ここでは,3つの独立変数をもつデータセットに対して,基本的な最小2乗(OLS)回帰を実行する例を示します.オプションの数が多すぎるためここではすべてを挙げることができませんが,以下の例でMathematica 関数のRegressの柔軟性や機能について少しはお分かりいただけるものと思います.


データのインポートと探究のためのプロット

まず,必要なパッケージをロードします.

[Graphics:Images/index_gr_1.gif]

Importは,ファイルC:\Windows\Desktop\example.datからデータを読み込むのに使われます.このデータセットは,以下のモデルを使って前もってMathematica で測定された100の結果から成っています.

        [Graphics:Images/index_gr_2.gif] = [Graphics:Images/index_gr_3.gif] + [Graphics:Images/index_gr_4.gif][Graphics:Images/index_gr_5.gif] + [Graphics:Images/index_gr_6.gif][Graphics:Images/index_gr_7.gif] + [Graphics:Images/index_gr_8.gif][Graphics:Images/index_gr_9.gif] + [Graphics:Images/index_gr_10.gif]
          where epsilon ~ N(0, 3)
               i  = 1, ... , 100
             [Graphics:Images/index_gr_11.gif]  = -8.0
             [Graphics:Images/index_gr_12.gif]  = 4.8
             [Graphics:Images/index_gr_13.gif]  = 0.35
             [Graphics:Images/index_gr_14.gif]  = 1.0          

[Graphics:Images/index_gr_15.gif]

TableFormを使って,最初のいくつかを見てみると,{{[Graphics:Images/index_gr_16.gif], [Graphics:Images/index_gr_17.gif], [Graphics:Images/index_gr_18.gif], [Graphics:Images/index_gr_19.gif]}, {[Graphics:Images/index_gr_20.gif], [Graphics:Images/index_gr_21.gif], [Graphics:Images/index_gr_22.gif], [Graphics:Images/index_gr_23.gif]}, ... }という形式になっていることが分かります.

[Graphics:Images/index_gr_24.gif]
[Graphics:Images/index_gr_25.gif]

個々のリストはデータを転置し,列を選択(この場合1から4)し,各リストへ変数を付加することで抽出できます.

[Graphics:Images/index_gr_26.gif]

これで,それぞれの独立変数と従属変数y との関係をビジュアル化することができます.

[Graphics:Images/index_gr_27.gif]

グラフを見る.


基本的な最小2乗(OLS)回帰を実行する

Regressは,RegressionReport -> { }のオプションのリストの他に,3つのパラメータを取ります.最初のパラメータは,{{[Graphics:Images/index_gr_29.gif], [Graphics:Images/index_gr_30.gif], [Graphics:Images/index_gr_31.gif], [Graphics:Images/index_gr_32.gif]}, {[Graphics:Images/index_gr_33.gif], [Graphics:Images/index_gr_34.gif], [Graphics:Images/index_gr_35.gif], [Graphics:Images/index_gr_36.gif]}, ... }の形式のデータです.2つ目のパラメータは,フィットさせるモデルの指定のリストです.これは,3つの独立変数とy 切片をもつ線形モデルなので,[Graphics:Images/index_gr_37.gif]と指定します.例えば,[Graphics:Images/index_gr_38.gif] = [Graphics:Images/index_gr_39.gif] + [Graphics:Images/index_gr_40.gif][Graphics:Images/index_gr_41.gif] + [Graphics:Images/index_gr_42.gif][Graphics:Images/index_gr_43.gif] + [Graphics:Images/index_gr_44.gif][Graphics:Images/index_gr_45.gif] + [Graphics:Images/index_gr_46.gif][Graphics:Images/index_gr_47.gif] + [Graphics:Images/index_gr_49.gif]のような変数間の相互作用があるモデルは,Regress[example,[Graphics:Images/index_gr_50.gif]]で指定されます.3つ目のパラメータは,行列exampleの列を指定するだけのリストです.

[Graphics:Images/index_gr_51.gif]
[Graphics:Images/index_gr_52.gif]
{ParameterTable ->
  Estimate SE TStat PValue
1 -6.66618 5.03394 -1.32425 0.188565
x1 4.66691 0.208694 22.3625 0.
x2 0.360126 0.0222749 16.1673 0.
x3 1.01063 0.0104686 96.539 0.      ,
RSquared -> 0.991364, AdjustedRSquared -> 0.991094, EstimatedVariance -> 7.84618,

ANOVATable ->
  DF SumOfSq MeanSq FRatio PValue
Model 3 86469. 28823. 3673.51 0.
Error 96 753.233 7.84618    
Total 99 87222.2           }

幸いなことに,パラメータの推定は,このデータを生成するのに使われた実際のパラメータに非常に近いものです.


回帰残差の探求

レポート方法にもいろいろなオプションがありますが,この例では,FitResidualsに焦点を当てます.

[Graphics:Images/index_gr_54.gif]
[Graphics:Images/index_gr_55.gif]

errorsは,先ほどの線形回帰からの残差のベクトルとして定義されました.以下の2つの定義は,残差のヒストグラムと正規分布の確率密度関数(PDF)のグラフィックスオブジェクトを生成します.誤差はおおよそ正常で,予想通り,標準偏差は生成されたデータセットのepsilon の標準偏差と等しくなっています.

[Graphics:Images/index_gr_56.gif]
[Graphics:Images/index_gr_57.gif]

[Graphics:Images/index_gr_58.gif]

最後に,残差に不均一分散あるいは他の問題がないかどうかを目で確かめます.もちろん数学的にテストすることも可能ですが,ここではすべて問題はないようです.

[Graphics:Images/index_gr_59.gif]

グラフを見る.


自分で行う

Mathematica の利点は,ユーザが自分のアルゴリズムを書き,テストを作成し,望みのレベルで自分のデータに対応することができるという点です.例えば,パラメータのベクトルbeta = [Graphics:Images/index_gr_61.gif]を返す最小2乗関数を書くことができます(関数の定義規則については,Mathematica を使って関数を作る方法を教えてください.」をご参照ください).

[Graphics:Images/index_gr_62.gif]

次に,データのリストが独立データポイントの行列になります.

[Graphics:Images/index_gr_63.gif]
[Graphics:Images/index_gr_64.gif]

計算しても,同じ結果となります.

[Graphics:Images/index_gr_65.gif]
[Graphics:Images/index_gr_66.gif]

これにより,最尤推定(MLE)や一般化積率法(GMM)のような新しく洗練された統計テストの実装につながります.このようなテストが必要な方にとっては,Mathematica は常に便利なツールとなります.これが必要ではない方も,MathematicaStatisticsパッケージの中のものに少しコマンドを加えるだけで必要なものになるはずです.