fMRI解析の基礎 (５)：GLMと検定

１. ANOVA

例として、２つの要因$A,B$によるANOVAを行います。

これも設計行列と回帰係数を用いて行列表示できるため、GLMのひとつと考えられます。

ノイズ同士の相関がほぼ無く、等分散の場合には、ガウス-マルコフ定理を用いてパラメータ推定ができます。

Gauss-Markov theorem：二乗和誤差を最小化する$\beta$が最適値（最良線形不偏推定量）となります (wiki)。

２. パラメータ推定

$\beta$を一意に定めるために、$X^T X$は正則である必要があります¹。そのためにも「Jittering」によって共線性を低くします。

ノイズの分散$\sigma_{\varepsilon}^2$の推定値は、以下の式で評価できます。

ノイズには時間相関がある場合もあります。これを修正する方法には、「hrfで時間軸に沿ったスムージング」や「Prewhitening」などがあります。Prewhitening では、ノイズに自己回帰過程（AR過程）を用いることで時間相関を除外します。

また、ノイズの時間相関をGLMに取り入れてしまう例もあります（一般化線形モデル）。具体的には、ノイズの共分散行列を対角要素以外も考慮してモデル化します。この場合、推論は難しくなります。

３. 仮説検定

帰無仮説検定に使用する帰無仮説は、コントラストベクトル $c'$ を用いて表現できます。

（例）$\mathrm{H_0}:\beta_1=0$ の場合、$c'=(1,0,0,0), \ \beta = (\beta_1,\beta_2, B_0, \Delta)^T$ とする。

$X^T X$が正則である限り全ての帰無仮説が検定できます。

t検定

帰無仮説 $\mathrm{H_0:}c^T\beta=c_0$ の下で、上の t 値は自由度 $N_{TR}-r$ の t 分布に従います。TRはふつう数百あるので、 t 分布は標準正規分布に近づきます。全てのボクセルに対して t 値を算出した結果を、Statistical Parametric Map と呼びます。

以下のような注意点があります。

• 第一の誤り（偽陽性）を減らすようにする（多重比較問題）。
• 設計行列の共線性を低くする（$c^T(X^T X)^{-1}c$が小さくなるように設計する）。

多重検定：one way ANOVA

３つのイベントに対して、帰無仮説は $\mathrm{H_0:}Q\beta=q_0$ と表現できます。

• 例）$Q\beta$ を次のように設計します。
  $$\begin{alignedat}{2} Q\beta &= \begin{bmatrix} 1 & -1 & 0 & 0 & 0 \\ 1 & 0 & -1 & 0 & 0 \\ 0 & 1 & -1 & 0 & 0 \\ \end{bmatrix} \begin{bmatrix} \beta_1 \\ \beta_2 \\ \beta_3 \\ B_0 \\ \Delta \end{bmatrix} \\ &= \begin{bmatrix} \beta_1 - \beta_2 \\ \beta_1 - \beta_3 \\ \beta_2 - \beta_3 \end{bmatrix} \end{alignedat}$$
  帰無仮説が、$\mathrm{H_0:}Q\beta=q_0=\mathbb{O}$ のとき、対立仮説 $\mathrm{H_1}$は「$\beta_1 \neq \beta_2$ または $\beta_1 \neq \beta_3$ または $\beta_2 \neq \beta_3$」となります。

帰無仮説検定の指標となる F 値は以下のように算出します。帰無仮説が正しい場合、F 値は m 次の F 分布に従います。

FBR model の場合：カイ二乗検定

１つのイベントに対して複数の$\beta$があるため、帰無仮説は $\mathrm{H_0:}\sum_{i=1}^M\beta_i^2=0$ となります。

この場合、以下のカイ二乗値を評価します。帰無仮設が正しい場合、カイ二乗値は M 次のカイ二乗分布に従います。

４. 多変量GLM

単変量解析では、それぞれのボクセルの独立性を仮定していました。ここでは、２つ以上のボクセルを同時に解析することで空間的な相関関係を考慮します。

今までのGLMを、D 個のボクセルに拡張したときの行列表示は以下のようになります。

多変量GLMでは、時間的な独立性のみを前提とします。誤差項$E$は全てのTRで同じものとします。

この最適解は単変量GLMと同様に正規方程式で与えられます。

誤差項は$E=Y-X\hat{B}$となります。誤差の分散共分散行列$\hat{\Sigma}_s$を推定しましょう。

５. ノンパラメトリック法

以上のGLMは、生のBOLD反応が正規分布に従い、神経活動と時系列相関があるという仮定の下で成立します。

分布の仮定を弱めるために、ノンパラメトリック手法であるブートストラップ法が使われることもあります。

1. 元データの統計量を算出する。
2. TRをシャッフルして、統計量を算出するのをM回繰り返す。
3. 求めた全統計量内における元データの統計量以上の割合を、帰無仮説検定の p 値とする。

７. Percent Signal Change

BOLD値は、以降の解析のために Percent Signal Change という統計量に変換することがあります。

Taskの効果量のみをパーセント指標で表現しましょう。Task と Rest がある場合、以下のように求めます。

Rapid Event-Related Desgin では Rest が定義しにくいので、実験全体の平均BOLDを$\bar{y}_R$とすることがあります。しかし、この方法は過大評価する傾向があります。

GLMを用いてNuisance 変数の回帰係数を計算し、それを用いて Rest の平均BOLDを推定することもできます。

この方法は Task-related BOLDの最大値を過小評価してしまうため、ボックスカーを元にした以下の推定値が利用されているようです。

---

Reference

• Ashby, F. G. (2019). Statistical analysis of fMRI data. MIT press. url