Fire Engine

化学の修士号→消防士→ITエンジニア(2016年11月〜)

GMOペパボ株式会社を退職しました

2019年7月24日が最終出社日でした。ペパボでは、ロリポップやへテムルといったレンタルサーバサービスのインフラエンジニアとして働いており、約1年半在籍していました。
今回は、ペパボでの一年半の振り返りと、転職に至った経緯などについて書いていきたいと思います。

ペパボでの一年半の振り返り

私はペパボカレッジの枠でペパボに採用されました。ペパボカレッジとは第二新卒エンジニア向け研修のことで、中途採用でも現場に配属される前にしっかり研修を受けることができます。私はエンジニア歴も浅いし、そもそもインフラの経験がなかったため、このペパボカレッジの応募はまたとないチャンスだと感じました。実際に研修では約一ヶ月間で広大なインフラの分野の知識を爆速で学んでいきます。具体的な研修内容などは以前のブログに書きました。頑張りたい人を応援してくれる本当にいい仕組みです。

研修が終わって、配属後も周りのエンジニアの方々に多大なるサポートをいただき、未経験の私でもなんとかインフラエンジニアとして仕事をすることができました。ペパボの「いるだけで成長できる環境」とはこのことだったのかと身に染みて感じました。 どれだけ成長できたかをざっと振り返ってみると、ペパボに入った時にはサーバーに入ってもlscdくらいしか叩けずにSSHって何?といった状態だった私が、入社後半年ちょいでSSHのプロキシサーバを書くまでに至ったので、おそらくかなり成長させてもらったんだろうなーと思います。

また、働く環境も私にとっては最高でした。私のチームは完全フレックスだったのですが、名ばかりのフレックスではなく、本当にフレックスでした。私は子どもがまだ小さいため、早めに出社して早めに退社したり、趣味の筋トレをしてから遅めに出社したりと本当に自分のリズムで働かせていただきました。残業時間も月に20時間を超えたことはありませんでした。スーパーホワイトな会社で心からおすすめできる会社です。実際に私の高校時代の友人をペパボに誘って、入社しました。

転職に至った経緯

以上のように、学ぶ環境も整っていて、生活リズムに合わせて働けるとても良い環境ではありましたが、一方で、ずっと心の中にもどかしい感情がありました。

私はもともとペパボの採用面接のときから、「いずれはインフラの領域に機械学習を活用する研究がしたい」と明言していました。配属後も常にそのことを目標に仕事をしてきました。インフラエンジニアになりたての時はインフラ関連の勉強でいっぱいいっぱいでしたが、最近ではコンテナのリソース割り当ての最適化の論文を読んだり、サーバの負荷状態の把握に機械学習を用いる構想について登壇をしたりなど、やりたいことの妄想は膨らんでいました。しかしながら、なかなか妄想を実行に移すことができませんでした。それは、私自身の力不足はもちろんですが、ロリポップのように長くやっているサービスの運用は非常に大変であることも要因の一つとしてあり、結局一年半で一度も機械学習を絡めた仕事をすることができませんでした。

また仕事では、私自身が自分の強みだと思っている「深く思考すること」がなかなかできず、リアクティブに急ぎの対応するような仕事が占める割合が非常に多いのが現状でした。もちろん、自分が優秀であれば、与えられた仕事をこなしながらも、やりたいことができるかもしれませんが、私にはできませんでした。そんな時に消防士からエンジニアに夢を持って転職した自分の感情が蘇り、本当にこのままでは何者にもなれないまま終えてしまいそうにだという焦りが抑えきれなくなり、転職を決意しました。

研究者の道を選んだ

2019年8月1日からさくらインターネット研究所の研究員として働かせていただくこととなりました。
正直、このタイミングでエンジニアから研究者になることに関してはかなり悩みました。それは、私が固定観念で「エンジニアとして優秀じゃないと研究者として活躍できない」と思っていたからです。
この固定観念による不安に対する一つの答えを与えてくれたのは、@matsumotoryさんが書いた下のブログでした。研究をいずれやりたいと思っているなら今すぐ研究に飛び込んだ方がいいという考えに変わるきっかけとなりました。

hb.matsumoto-r.jp

実は松本さんには、もう一年以上前に松本さんがペパボに在籍されていた時から、エンジニアとしての自分の理想と現実が乖離していることや、いずれ研究をやりたいけどまだやっていける自信が持てないことなど色々相談させていただいていました。会議室で二時間近く二人でお話させていただいたこともありました。
その時に松本さんがおっしゃった言葉のいくつかは時間をおいてからやっとその意味が自分の中で消化できて、私の人生に大きく影響を与えた言葉となりました。松本さんはペパボ時代から私の学生時代の研究(専攻は材料工学でした)の業績をしっかり見て下さっていて「過去にこれだけ研究で活躍できていた人間が分野が違うからできないなんてことはありえない」と言ってくれていました。その言葉が最終的に私に研究に飛び込む勇気を与えてくれて、研究への挑戦を心に決めました。

そんな中、現在さくらインターネット研究所の研究員をされている@yuuk1tさんや、@kumagalliumさんとお会いする機会に恵まれ、色々と研究のお話をさせていただく中で、研究に対する熱量や考え方に共感し「私も一緒に研究をやりたい」という気持ちが高まって、さくらインターネット研究所に応募させていただきました。

