to be continued...

Rや心理統計学の備忘録的な

There are three kinds of lies: lies, damned lies, and statistics.

- 嘘には三種類ある。嘘、大嘘、そして統計だ -
Benjamn Disraeli(19世紀のイギリス首相)

Windows RMeCabで NEologdを諦めるにはまだ早い(Dockerを使えばいいじゃない)

どうも。最近は「彼方のアストラ」を見ています。一話目でうるうるきてました。

Windowsを使っていると、RMeCabを使うときに辞書として mecab-ipadic-NEologd(neolog)を使うのは難しいですよね。出来なくはない(はず)ですし、いろいろと紹介してくれる記事はありますが、難しいんですよね。(私は挫折しました)

今日は、RMeCabにmecab-ipadic-NEologdを搭載する方法としてDockerを使った方法を紹介します。

※本記事はRMeCabを知っていることを前提としておりますので、ご了承ください。

結論

Windowsユーザーでmecab-ipadic-NEologdを辞書に指定してRMeCabを使いたい人は

$docker pull ando6oid/neo-mecab
$docker run -e PASSWORD=umr -p 8787:8787  -d  ando6oid/neo-mecab
//このときパスワードは「umr」となる

を打ち込んで、「http://localhost:8787/」を開くと使えます。

docker hubのページはこちら↓↓↓

hub.docker.com


mecab-ipadic-NEologdとは

MeCabの辞書です。形態素解析では、文章を分かち書き(を区切る)を行います。そのとき、どこからどこまでをとするのかの基準として辞書を使います。

こんな感じです↓

f:id:ando_Roid:20190821160958p:plain

雨宮天」が"雨宮"と"天"、闇営業が"闇"と"営業"になっています。

雨宮天さんと闇営業は一切関係ありません

(雨宮天が分からない人は義務教育からやり直してほしい)

本来は、「雨宮天」「闇営業」というように名詞として処理してほしいのに...ということが初期のRMeCabだと起こります。ここでneologを辞書にすると...↓↓

f:id:ando_Roid:20190821161009p:plain

望むべき処理をしてくれるようになります。

neologは週二で新しい単語を追加してくれるので、neologを辞書とすると流行語やドラマの名前、人の名前に対応できるようになります。

他にもNeologのすごさについてLINEさんも書いています↓↓

engineering.linecorp.com

形態素解析に必須レベルなのですが、Windowsだと入れるのが難しい。PATHがどうちゃらとかなんちゃらで...Linuxなら簡単なんだよなぁ

Linux?じゃあ、Dockerを使おう。


Dockerとは

Dockerは、コンテナと呼ばれるOSレベルの仮想化環境を提供するオープンソースソフトウェアである(Wikipedia)。私もDockerについては勉強中ですが、簡単に言うとDocker(具体的にはコンテナ)を使うことでWindowsでもLinux環境で作業することが可能になります。もっとほかのことも出来ます。

さらに、Dokcerhubを使うことで、先人たちが作ってくれたDocker imageを使ってコンテナを簡単に作ることができます(まさに、巨人の肩に乗る)。

(他のことはググってください)


Dockerのインストール方法

あたいのWindowsがproかhomeかわからない方

それ以外の方はググってください(私もよくわかりません、すみません)


Dockerがインストールできたら

  • power shellを起動(わからない人はcortanaにきいてね)

  • docker pull ando6oid/neo-mecabと打ち込む

    dockerhubからneologを持ったRstudioサーバーイメージを持ってくる

  • dockcer imagesでando6oid/neo-mecabイメージがあるか確認

  • docker run -e PASSWORD=すきなやつ -p 8787:8787 -d ando6oid/neo-mecabと打ち込む(passwordはお好きにどうぞ)

    イメージからコンテナを作成

  • chromeを立ち上げてアドレスバーにhttp://localhost:8787/とうつ

  • そこにRstudioがあれば完了

neologが初期辞書になっているので、いつも通りにRMeCabを使ってくれればok


ando6oid/neo-mecabについて

私が作ったDocker imageです。

ymattu/mecab-dをベースとしています(ymattuさん神)。ymattu/mecab-dもrocker/tidyverseというイメージをベースとしています。したがって、tidyverse系パッケージが基本的に入っている上、ymattuさんによって日本語に対応されつつRMeCabtidytextもあるというイメージです。私はそこにmecab-ipadic-NEologdを追加しただけ。

y-mattu.hatenablog.com

要するに、tidyverse,tidytext.RMeCabといったおおよそ必要になるだろうパッケージが入っていて、日本語にも対応しており辞書も最新なRstudio。RMeCabによるテキストマイニングの決定版イメージになりうる(自画自賛、ほとんど私は何もしていない)。

dockerfileはこんな感じ。

FROM ymattu/mecab-d:latest

LABEL mainaier = "ando_Roid"

RUN echo "now building..."

CMD echo "now running..."

RUN apt-get -qq update\
 && apt-get -qq -y install curl \
 && apt install -y xz-utils file \
 && apt install patch

WORKDIR /home/rstudio

RUN git clone --depth 1 https://github.com/neologd/mecab-ipadic-neologd.git


WORKDIR mecab-ipadic-neologd

RUN ./bin/install-mecab-ipadic-neologd --create_user_dic -n -y\
 && cat /usr/local/etc/mecabrc | sed "s/ipadic/mecab-ipadic-neologd/" > /usr/local/etc/mecabrc

EXPOSE 8787

CMD ["/init"]

注意

mecab-ipadic-NEologdはdocker imageをプルした時の最新のものを使用しており、辞書の自動更新機能はございません。
気が向いたら、追加するかも。

従って、最新の辞書を使いたい場合は、docker imageをプルしなおす必要があります。

なにか不備等ありましたらコメントください

to be continued...

Rで分位分け関数

Rで分位分け関数

どうも、最近は「女子高生の無駄遣い」見ています。

「輪! moon! dass! cry!」が頭から離れない。

今回は、ショート記事。

Rで分位分け関数を作ります

(分位分けという表現あってるかわからない)

というのも、機械学習や統計分析をしてるとデータの10のランクに分けたりして、データの上位40%はどうなのかを見たり、中央についてみたりします。

ググったけどなかったので、個人的なメモとして残しておきます。

(余裕でdplyrにありました)


使う関数は

  • ceiling():x未満ではない最小の整数を返す

    例:0→0, 1→1, 1.001→2, 2.999→3

  • rank():データの大きさの順番を返す関数(順序尺度変換)

  • length():ベクトルの要素数を返す関数

そして関数がこれ

sep_tile <- function(x,tile = 10){
  res <- ceiling(rank(x)/length(x)*tile)
  return(res)
}

ほんまに使えるの?


データ生成

> set.seed(42)
> umr <- runif(100,0,1)
> head(umr)
[1] 0.9148060 0.9370754 0.2861395 0.8304476 0.6417455 0.5190959

シードを固定して、[0,1]の一様乱数を100個生成したものをumrに格納

仮想データに自作関数を適用してみると

> sort(sep_tile(umr))
  [1]  1  1  1  1  1  1  1  1  1  1  2  2  2  2  2  2  2  2  2  2  3  3  3  3  3  3  3  3  3  3  4  4  4  4  4  4  4  4  4
 [40]  4  5  5  5  5  5  5  5  5  5  5  6  6  6  6  6  6  6  6  6  6  7  7  7  7  7  7  7  7  7  7  8  8  8  8  8  8  8  8
 [79]  8  8  9  9  9  9  9  9  9  9  9  9 10 10 10 10 10 10 10 10 10 10

