PL/CUDAでk-means法を実装する

前回のエントリでは、CBI学会で発表を行った、PL/CUDAによる類似化合物の検索について説明した。

今回は、コレとはまた別のワークロードに対する応用という事で、クラスタリング処理のIn-Database実装に挑戦してみた。
トライしてみたのは k-means法によるクラスタリング。非階層クラスタリングの領域では最も頻繁に使用される（と、言われている）アルゴリズムで、計算量もそこそこ大きい。

k-meansクラスタリングとは

教師なし学習方式の一つで、所与のデータ群を一定数(=k個)のグループに分類するためのアルゴリズムである。
以下のステップを一定回数、またはクラスタに属するデータ群に変化がなくなるまで繰り返す。

１．初期クラスタをランダムに設定する。

２．各クラスタの中心点を計算する。

３．各データ要素の属するクラスタを更新する。各データ要素はクラスタ中心点が最も近傍であるクラスタに割り当てられる。

４．新しいクラスタ定義に基づいて、各クラスタ中心点を再計算する。

５．以下繰り返し・・・。

これをPL/CUDAで実装するとどうなるか？

基本的には上記のアルゴリズムを愚直にCUDAで実装する事になるので、特別な事はしていない。
一回のPL/CUDA関数呼び出しの中で何度も繰り返し処理を行う形になるので、メインのGPU Kernelは1スレッドで起動し、Dynamic Parallelismを使って別のGPU Kernelを起動する形にする方が合理的である。

関数定義の中身に興味がある方は、こちらをどうぞ。
toybox/gpu_kmeans.sql at master · kaigai/toybox · GitHub

In-databaseでクラスタリングする

以下の図はお試し的に10000件のランダムデータを5つのクラスタに分割してみたもの。綺麗に分類されている。

ただ、ランダムデータだけだと面白くないので、来週のPGconf.SVでの発表に向け、そこそこ件数が多くてクラスタリングに適したデータは無いかと探してみたら・・・あった。

iot.ee.surrey.ac.uk

FP7ファンドの支援を受けたCityPulseというプロジェクトが収集しているデータで、デンマークのオーフス(Aarhus)市の自動車通行状況を調査したデータが公開されている。
これは、ある二地点間を通過した自動車の台数や平均速度がタイムスタンプと共に収集されているもので、2014年2月～6月の四ヵ月間のデータが全部で1350万レコード。
これをクラスタリングにかけると、どこの道路とどこの道路の通行状況が似ているのか分かるはず。

と、いう事でやってみた。

まず、PL/CUDAを用いて実装したgpu_kmeans()関数は以下のように定義されている。
第一引数にクラスタ化すべきデータとそのID値をreal[]行列として与え、第二引数にはクラスタ数を与える。

CREATE OR REPLACE FUNCTION
gpu_kmeans(real[],    -- ID + Data Matrix
           int,       -- k-value (number of clusters)
           int = 10,  -- max number of iteration
           int = 1)   -- seed of initial randomness
RETURNS int[]

これをクエリの中で使用すると以下のようになる。
少しGoogle Map Static APIを勉強して、SQLから直接URLを生成するようにしてみた。
ただ、URL長の制限が8KBという事で、線を何百本も引くような地図は作れないため、泣く泣く125件分だけに限定して地図を生成している。

WITH
summary AS (
SELECT report_id, k, c
  FROM (SELECT report_id, k, c,
               row_number() OVER (PARTITION BY report_id
                                  ORDER BY c DESC) rank
          FROM (SELECT report_id, k, count(*) c
                  FROM matrix_unnest(
                          (SELECT gpu_kmeans(array_matrix(
                                               int4_as_float4(report_id),
                                               avg_measured_time,
                                               avg_speed,
                                               vehicle_count),
                                             5)
                             FROM tr_rawdata
                          )
                       ) R(report_id int, k int)
                 GROUP BY report_id, k
               ) __summary_1
       ) __summary_2
   WHERE rank = 1
),
location AS (
SELECT point_1_lat, point_1_lng,
       point_2_lat, point_2_lng,
       CASE k WHEN 1 THEN 'red'
              WHEN 2 THEN 'blue'
              WHEN 3 THEN 'green'
              WHEN 4 THEN 'purple'
              ELSE 'orange'
       END col
  FROM summary s, tr_metadata m
 WHERE s.report_id = m.report_id
),
path_definition AS (
SELECT 'path=color:' || col || '|weight:3|' ||
       point_1_lat::text || ',' || point_1_lng::text || '|' ||
       point_2_lat::text || ',' || point_2_lng::text path_entry
  FROM location
 LIMIT 125 -- becuase of Goole Map API restriction
)
SELECT 'http://maps.google.com/maps/api/staticmap?' ||
       'zoom=11&' ||
       'size=640x480&' ||
       'scale=2&' ||
       string_agg(path_entry, '&') ||
       '&sensor=false'
  FROM path_definition;

$ wget -O map.png "`psql traffic -At -f ~/traffic.sql`"

その結果がこれ。

どうやら、青の区間は高速道路かバイパス道路で、比較的自動車の台数が多くスピードが出ているクラスタに見える。
一方、緑の区間はダウンタウンに向かう道路なので、慢性的な渋滞に悩まされているのだろうか？
赤は少し判断に迷うが、そのどちらでもない、といった所のようである。

また、これは中核となるkmeansクラスタリングをSQLの中に埋め込んでいるので、抽出条件を変えるだけで簡単に母集団を切替える事ができる。
以下の例は、個々のデータのタイムスタンプを用いて、ウィークディ（月曜～金曜）とウィークエンド（土曜、日曜）でそれぞれクラスタリング、図を作ってみたケース。

■平日

■週末

微妙な差ではあるが、週末のケースでは、ダウンタウンに向かう道路がハイウェイと同じ色に分類されている。
平日に比べると交通量が減る分、自動車がスイスイ進むという事だろうか。

閑話休題。今週のPGconf.SVでは、この辺のネタも交えてPL/CUDAメインで話をしてくる事になります。