Rで確率的ニューラルネットワーク
あっ, どーも僕です.
昨日Tokyo.Rがあったらしいですね. 一度は参加したり発表してみたいと思いつつも, 片道10時間の7000円なのでちょっと厳しいですね. 発表のネタとしてはこないだ自分の研究のために適合的棄却サンプリングを勉強したのでそれとかいいと思うですが.
まっ, それはおいといて今回は確率的ニューラルネットワークを簡単に試したいと思います.
Probabilistic Neural Network
といっても, 前に載せたものを知らぬ間に消していたので前載せたのを少し手直して載せます. 調べたところほかに紹介しているブログもなかったのでまぁいいのかな.
Rでニューラルネットワークといえば推奨パッケージでもあるnnet*1が有名です. nnetではシグモイド関数か線形関数を活性化関数にした最大3層のフィードフォーワード型のニューラルネットワークを提供しています. とても使い易いので, 多くのブログで取り上げられています.
そんななか, 今年に入ってから確率的ニューラルネットワークがpnn*2パッケージで使用できるようになりました. 作者の方はフランスの方のようです.
確率的ニューラルネットワークってなんだといわれると困るのですが, 簡単に説明すると中間層に活性化関数を使用しないでカーネル関数を使用するそうです. そして, 一般的なニューラルネットワークと同じように層は3層だそうです.....
まっ, 確率的ニューラルネットワークで検索するとフリーの論文がでてくるので詳しくはそちらを見てください.
数学的な部分はおいといて, リンク先を参考にirisを検証データと学習データの半々に分けてで試してみました. コードと結果を次に示します. コードですが, 予測がひとつずつしかできないのと, 結果の構造が扱いにくい形なので少し込み入った書き方をしてます.
結果うち各カテゴリーへ確率の部分をどう表すはちょっと考えましたがとりあえず棒グラフで. もっと良い手法あるかもですが, このくらいのデータ数なら棒グラフでもいいかなって. 縦軸が各カテゴリーに属するであろう確率で, 横軸がデータ番号です. やはりsetosaははっきりわかれますが, 残りは半々くらいになっているものもありますね.
# http://flow.chasset.net/category/software/r/pnn/ # pnnをirisで練習する library(pnn) # data myiris <- iris index <- seq.int(2, 150, 2) TrainDat <- myiris[ index, ] TestDat <- myiris[-index, ] # learn pnn <- learn (TrainDat, category.column = 5) pnn <- smooth (pnn, sigma = 1) # 最適化 # result result <- ldply(seq.int(nrow(TestDat)), function (i) unlist(guess(pnn, unlist(TestDat[i, -5])))) head(result) # category probabilities.setosa probabilities.versicolor probabilities.virginica # 1 setosa 0.985951140020672 0.0139619842011578 8.68757781698741e-05 # 2 setosa 0.987385666377889 0.0125660319928924 4.83016292180146e-05 # 3 setosa 0.987582582540735 0.0123435605294123 7.38569298525304e-05 # 4 setosa 0.986181870326146 0.0137534768220793 6.46528517746964e-05 # 5 setosa 0.982146716704304 0.017794089036588 5.91942591078885e-05 # 6 setosa 0.983339235205366 0.0165128482719832 0.000147916522650823 # 正当率 sum(diag(prop.table(table(result[,"category"], TestDat[,5])))) #[1] 0.9466667 # 結果のグラフ化 # とりあえず棒グラフで・・・ prob.result <- sapply (result[,-1], type.convert) # resultのクラスはみんなキャラクターなので数値に変換する col = brewer.pal(3, "Accent") def.par <- par(no.readonly = TRUE) # デフォルトの環境を保存 par(mai = c(1.5, 1.5, .5, 2.5)) tmp <- barplot (t(prob.result), space = 0, border = "lightgrey", las = 1, cex.axis = 1.5, cex.lab = 1.5, ylab = "prob", col = col, xlab = "index") abline(h = 0, col = "lightgrey") axis(side = 1, tcl = .5, cex.axis = 1.5, at = tmp, labels = 1:75) legend (x = tmp[75] + 1, y = .65 , legend = c(unique(result[,1])), xpd = TRUE, pt.cex = 3, pch = 15, title = "Species", cex = 1.5, col = col, bty = "n") par(def.par) # デフォルトの環境を復帰
sigmaの検討
学習の部分でいじれるパラメータにsigmaがあります. 詳しいことは勉強不足なのですが, ひとまずこの値をいじって上と同じシュミレーションをし正解率を比較してみました. 結果は次のとおりです. irisだと, sigmaを大きくしても正解率は0.92より小さくなりませんね. ってか縦軸名を入れ忘れていますが, 縦軸は正解率です.
# ********************************************************************* # sigmaを変えて正当率をみる sigma <- c(.1, .5, 1, 5, 10, 100) rtn <- numeric (length (sigma)) for (i in seq_along (sigma)) { pnn <- smooth (pnn, sigma = sigma[i]) result <- ldply(seq.int(nrow(TestDat)), function (i) unlist(guess(pnn, unlist(TestDat[i, -5])))) rtn [i] <- sum(diag(prop.table(table(result[,"category"], TestDat[,5])))) } #グラフ par (def.par) # デフォルトの環境を復帰 par (mai = c(1.5, 1.5, .5, .5)) plot (seq_along (sigma), rtn, type = "b", col = "#0068b7", cex = 2, lwd = 2, tcl = .5, las = 1, pch = 19, cex.axis = 1.5, cex.lab = 1.5, xlab = "sigma", ylab = "", ylim = c(.9, 1), xaxt = "n", lty = 3) axis(side = 1, tcl = .5, cex.axis = 1.5, at = seq_along (sigma), labels = sigma) def.par <- par(no.readonly = TRUE) # デフォルトの環境を保存
訂正
追記と, プログラムを一部直しました. (2013/7/13)