12.6 周期共分散分析

(1) 周期共分散分析の原理

時系列データに限らず関連性を持つ一連のデータに特定の関数を当てはめた時、その関数が表すグラフにはレベルを表すパラメーターとパターンを表すバラメーターが存在します。 例えば一連のデータに周期関数を当てはめた周期回帰曲線ではメサーがレベルを表し、振幅と位相がパターンを表します。 また一連のデータに直線y = α + βxを当てはめた回帰直線では切片a(定数)がレベルを表し、回帰係数β(傾き)がパターンを表します。

回帰直線を当てはめる群が複数ある時、共分散分析を利用して各群の回帰直線が平行かどうか、つまり回帰係数の値が同じかどうかを検討することができます。 そして各群の回帰直線が平行なら各群の目的変数の修正平均値——回帰直線を利用して説明変数の影響を取り除いた平均値——を比較することができます。 その場合、各群の修正平均値の差は各群の平行な回帰直線の切片の差と一致するので、修正平均値の比較は平行な回帰直線の切片の比較に相当します。 したがって共分散分析は複数の回帰直線のレベルとパターンを比較する手法と考えることができます。 (→第8章　共分散分析)

この共分散分析の原理を周期回帰曲線に応用すると、複数の群における周期回帰曲線のレベル(メサー)とパターン(振幅と位相)を比較することができます。 その手法を周期共分散分析(PERCOVA：periodic analysis of covariance)と呼ぶことにします。

(2) 単純な周期共分散分析

表12.1.1を疾患Aの2名の被験者のデータとし、これに疾患Bの3名の被験者のデータを追加したものが表12.6.1だとします。

表12.6.1のA群とB群の平均値について、基本周期を24時間とし、周期成分1と周期成分2だけの周期回帰式を当てはめると次のようになります。

図12.6.1において青色の折れ線と青色の曲線がA群の平均値の変動と周期回帰曲線であり、赤色の折れ線と赤色の曲線がB群の平均値の変動と周期回帰曲線です。 このグラフと上記の計算結果から、2群の周期回帰曲線のメサーはあまり変わらないものの、B群の振幅はA群の半分以下であり、位相が3時間ほど早いことがわかります。 つまり2群の周期回帰曲線はレベルはあまり変わらず、パターンが異なっていると解釈できます。

このデータに周期共分散分析を適用すると次のようになります。 _(注1)

「群差」は2群の24時間平均値の差の検定であり、周期回帰式を当てはめない時の2群のレベルの検定になります。 これは周期共分散分析の主目的ではなく、参考程度です。

「共通周期」は2群の周期回帰曲線が平行と仮定した時の周期回帰の検定、つまり平行な周期回帰曲線の振幅が0かどうかの検定です。 この検定結果が有意ではない時は2群のデータに周期回帰式を当てはめることが妥当ではない、つまり2群の周期回帰曲線のレベルとパターンを比較するのは無意味ということになります。 このデータの場合は有意水準5％で有意ですから、2群の周期回帰曲線のレベルとパターンを比較することが可能です。

「修正群差」は2群の周期回帰曲線が平行と仮定した時のメサーの差の検定つまりレベルの検定です。 原則として、この検定は2群の周期回帰曲線がほぼ平行の時だけ正確な結果になります。 このデータの場合は2群の測定時期が同じため群差の検定結果と同一になります。 これは2群の測定時期が同じ時は2群の24時間平均値の差とメサーの差が同一になるからです。

また測定間隔が等間隔の時は24時間平均値と周期回帰式のメサーが同じ値になります。 このデータの場合は2群の測定時期が同じで、しかも等間隔ですから、群差と修正群差の検定結果が一致し、しかも24時間平均値と周期回帰式のメサーが同じ値になります。 このような時は2群の周期回帰曲線が平行ではなくても検定結果が正確になります。

「全体周期」は2群を合わせて周期回帰式を計算した時の周期回帰の検定、つまり全体の周期回帰曲線の振幅が0かどうかの検定です。 これは周期共分散分析の主目的ではなく、参考程度です。 このデータの場合は2群の測定時期が同じため、全体の周期回帰曲線の振幅と位相が平行な周期回帰曲線の振幅と位相に一致し、共通周期の検定結果と一致します。

表12.6.2の周期共分散分析表において、共通周期と修正群差の間に隙間があり、全体周期と非平行性の間にも隙間があります。 これは群差と共通周期の平方和の合計と、修正群差と全体周期の平方和の合計が一致するからです。 このことはデータの変動のうちのある部分については2群の平均値の違いによる変動と平行な周期回帰曲線による変動を合わせたものと解釈できると同時に、修正平均値つまりメサーの違いによる変動と全体の周期回帰曲線による変動を合わせたものとも解釈できることを意味しています。 これらの事情は普通の共分散分析と同様です。