さいごに

お世話になったペパボの皆さんには「辞めてすみませんでした」ではなく「教えてもらったことを活かしてこれからさらに飛躍します」とポジティブな気持ちでいます!
一年半の間本当にありがとうございました!

遺伝的アルゴリズムの並列化とgoroutineによる実装

先日、「遺伝的アルゴリズムをGoで実装してみた」という記事を書きました。
この内容で2019年7月13日(土)に開催されるGo Conference'19 in Fukuokaに登壇させていただくことになったので、開発中のeagoというパッケージをもっと作り込んで行きたいと思います。
今回は、遺伝的アルゴリズム(以下、GA)の計算処理をgoroutineでサクッと並列化した話です。

概要

GAは複数の個体を用いた多点探索のアルゴリズムであるため、本質的に並列化と親和性が高い手法であると言えます。実用的な観点からも短い計算時間で良好な解を得ることが望まれているため、これまでに数多くの並列化手法が提案されています。

今回はGo言語の強みの一つであるgoroutineを使ってGAの処理の一部を並列化してみます。また、擬似問題に対してベンチマークを取り、実装した並列処理が全体の計算時間短縮に有効であることも検証しました。

代表的な並列化手法

マスタースレーブ方式

GAを並列化するための最も単純な手法は適応度計算を並列化することです。マスタースレーブ方式では下図のように、全体の制御および遺伝的操作を行うマスターと、適応度計算を行う複数のスレーブから構成されます。

f:id:hirotsuru314:20190619141143p:plain
マスタースレーブ方式

この方式は、同時に計算できる処理を並列化しただけの極めて単純かつ一般的な手法であり、GA固有の並列化処理ではありません。
今回はこちらをgoroutineを使って実装してみます。

島モデル

島モデルはGA自体を並列化するというアプローチをとっています。集団を複数の部分集団に分割することでGAを同時に実行し、並列処理を実現しています。
その概要を下図に示します。

f:id:hirotsuru314:20190619142435p:plain
島モデル

それぞれの部分集合を初期化して、評価・選択・交叉・突然変異などをそれぞれで行います。ここで本手法の重要な特徴となるのが、ある一定の条件を元に部分集団から一部の個体を選んで交換を行う点です。このプロセスにより、島モデルは並列化手法であるという側面のほかに、個体の多様性を維持するための手法としての側面も併せもちます。

goroutineによる並列化

今回はマスタースレーブ方式を実装してみます。前述の通り、マスタースレーブ方式では適応度計算を並列化するため、適応度関数(目的関数)が非常に複雑で計算処理に時間がかかる場合に有効な並列化であると言えます。
個体群に属するそれぞれの個体の適応度計算は独立して行えるため、この適応度計算をするためのスレーブとしてgoroutineを立ち上げて、並列で計算を行うようにします。マスターについては、それぞれのスレーブの計算結果を受け取るプロセスがいるというよりは、それぞれのスレーブが以下のように定義したIndividual構造体のFitnessに計算結果を書き込むことで、以降の処理で全ての結果を参照できるようにしています。

type Individuals []Individual

type Individual struct {
    Chromosome Genome
    Fitness    float64
    Evaluated  bool
}

実装にはsync.WaitGroupを用いました。
WaitGroupを用いると、独立した複数のタスクを平行で処理して、それら全ての終了を待ち合わせる処理が非常にシンプルな記述で実装できます。

func (indis Individuals) Evaluate(parallel bool) {
    var wg sync.WaitGroup
    if !parallel {
        for i := range indis {
            indis[i].Evaluate()
        }
    } else {
        for i := range indis {
            wg.Add(1)
            go func(i int) {
                defer wg.Done()
                indis[i].Evaluate()
            }(i)
        }
        wg.Wait()
    }
}

上記の Evaluate() の部分が適応度計算の関数です。
WaitGroupはカウンターのようになっていて、Addを呼び出すと引数に渡された整数だけカウンターを増やし、Doneを呼び出すとカウンターを1つ減らします。そしてWaitを呼び出すとカウンターがゼロになるまでブロックします。

パフォーマンス検証

実行時間の測定には、testingパッケージに備わっているベンチマークの仕組みを用いました。
並列化が有効かどうかを検証するために、実用的ではないですが、適応度関数の中で1 ms sleepするようにして、擬似的に計算処理に時間がかかる適応度関数を作りました。

ベンチマークのための全体コードの以下のようになります。

package main

import (
    "github.com/tsurubee/eago"
    "log"
    "testing"
    "time"
)

type Vector []float64

func (V Vector) Initialization() eago.Genome {
    return Vector(eago.InitFloatVector(2, 32, -32))
}

func (V Vector) Fitness() float64 {
    time.Sleep(1 * time.Millisecond) // Heavy process
    var s float64
    for _, x := range V {
        s += x * x
    }
    return s
}

func (V Vector) Mutation() {
    eago.AddNormalFloat(V, 0.5)
}

func (V Vector) Crossover(X eago.Genome) eago.Genome {
    return Vector(eago.BLXalpha(V, X.(Vector), 0.3))
}

func BenchmarkParallelGA(b *testing.B) {
    for i := 0; i < b.N; i++ {
        runGA(true)
    }
}

func BenchmarkNonParallelGA(b *testing.B) {
    for i := 0; i < b.N; i++ {
        runGA(false)
    }
}

