くりーむわーかー

プログラムとか。作ってて ・試しててハマった事など。誰かのお役に立てば幸いかと。 その他、いろいろエトセトラ。。。

マイニング

R SVMのあれこれ

RでSVMの分類器作るときの段取りまとめ。kernlabとe1071を使用。

いつものアイリスデータで。。。

##################################################################################################
#アイリスデータセットを使用
d <- iris

#学習データと検証データを作る
set.seed(256)
id.random <- sample(nrow(d),nrow(d)/2)
d.train <- d[id.random,]
d.test <- d[setdiff(1:nrow(d),id.random),]

##################################################################################################
#とりあえずデフォのパラメータでSVMをやる
library( kernlab )
#学習データでSVMの分類器を作る
d.ksvm <- ksvm(Species ~.,type="C-svc", kernel="rbfdot",data=d.train,cross=5,prob.model=TRUE)

#トレーニングエラー
d.ksvm@error
#交差検証エラー
d.ksvm@cross

#検証データで精度を確認
d.pred <- predict(d.ksvm, newdata=d.test)
rt <- table(d.pred, d.test$Species)
rt

#正答率
before_tune <- as.integer(rt[1]+rt[5]+rt[9]) / as.integer(length(d.pred))
before_tune

そしたら、SVMのパラメータの値をチューニングしてく。


##################################################################################################
#チューニング1回目(適当にコストとガンマの範囲を決めて、1刻みくらいでやる)
library(e1071)
t <- tune.svm(Species ~.,data=d.train,gamma=10^(seq(-10,10,1)),cost=10^(seq(-10,10,1)),tunecontrol=tune.control(sampling="cross",cross=5))
#コストの一番いい値
t$best.parameters$cost
#ガンマの一番いい値
t$best.parameters$gamma
#一番いいパフォーマンス
100-100*tuning$best.performance
#プロット x軸:gamma y軸:cost
plot(t,transform.x=log10,transform.y=log10)

プロットすると↓の感じのグラフが出る。 sample01
そしたら、青い部分が濃い辺りの範囲を再度指定してチューニングを実行する。(ちょっと時間かかる)

##################################################################################################
#チューニング2回目
#1回目のプロットで濃い青になってる辺りを中心に今度は0.1刻みで実行
library(e1071)
t <- tune.svm(Species ~.,data=d.train,gamma=10^(seq(-4,1,0.1)),cost=10^(seq(1,6,0.1)),tunecontrol=tune.control(sampling="cross",cross=5))
#コストの一番いい値
t$best.parameters$cost
#ガンマの一番いい値
t$best.parameters$gamma
#一番いいパフォーマンス
100-100*tuning$best.performance
#プロット x軸:gamma y軸:cost
plot(t,transform.x=log10,transform.y=log10)

プロットしたら↓の感じのが出る。 sample02
もういいでしょ。というこーとで、「best.parameters$cost」の値をSVMのパラメータに入れて実行する。

##################################################################################################
#チューニングしたコストを指定して実行。kparは"automatic"にしておけば勝手にやるらしい。
library( kernlab )
#学習データでSVMの分類器を作る
d.ksvm <- ksvm(Species ~.,type="C-svc", kernel="rbfdot",C=t$best.parameters$cost, kpar="automatic",data=d.train,cross=5,prob.model=TRUE)

#トレーニングエラー
d.ksvm@error
#交差検証エラー
d.ksvm@cross

#検証データで確認
d.pred <- predict(d.ksvm, newdata=d.test)
rt <- table(d.pred, d.test$Species)
rt
after_tune <- as.integer(rt[1]+rt[5]+rt[9]) / as.integer(length(d.pred))
after_tune

自分のPCでは正答率は↓の感じになった。

チューニング前:93.3%
チューニング後:96%

精度は上がった模様。

おまけ。”probabilities”の引数渡すとそれぞれの所属確率を出してくれる。↓の感じ。

d.allprob <- predict(d.ksvm,newdata=d,type="probabilities")
head(d.allprob)
        setosa versicolor  virginica