わー、どうやら出来てるっぽい。


関数の説明

tileの引数を10として、10分位に分けます

  • rank()でデータを順序尺度に変換します

    1~100のランクが与えられます

rank(umr)
  [1]  88  92  28  80  59  50  71  13  61  67  42  68  91  25  43  93  98  12  44  51  85  14 100  94   9  46  37  86  41
 [30]  82  72  78  36  65   3  81   4  20  87  55  35  40   6  97  39  95  84  58  96  56  31  32  38  77   7  74  64  17
 [59]  26  48  63  99  75  53  83  18  27  79  66  24   8  15  22  45  19  69   5  34  47   2  54  16  33  60  76  52  23
 [88]  11  10  29  62   1  21  90  89  70  30  49  73  57
  • length()で割ります

    0~1の値に変換されます

rank(umr)/length(umr)
  [1] 0.88 0.92 0.28 0.80 0.59 0.50 0.71 0.13 0.61 0.67 0.42 0.68 0.91 0.25 0.43 0.93 0.98 0.12 0.44 0.51 0.85 0.14 1.00
 [24] 0.94 0.09 0.46 0.37 0.86 0.41 0.82 0.72 0.78 0.36 0.65 0.03 0.81 0.04 0.20 0.87 0.55 0.35 0.40 0.06 0.97 0.39 0.95
 [47] 0.84 0.58 0.96 0.56 0.31 0.32 0.38 0.77 0.07 0.74 0.64 0.17 0.26 0.48 0.63 0.99 0.75 0.53 0.83 0.18 0.27 0.79 0.66
 [70] 0.24 0.08 0.15 0.22 0.45 0.19 0.69 0.05 0.34 0.47 0.02 0.54 0.16 0.33 0.60 0.76 0.52 0.23 0.11 0.10 0.29 0.62 0.01
 [93] 0.21 0.90 0.89 0.70 0.30 0.49 0.73 0.57
  • 10を掛けます

    ここで掛ける値が分けたい分位数になります(0~10の値に変換)

    今回は10分位に分けました

rank(umr)/length(umr)*10
  [1]  8.8  9.2  2.8  8.0  5.9  5.0  7.1  1.3  6.1  6.7  4.2  6.8  9.1  2.5  4.3  9.3  9.8  1.2  4.4  5.1  8.5  1.4 10.0
 [24]  9.4  0.9  4.6  3.7  8.6  4.1  8.2  7.2  7.8  3.6  6.5  0.3  8.1  0.4  2.0  8.7  5.5  3.5  4.0  0.6  9.7  3.9  9.5
 [47]  8.4  5.8  9.6  5.6  3.1  3.2  3.8  7.7  0.7  7.4  6.4  1.7  2.6  4.8  6.3  9.9  7.5  5.3  8.3  1.8  2.7  7.9  6.6
 [70]  2.4  0.8  1.5  2.2  4.5  1.9  6.9  0.5  3.4  4.7  0.2  5.4  1.6  3.3  6.0  7.6  5.2  2.3  1.1  1.0  2.9  6.2  0.1
 [93]  2.1  9.0  8.9  7.0  3.0  4.9  7.3  5.7
  • ceiling()でx未満ではない最小の整数を返します

    こんな感じ

> ceiling(0)
[1] 0
> ceiling(0.1)
[1] 1
> ceiling(1.03)
[1] 2

これで完成

ceiling(rank(umr)/length(umr)*10)
  [1]  9 10  3  8  6  5  8  2  7  7  5  7 10  3  5 10 10  2  5  6  9  2 10 10  1  5  4  9  5  9  8  8  4  7  1  9  1  2  9
 [40]  6  4  4  1 10  4 10  9  6 10  6  4  4  4  8  1  8  7  2  3  5  7 10  8  6  9  2  3  8  7  3  1  2  3  5  2  7  1  4
 [79]  5  1  6  2  4  6  8  6  3  2  1  3  7  1  3  9  9  7  3  5  8  6

ほんとに分位分け出来てるの?

> ceiling(rank(umr)/length(umr)*10)
  [1]  9 10  3  8  6  5  8  2  7  7  5  7 10  3  5 10 10  2  5  6  9  2 10 10  1  5  4  9  5  9  8  8  4  7  1  9  1  2  9
 [40]  6  4  4  1 10  4 10  9  6 10  6  4  4  4  8  1  8  7  2  3  5  7 10  8  6  9  2  3  8  7  3  1  2  3  5  2  7  1  4
 [79]  5  1  6  2  4  6  8  6  3  2  1  3  7  1  3  9  9  7  3  5  8  6
> sort(sep_tile(umr,tile = 5))
  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
 [59] 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5
> sort(sep_tile(umr,tile = 3))
  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 [59] 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
> sort(sep_tile(umr,tile = 100))
  [1]   1   2   3   4   5   6   8   8   9  10  11  12  13  15  15  16  17  18  19  20  21  22  23  24  25  26  27  29  29
 [30]  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  56  57  57  58
 [59]  59  60  61  62  63  64  65  66  67  68  69  70  71  72  73  74  75  76  77  78  79  80  81  82  83  84  85  86  87
 [88]  88  89  90  91  92  93  94  95  96  97  98  99 100
> sort(sep_tile(umr,tile = 20))
  [1]  1  1  1  1  1  2  2  2  2  2  3  3  3  3  3  4  4  4  4  4  5  5  5  5  5  6  6  6  6  6  7  7  7  7  7  8  8  8  8
 [40]  8  9  9  9  9  9 10 10 10 10 10 11 11 11 11 11 12 12 12 12 12 13 13 13 13 13 14 14 14 14 14 15 15 15 15 15 16 16 16
 [79] 16 16 17 17 17 17 17 18 18 18 18 18 19 19 19 19 19 20 20 20 20 20

こげな感じです。

今回のコードまとめ

#シードを固定
set.seed(42)
#[0,1]の一様乱数を100個→umrに格納
umr <- runif(100,0,1)
#10分位に分ける
ceiling(rank(umr)/length(umr)*10)
#分位分け自作関数
sep_tile <- function(x,tile = 10){
  res <- ceiling(rank(x)/length(x)*tile)
  return(res)
}
#関数使ってみた
sep_tile(umr)

わーい。


{dplyr}に搭載されていた...

> library(tidyverse)
> ntile(umr,10)
  [1]  9 10  3  8  6  5  8  2  7  7  5  7 10  3  5 10 10  2  5  6  9  2 10 10  1  5  4  9  5  9  8  8  4  7  1  9  1  2  9  6  4  4  1 10  4 10  9  6
 [49] 10  6  4  4  4  8  1  8  7  2  3  5  7 10  8  6  9  2  3  8  7  3  1  2  3  5  2  7  1  4  5  1  6  2  4  6  8  6  3  2  1  3  7  1  3  9  9  7
 [97]  3  5  8  
to be continued...

「おんやど恵」で開発合宿

どうも。最近は「ダンベル何キロ持てる?」を見ています。

先日おんやど恵さんにてゼミ合宿を行ったので共有を。

開発合宿って知ってる?

私のゼミでは夏に二泊三日の合宿を行い、そこでデータ分析や研究発表を行っています。

去年のゼミ合宿でも幹事をおこなっていた私は、今年のゼミ合宿を企画している時に開発合宿というワードを耳にしました。

開発合宿は、旅館などの宿泊施設を利用して泊りがけで開発(主に、プログラミング)を行うものです。

