December 6, 2015

Japan.R 2015を開催しました

R advent Calendarの12月6日担当です。更新日が12月7日になってしまいました。(寝落ちしてた) 12月5日にJapan.R 2015を開催しました。ATNDでの登録数は197人と昨年同様の盛り上がりでした。 スタッフの皆さん、参加者のみなさん、ありがとうございます。 懇親会の参加者は約70人で、二次会の参加者は約40人でした。眠かったので帰りましたが、三次会もあったようです。 また、会場提供をして頂いた株式会社リクルートコミュニケーションズさん、ありがとうございます。 当日の様子はtogetterにまとめています。 Japan.R 2015 前半セッション Japan.R 2015 後半セッション Japan.R 2015 LT大会以降 http://www.ustream.tv/channel/kbWLTTuAtLhに当日の録画もあります。 また、僕自身もLT大会に参加しました! gepuro task views from Atsushi Hayakawa 今年もR言語にお世話になりました。来年以降もR言語コミュニティに貢献していけたらと思います。来年もJapan.Rを開催しますので、よろしくお願いします。 Enjoy R!

December 2, 2015

ioslidesにおけるロゴのサイズを変える

スライドをRmdで書くgepuroです。 デフォルトのままでは、ロゴのサイズがピクセル単位で固定されているので、画像によっては歪んで表示されてしまう。 それでは困るので、 CSSを弄りました。 また、 タイトル以外のページに現れるロゴは邪魔なので、 background-size: 0px 0px;として、消しています。 slide.Rmd title: "タイトル" author: "@gepuro" output_format: ioslides_presentation: css: logo.css logo: logo.png mystyle.css /* トップのロゴの大きさを変える */ aside.gdbar img { width: 180px; height: 47.2px; position: absolute; right: -5px; top: 15px; margin: 8px 15px; } /* 他のページのロゴの大きさを変える */ slides > slide:not(.nobackground):before { font-size: 12pt; content: ""; position: absolute; bottom: 20px; left: 60px; background: url(logo.png) no-repeat 0 50%; -webkit-background-size: 30px 30px; -moz-background-size: 30px 30px; -o-background-size: 30px 30px; background-size: 0px 0px; padding-left: 40px; height: 30px; line-height: 1. Read more

November 26, 2015

Rmarkdownのioslidesで右下に画像を配置する方法

もう少し綺麗な方法があると思うけど、下で出来た。 <style> .bottom_right img { position: relative; left: 500px; top: 180px; }</style> <div class="bottom_right" markdown="1"> <img src="logo.png", width="300px"> </div> kazutanさんに教えてもらった方法 <aside class="gdbar right bottom"><img src="hoge.png"></aside> ありがとうございます!

September 10, 2015

library(ranger)とlibrary(randomForest)の速度を比較する

Predicting Titanic deaths on Kaggle V: Rangerを読んでいたら、 「ranger早いぜ」というのを見つけたので試してみた。 library("ranger") library("randomForest") ## randomForest 4.6-10 ## Type rfNews() to see new features/changes/bug fixes. ## ## Attaching package: 'randomForest' ## ## The following object is masked from 'package:ranger': ## ## importance data("iris") # http://d.hatena.ne.jp/a_bicky/20101225/1293258748 test <- function(code, ic = 1000) { expr <- substitute(code) print(system.time(for (i in 1:ic) { eval(expr) })) } test( ranger(Sepal.Length~., data=iris, mtry=2, num.trees = 500, write.forest=TRUE) ) ## user system elapsed ## 71. Read more

August 8, 2015

RでCSVファイルを一行ずつ読み込む

メモリを節約するんじゃ 学習データはサンプリングするからメモリに載せられるけど、 予測対象はメモリに乗らない時に一行ずつ読み込んで処理をさせたい。 f <- file("data.csv", "r") while(1){ line <- readLines(f, n=1) if(length(line) == 0){ break }else{ line.f <- textConnection(line) line.df <- read.csv(line.f, header = FALSE, sep=",") } write.table(line.df, stdout(), row.names=F, col.names=F, sep=",") } 標準入力から受け取る時は、 file("data.csv", "r") を file("stdin") にすれば良さそう 参考 http://stackoverflow.com/questions/9370609/piping-stdin-to-r 追記 奥村先生版です。 「RでCSVを1行ずつ」のいいかげん版: while (length((line = readLines(f,n=1))) > 0) { x = unlist(strsplit(line, ",")) https://t.co/fTpnA5CVIr — Haruhiko Okumura (@h_okumura) August 8, 2015 追記2 hoxo_mさんバージョンです。 Read more

June 17, 2015

kmeansとkmedoidsを比較してみる

やったこと kmeansとkmedoidsの比較をしてみました。 アルゴリズムの詳細は下の参考文献を見てください。 ざっくりとした説明をすると、 kmeansはクラスタの中心を重心として距離を2乗で評価するのに対して、 kmedoidsはクラスタの中心をデータ点として距離を1乗で評価します。 パッケージの読み込み library("ggplot2") library("cluster") サンプルデータ作成、可視化用関数 makeSampleData <- function(sd){ df1 <- data.frame(x = rnorm(20, mean=1, sd=sd), y = rnorm(20, mean=1, sd=sd), cluster = rep(1, 20)) df2 <- data.frame(x = rnorm(20, mean=1, sd=sd), y = rnorm(20, mean=3, sd=sd), cluster = rep(2, 20)) df3 <- data.frame(x = rnorm(20, mean=3, sd=sd), y = rnorm(20, mean=1, sd=sd), cluster = rep(3, 20)) df4 <- data.frame(x = rnorm(20, mean=3, sd=sd), y = rnorm(20, mean=3, sd=sd), cluster = rep(4, 20)) df5 <- data. Read more

February 24, 2015

dplyrの暗黒面

> library(dplyr) > iris %>% + summarise_(heikin = "mean(Sepal.Length)") heikin 1 5.843333 と実行できるが、 > myfunc <- function(x){ + mean(x) + } > iris %>% + summarise_(heikin = "myfunc(Sepal.Length)") Error in summarise_impl(.data, dots) : could not find function "myfunc" は実行出来ない。 > iris %>% + summarise_(heikin = as.formula("~myfunc(Sepal.Length)")) heikin 1 5.843333 とすれば、実行できた。

February 22, 2015

最強のハードディスクはどれだ?

第46回R勉強会@東京(#TokyoR)でLTをしてきました。 資料は、こちらです。 最強のハードディスクはどれだ? from Atsushi Hayakawa また、https://github.com/gepuro/gweibullplotでワイブルプロットのパッケージを公開しました。本資料を作成する際に用いた関数とは別のものです。(研究室秘伝のソースでライセンスが謎なため)

February 7, 2015

ハードディスクの寿命分布を比較

Hard Drive Data Setsに4万台以上のハードディスクに関するデータが公開されている。これらの寿命分布を推定してみようと思う。 公開されているデータをSQLiteに展開し、前処理 create table total_days as select serial_number, model, count(*) as total_days, max(failure) as failure, min(date(date)) as min_date, max(date(date)) as max_date from drive_stats group by serial_number ; Rで読み込み、解析準備 library(RSQLite) library(dplyr) drv<-dbDriver("SQLite") con<-dbConnect(drv,dbname="drive_stats.db") totaldays <- dbGetQuery(con,"select * from total_days") head(totaldays) totaldays %>% group_by(model) %>% summarise(model_sum = n()) %>% arrange(model_sum) hdd.df <- totaldays %>% filter(model %in% c("ST4000DM000", "HGST HMS5C4040ALE640")) weib.df <- data.frame( x = hdd.df$total_days, condition = ifelse(hdd.df$model == "ST4000DM000", 1, 2), status = hdd. Read more

January 26, 2015

指数分布の確率プロット

2年ぐらい前にお勉強用に書いたコードを再利用するために発掘したので、共有しておきます。 “exp”の箇所を好きな分布に書きかえれば流用できるかと思います。 exponential.plot <- function(x){ qq <- function(x){ n <- length(x) qexp((1:n - 1/2)/n, rate=mean(x))[order(order(x))] # ミラーランク } probs<- c(0.001, 0.01, 0.1, 1, 5, 10, 20, 30, 50, 70, 80, 90, 95, 99, 99.9, 99.99, 99.999) / 100 plot(c(min(x), max(x)), qexp(c(probs[1], probs[length(probs)])), type="n", xaxt="n", yaxt="n", ylab="Cumulative Percent", xlab="Observed Value") abline(h=qexp(probs), v=as.integer(seq(min(x), max(x))), col="grey") abline(h=qexp(0.5)) axis(side=2, at=qexp(probs[c(4,8:17)]), las=1, labels=as.character(probs[c(4,8:17)]*100), cex.axis=0.8) axis(side=1, as.integer(seq(min(x), max(x))), cex.axis=0.8) points(x, qq(x)) model <- lm(qq(x) ~ x - 1) abline(model) print(paste("r-squared",summary(model)$r. Read more

December 22, 2014

回帰木+ヒストグラムの図を描く

回帰木の下にヒストグラムがあると格好良いと思ったので、メモしておく。 library("party") data(cars) model <- ctree(dist~., data=cars) t.style <- node_hist(model, ymax=0.06, xscale=c(0,150), col="red", fill=hsv(0.6, 0.5, 1)) plot(model, terminal_panel = t.style) 参考 [Rによるデータサイエンス13「樹木モデル」 ](http://www.slideshare.net/takemikami/r13-9821987)

December 22, 2014

集計時の順序を曜日で並べる

weekdays()を使って、Date型から曜日を取得できる。 > (today <- Sys.Date()) [1] "2014-12-22" > str(weekdays(today)) chr "月曜日" 戻り値が文字列になってるので、集計時に曜日順にならない。 > x <- c("月曜日","日曜日", "火曜日", "日曜日", "土曜日") > table(x) x 火曜日 月曜日 土曜日 日曜日 1 1 1 2 factor型にしてlevelsを指定することによって、順番を制御できる。 > yobi <- factor(x, levels=paste(c("日", "月", "火", "水", "木", "金", "土"), "曜日", sep="")) > table(yobi) yobi 日曜日 月曜日 火曜日 水曜日 木曜日 金曜日 土曜日 2 1 1 0 0 0 1 嬉しい! twitterで教えてもらった書き方は、曜日順での集計に関するメモ書き です。

December 19, 2014

Japan.R 2014を支える技術

R Advent Calendar 2014の18日目担当です。 12月6日にJapan.R 2014が開催されました。 総勢250名を超える方に参加して頂いて非常に嬉しいです。@0kayuさんを始め、会場提供およびスポンサーになって頂いたフリークアウトさん、 当日に手伝いをしてくれた方々、ありがとうございます! Japan.R 2014を開催するまでの話を書きます。 主催の一人に選ばれる ヤフー株式会社で開催されたJapan.R 2013での二次会後に、 主催の@holidayworkingさんに無茶振りされました。その時にいた@0kayuさんも巻き添え?になってしまいました。 2010年のJapan.Rに参加して以来、頻繁にTokyo.Rに参加していてRコミュニティに貢献したかったので、受けることにしました。 開催準備 6月頃に、12月6日にフリークアウトで開催することを決定 7月に参加して勉強会の懇親会で、パネルディスカッションのアイディアを頂く パネラー探しを始める。最終的に決まったのは9月下旬 30分枠の発表をお願い出来た 7月に オプトの@shsaixさん フリークアウトの@yanaokiさん 後に RCOの@TJO_datasciさん @teramonagiさんや@dichikaさんの紹介で@AriLamsteinさん 8月のTokyo.RでLTの募集を開始する(google docsのフォームを利用) 開催内容を紹介する。http://www.slideshare.net/yurieoka37/japanr 10月にフリークアウトさんがスポンサーになってくれることが決まった。ありがとうございます! 10月下旬にATNDをたてる Bar doradoraのお願いをする 11月上旬に@0kayuさんとピザを食いながら、パネルディスカッションの内容を打ち合わせ。ピザ生地が旨かった。 11月下旬に某社の会議室を貸して頂いて、パネラーの皆さんに集まってもらってパネルディスカッションの内容を打ち合わせ 12月始め レガラートで食べ物を注文 パネルディスカッション用のマイクをレンタル 二次会のエイトを予約 当日 11:15に当日スタッフの集合 会場の準備 Read more

December 11, 2014

xargsまたはparalellを用いてRを並列処理させる

list.txtから入力を受け取り、hoge.Rで並列処理をさせる。 hoge.R argv <- commandArgs(TRUE) argv <- strsplit(argv, ",")[[1]] par1 <- as.numeric(argv[1]) par2 <- as.numeric(argv[2]) print(par1) print(par2) list.txt 1,1 1,2 2,1 2,2 xargsを用いる場合 cat list.txt | xargs -P8 -n1 Rscript hoge.R paralellを用いる場合 cat list.txt | parallel -j 8 Rscript hoge.R {}

December 10, 2014

Rstanで重回帰分析

ワイブルばかりだったので、メジャーな重回帰分析をしてみた。Rstanにハマりつつある。 linear_regression.stan data { int<lower=1> J; // number of data int<lower=1> K; // number of covariate real y[J]; matrix[J,K] x; } parameters { real<lower=0> sd0; vector[K] beta; real beta0; } model { for(j in 1:J){ increment_log_prob(normal_log(y[j], x[j] * beta + beta0, sd0)); } } linear_regression.R library(rstan) n <- 100 x1 <- rnorm(n, mean=1, sd=0.1) x2 <- rnorm(n, mean=1, sd=0.1) beta1 <- 1 beta2 <- 2 beta0 <- 3 y <- beta1 * x1 + beta2 * x2 + beta0 data. Read more

December 4, 2014

Rstanでワイブル回帰の推定(説明変数複数バージョン)

matrixやvectorの使い方を覚えた。 weibull_regression.stan data { int<lower=1> J; // number of data int<lower=1> K; // number of covariate int status[J]; real<lower=0> tt[J]; matrix[J,K] x; } parameters { real<lower=1, upper=3> m; vector[K] beta; real beta0; } transformed parameters { real<lower=0> eta[J]; for (j in 1:J) eta[j] <- exp(x[j]*beta + beta0); } model { for(j in 1:J){ if(status[j] == 1){ increment_log_prob(weibull_log(tt[j], m, eta[j])); } if(status[j] == 0){ increment_log_prob(weibull_ccdf_log(tt[j], m, eta[j])); } } } weibull_regression. Read more

December 4, 2014

Rstanでワイブル回帰の推定

weibull_regression.stan data { int<lower=1> J; // number of data int status[J]; real<lower=0> tt[J]; real x[J]; } parameters { real<lower=1, upper=3> m; real beta; real beta0; } transformed parameters { real<lower=0> eta[J]; for (j in 1:J) eta[j] <- exp(beta * x[j] + beta0); } model { for(j in 1:J){ if(status[j] == 1){ increment_log_prob(weibull_log(tt[j], m, eta[j])); } if(status[j] == 0){ increment_log_prob(weibull_ccdf_log(tt[j], m, eta[j])); } } } mを1から3ぐらいにしておかないと、推定値が遥か彼方に飛ばされることがある。ワイブル分布は厄介やね。 weibull_regression.R library(rstan) set.seed(0) n <- 1000 # 台数 m <- 2 # 形状パラメータ beta <- -1 # 係数 beta0 <- 10 # 切片 x <- rnorm(n) data. Read more

December 4, 2014

Rstanでワイブル分布のパラメータ推定(打切りのあるデータで)

weibull.stan data { int<lower=1> J; // number of data int status[J]; real<lower=0> tt[J]; } parameters { real<lower=1> m; real<lower=1> eta; } model { for(j in 1:J){ if(status[j] == 1){ increment_log_prob(weibull_log(tt[j], m, eta)); } if(status[j] == 0){ increment_log_prob(weibull_ccdf_log(tt[j], m, eta)); } } } weibull.R サンプルデータの生成 > set.seed(0) > data.df <- data.frame(tt = rweibull(1000, shape = 2, scale = 100), + status = rep(1,1000)) > data.df[data.df$tt > 30,]$status <- 0 > data. Read more

December 3, 2014

Rstanでワイブル分布のパラメータ推定

Rstanでワイブル分布のパラメータ推定をしてみる。打切りを考慮するにはどうしたら良いのだろうか? Stan Modeling Language User’s Guide and Reference Manual, v2.5.0の83ページ辺りを読めば出来そう。 weibull.stan data { int<lower=0> J; // number of data real<lower=0> tt[J]; // time } parameters { real<lower=0> m; real<lower=0> eta; } model { for(j in 1:J) tt[j] ~ weibull(m, eta); } weibull.R library(rstan) data.ls <- list(J = 1000, tt = rweibull(1000, shape = 2, scale = 100)) fit <- stan(file = 'weibull.stan', data = data.ls, iter = 1000, chains = 4) 結果 > fit Inference for Stan model: weibull. Read more

© gepuro 2013