全文 - terrapub

日本統計学会誌
第 44 巻, 第 1 号, 2014 年 9 月
217 頁 ∼ 240 頁
日本統計学会賞受賞者特別寄稿論文
多変量解析の感度分析：影響関数・Cook の局所影響分析・
固有値問題の摂動論
田中豊∗
Sensitivity Analysis in Multivariate Methods: Inﬂuence Function, Cook’s
Local Inﬂuence and Perturbation Theory of Eigenvalue Problems
Yutaka Tanaka∗
1980 年代の初め頃から統計モデルとくに多変量解析における感度分析に関心を持ち，数量化
法 (対応分析)，主成分分析，因子分析，共分散構造分析，正準相関分析など各種の感度分析法の
研究を行ってきた．本稿ではその間に出会い使用してきた主要な方法論である影響関数，法曲率
を用いた Cook の局所影響分析，固有値問題の摂動論について概説し，関連する話題について議
論する．
We have studied sensitivity analysis in statistical models, in particular, multivariate statistical models such as quantiﬁcation methods (correspondence analysis), principal component
analysis, canonical correlation analysis, factor analysis and covariance structure analysis, since
early 1980s. This article surveys major methodologies of sensitivity analysis, which we used
during our studies, including inﬂuence functions, Cook’s local inﬂuence analysis using normal
curvature and perturbation theory of eigenvalue problems, and discuss some topics related to
these methodologies.
キーワード: 影響関数，局所影響分析，影響グラフ，法曲率，固有値問題の摂動論
はじめに
1.
多変量解析における感度分析・影響分析法の研究にはじめて取り組んだのは 1980 年代の
初め頃のことであった．1970 年代に回帰分析における感度分析に関する研究が急速に進み，
1980 年代に入って Belsley et al. (1980)，Cook and Weisberg (1982)，Atkinson (1985)，
Chatterjee and Hadi (1988) などの専門書が次々に出版された．これらの回帰診断分野に
おける成功に刺激されて，その考え方や技法が，当時関心をもっていた数量化法の分析
結果の安定性や信頼性の評価に応用できないかと考えたことに始まる．その後，数量化法
∗
岡山大学名誉教授：〒 700-8530 岡山市北区津島中 3-1-1 (E-mail: [email protected]).
218
日本統計学会誌第44巻第1号 2014
(Tanaka (1984), Tanaka and Tarumi (1986))，主成分分析 (Tanaka (1988), Yamanishi and
Tanaka (2005))，因子分析 (Tanaka and Odaka (1989a, b))，共分散構造分析 (Tanaka et
al. (1991), Tanaka and Watadani (1992))，正準相関分析 (Castano-Tostado and Tanaka
(1989), Tanaka et al. (2002)) などにおける感度分析法について研究し，ソフトウェアを開
発した (Tarumi and Tanaka (1986), 森他 (1997), Mori et al. (1998))．また初期の頃には
田中 (1984, 1992)，Tanaka (1994) など，逐次，総合報告・解説のような形で報告したが，
その後，そういう形で発表する機会がなく現在に至っている．その間，多くの研究者がこ
の分野の研究に取り組み，研究対象も単純なモデルからより複雑なモデルへと拡大し，研
究のための方法論も充実してきている．本稿では，開発された個別統計モデルの感度分析
法ではなく，これまでの感度分析の研究の中で利用してきた主要な方法論である Hampel
や Belsley らの影響関数，Cook の幾何学的な概念に基づく局所影響分析，固有値問題の摂
動論に関連した話題について概観し議論する．
感度分析の方法としては，1 個ずつあるいは複数個の観測値を落として分析した結果と
元のすべての観測値を用いた分析結果を比較して特定の観測値の影響を測る global な感度
分析と，観測値の値や重みを微小に変化させ微分を用いて影響を測る local な感度分析が
ある．回帰分析など観測値を落として分析した結果が比較的容易に求まる場合には global
な感度分析が可能で，実際，回帰診断において global な診断統計量が主要な位置を占める
が，多変量解析などの複雑なモデルでは global な感度分析は計算的に容易でなく，微分を
用いた局所的な感度分析が中心となっている．
以下，2 節で影響関数と Cook の局所影響分析およびそれら相互の関係，3 節で記述的な
多変量解析において重要な役割を担う固有値問題の摂動論に関連した話題を取り上げる．
影響関数と Cook の局所影響分析
2.
2.1
Hampel の影響関数：分布関数への摂動
感度分析の主要なツールのひとつに影響曲線 (inﬂuence curve) あるいは影響関数 (inﬂuence
function, IF) がある (Hampel (1974))．p 変量観測ベクトルの分布関数を F ，分析結果を
表す分布関数の汎関数 (パラメータと呼んでおく) を θ = θ(F ) とし，分布関数 F に対して
微小な変化 (摂動)F → F˜ = (1 − ε)F + εδ(x) を加えてそれに対応する θ(F ) の変化を調べ
る．δ(x) は値 x を確率 1 で取る確率変数の分布関数である．このとき
IF (x, θ) = lim {θ(F˜ ) − θ(F )}/ε = θ(1) (x)
ε→+0
(2.1)
で定義される関数を (理論) 影響関数と呼ぶ．上付きの (1) は 1 次の微係数を表す．標本の
データに基づく影響関数として次の 3 種類が感度分析によく利用される (Cook and Weisberg
(1982), Critchley (1985))．1) 経験影響関数 (empirical inﬂuence function, EIF)：理論影響
多変量解析の感度分析
219
関数の定義式 (2.1) の中の分布関数 F を経験分布関数 Fˆn で置き換えた式で定義され，各
ˆ ≡ θˆ(1) がそれぞれの観測値
観測値 (x = xi , i = 1, 2, . . . , n) に対して計算した値 EIF (xi , θ)
i
の影響を測る量として使用される．2) 標本影響関数 (sample inﬂuence function, SIF)：1
つの観測値を落とすことに相当する摂動 (EIF の定義式において ε = −1/(n − 1) とおき極
ˆ = −(n − 1)−1 (θˆ(i) − θ)
ˆ により計算さ
限をとらない) に対する影響を測るもので SIF (xi , θ)
れる．ここに，下つきの (i) は i 番目の観測値を除いた標本に基づく推定量であることを表
す．3) deleted EIF：i 番目の観測値を除いた n − 1 個の観測値に基づく経験分布関数 Fˆn−1
を基準 (摂動前の分布関数) とし，i 番目の観測値の値を確率 1 でとる分布関数を微小な割
合 ε で加えたときの分布関数 F˜ = (1 − ε)Fˆn−1 + εδ(xi ) に対する汎関数 θ(F˜ ) の ε = 0 にお
ける微係数として定義される．i 番目の観測値が唯一の外れ値であるとき摂動前の分布関
数が正しく推定され，その観測値の影響が正確に測定される．
平均ベクトルおよび共分散行列の理論影響関数は µ(1) = x−µ, Σ(1) = (x−µ)(x−µ)T −Σ，
対応する経験影響関数は µ
ˆ(1) = x − µ, S (1) = (x − x
¯)(x − x
¯)T − S となる．ここに x
¯, S は
標本平均および標本共分散行列 (偏差平方和積和行列を n で割ったもの) である．一般に影
響関数 θ(1) (x) の期待値はゼロになる．
2.2
Belsley らの観測値への重みに関する微係数：重み摂動・分散摂動
観測値の影響を測る別の考え方として，特定の観測値への重み (case weight) を考え，そ
の重みに関する微係数を利用する方法がある (Belsley et al. (1980), Pregibon (1981) など)．
Belsley らは回帰分析において i 番目の観測値に対する重み wi を１からわずかに変化させ
ることと誤差分散に対して V (εi ) = σ 2 → V (εi ) = wi−1 σ 2 のような摂動を与えることを同
[
]
ˆ i )/∂wi
等とみなして，観測値の影響を微係数 ∂ θ(w
で測り，この微係数を影響関数と
w0
呼んだ．wi = ∞ は i 番目の観測値を落とすことに対応し，w0 = (1, 1, . . . , 1)T は非摂動状
態を表す．この摂動は多変量解析の場合
V (xi ) = wi−1 Σ
(2.2)
の形になる．
また，Pregibon (1981) はロジスティック回帰において共変量の値に関わらず個別の観測
値ごとに並べた unstructured data set と，同じ共変量をもつ観測値を集めてその個数を記
述する structured data set を対比して考え，unstructured の場合の各々の観測値に対し観
測値の個数に相当する重み wi を導入している．
通常の分析では分散摂動と observation count の意味での重み摂動の効果は同じになるの
で，ここでは特に断らない限り両方を含めて重み摂動と呼んでおくが，モデルが未知の変
数変換パラメータを含む場合には感度分析の結果が異なるので注意が必要である (Carroll
and Ruppert (1988), pp. 177–178; Tsai and Wu (1992))．
220
日本統計学会誌第44巻第1号 2014
分散摂動における重みの与え方としては，式 (2.2) の形が標準的であるが，Hampel の影
響関数との対応を考えると，全体の重みの和を一定に保つように修正した形

−1
/∑
V (xi ) = nwi
wj  Σ
(2.3)
j
で利用するのが望ましい (Tanaka (1994), Kwan and Fung (1998), Tanaka and Zhang
(1999))．式 (2.2) の形をタイプ 1 の重み，式 (2.3) の形をタイプ 2 の重みと呼ぶことにす
る．これら 2 つのタイプの重みを用いたときの平均ベクトルと共分散行列に対する最尤推
定量の重み wi に関する偏微分係数はそれぞれ
タイプ 1：∂ µ
ˆw /∂wi = n−1 (xi − x
¯),
タイプ２：∂ µ
ˆw /∂wi = n−1 (xi − x
¯),
∂Sw /∂wi = n−1 (xi − x
¯)(xi − x
¯)T
∂Sw /∂wi = n−1 {(xi − x
¯)(xi − x
¯)T − S}
となり，タイプ 2 の場合の微係数はタイプ 1 の場合の微係数からその平均を引いた形，す
なわち，当該観測値の影響から平均的な影響を引いた形になる．2.1 節の Hampel の影響関
数と比較すると，タイプ 2 の重みを用いたときの微係数 (Belsley らの影響関数) は Hampel
の影響関数の 1/n に等しく，観測値の影響の相対的な大きさは同じになる．平均ベクトル
と共分散行列だけでなく一般のパラメータの最尤推定量に関しても，Hampel の経験影響
関数 EIF とタイプ２の重みの微係数の間に同じ関係が成り立つ．
回帰分析においては回帰係数への影響が主要な関心の対象となるが，i 番目の観測値に対
する残差を ei と表すとき，回帰係数の推定量の重みに関する微係数はタイプ 1 の重みの場
合 βˆ(1) = (X T X)−1 xi ei でその平均がゼロになるため，タイプ 2 の重みに関する微係数と
一致し，両者とも Hampel の EIFβˆ(1) = n(X T X)−1 xi ei の 1/n となる．従って，回帰係数
に対する影響を議論するときにはタイプ 1，タイプ 2 のどちらの重みを用いてもよい．
2.3
尤度距離の法曲率を用いる Cook の局所影響分析
Cook (1986) は観測値の影響に限らず，さまざまな形でモデルやデータに導入される摂
動の影響を評価する，幾何学に基づく汎用的な方法を提案した．すなわち，任意に導入さ
れる摂動パラメータ・ベクトル (観測値の影響を評価する場合には観測値全体に重みベク
トルを付与) を一定方向に微小に変化させてその影響を評価し，局所的な影響が最大とな
る方向を探索する方法である．
n 個の p 変量観測ベクトル {x1 , . . . , xn } が得られ，これに対して p 次元パラメータ・ベ
クトル θ = θ(F ) を含む統計モデルが想定されたとする．θˆ を θ の最尤推定量とし，摂動を
表すため r 次元の摂動パラメータ・ベクトル w ∈ Ω (Ω : Rr の開部分集合) を導入，非摂動
状態 w0 (w0 ∈ Ω) から摂動状態 w ∈ Ω へパラメータ・ベクトルのすべての要素を同時に変
化させる．一般性をもたせるため，次元を r としているが，観測値の重みや値の影響を評
価する場合には r = n となる．本稿では観測値の影響評価を想定して w の次元は n として
多変量解析の感度分析
221
ˆ w0 )
おく．L(θ|w0 ), L(θ|w) をそれぞれ摂動前，摂動後の対数尤度関数とし，L(θ|w) は (θ,
の近傍で θ, w に関して 2 回連続微分可能と仮定，θˆw を摂動後の最尤推定量とする．パラ
メータに等式制約があってもよいが，制約式は摂動パラメータ w に依存せず，θ に関して
2 回微分可能とする．簡単のため，全パラメータに関心がある場合を考える．
摂動 w0 → w に対応してパラメータ推定量が θˆ → θˆw と変化するとき，その影響の大き
さを次の式で定義される尤度距離 (likelihood displacement) D(w) によって測る．
ˆ 0 ) − L(θˆw |w0 )]
D(w) = 2[L(θ|w
(2.4)
尤度距離を用いる理由は θ の漸近信頼領域が
ˆ 0 ) − L(θ|w0 )] ≤ χ2 (m)}
{θ| 2[L(θ|w
α
の形で与えられることによる．重みベクトル w の尤度距離 D(w) への影響の大きさを評価
するため影響グラフ α(w) = (wT , D(w))T の形状を調べる．D(w) を w = w0 の近傍でテ
イラー展開すると
D(w) ∼
= D(w0 ) + (w − w0 )T
∂ 2 D(w0 )
∂D(w0 ) 1
+ (w − w0 )T
(w − w0 )
∂w
2
∂w∂wT
と近似でき，最尤推定量の性質から 1 次の項は消える．任意方向の直線 w0 → w = w0 +
td, kdk = 1 に沿った摂動を導入して，この直線と w0 における曲面 α(w) の法線によって
決まる平面で切ったときの切り口の曲線 (lifted line) に注目し w = w0 における曲率 (曲面
の w0 における法曲率) Cd が大きいほど影響が大きいと考えて，法曲率 Cd を最大にする
方向 dmax を最大影響方向 (most inﬂuential direction)，ベクトル dmax の要素のうち絶対値
の大きい要素に対応する観測値の集合を影響の大きい部分集合 (inﬂuential subset) とみな
す．最大法曲率と最大影響方向は，次に示す固有値問題 (2.5) の解として得られる．
)
(
T
∂ θˆw
∂ θˆw
[−Lθθ ]
− τˆI d = 0
(2.5)
2
∂w
∂wT
ただし，−Lθθ = −∂ 2 L(θ)/∂θ∂θ T は観測情報行列を表し，式 (2.5) のカッコ内第 1 項は尤度
距離 D(w) の w に関する 2 次の偏微分係数行列 (ヘシアン行列) である．固有値問題 (2.5)
の最大固有値 τˆmax が最大法曲率 Cmax ，対応する固有ベクトル dmax が最大影響方向を与
える．
全パラメータ θ = (θ1T , θ2T )T でなくパラメータの一部 θ1 のみに関心がある場合には尤度
距離の代わりにプロファイル尤度距離を用いて同様に定式化される．
なお，方向を固定して d = Ei の方向 (Ei は i 番目の要素のみ 1，ほかの要素はゼロであ
T
るような n 次元ベクトル) の法曲率を考えると 2(∂ θˆw
/∂wi )Vˆ −1 (∂ θˆw /∂wi ) となり，i 番目
の観測値に対する影響関数を Vˆ −1 で規準化した 2 乗ノルムの 2 倍に等しい．
222
日本統計学会誌第44巻第1号 2014
影響グラフの目的関数が尤度距離以外の一般の関数の場合の f の w = w0 における法曲
率は以下のように表される (Cook (1986))．
Cd =
dT Hf d
1
(1 + ∇Tf ∇f ) 2 dT (In + ∇f ∇Tf )d
ここに ∇f = ∂f /∂w，Hf = ∂ 2 f /∂w∂wT は勾配ベクトルとヘシアン行列 (いずれも w = w0
における値) を表す．勾配がゼロでない場合，目的関数を k 倍すると法曲率は
Cd =
kdT Hf d
1
(1 + k 2 ∇Tf ∇f ) 2 dT (In + k 2 ∇f ∇Tf )d
となり，Cmax と dmax のどちらも変化し尺度不変でない (Fung and Kwan (1997))．
2.4
影響関数の主成分分析 (PCA)
関心のあるパラメータを θ，その次元を p，関心のある m 個の観測値の集合を S とする．
S に属する m 個の観測値に関連する方向の一般化影響関数を，影響関数の定義式 (2.1) で
∑
F˜ = (1 − ε)F + εG, G = m−1 j∈S δ(xj ) とおいた式で定義すると，観測値の集合 S を除
く標本に基づく θ の推定量は，ε の 1 次の項までのテイラー展開を用いて
θˆ(S) ∼
= θˆ − (n − m)−1
∑
(1)
θˆi
(2.6)
xi ∈S
により近似できる (Tanaka (1994))．すなわち，影響関数のように 1 次微分を用いて影響
を評価するとき，複数の観測値の影響は個々の観測値の影響の和によって決まる．Enguix∑
Gonzalez et al. (2005) は F˜ = (1 − mε)F + εG, G = j∈S δ(xj ) とおいて複数個の観測
値に対する一般化影響関数を定義しているが同じ考え方といえる．このことから影響の大
きい観測値の集合を見つけるためには “単独で比較的に影響が大きく，影響の方向が互い
に類似している観測値の集合”を探せばよい．ただし，方向の類似性を判断するときには
θˆ の要素の分散・共分散を考慮する必要がある．このことに注目してわれわれは影響関数
(1)
ˆ の逆行列を計量行列とする主成
{θˆi , i = 1, . . . , n} に対して漸近共分散行列 Vˆ = a cov(
ˆ θ)
分分析 (PCA) を行う方法を提案した (Tanaka et al. (1990), Tanaka (1994), Tanaka and
Zhang (1999))．基本的な考え方は次の通りである．
(1)
(1)
(1)
ˆ (1) = (θˆ , θˆ , . . . , θˆn ) の n 個の列を p 次元空間
標本に基づく影響関数の p × n 行列 Θ
1
2
(1)
内の n 個の点と考え，線形変換 zi = AT θˆi
, (A : p × q, q ≤ p) を用いて q 次元に縮約す
∑
(1)T
(1)
る．元の p 次元空間における点の散らばりは i Di , Di = θˆi Vˆ −1 θˆi により評価でき，
∑
(1)T
(1)
線形変換後の q 次元空間における散らばりは i Di∗ , D∗i = θˆi A(AT Vˆ A)−1 AT θˆi と表
される．散らばりの大きさをできるだけ保存しながら次元を縮約するためには，変換後の
散らばりを最大化すればよい．この最大化問題の解は一般固有値問題
ˆ Vˆ )a = 0
ˆ (1) Θ
ˆ (1) − λ
(n−1 Θ
T
(2.7)
多変量解析の感度分析
223
(1)
の解として得られる．影響関数 θˆi の平均はゼロであるから，式 (2.7) の左辺カッコ内第 1
項は影響関数の標本共分散行列に等しい．aj をこの一般固有値問題の j 番目に大きい固有
(1)
ˆ j に対応する固有ベクトル，zji = aT θˆ を第 j 主成分スコア (PC スコア) とすると
値λ
j i
(1)T
(1)
2
Di = θˆi Vˆ −1 θˆi = z1i
+ · · · + z 2pi
(1)
の関係が成り立ち，i 番目の観測値の影響 Di (影響関数ベクトル θˆi を θˆ の漸近共分散行
(1)
列の推定量で規準化した 2 乗ノルム，一般化 Cook 距離 Di と呼んでおく) は影響関数 θˆi
の PC スコアの 2 乗の和に分解される．この意味で影響関数の PCA は一般化 Cook 距離
Di として要約される影響の多次元的な構造に関する情報を与え，影響関数によってとらえ
られた観測値の影響の主要な次元を明らかにしてくれる．Vˆ が退化するときには Vˆ の逆行
列を一般逆行列で置き換えればよい (Tanaka and Zhang (1999))．
影響が大きいとみなしてある観測値の集合を落とすとき実際に影響が大きいかどうか
という意味での近似の良さは，経験影響関数 EIF と標本影響関数 SIF の近さに関係す
る．Critchley (1985) が述べているように，PCA (あるいは多くのほかの多変量解析法)
においては n → ∞ とするとき EIF と SIF は一致するが，回帰係数の影響関数の場合は
SIF i /EIF i → (1 − hii )−1 となって EIF と SIF は一致しない．とくに，てこ比 hii のばら
つきが大きいとき食い違いが大きくなり，式 (2.6) の近似の精度はよくない．
2.5
標準化影響関数の主成分分析
2.4 節の影響関数の PCA と類似の方法に Lu et al. (1997) の標準化影響関数の PCA があ
る．理論影響関数 IF (x; θ) の共分散行列を V ∗ ，標本空間内の n 個の点 {x1 , . . . , xn } における
影響関数の p×n 行列を IF (X; θ)，標準化した影響関数行列を SIF (X; θ) = V ∗−1/2 IF (X; θ)
とし，それに基づく標準化影響行列とその compliment を以下のように定義してそれぞれ
をスペクトル分解する．
SIM (X; θ) = n−1 SIF (X; θ)T SIF (X; θ) = LΛ∗ LT
SIM c (X; θ) = n−1 SIF (X; θ)SIF (X; θ)T = EΛ∗ E T
SIM は n 次，SIM c は p 次の正方行列である．この 2 つの行列の共通の固有値を Λ∗ =
diag(λ∗1 , . . . , λ∗r )，ノルム 1 に規準化された固有ベクトルを列ベクトルとする n × r，p × r 行
√
列をそれぞれ L = (l1 , ..., lr ), E = (e1 , ..., er ) とすると， λ∗i li は標準化影響関数の PCA
√
√
から得られる第 iPC スコアを表し， λ∗i li と λ∗i0 li0 の n 個の要素の散布図を描くことに
よって影響の多次元的な構造がわかる，というものである．ここに r は SIF (X; θ) の階数
である．SIM c (X; θ) は標準化影響関数の共分散行列，そのスペクトル分解は PCA に相当
するので，固有ベクトルは PCA の係数ベクトル，SIM c の双対形である SIM の固有ベ
224
日本統計学会誌第44巻第1号 2014
クトルは定数倍を除いて PC スコアとなり，2.4 節の PC スコアで影響の多次元構造を探る
方法と基本的に同じ方法ということができる．2.4 節の定式化では経験影響関数と推定量 θˆ
の漸近共分散行列の推定量 Vˆ ，2.5 節の定式化では理論影響関数とその漸近共分散行列 V ∗
ˆ に対する影響関数に基づく推定量である
ˆ (1) Θ
ˆ (1)T ∼
を用いているが，n−2 Θ
= n−1 Vˆ ∗ は V (θ)
ˆ∼
ˆ ∗ の関係が成り立つ．
ことを考慮すると，両者の固有値の間に λ
= nλ
2.6
Cook の局所影響分析と影響関数に基づく影響分析の関係
Cook の局所影響分析が全観測値に対して付与する重みベクトル (摂動パラメータ・ベク
トル) の要素全体を同時に変化させてその影響を評価するのに対して，影響関数に基づく方
法は個々の観測値の重みを 1 つずつ変化させて影響を評価するというように，感度分析へ
の接近の仕方は異なるが，分析結果からは同じ情報が得られる (Tanaka (1994), Shi (1997),
Tanaka and Zhang (1999) など)．
ˆ
タイプ 2 の重みを摂動パラメータ・ベクトル，影響関数を重みに関する微係数 ∂ θ/∂w
i
により定義し，簡単のため全パラメータに関心のある場合について議論する．Cook の法曲
率に基づく方法からは式 (2.5) の固有値問題が得られたが，観測情報行列の逆行列を用い
∂ θˆ
て θˆ の漸近共分散行列を推定し，式 (2.5) より d ∈ Span( ∂ww ) が成り立つため d =
T
∂ θˆT
w 0
∂w a
を満たす a0 が存在することより，式 (2.5) の固有値問題は
(
∂ θˆw ∂ θˆT
2 T w − τˆVˆ
∂w ∂w
)
a0 = 0
と変形できる．この式を影響関数の PCA の固有値問題の式 (2.7) と比較すると，左辺カッ
コ内の第 1 項の行列にかかる乗数以外は同じ形であるから，固有ベクトルは共通，一方の
固有値は他方の固有値の定数倍になり，2 つの固有値問題から同じ情報が得られる．タイ
プ 2 の重みに関する微係数で定義される影響関数は Hampel の影響関数の定数倍に等しい
ので，Hampel の影響関数に基づく分析と Cook の局所影響分析からも同値の情報が得ら
れることになる．一部のパラメータに関心のある場合やパラメータに制約のある場合にお
いても同様な関係が成り立つ (Tanaka and Zhang (1999), Tanaka et al. (2003))．
固有値の関係について詳しく述べると，Cook の法曲率 τˆ とタイプ 2 の重みの微係数で定
ˆ と間の関係は τˆ = 2nλ
ˆ ，Hampel の影響関数の PCA
義された影響関数の PCA の固有値 λ
ˆ となる．2.5 節の標準化影響関数 PCA の固有値 λ∗ は Hampel の影
の場合なら τˆ = 2λ/n
響関数に対して θˆ の共分散行列でなく影響関数 θ(1) の共分散行列に対する推定量を計量行
ˆ ∗ となる．
列として用いているので τˆ = 2λ
Shi (1997) は主成分分析における局所影響分析において w0 → w = w0 + td, kdk = 1 の
形の d 方向の摂動を影響関数の定義式に導入して t に関する微係数 (d 方向の方向微分) を
一般化影響関数 (generalized inﬂuence function) と呼び，Cook and Weisberg (1982) で用
多変量解析の感度分析
225
いられている一般的な計量行列 M/c (M は正定値または非負定値行列，c はスカラー) を
用いたノルムを一般化 Cook 統計量 (GC) と定義し，M/c として観測情報行列を用いると，
GC を最大化する dmax が Cook の最大影響方向，GCmax の 2 倍が Cook の最大法曲率 τˆmax
と等しくなる．
2.7
感度分析のロバスト化
簡単のためパラメータ推定量が標本共分散行列だけに基づいて定まる場合を考えると，
影響関数 (タイプ 2 の重みの微係数) に関して
ˆ
ˆ T
∂ θ/∂w
i = (∂ θ/∂s )(∂s/∂wi )
の関係が成り立つ．S は標本共分散行列，s = vech(S) であり，2.2 節より
∂s/∂wi = n−1 vech{(xi − x
¯)(xi − x
¯)T − S}
となる．もし，x
¯, S が母平均，母共分散行列を正しく反映していなければ，パラメータ
推定量の影響関数も影響を正しく反映せず，いわゆるマスク効果を受けることになる．そ
こで影響関数を計算する元になっている平均ベクトルと共分散行列などに対するロバス
ト推定を用いて外れ値や影響の大きい観測値の影響を受けにくい影響関数の式を計算し,
それに基づいて影響分析を行う方法がいろいろ提案されている．たとえば，Tanaka and
Watadani (1994) は共分散行列に対する minimum volume ellipsoid estimator を用い外れ
値を含まないデータを選択して影響関数をロバスト推定する形の因子分析のロバスト影響
分析法を提案し，Lu et al. (1997) は回帰分析において least median of squares (LMS) や
least trimmed squares (LTS) を用いて SIM のロバスト推定を求めた上で影響分析する方法
を提唱している．より丁寧な接近法に Atkinson らの前進探索法 (forward search method)
がある (Atkinson and Riani (2000), Atkinson et al. (2004))．それは，モデルにもっとも
よく当てはまる最小サイズ m0 の集合をまず探索し，現在のサイズ m の集合に対して残り
の n − m 個の中から 1 個ずつ追加してみて，あてはまりのよい m + 1 番目の観測値の集合
に拡大，同じ方法で 1 個ずつ集合を拡大しながら，関心のある統計量 (パラメータ推定値や
診断指標) の値をモニターしていく方法である．
2.8
Cook の局所影響分析に関連したいくつかの話題
1) 摂動パラメータ w の変数変換
最大法曲率や最大影響方向に対する摂動パラメータの変数変換の影響については，摂動
パラメータになめらかな 1 対 1 変換 wi → k(wi ) を施すとき，最大影響方向 dmax は不変，
最大法曲率 Cmax は (∂k(w0 )/∂wi )2 倍になる (Cook (1986), Rejoinder)．定数の臨界値を用
いるためには影響指標を一定の範囲に入るように規準化しておく必要があるという観点か
226
日本統計学会誌第44巻第1号 2014
ら，Poon and Poon (1999) は共形変換に対して不変で，[0, 1] の範囲の値をとるように規
準化した共形法曲率 (conformal normal curvature)
Bd = C d /{tr(Hf )2 }1/2
を提案した．共形変換とはヤコビアン M が正則で M M T = cIn (c は正の定数) を満たす
変換で，たとえば，wi → (1 + wi )/2 やなめらかな関数による変換 wi → k(wi ) などは非摂
動点において共形変換である．
2) 影響グラフの目的関数
Cook (1986) は影響グラフの任意の点における法曲率の式を提示したが，実際の応用で
は目的関数としてもっぱら尤度距離を用いられ，勾配がゼロでない非摂動点における法曲
率を用いる研究は少なかった．
Laurance (1988) は Box-Cox 変換 (変換パラメータ未知) された従属変数をもつ回帰モデ
ルの誤差分散に w0 → w = w0 + td, kdk = 1 の形の分散摂動を導入して，変換パラメー
[
]
ˆ を目的関数とする影響グラフを考え，勾配 ∂ λ(t)/∂t
ˆ
タの最尤推定量 λ(t)
を最大にす
t=0
る方向ベクトルに基づく分析法を提案した．Tsai and Wu (1992) は，Laurance (1988) の
Box-Cox 変換パラメータの最大勾配方向は，そのパラメータに対するプロファイル尤度距
離を目的関数とする Cook (1986) の最大法曲率方向と同じになること，また，分散摂動で
なく，観測値の数の意味での重み摂動を用いた場合も，同様なことがいえることを示した．
Wu and Luo (1993) は回帰モデルにおいて関心のあるスカラー・パラメータの最尤推定量
を目的関数とする場合の曲面を MLE 曲面と呼び，その非摂動点における最大勾配方向に
基づく分析 (プロファイル尤度距離を目的関数とする Cook (1986) の法曲率に基づく分析
と同値) を ﬁrst order approach と呼び，それに加えて global な効果をより精度よく推定す
るため MLE 曲面の 2 次の微係数に基づく分析 (second order approach) を行うことを提案
している．
3) 影響指標の臨界値
外れ値の検出のためにはしばしば検定が使用されるが，影響診断のための各種の影響
指標には検定はなじまない (Cook (1986), p. 168; Cook and Weisberg (1999), p. 357)．
Cook (1986) は最大法曲率に対して a useful general reference という程度の意味で臨界値
Cmax = 2 を提案し，その理由として，すべての観測値の影響が均一である場合に２となる
こと，具体的には，説明変数が定数項だけの回帰分析に分散摂動を導入し回帰係数ベクト
ルに関心がある場合の最大法曲率 Cmax が 2 になることをあげている．誤差分散を既知と
すると尤度距離の法曲率は
Cd = 2dT D(e)PX D(e)d/σ 2
多変量解析の感度分析
227
と表される．ただし，PX = X(X T X)−1 X T は説明変数 X の列ベクトルの張る部分空間
への射影行列，D(e) = diag(e1 , . . . , en ) は残差を対角要素に持つ対角行列，d はノルム 1
の方向ベクトルであり，誤差分散に対しては推定量を代入して利用する．説明変数が定
数項だけの場合を考えると，PX = n−1 1n 1Tn , D(e)PX D(e) = n−1 eeT (ただし，1Tn =
(1, . . . , 1), eT = (e1 , . . . , em ) は残差ベクトル) となり，これを最大化する固有値問題はゼロ
でない固有値 1 個と (n − 1) 個の固有値ゼロを持ち，ゼロでない固有値に対応する固有ベク
トルは dmax = e/ kek，対応する固有値は τmax = Cmax = 2 となるというものである．説明
変数が定数項のみの場合だけでなく，説明変数が複数個ある場合も，例えば，2n 完備型実
験計画 (L4 直交表の 3 つの列に 2 水準の因子を割り付け) の各水準組合せで r 回ずつの繰
り返しを行う線形モデルを考えると，上と同様な方法で τˆmax ∼
= 2 となることが導かれる．
Cmax = 2 は釣り合いのとれた計画では説明変数の個数によらず成り立つと予想される．
外れ値や異常に大きい影響を持つ観測値が存在せず影響関数が均一である場合の結果と
比較するという考え方は，2.4 節あるいは 2.5 節の影響関数・標準化影響関数の PCA に対
しても応用できる．2.4 節の影響関数の PCA における一般固有値問題 (2.7) の左辺カッコ
ˆ の推定量となり，n が大きくなるとカッコ内
内の第 1 項の 1/n は影響関数に基づく V (θ)
ˆ に近づき，固有値は λ
ˆ∼
第 2 項 Vˆ と同様に V (θ)
= n，換算すると τˆ ∼
= 2 になる．標準化影
響関数の PCA の固有値問題においても，同様に，n が大きくなると共分散行列 SIM c は
←∗
単位行列に近づき，その固有値は λ ∼
= 1 ，法曲率に換算すると τˆ ∼
= 2 に近づく．
Lu et al. (1997) は，さらに，固有値の漸近分布を用いた次のような臨界値を提唱してい
る．彼らは，Anderson (1984) にあるように多変量正規データの標本共分散行列の固有値
ˆ ∗ の漸近分布が N (λ∗ , 2λ∗2 /n) となるが，非正規の場合も同じ分布に近似的に従うものと
λ
√
√
仮定し，さらに分散安定化変換を用いて，SIM c の固有値 λ∗ の臨界値 1 + 2zα / n + zα2 /n
を導出している．ここに zα は標準正規分布の上側 α 点を表す．法曲率 τˆ の臨界値はこの
値の 2 倍になる．α = 0.50 とすると上の展開式の第 2 項以下が消えて，固有値 λ∗ の臨界
値は 1，曲率の臨界値は 2 となり，Cook の提唱する臨界値に等しい．
Poon and Poon (1999) の共形法曲率 conformal normal curvature
Bd = C d /{tr(Hf )2 }1/2
は，Cook の法曲率と比較すると，固有ベクトルは変わらず固有値 (規準化法曲率) だけが
/(
τ˜j = τj
n
∑
)1/2
τk2
k=1
√
に変化する．そこで，τ˜1 ≥ τ˜2 ≥ .. ≥ τ˜k ≥ q/ n > τ˜k+1 ≥ · · · ≥ τ˜n ≥ 0 を満たす k 個の固
有値を取り上げて，これらを “q-inﬂuential”と呼ぶことを提案している．すべての方向の
228
日本統計学会誌第44巻第1号 2014
√
曲率が一様なら，τ˜j = 1/ n, j = 1, . . . , n となるので，その q 倍 (q は任意に与える) を超
えるかどうかで影響の大きさを測ろうという考えである．
数値例に適用してみる．Critchley (1985) や Shi (1997) において分析された Kendall
(1975) の土壌成分データ (n = 20, p = 4) を取り上げて，張・田中 (2002) は相関行列に
基づく PCA に対して Cook の局所影響分析を適用した．重み摂動 (相関行列に対しては
タイプ 1 と 2 の重みの微係数は同じになる) を導入し，すべての固有値と固有ベクトル
に関心のあるものとして分析すると，固有値 (法曲率) は大きい方から，τˆ1 = 3.756, τˆ2 =
2.564, τˆ3 = 2.420, τˆ4 = 0.989 となる．Cook (1986) の臨界値 τˆ = 2 を用いると第 3 固有値
までが採択され，Lu et al. (1997) の臨界値を用いると，パーセント点の選び方により臨界
値は τ0.05 = 3.31, τ0.25 = 2.47, τ0.30 = 2.36 となるので，5%なら第 1 固有値まで，25%な
ら第 2 固有値まで，30%なら第 3 固有値までが臨界値を超える．張・田中 (2002) の p. 11
に第 1 軸と第 2 軸，第 1 軸と第 3 軸の固有ベクトルの要素の散布図が描かれ，第 2 軸まで
を採用すれば#1 と#14 が individually inﬂuential，第 3 軸までを考えれば，さらに#13 と
#17 が jointly inﬂuential と示唆されている．この例に関して言えば，第 3 軸まで取り上げ
て検討してもよいように思われる．
シミュレーション (一種のブートストラップ法) を用いて臨界値を設定する考えもある．
Atkinson (1985) は回帰分析における外れ値あるいは影響観測値を検出するため，Student
化残差や修正 Cook 統計量の Q-Q プロット上に simulation envelope を描いている．説明
変数 X を固定し誤差項として正規乱数を用いてサイズ n の標本データを K 組生成し，そ
の結果得られる診断統計量の順序統計量を用いて信頼区間に相当する範囲を示す方法であ
る．Atkinson et al. (2004) では構造を考えない多変量 (p 変量) データの平均から各点への
マハラノビス 2 乗距離と自由度 p のカイ 2 乗分布の分位点との Q-Q プロット (軸のスケー
ルは平方根) を描き，生成した多数組 (例では 99 組) のサイズ n の多変量正規標本に基づ
く 90%envelope をその上に重ねて描いている．また，Brooks (1994) は推定された共分散
行列をもつ人工データを生成して PCA の固有値，固有ベクトル，固有空間への射影行列
などの影響統計量の臨界値を求めている．
4) Cook の局所影響分析の拡張 (i)．EM アルゴリズムに基づく最尤推定への拡張
Cook (1986) の局所影響分析は，observed-data log-likelihood に基づく式 (2.4) の尤度
距離 D(w) を目的関数とする影響グラフの，非摂動点における法曲率によって影響の大き
さを評価する方法である．しかし，欠測値や潜在変数を含むモデルなど observed-data に
基づく尤度距離を扱うことが困難な場合がある．そのような場合に対して，Zhu and Lee
(2001) は尤度距離の代わりに EM アルゴリズムの推定プロセスで得られる complete-data
229
多変量解析の感度分析
log-likelihood の条件付期待値
[
]
ˆ θ)
ˆ − Q(θˆw |θ)
ˆ
fQ (w) = 2 Q(θ|
を目的関数として用いた影響グラフを考えて非摂動点における法曲率 Cd を最大化し，対応す
[
]
る最大影響方向 dmax を求める方法を提案した．ただし，
θˆw := arg maxθ E Lc (θ, w|Yc )|Yo , θˆ
である．この方法は各種の不完全データを含むモデルや observed-data log-likelihood を用
いると高次元の積分が必要となるため実施が困難であった一般線形混合モデルの感度分析
に応用できることを示した．Song and Lee (2004) は同じ方法を 2 レベル潜在変数モデル
の影響分析に応用している．
5) Cook の局所影響分析の拡張 (ii)．摂動パラメータの変換に対して不変な影響指標
目的関数の非摂動点での勾配がゼロでない場合の Cook (1986) の法曲率は尺度不変でな
い．Poon and Poon (1999) の共形法曲率も同じ欠点をもつ．この問題に対して Zhu et al.
(2007) は微分幾何学 (計量テンソル，歪みテンソル，アフィン接続) を応用して勾配がゼロ
でない点においても摂動パラメータの任意の変換 (再パラメータ化) に対して不変な影響指
標を提案した．Chen et al. (2009) と Chen et al. (2010) は，Zhu et al. (2007) の方法を
EM アルゴリズムに基づく最尤法に応用して，潜在変数を含む非線形構造方程式モデルと
一般線形混合モデルの局所影響分析法を構築している．
Zhu et al. (2007) の方法は次の通りである．n 次元摂動パラメータ・ベクトル w を含
む確率密度関数 (probability density function, pdf) を p(x|θ, w), w ∈ Ω とし，その pdf は
Amari (1985), p. 16 の 4 つの正則条件を満たすとする．パラメータ θ, w のうち，w0 の近
傍での w の挙動に関心があるので θ は既知 (あるいは一定) とみなし，統計的摂動モデル
M = {p(x|θ, w) : w ∈ Ω} を n 次元多様体と考えて，M 上の任意のなめらかな関数 f (w)
に対する影響グラフ α(w) = (wT , f (w))T の形状を検討する．この影響グラフの点 w にお
ける接空間 Tw に対して，統計モデルに対する自然な基底 ∂i l(w|x, θ) (ただし ∂i = ∂/∂wi ，
∑
(1)
l(w|x, θ) = log p(x|θ, w)) を用いたベクトル空間 Tw = {d(x) | d(x) = i di ∂i l(w|x, θ)}
(1)
を考えると，2 つのベクトル空間は同形であり，Tw は接空間 Tw の確率変数表現あるいは
(1)
1-表現と呼ばれる．E(∂i l(x|θ, w)) = 0 がなりたつので，Tw に属する任意の確率ベクトル
の期待値はゼロである．確率変数表現の 2 つの接ベクトルの内積を共分散の形で定義する
と，基底ベクトル間の内積は
gij (w) = h∂i , ∂j i = Ew [∂i l(x|θ, w)∂j l(x|θ, w)]
と表される．gij (w), i, j = 1, 2, . . . , n は計量テンソルと呼ばれ，行列 G(w) = (gij (w)) は
摂動パラメータ・ベクトル w に関する期待 Fisher 情報行列を表す．接ベクトル d の長さ
230
日本統計学会誌第44巻第1号 2014
2
は kdk = dT G(w)d により評価でき，C : w(t) = (w1 (t), . . . , wn (t)) を多様体 M 上の 2 つ
の点 w(t1 ) と w(t2 ) を結ぶなめらかな曲線とするとき，2 点間の経路に沿った距離は計量テ
ンソルを用いて計算した微小距離を積分することによって得られる．Zhu et al. (2007) は，
さらに，アフィン接続 (Amari の α 接続) の係数 Γα
ijk を用いて任意の d 方向に対する 1 次
および 2 次の影響指標
F If,d = dT ∇f ∇Tf d/dT Gd
˜ f0 d/dT Gd
SIf,d = dT H
を提案した．そのうち，2 次の指標 SIf,d は計量テンソルのレビ・チビタ接続 (α = 0 の α
˜ 0 は測地線に沿った 2 次偏微分を求めるための修正項を含んだ共
接続) に基づいており，H
f
変ヘシアン行列と呼ばれる対称行列で，その (i, j) 要素は
[
˜0
H
f (w)
]
(i,j)
= ∂i ∂j f (w) −
∑
g sr (w)Γ0ijs (w)∂r f (w)
s,r
により与えられる．任意の非摂動点において定義される影響指標 F If,d と SIf,d は摂動パ
ラメータの変換に対して不変で，目的関数 f を k 倍すると 1 次の影響指標は k 2 倍に，ま
た 2 次の影響指標は k 倍になる．
固有値問題の摂動論と布置の変化の測定
3.
主成分分析，正準相関分析，対応分析や数量化法など多くの多変量解析において対称行
列の固有値問題が中心的な役割を果たす．これらの多変量解析における感度分析法を構築
する際，固有値や固有ベクトルあるいは固有空間を特徴づける行列などの (2 節で取り上げ
た影響関数や一定方向の摂動に対する方向微分などの)1 次元の摂動パラメータ ε や t に関
するべき級数展開を求めることができると便利である．3.1 節は最も基本的な場合で，関心
のある固有値が単純固有値であるときの固有値と固有ベクトルの展開式，3.2 節は重複固有
値を含む場合も含め，関心のある固有ベクトルの張る部分空間を特徴づける行列の展開式
を扱う．3.3 節以下では多変量解析の結果として得られた係数行列や負荷行列，スコア行列
などの布置の類似度を測定する RV 係数などの指標について議論する．簡単のため 1 次の
項までの展開式を示す．詳細は原論文を参照されたい．
3.1
固有値・固有ベクトルのべき級数展開
次のような固有値問題において行列 H が微小な摂動により変化するとき，固有値や固有
ベクトルがどのように変化するかを知りたいとする．
(H − λs I)vs = 0
多変量解析の感度分析
231
ただし，H は p × p 対称行列，λs , vs は s 番目の固有値 (単純固有値と仮定) とそれに対応
するノルム１に規準化された固有ベクトルである．対称性を保ちながら行列 H が H(ε) に
微小に変化し H(ε) は ε = 0 の近傍で収束するべき級数
H(ε) = H + εH (1) + (ε2 /2)H (2) + O(ε3 )
に展開できるとする．このとき，その固有値・固有ベクトルも収束するべき級数
2
(2)
3
λs (ε) = λs + ελ(1)
s + (ε /2)λs + O(ε )
vs (ε) = vs + εvs(1) + (ε2 /2)vs(2) + O(ε3 )
に展開でき (Rellich (1969), Kato (1980, 1982))，それらの展開式が種々の多変量解析にお
ける感度分析に応用された (古典的 MDS：Sibson (1979)；数量化３類・対応分析：Tanaka
(1984), Pack and Jolliﬀe (1992)； PCA：Critchley (1985), Pack et al. (1987), Tanaka
(1988)；数量化 2 類：Tanaka and Tarumi (1986))．固有値と固有ベクトルの 1 次の項の
係数はそれぞれ以下のように与えられる．
(1)
λ(1)
s = ass
∑
(λs − λr )−1 a(1)
vs(1) =
rs vr
r6=s
(1)
ただし，ars = vrT H (1) vs ，δrs はクロネッカーのデルタ，上付きの (j) は j 次の微係数を
表す．
Hampel や Belsley らの影響関数に対応する摂動を導入するとき，固有値・固有ベクトル
のべき級数展開における 1 次の項の係数が固有値・固有ベクトルの影響関数となる．固有
ベクトルの展開式に，(λs − λr )−1 , r = 1, 2, . . . , p, r 6= s が含まれていることからも明ら
かなように，注目している s 番目の固有値は単純固有値でなければならない．正確に等し
くなくても近い値をもつ固有値が存在するとき，この展開式を用いて計算すると固有ベク
トルの影響関数の絶対値が大きくなり，計算誤差が混入しやすい．単純固有値を仮定した
展開式と重複固有値を仮定した展開式に対して，計算誤差の大きさと固有値の近さの関係
をシミュレーションによって調べた結果が Tanaka and Tarumi (1986, 1988) に示されてい
る．当然のことながら，固有値が近くなると単純固有値を仮定して求めた展開式の精度は
悪い．
3.2
部分空間を特徴づける行列のべき級数展開
対応分析や主成分分析の感度分析において，摂動前後の固有ベクトルやスコアの散布図
を描くと，点間の相対的距離はほとんど変わらないが，個々の点の位置が大きく変わる場
合がある (たとえば，田中 (1984))．この現象は，得られる部分空間は安定しているが，そ
232
日本統計学会誌第44巻第1号 2014
の中の座標が大きく変化していることを意味する．多変量解析においては，布置の相対的
な位置に重要な情報が含まれる場合が多いが，その場合には固有値・固有ベクトルの変化
を個別に測るのではなく，プロットしたときの相対的な位置の変化あるいは部分空間の変
化の大きさを評価することが望ましい．そのような目的のため，次のような行列のべき級
数展開式を求めた．
P =
∑
vs vsT
(3.1)
λs vs vsT
(3.2)
s∈K
T =
∑
s∈K
ただし，K は関心のある固有値番号の集合を表し，値の等しい固有値は集合 K とその補
集合 K c の各々には含まれてもよいが，両方にまたがって含まれていないものと仮定する．
式 (3.1) と (3.2) で定義された P と T は，それぞれ行列 H の関心のある固有ベクトルの張
る部分空間への射影行列および H のスペクトル分解における関心のある固有値に対応する
部分項である．これらの展開式の 1 次の項の係数，すなわち，摂動パラメータに関する 1
次の微係数は
P (1) =
∑∑
(λs − λk )−1 ask (v s vkT +v k vsT )
(1)
s∈K k∈K
/
T (1) =
∑∑
s∈K r∈K
T
a(1)
sr vs vr +
∑∑
(3.3)
λs (λs − λk )−1 ask (v s vkT +v k vsT )
(1)
(3.4)
s∈K k∈K
/
で与えられる．P (1) と T (1) は主成分分析の感度分析に関連して Tanaka (1988) は上の仮定
の下で，2 次の項は Tanaka and Castano-Tostado (1990) は関心のある固有値が互いに近
くてもよいが厳密には等しくないことを仮定して，それぞれ初等的な方法で導いている．
また，Tanaka (1989) では対称行列の一般固有値問題の固有値・固有ベクトル，部分空間
への射影行列などの展開式についても議論している．その後，Watson (1983) が球面上の
統計学に関連して，重複固有値を許す仮定の下で P と T に関して 1 次の展開式を求めてい
ることを発見した．導出のテクニックとして我々は行列代数を用いたが，Watson (1983)
は Kato (1980) による複素解析におけるコーシーの積分定理を応用する方法を用いて導い
ており，高次の項も同じ考え方で導くことができる．Tanaka (1992) は，Watson と同じ方
法を用いて Tanaka (1988) や Tanaka and Castano-Tostado (1990) の 2 次までの展開式の
確認を行っている．また，シャトラン (Chatelin (1988)) の訳本「行列の固有値」(伊理正
夫，伊理由美訳，1993) に射影行列の Rellich-Kato 級数展開に関する一般的な記述を見つ
けた．展開式 (3.3) と (3.4) は主成分分析 (Tanaka (1988)) や因子分析 (Tanaka and Odaka
(1989a, b)) などの感度分析に応用されている．
多変量解析の感度分析
3.3
233
行列の類似度の指標と摂動の影響の測定
各種の多変量解析において，係数行列あるいは負荷行列，スコア行列が主要な結果と
して利用され，微小な摂動のそれらの行列への影響を評価したい場合がある．2 つの行列
X : n × p, Y : n × q を Rp , Rq の中の n 個の点と考えるときの 2 組の点の布置の類似
度，あるいは 2 つの行列の間の類似度を測る量としていろいろ提案されているが，そのう
ち Escouﬁer の RV 係数 (Escouﬁer (1973), Robert and Escouﬁer (1976))，柳井の一般化決
定係数 (柳井 (1974))，Benasseni の係数 (Benasseni (1990)) を取り上げ，これらを用いた
摂動の影響の測定について議論する．
1) Escouﬁer の RV 係数 (RV coeﬃcient)
RV 係数とは，X と Y を p および q 次元空間の中の点として表したときの布置の類似度を
√
2(1 − RV (X, Y )) = kS(X)/ kS(X)k − S(Y )/ kS(Y )kk
で定義される RV (X, Y ) を用いて評価しようというものである．ただし，S(X) = XX T ,
S(Y ) = Y Y T であり，kAk = {tr(AT A)}1/2 は A の Frobenius ノルムを表す．布置の原点
の位置が重要でなく点の間の相対的な位置関係の類似度を測りたい場合には，列方向に中
心化した X, Y を用いる．RV 係数は 0 ≤ RV ≤ 1 の範囲の値をとり，行列全体への定数
倍および直交変換により一致するとき 1，X T Y = 0 のとき 0 という値をとる．実際の計算
においては
RV (X, Y ) = tr(XX T Y Y T )/{tr(XX T )2 tr(Y Y T )2 }1/2
(3.5)
が使用される．
2) 柳井の一般化決定係数 (Generalized Coeﬃcient of Determination, GCD)
X と Y の列ベクトルの張る部分空間への射影行列をそれぞれ ΠX , ΠY とするとき，一
般化決定係数は
GCD(X, Y ) = tr(ΠX ΠY )
で定義される (柳井 (1974))．RV 係数の場合と同様，原点の位置が重要でないときには，列
方向に中心化しておく．X または Y の一方が 1 変量のとき重相関係数の 2 乗，X と Y の
両方が 1 変量のとき相関係数の 2 乗に等しい．
3) Benasseni の係数
Benasseni (1990) は主成分の係数行列の変化の大きさを RV 係数で測るとき，摂動パラ
メータ ε のべき級数で展開すると 1 次の項が消えるため変化に対して敏感でないとして，1
次の項の残る次の形の係数を提案した．
234
日本統計学会誌第44巻第1号 2014
ノルム 1 に規準化された固有ベクトルを列ベクトルとする 2 組の p × q 行列を A, B とす
る．これらの固有ベクトルによって張られる部分空間の近さを
ρ(A, B) = 1 − q −1
∑
kak − ΠB ak k
k
を用いて測ろうというのが Benasseni の考え方である．A の各列を B の情報を用いて説明
するという非対称な指標であり，逆向きの説明をする場合は A と B を置き換える．
4) 摂動の影響の測定
係数行列，負荷行列などの摂動後の値を表す行列 Lε が摂動前の値 L を中心にしてべき
級数 Lε = L + εL(1) + (ε2 /2)L(2) + O(ε3 ) の形に展開されるとき，摂動前後の L の間の RV
係数は
RV (L, Lε ) = 1 − (ε2 /2)[tr(S ∗(1) )2 − {tr(S ∗ S ∗(1) )}2 /tr(S ∗2 )]/tr(S ∗2 ) + O(ε3 )
(3.6)
の形に展開される．ただし，S ∗ = LLTε , Sε∗ = Lε LTε である．ε の 1 次の項は消えて leading
term は ε2 の項となる．2 次の項の係数のかぎカッコ内は非負で 1 次の微係数 S ∗(1) のみを
含んでいて実質的に１次の影響を表すと考えられ，2 次の項の係数を影響指標として利用
することができる (Castano-Tostado and Tanaka (1991))．
PCA の場合について，p × q の係数行列 V ，負荷行列 V Λ1/2 ，n × q の PC スコア行列
XV への摂動の影響を考えよう．係数行列や負荷行列を q 次元空間内の p 個の点とみると，
原点からの距離は重要な情報を持つので，RV 係数や GCD を用いるとき中心化するのは適
当でない．一方，スコア行列についてはスコアの原点に特別な意味がないので列方向に中
心化しておくのがよい．PCA における係数行列と負荷行列の摂動前後の RV 係数の展開式
は以下のようになり，3.2 節で述べた射影行列 P およびスペクトル分解の部分項 T の 1 次
の微係数 (影響関数) の関数となる．
RV (V, Vε ) = tr(P Pε )/q = 1 − (ε2 /2)tr(P (1) )2 /q + O(ε3 )
RV (V Λ1/2 , (V Λ1/2 )ε ) = 1 − (ε2 /2){tr(T (1) )2 − (tr(T T (1) ))2 /tr(Λ2 )}/tr(Λ2 ) + O(ε3 )
これらの展開式においても一般の行列 L の場合の式 (3.6) と同様に 1 次の項は消えて leading
term である 2 次の項の係数に P や T の 1 次の微係数のみが現れ，この項が実質的に 1 次
の影響を表す．また，摂動前後のスコア行列間の RV 係数はつぎのように表される．
RV (XV, XVε ) = 1 − (ε2 /2)R(1) + O(ε3 )
ただし
R(1) = {tr(P (1) XX T )2 − (tr(P X T XP (1) X T X))2 /tr(P X T X)2 }/tr(P X T X)2
235
多変量解析の感度分析
は非負の行列である．
一般化決定係数を用いると，係数行列および負荷行列に対して
GCD(V, Vε ) = GCD(V Λ1/2 , (V Λ1/2 )ε ) = tr(P Pε )/q = 1 − (ε2 /2)tr(P (1) )2 /q + O(ε3 )
が成り立つ．
Benasseni (1990) の係数を用いて PC 係数行列の変化を測定すると
∑
∑
(1) ρ(V, Vε ) = 1 − q −1
kv k − Pε vk k = 1 − |ε| q −1
P vk + O(ε2 )
k
k
となる．Castano-Tostado and Tanaka (1990) は ρ(V, Vε ) の |ε| の係数の上限が kP (1) k に等
∑
しいこと，また，Prendergast (2008) は ρ(A, B) の定義式において右辺の k kak − ΠB ak k
の各項のノルムを 2 乗ノルムに置き換えると RV 係数や GCD と同じになることを指摘し
ている．
Prendergast (2008) は PC 係数行列の部分空間への影響指標として，1 − RV (V, Vε ) の ε2
の項の係数を用い，変数の数が大きいいわゆる高次元の場合の共分散行列の PCA におけ
る影響指標として PC スコアを用いた次の式を導いている．
lim (1 − RV (V, Vε ))/ε2 =
ε→0
∑∑
(λs − λk )−2 ys2 yk2 /q
s∈K k∈K
/
5) Smilde らの修正 RV 係数
Smilde et al. (2009) は変数の数の多い遺伝データの分析において RV 係数が大きくなり，
類似度を過大に評価しやすい現象に注目して RV 係数の性質を検討し，改良版として修正 RV
係数 (modiﬁed RV-coeﬃcient) を提案した．彼らは 2 組のデータ行列 X : n × p1 , Y : n × p2
の各要素が互いに独立に標準正規分布に従う場合について検討し，n を大きくすると，RV
係数 RV (X, Y ) の値は
RV (X, Y ) ∼
= p1 p2 /[{p21 + (n + 1)p1 }1/2 {p22 + (n + 1)p2 }1/2 ]
と近似でき，その値はゼロでなく正の値をもつこと，n を固定して p1 , p2 を大きくすると
RV (X, Y ) → 1 となることを示して，高次元データの場合に RV 係数が大きくなりやすい
現象を説明した．RV 係数のこの欠点を解消するため，彼らは XX T と Y Y T の対角要素を
ゼロで置き換えた行列を式 (3.5) に代入して RV 係数を計算し，その値 RV2 を修正 RV 係
数 (modiﬁed RV-coeﬃcient) と呼んで RV 係数の代わりに用いることを提案した．シミュ
レーションにより互いに独立なデータ行列の間の修正 RV 係数 RV2 は，変数の数が多い場
合でも n を大きくするとゼロに近づくことを確認している．
p × q 行列 L の摂動前後の変化の測定に修正 RV 係数 RV2 を応用することは容易で，式
(3.6) の ε2 の係数に現れる行列 S ∗ と S ∗(1) の対角要素をゼロで置き換えて計算すればよい．
236
4.
日本統計学会誌第44巻第1号 2014
おわりに
2 節において，感度分析・影響分析の主要なツールである影響関数法 (IF 法と略記) と
Cook の局所影響分析法 (LI 法と略記) について概説し，両者の関係について議論した．IF
法は観測値の影響の評価に焦点を絞っており影響指標の統計的な意味が明確である．これ
に対して，LI 法は任意の摂動パラメータに対する目的関数の曲面の法曲率という幾何学的
な量に基礎をおいているため，影響指標の統計的な意味は直ちには明らかでないが，その
反面，摂動パラメータが自由に導入できて柔軟に感度分析に行うことができるという利点
をもち，より広い目的の感度分析に適用できる．目的を分散摂動あるいは重み摂動に対す
る観測値の影響の評価に絞ると，2.6 節で説明したように IF 法と LI 法は双対的な関係に
あって，局所影響分析に関して同じ情報が提供してくれる．この範囲では，一方の方法か
ら得られた結果は他方の観点からも解釈が可能であり，その意味で両者は補完的な役割を
果たすということができる．
IF 法の影響関数の PCA は “経験影響関数の加法性”に基礎をおいているので，(微分に基
づいているので当然であるが)LI 法も影響の加法性に基づいていることになる．2 節で述べ
たように，影響関数の標本版のうち微分の形で定義される経験影響関数 EIF と差分 (標本全
体を用いた場合と観測値を落とした場合の差) から定義される標本影響関数 SIF に関して，
サンプルサイズを大きくするとき，PCA をはじめ多くの多変量解析において SIF ∼
= EIF
が成り立つのに対して，回帰係数に対しては SIF ∼
= EIF/(1 − hii )，ただし，hii はてこ
比，となり，ほかの多変量解析の場合と比較して，回帰分析においては加法性が成り立ち
にくい場面が多いと考えられ，2.7 節のロバスト化の考え方が LI 法においても実用上重要
になってくると思われる．
LI 法は上にあげた利点のため，経時測定データや潜在変数モデルなど複雑な問題を含め
て幅広い分野に利用されている．Verbeke et al. (2001) は脱落を含む動物実験データの分
析において，missing-at-random (MAR) model を非摂動モデルとし nonrandom model の
方向の摂動パラメータ・ベクトルを導入して LI 法による脱落モデルの感度分析法を提案し
ているが，LI 法の特徴を活かした感度分析といえよう．2.8 節で取り上げた LI 法のいくつ
かの拡張は LI 法の適用範囲を広げたものである．EM アルゴリズムへの拡張のほうはすで
にいろいろなモデルに応用がされているが，Zhu et al. (2007) の目的関数を尤度距離から
任意の関数に拡げた拡張については，まだ応用例は少なく今後の活用が期待される．
2.2 節では重み (摂動パラメータ) の定義の仕方として，タイプ 1 とタイプ 2 の 2 つのタ
イプについて述べた．そのうち，タイプ 2 の重みに関する微係数は Hampel の影響関数の
場合と同様，微係数で測る影響の平均がゼロになるのに対し，タイプ 1 の重みに関する微
係数の平均はゼロでなく，サンプルサイズが 1 増加することに対応する影響を受けること
多変量解析の感度分析
237
になる．LI 法においても，観測値の影響の分析という観点からは摂動パラメータとして平
均的な影響がゼロになるように設定することが望ましいが，Critchley and Marriott (2004)
の 1 次および 2 次の微係数を w
¯ = c を満たすアフィン部分空間に射影した上で議論する試
みを除き，あまりされていない．タイプ 2 の重みを用いると，重みベクトルの次元は n − 1
となる．Poon and Poon (1999) の共形法曲率や Zhu et al. (2007) の拡張 LI 法では重みベ
クトルの要素が 1 次独立であることを前提としているので，タイプ 2 の重みを用いてこれ
らの方法を利用する際には個別に工夫が必要となる．一般論として重みベクトルが 1 次独
立でない場合への LI 法の拡張が望まれる．
3 節の固有値問題の摂動論については，高次元データの主成分の部分空間への影響指標の
計算を p 次元の P (1) でなく n 次元の PC スコアを用いて計算する工夫 (Prendergast (2008))
が注目される．これは遺伝子データのように p が大きくなる状況に対応した改良といえる．
現段階では，共分散行列の PCA の場合に限られるが，今後そういった観点からの研究が
望まれる．
近年，分析対象のデータの性質や量が大きく変わり，それに伴って新しい分析法の必要
性が高まってきている．データの量や性質が変わっても分析結果の安定性・信頼性の評価
が重要な問題であることには変わりはない．摂動パラメータを柔軟に定義できる LI 法の適
用場面が多いと思われる．近い将来，新しいタイプのデータに対する感度分析法が開発さ
れることを期待したい．
参考文献
Amari, S. (1985). Diﬀerential-Geometrical Methods in Statistics, Springer-Verlag.
Anderson, T. W. (1984). An Introduction to Multivariate Statistical Analysis, 2nd Edition, Wiley.
Atkinson, A. C. (1985). Plots, Transformations, and Regression: An Introduction to Graphical Methods of
Diagnostic Regression Analysis, Oxford University Press.
Atkinson, A. C. and Riani, M. (2000). Robust Diagnostics Regression Analysis, Springer-Verlag.
Atkinson, A. C., Riani, M. and Cerioli, A. (2004). Exploring Multivariate Data with Forward Search, SpringerVerlag.
Belsley, D. A., Kuh, E. and Welsch, R. E. (1980). Regression Diagnostics: Identifying Inﬂuential Data and
Sources of Collinearlity, Wiley.
Benasseni, J. (1990). Sensitivity coeﬃcients for the subspaces spanned by principal components, Commun. Stat.
Theory Methods, 19, 2021–2034.
Brooks, S. P. (1994). Diagnostics for principal components: Inﬂuence functions as diagnostic tool, The Statistician, 43, 483–494.
Carroll, R. J. and Ruppert, D. (1988). Transformation and Weighting in Regression, Chapman and Hall.
Castano-Tostado, E. and Tanaka, Y. (1989). Sensitivity analysis in canonical correlation models for contingency
tables, J. Jpn. Soc. Comput. Stat., 2, 31–54.
Castano-Tostado, E. and Tanaka, Y. (1990). Some comments on Escouﬁer’s RV-coeﬃcient as a sensitivity
measure in principal component analysis, Commun. Stat. Theory Methods, 19, 4619–4626.
Castano-Tostado, E. and Tanaka, Y. (1991). Sensitivity measures of inﬂuence on the loading matrix in exploratory factor analysis, Commun. Stat. Theory Methods, 20, 1329–1343.
Chatterjee, S. and Hadi, A. S. (1988). Sensitivity Analysis in Linear Regression, Wiley.
238
日本統計学会誌第44巻第1号 2014
Chatelin, F. (1988). Valeurs propres de matrices, Masson. (シャテラン著, 伊理正夫, 伊理由美訳 (1993).『行列
の固有値』シュプリンガー・フェアラーク東京).
Chen, F., Zhu, H. T. and Lee, S. Y. (2009). Perturbation selection and local inﬂuence analysis for nonlinear
structural equation model, Psychometrika, 74(3), 493–516.
Chen, F., Zhu, H. T., Song, X. Y. and Lee, S. Y. (2010). Purturbation selection and local inﬂuence analysis of
generalized linear mixed models, J. Comput. Graph. Stat., 19(4), 826–842.
Cook, R. D. (1986). Assessing local inﬂuence, J. R. Stat. Soc. B , 48, 133–169.
Cook, R. D. and Weisberg, S. (1982). Residuals and Inﬂuence in Regression, Chapman and Hall.
Cook, R. D. and Weisberg, S. (1999). Applied Regression including Computing and Graphics, Wiley.
Critchley, F. (1985). Inﬂuence in principal component analysis, Biometrika, 72, 627–636.
Critchley, F. and Marriott, P. (2004). Data-informed inﬂuence analysis, Biometrika, 91, 125–140.
Enguix-Gonzalez, A. et al. (2005). Inﬂuence analysis in principal component analysis through power series
expansions, Commun. Stat. Theory Methods, 34, 2025–2046.
Escouﬁer, Y. (1973). Le traitment des variables vectorierlles, Biometrics, 29, 751–760.
Fung, W. K. and Kwan, C. K. (1997). A note on local inﬂuence based on normal curvature, J. R. Stat. Soc.
B , 59, 839–843.
Hampel, F. R. (1974). The inﬂuence curve and its role in robust estimation, J. Am. Stat. Assoc., 69, 383–393.
Kato, T. (1980). Perturbation Theory of Linear Operators, Springer-Verlag.
Kato, T. (1982). A Short Introduction to Perturbation Theory of Linear Operators, Springer-Verlag. (加藤敏夫
著, 丸山徹訳 (1999). 『行列の摂動』シュプリンガー・フェアラーク東京)
Kendall, M. G. (1975). Multivariate Analysis, Griﬃn.
Kwan, C. W. and Fung, W. K. (1998). Assessing local inﬂuence for speciﬁc restricted likelihood: application to
factor analysis, Psychometrika, 63, 35–46.
Laurance, A. J. (1988). Regression transformation diagnostics using local inﬂuence, J. Am. Stat. Assoc., 83,
1067–1072.
Lu, J., Ko, D. and Chang, T. (1997). The standardized inﬂuence matrix and its applications, J. Am. Stat.
Assoc., 92(440), 1572–1580.
Mori, Y., Watadani, S., Tarumi, T. and Tanaka, Y. (1998). Development of statistical software SAMMIF for
sensitivity analysis in multivariate methods, in COMPSTAT 98 (eds. R. Payne and P. Green), Physica-Verlag,
pp. 395–400.
森裕一，山本義郎，渡谷真吾，尾高好政，垂水共之，田中豊 (1997). 「感度分析プログラム SAMMIF」
『計算機統
計学』10(2), 145–153.
Pack, P. and Jolliﬀe, I. T. (1992). Inﬂuence in correspondence analysis, Appl. Stat., 41, 365–380.
Pack, P., Jolliﬀe, I. T. and Morgan, B. J. T. (1987). Inﬂuential observations in principal component analysis:
A case study, J. Appl. Stat., 15, 39–50.
Poon, W. Y. and Poon, Y. S. (1999). Conformal normal curvature and assessment of local inﬂuence, J. R. Stat.
Soc. B , 61, 51–61.
Pregibon, D. (1981). Logistic regression diagnostics, Ann. Stat., 9, 705–724.
Prendergast, L. A. (2008). A note on sensitivity of principal component subspaces and the eﬃcient detection
of inﬂuential observations in high dimensions, Electron. J. Stat., 2, 454–467.
Rellich, F. (1969). Perturbation Theory of Eigenvalue Problems, Gordon and Breach.
Robert, P. and Escouﬁer, Y. (1976). A unifying tool for linear multivariate statistical methods: The RVcoeﬃcient, Appl. Stat., 25, 257–265.
Shi, L. (1997). Local inﬂuence in principal components analysis, Biometrika, 84, 175–186.
Sibson, R. (1979). Studies in the robustness of multidimensional scaling: perturbation analysis of classical
scaling, J. R. Stat. Soc. B , 41, 217–229.
Smilde, A. K., Kiers, H. A. L., Bijlsma, S., Rubingh, C. M. and van Erk, M. J. (2009). Matrix correlations for
high dimensional data: The modiﬁed RV-coeﬃcient, Bioinformatics, 25(3), 401–405.
多変量解析の感度分析
239
Song, X. Y. and Lee, S. Y. (2004). Local inﬂuence analysis of two-level latent variable models with continuous
and polytomous data, Statistica Sinica, 14, 317–332.
Tanaka, Y. (1984). Sensitivity analysis in Hayashi’s third method of quantiﬁcation, Behaviormetrika, 16, 31–44.
田中豊 (1984). 「数量化理論における感度分析とその応用」『品質』14, 21–29.
Tanaka, Y. (1988). Sensitivity analysis in principal component analysis: Inﬂuence on the subspace spanned
by principal components, Commun. Stat. Theory Methods, 17, 3157–3175. (Corrections , Commun. Stat.
Theory Methods, 18, 4305 (1989))
Tanaka, Y. (1989). Inﬂuence functions related to eigenvalue problems which appear in multivariate analysis,
Commun. Stat. Theory Methods, 18, 4067–4084.
Tanaka, Y. (1992). Sensitivity analysis in multivariate methods: Principles, methods and software, Technical
Report No. 52, Okayama Statisticians Group, 1–32.
田中豊 (1992). 「多変量解析における感度分析」『行動計量学』19, 3–17.
Tanaka, Y. (1994). Recent advance in sensitivity analysis in multivariate methods, J. Jpn. Soc. Comput. Stat.,
7, 1–25.
Tanaka, Y. and Castano-Tostado, E. (1990). Quadratic perturbation expansions of certain functions of eigenvalues and eigenvectors and their application to sensitivity analysis in multivariate methods, Commun. Stat.
Theory Methods, 19, 943–965.
Tanaka, Y. and Odaka, Y. (1989a). Inﬂuential observations in principal factor analysis, Psychometrika, 54,
475–485.
Tanaka, Y. and Odaka, Y. (1989b). Sensitivity analysis in maximum likelihood factor analysis, Commun. Stat.
Theory Methods, 18, 3991–4010.
Tanaka, Y. and Tarumi, T. (1986). Sensitivity analysis in Hayashi’s second method of quantiﬁcation, J. Jpn.
Stat. Soc., 16, 44–60.
Tanaka, Y. and Tarumi, T. (1988). A numerical investigation on sensitivity analysis in multivariate methods,
in Data Analysis and Informatics V (ed. E. Diday), Elsevier Science Publisher, pp. 291–301.
Tanaka, Y. and Watadani, S. (1992). Sensitivity analysis in covariance structure analysis with equality constraints, Commun. Stat. Theory Methods, 21, 1501–1515.
Tanaka, Y. and Watadani, S. (1994). Unmasking inﬂuential observations in multivariate methods, in COMPSTAT ’94 (eds. R. Dutter and W. Grossman), Physica-Verlag, pp. 292–297.
Tanaka, Y. and Zhang, F. (1999). R-mode and Q-mode inﬂuence analyses in statistical modelling: Relationship
between inﬂuence function approach and local inﬂuence approach, Comput. Stat. Data Anal., 32, 197–218.
Tanaka, Y., Castano-Tostado, E. and Odaka, Y. (1990). Sensitivity analysis in factor analysis: Methods and
software, in COMPSTAT 1990 (eds. K. Momirovic and V. Mildner), Physica-Verlag, pp. 205–210.
Tanaka, Y., Watadani, S. and Moon, S. H. (1991). Inﬂuence in covariance structure analysis: With an application
to conﬁrmatory factor analysis, Commun. Stat. Theory Methods, 20, 3805–3821.
Tanaka, Y., Zhang, F. and Yang, W. S. (2002). On local inﬂuence in canonical correlation analysis, Commun.
Stat. Theory Methods, 31, 2325–2347.
Tanaka, Y., Zhang, F. and Mori, Y. (2003). Local inﬂuence in principal component analysis: Relationship
between the local inﬂuence and inﬂuence function approaches revisited, Comput. Stat. Data Anal., 44, 143–
160.
Tarumi, T. and Tanaka, Y. (1986). Statistical software SAM—sensitivity analysis in multivariate methods, in
COMPSTAT 1986 (eds. F. De Antoni, N. Lauro and A. Rizzi), Physica-Verlag, pp. 351–356.
Tsai, C. and Wu, X. (1992). Transformation-model diagnostics, Technometrics, 34, 197–202.
Watson, G. S. (1983). Statistics on Spheres, John Wiley & Sons.
Verbeke, G., Molenberghs, G., Thijs, H., Lasaﬀre, E. and Kenward, M. G. (2001). Sensitivity analysis for
nonrandom dropout: A local inﬂuence approach, Biometrics, 57(1), 7–14.
Wu, X. and Luo, Z. (1993). Second-order approach to local inﬂuence, J. R. Stat. Soc. B , 55, 929–936.
Yamanishi, Y. and Tanaka, Y. (2005). Sensitivity analysis in functional principal component analysis, Comput.
Stat., 20, 311–326.
240
日本統計学会誌第44巻第1号 2014
柳井晴夫 (1974). 「一般化決定係数による多変量解析の各種技法の統一的試み」『行動計量学』1, 46–54.
張方紅，田中豊 (2002). 「相関行列に基づく主成分分析における Cook の局所影響分析」
『計算機統計学』15, 1–17.
Zhu, H. T. and Lee, S. H. (2001). Local inﬂuence for incomplete-data analysis, J. R. Stat. Soc. B , 116–126.
Zhu, H. T., Ibrahim, J. G., Lee, S. H. and Zhang, H. P. (2007). Perturbation selection and inﬂuence measures
in local inﬂuence analysis, Ann. Stat., 36(6), 2565–2588.

Download Report