func runGA(parallel bool) {
    var v Vector
    ga := eago.NewGA(eago.GAConfig{
        PopulationSize: 30,
        NGenerations:   20,
        CrossoverRate:  0.8,
        MutationRate:   0.01,
        ParallelEval:   parallel,
    })
    ga.PrintCallBack = func() {} // Do not print messages while running
    if err := ga.Minimize(v); err != nil {
        log.Fatal(err)
    }
}

このファイルを *_test.goという形式で保存して、ファイルが存在するディレクトリで以下のようなコマンドを実行することでベンチマークが取得できます。

$ go test -bench .
goos: darwin
goarch: amd64
pkg: github.com/tsurubee/eago/examples/ga
BenchmarkParallelGA-4             50      37193820 ns/op
BenchmarkNonParallelGA-4           2     821147753 ns/op
PASS
ok      github.com/tsurubee/eago/examples/ga    4.365s

適応度関数内でのsleepの時間を1・5・10 msにして、並列化あり(parallel)並列化なし(non-parallel)における処理時間をグラフにすると以下のようになりました。

f:id:hirotsuru314:20190623114700p:plain

今回のベンチマークの条件では、個体数30、世代数20であるため、適応度計算は600回程度行われ、sleepが10 msのときに処理を並列化しない場合は、適応度計算だけでも6000 msかかることになります。
このことを考えると、適応度計算が重い場合に、今回のようなマスタースレーブ方式の並列化が、全体の処理時間の短縮に有効であることがわかります。

今後

Go Conference'19 in Fukuokaまでに以下の2点について取り組んで行きたいと思います。

  • 他の並列化手法の検証(島モデルなど)

  • 実用的な問題への適用

参考

www.morikita.co.jp

www.oreilly.co.jp

遺伝的アルゴリズムをGoで実装してみた

こんにちは!つるべーです!
最近は、進化計算と呼ばれるバイオミメティックな計算技法に興味を持っており、実装しながら勉強しています。
前回の記事では粒子群最適化(Particle Swarm Optimization: PSO)を実装しました
今回は、遺伝的アルゴリズムという生物の進化の過程を模倣して作られた最適解探索アルゴリズムをGoでスクラッチで書いてみました。

遺伝的アルゴリズム(Genetic Algorithm: GA)

GAとは何か?という問いに対する答えとして、以下の一文が端的にその特徴を表していると思います。

生物進化における遺伝と適者生存による自然淘汰の仕組みをソフトウェア的に模すことで複雑な問題に対する最適解を探索する手法
引用:「遺伝的アルゴリズム(Genetic Algorithm)を始めよう!」

GAでは、問題に対する解を個体の染色体、解の構成要素を遺伝子として表現し、複数の個体を進化させて、より環境に適応できる個体を見つけることで最適解を導いていきます。その過程をフローチャートとしてまとめたものが以下になります。

f:id:hirotsuru314:20190503224038p:plain

コード

github.com

eago(イーゴ)という進化計算パッケージを鋭意開発中です!
ちなみに、Pythonで書かれたdeapという素晴らしいライブラリがあり、遺伝的アルゴリズムを始めとした様々な進化計算の手法が実装されています。

最適化問題を解いてみる

今回はAckley functionと呼ばれる関数の最小化問題を解いてみたいと思います。数式およびグラフは以下の通りです。グラフは左側が3Dに描いたもので、右側が2次元に写像したもの(左のグラフを上から見たもの)です。

f:id:hirotsuru314:20190504145941p:plain

今回は n=2の条件を採用しています。
Ackley functionは多峰性関数であり、多くの局所解をもつため、最適化アルゴリズムベンチマーク関数として用いられることもあるようです。下のサイトに様々なベンチマーク関数が載っていて眺めるだけでも面白かったです。

Test functions for optimization - Wikipedia

グラフからわかるように、この関数は (x_1,x_2)=(0,0)で最小値 0を取るため、今回の探索の目的は、 (x_1,x_2)=(0,0)を見つけることになります。

サンプルコードは以下の通りです。

package main

import (
    "github.com/tsurubee/eago"
    "log"
    "math"
)

type Variables []float64

func (V Variables) Initialization() eago.Genome {
    return Variables(eago.InitFloatVector(2, 32, -32))
}

func (V Variables) Fitness() float64 {
    return Ackley(V)
}

func (V Variables) Mutation() {
    eago.AddNormalFloat(V, 0.5)
}

func (V Variables) Crossover(X eago.Genome) eago.Genome {
    return Variables(eago.BLXalpha(V, X.(Variables), 0.3))
}

func Ackley(X []float64) float64 {
    a, b, c, d := 20.0, 0.2, 2*math.Pi, float64(len(X))
    var s1, s2 float64
    for _, x := range X {
        s1 += x * x
        s2 += math.Cos(c * x)
    }
    return -a*math.Exp(-b*math.Sqrt(s1/d)) - math.Exp(s2/d) + a + math.Exp(1)
}

func main() {
    var v Variables
    ga := eago.NewGA(eago.GAConfig{
        PopulationSize: 30,
        NGenerations:   20,
        CrossoverRate:  0.8,
        MutationRate:   0.01,
    })
    if err := ga.Minimize(v); err != nil {
        log.Fatal(err)
    }
}

