読書要約 道徳の系譜学 第一論文

思想 読書

本記事は、ニーチェの道徳の系譜学のうち、第一論文に限って要約したものです。

問いの設定

本論文では、今日（19世紀末当時のヨーロッパ）に行われている「善いとか悪いとかの価値判断」ということが慣行として根付くためには、どのような経緯が必要であったのかを明らかにしようとします。

研究のアプローチは本のタイトルにもある通り「系譜学」の方法をとります。ニーチェが本論文で取ったアプローチは文献の蓄積から解を発見しようとする点で歴史研究の一種です。どのような点にあえて系譜学と呼ぶべき視点が含まれているのかというと、「善という価値が、一体誰によって、どのような歴史的経緯を経て設定されたのか？」という政治学的な切り口によって、価値判断の仕方の変遷に着目した点です（中山元 2005, 第4章 ニーチェの系譜学）。

ニーチェの主張

今日に行われているキリスト教的な価値判断は、強き支配者に対する否定的な感情にルーツがあり、妬みの反動から弱者としての自分を「善人」と自称することで「善さ」が成立したというのがニーチェの主張です。

この結論の論拠として、ニーチェは「よい」という言葉が史上で登場し始めてから今日に至るまでの系譜を紐解きます。 はじめ、「良い」という言葉は貴族が自らの支配者としての高潔さを自画自賛して使った言葉でした。ところが、支配者に対する妬みから貴族的な支配者を「悪人」とし、そのような悪から対置されるところの弱者性を「善さ」として尊ぶ価値基準が登場しました。このような価値の転換をやってのけたのがキリスト教であり、今日の価値判断の主流を占めているーー以上がニーチェが文献から導き、論拠とするところの洞察です。

疑問点

今後更に本書を読み深めたり、他の文書にあたるときに備えて、ここでは疑問に思った箇所を書き記しておきます。

本論文の後半でニーチェはキリスト教的な道徳に基づく価値判断を辛辣に批判します。例えばフランス革命を例に取って次のように述べています。

ユダヤはさらにフランス革命において、ふたたび古典的な理想に勝利を収めたが、これはさらに深く、決定的な意味を持つ出来事だった。ヨーロッパに存在していた最後の政治的な高貴さが、すなわち十六世紀と十七世紀というフランスの世紀の政治的な高貴さが、民衆的なルサンチマンの本能のもとで崩壊したのである

私が疑問に思ったのは、ニーチェによる道徳批判の論拠です。例えば貴族的な独特を「善さA」と呼び、奴隷的な道徳を「善さB」と呼ぶとしたなら、「AとBを比較する議論」はどのような観点に基づいて行うべきか？「『善さ』のよさ」をどうやったら評価できるのか？

参考文献

CtoCで協調フィルタリングした記録〜訓練データ作り・評価指標・実装〜

python 論文

以前、協調フィルタリングをお仕事で使う機会があったので、その時の記録を書きます。

訓練データづくりで注意したこと

「user-item行列においてuser数 が item数に対して十分に多くなるようにする」ということ*1に注意しました。

f:id:sakatoken:20190627230134p:plain
左側のような形状の行列の方がユーザ選好を発見しやすい

userとかitemとか言っておりますが、これらの単語は次のような意味で使っています。

  • item：購入履歴や星5段階ratingなど、評価を下されるところの商品などを指すものです。
  • user：itemに対して評価を下すところの主体に相当するユーザを指すものです。
  • user-item行列：この行列の i ,j成分は、user_i が item_j に下した評価を表します。

例えば「10人のユーザいて、10種類の商品をratingした。ただし１人につきratingは１件しか下さず、なおかつ全てのユーザが相互に異なる商品に対してratingをつけた。」のようなuser-item行列を想定してみます。協調フィルタリングではuser-item行列に対して種々のmatrix-factorizationを用いますが、このように「user同士がitemへの選好を分かち合っていない」ような行列からは、潜在的なuserベクトルとitemベクトルを発見することはできません。

より具体的には、訓練データを準備する際には以下の要件を満たしたい。

  • user数が十分にitem数よりも多い
  • userは十分に多くのitemに評価を下している

特に商品在庫が「一点もの」になりやすいCtoCのwebサービスでは、この問題が起こりやすいと思います。（私はそうなりました）

