SPM-MLに勉強になる話題が流れていたので、共有します。
出典はこちら。
https://www.jiscmail.ac.uk/cgi-bin/webadmin?A2=SPM;a65dd354.1307
コントロールAと疾患B、疾患Cという3群がある時、Design Matrixに全部を入れてしまった方がいいのか、
もしくは、コントロールと疾患A、コントロールと疾患Bというように、別々にDesign Matrixを作った方がいいのかという質問です。
これに対し、Cyril Pernetが非常にわかりやすいたとえを使って説明をしています。
ただ、少しだけtypoがあったのでそこを補足して説明します。
以下はMatlabを実際に動かしてみながらやってみるといいでしょう。
非常にわかりやすくするために、各群3人だとします。
そして、あるボクセルのデータが以下のようになっているとします。
>> A=[9; 10; 11];
>> B=[19; 20; 21];
>> C=[29; 30; 31];
3群を全部ひとつのモデルに入れることにしましょう。
この時、GLMで考えると、Yはデータ行列になりますので、以下のようになります。
>> Y=[A; B; C]
Y =
9
10
11
19
20
21
29
30
31
この時、Design Matrixは次のようにあらわされます。(今はわかりやすく共変量はすべて1としています。)
>> X=[kron(eye(3), ones(3,1)) ones(9,1)]
X =
1 0 0 1
1 0 0 1
1 0 0 1
0 1 0 1
0 1 0 1
0 1 0 1
0 0 1 1
0 0 1 1
0 0 1 1
GLMでは、Y=XB+Eであらわされます。
そして、Bを求めたい時は、Matlabで以下のコマンドを叩くことで簡単に求まります。
ただ、今、上でBという変数を使ってしまっているので、わかりやすくBETAとします。
>> BETA=pinv(X)*Y
そうすると、Bは以下のように表示されます。
BETA =
-5.0000
5.0000
15.0000
15.0000
今、BETAは上からAの平均値、Bの平均値、Cの平均値、そして定数となります。
なので、AとBの差の絶対値は10となりますし、AとCの差の絶対値は20となります。
それでは、次にAとBの2群をモデルしてみましょう。
今の場合、Yは以下のようになります。
>> Y=[A;B]
Y =
9
10
11
19
20
21
このときのDesign Matrixは以下のようになります。
>> X = [kron(eye(2), ones(3,1)) ones(6,1)]
X =
1 0 1
1 0 1
1 0 1
0 1 1
0 1 1
0 1 1
同じようにBETAを求めましょう。
>> BETA=pinv(X)*Y
BETA =
-0.0000
10.0000
10.0000
今の場合、BETAは上からAの平均値、Bの平均値、定数となります。
Aの平均値とBの平均値の差の絶対値は10となります。上のモデルと同じです。
では、AとCの場合はどうでしょうか。
>> Y=[A;C]
Y =
9
10
11
29
30
31
Xは同じものが使えますので、BETAを求めます。
>> BETA=pinv(X)*Y
BETA =
-3.3333
16.6667
13.3333
今の場合、BETAは上からAの平均値、Cの平均値、定数となりますので、
Aの平均値とCの平均値の差の絶対値は(約ではありますが)20となります。
つまり、このような場合、Design Matrixを複数作るよりも、1つのDesign matrixを作ってしまった方がエレガントだということになります。
今までいつも迷っていたことだったので、すっきりしました。