テューキー回帰を実装して実験してみた-ぐぐりらにっき

今回は理論編ではなく実装編ということで、テューキー回帰を実装して、外れ値を含むような人工データを生成して実験してみます。

ベースラインとして、ただのlinear regressionと比較しています。

## 実装

pythonで実装しました。
Tukey回帰には、

$$
f_\theta(x)=\theta_1x+\theta_2
$$

を用いました。直線モデルです。

```python
class TukeyRegression:
    def __init__(self, eta, epsilon=1e-10, iter_num=10000, seed=None):
        self.eta = eta
        self.epsilon = epsilon
        self.seed = seed
        self.theta = np.random.randn(2)
        self.iter_num = iter_num

def fit(self, x, y):
        np.random.seed(self.seed)
        const = np.ones(len(x))
        Phi = np.stack([const, x], axis=1)

for i in range(self.iter_num):
            W = self.get_W(Phi, y)
            theta = np.linalg.solve(
                np.dot(np.dot(Phi.T, W), Phi), np.dot(np.dot(Phi.T, W), y))
            if np.linalg.norm(theta - self.theta) < self.epsilon:
                print("Converged!")
                self.theta = theta
                break
            self.theta = theta
            print(f"iter {i+1} loss {self.loss(Phi,y)}")

def predict(self, x):
        const = np.ones(len(x))
        Phi = np.stack([const, x], axis=1)
        return self.f(Phi)

def f(self, Phi):
        return np.dot(Phi, self.theta)

def W(self, r):
        w = 1 - r / self.eta**2
        w[abs(r) >= self.eta] = 0
        return np.diag(w)

def get_W(self, Phi, y):
        y_hat = self.f(Phi)
        return self.W(y_hat - y)

def loss(self, Phi, y):
        return np.linalg.norm(self.f(Phi) - y) / 2

```

ついでに、人工データの生成を行う`AnomalyDataGenerator`の実装も載せます。

```python
class AnomalyDataGenerator:
    def __init__(self, theta1, theta2, outlier_num=5, outlier=15, start=-3, stop=3, seed=None):
        self.seed = seed
        self.stop = stop
        self.start = start
        self.theta1 = theta1
        self.theta2 = theta2
        self.outlier = outlier
        self.outlier_num = outlier_num

def run(self, n):
        """データ生成.

:params n:data num
        """
        np.random.seed(self.seed)
        x = np.linspace(self.start, self.stop, n)
        epsilon = 2 * np.random.randn(n)
        y = self.theta1 * x + self.theta2 + epsilon
        for i in np.random.choice(range(n), self.outlier_num, replace=False):
            y[i] *= self.outlier
        return x, y
```

$$
y=\theta_1x+\theta_2+2\epsilon
$$
となるように生成します。
そのうち`outlier_num`個のデータをランダムに選んで`outlier`倍して、外れ値とします。

## 実験結果

まず生成した人工データを見てみます。

赤線が元となる直線で、だいたいこの直線の近くにデータが生成されています。
今回は5つ、外れ値を含むように生成しました。全部で50点生成しているので10%が外れ値、というかんじです。

ロバストでない、普通の線形回帰の結果がこちら。

すこし外れ値の影響を受けています。
ちなみに、学習の結果得られた重みは

$$
\theta_1= 2.23 ,\theta_2=-0.16 
$$
でした。真の値は$\theta_1=1.0,\theta_2=0.2$なので、結構ずれています。

次に、Tukey回帰の結果です。

人工データのグラフの直線と見比べても、結構近いです。
得られた$\theta$は

$$
\theta_1=0.96 ,\theta_2= 0.76 
$$

切片がちょっとずれていますが、$\theta_1$はかなり近づきました。
ロバストですね〜

## まとめ

前回理論編で紹介したTukey回帰を実装して、簡単に比較実験をしてみました。
最近は技術書や講義で勉強したモデルを実装してみるのが楽しいので、そういう系の記事が増えるかもしれません。

ゆくゆくは論文の再現実装とかもアップできたらいいなと思っているのですが、とりあえずはこんなかんじでやっていきたいです。

ありがとうございました。

テューキー回帰を実装して実験してみた

0 件のコメント:

コメントを投稿

このブログを検索

POPULAR POSTS

SPONSORED LINK

RECENT POSTS

TWEETS

ARCHIVES

LABEL

CONTACT