【機械学習】Goでkmeans実装

Goで機械学習フルスクラッチシリーズ

* [【機械学習】Goでパーセプトロン実装](http://qiita.com/shu_marubo/items/5081ce10c5490eaac207)
* [【機械学習】GoでKNN(K近傍法)実装](http://qiita.com/shu_marubo/items/25962cfabebf483af4ee)


この記事は、自身のブログ、[Data Science Struggle](http://marubon-ds.blogspot.jp)でも掲載予定。許可なき掲載とかではない。

#概要
機械学習アルゴリズムのうちの一つ、kmeansをGoで書く。最近Goを使うようになったので、コーディング練習がてら、機械学習アルゴリズムを書いていく。 
Goでも機械学習用のパッケージはあるが、それぞれの機械学習アルゴリズムの書き方に焦点を当てた情報はなかったので少しづつ公開していく。 
アルゴリズム理解のための書き下しなので、導入することで見通しの悪くなるような精度向上要素はいれず、多少の問題が生じることは気にせずにシンプルに書いていく。 

#kmeansとは？
kmeansとは教師なし機械学習アルゴリズムの一つで、データのクラスタリングを行う。クラスタの個数を指定することができるため、個人的にはクラスタリング系のアルゴリズムの中では使い勝手が良いように感じる。 

#アルゴリズム
kmeansのアルゴリズムは、適当な代表点をk個定めて、データを最も距離的に近い代表点に帰属させる。それぞれの代表点について、その代表点に属するデータの平均求め、それをその代表点の更新先とし、再びデータを最も距離的に近い代表点に帰属させる。これを状態の変化がなくなるまで繰り返す。 
まとめると以下ステップを踏むことになる。 

適当な代表点を定めた後
1. データを代表点に帰属させる
2. 帰属させたデータのセントロイドを新たな代表点とする
3. 上記の二つを、状態の変化がなくなるまで繰り返す

#今回書いたアルゴリズム
今回書いたアルゴリズムでは、具体的に以下の方法を取っている。 

* 代表点の初期化は、範囲を定めた乱数で行なっている
* 使用している距離関数はユークリッド距離
* 打ち切りは『データのラベルが更新により変化しなかった』ときに行う

実際には、代表点の初期化は、データをランダムにk個のクラスタに分けてそのセントロイドを取るのが一般的であり、打ち切りも様々な条件をつけることができるが、今回はなんとなく上記の方法を採択した。特に理由はない。（書いて見てわかったが問題はあったが、それはラストで。）

#使用コード

```go

package main

import (
	"os"
	"encoding/csv"
	"io"
	"strconv"
	"math"
	"math/rand"
	"reflect"
	"fmt"
	"time"
)

func main() {
	//データ読み込み
	irisMatrix := [][]string{}
	iris, err := os.Open("iris.csv")
	if err != nil {
		panic(err)
	}
	defer iris.Close()

	reader := csv.NewReader(iris)
	reader.Comma = ','
	reader.LazyQuotes = true
	for {
		record, err := reader.Read()
		if err == io.EOF {
			break
		} else if err != nil {
			panic(err)
		}
		irisMatrix = append(irisMatrix, record)

	}

	X := [][]float64{}
	Y := []string{}
	for _, data := range irisMatrix {
		//strスライスデータをfloatスライスデータに変換
		temp := []float64{}
		for _, i := range data[:4] {
			parsedValue, err := strconv.ParseFloat(i, 64)
			if err != nil {
				panic(err)
			}
			temp = append(temp, parsedValue)
		}
		//説明変数へ
		X = append(X, temp)
		//被説明変数
		Y = append(Y, data[4])
	}
	km := Kmeans{}
	fmt.Printf("predict:%v\n", km.fit(X, 3))
	fmt.Printf("teacher:%v\n", Y)
}

type Kmeans struct {
	data            [][]float64
	labels          []int
	representatives [][]float64
}

func Transpose(source [][]float64) [][]float64 {
	var dest [][]float64
	for i := 0; i < len(source[0]); i++ {
		var temp []float64
		for j := 0; j < len(source); j++ {
			temp = append(temp, 0.0)
		}
		dest = append(dest, temp)
	}

	for i := 0; i < len(source); i++ {
		for j := 0; j < len(source[0]); j++ {
			dest[j][i] = source[i][j]
		}
	}
	return dest
}

//二点の距離を計算。ユークリッド距離。
func Dist(source, dest []float64) float64 {
	var dist float64
	for i, _ := range source {
		dist += math.Pow(source[i]-dest[i], 2)
	}
	return math.Sqrt(dist)
}

//スライスに格納されている値が最も小さくなるときにそのインデックスを返す
func ArgMin(target []float64) int {
	var (
		index int
		base  float64
	)
	for i, d := range target {
		if i == 0 {
			index = i
			base = d
		} else {
			if d < base {
				index = i
				base = d
			}
		}

	}
	return index
}

func (km *Kmeans) fit(X [][]float64, k int) []int {
	//乱数のタネ
    rand.Seed(time.Now().UnixNano())
	//データを構造体に格納
	km.data = X

	//代表ベクトルの初期化
	//乱数で与える初期値の範囲を定めるために、データを特徴量ごとに走査して、各特徴量の最小値と最大値を獲得する
	transposedData := Transpose(km.data)
	var minMax [][]float64
	for _, d := range transposedData {
		var (
			min float64
			max float64
		)
		for _, n := range d {
			min = math.Min(min, n)
			max = math.Max(max, n)
		}
		minMax = append(minMax, []float64{min, max})
	}
	//各特徴量の最小値から最大値までの範囲に収まる乱数を生成し初期値とする
	for i := 0; i < k; i++ {
        km.representatives = append(km.representatives, make([]float64, len(minMax)))
	}
	for i := 0; i < k; i++ {
		for j := 0; j < len(minMax); j++ {
			km.representatives[i][j] = rand.Float64()*(minMax[j][1]-minMax[j][0]) + minMax[j][0]
		}
	}
	//ラベルの初期化
	//各データ点と代表ベクトルとの距離を計算し、最も近い代表ベクトルのインデックスをデータ点のラベルとする
	for _, d := range km.data {
		var distance []float64
		for _, r := range km.representatives {
			distance = append(distance, Dist(d, r))
		}
		km.labels = append(km.labels, ArgMin(distance))
	}
	for {
		//代表ベクトルの更新
		//ラベルに属するデータ点の平均値を代表ベクトルの更新値とする
		//インデックスiはラベルを表す
		var tempRepresentatives [][]float64
		for i, _ := range km.representatives {
			var grouped [][]float64
			for j, d := range km.data {
				if km.labels[j] == i {
					grouped = append(grouped, d)
				}
			}
			if len(grouped) != 0 {

				transposedGroup := Transpose(grouped)
				updated := []float64{}
				for _, vectors := range transposedGroup {

					value := 0.0
					for _, v := range vectors {
						value += v
					}
					//特徴量ごとの平均を格納
					updated = append(updated, value/float64(len(vectors)))
				}
				tempRepresentatives = append(tempRepresentatives, updated)
			}
		}
		km.representatives = tempRepresentatives

		//ラベル更新
		tempLabel := []int{}
		for _, d := range km.data {
			var distance []float64
			for _, r := range km.representatives {
				distance = append(distance, Dist(d, r))
			}
			tempLabel = append(tempLabel, ArgMin(distance))
		}
		if reflect.DeepEqual(km.labels, tempLabel) {
			break
		} else {
			km.labels = tempLabel
		}
	}
	return km.labels
}


```

#順を追って見ていく
コードの大半は実行部かfit()メソッドとなる。fit()メソッドを分割して見ていく。 


```go

func (km *Kmeans) fit(X [][]float64, k int) []int {
	//乱数のタネ
    rand.Seed(time.Now().UnixNano())
	//データを構造体に格納
	km.data = X

```
fit()メソッドはクラスタリング対象データとクラス数を引数に取る。まず、対象データを構造体に格納する。 

```go
	//代表ベクトルの初期化
	//乱数で与える初期値の範囲を定めるために、データを特徴量ごとに走査して、各特徴量の最小値と最大値を獲得する
	transposedData := Transpose(km.data)
	var minMax [][]float64
	for _, d := range transposedData {
		var (
			min float64
			max float64
		)
		for _, n := range d {
			min = math.Min(min, n)
			max = math.Max(max, n)
		}
		minMax = append(minMax, []float64{min, max})
	}
```
kmeansでは指定したクラスの数だけ、代表点を用意し、それをデータで更新していくが、この代表点の初期値がかなり重要。今回は、データの説明変数ごとに最小値と最大値を調べ、その範囲内で乱数を取って初期値を与えている。そのため、ここでは説明変数ごとの最小値と最大値を調べ、minMaxに格納している。 

```go

	//各特徴量の最小値から最大値までの範囲に収まる乱数を生成し初期値とする
	for i := 0; i < k; i++ {
        km.representatives = append(km.representatives, make([]float64, len(minMax)))
	}
	for i := 0; i < k; i++ {
		for j := 0; j < len(minMax); j++ {
			km.representatives[i][j] = rand.Float64()*(minMax[j][1]-minMax[j][0]) + minMax[j][0]
		}
	}

```
分類対象データの持つ説明変数の数に合わせて乱数によって代表点の初期値を与えている。 

```go
	//ラベルの初期化
	//各データ点と代表ベクトルとの距離を計算し、最も近い代表ベクトルのインデックスをデータ点のラベルとする
	for _, d := range km.data {
		var distance []float64
		for _, r := range km.representatives {
			distance = append(distance, Dist(d, r))
		}
		km.labels = append(km.labels, ArgMin(distance))
	}
```
分類対象データの属する代表点を定める。データとそれぞれの代表点との距離を計算し、最も距離が小さい代表点をデータのラベルとして与えている。今回は代表点のインデックスをそのままラベルとして扱っている。 

```go

	for {
		//代表ベクトルの更新
		//ラベルに属するデータ点の平均値を代表ベクトルの更新値とする
		//インデックスiはラベルを表す
		var tempRepresentatives [][]float64
		for i, _ := range km.representatives {
			var grouped [][]float64
			for j, d := range km.data {
				if km.labels[j] == i {
					grouped = append(grouped, d)
				}
			}
			if len(grouped) != 0 {

				transposedGroup := Transpose(grouped)
				updated := []float64{}
				for _, vectors := range transposedGroup {

					value := 0.0
					for _, v := range vectors {
						value += v
					}
					//特徴量ごとの平均を格納
					updated = append(updated, value/float64(len(vectors)))
				}
				tempRepresentatives = append(tempRepresentatives, updated)
			}
		}
		km.representatives = tempRepresentatives

		//ラベル更新
		tempLabel := []int{}
		for _, d := range km.data {
			var distance []float64
			for _, r := range km.representatives {
				distance = append(distance, Dist(d, r))
			}
			tempLabel = append(tempLabel, ArgMin(distance))
		}
		if reflect.DeepEqual(km.labels, tempLabel) {
			break
		} else {
			km.labels = tempLabel
		}
	}
```

ここが学習部。前半で代表点の更新、後半でデータのラベルの更新を行なっている。代表点の更新は、その代表点のラベルに属するデータの平均点を新たな代表点とする。データラベルの更新は、更新された代表点のうち、最も近いものをそのデータの属する代表点としてそのデータのラベルとする。 
更新の打ち切りは***状態が変化しなくなった***ら行う。代表点の変化によるラベルの変化や、距離の変化などで閾値を設けることができるが、今回は簡便で、『ラベルの更新前後でラベルに変化がなくなったら打ち切り』としている。実はこれだと、初期値によっては全てのデータが一つのラベルに割り振られ、一回の更新では状態が変わらずに打ち切りがなされ、全てのデータが一つのラベルに割り振られるのが最終アウトプットになってしまったりするのであまりよくはない。実際、このコードも何回かに一回は全てのデータが一つのラベルに割り振られてしまうため、打ち切り方法、もしくは代表点の初期値の与え方は、実際に使用していくためにはもう少ししっかりと作る必要がある。このまま使用しては絶対にいけない。 
たぶんあとで一般的な方法で代表点の初期化を行なったコード載せます。