普段、学校やオフィスにて開発をしていると可処分時間に制限あったり、開発以外の作業があったりと障害となるものがあると思います。開発合宿では、宿泊施設を利用しているので食事や風呂など身の回りの世話はしてくれるし、学校でも会社でもないので作業に没頭できるというわけです。更に、自然に囲まれた非日常的空間やおいしいご飯に温泉は、作業効率をあげるのはずです。

開発合宿の条件

そんなわけで、今年は開発合宿に適してる宿泊施設を探しておんやど恵さんに決めました。

おんやど恵さんは、主が昔プログラミングをしていたことから、「集中して開発できる環境を提供したい」ということでリーズナブルな価格で開発合宿プランを提供しています。

(ありがてぇ...)

都心から近いことからも様々な企業が開発合宿の場として利用しているようです。

詳しくは、以下を参照↓↓↓

http://www.onyadomegumi.co.jp/plan3_lp/

開発合宿に適している条件としては...

Must

  • Wi-Fiがある。しかも、早い
  • 開発用の部屋がある。大きな会議室的な
  • プロジェクターを使って大きな画面を映し出せる
  • 一日中、開発しててもお咎めなし。夜中ににぎやかにしてても問題なし

Want

といったものがあるかなと思います。

おんやど恵では、開発合宿プランを提供しており上記の条件をほぼすべて満たしていると言っても過言ではありません。旅館自体、なかなかにランクの高い旅館なのですが、そこが開発環境を提供してくれているので正直、最強です。

以下では、おんやど恵利用したゼミ合宿を振り返っていきます。なにか旅館、ゼミ合宿、開発合宿探しの一助になれば幸いです。


合宿の内容

私の所属するゼミでは、二つのデータ分析コンペの出場を視野に入れており合宿では二つのチームに分かれて分析することになっていました。以下では、合宿の行程を書いていきます。


一日目

12:45 JR湯河原駅集合

旅館はJR湯河原駅から出ているバスで約10分程度いったところにある「理想郷」バス停の近くにあります。

バス自体は15分間隔で出ているので不便はないです。

13:10 「おんやど恵」到着

チェックインです。客室の利用は15:00からですが、今回利用した開発合宿プランでは会議室を13:00から利用することが出来ます。これは嬉しい。会議室は広々としており参加者9人では持て余すほどです。到着してすぐに部屋のセッティングと買い出しをします。買い出しは旅館のすぐ隣にローソンがあるので非常に便利!!

f:id:ando_Roid:20190809213057j:plain

(おい)

13:45 開発

買い出し、セッティングも終わって作業に着手し始めます。

最初に、各チームの方針と今日の目標を発表して、あとは作業。

部屋が広いので、机で二つの島を作ってホワイトボードでパーテーション代わりにします。

チームA

データを可視化して傾向や特徴を見出しています。

こっちのチームは教授と3年生がメインのチームなので、学部生のアイデアを教授が実装するという感じ。3年生にとって初めてのデータ分析コンペとなりますが闊達な議論が繰り広げられています。

f:id:ando_Roid:20190809213128j:plain

チームB

こっちは主として院生のチーム(うp主はこっち)

チームAと同様にデータ分析のストーリーを作る段階ですが、可視化よりも探索的な分析をかけることを軸としています。ホワイトボードに書かれているのはLasso回帰について

f:id:ando_Roid:20190809213156j:plain

こんな感じで、各チームとも適宜休憩を取りながら開発を進めていきます。無駄話も交じえつつ作業は進む

18:00 夕食&温泉

一日目の作業を終えて夕食です。

写真撮り忘れましたがおいしかった。。。

温泉もいいぞ~

21:00 自由時間

夕食後は風呂に入った後に会議室に集合。

開発を再開してもよいし、遊んでてもいい。素晴らしい時間

麻雀してる人(の後ろで事務作業してる人)

f:id:ando_Roid:20190809213249j:plain

ビール片手に開発してる人(とスマホ見てる人)

f:id:ando_Roid:20190809213330j:plain

プロジェクターと旅館にあるBlue-toothスピーカーを使ってガルパン上映会(ガルパンはいいぞ~)

低音が響くので戦車のSEがすごいことに...!!!

f:id:ando_Roid:20190809213418j:plain

こんな感じで、個人がやりたいことを過ごして夜は更けていきます。

みんな明日も頑張ろうな


二日目

f:id:ando_Roid:20190809213457j:plain

8:00 朝食

昨晩はいつ寝たのかわからないけど朝食。

結構がっつり朝ご飯を食べた気がする。

10:00 開発

食休みを挟んで開発作業。

昨日の進捗をお互い報告しあいます。

各チームから意見をしあい、今日の方針を発表して作業に着手。

f:id:ando_Roid:20190809213524j:plain

13:00 昼食&観光

昼食は旅館では出ないので食事は外で。

ここで教授から、「せっかくなら観光もしたいよね」と鶴の一声。

急遽、レンタカーを借りて箱根は芦ノ湖を目指します。

車で40分くらい。旅館周辺は飲食店が少し少ないので、良い選択でした。

f:id:ando_Roid:20190809213602j:plain

芦ノ湖のほとりでランチ

f:id:ando_Roid:20190809213635j:plain

炎天下の中、アヒルボート

f:id:ando_Roid:20190809213709j:plain

17:00 開発

がっつり昼食&観光をして既に疲れていますが、作業はします。

眠いですが、作業をします。

18:00 夕食

ステーキが出ました!!

f:id:ando_Roid:20190809213741j:plain

21:00 自由時間

昨日同様の自由時間。

大富豪したり、開発したり、ロケットリーグしたり、ガルパンのゲームしたりします。

(どんだけガルパン好きやねん)

この日は何時に寝たのか覚えていません。


三日目

8:00  朝食

寝かしてくれと思いながら朝食。

やっぱ朝食はかなりがっつり出るので、おなか一杯。

10:00 進捗・研究発表

作業に終始してもよかったのですが、

せっかく院生もいるので研究発表をしてもらいます。

f:id:ando_Roid:20190809213810j:plain

学部生も院生の研究を聞いて刺激をもらいます。

その後、最終的なチームの進捗を報告し合います。

12:00 合宿終了

12時には退館。

チェックアウト自体は9:30ですが会議室は12:00まで利用できます。

結構、うれしいですよね。


終わりに

ということで、おんやど恵さんの開発合宿プランを利用してゼミ合宿を行ってきました。

開発合宿の醍醐味は開発に集中できる環境が用意されつつ、それ以外の観光や旅館のサービスを享受できることだと思います。開発に没頭するもよし、プライベートな時間を大切にするもよし。

いつも同じ環境での作業がマンネリ化しているあなた

数日間チームで集中して作業を行いたいあなた

是非、開発合宿...もといおんやど恵で開発合宿をしてみては?

to be continued...

ごん(R)、お前だったのか...ベクトルをリサイクルしてくれていたのは...

初めに

どうもです。最近はガルパンを見ています。

(カチューシャが頭から離れない)

(最終章第二話面白かったです)

生きていると、R上でのベクトルのリサイクルを感じたいですよね。

今日は、Rで起こるベクトルリサイクルの性質と留意すべきポイントについてのお話です

こんなミス知らないうちにやっていそう....

結論

  • Rは論理式でベクトル評価をするとき、両ベクトルの長さが等しくないと短い方のベクトルを長いベクトルの要素数と等しくなるまでリサイクルする(論理式だけじゃないかも...)
  • 上述のため、エラーは返されないが予期せぬアウトプットになっているかもしれないことを留意せよ