[1,] 0.9770838 0.01094580 0.01197041
[2,] 0.9563003 0.01973656 0.02396311
[3,] 0.9708821 0.01357799 0.01553994

そしたら、元データにくっつけてCSVあたりに吐けば、個別データを吟味しやすい。

#元データに結合(列)
d.withprob <- cbind(d,d.allprob)
#CSVに出力
write.csv(d.withprob, "iris_withprob.csv")

CSVに出すと↓の感じ。 sample03

R アソシエーション分析用のデータを作るあれこれ

すぐ忘れるのでひとまずメモ。

データ分析にアソシエーション分析(バスケット分析)っていうのがある。Rを使ってやるんですが、これをやるとき、元ネタデータの整形をどうするかいっつも忘れるので特筆。

基本的には↓感じのデータを使う想定。

1課	2課	3課	4課	5課	6課	7課	8課	9課
1	1	1	0	0	0	0	0	0
0	0	0	1	1	0	0	0	0
0	0	0	0	0	1	1	1	1

で、これをRのアソシエーション分析(arules)で使う場合、最終的に↓の感じに変換する。

1課,2課,3課
4課,5課
6課,7課,8課,9課

こんな感じのデータをトランザクションデータというらしい。

Rでこれをやる場合、トランザクションデータをそのまま読む場合は↓の感じ。

hoge.tran <- read.transactions("test.csv",sep=",",format="basket")

トランザクションデータの場合はこのままaruleに食わせられる。

hoge.ar <- apriori(hoge.tran)
inspect(hoge.ar)

↓のタイプのデータを読む場合は

1課	2課	3課	4課	5課	6課	7課	8課	9課
1	1	1	0	0	0	0	0	0
0	0	0	1	1	0	0	0	0
0	0	0	0	0	1	1	1	1

R上は↓の感じで読む

#CSVから読み込み
hoge <- read.table("test.csv",header=TRUE,sep=",",row.names=1)
#行列に変換
hoge.mat <- as(hoge,"matrix")
#トランザクションに変換
hoge.tran <- as (hoge.mat,"transactions")
#あとは一緒
hoge.ar <- apriori(hoge.tran)
inspect(hoge.ar)

データフレームからは直で変換できないっぽいので行列に一回変換する。

まー結局、Rに食わせる用のデータ形式を色々準備しないとね。基本的にはアンケートの選択肢とか、商品のマスタとかマスタを全部横に並べて、0か1のON・OFFデータを作るのが一番楽かな。

数値系のデータは値の範囲でコード化してON/OFFのデータを作らないとダメかしらね。

C#でEntityFrameworkはやっぱ遅い?

EF使ってちょろちょろ作ってるけど、発行してるSQLが微妙な感じがたまにする。

特にidentity使ってるテーブルにインサートするとき、インサート文に↓の奴がくっつく。

SELECT ID
FROM [dbo].[TABLENAME]
WHERE @@ROWCOUNT > 0 AND [TABLENAME] = scope_identity()

インサートした直後のidの値をとってるっぽい。インサート後にそのオブジェクトのIDに勝手にID振ってくれてるけど、それようのクエリでしょうね。ただ、クエリとしては遅くないんだけど、プロファイラで見ると、Readsが大きいし、CPUを割とくってる。大量に発行すると致命的な気がするのですが・・・。 インサート直後にそのオブジェクトのID見たいとかあんまりないからなぁ。

んーインサート系は自分でオブジェクトから動的にクエリ作るの作らないとダメかな。

Rでワードクラウド その2

ワードクラウドを作るおおまか手順。
  1. テキストを単語にバラし、単語の出現数(率)を集計する
  2. Rで集計結果を読み込んでwordcloudでプロット
単語の出現数を求められればあとはそれをwordcloudに喰わせるだけ。

単語をバラすのにMeCabを使う。MeCabをRから呼び出すRMeCabのライブラリを使ってもできる。ただ、簡単な文章はよいのだけれども、多少量があるとなんかの制限にかかって動かない。なので、Rubyとかテキトーなスクリプトなりプログラムなりで自分でカウントする方がてっとりばやい。

