名義尺度の因子は、計画行列においてダミー変数に変換されます。SASシステムのGLMプロシジャでは、すべての水準それぞれに、ダミー変数が1つずつ作成されます。一方、JMPでは、最後の水準を除くすべての水準それぞれに、ダミー変数が1つずつ作成されます。JMPでは、最後の水準に対しては、すべてのダミー変数において値が-1になっています。たとえば、名義尺度の因子であるAに3つの水準があった場合、GLMとJMPでは、それぞれ次のようにコード変換されます。
|
GLM |
JMP |
|||
A |
A1 |
A2 |
A3 |
A13 |
A23 |
A1 |
1 |
0 |
0 |
1 |
0 |
A2 |
0 |
1 |
0 |
0 |
1 |
A3 |
0 |
0 |
1 |
-1 |
-1 |
GLMの線形モデルでは、計画行列の列の間に一次従属性があり、最小2乗法の解は一般逆行列を使用します。この時、A3に対するパラメータがゼロになるような解が選択されます。
JMPにおける線形モデルの計画行列は、欠測セルや偶然などによって一次従属性が生じていない限り、フルランクになります。JMPのパラメータ表現では、最後の水準に対するパラメータは、それ以外の水準におけるパラメータを合計してマイナスをつけると計算することができます。よって、すべての水準にわたってパラメータを合計するとゼロになります。
注: JMPにおいて、名義水準に対するパラメータは、その水準における応答の予測値と、全水準における応答の予測値を平均したものとの差であると解釈できます。
JMPにおいて、各水準に対する計画列は、その水準の指示変数(0/1を含んだインジケータの変数)から最後の水準の指示変数を引けば作成できます。パラメータが上記のように解釈されるのは、このようにコード変換しているためです。
JMPパラメータレポート |
解釈 |
計画行列における値 |
切片 |
全水準の平均 |
1´ |
A[1] |
(A==1) - (A==3) |
|
A[2] |
(A==2) - (A==3) |
GLMとJMPのどちらでも、交互作用の効果は、交差している因子の計画列に属する行すべてにわたって直積を取ったものです。たとえば、GLMコード
PROC GLM;
CLASS A B;
MODEL A B A*B;
からは、次の計画行列が作成されます。
|
|
A |
B |
AB |
||||||||||||
A |
B |
1 |
2 |
3 |
1 |
2 |
3 |
11 |
12 |
13 |
21 |
22 |
23 |
31 |
32 |
33 |
A1 |
B1 |
1 |
0 |
0 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
A1 |
B2 |
1 |
0 |
0 |
0 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
A1 |
B3 |
1 |
0 |
0 |
0 |
0 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
0 |
A2 |
B1 |
0 |
1 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
A2 |
B2 |
0 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
A2 |
B3 |
0 |
1 |
0 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
1 |
0 |
0 |
0 |
A3 |
B1 |
0 |
0 |
1 |
1 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
1 |
0 |
0 |
A3 |
B2 |
0 |
0 |
1 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
1 |
0 |
A3 |
B3 |
0 |
0 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
0 |
1 |
JMPの[モデルのあてはめ]コマンドで列Aと列Bの因子モデルを作成すると、次のような計画行列が作成されます。この行列のA13は前出の行列でのA1-A3に該当し、また、A13B13はこの行列でのA13にB13を掛けたものです。
|
|
A |
B |
|
|
|
|
||
A |
B |
13 |
23 |
13 |
23 |
A13 B13 |
A13 B23 |
A23 B13 |
A23 B23 |
A1 |
B1 |
1 |
0 |
1 |
0 |
1 |
0 |
0 |
0 |
A1 |
B2 |
1 |
0 |
0 |
1 |
0 |
1 |
0 |
0 |
A1 |
B3 |
1 |
0 |
-1 |
-1 |
-1 |
-1 |
0 |
0 |
A2 |
B1 |
0 |
1 |
1 |
0 |
0 |
0 |
1 |
0 |
A2 |
B2 |
0 |
1 |
0 |
1 |
0 |
0 |
0 |
1 |
A2 |
B3 |
0 |
1 |
-1 |
-1 |
0 |
0 |
-1 |
-1 |
A3 |
B1 |
-1 |
-1 |
1 |
0 |
-1 |
0 |
-1 |
0 |
A3 |
B2 |
-1 |
-1 |
0 |
1 |
0 |
-1 |
0 |
-1 |
A3 |
B3 |
-1 |
-1 |
-1 |
-1 |
1 |
1 |
1 |
1 |
JMPではこのようなコード変換が行われるため、交互作用がより少ない水準になっており、必要なメモリが少なく、計算時間も短くなります。
このモデルにおいて、セルの期待値をパラメータで表すと、次のようになります。
|
B1 |
B2 |
B3 |
---|---|---|---|
A1 |
|||
A2 |
|||
A3 |
GLMでは、枝分かれ効果に対しても、交互作用効果と同じコード変換が使われています。そして、枝分かれ効果の検定は、モデルに含まれていない効果を考慮して導出されています。GLMでは、枝分かれ効果は、モデルに含まれていない効果を内在しています。
JMPでは、枝分かれ効果がGLMとは異なる方法でコード変換されます。JMPでは、括弧内に指定された効果が、各グループをグループ分けするために使用されます。括弧内に指定された効果の水準ごとに、括弧外に指定された効果の水準の列が計画行列に作成されます。括弧外に指定された効果の全水準が、括弧内に指定された各水準においてすべて存在している必要はなく、いくつかの水準は抜けていてもかまいません。各水準ごとに考慮されて、計画行列の列は作成されます。
|
|
|
B(A) |
||||||
|
|
|
A1 |
A1 |
A2 |
A2 |
A3 |
A3 |
|
A |
B |
A13 |
A23 |
B13 |
B23 |
B13 |
B23 |
B13 |
B23 |
A1 |
B1 |
1 |
0 |
1 |
0 |
0 |
0 |
0 |
0 |
A1 |
B2 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
0 |
A1 |
B3 |
1 |
0 |
-1 |
-1 |
0 |
0 |
0 |
0 |
A2 |
B1 |
0 |
1 |
0 |
0 |
1 |
0 |
0 |
0 |
A2 |
B2 |
0 |
1 |
0 |
0 |
0 |
1 |
0 |
0 |
A2 |
B3 |
0 |
1 |
0 |
0 |
-1 |
-1 |
0 |
0 |
A3 |
B1 |
-1 |
-1 |
0 |
0 |
0 |
0 |
1 |
0 |
A3 |
B2 |
-1 |
-1 |
0 |
0 |
0 |
0 |
0 |
1 |
A3 |
B3 |
-1 |
-1 |
0 |
0 |
0 |
0 |
-1 |
-1 |
最小2乗平均は、他のすべての因子を何らかの中立な値に固定した状態で、水準の特定の組み合わせに対して予測した値です。連続尺度の説明変数に対しては、標本平均が中立な値として使用されます。名義尺度の因子に対しては、全水準にわたる平均が中立な値として使用されますが、これはJMPのコード変換では常にゼロです。順序尺度の因子には、別の中立な値が使用されます(順序尺度の因子における最小2乗平均)。最小2乗平均は推定できないことがあり、その場合は、「推定不可能」と表示されます。JMPの最小2乗平均は、重み付けされている場合を除き、SAS PROC GLM(Goodnight and Harvey 1978)と一致しています。重み変数が使用されている場合、JMPでは加重平均、SAS PROC GLMでは非加重平均を中立平均としています。
一般にJMPの仮説検定は、SASのプロシジャGLM(タイプIIIとIVの仮説)を始めとする、信頼性の高い多くのプログラムと一致します。この節と次の節で、一致しない点について説明します。
SAS PROC GLMでは、タイプIIIとIVの仮説検定を作成する際、推定可能な関数の一般式を使用し、分析対象である効果、およびそれらの効果に含まれる効果だけに関連する関数が見つけ出されます(Goodnight 1978)。
同様の検定がJMPで作成されます。ただし、GLMとはパラメータ化の方法が異なるため、ある効果を検定するとき、(とりあえずフルランクと仮定して)その効果のパラメータすべてをまとめて検定します。JMPにおける主効果の検定には、交互作用パラメータは含まれません。なぜなら、JMPのコード変換においては、交互作用パラメータは、その主効果と無関係になっているからです。
欠測セルなどの特異性があると、JMPとGLMの検定が違ってきます。これは、いくつかの観点から説明できます。
• JMPの検定は、少なくとも主効果について言えば、最小2乗平均が互いに異なることを検定するのと同じです。いくつかの最小2乗平均が推定できないときは、検定に含まれない比較が出てくるので、自由度が足りなくなります。交互作用の検定では、交互作用から計算される最小2乗平均の値が、該当する効果によって表される周辺パターンと異なるかどうかをJMPは検定しています。
• JMPでは、効果を検定するために、その効果を入れたモデルと、その効果を除いたモデルの誤差平方和が比較されます。JMPでは、そのような解釈ができるようなパラメータ化が採用されています。
• JMPは、Hocking(1985, pp. 80–89, 163–166)に説明のある有効仮説検定を採用しています。ただし、セル平均のパラメータ化ではなく、構造パラメータ化を採用しています。有効仮説検定は、その効果に望ましい仮説から始まり、検定の「できるだけ多く」の部分を含めていきます。もちろん、効果に欠測セルがある場合は、仮説全体を推定することができないので、一部を除外しなければなりません。有効仮説検定では、仮説全体のうち、除外する部分ができるだけ少なく抑えられます。
• 欠損セルがあると、JMPとGLM(および他のプログラム)の仮説検定が異なることがありますが、そのような場合はどちらにしても重要な検定ではありません。交互作用が有意であれば、それに含まれる主効果を検定する必要はなくなります。交互作用が有意でないなら、その交互作用はモデルから除外できます。GLMで採用されている検定には、一意に決めることができない検定もあります。欠損セルのある計画で、水準のラベルを変更すると、GLMのタイプIV検定が変化することがあります。
欠損セルがある場合については、次の節でさらに詳しく解説します。
計画列の間に一次従属性が見られるケースについて考えてみましょう。一次従属性は、JMPのコード変換では、推定に必要なセルに対するデータがないとき、または効果に何らかの交絡または共線性があるときにしか生じません。
一次従属性がある場合、最小2乗法によるパラメータの解は一意でなくなり、一部の仮説が検定不可能になります。JMPでは、計画行列のある列が、モデルにおいて前に指定されている効果の列に一次従属しているとわかったとき、そのパラメータ推定値をゼロに固定します。この場合、レポートには、どのパラメータ推定値がゼロに固定されて、どのパラメータ推定値が推定可能であるかが表示されます。さらに、一次従属性について特異性レポートが作成されます。
特異性が見られる場合は、JMPとGLMとにおいて、仮説検定の結果が異なることがあります。JMPの検定では、周辺平均に関してより厳しい基準を採用しているため、一般にGLMよりも自由度が少なくなります。つまり、JMPの検定では、常にその効果の最小2乗平均ごとの比較が行われますが、GLMの場合は必ずしもそうではありません。
たとえば、交互作用がある2元モデルにおいて、1つの欠損セルがある場合を考えてみます。Aには3水準、Bには2水準があり、A3B2が欠損セルです。
A B |
A1 |
A2 |
B1 |
A1B1 |
A2B1 |
|
---|---|---|---|---|---|---|
A1 B1 |
1 |
0 |
1 |
1 |
0 |
|
A2 B1 |
0 |
1 |
1 |
0 |
1 |
|
A3 B1 |
-1 |
-1 |
1 |
-1 |
-1 |
|
A1 B2 |
1 |
0 |
-1 |
-1 |
0 |
|
A2 B2 |
0 |
1 |
-1 |
0 |
-1 |
|
A3 B2 |
-1 |
-1 |
-1 |
1 |
1 |
欠損 |
各セルの期待値は次のとおりです。
|
B1 |
B2 |
A1 |
||
A2 |
||
A3 |
データのあるセルはすべて期待値が推定可能で、欠損セルの期待値は推定不可能であることは明白です。実際、この例において、欠損セルの期待値は、次のような計画列の線形結合になっています(この情報は、特異性レポートに表示されます)。
B1とB2の最小2乗平均を比較する検定を作成するとしましょう。この例では、上記の表の行平均が最小2乗平均になります。
LSM(B1) =
LSM(B2) =
LSM(B1) - LSM(B2) =
これを見ると、b1パラメータに対する検定は、最小2乗平均が等しいことを検定するのと同じであることがわかります。ただし、b1は推定できないため検定ができず、検定に対する自由度はありません。
今度は、Aの水準ごとの最小2乗平均に対する検定を作成します。
LSM(A1) =
LSM(A2) =
LSM(A3) =
LSM(A1) – LSM(A3) =
LSM(A2) – LSM(A3) =
これらの2つの差は推定不可能です。しかし、もう1つの、欠損セルがないAの2列を比較した差は推定可能です。
LSM(A1) – LSM(A2) =
この差にはパラメータが2つありますが、JMPでこの差を検定するときは、自由度が1である検定が使用されます。
ある差が推定が可能かどうかは、特異性を示す線形結合との内積を取り、それがゼロであることを調べればわかります。
|
特異性を示す |
仮説検定の |
パラメータ |
線形結合 |
線形結合 |
m |
1 |
0 |
a1 |
-1 |
1 |
a2 |
-1 |
-1 |
b1 |
-1 |
0 |
ab11 |
1 |
0 |
ab21 |
1 |
0 |
交互作用があるモデルにおいて、欠損セルが存在すると、主効果の自由度が少なくなります(因子が名義尺度である場合)。言い換えると、最小2乗平均が推定できないと、主効果の検定における自由度が減少します。
それでは、GLMではどうでしょうか。欠損セルがない場合は、GLMとJMPの検定は同じです。欠損セルがない場合は、両方とも、最小2乗平均に対する有効仮説検定になっています。ただし、特異性があると、GLMでは、これらのセルがタイプIIIまたはIVのどちらかによって異なる方法で除外されます。タイプIVでは、見つけられる推定可能な線形結合が検索されます。これらの検定は必ずしも一意ではなく、水準の順序を変えると異なる結果が出る可能性があります。タイプIIIでは、一意の検定にするため、推定可能な関数が正規直交化されます。ただし、この検定は、セル平均の解釈にはあまり役立たないことがあります。
JMPの方法は、GLMより自由度が少なくなることが最初は欠点のように思われるかも知れませんが、次のような長所があります。
1. JMPの検定は、本質的に最小2乗平均と関連しています。
2. JMPの検定は、JMPの方法でパラメータ化されたモデルでの平方和の差から計算できます。そのため、検定の計算が簡単です。
3. JMPの検定は、Hockingの「有効仮説検定」と一致します。
4. JMPの検定は完全な周辺検定(whole mariginal test)なので、交互作用内において比較できる効果がすべて検討されます。
最後の点について、もう少し詳しく説明しましょう。前出のA3B2が欠損セルである例で、セル平均の期待値をプロットしたグラフを検討してみます。
図A.1 セル平均の期待値
グラフには、セル平均の期待値と欠損セルが表示されています。A1セルとA2セルの平均は、Bの水準ごとにプロットされています。JMPの方法では、A3B2セルが欠損しているときは主効果Bが検定できません。欠損セルの平均は、交互作用項が許す範囲でどのような値でも取り得るからです。欠損セルの平均が、グラフで示す高い方の値であった場合には、主効果Bは有意になるかもしれませんし、低い方の値であった場合には、有意でなくなるかもしれません。つまり、有意か有意でないか判断できないのです。それが、最小2乗平均が推定できない場合に最小2乗平均から結論できることであり、効果の仮説からも同じこと(つまり、判断できないということ)が結論できなければなりません。
セルの一部分から計算される周辺平均に基づき仮説を検定するのは、GLMのなかではタイプIV検定です。JMPでタイプIV検定を行うには、効果を分けてたくさん計算を行うか、対比を使って自分で検定を作成しなければなりません。
欠測セルがある場合には、GLMでは、JMPとは異なる方法で処理が行われ、違う仮説検定が実行されます。特に、JMPでは区別される枝分かれ効果と交差効果が、GLMでは区別されません。ある計画が、A、B(A)、C(A B)の3層に枝分かれしていて、それぞれに異なる数の水準があるとします。
図A.10は、主効果AをGLMにおけるパラメータ化で検定したものです。第1の列セットはJMPの検定をGLMにおけるパラメータ化で検定したもの、第2の列セットはGLMのタイプIVで行われる検定です。第3の列セットはJMPの検定のものですが、GLMの検定と比較できる形にするため、次の行列を掛け合わせてあります。
最後の列セットは、GLMのタイプIII検定です。これらの検定の違いは、効果の水準に対してどのように係数を配分しているかにあります。JMPでは階層構造がトップダウン型で、GLMのタイプIV検定では、ボトムアップ型のように見えます。ほとんどの場合、実際の検定統計量にはあまり差がありません。
パラメータ | JMP 検定A | GLMタイプIV検定 A | JMP 回転 検定 | GLMタイプIII検定 A | ||||
u | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
a1 | 0.6667 | -0.3333 | 1 | 0 | 1 | 0 | 1 | 0 |
a2 | -0.3333 | 0.6667 | 0 | 1 | 0 | 1 | 0 | 1 |
a3 | -0.3333 | -0.3333 | -1 | -1 | -1 | -1 | -1 | -1 |
|
|
|
|
|
|
|
|
|
a1b1 | 0.1667 | -0.0833 | 0.2222 | 0 | 0.25 | 0 | 0.2424 | 0 |
a1b2 | 0.1667 | -0.0833 | 0.3333 | 0 | 0.25 | 0 | 0.2727 | 0 |
a1b3 | 0.1667 | -0.0833 | 0.2222 | 0 | 0.25 | 0 | 0.2424 | 0 |
a1b4 | 0.1667 | -0.0833 | 0.2222 | 0 | 0.25 | 0 | 0.2424 | 0 |
|
|
|
|
|
|
|
|
|
a2b1 | -0.1667 | 0.3333 | 0 | 0.5 | 0 | 0.5 | 0 | .5 |
a2b2 | -0.1667 | 0.3333 | 0 | 0.5 | 0 | 0.5 | 0 | .5 |
|
|
|
|
|
|
|
|
|
a3b1 | -0.1111 | -0.1111 | -0.3333 | -0.3333 | -0.3333 | -0.3333 | -0.3333 | -0.3333 |
a3b2 | -0.1111 | -0.1111 | -0.3333 | -0.3333 | -0.3333 | -0.3333 | -0.3333 | -0.3333 |
a3b3 | -0.1111 | -0.1111 | -0.3333 | -0.3333 | -0.3333 | -0.3333 | -0.3333 | -0.3333 |
|
|
|
|
|
|
|
|
|
a1b1c1 | 0.0833 | -0.0417 | 0.1111 | 0 | 0.125 | 0 | 0.1212 | 0 |
a1b1c2 | 0.0833 | -0.0417 | 0.1111 | 0 | 0.125 | 0 | 0.1212 | 0 |
a1b2c1 | 0.0556 | -0.0278 | 0.1111 | 0 | 0.0833 | 0 | 0.0909 | 0 |
a1b2c2 | 0.0556 | -0.0278 | 0.1111 | 0 | 0.0833 | 0 | 0.0909 | 0 |
a1b2c3 | 0.0556 | -0.0278 | 0.1111 | 0 | 0.0833 | 0 | 0.0909 | 0 |
a1b3c1 | 0.0833 | -0.0417 | 0.1111 | 0 | 0.125 | 0 | 0.1212 | 0 |
a1b3c2 | 0.0833 | -0.0417 | 0.1111 | 0 | 0.125 | 0 | 0.1212 | 0 |
a1b4c1 | 0.0833 | -0.0417 | 0.1111 | 0 | 0.125 | 0 | 0.1212 | 0 |
a1b4c2 | 0.0833 | -0.0417 | 0.1111 | 0 | 0.125 | 0 | 0.1212 | 0 |
|
|
|
|
|
|
|
|
|
a2b1c1 | -0.0833 | 0.1667 | 0 | 0.25 | 0 | 0.25 | 0 | 0.25 |
a2b1c2 | -0.0833 | 0.1667 | 0 | 0.25 | 0 | 0.25 | 0 | 0.25 |
a2b2c1 | -0.0833 | 0.1667 | 0 | 0.25 | 0 | 0.25 | 0 | 0.25 |
a2b2c2 | -0.0833 | 0.1667 | 0 | 0.25 | 0 | 0.25 | 0 | 0.25 |
|
|
|
|
|
|
|
|
|
a3b1c1 | -0.0556 | -0.0556 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 |
a3b1c2 | -0.0556 | -0.0556 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 |
a3b2c1 | -0.0556 | -0.0556 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 |
a3b2c2 | -0.0556 | -0.0556 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 |
a3b3c1 | -0.0556 | -0.0556 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 |
a3b3c2 | -0.0556 | -0.0556 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 | -0.1667 |
JMPのパラメータ化では、Aに対する検定は次のようになります。
パラメータ | GLM-IV検定 | JMP検定 | ||
m | 0 | 0 | 0 | 0 |
a13 | 2 | 1 | 1 | 0 |
a23 | 1 | 2 | 0 | 1 |
a1:b14 | 0 | 0 | 0 | 0 |
a1:b24 | 0.11111 | 0 | 0 | 0 |
a1:b34 | 0 | 0 | 0 | 0 |
a2:b12 | 0 | 0 | 0 | 0 |
a3:b13 | 0 | 0 | 0 | 0 |
a3:b23 | 0 | 0 | 0 | 0 |
a1b1:c12 | 0 | 0 | 0 | 0 |
a1b2:c13 | 0 | 0 | 0 | 0 |
a1b2:c23 | 0 | 0 | 0 | 0 |
a1b3:c12 | 0 | 0 | 0 | 0 |
a1b4:c12 | 0 | 0 | 0 | 0 |
a2b1:c13 | 0 | 0 | 0 | 0 |
a2b2:c12 | 0 | 0 | 0 | 0 |
a3b1:c12 | 0 | 0 | 0 | 0 |
a3b2:c12 | 0 | 0 | 0 | 0 |
a3b3:c12 | 0 | 0 | 0 | 0 |
GLMの検定では、a1b24パラメータに係数が使われ、JMPから見ると少し奇妙な感じがします。