LOESS curve fitted to a population sampled from a sine wave with uniform noise added. The LOESS curve approximates the original sine wave.
局所回帰 (きょくしょかいき、英語 : local regression )または局所多項式回帰 (きょくしょたこうしきかいい、英語 : local polynomial regression )は、移動回帰 (いどうかいき、英語 : moving regression )とも呼ばれ、移動平均 や多項式回帰 を一般化したものである。
概要
局所回帰の最も一般的な方法がLOESS (locally estimated scatterplot smoothing) およびLOWESS (locally weighted scatterplot smoothing) であり、いずれも と発音される。いずれもノンパラメトリック回帰の手法であり、多数の回帰モデルをK近傍法 に基づくメタモデルで組み合わせる。
LOESSは、線形最小二乗回帰 の単純さと非線形回帰 の柔軟性の多くを兼ね備えている。グローバルな関数を指定する必要はなく、データの局所的部分集合localized subsets of dataに単純なモデルを当てはめればよい。最小二乗回帰と比較すると計算量は膨大である。この統計手法で得られる滑らかな曲線は、loess曲線ないしlowess曲線と呼ばれる。
モデルの定義
1964年、SavitskyとGolayがLOESSと等価な手法を提案し、Savitzky–Golayフィルタと呼ばれるようになった。1979年、William S. Clevelandがこの手法を再発見し、別の名前を付けた。1988年、ClevelandとSusan J. Devlinが、この手法をさらに発展させた。
データセット の範囲内の各ポイントで、低次の局所多項式 がデータの部分集合にフィットされる。多項式は、応答が推定されるポイントに近くのポイントに重みを与える、重み付き最小二乗法を用いてフィットされる。各データ点のそれぞれについて回帰関数の値が計算されたところで、LOESSのフィットが完了する。
データの局所的部分集合
LOESSの重み付き最小二乗法によるフィッティングに使用されるデータの部分集合 は最近傍アルゴリズムによって決定される。平滑化パラメータ
α
{\displaystyle \alpha }
は、各局所多項式を適合させるために、どのくらいの割合のデータを使用するかを決定する。
k 次の多項式のフィッティングにはk + 1 以上のポイントが必要であるため、平滑化パラメータ
α
{\displaystyle \alpha }
は、
(
λ
+
1
)
/
n
{\displaystyle \left(\lambda +1\right)/n}
と1の間にある必要がある。ここで、
λ
{\displaystyle \lambda }
は局所多項式の次数を示す。
α
{\displaystyle \alpha }
が小さいほど回帰関数がデータに近くなるが、データの変動に伴うブレが大きくなる。
局所多項式の次数
局所多項式は、ほとんどの場合、1次か2次である。より高次の多項式は理論的には有効であるが、「どんな関数も局所では低次の多項式で近似できる」というLOESS の精神にはそぐわないし、過剰適合 のリスクがある。
重み関数
重み関数は、説明変数空間において推定点に近いデータ点に最も大きな重みを与え、最も遠いデータ点に最も小さな重みを与える。
伝統的には三次元重み関数 が用いられる。
w
(
x
)
=
(
1
−
|
d
|
3
)
3
{\displaystyle w(x)=(1-|d|^{3})^{3}}
ここで、d はデータポイント間の距離で、0から1の範囲にスケーリングされる。
ターゲット空間
R
m
{\displaystyle \mathbb {R} ^{m}}
上の計量
w
(
x
,
z
)
{\displaystyle w(x,z)}
x
,
z
∈
R
n
{\displaystyle x,z\in \mathbb {R} ^{n}}
による線形回帰の一般化を考える。
x
↦
x
^
:=
(
1
,
x
)
{\displaystyle x\mapsto {\hat {x}}:=(1,x)}
により
n
{\displaystyle n}
個の入力パラメータを
R
n
+
1
{\displaystyle \mathbb {R} ^{n+1}}
に埋め込んで、次の損失関数を考える。
RSS
x
(
A
)
=
∑
i
=
1
N
(
y
i
−
A
x
^
i
)
T
w
i
(
x
)
(
y
i
−
A
x
^
i
)
.
{\displaystyle \operatorname {RSS} _{x}(A)=\sum _{i=1}^{N}(y_{i}-A{\hat {x}}_{i})^{T}w_{i}(x)(y_{i}-A{\hat {x}}_{i}).}
ここで、
A
{\displaystyle A}
は
m
×
(
n
+
1
)
{\displaystyle m\times (n+1)}
の実行列であり、
w
i
(
x
)
:=
w
(
x
i
,
x
)
{\displaystyle w_{i}(x):=w(x_{i},x)}
と定義される。添え字i は訓練データの入出力のベクトルを示す。
w
{\displaystyle w}
は計量なので対称な正定値行列であり、
w
=
h
2
{\displaystyle w=h^{2}}
を満たす対称行列
h
{\displaystyle h}
が存在する。損失関数は次のように変形できる。
y
T
w
y
=
(
h
y
)
T
(
h
y
)
=
Tr
(
h
y
y
T
h
)
=
Tr
(
w
y
y
T
)
{\displaystyle y^{T}wy=(hy)^{T}(hy)=\operatorname {Tr} (hyy^{T}h)=\operatorname {Tr} (wyy^{T})}
ベクトル
y
i
{\displaystyle y_{i}}
を
m
×
N
{\displaystyle m\times N}
型の行列
Y
{\displaystyle Y}
、ベクトル
x
^
i
{\displaystyle {\hat {x}}_{i}}
(
n
+
1
)
×
N
{\displaystyle (n+1)\times N}
型の行列
X
^
{\displaystyle {\hat {X}}}
とすることで、損失関数は次のように変形できる。
Tr
(
W
(
x
)
(
Y
−
A
X
^
)
T
(
Y
−
A
X
^
)
)
{\displaystyle \operatorname {Tr} (W(x)(Y-A{\hat {X}})^{T}(Y-A{\hat {X}}))}
ここで、
W
{\displaystyle W}
は
N
×
N
{\displaystyle N\times N}
の対角行列 であり、 その成分は
w
i
(
x
)
{\displaystyle w_{i}(x)}
である。
A
{\displaystyle A}
に関して微分した値をゼロとすることで
A
X
^
W
(
x
)
X
^
T
=
Y
W
(
x
)
X
^
T
.
{\displaystyle A{\hat {X}}W(x){\hat {X}}^{T}=YW(x){\hat {X}}^{T}.}
さらに、正方行列
X
^
W
(
x
)
X
^
T
{\displaystyle {\hat {X}}W(x){\hat {X}}^{T}}
が可逆行列 であるとき、損失関数
RSS
x
(
A
)
{\displaystyle \operatorname {RSS} _{x}(A)}
は下記で最小値を取る。
A
(
x
)
=
Y
W
(
x
)
X
^
T
(
X
^
W
(
x
)
X
^
T
)
−
1
.
{\displaystyle A(x)=YW(x){\hat {X}}^{T}({\hat {X}}W(x){\hat {X}}^{T})^{-1}.}
w
(
x
,
z
)
{\displaystyle w(x,z)}
として、主にガウス関数 が選択される。
w
(
x
,
z
)
=
exp
(
−
‖
x
−
z
‖
2
2
α
2
)
{\displaystyle w(x,z)=\exp \left(-{\frac {\|x-z\|^{2}}{2\alpha ^{2}}}\right)}
利点
平滑化パラメータと局所多項式の次数を与えるだけで、モデルをサンプルデータに柔軟に適合させることができる。
欠点
局所的なデータ構造に基づいてフィッティングするため、充分な大きさの標本が必要である。また、数式で簡単に表現できる回帰関数を生成しないため、分析結果を他者に伝えることが困難である。また、他の最小二乗法と同様に、外れ値 の影響を受けやすい。
脚注
出典
関連項目
外部リンク
実装
この記事にはパブリックドメイン である、アメリカ合衆国連邦政府 が作成した次の文書本文を含む。アメリカ国立標準技術研究所 .