私はMeCabに解析させたいテキストファイルを食わせて結果を取得⇒カウントしてCSVに保存⇒Rで読み込んでプロットしてる。あとで色々使い回したいから、MySqlにMeCabの解析結果をつっこんで、Rubyで拾って集計してる。集計部分はまた今度載せるとして、Rで読み込みプロットする部分。

↓の感じで「単語,出現数」の集計結果のCSVを準備
・・・
項,2
通り,3
否定,2
対,61
懐疑,3
ECB,17
国債,11
購入,5
・・・
そしたらRでCSVを読み込む
data1 <- read.csv("C:\\R\\wordcount.csv",header=FALSE,row.names=1)
このまま使ってもいいけど、出現数が10回以上の単語だけ出すみたいな事をしたい場合は↓の感じで絞る
data1sub <- subset(data1, data1[, 1] >= 10)
プロットする際の色を決める
Color <- brewer.pal(9, "Blues") #文字色
Color <- Color[-(1:3)] #薄い色を消す
プロット
wordcloud(row.names(data1sub), data1sub[, 1], scale=c(6,.2),random.order = T, rot.per = .15, colors = Color)
ライブラリの読み込みからまとめて書くと↓の感じ
library("RColorBrewer")
library("wordcloud")
data1 <- read.csv("C:\\R\\wordcount.csv",header=FALSE,row.names=1)
data1sub <- subset(data1, data1[, 1] >= 10)
Color <- brewer.pal(9, "Blues") #文字色
Color <- Color[-(1:3)] #薄い色を消す
wordcloud(row.names(data1sub), data1sub[, 1], scale=c(6,.2),random.order = T, rot.per = .15, colors = Color)

Rでワードクラウド作るのはいい感じなんだけど画像なのがまれに傷。FlashとかJQueryのプラグインとかでなんかかっこいいのないかなーって探してみたけど、グッとくるものがない。さらにJSONとかでデータ渡せると胸が熱くなると思う。じゃ作れとか言われそうですが。 見た感じはRのこのwordcloudっぽく。ホイールクリクリしたら中まで入っていけたり、文字をクリックでリンクしてみたりとかを夢見てみる。

Rでワードクラウド

テキストマイニングで単語の出現率を出すけど、いまいちエクセル的な一覧だとぱっとしない。なんかぱっと見、ふーんで終わってしまう。なんかいい表現方法がないものかと思ってワードクラウドに行きついた。

↓の感じ。

plot

↑はとある為替関係のサイトのアナリストさんのRSSを単語の出現率でワードクラウド表示してみたもの。

こういう見せ方をするとなーんか傾向が見えてきそうなモヤモヤ感が出てくるじゃない。。。

Rのwordcloudを使って作る。

具体例はまた今度。

MecabをRubyから呼び出すお手軽方法

Rubyから呼び出す方法はいくつかあるみたいだけど、設定とかハマりまくった覚えがあるので、お手軽版。ただし、一番遅いと思う。まぁ、膨大なデータを解析するわけじゃなければ気にならない遅さなので十分かと。Mecabをコマンドプロンプトで実行するのと同じ事をRubyで書くだけです。コマンドの文字列作ってopen3で投げる感じ。後はその結果を受け取ってるだけです。

require "open3"

mpath       = "\"C:\\Program Files\\MeCab\\bin\\mecab.exe\"" # MeCabへのパス
arg         = ""#Mecabにオプション渡したい場合
target_file = "hoge.txt"#解析対象のファイル
cmd_string = [mpath, arg, target_file].join(" ")

out, err, status = Open3.capture3(cmd_string)
out.split("\n").each  do |wkStr|
  p wkStr.split("\t")
end

↓実行結果サンプル