「非平行性」は2群の周期回帰曲線の振幅と位相の差の検定つまりパターンの検定です。 この検定結果が有意の時は2群の周期回帰曲線が平行ではなく、曲線のパターンが異なっていることになります。 このデータの場合は有意水準5％で有意であり、図12.6.1からわかるように2群の周期回帰曲線のパターンが少し異なっています。

表12.6.1のデータについて、群ごとの平均値ではなく個々の被験者のデータに周期共分散分析を適用することもできます。 しかし個々の被験者は周期的な活動以外に突発的な活動もしていることが多く、データに周期変動だけでなく突発的な変動も含まれている可能性があります。 そのため周期回帰分析と同様に、個々の被験者に周期共分散分析を適用すると往々にして周期回帰曲線がうまく適合しないことがあります。

個々の被験者のデータを平均すると突発的な変動が小さくなり、群に共通する変動——例えば日内変動のような変動が浮かび上がります。 そのため個々の被験者のデータに周期共分散分析を適用するよりも、平均値に適用する方が周期回帰曲線がうまく適合する可能性が高く、結果の解釈も妥当なものになる可能性が高くなります。

(3) 周期回帰曲線の一致係数

普通の共分散分析と同様に、周期共分散分析における検定はほとんど有意性検定になります。 そのため検定結果よりも、周期回帰式そのものや重寄与率を科学的に検討する方が有意義です。 しかし周期回帰曲線の非平行性つまりパターンがどの程度異なっているかについて、周期回帰式から判断するには色々と問題があります。

回帰直線のパターンは回帰係数つまり傾きだけなので、「2群の回帰直線のパターンが異なっている」ということは「2群の回帰直線の傾きが異なっている」という意味になります。 しかし周期回帰曲線のパターンには振幅と位相という2つのパラメーターがあるので、「2群の周期回帰曲線のパターンが異なっている」といっても振幅と位相のどちらが異なっているのかわかりません。 そこで、2群の周期回帰曲線の振幅と位相が一致している程度を表す指標があると便利です。 そしてそれらの指標は関数の相関係数と一致係数の考え方に基づいて導くことができます。

詳しい説明は(注2)を見ていただくとして、関数の相関係数は2つの関数の大小関係が同じかどうかの指標であり、周期回帰曲線の場合は2つの周期回帰曲線の位相が一致しているかどうかを表す位相の一致係数になります。 一方、関数の一致係数は2つの関数の値が同じかどうかの指標であり、周期回帰曲線の場合は2つの周期回帰曲線のメサーが一致するように平行移動した時に、周期回帰曲線が一致するかどうかを表す周期回帰曲線全体の一致係数になります。 そしてそれら2つの一致係数から、2つの周期回帰曲線の振幅が一致しているかどうかを表す振幅の一致係数を導くことができます。

表12.6.1のデータについて、それらの一致係数を計算すると次のようになります。 _(注2)

これらの一致係数から、A群とB群の周期回帰曲線は位相と振幅が同じ程度にずれていて、全体としてはかなりずれていることがわかります。

(4) 二元配置型周期共分散分析

表12.6.1は繰り返し測定した二元配置型のデータですから、繰り返し測定型二元配置分散分析と周期共分散分析を組み合わせて、より詳細な分析をすることができます。 その手法を二元配置型周期共分散分析(two-way layout periodic analysis of covariance)と呼ぶことにします。 それを表12.6.1のデータに適用すると次のようになります。 _(注3)

「群差」は2群の24時間平均値の差の検定であると同時に、2群の周期回帰曲線のメサーの差の検定つまりレベルの検定でもあります。 繰り返し測定した二元配置型のデータは2群の測定時期が必ず同じになるので、2群の24時間平均値の差とメサーの差が同一になります。 そのため群差の検定結果とメサーの差つまり修正群差の検定結果が必ず一致します。

「個体残差」は群差の検定のための誤差であり、群を要因として24時間平均値を一元配置分散分析法で比較する時の残差に相当します。 「個体」は被験者ごとの24時間平均値がばらついているかどうかの検定です。 これは個体差を取り除いて効率的な分析を行うためのものであり、周期共分散分析の主目的ではありません。 「個体残差」と「個体」の間に隙間があるのは、群差と個体残差の平方和を合計したものが個体の平方和になるからです。 そしてこの部分は個体の変動を全変動とし、群を要因Aとした一元配置分散分析に相当します。 これは繰り返し測定型二元配置分散分析と同様です。

「全体周期」は2群を合わせて周期回帰式を計算した時の周期回帰の検定つまり全体の周期回帰曲線の振幅が0かどうかの検定であると同時に、2群の周期回帰曲線が平行と仮定した時の周期回帰の検定つまり共通周期の検定でもあります。 繰り返し測定した二元配置型のデータは2群の測定時期が必ず同じになるので、全体の周期回帰曲線の振幅と位相が平行な周期回帰曲線の振幅と位相に一致し、共通周期の検定結果と一致します。