コードの詳細については最後に補足します。
上のコードをmain.goとして以下のように実行します。

$ go run main.go
Generation   1: Fitness=11.208 Solution=[-4.501 0.087]
Generation   2: Fitness=10.001 Solution=[-3.592 -0.833]
Generation   3: Fitness=9.165 Solution=[-3.841 -0.960]
Generation   4: Fitness=5.747 Solution=[-1.426 -0.235]
Generation   5: Fitness=4.258 Solution=[-0.824 0.380]
Generation   6: Fitness=2.797 Solution=[-0.843 -0.011]
Generation   7: Fitness=2.605 Solution=[-0.955 -0.030]
Generation   8: Fitness=2.605 Solution=[-0.955 -0.030]
Generation   9: Fitness=2.604 Solution=[-0.925 0.016]
Generation  10: Fitness=2.604 Solution=[-0.925 0.016]
Generation  11: Fitness=2.479 Solution=[-0.346 -0.054]
Generation  12: Fitness=1.216 Solution=[-0.184 -0.006]
Generation  13: Fitness=0.615 Solution=[-0.111 -0.001]
Generation  14: Fitness=0.501 Solution=[-0.096 -0.007]
Generation  15: Fitness=0.267 Solution=[-0.055 0.024]
Generation  16: Fitness=0.056 Solution=[0.010 -0.014]
Generation  17: Fitness=0.056 Solution=[0.010 -0.014]
Generation  18: Fitness=0.053 Solution=[-0.005 -0.015]
Generation  19: Fitness=0.031 Solution=[0.010 0.003]
Generation  20: Fitness=0.013 Solution=[-0.000 0.004]

出力結果を見てみると、世代数が大きくなるにつれて、Fitnessが 0、パラメータの解が (x_1,x_2)=(0,0)に収束しようとしているため、うまく探索できていることがわかります。

最後に、今回用いたコードとアルゴリズムの補足をしておきます。

実装とアルゴリズムの補足

Genomeインターフェース

実装では、下記のようなGenomeインターフェースを定義しており、上記のサンプルコードのVariablesはGenomeインターフェースを満たすように各メソッドを実装しています。

type Genome interface {
    Initialization() Genome
    Fitness() float64
    Mutation()
    Crossover(Genome) Genome
}

このように実装した理由は、eagoを使う側で自由に独自に実装した選択や交叉などの手法を組み込めるようにするためです。できればメジャーな手法については全てeago側に実装して、特殊なケース以外は全てeagoに実装されたメソッドを呼び出すだけで事足りるようにしていきたいと思っています。

実数値遺伝的アルゴリズム(Real coded GA)

今回の最適化では、解表現に浮動小数点数を用いており、実数値GAと呼ばれる手法を実装しました。
GAにおいては、各遺伝子を0か1で表現したバイナリエンコーディングが最も一般的な解表現だと思います。このような解の表現の仕方により、探索過程の各ステップ(選択・交叉・突然変異)で用いる手法も大きく異なってきます。
それぞれのステップにおいて、今回用いた手法を説明します。

選択

今回用いている選択の手法は、トーナメント法と呼ばれるものです。
任意のn個の個体を組にして、それらの中で最も適応度が高いものを選択するという手法になります。
他の手法として、ルーレット法やランキング法などがあります。

交叉

交叉にはブレンド交叉(BLX-α)という実数値GA特有の手法を用いました。
BLX-αについては下のサイトがわかりやすく解説されていました。

実数型GAに於ける交叉法の改良

親となる2個体によってできる長方形を一定数倍引き伸ばした長方形で囲まれる範囲内でランダムに次世代の個体が生成されます。

突然変異

突然変異として、個体の各遺伝子に任意の確率で一様乱数を加算するようなものを実装しました。一般的に突然変異は、局所解に陥ることを防ぐために取り入れられます。

さいごに

GAは、解の表現の仕方や、選択・交叉・突然変異など各遺伝子操作に様々な手法が提案されていて、何を使うかは問題を解く側に委ねられており、以前実装したPSOと比べて、非常に自由度が高いという印象を受けました。なので、GAの理論や実装をまだまだ深掘りしながら楽しめそうです!
今後は、処理の並列化や多目的最適化への拡張などを中心に勉強していきたいと思っています。

粒子群最適化(Particle Swarm Optimization: PSO)をGoで実装してみた

粒子群最適化とは、群知能による最適化手法の一種です。この手のバイオミメティクス(生物模倣)によるアプローチは、私の学生時代の専門である材料工学でも非常に盛んに研究されていましたが、データサイエンスでも応用されているのを知って興味を持ったので、実装しながら勉強しました。