この問題を解決するために、私が行なった対策は以下です。

  1. itemの単位を「商品id」ではなく「商品細分類」に変える。つまり特徴(テキストや画像など)に基づくクラスタリングを事前に行うということ。
  2. user1人あたりの評価件数に閾値を設けて、ほんの少ししか評価していない人は訓練データから除外する。（閾値の適切さはケースバイケースだと考えています）

私の事例では、協調フィルタリングの訓練段階よりも、むしろ上記の商品細分類を発見する作業の方により多くの労力を使いました。こういう風にitemの本質を量的に表現することはレコメンド以外のタスクにも活きると思うので、おそらく色々な局面で「非構造データからitem特徴を入手する」ようなアプローチが重要なんだと思います。

下記のwantedlyのrei kubonaga さんの事例（すごくイケてる！）では、DNNを使って「userとitemの相性の良さを予測する」というアプローチを取られていますが、これも特徴抽出ありきの話ですね。

レコメンデーションのターゲットメトリックスの設定 / ML loft 3 - Speaker Deck

性能評価で注意したこと

オフライン評価（サービスにリリースしたときのビジネス成果による評価ではなく、開発時のテストスコア）におけるレコメンドエンジンの評価方法について、Collaborative Filtering for Implicit Feedback Datasets(2008)という論文が大変示唆に富んでいると思います。（少なくとも私にとっては。もしかしたら2008年よりずっと前から指摘され続けてる話かもしれませんが...）

これを参考に以下の２点に注意しました。

1. テストデータの用意の仕方を実際の問題設定に近づける