「非平行性」は2群の周期回帰曲線の振幅と位相の差の検定つまりパターンの検定です。 この検定結果が有意の時は2群の周期回帰曲線が平行ではなく、パターンが異なっていることになります。

「ズレ合計」は2群の時間ごとの平均値の変動と周期回帰曲線とのズレ(LOF：Lack Of Fit)を合計した値が0かどうかの検定です。 この検定結果が有意の時は時間ごとの平均値の変動と周期回帰曲線のズレが大きいことになります。

表12.6.3より平均値の変動と周期回帰曲線のズレはあまり大きくなく、平均値の変動を周期回帰曲線でうまく近似できることがわかります。 そして平均値だけを用いた単純な周期共分散分析結果と同様に、2群の周期回帰曲線のメサーはあまり変わらず、パターンだけが異なっていることもわかります。

このように二元配置型周期共分散分析ではズレの検討を行うことができるので、できるだけ少ない周期成分で、できるだけ信頼性の高い合理的な周期回帰曲線を用いてレベルとパターンの検討することが可能です。 これがこの手法の特徴であり、平均値を用いた単純な周期共分散分析よりも優れた点です。 そして実際の研究現場で得られる時系列データは繰り返し測定した二元配置型データが多いので、この手法が最も有用だと思います。

(5) 一元配置型周期共分散分析

繰り返し測定型二元配置分散分析と同様に、二元配置共分散分析では全ての測定時点のデータが揃っている被験者だけが解析対象になります。 しかし測定時点が多いとどうしても欠測値が生じやすく、解析対象から除外される被験者が出てくる可能性が高くなります。 また群によって測定間隔が異なり、測定時点数が異なっている場合も有り得ます。

そのような時は第3節の分散分析を応用した周期回帰分析と同様に、データを同一の被験者で連続測定されたものと考えず、測定時点ごとに独立した被験者で測定された一元配置型のデータと考え、一元配置分散分析を応用した周期共分散分析を適用することができます。 その手法を一元配置型周期共分散分析(one-way layout periodic analysis of covariance)と呼ぶことにします。

表12.6.1のデータを測定時点ごとに独立した被験者から得られた一元配置型のデータと考え、周期成分1と周期成分2を用いた一元配置型周期共分散分析を適用すると次のようになります。 _(注4)

「群差」、「共通周期」、「修正群差」、「全体周期」、「非平行性」の解釈は平均値を用いた単純な周期回帰分析と同様であり、「ズレ合計」の解釈は二元配置型周期回帰分析と同様です。

共通周期と修正群差の間に隙間があるのは、群差と共通周期の平方和の合計と、修正群差と全体周期の平方和の合計が一致するからです。 これは単純な周期共分散分析と同様です。 また全体周期と非平行性の間と、ズレ合計と時期の間にも隙間があります。 これは群差＋共通周期＋非平行性+ズレ合計、または修正群差＋全体周期＋非平行性+ズレ合計の平方和の合計が時期の平方和に一致するからです。 そして「時期」以下は時期を要因Aにした一元配置分散分析に相当します。

表12.6.1のデータは例数が少なく、しかも個人差が比較的小さいので、二元配置型周期共分散分析と一元配置型周期共分散分析の結果はあまり変わりません。 しかし普通はもっと例数が多く、しかも個人差が大きいので多少の欠測値があったとしても二元配置型周期共分散分析の方が効率的です。 一元配置型周期共分散分析を適用した方が良い場合は、本当に各時点ごとに独立した被験者から得られたデータか、あるいは測定間隔が被験者ごとに異なるので二元配置型にすると欠測値が非常に多くなってしまうデータを解析する時だけです。

ただし各時点ごとに独立した被検者から得られたデータの時または欠測値が多い時は、平均値の変動が時間によるものか被検者の違いによるものか厳密には区別できません。 そのため結果の解釈は慎重に行う必要があります。 (→4.3 繰り返しのある多標本・多時期の計量値　(6) 欠測値の処理方法)

以上のような3種類の周期共分散分析のどれを適用するにしても、次のような周期回帰分析の一般的条件を満たしている必要があります。

1. 周期回帰式に組み込んだ周期成分が科学的に解釈可能か？
2. 周期回帰式のメサー、振幅、位相が科学的に納得のいくものか？
3. 周期回帰式の重寄与率が大きな値か？

周期回帰分析と同様に、この場合も時系列データの変動を一般的な関数――例えば1次時間数等――で回帰し、共分散分析を行うことができます。 そこで時系列データの変動を一般的な関数で回帰して共分散分析を行う手法を時系列共分散分析と呼ぶことにすると、周期共分散分析はその一種ということになります。 (→第8章　共分散分析)

(注1)　平均値だけを用いた単純な周期回帰分析は、群ごとの平均値ではなく被験者ごとのデータに対して適用することもできます。 被験者ごとのデータに対して適用する時は、以下の説明の「群」を「被験者」と読み替えてください。