粒子群最適化(Particle Swarm Optimazation: PSO

粒子群最適化とは、鳥や魚などの群れに見られる社会的行動のシミュレーションを基にモデル化されたヒューリスティックな最適解探索アルゴリズムです。

全体的な最適化の流れを説明します。 まず最初に探索空間内で、決められた個数の粒子を全て初期化します。初期化とは空間内の位置をランダムに決めることを指します。

それぞれの粒子は、目的関数により適合度(どれだけその解が優れているか)を算出できます。全ての粒子をある移動ルールを元に動かしながら、適合度が高い位置を探していくのがPSOの探索方法です。

ここでPSOにおいて重要な移動ルールについて説明します。

粒子 P_{i}の時刻 tから t+1の移動速度 \vec{v}_{i}(t+1)は以下の式で表されます。

f:id:hirotsuru314:20190409205526p:plain:w500

ここで重要なのが、グローバルベスト \vec{g}(t)とパーソナルベスト \vec{p}_{i}(t)です。
グローバルベストとは全ての粒子群で過去の探索も含めて最も適合度が高い位置、パーソナルベストとは特定の1粒子におけるこれまでの探索の中で最も適合度が高い位置を示します。これは繰り返しの探索の中で随時更新されていく変数です。
 Iは慣性係数と呼ばれる定数で、前の移動速度がどれだけ維持されるかを表します。 A_{g}および A_{p}は加速係数と呼ばれる定数です。

この数式からわかるように粒子群の中ではグローバルベストの情報が共有されていて、その方向に向かうように移動速度を調整しながら、より良い解の探索を繰り返します。これはあたかも餌を見つけた鳥が仲間にその情報を伝えて、仲間が餌のある位置に集まってくるようです。

作ったもの

github.com

『Evolutionary Algorithm implemented in Go』の略で名前を付けました。 進化計算全般に興味があるので、ここに実装しながら勉強したいと思っています。

最適化問題を解いてみる

今回は下のシンプルな目的関数の最小値を最適解として、PSOでパラメータを探索してみました。

f:id:hirotsuru314:20190409212421p:plain:w120

 N=2としたときの関数のグラフを描くと以下のようになります。

f:id:hirotsuru314:20190409215708p:plain

左側が3Dに描いたもので、右側が2次元に写像したもの(左のグラフを上から見たもの)です。
グラフや数式から直感的にわかるように、この関数は (x_1,x_2)=(0,0)で最小値 0を取るため、今回の探索の目的は、 (x_1,x_2)=(0,0)を見つけることになります。

サンプルコードは以下の通りです。

package main

import (
    "github.com/tsurubee/eago"
    "log"
)

func objectiveFunc(x []float64) float64 {
    return x[0] * x[0] + x[1] * x[1]
}

func main() {
    pso := eago.NewDefaultPSO()
    pso.NParticle =  5
    pso.NStep = 20
    pso.Min = -20
    pso.Max = 10

    if err := pso.Minimize(objectiveFunc, 2); err != nil {
        log.Fatal(err)
    }
}

上のコードをmain.goとして以下のように実行します。

$ go run main.go
Step   0: Fitness=8.124 Position=[-2.168 1.851]
Step   1: Fitness=8.124 Position=[-2.168 1.851]
Step   2: Fitness=7.394 Position=[-1.618 2.186]
Step   3: Fitness=6.958 Position=[1.068 2.412]
Step   4: Fitness=6.958 Position=[1.068 2.412]
Step   5: Fitness=4.809 Position=[-2.106 0.611]
Step   6: Fitness=2.215 Position=[-1.052 1.053]
Step   7: Fitness=1.535 Position=[0.580 1.095]
Step   8: Fitness=0.163 Position=[0.023 0.403]
Step   9: Fitness=0.163 Position=[0.023 0.403]
Step  10: Fitness=0.055 Position=[-0.175 0.156]
Step  11: Fitness=0.055 Position=[-0.175 0.156]
Step  12: Fitness=0.055 Position=[-0.175 0.156]
Step  13: Fitness=0.055 Position=[-0.175 0.156]
Step  14: Fitness=0.055 Position=[-0.175 0.156]
Step  15: Fitness=0.055 Position=[-0.175 0.156]
Step  16: Fitness=0.055 Position=[-0.180 0.151]
Step  17: Fitness=0.044 Position=[0.019 0.210]
Step  18: Fitness=0.013 Position=[-0.111 0.023]
Step  19: Fitness=0.008 Position=[-0.083 -0.031]

出力結果を見てみると、Stepが進むにつれて、Fitnessが 0、パラメータの解が (x_1,x_2)=(0,0)に収束しようとしていることがわかります。
粒子が速度を調整しながら最適解に収束して行く様子は下のサイトのアニメーションが非常にわかりやすかったです。

Particle swarm optimization - Wikiwand

さいごに

群知能とか生物模倣的なアプローチはもう少し勉強を続けて深掘りしていきたい。
あと、もっと実践的なデータに対して適用してみて、自分の今の領域における応用についても考えていきたい。

参考

shop.ohmsha.co.jp

遺伝的アルゴリズムを用いてコンテナの配置を最適化する論文:「Genetic Algorithm for Multi-Objective Optimization of Container Allocation in Cloud Architecture」

ふとタイトルが目に入ってきて気になって読んでみた。私自身そんなにコンテナ技術を触ってないけど、前から遺伝的アルゴリズムが気になっていた。 ちょいちょい知識が足りずに理解しきれないところ出てきたけど、とりあえず最後まで読んだのでざっくりまとめを書いてみる。

読んだ論文

pdfはダウンロードできなかったけど、ここからOnlineで読めた。

所感

先に論文を読んだ所感を書いておく。

  • コンテナのリソース割り当ての最適化問題遺伝的アルゴリズム(GA)を使って解くというコンセプトが面白かった。Related Workとして多くの論文が挙げられていた(Table 1およびTable2)ので、この辺りは追っていくと面白そう。

  • GAはパラメータの全検索が不可能と思われるくらい膨大な組み合わせがあって、かつ有効な探索方法が定かでないパターンに有効っぽい。さらに、論文内で用いているNon-dominated Sorting Genetic Algorithm Ⅱ(NSGA-Ⅱ)などの手法は、多目的最適化に有効な手法だから、けっこう応用の幅が広そうだなぁと思った。

  • いかに最適化したいものを目的関数として定式化するかっていうところが大変で重要なところなんだろうけど、それができるならば、目的は複数あっても、目的同士が相関関係にあっても、近似解をヒューリスティックに探索してくれるから、GAという手法は学んでいくといろんな問題に適用できそう。

  • 遺伝的アルゴリズムを始めとした進化的計算という分野に興味が湧いたのでやっていきたい。

概要

クラスター内の物理ノード群に対するコンテナの配置は、システム全体のパフォーマンスや信頼性、スケーラビリティに大きな影響を与えるため、システム設計において重要な要素である。Kubernetesのような既存のコンテナクラスター管理は、コンテナのオートスケールや物理マシンへのコンテナの配置に対して単純な実装しかされておらず、それらは物理的なリソースの使用率や閾値にのみ焦点を当てている。

このようなリソース最適化はNP完全問題であり、メタヒューリスティックなアプローチが有効である。中でも、GAのような進化的アプローチはリソース最適化問題に非常に有効な解法の一つである。本論文では、コンテナのリソース割り当てに対してNon-dominated Sorting Genetic Algorithm Ⅱ(NSGA-Ⅱ)というアルゴリズムで最適化している。NSGA-Ⅱという手法はMulti-Objective(多目的)な最適化に有効であるため、今回のように最適化したい目的関数が複数あって、それぞれの目的関数同士が独立でない場合でも適用可能である。

提案手法を用いて最適化したコンテナ配置は、Kubernetesによるコンテナ管理のアプローチよりも、論文中に定義した目的(後述)において、優位な結果が出ることが実験的に示されている。

遺伝的アルゴリズム(GA)

GAとは何か?という問いに対する答えとして「遺伝的アルゴリズム(Genetic Algorithm)を始めよう!」の7枚目のスライドがわかりやすかった。

生物進化における遺伝と適者生存による自然淘汰の仕組みをソフトウェア的に模すことで複雑な問題に対する最適解を探索する手法

私の印象としては、まじめに解こうとしても解けないくらい解の組み合わせの数が膨大なときに、実用的な時間スケールで近似的な解を探索できる点が大きな特徴だと思った。
また、実際に論文内の最適化で用いているNSGA-Ⅱについてもいくつか解析記事があった。

http://mikilab.doshisha.ac.jp/dia/research/mop_ga/moga/3/3-5-5.html

この辺りはこれから学んで数学的な背景からしっかり押さえていくので、また詳しくブログに書こうと思う。

システムのモデル化

最適化問題として解くためには、当然ながら何が最適化を評価できなければならない。さらに出てきた解(GA的には染色体)の候補同士に優劣をつけるために、定量的に評価できる形、つまり定式化しておく必要がある。

本論文ではコンテナのリソース割り当てに対して、以下の4つの目的を置いて定式化している。式中の変数についてはTable3を、定式化の詳細については「3.1 Optimization Objectives」を参照。

これらの定式化したものは、個体(解の候補)の適応度(fitness)を定量的に評価するために用いられる。

1. ボトルネックを回避するためにクラスタ全体でワークロードを均等に分散させる(Threshold Distance)

f:id:hirotsuru314:20190328222537p:plain

2. 新しいアプリケーションの追加を考慮した上で、クラスタのリソースをバランスよく利用する(Balanced Cluster)

f:id:hirotsuru314:20190328222715p:plain

3. クラスタ内のノードに適切にコンテナを分散させることによってマイクロサービスの信頼性を高める(System Failure)

f:id:hirotsuru314:20190328222655p:plain

4. 関連するマイクロサービスのコンテナを短いネットワーク距離の物理マシンに配置することによって、マイクロサービス間の通信のオーバーヘッドを小さくする(Network Distance)

f:id:hirotsuru314:20190328222846p:plain

実験と結果

実験については、GAで最適化したコンテナの配置とKubernetesアロケーションポリシーによるアプローチのそれぞれをサンプルアプリケーションを使ってベンチマークを取って比較している。

比較対象のKubernetesについては以下のように説明されている(私はKubernetesを触ったことないのでピンときていない)
KubernetesはコンテナやPodの配置を管理するために2つのポリシーを採用している。

  1. PodFitResources:物理マシンへのリソース要求の合計がマシンのキャパシティを超えないようにする

  2. LeastRequestPriority:コンテナは一番リソース消費が少ないマシンに配置される

Kubernetesの実装についてはよく知らないが、これを見るだけでも、論文中の目的関数4つの方が多くのことを考慮して最適化しようとしているのがわかる。
論文中では、GAの世代ごとの成長過程についても色々と議論されている。結論としてはGenerations 100、Population size 200というリーズナブルな値で最適化された解が得られたと述べられている。

Kubenetesとの比較の結果については、前述の4つの目的関数全てにおいて、提案手法のアプローチの方が優位な結果が出たと述べられている。例えば、一つの目的関数をピックアップすると以下のような感じ

f:id:hirotsuru314:20190328225640p:plain

変数として、物理マシンの台数を250・300・350・400台と降っている。こんな規模で実験やるんだなぁ。

さいごに

GAのこともKubernetesのこともよく知らないから、へーって感じだったけど、こういったインフラ周りで機械学習を適用している事例は面白いなと思った。GAについては数学的な理論からしっかり学んでいく。

Dynamic Time Warping(動的時間伸縮法)で時系列データをクラスタリングする

最近時系列データのクラスタリングに興味を持ち始めて、いくつか論文読んだり、アルゴリズムについて調べていたら、実装してみたくなったので勉強のために作ってみました。
実装の言語にはGolangを用いていて、クラスタリングアルゴリズムは、Dynamic Time Warping(以下、DTW)とk-medoids法を組み合わせたものです。

作ったもの

github.com

tsclusterはtime series clusteringの略です。 今回は作ったのは、特定のアルゴリズムのみですが、今後興味があるアルゴリズムがあれば、ここに実装していきます。

使い方

func main() {
    var dataset [][]float64
    dataset = append(dataset, []float64{1, 1, 1, 1, 1, 1})
    dataset = append(dataset, []float64{1, 1, 1, 1, 1, 1, 1})
    dataset = append(dataset, []float64{2, 2, 2, 2, 2})
    dataset = append(dataset, []float64{2, 2, 2, 2, 2, 2, 2})

    tc := tscluster.NewTscluster(tscluster.DTW) // 引数:距離関数
    labels, err := tc.Kmedoids(dataset, 2, 20)  // 引数:データ, クラスタ数, 最大反復回数
    if err != nil {
    log.Fatal(err)
    }
    fmt.Println(labels)
}

#=>
[0 0 2 2]

いくつか補足すると、与えるデータセットの型は[][]float64であり、1次元時系列データ([]float64)のスライスです。
また、出力結果labelsは、与えたデータセットのスライスに対してクラスタリングした結果のラベルを返しています。ラベルの値として、データが所属するクラスタのmedoids(後述するがクラスタの代表点のようなもの)のindex番号を返しています。
すなわち、[0 0 2 2] とは与えたデータセットの0番目と1番目は同じクラスタに属しており、datasetの0番目が代表点であるということです。

検証

それでは、もう少し実践的なデータで試してみます。
用いたデータとプログラムは、tsclusterのリポジトリexamplesディレクトにおいています。
1次元の時系列データを30個用意してグラフを描いたものが、以下になります。このデータは意図的に、上昇トレンド・停滞トレンド・下降トレンドそれぞれを持つものを10個ずつプログラムで用意しました。

f:id:hirotsuru314:20190217160343p:plain

これらのデータセットをDTW + k-medoidsでクラスタリングして、クラスタごとに色分けして描いたグラフが以下になります。

f:id:hirotsuru314:20190217160405p:plain

グラフを見るとそれぞれのトレンドごとにきれいにクラスタリングできているように見えます。
今回のような単純な線形のトレンドごとにクラスターに分けるにはもっと単純な手法もありそうです。本来DTWという距離尺度を用いていると、形状は似ているけど時間軸に対してずれているといったデータの類似性が高くなるような特徴を持つため、そういったデータをクラスタリングしてみるのも良いかもしれません。

実装したアルゴリズム

時系列クラスタリングの概要を把握するには、下のサーベイ論文が非常に情報がまとまっていました。

Time-series clustering – A decade review

Dynamic Time Warping(DTW)

DTWの特徴としては、

  • 時系列のデータ長がそろっていなくても計算可能
  • 時間・位相のずれていても形状が似ていれば類似度が高くなる

などが挙げられます。
DTWについてはこちらのブログに直感的に大変わかりやすいアニメーションがあったので、見てみるとイメージが湧きやすいと思います。

データ長がmとnの時系列データのDTW距離を求めるアルゴリズムは以下のステップで構成されます。

  1. 各データ点に対して総当たりで距離を計算し、要素がデータ点同士の距離となる行列D(m × n)を作成
  2. D[0][0] からD[m][n]までに通るパスのうち、行列の要素の合計が最小となるパスを探索する。その際のパスの要素の合計がDTW距離となる

k-medoids

k-medoidsはクラスタリングの手法で大変メジャーなk-meansと非常に似た手法です。
両者の大きな違いはクラスタの代表点をcentroid(重心)ではなくmedoidを選択するという点です。
medoidとは、クラスタ内のデータ点で、同一クラスタ内の他の全ての点までの距離の総和が最小になる点のことを言います。
Wikipediaには、k-meansと比較したk-medoidsの特徴が以下のように表現されていました。

It is more robust to noise and outliers as compared to k-means because it minimizes a sum of pairwise dissimilarities instead of a sum of squared Euclidean distances. (https://en.wikipedia.org/wiki/K-medoids より引用)

k-meansよりノイズや外れ値に強いとのことでした。また、k-medoidsはクラスタの代表点を実際のデータの中から選ぶため、一度データ同士の距離を全て計算してしまえば(距離行列を作成すれば)、k-meansのような重心の計算や重心からのデータ点までの距離の計算などといった処理がないため、計算量も少なくなると思います。

アルゴリズムとしては以下の通りです。

  1. データセットの中からk個のデータをランダムに選択する(modoid)
  2. 全てのデータを一番近いmedoidのクラスタに割り当てる
  3. それぞれのクラスタ内において、クラスタ内の他のすべての点との距離の合計が最小になる点を新たにmedoidとする
  4. 2と3を繰り返す。ただし、3においてmedoidに変化がない、もしくは指定した最大の反復回数に到達すれば終了する

さいごに

今回のようなデータのクラスタリング教師なし学習で、データがあればすぐに適用できるというメリットはあるものの、当然ながら出てきたクラスター群に対する解釈は人間が与えないといけないですし、そもそも何個のクラスターに分けるかはある程度事前知識が必要だと思います(クラスタ数を最適化する方法はあるものの)
また何か気になるものがあれば実装してみたいと思います!

k-Shapeによる時系列クラスタリングの論文:「k-Shape: Efficient and Accurate Clustering of Time Series」を読んだ

最近、時系列データのクラスタリングに興味を持っているので、k-Shapeというクラスタリング手法に関する論文を読みました。
なぜ興味を持っているかというと、サーバの各種メトリクス(CPU使用率・メモリ使用率など)を使って、似たような特徴を持っているサーバ群をクラスタリングできないかと考えいるためです。例えば、負荷の高いサーバ群と負荷の低いサーバ群などにグルーピングできると、面白いのではないかと考えています。

元論文はこちらです。

k-Shape: Efficient and Accurate Clustering of Time Series

この論文では、提案するk-Shapeのアルゴリズムの説明だけでなく、時系列データをクラスタリングする上での理論的な背景から説明されていて、時系列クラスタリング手法の全体像を把握するのに大変良い論文でした。

概要

本論文では、k-Shapeという時系列クラスタリング手法を提案しています。

k-Shapeの特徴は、

  • 時系列データの形状に着目したshape-basedクラスタリングである
  • データ間の距離尺度として、規格化した相互相関を用いている
  • 高効率かつ高精度で、幅広いデータに適用できる

などが挙げられます。
以下、本論文の内容をまとめていきますが、理論の詳細に関しては、論文をご参照ください。

Time-Series Invariances

時系列データをクラスタリングする際には、何に着目してクラスタリングするか(何をもって似ているとするのか)というのが重要になります。
例えば、形状は似ているけど、時間軸に対してずれている2つの時系列を同一クラスと捉えるかどうかは、そのときの用途によって異なると思います。

k-Shapeの場合、ScalingShiftingに対する不変性 (Invariances)に着目しています。
簡単に言うと、Scalingは軸に対してスケールした際にデータ同士の性質が似ているかどうか、Shiftingは時間軸(位相)に対してずらした場合にデータ同士の性質が似ているかどうかに着目しており、この特徴からk-Shapeはデータの形状に着目したshape-basedなクラスタリングであるとされています。

Time-Series Distance Measures

クラスタリングを行う際に、データ間の類似性は距離として表現されます。すなわち、距離が近いほど類似性が高く、距離が遠いほど類似性は低いということです。

このことから、前述の「何を持って似ているとするか」は、データ間の距離をどのように定義するかに強く依存するため、距離尺度の定義が極めて重要となってきます。論文中でも距離をどうやって決めるかが、クラスタリングアルゴリズム自体より重要だと述べらていました。

時系列データの「形状」に着目したクラスタリングをするには、振幅と位相のずれをうまく扱う距離尺度が必要となります。
k-Shapeでは、一般的によく用いられるユークリッド距離や動的時間伸縮法ではなく、独自のShape-based distance(SBD)という距離尺度を用いています。

Shape-based distance(SBD)

SBDは、規格化した相互相関を用いています。
相互相関とは、信号処理や画像処理に広く使われている信号間の類似性の尺度です。
ざっくりいうとデータをずらしながら内積を計算していき、内積が一番大きくなった位置を探す感じです。

2つの時系列データxyにおけるSBDは以下の式で表されます。

f:id:hirotsuru314:20190206220202p:plain

ただし、

f:id:hirotsuru314:20190206220243p:plain

Shape-based Time-Series Clustering

距離尺度が決まってしまえば、クラスタの重心が計算できます。(重心の計算方法は論文中の「3.2 Time-Series Shape Extraction」を参照)
あとのクラスタリングアルゴリズムはk-meansとほぼ同じような感じで、以下のように反復的な処理を行います。

  1. 時系列データを各クラスタの重心と比較し、最も近い重心をもつクラスタにデータを割り当てる
  2. クラスタに属するメンバーから重心を計算して更新する

これを繰り返したのち、クラスタのメンバーに変更がないか、反復回数の最大値に達するまで上の手順を繰り返します。

感想

手持ちのデータをクラスタリングしたいときには、どの不変性に着目して、何を似ているとしたいのかを考えて、それに適したアルゴリズムを選定するかが重要だと思いました。
今回のアルゴリズムは、時間軸のずれや、伸縮に対して頑強なアルゴリズムである印象を受けました。サーバメトリクスの場合、scalingは重要な違いを表す要素だと思うので、サーバメトリクスのクラスタリングに対してはもっと適したアルゴリズムがあるかもなーと思いました。もっといろんな手法をサーベイしていきます。

k-Shapeに関しては、論文中に擬似コードが載っていたのでスクラッチで書けそうだし、tslearnというPythonライブラリもあるようなので試したい。