現象

> sessionInfo()
R version 3.6.0 (2019-04-26)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 18.04.2 LTS

> library(tidylog)
# library(tidyverse)でも大丈夫

data("diamonds")
# A tibble: 53,940 x 10
   carat cut       color clarity depth table price     x     y     z
   <dbl> <ord>     <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1 0.23  Ideal     E     SI2      61.5    55   326  3.95  3.98  2.43
 2 0.21  Premium   E     SI1      59.8    61   326  3.89  3.84  2.31
 3 0.23  Good      E     VS1      56.9    65   327  4.05  4.07  2.31
 4 0.290 Premium   I     VS2      62.4    58   334  4.2   4.23  2.63
 5 0.31  Good      J     SI2      63.3    58   335  4.34  4.35  2.75
 6 0.24  Very Good J     VVS2     62.8    57   336  3.94  3.96  2.48
 7 0.24  Very Good I     VVS1     62.3    57   336  3.95  3.98  2.47
 8 0.26  Very Good H     SI1      61.9    55   337  4.07  4.11  2.53
 9 0.22  Fair      E     VS2      65.1    61   337  3.87  3.78  2.49
10 0.23  Very Good H     VS1      59.4    61   338  4     4.05  2.39
# … with 53,930 more rows

今回は53940行あるdiamondsデータを使います

みんなのR[第二版]の第12章「dplyrによる高速なグルーピング処理」を読んでいたときのこと。

filterについて、条件式を複数用いて抽出する際の書き方として論理和(積)を使う方法と%in%(in演算子)を使う方法が書かれています。

参考までに↓↓↓

(dplyr::filterの再確認)

(%in% operator in R)

2つをそれぞれ書いてみると

> diamonds %>% filter(cut %in% c("Ideal", "Good"))
filter: removed 27483 out of 53940 rows (51%)
# A tibble: 26,457 x 10
   carat cut   color clarity depth table price     x     y     z
   <dbl> <ord> <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1  0.23 Ideal E     SI2      61.5    55   326  3.95  3.98  2.43
 2  0.23 Good  E     VS1      56.9    65   327  4.05  4.07  2.31
# … with 26,447 more rows

> diamonds %>% filter(cut == "Ideal" | cut == "Good")
filter: removed 27483 out of 53940 rows (51%)
# A tibble: 26,457 x 10
   carat cut   color clarity depth table price     x     y     z
   <dbl> <ord> <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1  0.23 Ideal E     SI2      61.5    55   326  3.95  3.98  2.43
 2  0.23 Good  E     VS1      56.9    65   327  4.05  4.07  2.31
# … with 26,447 more rows

という感じになります。

どちらも残っている要素は26457となっています。

ここでin演算子について「え?いらなくね?」と考えた愚かな私。

試しに書いてみると

> diamonds %>% filter(cut == c("Ideal", "Good"))
filter: removed 40711 out of 53940 rows (75%)
# A tibble: 13,229 x 10
   carat cut   color clarity depth table price     x     y     z
   <dbl> <ord> <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1  0.23 Ideal E     SI2      61.5    55   326  3.95  3.98  2.43
# … with 13,219 more rows

行けたわ...

あり?要素が13229になってる。大体半分くらい?

どういうことだってばよ?

説明

ある塾で聞いてみると、

「基本的に、2つのベクトルオブジェクト(A , B)があって、それをマッチさせるような処理( 例えば A + B とか A == B など)をする場合、かつその二つの長さが違う場合には、R言語は短い方のベクトルをリサイクルして長さを合わせようとします。」

という解を得られました。

記号の書き方と集合の問題ですね。

in演算子は、要素がベクトルに属するか否かを返すものです。

対して、filter(cut == c("Ideal", "Good"))という書き方は、

元のdaiamondsデータが53940行要素あるのに対して、c("Ideal", "Good")の2つの要素しかありません。

他の言語ではわかりませんが、

Rでは論理式を評価するときにどちらかのベクトルが短いと短い方をリサイクルして長さを合わせようとするのです。

つまり、"Ideal", "Good"を26970回リサイクルして53940要素に増やすのです。

感動....!!!

例えば、diamonds %>% filter(cut == "Ideal")と書いたとき、R上では

> diamonds %>% slice(1:10) %>% select(carat,cut,color) %>% 
+   mutate(cut_check = rep("Ideal",10),
+          answer = cut_check == cut)

# A tibble: 10 x 5
   carat cut       color cut_check answer
   <dbl> <ord>     <ord> <chr>     <lgl> 
 1 0.23  Ideal     E     Ideal     TRUE  
 2 0.21  Premium   E     Ideal     FALSE 
 3 0.23  Good      E     Ideal     FALSE 
 4 0.290 Premium   I     Ideal     FALSE 
 5 0.31  Good      J     Ideal     FALSE 
 6 0.24  Very Good J     Ideal     FALSE 
 7 0.24  Very Good I     Ideal     FALSE 
 8 0.26  Very Good H     Ideal     FALSE 
 9 0.22  Fair      E     Ideal     FALSE 
10 0.23  Very Good H     Ideal     FALSE 

こんな感じで、Idealという単一要素をデータ行分リサイクルしています。

ごん、お前だったのか...ベクトルをリサイクルしてくれていたのは...

そして、TRUEだけを抽出しています。

(cut_checkanswerを擬似的にたした)

ということはdiamonds %>% filter(cut == c("Ideal", "Good"))では、

> diamonds %>% slice(1:10) %>% select(carat,cut,color) %>% 
+   mutate(cut_check = rep(c("Ideal","Good"),5),
+          answer = cut_check == cut)
# A tibble: 10 x 5
   carat cut       color cut_check answer
   <dbl> <ord>     <ord> <chr>     <lgl> 
 1 0.23  Ideal     E     Ideal     TRUE  
 2 0.21  Premium   E     Good      FALSE 
 3 0.23  Good      E     Ideal     FALSE 
 4 0.290 Premium   I     Good      FALSE 
 5 0.31  Good      J     Ideal     FALSE 
 6 0.24  Very Good J     Good      FALSE 
 7 0.24  Very Good I     Ideal     FALSE 
 8 0.26  Very Good H     Good      FALSE 
 9 0.22  Fair      E     Ideal     FALSE 
10 0.23  Very Good H     Good      FALSE 

c("Ideal", "Good")を元データに足りるまで繰り返し続けるわけです(暁美ほむら?)

53940は3の倍数でもあるので

> diamonds %>% filter(cut == c("Ideal", "Good", "premium"))
filter: removed 45079 out of 53940 rows (84%)
# A tibble: 8,861 x 10
   carat cut   color clarity depth table price     x     y     z
   <dbl> <ord> <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1  0.23 Ideal E     SI2      61.5  55     326  3.95  3.98  2.43
 2  0.31 Good  J     SI2      63.3  58     335  4.34  4.35  2.75
 3  0.3  Good  J     SI1      64    55     339  4.25  4.28  2.73
 4  0.31 Good  H     SI1      64    54     402  4.29  4.31  2.75
 5  0.33 Ideal I     SI2      61.8  55     403  4.49  4.51  2.78
 6  0.26 Good  D     VS1      58.4  63     403  4.19  4.24  2.46
 7  0.23 Ideal G     VS1      61.9  54     404  3.93  3.95  2.44
 8  0.35 Ideal I     VS1      60.9  57     552  4.54  4.59  2.78
 9  0.3  Ideal D     SI1      62.1  56     552  4.3   4.33  2.68