この論文では「4週間分のテレビ視聴ログからユーザの選好を学習して、推薦によって次の１週間の視聴内容を言い当てる」という実験を行なっています。(下記は該当箇所の引用）

we use a similarly constructed test set, which is based on all channel tune events during the single week following a 4-week training period. Our system is trained using the recent 4 weeks of data in order to generate pre-dictions of what users will watch in the ensuing week.

単にデータの一部を（色々と条件づけた上で）ホールドアウトすることでテストデータとすることも出来るわけですが、テレビ番組の視聴データで言えば「データ収集期間の途中で新番組が始まった」のようなケースを考えると、未来のデータを使って起こりえない過去を予測するような学習が起こり得ます。テストデータにはあくまで未来のデータを用いた方が、より本番に近い状態で実験が出来ると思います。

2. precisionあまり気にしない(問題設定によるけど)

この実験では「推薦の順位と、実際のユーザ選好の順位の一致度」のようなものを評価指標とするわけですが、この評価指標に関して次のように言及しています。

precision based metrics are not very appropriate, as they require knowing which programs are undesired to a user. However, watching a program is an indication of liking it, making recall-oriented measures applicable.

「テレビ視聴データではユーザが何を好きか(いわば陽性)は分かるが、何を嫌いか(陰性)は分からない。よってprecisionよりもrecall志向で評価する方が適切。」という旨のことを言っています。

例えばこの論文では「あるテレビ番組を視聴していないということは、単に『好きな番組の裏番組だから見れなかった』だけかもしれないし、見ていないというだけで積極的に『嫌い』の烙印を押すのは誤りだ」ということに言及していて、これは言い換えると「本当はレコメンドが正しいのに、ユーザが視聴ログを残してくれなかった」というケースがありえることを意味します。上記で触れたrecall志向とは「推薦されたitemが正しかろうと、実際のuserの視聴ログは歯抜けになっちゃうのだから、推薦上位にやってくるitemが十分なrecallを持つのかどうかに着目しましょう」ということです。

前述のテストセットの用意の仕方と同じく、レコメンドエンジンの使途目的に応じて評価指標を選ぶと、納得のいく実験をデザインできますね。
その上で「具体的に評価指標にどんな尺度を使うのか？」についてはbrainpadさんの↓のブログが大変参考になります。

blog.brainpad.co.jp

実装時に使ったライブラリ・参考にした記事など

私が扱ったデータは「クリック数をuserによるitem選好として解釈する」ようなものでした。これを協調フィルタリングの界隈ではimplicitなデータと呼ぶようです。

そこで既に本記事でたびたび引用しているCollaborative Filtering for Implicit Feedback Datasets(2008)で紹介されている手法を扱うことのできる"implicit"というそのものズバリな名前のpythonライブラリを使うことにしました。なおこのライブラリを使う際には、元の論文で登場するαパラメータは手動でチューニングする必要がある点に注意です。

github.com

さらにimplicitを使ったサンプル実装としてこちらの記事が大変参考になりました。この記事が素晴らしいのは、グリッドサーチの実装を公開してくれている点です(下記)。αパラメータ(下記実装ではconfidenceという変数) のチューニングも含めた内容となっています。

Recommending GitHub repositories with Google Big Query and implicit library: https://medium.com/@jbochi/recommending-github-repositories-with-google-bigquery-and-the-implicit-library-e6cce666c77 · GitHub
※implicitのAPI仕様書でcsr_matrixが指定されている箇所で、何故かcoo_matrixを使っている点だけひっかかります。

仕上げにテストスコアを「単純ランキング」と比較する

前節でのグリッドサーチを行うことでベストなパラメータは求まりますが、これはあくまで協調フィルタリングという枠の中における相対的なベストに過ぎません。つまり、ここで言うベストとは推薦内容が受け取り手のユーザにとって優れていることを意味しません。よって、推薦手法の対立候補として「単純な人気ランキングに基づく推薦」などを用意し、テストスコアを比較することで「パーソナライズに意義がある！」という主張のウラを取る必要があります。

（さらに、オフライン評価での成功はただちにKPIへの寄与を意味しないので、この後にはオンライン評価が控えています。）

*1:これはquoraの投稿：https://www.quora.com/How-do-I-speed-up-matrix-factorization-by-sampling-users-without-losing-precision を参考にしました。元ネタの文献が何なのかがわからず。。。

時系列クラスタリング手法のk-shapeの論文を一部読んでみた

統計 時系列解析 論文

以前の記事では適当にダミーデータを用意して、愚直にユークリッド距離 & k-meansによる時系列データのクラスタリングを試みましたが、やり残した感があったのでもう少しインプットを深めてみました。今回は*1k-shapeというクラスタリング手法の論文の一部に目を通してみました。k-shapeの基礎を知るとともに、時系列データクラスタリング一般での注意点に気づきがあったので、これを記録します。

元論文はここ

k-shapeのざっくりとしたイメージ

下図の矢印のような「位相と振幅を揃えてあげれば、似ている波だよね」という直感を反映した距離尺度(後述するShape-based distance)を使ってクラスタリングする手法です。

f:id:sakatoken:20190611232753p:plain

原論文2.6節では次のように述べられています。

we focus on distance measures that offer invariances to scaling and shifting, which are generally sufficient

距離尺度の重要性

元論文のintroductionでは次のような記述があります。

more attention has been given to the creation of new distance measures rather than to the creation of new clustering algorithms

「時系列データのクラスタリングの研究においては、クラスタリング手法よりもその前段の距離尺度の定義に労力が割かれてきた」ということを言っています。前述の「位相と振幅を揃えてあげれば、似ている波だよね」のような直感をデータで上手いこと表現するには、位相と振幅に影響されない距離尺度の定義が重要だということですね。

言い換えると、位相と振幅のズレを無視したくない場合には、k-shapeは不適だということです。「どのような類似性に基づいてクラスタを抽出したいのか？」という問題設定に応じて距離尺度を選ぶべき。

「似ている」とは何か？ーーinvarianceの観点

現実の問題で時系列データに類似性を見出すときには、「位相と振幅を揃えてあげれば、似ている」のように「もし（ 因子 ）さえ無ければ、似ているのになぁ」という考え方をすることになります。つまりノイズや歪みの不在を仮想する洞察が要求されます。

では、どのような因子に対してその不在を仮想すべきか？これは問題設定に依存しますが、本論文ではいくつかの選択肢を提示してくれています。

  • scaling and translation(amplitude and offset)：振幅とオフセットの差異
  • shift：位相の差異
  • uniform scaling：一つ一つの時系列データの長さ

などなど。*2

これらの「クラスタリングの障害となりうるもの」を指して本論文ではdistortionと呼び、またそれらに対して適切な前処理なり距離尺度なりを施した状態を指してinvarianceと呼んでいます。（こういった用語が時系列分析の界隈で一般的なのかは存じません。）

まず問題設定がありきで「何をdistortionとみなすか？」を決めて「ゆえにどんなinvarianceが必要か？」を検討するーーそのように考えると時系列クラスタリングの様々な手法の選択について、スッキリした視点を持つことが出来るのではと思います。

Shape-based distanceという尺度

さて、冒頭で述べたk-shapeの元にある「位相と振幅を揃えてあげれば、似ている波だよね」という直感を前節での用語で言い直すと、

  • 位相と振幅をdistortionとして扱いたい

となります。この観点からinvarianceを実現するための距離尺度とは何か？

k-shapeの回答はShape-based distance(SBD)という距離尺度です。*3二本の時系列データのベクトルをxとyとおくと、これらのSBD距離は以下の通りです。

\displaystyle{ SBD(\boldsymbol{x}, \boldsymbol{y}) = 1 - \underset{w}{\max}\left( \frac{CC_w(\boldsymbol{x}, \boldsymbol{y})}{\sqrt{R_0(\boldsymbol{x}, \boldsymbol{x})\cdot{}R_0(\boldsymbol{y}, \boldsymbol{y})}} \right) }

最大化されるところのCC_wとは、xとyの相互相関関数を意味しています。タイムラグwについてCC_wの最大化を検討することは「位相のズレに引っ張られずにxとyの相関を測る」ことに相当するわけです。

前回の記事を書くときに参考にした以下のブログではDTWに触れていらっしゃるのですが、タイムラグへの配慮というのは時系列クラスタリングで共通するテーマなんですね。

sinhrks.hatenablog.com

おわりに

以前の記事に関連して、新たに実装して記事にしてみようと思います。

*1:最初にk-shapeを選んだ動機として、こちらのブログをを見て論文が面白そうだと思ったというのがありますhttps://blog.tsurubee.tech/entry/2019/02/06/223555

*2:他にもocclusionとcomplexityという概念が並んでいましたが、私自身が波の解析に関する用語に疎いので訳をあてるのは避けます。これらを省くことはk-shapeの概要を述べるのには差し支えありません。

*3:振幅に関しては、k-shapeにおいては距離尺度の計算に先立って単に標準化スケーリングを行うことで各時系列を揃えます。

裾の長いデータに対数正規分布を仮定してベイズモデリング

R stan ベイズ 統計

顧客１人あたりの売上やら何やら、本当に様々なデータが裾の長い度数分布の形状を持ちます。
そういったデータに対して何か説明変数を設定してモデリングをする必要に迫られるシーンを想定して、下図のごときデータにフィットするモデリングを考えてみます。本記事では目的変数をyとおき、説明変数(単変量)をxとおきます。

f:id:sakatoken:20190611000403p:plain

ダミーデータ生成など

ダミーデータ生成とライブラリ読み込みをしておきます。

library(tidyverse)
library(GGally)
library(rstan)
library(tidybayes)

d <- tibble(
    x = abs(rnorm(2000, mean =1.5))
) %>%
    mutate(x = x / max(x)) %>%
    mutate(
        mu = -4 * x,
        sigma = exp(-0.5*x -0.5)
    ) %>%
    mutate( 
        y = as.numeric(map2(mu, sigma, function(m, s){rlnorm(1, meanlog = m, sdlog = s)}))
    )

分布を仮定する

この問題を考えるには「裾の長い非負の度数分布を見たときにどのような分布を仮定するべきか？」に明快に答えを出したいところですが、今回は大変安直ながら指数分布と対数正規分布を見比べて

  1. xに応じて平均と分散の両方をモデリングできる
  2. xに応じて分布の峰が0から離れる分布のほうが、ビジネスの実問題にはフィッティングさせやすいと思われる

という理由で対数正規分布を選びました。

※ダミーデータでは対数正規分布が真のモデルになっています。

モデルの記述

xと対数正規分布のパラメータの間には次の仮定を設定しました。

  • 対数正規分布のパラメータには平均と分散があります。今回はそれらの両方に対し、xによる線形予測子を仮定します。
  • 問題によりますが、ここではxは実数全体を取るため、分散（非負）と結びつけるにあたって指数を使います。
  • 線形予測子に含む係数部分には無情報事前分布を仮定します。

stanコードは次の通りです。単変量であるはずのxにmatrixを指定しているのは、Rコード側で切片項を追加するためです。

data {
    int N;
    int D;
    vector<lower=0>[N] Y;
    matrix<lower=0, upper=1>[N, D] X;
}

parameters {
    vector[D] b;
    vector[D] a;
}

transformed parameters {
    vector<lower=0>[N] s;
    vector[N] mu;
    
    s = exp(X * b);
    mu = X * a;
}

model {
    Y ~ lognormal(mu, s);
}

前処理とMCMC実行

以上でモデルを定義できたので、前処理してMCMCを実行する部分のRコードを記述していきます。

# 外れ値 除去
d.f <- d %>%
    filter(
        x < quantile(x, .98), 
        y < quantile(y, .98)
    )

# 最大値を使ってスケーリング
d.f.max <- d.f %>%
    summarize_all(max) %>%
    select(x,y)
d.f <- d.f %>%
    mutate(
        x = x / d.f.max$x,
        y = y / d.f.max$y
    )

# stanにわたす用のlist型データを作成
d.stan <- list(
    N = nrow(d.f),
    Y = d.f$y,
    X = d.f %>%
        select(x) %>%
        mutate(intercept = 1) %>% # 切片項を追加
        as.matrix(),
    D = 2
)

# MCMCを実行
fit <- stan(file = './modeling.stan', data = d.stan, seed = 0,
            iter = 700, warmup = 200, chains = 3, pars=c('b', 'a'))

# 収束判定にはtraceplot とRhat値のチェックをするけどコード略。

モデルの実データへのフィット具合をグラフで確認

以上で推定したモデルを使ってyの予測値を新たに生成し、これが実測のyを言い当てるかどうかを確認します。ただし冒頭の図で確認した通り、今回は特定のxの値に対してyは広くばらつく分布を持っているため、予測値と実測値の誤差の小ささという基準でフィットを評価しようとすると「まるでデタラメなモデルである」という結果になります。

そこで「モデルは個別のyを言い当てることはできないが、xがある値に固定されたときのyの分布を言い当てることが出来る」という立場を取り、グラフを作成することとします。

具体的には以下の手順を取ります。

  1. 推定されたパラメータを使って、xからμとσを算出する。
  2. μとσにしたがって、新しいy_predを生成する。
  3. μの値でデータを等間隔の階級に区切る。
  4. 階級別にyとy_predのヒストグラムを重ねてプロットする。
  5. グラフでy_predが十分にyに追従できているのかどうかを確認する。

以下、コードなどです。

まずは先程のMCMCによる推定結果から、パラメータの値を控えておきます。

Inference for Stan model: modeling.
3 chains, each with iter=700; warmup=200; thin=1; 
post-warmup draws per chain=500, total post-warmup draws=1500.

       mean se_mean   sd   2.5%    25%    50%    75%  97.5% n_eff Rhat
b[1]  -0.30    0.00 0.07  -0.43  -0.34  -0.30  -0.25  -0.15   930    1
b[2]  -0.53    0.00 0.03  -0.59  -0.55  -0.53  -0.50  -0.46   906    1
a[1]  -2.71    0.00 0.05  -2.81  -2.75  -2.71  -2.68  -2.61   592    1
a[2]  -0.65    0.00 0.03  -0.70  -0.66  -0.65  -0.63  -0.59   657    1
lp__ 298.37    0.06 1.40 294.74 297.71 298.69 299.39 300.09   619    1

Samples were drawn using NUTS(diag_e) at Thu Jun  6 18:01:56 2019.
For each parameter, n_eff is a crude measure of effective sample size,
and Rhat is the potential scale reduction factor on split chains (at 
convergence, Rhat=1).
beta = c(-0.3, -0.53)
alpha = c(-2.71, -0.65)

続いて、グラフをプロットするまでのRコードです。

# μとσにしたがって、新しいy_predを生成する。
d.f.pred <- d.f %>%
    mutate(
        sigma = exp(x * beta[1] + beta[2]),
        mu = x*alpha[1] + alpha[2]
    ) %>%
    mutate(
        y_pred = rlnorm(n(), meanlog = mu, sdlog = sigma)
    ) %>%
    filter(y_pred < quantile(y_pred, .98)) # 外れ値 除去

# 階級別にyとy_predのヒストグラムを重ねてプロットする。
d.f.pred %>%
    mutate(class = cut_interval(mu, n = 8)) %>% # μの値でデータを等間隔の階級に区切る。
    select(y, y_pred, class) %>%
    gather("name", "value", y, y_pred) %>% 
    arrange(class) %>%
    mutate(name = if_else(name == "y", "実測", "予測")) %>%
    ggplot(
        aes(value, fill = name)
    ) +
    geom_histogram(bins = 15, position = "identity", alpha = 0.5) +
    facet_wrap(~class, ncol = 2, scales="free") +
    labs(
        x = "推定μ(対数正規分布)",
        fill = "予測と実測の色分け",
        title="推定μの階級別ヒストグラム"
    )

以上によって生成されたグラフが以下です。狙いどおりに「分布を言い当てる」ことに成功していることが見て取れます。私は実務で応用するに際しては、同じモデルでテストデータに対する予測の良さを見ておくなどしました。

f:id:sakatoken:20190611005635p:plain

参考文献

StanとRでベイズ統計モデリング (Wonderful R)

StanとRでベイズ統計モデリング (Wonderful R)

word2vecで専門用語バリバリのコーパスを扱うため、語彙をMeCabに覚えさせる

f:id:sakatoken:20190606165702p:plain
《渦巻く知識》というカード

mecabに流行語や固有名詞を覚えさせて分かち書きの精度を上げようと考えた時に、いつも大変お世話なっている辞書はipadic neologdです。
しかしneologdといえどもマニアックすぎる語彙は分かち書きしてくれません。

本記事ではmecab公式にて解説のあるユーザ辞書機能にマニアックな語彙を登録した記録を残します。

例：某カードゲームの戦術についてのテキストを分かち書きしてみる。

Magic: the Gatheringというカードゲームがあります。
題材として、このカードゲームの戦術に関する例文を分かち書きしてみることとします。

以下例文。

初手にフェッチランドと島と《渦巻く知識》を引いていた場合に、ひとまずワンマリの憂き目は逃れるとして、どの土地をセットするかの決断を迫られる。

専門用語満載です。これらのうち、固有名詞として捉えておきたい語彙は以下です。

  • フェッチランド
  • ワンマリ
  • 渦巻く知識

さて、これを分かち書きすると次のようになります。(名詞のみ列挙）
当然ながらマニアックな語彙を拾うことはできません。

初手   名詞,一般,*,*,*,*,初手,ショテ,ショテ
フェッチ    名詞,サ変接続,*,*,*,*,フェッチ,フェッチ,フェッチ
ランド   名詞,一般,*,*,*,*,ランド,ランド,ランド
島 名詞,一般,*,*,*,*,島,シマ,シマ
知識  名詞,一般,*,*,*,*,知識,チシキ,チシキ
場合  名詞,副詞可能,*,*,*,*,場合,バアイ,バアイ
ワン  名詞,一般,*,*,*,*,ワン,ワン,ワン
マリ  名詞,固有名詞,地域,国,*,*,マリ,マリ,マリ
憂き目   名詞,一般,*,*,*,*,憂き目,ウキメ,ウキメ
土地  名詞,一般,*,*,*,*,土地,トチ,トチ
セット   名詞,一般,*,*,*,*,セット,セット,セット
決断  名詞,サ変接続,*,*,*,*,決断,ケツダン,ケツダン

MeCabにマニアックな語彙を覚えてもらおう

MeCabのユーザ辞書という機能を使ってマニアックな語彙を覚えさせましょう。

先ほど上げた固有名詞をユーザ辞書に登録するため、以下の内容でcsvファイルを作ります。

フェッチランド,,,,名詞,固有名詞,一般,*,*,*,フェッチランド,フェッチランド,フェッチランド
ワンマリ,,,,名詞,固有名詞,一般,*,*,*,ワンマリ,ワンマリ,ワンマリ
渦巻く知識,,,,名詞,固有名詞,一般,*,*,*,渦巻く知識,ウズマクチシキ,ウズマクチシキ

そして、下記コマンドでユーザ辞書を作成します。

/usr/lib/mecab/mecab-dict-index \
-m /usr/share/mecab/dic/mecab-ipadic-2.7.0-20070801.model \ # 公式によると「モデルファイルの指定」をしている。よくわかっていません。
-d /usr/share/mecab/dic/mecab-ipadic-2.7.0-20070801 \ # システム辞書のあるディレクトリ
-u /お好きなディレクトリ/mtg_dic.dic \ # ユーザ辞書ファイルの出力先
-f utf-8 \
-t utf-8 \
mtg_dic.csv # 先ほどのcsvファイル

このステップで参考にしたのは以下の記事です。

ユーザ辞書を使って分かち書きリトライ

というわけで、先の例文を改めて分かち書きします（名詞のみ抽出）。期待通り、マニアックな固有名詞をつかまえることができています。

初手にフェッチランドと島と《渦巻く知識》を引いていた場合に、ひとまずワンマリの憂き目は逃れるとして、どの土地をセットするかの決断を迫られる。

$ mecab -u /root/work/pstack/ppp/utils/mtg_dic.dic

初手  名詞,一般,*,*,*,*,初手,ショテ,ショテ
フェッチランド   名詞,固有名詞,一般,*,*,*,フェッチランド,フェッチランド,フェッチランド
島 名詞,一般,*,*,*,*,島,シマ,シマ
渦巻く知識 名詞,固有名詞,一般,*,*,*,渦巻く知識,ウズマクチシキ,ウズマクチシキ
場合  名詞,副詞可能,*,*,*,*,場合,バアイ,バアイ
ワンマリ    名詞,固有名詞,一般,*,*,*,ワンマリ,ワンマリ,ワンマリ
憂き目   名詞,一般,*,*,*,*,憂き目,ウキメ,ウキメ
土地  名詞,一般,*,*,*,*,土地,トチ,トチ
セット   名詞,一般,*,*,*,*,セット,セット,セット
決断  名詞,サ変接続,*,*,*,*,決断,ケツダン,ケツダン

おわりに

このようにマニアックな語彙を覚えさせてから単語埋込表現を作るなどすることで、文書分類タスクの性能向上を図ることができます。
僭越ながら私がgensimでword2vecを訓練した実装はgithubで公開しています。

github.com

ユーザがじわじわと課金を積んでいく時系列データをクラスタリングしたい

R 統計 時系列解析

下図ではサービス利用開始の0日目から30日目までの、個々のユーザが課金を積み上げていく様子(のダミーデータ)を示しました。

f:id:sakatoken:20190605220208p:plain

この記事では
「このデータから『早熟ユーザ』『毎日じっくり課金するユーザ』『遅咲きユーザ』といった課金の積み上げ方のスタイルを抽出できないか？」という観点でクラスタリングを行っていきます。

前準備としてのスケーリング

とはいっても、上述のようなデータをいきなりクラスタリングしても「ものすごい課金してる人&その他大勢」のような結果が出てしまい、今回の関心の対象であるところの課金の積み上げ方のスタイルを表現することができなかったりします。

そこで「個々のユーザの最終的な合計課金金額」を最大値とするスケーリングを施し、全てのユーザのデータが0〜1の範囲におさまるように調節します。金額の多寡についての情報を捨てて、「いつ課金してくれたか？」を際立たせます。

すると、グラフとコードは以下のようになります。

f:id:sakatoken:20190605222245p:plain

序盤に課金したっきり課金がストップする人、レイトスターターな人などが見えるようになりました。

d <- d %>%
    arrange(user_id, time) %>%
    group_by(user_id) %>%
    mutate(sum_pay = sum(pay)) %>%
    mutate(
        cum_pay = cumsum(pay),
        cum_pay_scaled = if_else(sum_pay == 0, 0, cumsum(pay) / sum(pay))
    ) %>%
    ungroup()

d %>%
    ggplot(
        aes(time, cum_pay_scaled, color=as.factor(user_id))
    ) +
    geom_line() +
    guides(color=F)

どのような距離尺度でクラスタリングするか？

時系列データのクラスタリングというテーマについて、今回は以下の記事を参考にさせていただきました。

例えば時系列データ向けの距離尺度の一例としてDTWというものがあります。次に引用するのは、後者のNHNさんのブログからのものです。

音声認識タスクでは同じ発話でも話者によって(または同じ話者でも状況によって)そのスピードが違います。このようなズレを補完するためにDynamic Time Wrapping (DTW)のような距離関数が用いられます。

今回の解析では「早熟か？遅咲きか？」のような時系列上の隔たりを距離としてキチンと評価したいので、選択肢は色々とありますが、まずはシンプルにユークリッド距離を使うこととしました。また、クラスタリング手法にはkmeans法を使うこととしました。

クラスタリングしてみた結果

以下、クラスタリングした結果のグラフとコードです。
ここでは簡単のため、クラスタ数の検討をスキップして適当に k=6としました。

クラスタ2では全く課金しないクラスタ、クラスタ3では遅咲きなクラスタを表現できているように見えます。一方、クラスタ4と6は何度も課金している人と一回しか課金しない人の両者を含むように見えます。このデータに関してはクラスタリングが上手くいったとは言えません。

※今回はダミーデータの生成の過程で正規分布を取り入れているので「毎日課金してるクラスタの人数が多い」のような結果になっていますが、実際はこうなりません。おそらくほとんどのサービスでは。

f:id:sakatoken:20190605223731p:plain

# 時系列を横持ちにして余計なデータを削ぎ落とす
d.f <- d %>%
    select(user_id, time, cum_pay_scaled) %>%
    spread(time, cum_pay_scaled)

# クラスタリング
d.f.km <- kmeans(d.f %>% select(-user_id), centers = 6)
d.f$cluster <- d.f.km$cluster

# グラフ表示
d.f %>%
    gather("time", "cum_pay_scaled", `0`:`29`) %>%
    mutate(time = as.numeric(time)) %>%
    ggplot(aes(
        time, cum_pay_scaled, group=as.factor(user_id)
    )) +
    geom_line() +
    facet_wrap(~cluster)

ダミーデータ生成のコード

d <- tibble(
    user_id = rep(1:500, 30), # ユーザID 500名分 * 30日分
    q = rep(abs(rnorm(500, mean = 0, sd = 10)), 30) # 課金確率 500名分 * 30日分
) %>%
    mutate(
        time = trunc((row_number()-1) / 500), # 日付番号の付与
        q = q/(max(q) + 0.2) # 0< q <1 になるように調整
    ) %>%
    mutate(
        pay = map(q, function(q){ # 売上の確率的生成
            rbernoulli(1, p=q) * as.numeric(rmultinom(n=1, size=1, prob=c(1/3,1/3,1/3))[1:3,] %*% c(1000, 2000, 3000))
        }) %>% as.numeric()
    )

おわりに

太く短く課金してくれる人と、細く長く課金してくれる人と、どちらが利益への貢献度が高いか？ーーのようなテーマについて考えたとき、今回の分析を思いつきました。

アドホック案件をtidyにグラフ作ってpdfでレポート出すまでの流れ

R 統計

機械学習や統計モデリングを使った格好のよろしいデータ分析で成果を出したいものですが、「SQLで引いてきたデータからインサイトを得て、適切にグラフ作ってレポートにまとめる」のような作業も私は日々やっています。この手の分析はR言語のtidyverseを使うと効率的にさばけます。DWH（ここではBigQuery）からデータを引いてきて、集計結果をpdfに出力するまでの一連の流れをここに示します。

余談ですが、単なる可視化を行うときも、因果推論の見地から鳥瞰視点で「交絡ない？大丈夫？」を意識することで、より実り多き分析が出来ると思っています。

以下、もくじです。

R markdownじゃなくてJupyterを使う

私は今回紹介するような軽いアドホック分析においては、Jupyterの方が使いやすいと考えています。以下の動画では「分析の最後にpdfにアウトプットする様子」を示しました。

f:id:sakatoken:20190603203252g:plain

動画中で行なっている操作は以下の２つです。

  • 「Hide input all」ボタンを押してコードブロックを非表示にする
  • ブラウザの印刷機能でページを印刷する

最終的にpdfを出力するためのステップがこれだけで済むので、分析作業を行なっている最中は「ただただ必要なコードをコードブロックに書く」「マークダウンのブロックを挿入して説明文を書く」ことだけに専念できます。

なお、前者のコードブロックを非表示にするという操作は、nbextensionsというjupyter拡張機能群(非公式)の中の「Hide input all」で実現しています。ちなみに同拡張機能群の１つである「Table of Contents」を使うと、markdownで定義した段落に対して自動で段落番号を補完してくれます。（動画中でも適用しています）

DWHからデータを引いてくる

さて、pdf出力をキレイに行えることがわかったので、安心してnotebook上に必要なコードを書きなぐっていきます。まずはデータの取得です。
私の場合はbigqueryからデータを引いてくるのでこんな感じです。

library(bigrquery)
set_service_token("your_project_key")
project <- "project_name"

sql = 'SELECT ....'
res <- bq_project_query(x = project, query = sql)

tidyに前処理してggplotでグラフを書く

これについては偉大なる先人が編纂した書物（通称宇宙本）を中心に、様々に優れた解説があるのでここでは触れません。

RユーザのためのRStudio[実践]入門−tidyverseによるモダンな分析フローの世界−

RユーザのためのRStudio[実践]入門−tidyverseによるモダンな分析フローの世界−

markdownで説明文を書き加えてpdf出力

markdownブロックをnotebookのあちらこちらに挿入して全体の体裁を整えたら、上記の動画の手順でpdfを出力します。

参考にしたもの