はじめに
東京大学/株式会社Nospare リサーチャーの栗栖です.
この記事から数回にわたり,ノンパラメトリックな回帰分析の統計手法の一つである局所多項式回帰について紹介します.今回は局所多項式回帰の特別な場合として,局所定数回帰について解説します.
モデル
ノンパラメトリックな回帰モデル
Y_i = m(\boldsymbol{X}_i) + U_i,\ E[U_i|\boldsymbol{X}_i] = 0
から独立同分布 (independent and identically distributed, i.i.d.) なデータ $(Y_i, \boldsymbol{X}_i)$, $i = 1,\dots,n$ が得られたとします.ただし,$Y_i \in \mathbb{R}$, $\boldsymbol{X}_i=(X_{1,i},\dots,X_{d,i})' \in \mathbb{R}^d$とします.
以下ではデータから回帰関数 $m(x) = E[Y_i|\boldsymbol{X}_i = x]$, $x = (x_1,\dots,x_d) \in \mathbb{R}^d$ を推定するという問題を考えましょう.推定にあたり,上記のモデルに対して以下を仮定します.
(A1:LC)
[1] $X_i$は密度関数 $f$ をもち,$f$ は $x$ のある近傍 $N_x$ で2回連続偏微分可能.
[2] $m$ は $N_x$ で2回連続偏微分可能.
[3] $\sigma^2(x) = E[U_i^2|\boldsymbol{X}_i = x]$ は $x$ で連続.
また $x \in \mathbb{R}^d$ で2回連続偏微分可能な関数 $g: \mathbb{R}^d \to \mathbb{R}$ に対して,
\partial_j g(x) = {\partial \over \partial x_j}g(x),\ \partial_{jk}g(x) = {\partial^2 \over \partial x_j \partial x_k}g(x).
とします.上記の仮定の下で,平均関数 $m$ は $x$ の周りで以下のように近似できることが分かります.
\begin{align*}
m(\boldsymbol{X}_i) &= m(x) + \sum_{j=1}^d \partial_j m(x)(\boldsymbol{X}_i - x) + {1 \over 2}\sum_{j,k=1}^d \partial _{jk} m(\tilde{\boldsymbol{X}}_i)(\boldsymbol{X}_i - x)^2,\\ \tilde{\boldsymbol{X}}_i &= \theta_i x + (1-\theta_i)\boldsymbol{X}_i,\ \theta_i \in (0,1].
\end{align*}
局所定数回帰
回帰関数 $m$ の推定にあたり,以下の最適化問題を考えましょう.
\begin{align*}
\hat{m}(x) &= \underset{\beta_0 \in \mathbb{R}}{\text{argmin}}\sum_{i=1}^n (Y_i - \beta_0)^2 K\left({\boldsymbol{X}_i - x \over h}\right),\\
K\left({\boldsymbol{X}_i - x \over h}\right) &= K\left({X_{1,i} - x_1 \over h},\dots, {X_{d,i} - x_d \over h}\right)
\end{align*}
ここで,$h = h_{n} \to 0$, $n \to \infty$ はバンド幅,$K:\mathbb{R}^d \to \mathbb{R}$ はカーネル関数と呼ばれ,以下を満たすとします.
(A2:LC)
[1] $\int K(u)du=1$
[2] $K(u) = K(-u)$
[3] $K(u)$ は有界かつ $[-1,1]^d$ に台をもつ ($K(u) = 0$, $\|u\| >1$)
上記の最適化問題は $\boldsymbol{X}_i = x$ の近傍で「局所」的に $Y_i$ を「定数」に「回帰」していると解釈できるので,「局所定数回帰」と呼ばれます.また得られる解は「局所定数推定量」(local constant estimator, LC 推定量) と呼ばれることもあります.特にこの方法では,$\boldsymbol{X}_i = x$ の近傍で
Y_i \approx m(\boldsymbol{X}_i) \approx m(x)
と近似できることを利用して $m(x)$ を推定しています.上記の最適化問題を解くと,$m(x)$ の LC推定量が以下で与えられることが分かります.
\hat{m}(x) = {\sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right)Y_i \over \sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right)}.
ここで,
W(\boldsymbol{X}_i) = {K\left({\boldsymbol{X}_i - x \over h}\right) \over \sum_{j=1}^{n}K\left({\boldsymbol{X}_j - x \over h}\right) }
とおくと,
\hat{m}(x) = \sum_{i=1}^{n}W(\boldsymbol{X}_i)Y_i,\ \sum_{i=1}^{n}W(\boldsymbol{X}_i) = 1
と書けることから,$W(\boldsymbol{X}_i)$ は $m(x)$ を推定する際の $Y_i$ の $\boldsymbol{X}_i$ 貢献度と解釈することができ, $\hat{m}(x)$ は $\boldsymbol{X}_i$ の貢献度 $W(\boldsymbol{X}_i)$ による $Y_i$ の加重平均の形で表現されることが分かります.
LC推定量のその他の構成方法
LC推定量は局所的な重み付き回帰の考え方以外の方法でも導出することができます.ここでは条件付き期待値の推定量としての導出方法について紹介します.まずいくつか記号を導入しておきます.
- $f(y,x)$: $(Y_i,\boldsymbol{X}_i)$ の同時密度関数,
- $f(y|x)$: $\boldsymbol{X}_i = x$ を与えたときの $Y_i$ の条件付密度関数.$f(y|x) = f(y,x)/f(x)$ であることに注意しておきます.
カーネル密度推定
ここで$X_i$ の密度関数 $f(x)$ の推定方法についても紹介しておきます.以下の推定量を考えましょう.
\hat{f}(x) = {1 \over nh^d}\sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right).
カーネル関数 $K$ が $K \geq 0$ を満たすとき,$\hat{f}(x)$ は密度関数になることが分かります.実際,$\hat{f}(x) \geq 0$ であり,
\begin{align*}
\int \hat{f}(x)dx &= {1 \over nh^d}\sum_{i=1}^{n}\int K\left({\boldsymbol{X}_i - x \over h}\right)dx = {1 \over n}\sum_{i=1}^{n}\int K(u)du = 1
\end{align*}
となります.
カーネル密度推定量の漸近正規性
$\hat{f}(x)$ の漸近的性質について紹介しておきます.(A1:LC), (A2:LC) に加えていくつか仮定を置きます.
- $n \to \infty$ として,$h \to 0$, $nh^d \to \infty$, $nh^{d+4} \to c \in [0,\infty)$ .
- $x \in \mathbb{R}^d$ に対して $f(x)>0$.
このとき,$n \to \infty$ として,(A1:LC)[1], (A2:LC) の下で以下が成り立ちます.
\begin{align*}
\sqrt{nh^d} &\left(\hat{f}(x) - f(x) - b_{f,n}(x)\right) \stackrel{d}{\to} N\left(0, \quad f(x)\int K^2(u)du\right),\\
b_{f,n}(x) &= h^2\sum_{j,k=1}^d {\partial_{jk} f(x) \over 2}\int u_ju_k K(u)du.
\end{align*}
LC推定量の導出
$x \in \mathbb{R}^d$ に対して $f(x)>0$ とすると,
m(x) = E[Y_i|\boldsymbol{X}_i = x] = \int yf(y|x)dy = {1 \over f(x)}\int yf(y,x)dy
が成り立ちます.ここで,$f(x)$, $f(y,x)$ は以下のカーネル密度推定量でそれぞれ推定できます.
\begin{align*}
\hat{f}(x) &= {1 \over nh^d}\sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right),\\
\hat{f}(y,x) &= {1 \over nh^d b}\sum_{i=1}^{n}W\left({Y_i - y \over b}\right)K\left({\boldsymbol{X}_i - x \over h}\right).
\end{align*}
ここで,$b = b_n \to \infty$, $n \to \infty$ は $h$ とは別のバンド幅, $W:\mathbb{R} \to \mathbb{R}$ は $d=1$ として (A2:LC) を満たすカーネル関数です.これらを用いると,$m(x)$ の自然な推定量の一つとして以下の推定量が考えられます.
\hat{m}(x) = {1 \over \hat{f}(x)}\int y\hat{f}(y,x)dy.
ここで,
\begin{align*}
\int y\hat{f}(y,x)dy &= {1 \over nh^d b}\sum_{i=1}^{n}\int yW\left({Y_i - y \over b}\right)K\left({\boldsymbol{X}_i - x \over h}\right)dy\\
&= {1 \over nh^d}\sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right)\int (Y_i - bz)W(z)dz\\
&= {1 \over nh^d}\sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right)Y_i
\end{align*}
が成り立つので,結局
\hat{m}(x) = {\sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right)Y_i \over \sum_{i=1}^{n}K\left({\boldsymbol{X}_i - x \over h}\right)}
となり,LC推定量が得られます.
LC推定量の漸近正規性
最後にLC推定量の漸近的性質について紹介しておきます.まず(A1:LC), (A2:C) に加えて以下を仮定します.
(A3:LC)
[1] $n \to \infty$ として,$h \to 0$, $nh^d \to \infty$, $nh^{d+4} \to c \in (0,\infty]$.
[2] また $x \in \mathbb{R}^d$ に対して $f(x)>0$,
[3] ある $\delta>0$ と $x$ の近傍 $N_x$ に対して $\sup_{z \in N_x}E[|U_i|^{2+\delta}|\boldsymbol{X}_i = z] \leq U(x)<\infty$ .
仮定(A1:LC)-(A3:LC)の下で,$n \to \infty$ として以下が成り立ちます.
\begin{align*}
\sqrt{nh^d} &\left(\hat{m}(x) - m(x) - b_{m,n}(x)\right) \stackrel{d}{\to} N\left(0, {\sigma^2(x) \over f(x)}\int K^2(u)du\right),\\
b_{m,n}(x) &= h^2\sum_{j,k=1}^d \left({\partial_{jk} m(x) \over 2} + {\partial_j m(x)\partial_k f(x) \over f(x)}\right)\int u_ju_k K(u)du.
\end{align*}
まとめ
この記事ではノンパラメトリック回帰分析の1つとして局所定数回帰を紹介しました.次の記事では局所定数回帰の拡張として,局所線形回帰について紹介する予定です.
株式会社Nospareには統計学の様々な分野を専門とする研究者が所属しております.統計アドバイザリーやビジネスデータ分析につきましては株式会社Nospareまでお問い合わせください.