10  0.32 Ideal I     VVS1     62    55.3   553  4.39  4.42  2.73
# … with 8,851 more rows

いけます。更に要素数が減っています。

53940では割り切れない7にしてみると

> diamonds %>% filter(cut == c("Ideal", "Good", "premium","umr","kirie","ebina","Shilfinford"))
filter: removed 50186 out of 53940 rows (93%)
# A tibble: 3,754 x 10
   carat cut   color clarity depth table price     x     y     z
   <dbl> <ord> <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1 0.23  Ideal E     SI2      61.5    55   326  3.95  3.98  2.43
 2 0.23  Good  E     VS1      64.1    59   402  3.83  3.85  2.46
 3 0.26  Good  D     VS1      58.4    63   403  4.19  4.24  2.46
 4 0.3   Ideal D     SI1      62.1    56   552  4.3   4.33  2.68
 5 0.75  Ideal G     SI1      62.2    55  2760  5.87  5.8   3.63
 6 0.8   Ideal F     SI2      59.9    59  2762  6.01  6.07  3.62
 7 0.580 Ideal F     VVS1     61.7    56  2772  5.33  5.37  3.3 
 8 0.71  Good  E     VS2      59.2    61  2772  5.8   5.88  3.46
 9 0.72  Ideal G     SI1      61.8    56  2776  5.72  5.75  3.55
10 0.71  Good  F     VS2      57.8    60  2777  5.87  5.9   3.4 
# … with 3,744 more rows
 警告メッセージ: 
1:  `==.default`(cut, c("Ideal", "Good", "premium", "umr", "kirie",: 
   長いオブジェクトの長さが短いオブジェクトの長さの倍数になっていません 
2:  is.na(e1) | is.na(e2): 
   長いオブジェクトの長さが短いオブジェクトの長さの倍数になっていません 

いや、いけるのかい!!!

でも、警告メッセージがしっかり出ていますね。

今回は、filter()で行いましたが、subset()でも当然起きます。

このように、今回はRの親切心について触れることができました。

間違った書き方でも評価はされてしますので、%in%を使うことをお忘れないように(自戒)

to be continued...

Rで取得する(疑似)セッション

※ここでのセッションは厳密にはセッションとは言えないかもしれないので、{lubridate}をつかったパズル的な記事として読んでいただけると幸いです。

最近は、「冴えカノ」みてます。

生きてると、

「ログデータをセッション単位で見たい」

ってことがあると思うんです。

今日は、Rでの(疑似)セッション情報の取得をします。

なぜ、疑似なのかは後で説明します。

(SQLからデータ落とすときに、セッションも取得できるよ)

(でも、ワイはRをいじりたいんや)

セッションとは

ドラマーと怖い指揮者の映画ではありません。

コンピュータ用語で、「一連のインタラクティブな操作のこと。典型的にはログイン(ログイン)してからログアウト(ログオフ)するまでが一つのログインセッション英語版)。」(Wiki引用)

「典型的には」と表現されているようにセッションには、言語やシステムごとにことなった定義があるそうです。

セッションとは通信中に使うトンネルのようなものです。Oracleはログインするとクライアントとサーバの間にトンネルを作り、トンネルの中にデータを流して、データベースを検索したり更新します。

セッションとは、「ユーザーがサイトに流入し離脱するまでの一連の行動」のことを表します。

セッションとは、一連の処理の始まりから終わりまでを表す概念のことです。

要するに、ネット操作における時間的単位のようなもの

例えば、特定のサイトにログインしてからログオフするまで、どのページどのくらい滞在していたかとか、amazonでどうページを見てから購入に至ったかなど。

これを使って、ユーザーがどんな風にネットサーフィンしたのかを時間的に把握することが可能です

ここでは、Googleアナリティクスさんの定義を使って実装していきたいと思います。

定義(Googleアナリティクスさんの)

では、どうなったらセッションが切れて、次のセッションとみなされるのかを説明します

  1. 日付が変わるタイミング

    →23:59に操作、次の日の0:00に操作したとき次のセッションになるよ

  2. セッションがスタートしてから何も操作されずに30分が経過する

    →セッションがスタートして、次の操作が30分後なら次のセッションになるよ

  3. 参照元が変わるタイミング

    検索エンジンで「冴えカノ」と入れて調べてから、同日の30分以内に「冴えカノ♭」で調べたときでセッションは変わるよ

上記が定義になりますが、今回は3番目の参照元による条件を除外します。

というのも、私がいじったデータに検索エンジンについての情報がなかったからです。

(というか「冴えカノ」も「冴えカノ♭」も「冴えカノ」についてググりたいから、同一の流れとして認識しても良いのでは)

そういう点で、今回は上2つの条件を用いた(疑似)セッションの取得ということになります。

大抵のセッションが特定のサイトやサービスでの周遊にフォーカスしているのに対して、

今回やるのは、Web周遊全体にフォーカスしたセッションということになるのでしょうか

前置きが長くなりましたが、Rで実装していきたいと思います

実装

環境と今回使うパッケージ

> sessionInfo()
R version 3.6.0 (2019-04-26)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 18.04.2 LTS
# Package
> # データ整形に必要
> library(tidyverse)
> # 日付データを扱うパッケージ(今回の主人公)
> library(lubridate)
> # なくてもいい(俺が好き)
> library(tidylog)

今回の主役は{lubridate}です。こと日付データを扱う上で右に出るやつはいないらしい。

(いつか使ってみたかった)

{lubridate}については以下が詳しく説明しています。

ログデータ処理で始めるlubridate入門

lubridateパッケージ入門

lubridate::date()base::date()がぶつかることがあるので注意

仮想データセット

本来、ログデータはDBに蓄積されたものを使うものですが、

今回は以下の仮想データフレーム(umr)を使います。

> umr <- data.frame(
+   #id:web周遊した人の記号
+   id = 1,
+   #action_log:サイト接触のタイムスタンプ
+   action_log = c("2018-2-24 15:30:15",
+                  "2018-2-24 15:40:34",
+                  "2018-2-24 16:30:19",
+                  "2018-2-24 16:32:49",
+                  "2018-2-25 09:10:37",
+                  "2018-2-25 09:30:26",
+                  "2018-2-25 23:47:55",
+                  "2018-2-25 23:59:44",
+                  "2018-2-26 00:02:28",
+                  "2018-2-26 00:04:26")  ,
+   #URL:idが接触したサイトのURL
+   URL = c("aaa", "bbb","ccc", "ddd","eee", "fff","ggg", "hhh","iii", "jjj")
+ ) %>% mutate(action_log = ymd_hms(action_log)) %>%
+   as_tibble() %>% print()
mutate: converted 'action_log' from factor to double (0 new NA)
# A tibble: 10 x 3
      id action_log          URL  
   <dbl> <dttm>              <fct>
 1     1 2018-02-24 15:30:15 aaa  
 2     1 2018-02-24 15:40:34 bbb  
 3     1 2018-02-24 16:30:19 ccc  
 4     1 2018-02-24 16:32:49 ddd  
 5     1 2018-02-25 09:10:37 eee  
 6     1 2018-02-25 09:30:26 fff  
 7     1 2018-02-25 23:47:55 ggg  
 8     1 2018-02-25 23:59:44 hhh  
 9     1 2018-02-26 00:02:28 iii  
10     1 2018-02-26 00:04:26 jjj  

言うなれば、

ある特定の人(id = 1)が2018/02/24 ~ 2018/02/26期間内にネットサーフィンしたログデータ

本来は、他のidがついてるもんですが、それはgroup_by()して処理すれば問題ないっす

ではこのデータに、セッションが付与すると...

> session <- c(1,1,2,2,3,3,4,4,5,5)
> umr %>% cbind(session) %>% select(action_log, session) %>% print()
select: dropped 2 variables (id, URL)
            action_log session
1  2018-02-24 15:30:15       1 # 始めのログは1
2  2018-02-24 15:40:34       1 # 変動なし
3  2018-02-24 16:30:19       2 # 30分以上経過しいるので次のセッション
4  2018-02-24 16:32:49       2 # 変動なし
5  2018-02-25 09:10:37       3 # 日付が変わったので次のセッション
6  2018-02-25 09:30:26       3 # 変動なし
7  2018-02-25 23:47:55       4 # 30分以上経過したので次のセッション
8  2018-02-25 23:59:44       4 # 変動なし
9  2018-02-26 00:02:28       5 # 日付が変わったので次のセッション
10 2018-02-26 00:04:26       5 # 変動なし

といった、感じで数字が割り振られるように自作関数を作れれば成功。

Let's go!!

セッション取得関数の作成

まず、1つめのセッションは1とします。

は?当然だろと思いますが、私はここで死んでました

function(action_log){
  res <- vector()
  a <- 1
  for (i in 1:NROW(action_log)) {
    if(i == 1){
      res[i] <- a
    }else {
        #二個目のログデータから処理するコード
        }
      }
  return(res)
}
条件1. 日付が変わると次のセッションへ

Rで書くとこんな感じ

function(action_log){
  res <- vector()
  a <- 1
  for (i in 1:NROW(action_log)) {
    if(i == 1){
      res[i] <- a
    }else {
##################
        #時間まである日付データをymdに変換
        #対象をymdに変換
      date <- action_log[i] %>% lubridate::date()
        #対象の一個前をymd変換
      pre_date <- action_log[(i-1)] %>% lubridate::date()
      if(pre_date != date){
         #対象と対象の一個前の日付が違うなら次のセッションへ
        a <- a+1
        res[i] <- a
      }else{
        #対象と対象の一個前の日付が同じとき
        #対象のログが前回のログから30分以上経過していたら次のセッションへ処理するコード
        }
      }
    }
  }
  return(res)
}
  • lubridate::date() ... 日付データをymd(年, 月, 日)に変換する関数(baseとのコンフリクトに注意)