表12.6.1の群ごとの平均値を一般化すると次のようになります。 時期t₁〜t_(pi)は全ての群で同一でもかまいませんし、群ごとに異なっていてもかまいません。

表12.6.5 周期共分散分析の
一般的データ

群	時期					平均
	t1	…	tj	…	t(pi)
1	y11	…	y1j	…	y1(p1)	m1.
:	:	…	:	…	:	:
i	yi1	…	yij	…	yi(pi)	mi.
:	:	…	:	…	:	:
a	ya1	…	yaj	…	ya(pa)	ma.
全体	m.1	…	m.j	…	m.(pi)	mT

周期共分散分析ではデータy_ijに対して次のような周期回帰式を当てはめて考えます。

○群別周期回帰モデル：群ごとに計算した周期回帰式

k = 1,…,m ≦ p/2　　j = 1,…,p_i (群ごとに異なっていても良い)
：群別周期回帰式による推定値 (i = 1,…,a、j = 1,…,p_i) 　　A_0i：A_i群のメサー、群によって異なる (i = 1,…,a)
A_ki：A_i群の振幅、群によって異なる (k = 1,…,m、i = 1,…,a) 　　θ_ki：A_i群の位相、群によって異なる (k = 1,…,m、i = 1,…,a)
○共通周期回帰モデル：各群の周期回帰曲線が平行と仮定して計算した周期回帰式

：共通周期回帰式による推定値 (i = 1,…,a、j = 1,…,p_i) 　　A_c0i：A_i群のメサー、群によって異なる (i = 1,…,a)
A_ck：全群共通の振幅 (k = 1,…,m) 　　θ_ck：全群共通の位相 (k = 1,…,m)
○全体周期回帰モデル：全群を合わせて計算した周期回帰式

：全体周期回帰式による推定値 (j = 1,…,p_i) 　　A_T0：全体のメサー
A_Tk：全体の振幅 (k = 1,…,m)　　θ_Tk：全体の位相 (k = 1,…,m)

以上の周期回帰式に共分散分析の原理を適用します。 まず3通りの周期回帰式による推定値を用いてデータy_ijを3通りに分解し、その基本式に対応する平方和と自由度を求めると次のようになります。

群別周期回帰：
全体周期回帰：
共通周期回帰：
○全体変動：(y_ij - m_T)
総例数： 　　平方和： 　　自由度：φ_T = N - 1
○群別周期回帰変動：
平方和：
　　β_1i = a_1i 　　β_2i = b_1i　…　β_(2k-1)i = a_ki 　　β_(2k)i = b_ki　…　β_(2m-1)i = a_mi　β_(2m)i = b_mi
x_j0 = 1
　　
　　：
　　
　　：
　　

　　 　　 　　

　　 　　 　　 　　
tan^-1(x)関数は-π/2〜π/2の値を返すので、次の図のようにしてθ_kを求める。

　　 　　
　　


自由度： (2m = p_iの時は(2m-1)a)
○共通周期回帰変動：
平方和：



　　自由度：φ_βc = 2m (2m = pの時は(2m - 1))
○全体周期回帰変動：
平方和：
　　
　　 　　自由度：φ_βT = 2m (2m = pの時は(2m - 1))
○群差変動：(m_i. - m_T)
平方和： 　　自由度：φ_A = a - 1
○修正群差変動：
平方和： 　　自由度：φ_AA = a - 1
○非平行性変動：
平方和：
自由度： (2m = pの時は(2m - 1)(a - 1))
○残差変動：
平方和： 　　自由度：

これらの変動の関係をグラフで表すと図12.6.2のようになり、平方和間の関係を模式図にすると図12.6.3のようになります。 そして平方和と自由度を表12.6.6のような周期共分散分析表にまとめます。

表12.6.6 周期共分散分析表(PERCOVA table)

要因	平方和	自由度	平均平方和(分散)	F値
群差	SA	φA	VA=SA/φA	FA=VA/VR
共通周期	Sβc	φβc	Vβc=Sβc/φβc	Fβc=Vβc/VR
修正群差	SAA	φAA	VAA=SAA/φAA	FAA=VAA/VR
全体周期	SβT	φβT	VβT=SβT/φβT	FβT=VβT/VR
非平行性	SD	φD	VD=SD/φD	FD=VD/VR
残差	SR	φR	VR=SR/φR
全体	ST	φT

