multicore (R)

Rのライブラリ multicore を試した結果です．

mclappy()という関数を使いました．

この関数は lapply() を並列でできるようにしたもので，
独立に繰り返し計算する必要がある場合（乱数や初期値を変えた結果が欲しいときなど）
に簡単に並列で計算できます．

例えば k-means 法を R で行う場合に kmeans() という関数があります．
この関数は

kmeans(x, centers, iter.max = 10, nstart = 1,
       algorithm = c("Hartigan-Wong", "Lloyd", "Forgy",
                     "MacQueen"))

という形になっていて，x としてクラスタリングしたいデータ，centers に
クラスタ数をとります．このときクラスターの初期値を nstart 回ランダムに選び，
その結果もっともクラスター内平方和（tot.withinss）を小さくした結果を返します．

k-meansの各初期値からの反復は，ほかの初期値を選んだ場合の結果に依存しないので
簡単に並列化できると考えられます．

そこで並列で計算させた場合にどれくらい早くなるか実験しました．
まず mclapply() は

mclapply(X, FUN, ..., mc.preschedule = TRUE, mc.set.seed = TRUE,
         mc.silent = FALSE, mc.cores = 1L,
         mc.cleanup = TRUE, mc.allow.recursive = TRUE)

となっています．

lapply() と同じで，X の要素に関数 FUN を作用させた結果をリストで返します．
mc.cores が（最大限）使うコアの数です．

これを踏まえ，ライブラリMASSにあるデータ Boston を4つのクラスタに k-means 法を使ってクラスタリングしました．
初期値は 5000*n 個 (n=1,2,...,10) 与え，並列したときとかかった時間の差を調べました．

並列させないときは

t <- proc.time()
invisible(kmeans(Boston, 4, nstart=5000*n))
duration <- (proc.time() - t)[3]

で時間を計り，並列させたときは

t <- proc.time()
results <- mclapply(rep(5000*n/8, 8), function(nstart) kmeans(Boston, 4, nstart=nstart))
#collect results
objective.values <- sapply(results, function(result) result$tot.withinss)
result.best <- results[[which.min(objective.values)]]
duration <- (proc.time() - t)[3]

で時間を計りました．
#collect results の後の二行はクラスター内平方和（tot.withinss）を小さくした結果を返すための部分です．

mclapply() は mc.cores を指定しないと，自動で使えるコアの数を調べます．
自分が実験した環境だとコアの数が8だったので，5000*n個の初期値を
8個に分割しています．分割数に差があるかどうかを調べるために
9個に分割した場合も調べました:

t <- proc.time()
results <- mclapply(rep(5000*n/9, 9), function(nstart) kmeans(Boston, 4, nstart=nstart))
objective.values <- sapply(results, function(result) result$tot.withinss)
result.best <- results[[which.min(objective.values)]]
duration <- (proc.time() - t)[3]

細かいですが kmeans() は nstart の値が小数でもエラーがでず，nstartの小数以下を切り捨てた回数だけ初期値をランダムに発生させて計算しています（1:4.5と1:4が同じため）．

結果としては次のようになりました:

横軸は初期値の数(5000*n)，縦軸はかかった時間です．
凡例の cores の値は使ったコアの数で，list の値は 5000*n 個の初期値の分割数です．

コア数が多い方が当然早く計算ができることが分かりました．
また初期値の分割数はコア数の倍数にした方がいいと思われます．
実際 multicore のマニュアルを読むと

"By default (mc.preschedule=TRUE) the input vector/list X is split into as many parts as there are cores (currently the values are spread across the cores sequentially, i.e. first value to core 1, second to core 2, ... (core + 1)-th value to core 1 etc.) and then one process is spawned to each core and the results are collected."

とあります．そこでコア数の倍数である16, 24分割した場合は8分割の場合とほとんど同じ結果になりました（凡例の一番上はcores=8の間違いです）:

以上をまとめると

・並列計算したほうが早く，multicoreのmclapply()で簡単にできる
・初期値や乱数の分割はコア数の倍数にした方が良い（一つ一つの計算量が同じなら）

snowとかも調べてみたいと思います．