条件2. 対象のログが前回のログから30分以上経過していたら次のセッションへ

この条件について書き下すと、

(対象の一個前のログデータ+30m) - (対象ログデータ) < 0 → 次のセッション

(対象の一個前のログデータ+30m) - (対象ログデータ) > 0 → 同一セッション

Rで書くとこんな感じ

function(action_log){
  res <- vector()
  a <- 1
  for (i in 1:NROW(action_log)) {
    if(i == 1){
      res[i] <- a
    }else {
      date <- action_log[i] %>% lubridate::date()
      pre_date <- action_log[(i-1)] %>% lubridate::date()
      if(pre_date != date){
        a <- a+1
        res[i] <- a
      }else{
##################
        #対象の一個前ログデータ+30mと対象ログデータの差分を取る
        dif <- (action_log[i-1] + dminutes(30))- action_log[i]
        #time_length()で差を分に変換。負の値ならTRUE
        judge <-time_length(dif, "minute") < 0
        if(judge == FALSE){
          #差が30分未満なら同一セッション
          res[i] <- a
        }else{
          #差が負の値なら次のセッションへ
          a <- a + 1
          res[i] <- a
        } 
      }
    }
  }
  return(res)
}
  • dminutes() ... 引数を分に変換する関数
  • time_length() ... 2つの日付データの間隔を割り出す関数。第二引数(unit)に単位(今回は"minute")を入れる

(あとから、気づきましたが「judge = o」は「FALSE = FALSE」なんでちょっと手間な処理かも)

完成した関数(get_session())はこちら

# セッション情報取得関数
get_session <- function(action_log){
  res <- vector()
  a <- 1
  for (i in 1:NROW(action_log)) {
    if(i == 1){
      res[i] <- a
    }else {
      date <- action_log[i] %>% lubridate::date()
      pre_date <- action_log[(i-1)] %>% lubridate::date()
      if(pre_date != date){
        a <- a+1
        res[i] <- a
      }else{
        dif <- (action_log[i-1] + dminutes(30))- action_log[i]
        judge <-time_length(dif, "minute") < 0    
        if(judge == 0){
          res[i] <- a
        }else{
          a <- a + 1
          res[i] <- a
        }
      }
    }
  }
  return(res)
}

では実行してみましょう。

sessionに(1,1,2,2,3,3,4,4,5,5)と入っていれば成功です

> umr %>% mutate(session = get_session(action_log)) %>% print()
mutate: new variable 'session' with 5 unique values and 0% NA
# A tibble: 10 x 4
      id action_log          URL   session
   <dbl> <dttm>              <fct>   <dbl>
 1     1 2018-02-24 15:30:15 aaa         1
 2     1 2018-02-24 15:40:34 bbb         1
 3     1 2018-02-24 16:30:19 ccc         2
 4     1 2018-02-24 16:32:49 ddd         2
 5     1 2018-02-25 09:10:37 eee         3
 6     1 2018-02-25 09:30:26 fff         3
 7     1 2018-02-25 23:47:55 ggg         4
 8     1 2018-02-25 23:59:44 hhh         4
 9     1 2018-02-26 00:02:28 iii         5
10     1 2018-02-26 00:04:26 jjj         5

やったぜ!!!!!!!!!

こんな感じで、Rでセッション情報を取得することが出来ました。

めでたし、めでたし

今回は、Googleアナリティクスさんの定義の「3. 参照元~」以外の条件を使って実装しましたが、もし参照元の情報があれば、対象と対象の直前を{stringr}とかで一致判定させれば、容易に厳密なセッション情報取得ができると思います。


別解

追記

ブログを投稿したところ、ネット上のとある猛者殿からよりスマートなコードを頂きました!

許可を頂いたので載せます。

> umr %>%
+   mutate(
+     d = as.numeric(
+       as.Date(action_log) - lag(as.Date(action_log), default = 0)),
+     m = as.numeric(
+       action_log - lag(action_log, default = 0)),
+     flag = d | (m > 30), session = cumsum(flag)) %>%
+   select(-d, -m, -flag)
mutate: new variable 'd' with 3 unique values and 0% NA
mutate: new variable 'm' with 10 unique values and 0% NA
mutate: new variable 'flag' with 2 unique values and 0% NA
mutate: new variable 'session' with 5 unique values and 0% NA
select: dropped 3 variables (d, m, flag)
# A tibble: 10 x 4
      id action_log          URL   session
   <dbl> <dttm>              <fct>   <int>
 1     1 2018-02-24 15:30:15 aaa         1
 2     1 2018-02-24 15:40:34 bbb         1
 3     1 2018-02-24 16:30:19 ccc         2
 4     1 2018-02-24 16:32:49 ddd         2
 5     1 2018-02-25 09:10:37 eee         3
 6     1 2018-02-25 09:30:26 fff         3
 7     1 2018-02-25 23:47:55 ggg         4
 8     1 2018-02-25 23:59:44 hhh         4
 9     1 2018-02-26 00:02:28 iii         5