○群差の検定
帰無仮説 H₀：m₁. = … = m_a. = m_T
F_A ＞ F(φ_A,φ_R,α)の時、有意水準100α％で有意
○共通周期の検定
帰無仮説 H₀：β_c1 = … = β_c(2m) = 0
F_βc ＞ F(φ_βc,φ_R,α)の時、有意水準100α％で有意
○修正群差の検定
帰無仮説 H₀：β_c01 = … = β_c0a
F_AA ＞ F(φ_AA,φ_R,α)の時、有意水準100α％で有意
○全体周期の検定
帰無仮説 H₀：β_T1 = … β_T(2m) = 0
F_βT ＞ F(φ_βT,φ_R,α)の時、有意水準100α％で有意
○非平行性の検定
帰無仮説 H₀：β₁₁ = … = β_1a = β_c1、… 、β_(2m)1 = … = β_(2m)a = β_c(2m)
F_D ＞ F(φ_D,φ_R,α)の時、有意水準100α％で有意

表12.6.1のA群とB群の平均値について、基本周期を24時間とし、周期成分1と周期成分2を含んだ周期回帰モデルを用いて実際に計算してみましょう。

○A群の群別周期回帰式 → 第2節の(注2)参照
○B群の群別周期回帰式

A₀ = 130.542　　A₁ = 9.97588　　A₂ = 7.76766 　　θ₁ = 2.98392(170.966°) 　　θ₂ = 4.18886(240.004°)
　　
S_T = 821533 - 48×129.74² ≒ 13579.9 　　φ_T = 48 - 1 = 47

S_βc = (-15.8875)×(-381.299) + (-5.67313)×(-136.155) + (-2.61483)×(-62.7559) + (-8.5352)×(-204.848) ≒ 8742.84
φ_βc = 4

S_βT = (-15.8875)×(-381.299) + (-5.67313)×(-136.155) + (-2.61483)×(-62.7559) + (-8.5352)×(-204.848) ≒ 8742.84
φ_βT = 4
S_A = 24×128.938² + 24×130.542² - 48×129.74² ≒ 30.8802 　　φ_A = 1
S_AA = 30.8802 + 8742.84 - 8742.84 = 30.8802 　　φ_AA = 1
S_D = 10991.9 - 8742.84 = 2249.13 　　φ_D = 8 - 4 = 4
S_R = 13579.87 - 30.8802 - 10991.9 = 2557.09 　　φ_R = 48 - 2×(1+2×2) = 38

これらの統計量を用いて表12.6.2の周期共分散分析表を作成することができます。

(注2)　第6章第2節の(注2)で説明したように、2つの変数y₁とy₂の相関係数をベクトル表示すると次のようになります。

y₁の偏差ベクトル： 　　y₂の偏差ベクトル：
y₁とy₂の相関係数：

ここでy₁とy₂を周期関数とすると、その相関係数は次のようになります。

　　
ｄ₁ = y₁ - m₁ = Ｄ_xβ₁^* 　　ｄ₂ = y₂ - m₂ = Ｄ_xβ₂^*