["この", "連体詞,*,*,*,*,*,この,コノ,コノ"]
["サイト", "名詞,一般,*,*,*,*,サイト,サイト,サイト"]
["について", "助詞,格助詞,連語,*,*,*,について,ニツイテ,ニツイテ"]
["EOS"]
["個人", "名詞,一般,*,*,*,*,個人,コジン,コジン"]
["情報", "名詞,一般,*,*,*,*,情報,ジョウホウ,ジョーホー"]
["保護", "名詞,サ変接続,*,*,*,*,保護,ホゴ,ホゴ"]

わかち書きにしたい場合はオプションに「-O wakati」をセットしてあげる。

・・・
arg         = "-O wakati"#Mecabにオプション渡したい場合
・・・

実行結果サンプル↓

["この サイト について "]
["個人 情報 保護 "]

わかち書きをさらに単語にばらしたい場合はwkStr.split(' ')みたいにして空白でsplitしてあげる。

ちなみに、参考にしたサイトではopen3じゃなくて、「IO.popen」を使ってた。↓の感じ。

io = IO.popen(cmd_string, "r")
until io.eof?
  word_list.concat io.gets.split(' ')
end

これで問題なく動くように見えるけど、数百ファイル以上を処理しようとすると、何故かたまに処理が止まる。理由は不明。ふらっと調べた感じではバグ?らしい。なので今は「open3」使えという事らしいので、素直に従って書いてます。

あとは、上記の結果の戻りをHashなり配列なりで戻すClassにしておけば、色々使い回せます。

Rでplot結果を画像に出力する方法

wikiにも載ってるけど、Rでplotを結果を画像に出力する方法。

png("C:/R/plot.png", width = 480, height = 480, pointsize = 12, bg = "white", res = NA, restoreConsole = TRUE)
plot(1)
dev.off()

dev.off()しないと出力が完了してくれない模様。

Rのスクリプトをコマンドプロンプトから実行して、plotの結果を画像として出力させたい場合に使ってる。とある統計を毎日実行でフォルダ指定して画像を出力⇒Webアプリ側から表示みたいな事に使用中。

Rのライブラリをインストールするあれこれ

最近、テキストマイニングを試してたり。

Rのライブラリのインストールで時間くったのでメモ。

Rでライブラリをインストールする場合、普通は↓。
install.packages("ggplot2")

自宅とかだとこれでいいんだけど、会社とか外部への接続が
わんさか切られてるとことかだとムリ。

プロキシの設定をちゃんと通せばいけるらしく、試してみたけどムリ。
ちなみにやってみたのは↓。

Rのコンソールを起動後に
Sys.setenv(http_proxy="http://プロキシサーバ:ポート")
認証がダメかと思い↓でやってみる
Sys.setenv(http_proxy="http://user:password@プロキシサーバ:ポート")

これも失敗。

なので、パッケージを直接手でダウンロードしてそれからインストールした。

欲しいパッケージの名前で「R CRAN パッケージ名」でググればCRANのページが出てくるので、
そこの「Windows binaries」のZipのどれかをダウンロード。たぶん「r-devel」か「r-release」あたり。
ぁ、OSはWindowsの話です。

Zipをダウンロードしたら”解凍しないで”適当なフォルダに置いて↓を実行

install.packages("[ダウンロードしたファイルのフルパス]", repos = NULL)

↓の感じ

install.packages("C:/R/ggplot2_1.0.0.zip", repos = NULL)

ただ、これだとライブラリの中でさらに使っているライブラリまではインストールできない。
なので、「library(ggplot)」みたいにすれば「このライブラリがねー」みたいなエラーが出るのでそれを頼りに
ダウンロード⇒インストール⇒エラー⇒ダウンロード⇒インストール・・・・をエラーで無くなるまで繰り返す。

めんどくさい。。。
まぁ、ネットワーク周りはハマると時間がかかるので、手っ取り早く使い始める方法。。。

あと、Windowsのアレなところで、ファイルを置いておくフォルダは下手に深い階層だったり、
日本語が含まれてると動かないまくりになるので、漢らしくCドラ直下くらいのフォルダに
置いておくのが吉みたいよ。。。

問合せ