10     1 2018-02-26 00:04:26 jjj         5

おいおい瞬殺だよ...

強すぎる!スマートすぎる!!

僭越ながら、説明させていただきますと

umr %>%
  mutate(
      #dを数値型に変換
    d = as.numeric(
        #action_logをdttm型からdate型に変換
        #dplyr::lag()で一つ後のログにずらしたものと差をとる
      as.Date(action_log) - lag(as.Date(action_log), default = 0)),
    m = as.numeric(
        #dttm型で一つずらしたものを差をとる
      action_log - lag(action_log, default = 0)),
      #論理式(dが1以上または、mが30以上の時TRUE(1))
      #session列にcumsum()でflagの累積和をとる
    flag = d | (m > 30),
    session = cumsum(flag)) %>%
  select(-d, -m, -flag)
  • dplyr::lag()...データを後ろにずらす。defaultでずらして空いた要素に何を入れるか指定(ここでは0)
  • base::as.Date()...データを日付(date)型に変換
  • dplyr::cumsum()...引数に指定したベクトルの累積和をとる

dplyr::lag(default = 0)にしてるのがみそだと思います。

引数に何も指定しないとNAがはいります。

flagで論理式を適用するのに、d,mが両方ともNAだとエラーが発生します。

それに、0を入れることで一つ目の要素とlagとの差がないので、

sessionに一行目のflagTRUE(1)にすることが出来ます。

そして、累積和(cumsum())を使うことでif文無しでsessionを作り出せます。

(累積和は、処理を早めるのに有効らしい)

一応、分かりやすいように先ほどのコードを書くと

> umr %>%
+   mutate(
+     lag = lag(action_log, default = 0),
+     d = as.numeric(
+       as.Date(action_log) - lag(as.Date(action_log), default = 0)),
+     m = as.numeric(
+       action_log - lag(action_log, default = 0)),
+     flag = d | (m >30),
+     session = cumsum(flag)) %>%
+   select(action_log,lag,d,m,flag,session)
# A tibble: 10 x 6
   action_log          lag                     d           m flag  session
   <dttm>              <dttm>              <dbl>       <dbl> <lgl>   <int>
 1 2018-02-24 15:30:15 1970-01-01 00:00:00 17586 25324770.   TRUE        1
 2 2018-02-24 15:40:34 2018-02-24 15:30:15     0       10.3  FALSE       1
 3 2018-02-24 16:30:19 2018-02-24 15:40:34     0       49.8  TRUE        2
 4 2018-02-24 16:32:49 2018-02-24 16:30:19     0        2.5  FALSE       2
 5 2018-02-25 09:10:37 2018-02-24 16:32:49     1      998.   TRUE        3
 6 2018-02-25 09:30:26 2018-02-25 09:10:37     0       19.8  FALSE       3
 7 2018-02-25 23:47:55 2018-02-25 09:30:26     0      857.   TRUE        4
 8 2018-02-25 23:59:44 2018-02-25 23:47:55     0       11.8  FALSE       4
 9 2018-02-26 00:02:28 2018-02-25 23:59:44     1        2.73 TRUE        5
10 2018-02-26 00:04:26 2018-02-26 00:02:28     0        1.97 FALSE       5

lag(default = 0)にすることで、date型の最初の値である1970/01/01になります。

ほんで、flagdになにか要素があるか、mが30以上ならTRUEになります

dmがdbl型なのはas.numeric()のおかげ)

前の記事でも挙げたように、TRUEは1、FALSEは0を示します。

そのため、累積和をとるとsession列が完成します

美しすぎる~!!!!!!!!

最後に私の自作関数と猛者どのコードどちらが強い(速い)かバトルします。

> library(tictoc)
> tic()
> library(tidyverse)
-- Attaching packages --------------------------------------- tidyverse 1.2.1 --
√ ggplot2 3.1.1       √ purrr   0.3.2  
√ tibble  2.1.1       √ dplyr   0.8.0.1
√ tidyr   0.8.3       √ stringr 1.4.0  
√ readr   1.3.1       √ forcats 0.4.0  
-- Conflicts ------------------------------------------ tidyverse_conflicts() --
x dplyr::filter() masks stats::filter()
x dplyr::lag()    masks stats::lag()
> toc()
8.28 sec elapsed

{tictoc}を使って、処理速度を測ります。

ちなみに、僕のPCでlibrary(tidyverse)をすると8.28秒かかります。

(おじいちゃん...)

行くぜ!!!!デゥエル(決闘)!!!!!!

決闘

ルール

  • kazutan先生の1万行の仮想ログデータ(df_log)を使う
  • df_log %>% (session処理)で行う
  • 上記のコードをtic()toc()で処理時間測定
  • 個別idがあるけど無視
  • tail()で最後のsessionナンバーが同じか確認

俺のターン!get_session()を召喚!!!

 > tic()
> df_log %>%mutate(session = get_session(action_log)) %>% tail()
mutate: new variable 'session' with 201 unique values and 0% NA
               action_log      id   item session
9995  2018-02-20 23:24:34 1000032 item_4     201
9996  2018-02-20 23:37:01 1000068 item_1     201
9997  2018-02-20 23:40:34 1000121 item_2     201
9998  2018-02-20 23:55:01 1000152 item_1     201
9999  2018-02-20 23:55:14 1000144 item_1     201
10000 2018-02-20 23:59:39 1000267 item_2     201
> toc()
51.89 sec elapsed

処理時間51.89秒、最終セッション201

これを実際のデータベースにあるログデータに使うのは結構きついか...

猛者殿のターン!7行の処理コードを召喚!!!

> toc()
1.69 sec elapsed
> tic()
> df_log %>%
+   mutate(
+     action_log = ymd_hms(action_log),
+     d = as.numeric(
+       as.Date(action_log) - lag(as.Date(action_log), default = 0)),
+     m = as.numeric(
+       action_log - lag(action_log, default = 0)),
      #なぜか秒換算になったので30分を秒変換
+     flag = d | (m > 1800), session = cumsum(flag)) %>%
+    tail()
               action_log      id   item d   m  flag session
9995  2018-02-20 23:24:34 1000032 item_4 0 342 FALSE     201
9996  2018-02-20 23:37:01 1000068 item_1 0 747 FALSE     201
9997  2018-02-20 23:40:34 1000121 item_2 0 213 FALSE     201
9998  2018-02-20 23:55:01 1000152 item_1 0 867 FALSE     201
9999  2018-02-20 23:55:14 1000144 item_1 0  13 FALSE     201
10000 2018-02-20 23:59:39 1000267 item_2 0 265 FALSE     201
> toc()
1.56 sec elapsed

処理時間1.56秒!!!!!!!!!!!!!!!!!!!!!!!!!!!!(セッションは同じ)

俺の人生...

いや、すごいっすね。コードの長さも処理時間も段違いすぎる...

{lubridate}だけでなく、lag(),cumsum()を知れたのは非常にためになりました!

ありがとうございます!!

to be continued...

R論理演算子の数とmutate()での条件分岐

初めましてのブログです。最近はダーリン・イン・ザ・フランキス見てます。

今回は、

Rにおける理論演算子(&,|)の数によるエラーtidyverse::mutate()を用いた条件分岐 をいくらか紹介するコーナー

昨日まで、私はRにおける論理演算子である& ,|は二つ並べて使うものだと思っていました。そこでこんなエラーと出会います。

まず、こんなデータセットがあります。