測定間隔が等間隔の時は[Ｄ_x'Ｄ_x]が対角化し、その対角成分が(n/2)になります。 そのため相関係数は次のように2つのフーリエ係数ベクトルβ₁^*とβ₂^*を標準化してお互いに正射影した時の影の射影、つまりβ₁^*とβ₂^*のなす角θの余弦になります。

第4節の(注1)で説明したようにフーリエ係数ベクトルを極座標表示した時の偏角が位相になり、動径つまりベクトルの大きさが振幅になります。 そのため2つの周期関数の相関係数は位相が一致している程度を表す位相一致係数と解釈することができます。

また第5節で説明したように振幅を平方した値をパワーと呼び、周期成分ごとにパワーをプロットしたものをパワースペクトルといいます。 さらに2つの周期関数のフーリエ係数同士を掛け合わせ、それを周期成分ごとにプロットしたものをクロススペクトル(cross spectrum)といいます。

そして2つの周期関数の周期成分kについて、そのクロススペクトルの平方をそれぞれのパワースペクトルの積で割った値のことをコレーレンス(Coherence、関連度)と呼び、それを周期成分の関数としたものをコヒーレンス関数(Coherence function)と呼びます。 これらは2つの周期関数の位相の類似性を表す指標であり、2つの波の干渉しやすさを表す値になります。

周期成分kのコヒーレンス：

フーリエ係数ベクトルの大きさの平方つまりフーリエ係数ベクトル自身の内積はパワースペクトルの合計を表すことになり、2つのフーリエ係数ベクトルの内積はクロススペクトルの合計を表すことになります。 このことから2つの周期関数の相関係数(位相一致係数)を平方した値つまり寄与率は、2つの周期関数のクロススペクトル合計の平方をそれぞれの周期関数のパワースペクトル合計の積で割った値と解釈することができ、周期関数全体のコヒーレンスを表す指標になることがわかります。

一方、第6章第2節の(注2)で説明したように、2つの周期関数y₁とy₂の一致度を表すエーベルの級内相関係数は次のようになります。

最後の式のrの後ろの項は、2つの周期関数のフーリエ係数ベクトルの大きさの平方の幾何平均値と算術平均値の比になっています。 この値は2つの周期関数のフーリエ係数ベクトルの大きさが一致している時は1になり、一致していない時は1よりも小さくなり、どちらかのフーリエ係数ベクトルの大きさが0の時は0になります。 したがってこれは2つの周期関数の振幅が一致している程度を表す振幅一致係数と解釈することができます。

そこで周期関数の相関係数つまり位相一致係数をr_θで表し、振幅一致係数をr_Aで表し、エーベルの級内相関係数つまり周期関数全体の一致係数をr_pで表すと、これらの指標の間には次のような関係があります。

位相一致係数：
※測定間隔が等間隔ではない時：
振幅一致係数：
※測定間隔が等間隔ではない時：
周期関数全体の一致係数：r_p = r_θ・r_A

ちなみに周期共分散分析において、群数が2つの時の非平行性の平方和は次のようになります。

これは2つのフーリエ係数ベクトルの差ベクトルに関するマハラノビスの平方距離的な値になります。 そして測定間隔が等間隔の時は[Ｓ_xy1^-1 + Ｓ_xy2^-1]^-1が対角化し、その対角成分が(n/4)になります。 すると非平行性の平方和は2つのフーリエ係数ベクトルの差ベクトルの大きさ||β₁^* - β₂^*||を平方して(n/4)をかけた値になります。 このことから群が2つの時の非平行性の検定は本質的に2つのフーリエ係数ベクトルの差ベクトルの大きさが0かどうか、つまり2つのフーリエ係数ベクトルが一致しているかどうかの検定になり、結局のところ周期関数全体の一致係数が0かどうかの検定になることがわかると思います。

(注1)で求めたA群とB群の周期回帰曲線について、3種類の一致係数を実際に計算してみましょう。

β₁^*'β₁^* = (-21.9228)² + (-12.9126)² + (-1.34629)² + (-10.3434)² = 756.1428
β₂^*'β₂^* = (-9.85214)² + 1.56637² + (-3.88337)² + (-6.72726)² = 159.8548
β₁^*'β₂^* = (-21.9228)×(-9.85214) + (-12.9126)×1.56637 + (-1.34629)×(-3.88337) + (-10.3434)×(-6.72726) = 270.5715
　　
r_p = 0.778×0.759 = 0.591

(注3)　表12.6.1を一般化すると次のようになります。

表12.6.7 二元配置型周期共分散分析の一般的データ

群	被験者	時期					平均
		t1	…	tj	…	tp
1	1	y111	…	y1j1	…	y1p1	m1.1
	:	:	…	:	…	:	:
	n1	y11(n1)	…	y1j(n1)	…	y1p(n1)	m1.(n1)
	平均	m11.	…	m1j.	…	m1p.	m1..
:	:	:	…	:	…	:	:
i	1	yi11	…	yij1	…	yip1	mi.1
	:	:	…	:	…	:	:
	ni	yi1(n1)	…	yij(n1)	…	yip(n1)	mi.(n1)
	平均	mi1.	…	mij.	…	mip.	mi..
:	:	:	…	:	…	:	:
a	1	ya11	…	yaj1	…	yap1	ma.1
	:	:	…	:	…	:	:
	na	ya1(na)	…	yaj(na)	…	yap(na)	ma.(na)
	平均	ma1.	…	maj.	…	map.	ma..
全体	N	m.1.	…	m.j.	…	m.p.	mT

二元配置型周期共分散分析では、群ごとの時期平均値m_ij.に3種類の周期回帰モデルを当てはめます。 その式は(注1)の群別周期回帰モデル、共通周期回帰モデル、そして全体周期回帰モデルのyをmに変えたものになるので(注1)を参照してください。

ただし二元配置型周期共分散分析では全ての群の測定時期が同じため、全体周期回帰モデルと共通周期回帰モデルのメサー以外のフーリエ係数が同じになり、全体周期回帰変動と共通周期回帰変動が同一になります。 また各群の平均値の差は共通周期回帰モデルのメサーの差と同じになり、群差と修正群差が同一になります。 そのため二元配置型周期共分散分析ではデータを次のように分解し、その基本式に対する平方和と自由度を求めます。

被験者と時期の二元配置：
群と時期の二元配置：
群別周期回帰：
全体周期回帰：
○全体変動：(y_ijl - m_T)
総例数： 　　平方和： 　　自由度：φ_T = Np - 1
○被験者変動：(m_i._l - m_T)
平方和： 　　自由度：φ_SUB = N - 1
○時期変動：(m._j. - m_T)
平方和： 　　自由度：φ_P = p - 1
○被験者×時期変動：{y_ijl - (m_i._l + m._j. - m_T)}
平方和：
自由度：φ_S×P = φ_T - φ_SUB - φ_P=φ_SUB×φ_P = (N - 1)(p - 1)
○群差変動：(m_i.. - m_T)
平方和： 　　自由度：φ_A = a - 1
○被験者残差変動：(m_i._l - m_i..)
平方和： 　　自由度：φ_SR = φ_SUB - φ_A = N - a
○群−時期変動：(m_ij. - m_T)
平方和： 　　自由度：φ_AP = ap - 1
○群×時期変動：{m_ij. - (m_i.. + m._j. - m_T)}
平方和：
自由度：φ_A×P = φ_AP - φ_A - φ_P=φ_A×φ_P = (a - 1)(p - 1)
○残差=被験者残差×時期変動：{y_ijl - (m_i._l + m_ij. - m_i..)}
平方和：
自由度：φ_R = φ_SR×P = φ_T - φ_SUB - φ_P-φ_A×P = φ_SR×φ_P = (N - a)(p - 1)
○群別周期回帰変動：
平方和：
※フーリエ係数ベクトルと積和ベクトルの内容は(注1)に準ずる
自由度： (2m = pの時は(2m - 1)a)
○群別時期変動：(m_ij. - m_i..)
平方和： 　　自由度：
○群別ズレの変動：
平方和：
自由度： (2m = pの時は0)
○全体周期回帰変動＝共通周期回帰変動：
平方和：
　　 　　
　　
自由度：φ_βT = φ_βc = 2m (2m = pの時は(2m - 1))
○非平行性変動：
平方和：
自由度： (2m = pの時は(2m - 1)(a - 1))

表12.6.8 二元配置型周期共分散分析表(PERCOVA table)

要因	平方和	自由度	平均平方和(分散)	F値
群差	SA	φA	VA=SA/φA	FA=VA/VSR
個体残差	SSR	φSR	VSR=SSR/φSR
個体	SSUB	φSUB	VSUB=SSUB/φSUB	FSUB=VSUB/VSR×P
全体周期	SβT	φβT	VβT=SβT/φβT	FβT=VβT/VSR×P
非平行性	SD	φD	VD=SD/φD	FD=VD/VSR×P
ズレ合計	∑SLOFi	∑φLOFi	VLOF=∑SLOFi/∑φLOFi	FLOF=VLOF/VSR×P
残差	SSR×P	φSR×P	VSR×P=SSR×P/φSR×P
全体	ST	φT

○群差の検定
帰無仮説 H₀：m₁.. = … = m_a.. = m_T
F_A ＞ F(φ_A,φ_SR,α)の時、有意水準100α％で有意
○個体の検定
帰無仮説 H₀：m₁.₁ = … = m_a._na = m_T
F_SUB ＞ F(φ_SUB,φ_SR×P,α)の時、有意水準100α％で有意
○全体周期の検定
帰無仮説 H₀：β_T1 = … β_T(2m) = 0
F_βT ＞ F(φ_βT,φ_SR×P,α)の時、有意水準100α％で有意
○非平行性の検定
帰無仮説 H₀：β₁₁ = … = β_1a = β_c1、… 、β_(2m)1 = … = β_(2m)a = β_c(2m)
F_D ＞ F(φ_D,φ_SR×P,α)の時、有意水準100α％で有意
○ズレ合計の検定
帰無仮説 H₀：
F_LOF ＞ F(φ_LOF,φ_SR×P,α)の時、有意水準100α％で有意

表12.6.1のデータについて、基本周期を24時間とし、周期成分1と周期成分2を含んだ周期回帰モデルを用いて実際に計算してみましょう。

S_T = 2062810 - 120×129.9² ≒ 37930.8 　　φ_T = 120 - 1 = 119
S_SUB = 2025532 - 120×129.9² ≒ 652.883 　　φ_SUB = 5 - 1 = 4
S_P = 5×409340 - 120×129.9² ≒ 21821.2 　　φ_P = 24 - 1 = 23
S_S×P = 37930.8 - 652.883 - 21821.2 = 15456.72 　　φ_S×P = 4×(24 - 1) = 92
S_A = 2024953 - 120×129.9² ≒ 74.1125 　　φ_A = 2 - 1 = 1
S_SR = 652.883 - 74.1125 = 578.771 　　φ_SR = 5 - 2 = 3
S_AP = 2054392 - 120×129.9² ≒ 29512.6 　　φ_AP = 2×24 - 1 = 47
S_A×P = 29512.6 - 74.1125 - 21821.2 = 7617.287 　　φ_A×P = (2 - 1)×(24 - 1) = 23
S_R=S_SR×P = 37930.8-652.883 - 21821.2 - 7617.287 = 7839.4 　　φ_R = φ_SR×P = (5 - 2)×(24 - 1) = 69
※群別周期回帰式は(注1)と同じ
　　
　　
　　
S_βT = 5×3700.856 = 18504.28　　φ_βT = 2×2 = 4
S_D = 23902.19 - 18504.28 = 5397.91 　　φ_D = 2×2×(2 - 1) = 4

これらの統計量を用いて表12.6.3の周期共分散分析表を作成することができます。

(注4)　一元配置型周期共分散分析の一般的データは、表12.6.7において各群の被験者が測定時期ごとに独立したものになります。 そのため群ごとに時期t₁〜t_pが必ずしも同一ではなく、被験者数n_iも測定時期ごとに必ずしも同一ではなくなります。 したがって単純な周期共分散分析と同様に、全体周期回帰変動と共通周期回帰変動が同一になるとは限らず、群差と修正群差も同一になるとは限りません。

しかし一元配置型周期共分散分析でも群ごとの時期平均値m_ij.に群別周期回帰モデル、共通周期回帰モデル、そして全体周期回帰モデルを当てはめます。 そしてその式は(注1)の群別周期回帰モデル、共通周期回帰モデル、そして全体周期回帰モデルのyをmに変えたものになります。

この手法ではデータを次のように分解し、その基本式に対する平方和と自由度を求めます。

時期の一元配置：(y_ijl - m_T) = (m_ij. - m_T) + (y_ijl - m_ij.)
群別周期回帰：
全体周期回帰：
共通周期回帰：
○全体変動：(y_ijl - m_T)
総例数： 　　平方和：
自由度：φ_T = N - 1
○時期変動：(m_ij. - m_T)
平方和：
自由度：φ_P = P - 1 　　
○残差変動：(y_ijl - m_ij.)
平方和： 　　自由度：φ_R = φ_T - φ_P = N - P
○群差変動：(m_i.. - m_T)
平方和： 　　自由度：φ_A = a - 1
○群別周期回帰変動：
平方和：
　　
　　

　　
d_xjki = n_ij(x_jk-m_x._ki) 　　


自由度：
○群別ズレの変動：
平方和：
自由度：
○全体周期回帰変動：
平方和：
　　
　　自由度：φ_βT = 2m (2m = pの時は(2m - 1))
○共通周期回帰変動：
平方和：


　　自由度：φ_βc = 2m (2m = pの時は(2m - 1))
○修正群差変動：
平方和： 　　自由度：φ_AA = a - 1
○非平行性変動：
平方和： 　　自由度：

表12.6.9 一元配置型周期共分散分析表(PERCOVA table)

要因	平方和	自由度	平均平方和(分散)	F値
群差	SA	φA	VA=SA/φA	FA=VA/VR
共通周期	Sβc	φβc	Vβc=Sβc/φβc	Fβc=Vβc/VR
修正群差	SAA	φAA	VAA=SAA/φAA	FAA=VAA/VR
全体周期	SβT	φβT	VβT=SβT/φβT	FβT=VβT/VR
非平行性	SD	φD	VD=SD/φD	FD=VD/VR
ズレ合計	∑SLOFi	∑φLOFi	VLOF=∑SLOFi/∑φLOFi	FLOF=VLOF/VR
時期	SP	φP	VP=SP/φP	FP=VP/VR
残差	SR	φR	VR=SR/φR
全体	ST	φT

○群差の検定
帰無仮説 H₀：m₁.. = … = m_a.. = m_T
F_A ＞ F(φ_A,φ_R,α)の時、有意水準100α％で有意
○共通周期の検定
帰無仮説 H₀：β_c1 = … = β_c(2m) = 0
F_βc ＞ F(φ_βc,φ_R,α)の時、有意水準100α％で有意
○修正群差の検定
帰無仮説 H₀：β_c01 = … = β_c0a
F_AA ＞ F(φ_AA,φ_R,α)の時、有意水準100α％で有意
○全体周期の検定
帰無仮説 H₀：β_T1 = … β_T(2m) = 0
F_βT ＞ F(φ_βT,φ_R,α)の時、有意水準100α％で有意
○非平行性の検定
帰無仮説 H₀：β₁₁ = … = β_1a = β_c1、… 、β_(2m)1 = … = β_(2m)a = β_c(2m)
F_D ＞ F(φ_D,φ_R,α)の時、有意水準100α％で有意
○ズレ合計の検定
帰無仮説 H₀：
F_LOF ＞ F(φ_LOF,φ_R,α)の時、有意水準100α％で有意
○時期の検定
帰無仮説 H₀：m₁₁. = … = m_ap.
F_P ＞ F(φ_P,φ_R,α)の時、有意水準100α％で有意

表12.6.1のデータについて、基本周期を24時間とし、周期成分1と周期成分2を含んだ周期回帰モデルを用いて実際に計算すると残差以外は(注3)の計算結果と同じ値になります。 そして残差の平方和S_Rと自由度φ_Rは(注3)の個体残差の平方和S_SRと自由度φ_SR、残差の平方和S_SR×Pと自由度φ_SR×Pをそれぞれ合わせた値になります。 それらの統計量を用いて表12.6.4の周期共分散分析表を作成することができます。

最終更新日：2016年3月8日 　第5節へ 第7節へ

　webmaster@snap-tck.com 　Copyleft (C) 2000 SNAP(Sugimoto Norio Art Production)