R version 3.5.3 (2019-03-11) -- "Great Truth"
Copyright (C) 2019 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64 (64-bit)

>> library(tidyverse)
> umr <- data.frame(
+   ID = as.factor(1:18),
+   method = as.factor(c("A","B","B", "A","A","A","A","A",
+                        "B","B","B","B","A", "B","B","B","B","B")),
+   critelia = c(1,0,0,1,1,1,1,1,1,0,0,0,1,1,1,1,1,1)
+ ) %>% 
+   as.tibble %>% print()
# A tibble: 18 x 3
   ID    method critelia
   <fct> <fct>     <dbl>
 1 1     A             1
 2 2     B             0
 3 3     B             0
 4 4     A             1
 5 5     A             1
 6 6     A             1
 7 7     A             1
 8 8     A             1
 9 9     B             1
10 10    B             0
11 11    B             0
12 12    B             0
13 13    A             1
14 14    B             1
15 15    B             1
16 16    B             1
17 17    B             1
18 18    B             1

ID : 参加者ID

method : 試した方法(例えば、2種類のダイエット方法)

critelia : 成功したか(例えば、1は成功、0は失敗)

※架空のデータセット


このデータに対して、以下のような条件分岐を行いたい

method == A かつcritelia == 1なら新たな列に`1

または

method == Bかつcritelia == 0なら新たな列に1

それ以外は0 を入れる

新たな列をclear_orderと名付け、で私はこうした。

> umr %>%
+   dplyr::mutate(
+     clear_order = ifelse((method == "A" && critelia == 1) ||
+                            (method == "B" && critelia == 0), 1, 0)
+   ) %>%print()
# A tibble: 18 x 4
   ID    method critelia clear_order
   <fct> <fct>     <dbl>       <dbl>
 1 1     A             1           1
 2 2     B             0           1
 3 3     B             0           1
 4 4     A             1           1
 5 5     A             1           1
 6 6     A             1           1
 7 7     A             1           1
 8 8     A             1           1
 9 9     B             1           1
10 10    B             0           1
11 11    B             0           1
12 12    B             0           1
13 13    A             1           1
14 14    B             1           1
15 15    B             1           1
16 16    B             1           1
17 17    B             1           1
18 18    B             1           1

全部1やん...

どうやらこれは、理論演算子が関係しているようだ。

詳しくは、以下の記事を参照↓↓

Rの理論演算子&と&&、|と||の違い

要するに
理論演算子1つで要素ごとに演算を行い、2つだとベクトルの1つ目の要素のみを演算する

つまり、今回のエラーは理論演算子を1つにすれば解決する。

Let's go!!

> umr %>%
+   dplyr::mutate(
+     clear_order = ifelse((method == "A" & critelia == 1) |
+                            (method == "B" & critelia == 0), 1, 0)
+   ) %>% print()
# A tibble: 18 x 4
   ID    method critelia clear_order
   <fct> <fct>     <dbl>       <dbl>
 1 1     A             1           1
 2 2     B             0           1
 3 3     B             0           1
 4 4     A             1           1
 5 5     A             1           1
 6 6     A             1           1
 7 7     A             1           1
 8 8     A             1           1
 9 9     B             1           0
10 10    B             0           1
11 11    B             0           1
12 12    B             0           1
13 13    A             1           1
14 14    B             1           0
15 15    B             1           0
16 16    B             1           0
17 17    B             1           0
18 18    B             1           0

tidyverse::mutate()を使って要素を一つずつ評価したいときは&|の理論演算子は1つで行うべきなんですね。学びが深いなぁ...

mutate(new = if_else(条件A │ 条件B, 1, 0))

ちなみに、今回このエラーをとある秘密集団に投げかけたところ、別解を得られたので紹介。


別解1. case_when()を使ったパターン

> umr %>% 
+   dplyr::mutate(
+     clear_order = case_when(
+       critelia == 0 & method == "A" ~ 1,
+       critelia == 0 & method == "B" ~ 1,
+       TRUE ~ 0
+     )
+   ) %>% print()
# A tibble: 18 x 4
   ID    method critelia clear_order
   <fct> <fct>     <dbl>       <dbl>
 1 1     A             1           0
 2 2     B             0           1
 3 3     B             0           1
 4 4     A             1           0
 5 5     A             1           0
 6 6     A             1           0
 7 7     A             1           0
 8 8     A             1           0
 9 9     B             1           0
10 10    B             0           1
11 11    B             0           1
12 12    B             0           1
13 13    A             1           0
14 14    B             1           0
15 15    B             1           0
16 16    B             1           0
17 17    B             1           0
18 18    B             1           0

case_when()も複数の条件を指定することが出来きる

~の左辺が論理式(条件)

~の右辺に置き換える値

って感じ。今回のケースは比較的論理式のコードが短いから良いが、コードが長くなるようであれば、こっちのほうが可読性が高い気がする。

case_when()およびmutate()については以下を参照↓↓

列の変換 - mutate関数


別解2. めっちゃ短いパターン

> umr %>% mutate(clear_order = ((method == "A") == critelia) %>% 
+                as.integer()) %>% print()
# A tibble: 18 x 4
   ID    method critelia clear_order
   <fct> <fct>     <dbl>       <int>
 1 1     A             1           1
 2 2     B             0           1
 3 3     B             0           1
 4 4     A             1           1
 5 5     A             1           1
 6 6     A             1           1
 7 7     A             1           1
 8 8     A             1           1
 9 9     B             1           0
10 10    B             0           1
11 11    B             0           1
12 12    B             0           1
13 13    A             1           1
14 14    B             1           0
15 15    B             1           0
16 16    B             1           0
17 17    B             1           0
18 18    B             1           0

おいおい、瞬殺だよ...

これは短い。知らなかったけど、RはTRUEを1、FALSEを0の整数として扱うんですね。

  1. (method == "A")"A"ならTRUE(1), "B"ならFALSE(0)に変換
  2. criteliaは元から0,1データなので合致か否かでclear_order行にTRUE,FALSEを返す。
  3. as.integer()TRUE,FALSEを1,0に変換

as.integer()を使わないとTRUE,FALSEで返ってくるよ

ナイス、スマート


別解3. 自作関数を使ったパターン

> judge_order <- function(method, critelia) {
+   #空のベクトルを生成
+   res <- vector()
+   #要素ごとに評価するためにfor文
+   for (i in 1:length(method)) {
+     if (method[i] == "A" & critelia[i] == 1) {
+       res[i] <- 1
+     } else if (method[i] == "B" & critelia[i] == 0) {
+       res[i] <- 1
+     } else {
+       res[i] <- 0
+     }
+   }
+   return(res)
+ }
> umr %>% mutate(clear_order = judge_order(method, critelia)) %>% print()
# A tibble: 18 x 4
   ID    method critelia clear_order
   <fct> <fct>     <dbl>       <dbl>
 1 1     A             1           1
 2 2     B             0           1
 3 3     B             0           1
 4 4     A             1           1
 5 5     A             1           1
 6 6     A             1           1
 7 7     A             1           1
 8 8     A             1           1
 9 9     B             1           0
10 10    B             0           1
11 11    B             0           1
12 12    B             0           1
13 13    A             1           1
14 14    B             1           0
15 15    B             1           0
16 16    B             1           0
17 17    B             1           0
18 18    B             1           0

実は、最初は自作関数を使ってトライして失敗しました...

詳しくは分からないがtidyval案件で、難しいみたい

to be continued...