警察庁オープンデータの前処理と死亡事故発生予測のための機械学習について

はじめに

先日、近所の幹線道路の信号付き横断歩道を横断しようとした際に、2メートルほど先で回避不能な速度で車が赤信号を無視して横断歩道を突っ切きるという非常に恐ろしいイベントに遭遇しました。
そこで、私は「このイベント、何回かあったら死亡事故が起きるぞ」と思いました。警察署に近隣道路の不満について電話をしてみたのですが、「警察もちゃんと考えて取り締まっていたり対策してますので」と言われ、何の成果も獲れない電話となりました。
そのようなモヤモヤの中、私は民間でできることは「客観的なデータでもって、エリア、曜日、天候、時間帯ごとの、死亡事故のリスクを見積もり、それを公開して注意喚起していくこと」なのではないかと思うに至りました。
よって、警察庁が公開している交通事故のオープンデータについて、その前処理と機械学習によるエリアの死亡事故の発生予測などについて記していきます。

前処理のためのデータやコードはこちらのGitHubにあります。
https://github.com/KamonohashiPerry/traffic_accident_analysis
スターを10個以上いただけて、データ作って良かったなと思いました。

目次

・データについて
・先行研究について
・前処理について
・EDAからわかること
・機械学習
・さいごに

データについて

2019年と2020年の交通事故に関するデータが警察庁よりオープンデータとして公開されています。
オープンデータ 2020年(令和2年)
オープンデータ 2019年(平成31年/令和元年)

どのようなデータかというと、
それぞれ2019年と2020年の間に発生した全国の事故について、発生した状況(時、場所、天候)やその結果に関するデータからなり、項目の数としては58個とオープンデータとしてはかなり豊富な情報量となっています。

データは本票、補充票、高速票の3つからなるのですが、今回は本票だけを用います。
高速道路の精緻な分析をしようと思うと、高速票を繋ぐ必要がありますが、私の目的は一般道にあるので、繋がないものとします。

先行研究について

事故の予測に関する先行研究

交通事故の分析に関する研究は以下のものがあり、参考にしました。
死亡事故を予測するという問題設定として、その特徴量としては場所、時間帯、曜日、場所かつ時間、場所かつ曜日かつ時間、天候、季節、風速・風向き、車線の数、歩道の有無、制限速度などを用いていることがわかりました。
そして、データをクラスタリングするためにDBSCANを用いたり、不均衡データであることから、アンダーサンプリングを行ったりしているようです。
最後に、アルゴリズムはSVMやロジスティック回帰、深層学習と色々と使われているようです。

・Live Prediction of Traffic Accident Risks Using Machine Learning and Google Maps

・Data-Driven Urban TrafficAccidentAnalysis and Prediction Using
Logit and Machine Learning-Based Pattern Recognition Models

・A Deep Learning Approach for Detecting Traffic
Accidents from Social Media Data

・A BASIC STUDY ON TRAFFIC ACCIDENT DATA ANALYSIS USING
SUPPORT VECTOR MACHINE

・Predicting Crash Injury Severity with Machine Learning Algorithm Synergized with Clustering Technique: A Promising Protocol

・Comparison Analysis of Tree Based and Ensembled Regression Algorithms for Traffic
Accident Severity Prediction

その他

・Taxi Demand Prediction System

・2019年04月号トピックス 警察庁のオープンデータを活用した独自AIによる犯罪発生予測

前処理について

ソースコードに関しては、GitHubに上げておりますので、
.pyなら
https://github.com/KamonohashiPerry/traffic_accident_analysis/blob/main/src/preprocess.py
https://github.com/KamonohashiPerry/traffic_accident_analysis/blob/main/src/preprocess_second.py

.ipynbなら
https://github.com/KamonohashiPerry/traffic_accident_analysis/blob/main/Traffic_Accident_Analysis.ipynb
を参照ください。

基本方針

・天候、時間帯、曜日、平日、道路の形状、スピード制限など、事故の結果(死亡の有無)を明らかに予測できそうなもの以外を特徴量にします。
・メッシュ単位での過去の事故に関する情報に関する特徴量の作成(学習を行う時点よりも過去の時点のデータのみ利用)

詳細

・事故の発生日時のdatetime化(列で年・月・日・時・分が分かれているため)
・マスタ情報(PDF)をもとにマスタのIDを名前に変換
・緯度経度のデータについて、度分秒から度への変換を行い、3次メッシュ(1km四方)のデータを作成
・メッシュ単位での特徴量として、過去の事故発生件数、死亡数などを集計
・メッシュ単位での特徴量として、過去に発生した事故に関するもろもろのデータの集計(そのメッシュで事故を起こした平均年齢や事故が発生してきた時間帯など)
・もろもろの変数のダミー変数化

EDAからわかること

EDAに関しては.ipynbのみで行っております。
https://github.com/KamonohashiPerry/traffic_accident_analysis/blob/main/Traffic_Accident_Analysis.ipynb

私がよく分析で使うのは、データの該当件数とカテゴリごとの率の可視化なのですが、以下の図の見方としては、青の棒グラフが各カテゴリの件数、赤の折れ線グラフが死亡率となっています。

都道府県


福井県、鳥取県、秋田県、島根県、岩手県などは死亡率が高そうです。

道路種別


農道(815件)や林道(88件)で起きる事故非常に死亡率が高くなるようです。

時間帯


夜間や明け方は死亡率が高くなるようです。

天候


霧の天候(294件)では死亡率が高くなるようです。

人口密集地かどうか


市街地ではないところでは死亡率が高くなるようです。

道路の形状


踏切-第一種(254件)はかなり死亡率が高くなるようです。カーブも死亡率が高くなります。そういえば実家の近くにも魔のカーブと呼ばれるところはありますし、事故も起きていました。

信号機のタイプ


点滅-3灯式(2566件)や点灯-押ボタン式(7703件)は死亡率が高いようです。

道路の幅


大きな交差点では死亡率が高い傾向にあります。単路に関しては一様な傾向ではないですが、小さい方が死亡率が高そうです。

道路の区切り


道路の区切りがペイントだけのケースは、中央分離帯がある場合よりも死亡率が高い傾向にあります。あと、中央線がポストコーン(チンアナゴみたいなやつ)だと死亡率が上がるようです。

歩道との区別


道路について歩道との区分がない場合に死亡率が高くなるようです。

速度制限


速度域が高くなるにつれ、死亡率が高くなるようです。あと、速度規制がされていないところも高くなっています。

曜日


土日の死亡率が高い。

祝日・その前日


祝日とその前日の死亡率は高い傾向にある。

以上のことから、得られた考察としては、
・土日や祝前日の道路はいつもより危ない
・夜間や明け方の道路は危ない
・制限速度の高い道路は危ない
・大きな交差点は危ない
・点滅信号は危ない
・カーブは危ない
・霧は危ない
でした。普段、出歩く際に意識しておくといいと思います。

機械学習

教師データを事故発生からの死亡有無(1or0)、特徴量を先行研究を参考にオープンデータから作れるものとして、
基本的には、LightGBMでの学習、ハイパーパラメータ探索、モデルの訓練(メッシュ番号でグループ化したクロスバリデーション)、テストデータ(学習は2020年11月まで、テストは2020年12月以降としている)での予測、日本地図への可視化を行っています。

ソースコードはGitHubに上げてあるので、ここでは記しません。

.pyなら(学習のみ)
https://github.com/KamonohashiPerry/traffic_accident_analysis/blob/main/src/train_model.py

.ipynbなら(学習だけでなく、予測と日本地図への可視化もしている)
https://github.com/KamonohashiPerry/traffic_accident_analysis/blob/main/Traffic_Accident_Analysis.ipynb

10-foldのCVをメッシュ番号でグループ化して学習させた際のAUCです。高い部類に入りますが、今回は不均衡データなので、抜群に当たっている感じではないです。


事故の死亡率は、0.9%なので、この予測の場合、2.2%が死亡すると当てれているということになります。また、予測から、全体の72%の死亡事故を当てているという結果です。
この性能に関しては今後、もっと高めていく余地はいくらでもあると思います。

こちらは、テストデータでの予測スコアの階級値(青い棒グラフ)と実際の死亡率(赤い線グラフ)を示したものです。スコアが高いものは信じてもいいかなぁと思いました。

特徴量重要度の上位10件は以下の通りです。

時間帯とか、防護策とか、信号機のタイプとか、市街地かどうかなどが予測に役立っているようです。

死亡事故発生を予測した確率を日本地図の関東エリアにプロットしたものが以下の図になります。

市街地のほうが死亡事故が起きやすいように見えますね。

東京都だけに絞ってみましたが、危険なエリアは都心に点在はしていますが、郊外のほうが死亡率が高そうです。

さいごに

交通事故にあいかけて、その怒りの感情を警察庁のデータにぶつけてみましたが、死亡事故が起きやすい条件などがある程度見えてきたと思いました。もちろん、なんとなく思っていた傾向もあったのですが、数字でハッキリと出せるのは良いですね。今後、旅行などをする際に、確率の高いエリアを警戒して行動するようにしたいです。
そして、地図への可視化は普段の仕事でほとんど行わないのですが、こうしてデータに触れて可視化の練習ができたのも有意義に思いました。
この警察庁のデータは分析するまでに工数がかかるデータとも言えるので、この記事をとっかかりとして、オープンデータの分析に挑戦して安全な生活を送る上での知見をシェアしていけると良いな、というか民間としてはこれくらいしかできないなと思っています。

参考情報

・[第3版]Python機械学習プログラミング 達人データサイエンティストによる理論と実践 (impress top gear)
・Pythonでプロットした地図を出力する
・Pythonで度と度分秒を相互に変換する

Bias Correction For Paid Search In Media Mix Modeling[A4一枚まで備忘録]

A4用紙1枚にまとめるイメージでメモを残そうという取り組みです。
今回は2018年のGoogleのMedia Mix Modeling(マーケティング・ミックス・モデリング)に関する論文です。Bias Correction For Paid Search In Media Mix Modeling前回と違い、因果推論の観点で考察などが書かれています。前回のものがモデリングのテクニカルな側面を描いたものとすると、今回はテクニカルな話は置いておいて、因果関係について様々なケーススタディを扱うというものとなります。

目次

Abstract
1.Introduction and problem description
2.Related work
3.Preliminary to Pearl’s causal theory
4.Methodology
5.Implementation
6.Case studies in simple scenarios
7.Case study with complex scenario
8.Discussion

Abstract

バイアスを回避するようなMedia Mix Modeling(以下、MMM)のモデルの推定について、バックドア基準に基づき、因果ダイアグラムを用いたアプローチをする。
ランダム化実験の結果との比較をすることで推定のバイアスについて考察を行なっている。
シンプルなモデルなら大丈夫だが、複雑なMMMのモデルに関しては不偏推定量は得られていないらしい。

1.Introduction and problem description

・MMMの歴史は長く、1964年のBordenの研究に始まるらしい。
・Fortune500にいるような色々な企業がMMMを試してきたが、課題が色々(データ収集、セレクションバイアス、広告の長期効果、季節性、ファネル効果など)ある。
・ターゲットとしているユーザー群の需要や関心が高まっているとき、広告支出と売上は共に伸びるので、セレクションバイアスが発生する。(実際のところ、儲かっているから広告を打つという場合がありうる。)
・この論文ではMMMにおける検索広告のセレクションバイアスについて扱う。Pearl流の因果推論のアプローチに従い、バックドア基準を満たすようなMMMの推定を行う。
・実験結果と照らし合わせて、不偏推定量がもとまることがわかった。(検索広告に関してだが。)

2.Related work

広告効果の検証などは以下の3つの研究がなされているとのこと。
・1.ユーザー単位のデータを利用して広告を触れたか触れていないかを比較する。傾向スコアなどを計算してマッチングや共変量の調整などを行う。
・2.キャンペーン単位で集計したデータを用いたアプローチで、KPIに関してキャンペーンを打たなかった場合の反事実との比較を行う。
・3.クエリ単位でオーガニック検索とペイド検索をランダムに出し分けてその増分を推定する方法。

この論文で扱うMMMは2つ目のアプローチに近い。違いとしては複数のメディアのキャンペーンを扱うこと。

なお、回帰不連続デザインを用いたケースなどの研究もあるそうな。

3.Preliminary to Pearl’s causal theory

以下のようなDAGを考える。

・Aはクエリに関するGoogle広告のオークション
・Qはユーザーの検索するクエリ
・Pはペイド広告でのインプレッション
・Oはオーガニック検索結果
・Yは売上

ようは、 Qが与えられたもとでOが決まり、QとAに影響を受けPが決まり、OとPによってYが決まると言う構図となります。
MMMにおいて広告効果を測るために因果ダイアグラムを理解することは大事と書かれています。
ここではPearl流の因果推論のフレームワークが紹介されています。
$$ X_i = f_i(pa_{i} , \epsilon_i) $$
\( X_i \)は子を表し、\( pa_{i} \)は親を表します。なお、\( f_i \)は関数、\( \epsilon_i \)は任意に決まるランダムな項で他の変数と独立するものとしています。

「因果効果の定義」
二つの変数XとYについて、XのYにおける因果効果は\( Pr(y | \check{x} ) \)で表現できる。Xがxとして与えられたもとで、\( X_i = f_i(pa_{i} , \epsilon_i) \)を解くことで\( Y=y \)の確率を計算できる。ここでの\( \check{x} \)はXをxに設定するという介入を意味している。

「識別性」
ダイアグラムで互換性のある、観測された変数の正の確率から、\( Pr(y | \check{x} ) \)を一意に計算することができるのであれば、XのYにおける因果効果は識別可能となる。

「d-分離」
因果ダイアグラムの二つのノードの間のパスについて、以下の二つの条件のいずれかが満たされる場合についてZノードによって分離されている、あるいはブロックされていると表現する。
・1.\( m \in Z \)として、そのパスが「\( i \to m \to j \)」あるいは、「\( i \gets m \to j \)」を含んでいる。
・2.\( m \notin Z \)として、mがZに依存しないものとして、そのパスが「\( i \to m \gets j \)」を含んでいる。

「バックドア基準」
ダイアグラムにおいて以下の二つを満たす場合、変数Zはバックドア基準を満たす。
・1.変数ZのノードがXの子孫でない。(Xの共変量を除外することに関する基準)
・2.変数ZがXとY(Xへの矢印を含んだもの)の間の全てのパスをブロックしている。(交絡因子の適切な集合をZが含んでいるための基準)

もし、変数Zのバックドア基準が満たされるのであれば、XのYにおける因果効果は以下のように表現することができる。つまりZでXのYにおける因果効果を計算できる。
$$ Pr(Y | \check{x}) = \sum_{z} Pr(Y | x, z) Pr(z) $$

図を再掲しますが、ここでペイド広告の売上への因果効果を知りたいとします。

その場合、QはPの子孫のYを持たないため、バックドア基準の1つ目を満たし、「\( P \gets Q \to O \to Y \)」であることからd分離の条件の1つ目を満たすことからバックドア基準の2つ目を満たす。
つまり、バックドア基準を満たすことから、検索クエリを使うことでペイド広告の売上への因果効果を測れることになります。

Pearl流の因果推論のフレームワークの3つの課題
・1.因果ダイアグラムをどう構築するのか?
・2.全ての必要な変数を正確に観測、計測することができるのか?
・3.識別が可能だとして、有限であるサンプルサイズにおいて、\( Pr(Y | X, Z) \)の関数系はどうなるか?

広告まわりのデータは潤沢にデータがあるわけではないし計測できないものもあるため、どの課題も無視できない。

4.Methodology

「シンプルなシナリオ」
売上への他のメディアの貢献は無視できるものとして、広告のみをチャネルとして想定する。

$$Y_{t}=\beta_0 + \beta_1 X_t + \epsilon_t$$

\( X_t \)はあるプロダクトの検索広告の支出、\( Y_t \)はあるプロダクトの売上。\( \beta_1 \)は広告の売り上げにもたらす効果でROASと呼ばれる、\( \epsilon_t \)は\( X_t \)では説明づけることができなかった売上への影響を指す。計量経済学などで問題として扱われる内生性はこの\( \epsilon_t \)と\( X_t \)の相関によって起きる。

実際に、\( \epsilon \)がXと相関するとして式を変形すれば、不偏推定ができないことがわかります。

ここで、潜在的にプロダクトの売上に影響を与える検索クエリを集約した\(V\)という統計量を考えます。
誤差項\( \epsilon \)が以下のように、オーガニックサーチによる影響(\( \epsilon_1 \))と経済要因(\( \epsilon_0 \))によるものからなるとします。

$$ \epsilon = \epsilon_0 + \epsilon_1 $$

ここで、
$$ \epsilon_1 \perp X | V $$

$$ \epsilon_0 \perp X | V $$
を満たすようなVを見つけることができれば、因果効果を測れることになります。


※想定1:検索広告の予算はないとする。
※想定2:検索クエリで条件付けした際に、広告主の入札や競合の行動の影響は無視できる。

理論1
ペイド広告ついて上の図を想定した際に、XとY(sales)が完全に相関していないならば、正則化条件のもとで、検索広告のROAS(\( \beta_1 \))を推定することができる。

$$Y = \beta_0 + \beta_1 X + f(V) + \eta$$

ここで\( f(\cdot) \)は未知の関数、\( \eta\)はXや\(f(V)\)と相関しない残差とします。

上の図において、変数ZのノードがXの子孫でないというバックドア基準の一つ目はクリアしていることは明らかで、d分離の観点に関しても一つ目を満たしているためブロックしているため、バックドア基準の二つ目もクリアしている。

$$ E(Y | X, V) = \beta_0 + \beta_1 X + E(\epsilon | X, V)$$

Vで条件付けるとXは誤差項と直交することから、
$$E(\epsilon|X,V)=E(\epsilon|V)$$
となり、
$$ E(Y | X, V) = \beta_0 + \beta_1 X + f(V) $$
と表すことができる。
この式で回帰分析をすることでROASである\(\beta_1\)を推定することができる。

「複雑なシナリオ」
ここでは検索広告以外(\( X_2 \))も売上に影響を与えるというケースを扱っています。

まず、予算に関して検索広告と検索広告以外の広告(\( X_2 \))の両方に制約がある場合を考えます。以下の図が、今回のダイアグラムですが、先程の図に予算が追加されていることがわかります。

なお、ここでは広告のヒストリカルな影響などの複雑さは考慮していません。

一方で、予算が検索広告に関しては制約となっていないケースだと以下のダイアグラムになります。

理論2
1.広告予算が検索広告に制約を設けているケースで、広告のラグ効果が無視できるとして、もし\(X_1\)が\(V\)や\(X_2\)と完全に相関しておらず、正則化条件を満たすなら検索広告のROAS(\( \beta_1 \))は以下の回帰モデルで推定することができる。

$$ Y = \beta_0 + \beta_1 X_1 + f(V, X_2) + \eta$$
ここで

$$f(v, x_2) = E(\epsilon_0 | V = v, X_2 = x_2) + E(\epsilon_1 | V = v) + E(\epsilon_2 | X_2 = x_2)$$

となり、\( \eta\)は残差で\(X_1\)と\(f(V,X_2)\)と相関しない。

ダイアグラムを再掲します。

2.広告予算が検索広告に制約を設けていないケースで、正則化条件を満たすなら検索広告のROAS(\( \beta_1 \))は以下の回帰モデルで推定することができる。

$$Y = \beta_0 + \beta_1 X_1 + f(V) + \eta$$
ここでfは未知の関数とする。

ダイアグラムを再掲します。

以上、因果ダイアグラムを描き、バックドア基準を満たしバイアスのないパラメータ推定をするということが書かれていました。

ここではシンプルに特定の広告の効果を見ていましたが、MMMにおいては複数のメディアの効果を求めなければならないです。

考えられるアプローチとしては、最初にバイアス補正を行い\( X_1 \)の影響を推定し、それを使って\( X_2 \)を推定するというもの。

5.Implementation

・検索クエリのデータの集計の仕方
・STEP1:広告主と競合のサイトを識別する。
・STEP2:特定のエリアで自然検索でURLが見つかる回数をカウントする。
・STEP3:クエリを3つのグループに分類。おのおの閾値を設ける。
・STEP4:クエリをグループごとに分けたものを推定の際に誤差項と直交するものとして扱う。(閾値は50%がいいらしい。)
・TBA

6.Case studies in simple scenarios

・ナイーブな推定手法では過大に効果を推定していることがわかった。
・需要による調整だけではバイアスを十分にとりきれないことがわかった。
・SBC(search bias correction)の手法によりROASがランダム化実験の結果に近くなった。
・TBA

7.Case study with complex scenario

・TBA
・TBA
・TBA

8.Discussion

・TBA
・TBA
・TBA

DAG描くのに使ったコード

DAG_draw.ipynb

[PRの巻]RユーザのためのRStudio[実践]入門(改訂2版)

今回はご献本PRの巻です。お題はこちらの本、『改訂2版 RユーザのためのRStudio[実践]入門〜tidyverseによるモダンな分析フローの世界』です。

リアルでもバーチャルでも大変お世話になっているy__mattuさんより、ありがたく献本いただけました。どうやら、今回の改訂では56ページほどページ数が増えているようです。

目次やら詳細は技術評論者のページをご覧ください。

端的に言うと、データ分析でR言語を使ったことがあるユーザーで、Rの比較的新しいライブラリの使い方を知りたい方、RStudioを使って生産性高く分析をしたい方にピッタリの本です。データ収集・前処理・可視化・レポーティングを学ぶことができます。

元よりこの本のユーザーであったことから、これまで色々なマーケター・エンジニアにこの本を勧めてきました。

社内にデータ分析に関心のある人が一定数いて、そのような人にR言語を勧める際にこの本はちょうど良いと思っています。
この本は、統計学や機械学習などについて触れてはいないですが、データ分析のほとんどの工程は前処理や可視化ですし、それを行えずして統計学も機械学習も成し遂げることはできないので、現実を知ってもらう上でも良い本です。
この本で自由自在に前処理・EDAができるようになったら、統計学や機械学習について書かれたR言語の本をやるように伝えています。

ここでは、初めてこの本を手に取る方、1版をすでにお持ちの方のそれぞれのお気持ちになって見どころを述べていきます。

  • 初めてこの本を手に取る方
    • 見どころ1
      RStudioの画面の説明やショートカットの説明などが丁寧で、最初の段階で生産性の高いコーディングのための準備が捗る。
    • 見どころ2
      Webスクレイピングにしっかりと章が割かれている。
      データ分析を学ぶに際して、多くのケースで手元にデータがなかったり、関心のあるデータが表形式で提供されていなかったりします。そのため、Webスクレイピングをしっかりと学べるのはデータ分析を自分でやっていくために不可欠な技術だろうと思います。
    • 見どころ3
    • モダンなR言語のツール、dplyr・tidyr・ggplot2・rmarkdownなどについて学ぶことができます。前処理して可視化してレポートするという一連の操作をこの本で学ぶことができます。

  • 1版をすでにお持ちの方
    • 見どころ1
      stringrやlubridateなどのライブラリの使い方についてAppendixで40ページ以上が割かれています。RやPythonを行ったり来たりする生活をしていると忘れがちなので、充実のAppendixは嬉しいと思います。stringrは何が出来たっけ?てなりますし、lubridateも非常にしばしば忘れる気がします。
    • 見どころ2
      tidyrの記述が増えていました。Pythonに浮気をしていると、Rのモダンな書き方に付いていけなくなるかもしれないので、こういった本でキャッチアップできるのは良いことに思います。

1版を会社のデスクに置いたまんまにしているので、今度しっかり比較して追記しようと思います。

Bayesian Methods for Media Mix Modeling with Carryover and Shape Effects[A4一枚まで備忘録]

A4用紙1枚にまとめるイメージでメモを残そうという取り組みです。

今回はメディアミックスモデリングに関して仕事があった際に参照できるメモの一つとして残します。
先日、Python/STAN Implementation of Multiplicative Marketing Mix Modelに参考文献としてBayesian Methods for Media Mix Modeling with Carryover and Shape Effectsの存在を知りました。時系列の重回帰の一種じゃないかと思っていたんですが、天下のGoogle様がどのようにマーケティングの成果を測ろうとしているのか、この論文から知ってみようと思います。

ちなみに、Python/STAN Implementation of Multiplicative Marketing Mix Modelの写経したものはMarketingMixModelingに残しています。

このコードでマーケティングにおけるチャネルごとのROAS(Return On Advertising Spend)を推定することができます。

以下の図はコードを回すことで得られるものですが、チャネルごとのROASの事後分布となります。

この図もコードを回すことで得られるものですが、広告などの支出に関する売上の関数をフィットしたものになります。精度高くこれらの関数を求めれば、過剰に広告に投資している可能性のあるチャネルを見つけることが可能になります。

目次

Abstract
Introduction
Model Specification
Estimating the Bayesian Model
Attribution Metrics and Optimal Media Mix
Application to a Simulated Data Set
Impact of Sample Size on Estimation Accuracy
Choice of Priors
Application to Real Data and Model Selection
Conclusion

Abstract

ROASやmROASに関するアトリビューション分析をベイズモデリングで行い、事後分布から効果を推定している。
広告関連のデータが、売上などに対して、数ヶ月後に遅れて効果が現れるキャリーオーバー効果の想定、支出を増やしても効果がなくなる収獲低減の想定などをモデルに置いている。
モデルの識別に関してはBICを評価基準として用いている。

Introduction

  • Media mix models(MMM)はメディア支出が売上にどのように影響を与えるのかを理解するために、最適なメディア投資を行うための支出の配分を決めるために使われる。
  • 扱うデータとしては売上、物価、プロダクトの分布、様々なメディアの支出、マクロ経済や天候、季節性、競合の情報などの外部要因が含まれる。
  • RCTが難しいため、あるいは傾向スコアなどのモデリングもデータが少なく難しいため、回帰分析が行われることが多い。
  • 線形回帰だと、広告がサチっている効果(広告が飽和している状態)や収獲低減の効果も表現できない。
  • 広告にはキャリーオーバー効果と言って、売上への効果が遅れて現れるという可能性が広く信じられている。
  • 事前分布に広告に関する様々なノウハウを反映させることで、データ数に比して推定するパラメータの数が多い分析に向きあう。

Model Specification

  • モデル構築には週次データが使われるケースが多い。
  • 売上のデータ、メディア支出のデータ、コントロール変数として売上と関係してそうなデータが扱われる。コントロール変数は業界によって異なる。
  • キャリーオーバー効果
    $$ adstock(x_{t-L+1, m}, \dots, x_{t, m} ; w_{m} ,L) = \frac{\sum_{t=0}^{L-1} w_m(l) x_{t-l, m} }{\sum_{l=0}^{L-1} w_m(l)} $$

    \( w_m(l) \)は非負のウェイトでラグ期間の関数になっている。
    \( L \)はキャリーオーバー効果の最大期間。
    \( x_{t, m} \)はある期間のあるメディアへの支出。
    論文ではL=13とされている。
    ウェイトの関数について、メディアの売上に対する効果は徐々に減っていくが、効果が遅れて生じることもあるので、効果のピークをずらした表現も扱われている。
    ウェイトの関数で最初にピークがくる場合、

    $$ w_{m}^{g} (l; \alpha_{m}) = \alpha_{m}^{l}, \\ l = 0, \dots, L-1 , \\ 0 < \alpha_m < 1 $$

    ウェイトの関数で遅れてピークがくる場合、

    $$ w_{m}^{d} (l; \alpha_{m}, \theta_m) = \alpha_{m}^{(l – \theta)^2}, \\ l = 0, \dots, L-1 , \\ 0 < \alpha_m < 1 ,\\ 0 \leq \theta_m \leq L -1 $$

    となる。
    \( \alpha_m \)はメディアの効果の、ある期から次の期への維持率。
    \( \theta_m \)は遅れてピークが現れる程度に関するパラメータ。

  • 形状効果
    $$Hill(x_{t,m};\mathcal{K}_m, \mathcal{S}_m) = \frac{1}{1+ \left (
    \frac{x_{t, m}}{\mathcal{K}_m} \right )^{- \mathcal{S}_m} }, \\ x_{t, m} \geq 0 $$

    \( x_{t, m} \)はある期間のあるメディアへの支出。
    \( \mathcal{S}_m \)は形状に関するパラメータで傾きに関わるもの。
    \( \mathcal{K}_m \)は形状に関するパラメータでサチるポイントに関わるもの。
    分子に\( \mathcal{K}_m ^{\mathcal{S}_m} \)を足して引くことで、

    $$ \beta_{m} Hill_{m} (x_{t, m}) = \beta_{m} – \frac{\mathcal{K}_m^{\mathcal{S}_m}\beta_{m}}{x_{t,m}^{\mathcal{S}_m}+\mathcal{K}_m^{\mathcal{S}_m}} $$

    に変換できる。
    なお、\( \beta_m \)は各メディアの回帰係数となる。これでもって各メディア支出の売上に対する効果を見れる。パラメータの大きさによって様々な形状をとりうる。

  • キャリーオーバー効果と形状効果をモデルで一緒に表現する。
    シンプルにするためにメディア間のシナジー効果を無視している。

    $$ y_t = \tau + \sum_{m=1}^{M} \beta_m Hill(x_{t, m }^{*}; \mathcal{K}_m,\mathcal{S}_m) + \sum_{c=1}^{C} \gamma_{c}z_{t, c} + \epsilon_t $$

    \( y_t \)は売上。
    \( x_{t, m }^{*} = adstock(x_{t-L+1, m}, \dots, x_{t, m} ; w_{m}, L) \)は前述の広告支出のキャリーオーバー効果を表現したもの。
    \( \tau \)は広告支出によらないベースラインの売上。
    \( z_{t,c} \)はコントロール変数。
    \( \gamma_c \)はコントロール変数の係数。
    \( \epsilon_t \)はホワイトノイズ。

Estimating the Bayesian Model

  • 階層ベイズに拡張することもできる。
  • ハイパーパラメータはメディアの効果が非負であるとか、割合は0から1だとか色々と事前情報として設定する。

Attribution Metrics and Optimal Media Mix

  • 経路ごとのROSAやmROASを計算してメディアの最適化をしたい。
  • メディアに投資することでの予測売上と、投資なかった場合の売上の差分でもって成果を測る。
  • ROASの定義式
    $$ ROAS_{m} = \frac{\sum_{t_{0} \leq t \leq t_{1} + L – 1 } \hat Y_{t}^{m} (x_{t-L+1, m}, \dots, x_{t, m}; \Phi) – \hat Y_{t}^{m} (\tilde{x}_{t-L+1, m}, \dots, \tilde{x}_{t, m}; \Phi) }{\sum_{t_{0} \leq t \leq t_{1}} x_{t, m} } $$

    \( \hat{Y} \)は売上の予測。
    \( \tilde{x} \)はメディアの支出の変化。

  • mROASの定義式
    $$ mROAS_{m} = \frac{\sum_{t_{0} \leq t \leq t_{1} + L – 1 } \hat Y_{t}^{m} (\tilde{\tilde{x}}_{t-L+1, m}, \dots, \tilde{\tilde{x}}_{t, m}; \Phi) – \hat Y_{t}^{m} (x_{t-L+1, m}, \dots, x_{t, m}; \Phi) }{ 0.01 \times \sum_{t_{0} \leq t \leq t_{1}} x_{t, m} } $$

    \( \hat{Y} \)は売上の予測。
    \( \tilde{\tilde{x}} \)はメディアの支出の1%の変化。

  • ROASの分布は事後分布と上述の定義式から得られる。
  • メディアミックスの最適化
    予算という制約付きの極値問題を解く。ラグランジュ未定乗数法を用いる。

Application to a Simulated Data Set

  • シミュレーションしたデータを使っている。2年間の週次のデータで、売上、3つのメディア、一つのコントロール変数からなる。

Impact of Sample Size on Estimation Accuracy

  • 単一のデータセットのみで評価することが難しいので、それぞれが独立した500個のデータセットを生成し、それぞれ推定をしている。

Choice of Priors

  • サンプルサイズが小さい場合、事後分布が事前分布の影響を大きく受けるため、事前分布に関して色々と考察している。
  • \( \beta \)に関する事前分布。設定次第で結構変わる。メディアの支出にかかってくるので分析に与える影響は大きい。
  • \( \mathcal{K}_m \)に関する事前分布。要はサチるポイントに関するもの。設定次第でグラフの挙動が結構変わる。

Application to Real Data and Model Selection

  • シャンプーの広告主の実際のデータを使って分析。
  • 2.5年の週次データで、反応変数として売上、メディアデータとしてTV・雑誌・ディスプレイ・YouTube・検索、コントロール変数として小売のデータ(平均オンスあたり価格)・重み付けされた全商品の流通量や広告量が扱われている。
  • BICを頼りにモデルを選択した。BICは事後分布のデータから計算できる。

Conclusion

  • メディアミックスモデルはモデルの識別によって結果が大きく変わる。
  • ベイズモデルでキャリーオーバー効果や形状効果を考慮したモデルを求め、それらの事後分布からROASやmROASを計算した。
  • データ数が少ないとバイアスの伴った推定となる。事前情報の選定にも慎重になるべき。
  • BICでモデルを識別した結果、キャリーオーバー効果や形状効果を想定した設定が選ばれた。
  • 残差に自己相関があったため、今回のモデルはまだ不完全で改善の余地がある。

Recommender Systems: The Textbookの要点まとめ(随時更新)

はじめに

レコメンド関連の書籍を探していた際に、
Recommender Systems: The Textbook (English Edition)の無料公開されているPDFを見つけたので、それについて読んでは追記するスタイルでメモを残していこうと思います。(すごく長くなる予感)

これまで読んできたレコメンド関連の本の中では、説明が丁寧だったり事例が豊富に思います。数式はあまり出てこないですが、言葉でちゃんと説明しようとしているのが感じられます。『AIアルゴリズムマーケティング 自動化のための機械学習/経済モデル、ベストプラクティス、アーキテクチャ (impress top gear)』のレコメンドの章もわかりやすく幅広いトピックが扱われていますが、それに匹敵する本にも思います。

第1章(An Introduction to Recommender Systems)

Goals of Recommender Systems

  • 1.ユーザーとアイテムの組み合わせの評価値の予測問題
    過去のユーザーのプリファレンスデータを使って学習する。未知のアイテムに関しては行列補完問題とも見なせる。
  • 2.上位k個の推薦問題
    あるアイテムに対して、あるユーザーの上位k個のアイテムを推薦する。

推薦システムの第一の目的はプロダクトの売上向上。

  • 推薦システムの運営上ないし技術的な目的
    1.関連性(Relevance)
    2.新規性(Novelty)
    3.意外性(Serendipity)
    4.推薦の多様性(Diversity)
  • 推薦システムの例
    1.GroupLensの推薦
    2.Amazon.comの推薦
    3.Netflixの映画推薦
    4.Googleニュース推薦
    5.Facebookの友達推薦

Basic Models of Recommender Systems

  • 協調フィルタリング(Collaborative Filtering)
    1.メモリーベース
    ・ユーザーベース協調フィルタリング
    ・アイテムベース協調フィルタリング
    2.モデルベース
  • コンテンツベース(Content-Based)
  • ナレッジベース(Knowledge-Based)
    滅多に買われないようなアイテム(自動車、不動産)の際に考える。評価は使わない。その代わり、ユーザーの要望とアイテムの近さや、ユーザーの制約などを用いる。
    1.制約ベースの推薦システム
    2.条件ベースの推薦システム
  • デモグラフィック(Demographic)
  • ハイブリッド、アンサンブルベース(Hybrid and Ensemble-Based)
  • 評価の種類
  • 明示的な評価と暗黙的な評価
  • 欠損値分析との関係性

Domain-Specific Challenges in Recommender System

  • コンテキストベース(Context-Based)
    時間や場所、ソーシャルデータを使う。季節性を考慮したりする。
  • 時間依存(Time-Sensitive)
    1.時間と共に評価は変わる。流行りすたりがある。
    2.特定の時期だけ、特定のアイテムの評価があがることはある。(レインコート、ダウン)
  • 位置情報(Location-Based)
  • ソーシャル(Social Recommender System)
    ネットワーク構造やソーシャルでの繋がり、タグ、組み合わせに基づいている。

Advanced Topics and Applications

  • コールドスタート問題(Cold-Start Problem)
  • 攻撃耐性(Attack-Resistant)
    アルゴリズムをハックされると情報推薦の効率性が阻害される。本来、上がってくるべきものが上がらないので。
  • 集団への影響(Group)
    推薦アルゴリズムは個々人よりも集団のプリファレンスに影響を与える。
  • 複数の基準(Multi-Criteria)
    ユーザーは単一の基準でアイテムを評価はしない。
  • アクティブラーニング(Active Learning)
    ユーザーが評価をしやすいようにするための取り組み、たくさん過去に評価したユーザーがいた場合は、新しい映画に関して評価をしやすいように予測するとか色々。
  • プライバシー(Privacy)
  • 様々な適用領域
    小売、音楽、コンテンツ、web検索、クエリ、デジタル広告

第2章

TBD

第3章

TBD

第4章

TBD

第5章(Knowledge-Based Recommender Systems)

  • コンテンツベースは新しいアイテムの推薦には良いが、新しいユーザーには推薦できない。協調フィルタリングやコンテンツベースは過去の履歴といったデータをベースに扱うが、ナレッジベースはユーザーが何を求めているのかに基づいている。
  • 不動産や自動車、旅行商品などは滅多に購入しないため、評価のデータは十分に得ることができない。
  • 映画はレコメンドされた際に、その理由が明示されなくとも受け入れるが、不動産や自動車は詳細な情報抜きには、そのレコメンドを受け入れることは難しい。
  • ナレッジベースのレコメンドが有効なケース
    1.ユーザーの要望が明確なケース。協調フィルタリングやコンテンツベースでは満たせない。
    2.選択肢が幅広く商品の複雑さが大きい商品のケース。詳細な評価を得にくい。
    3.商品の評価が時間と共に大きく変わるケース。パソコンとか車はすぐ新しいものが出てきたりして評価が変わってしまう。
  • ナレッジベースのレコメンドにおけるユーザーとのコミュニケーション
    1.対話システム
    対話でのフィードバックループのコンテキストからユーザーのプリファレンスを把握する。
    2.検索ベースシステム
    質問に答えることでプリファレンスを把握していく。
    3.ナビゲーションベース
    ユーザーがアイテムに対して反応をしてしていき、その結果に応じて望んでいるアイテムに到達できる。

Constraint-Based Recommender Systems

  • アイテムの属性に関する強い要望や制約を識別する。
  • 制約ベースのレコメンドシステムにおける3つのインプット
    1.ユーザー固有の情報(デモグラ、リスク嗜好)、必須要件(バストイレ別とか)
    2.ナレッジベースの属性情報。アイテムに応じた属性を顧客の属性にマッピングする。
    直接:都心か郊外か→郊外→近隣のものを推薦
    間接:家族の人数→5人以上→寝室が2つ以上のものを推薦
    3.アイテムに関するデータを全部使う
  • 関連した結果を返す際の3つのアプローチ
    1.ユーザーの回答をもとに条件を類推し、ユーザーの要望に加える。(例えば、家族の人数が5以上と回答された場合に、寝室の数は3以上必要だろうとして要望に加える。)
    2.要望をもとにデータベースのクエリを生成
    3.ユーザーの要望に関連したカタログを抽出する
  • 対話のアプローチ
    1.インターフェースを通してユーザーのプリファレンスの初期値を得る。
    2.マッチさせたアイテムを並べてユーザーに示す。なぜそれを表示したのかの理由も付ける。
    3.条件緩和や追加の要望を受け付ける。レコメンドの件数が多すぎたら要件を増やして、少なすぎたら要件を減らす。
  • マッチしたアイテムの並べ方
    効用関数を使ったりする。
    $$ U( \overline V ) = \sum_{j=1}^{d}w_j \cdot f_j(v_j) $$
    \( w_j \)は評価のウェイト、\( v_j \)はマッチした属性の値。
  • 受け入れられない結果や空集合への対応
  • 制約の追加

Case-Based Recommenders

  • 制約ベースのレコメンドとは違い、強い制約がない。良いレコメンド結果になるまでユーザーの要望を繰り返し修正していくアプローチとなる。
  • ケースベースのレコメンドの場合の重要な要素
    1.類似度の指標

    $$f(\overline{T}, \overline{X}) = \frac{\sum_{i \in S}w_i \cdot Sim(t_i, x_i)}{\sum_{i \in S} w_i} $$

    Sは属性の集合で、属性ごとに重みが違う。tはターゲット。xはあるアイテム。

    ・対称
    $$ Sim(t_i, x_i) = 1 – \frac{|t_i – x_i|}{max_i – min_i} $$
    ・非対称

    $$ Sim(t_i, x_i) = 1 – \frac{|t_i – x_i|}{max_i – min_i} + \alpha_i \cdot I( x_i > t_i ) \cdot \frac{|t_i – x_i|}{max_i – min_i} $$

    ・多様性
    1から類似度を差し引いたものが多様性の尺度となる。
    $$D(\overline{X}, \overline{Y}) = 1 – f(\overline{X}, \overline{Y})$$
    2.批評
    ・推薦結果がユーザーに示されてから、批評を通じてフィードバックを受ける。
    ・単純な批評:特徴の一つを変更してもらう
    ・複合的な批評:複数の特徴量を変更してもらう
    ・動的な批評:データマイニングからフィードバックを掴み取る(アソシエーション分析とか)

Persistent Personalization in Knowledge-Based Systems

  • ナレッジベースのレコメンドでの過去ログの蓄積があれば、効用関数や類似度関数、制約の提案、動的な批評に関してパーソナライズすることが可能。

第6章

TBD

第7章

評価の目標に正確度がよく使われるが、ユーザー体験において重要な指標は目新しさや信頼度やカバー範囲やセレンディピティなど、それ以外の指標によるところが多い。正確度は大事な指標ではあるが、それだけでは良いレコメンドシステムを作ることは難しい。

Evaluation Paradigms

レコメンドシステムの評価方法にはユーザースタディ、オンライン評価、オフライン評価の3つがある。

  • ユーザースタディ
    特定のタスクを行い、レコメンドシステムを触ってもらい質問をするなどして、ユーザーからのフィードバックをもらう。ただし大規模に行うとコストもかかるし、参加する時点でユーザーの性質が違う可能性もあり、バイアスが含まれる可能性がある。
  • オンライン評価
    レコメンドシステムを実装したサイトなどで、ユーザーのグループをAとBに分け、アルゴリズムAとアルゴリズムBをそれぞれ用いて、その結果としてコンバージョンレートなどを比較する。バンディットアルゴリズムを使うのも良い。ただ、オンライン評価は実装含めコストがかかる。
  • オフライン評価
    レコメンドアルゴリズムにおいて最もポピュラーな評価方法で、過去のデータにおける評価を行う。当然、オンライン評価のほうが望ましいが、統計的にロバストであろうと観点で非常によく使われる。ただし、セレンディピティや目新しさを評価することは難しい。

General Goals of Evaluation Design

  • Accuracy
    モデルの訓練にも、評価にも使われる。レコメンドシステムにおいても、機械学習の分類や回帰タスクと同様に、ホールドアウトや交差検証などを行って評価する。
    あるユーザーのあるアイテムに対する評価の誤差を以下の式で表す。
    \(e_{uj} = \hat r_{uj} – r_{}uj \)
    なお、Eは評価対象の集合とする。その際、正確度の評価は以下のようになる。

    $$ MSE = \frac{\sum_{(u, j) \in E} e_{uj}^2 }{|E|} $$

    $$ RMSE = \sqrt{ \frac{\sum_{(u, j) \in E} e_{uj}^2 }{|E|} }$$

    なお、他にもランキングベースの評価方法もある。
    ただ、精度が高くても、ユーザーが絶対に買うであろうアイテムしかレコメンドできない場合は、レコメンドシステムの意味合いとしては薄れる。

  • Coverage
    レコメンドシステムが精度が高すぎる場合、一部のアイテムは一部のユーザーに全く推薦されなくなる。
    k点以上の評価を付けるであろうと予測したユーザーの割合は、ユーザー空間でのカバレッジと呼ばれる。
    一般的に、近傍ベースの手法を考えた際に、近傍を広げた場合、精度は落ちてしまうが、カバレッジは高くなる。ようは、精度とカバレッジはトレードオフの関係にある。
  • Confidence and Trust
    予測した評価の範囲として信頼区間が与えられる。予測した値の誤差が信頼区間にマッチしているかどうかで測る。同じアルゴリズムにおいて、信頼区間が狭いものほど良いと判断する。ただし、予測が正確だとしても、評価をユーザーに信じてもらえないと役に立たない。推薦理由がロジカルかどうかが重要となる。
  • Novelty
    ユーザーがこれまで気づかなかったような推薦ができているかどうかを表す指標。オンライン評価として画面上でユーザーに回答してもらうという方法で測るケースが多いが、いつでもその評価を実行することはできない。
    一方、オフラインであれば、それが可能となる。現時点よりも将来に選びそうなアイテムを推薦できるようになれば、それがノベルティに繋がる。
    過去に評価されているアイテムを、正しく推薦されてしまうとノベルティに関するペナルティが与えられるとする。将来時点で評価されるアイテムを正しく推薦できたらノベルティに関してリワードが与えられると考える。ようは将来と過去の予測の正確度の差分でもってノベルティを測ることとなる。
  • Serendipity
    セレンディピティは幸運な発見を意味する。レコメンドの成功における驚きのレベルを測る。セレンディピティはノベルティよりも強い条件で、セレンディピティがある場合、そのアイテムは目新しいが、逆は真ではない。
    オンライン評価の場合、ユーザーからのフィードバックで意外性があったかどうかを測ればいい。
    オフライン評価の場合、コンテンツベースで選んだものではないアイテムがユーザーに選ばれた場合に、セレンディピティがあると見なすことができる。コンテンツベースのレコメンドとそうでないレコメンドの結果の二つを用意する必要がある。
  • Diversity
    推薦リストの中でのアイテムの多様性を測る。多様性が高まれば、ユーザーのニーズを捉える可能性が高まると考えられる。
    多様性はコンテンツを元にしたアイテム間の類似度で測る。全てのアイテムのペアに対しての平均的な類似度が多様性を表す。
  • Robustness and Stability
    偽の評価や時間とともに急激にデータが変化した際に、レコメンドシステムが影響をあまり受けない場合、そのレコメンドシステムはロバストで安定していると見なすことができる。レコメンドシステムの利用者の中には、ランキングを不正に高めることが自身の利益に繋がることがある。
  • Scalability
    近年、暗黙フィードバックも考慮すれば、レコメンドシステムを構築するためのデータは比較的簡単に入手できるようになった。そこで、データの大きさがもたらす問題に関して注意を払う必要がある。
    1.機械学習モデルの訓練時間
    2.予測(推論)時間
    3.計算の際に要求するメモリの量

Design Issues in Offline Recommender Evaluation

  • レコメンドシステムを開発するに際して、データは3つのパーツに分割される。
    1.訓練データ(モデルの構築)
    2.検証データ(モデル選択・パラメータチューニング)
    3.テストデータ(二度漬け禁止)
  • ※豆知識:Netflix Prizeのデータは学習データとテストデータがかなり違っていたらしい。学習データに古いレーティングが入っていて、テストデータに新しいレーティングがあるという状態だったとのこと。
  • Segmenting the Rating for Training and Testing
    ・Hold-Out
    ・Cross-Validation
  • Comparison with Classification Design
    回帰や分類タスクと同様にサンプルセレクションバイアスの問題がはらんでいる。

Accuracy Metrics in Offline Evaluation

  • レコメンドシステムの評価方法としてはランキング系のものが良いが、わかりやすさを重視するならばRMSEが良い。
  • Measuring the Accuracy of Ratings Prediction
    ・MSE
    ・RMSE
    ・MAE
    ・NRMSE
    ・NMAE
  • RMSE versus MAE
    RMSEは異常の影響を受けやすい。一方でMAEは影響を受けにくい。
  • Impacr of the Long Tail
    主要なアイテムによって評価は影響を受けるが、売り上げの大半はロングテールなアイテムであり、主要なアイテムがロングテールなアイテムに悪影響を与えることがありうる。
  • Evaluating Ranking via Correlation
    レコメンドシステムはユーザーに対しアイテムのランキングを提供し、上位k個のアイテムが推薦されることとなる。
    ランキングのための評価指標について。
    1.スピアマン順位相関係数
    2.ケンドール順位相関係数
  • Evaluating Ranking via Utility
    効用関数ベースのランキング評価指標。
    $$ F(u, j) = \frac{max\{ r_{uj} – C_{u} , 0 \}}{2^{\frac{v_j – 1}{\alpha}}} $$
    ここで、\( C_{u} \)はユーザーの平均的な評価とされる。

    $$ R-score(u) = \sum_{j \in I_{u} F(u, j) } $$

    ・DCG(discounted cumulatiove gain)
    $$ DCG = \frac{1}{m} \sum_{u=1}^{m} \sum_{j \in I_u, v_j \leq L } \frac{g_{uj}}{log_2 (v_j + 1)} $$
    \( v_j \)はアイテムjの順位を表す。
    $$ g_{uj} = 2^{rel_{uj} – 1 } $$
    ここで\( rel_{uj} \)は実際の評価が使われることが多いとのこと。

    ・NDCG(normalized discounted cumulatiove gain)
    $$ NDCG = \frac{DCG}{IDCG} $$

Limitations of Evaluation Measures

第8章

TBD

第9章

TBD

第10章

TBD

第11章

TBD

第12章

TBD

第13章

TBD

[Python]Music×AnalyticsというイベントでLTをやってきましたレポート

今回はゆるめで、あるイベントのLTに参加してきた感想やら、資料やら動画やらソースコードを共有するだけの記事となっております。

ことの成り行き

Tokyo.Rでお世話になっているy__mattuさんから、「Music×Analytics Meetup Vol.4」のイベントのことを教えていただき、せっかくなのでLTをさせてもらえることになりました。

音楽の分析に関して入門レベルなので果たして参加者が盛り上がってくれる何かを提供できるか不安でした。

納期的にも、今の実力的にも過去に書いたブログ([Python]スペクトグラムで楽器の試奏データを比較(プロ奏者とも))の延長線上のことでいいなと思ったんですが、
先行研究は欲しいし、分析対象の音楽もスケール(音階)ではなく何かの曲が良いと思いました。あと、意外な情報も載せたい、分析結果を活かして便益を得たいとも思いました。

なので、LTまでに先行研究を見つけ、意外な情報としてjupyter notebookで録音する機能を見つけ、分析対象の曲にギャツビーのCMのものを採用しました。
そして、この研究から自分の楽器の演奏の腕前を向上させることができました。
これならLTとして発表しても入門レベルでも楽しいのではないかなと思いました。

LT内容

Google Presentationと事前に録音したYoutubeをこちらに載せておきます。

発表音声の文字起こし

LTを7分以内に抑えるリハーサルをした際の音声があったので、それをGoogle Speech To Text APIを使って文字起こしをしておきます。
「それではその音に近づくための研究と練習ということで江口さん自己紹介ですが名前は救えバックグランドは経済学や統計学仕事では同席部のマネージャー兼データ分析さんはこんな感じで趣味はトランペットブログ作成をやってるので良かったみてください結局どうしたんですけどその他データ分析関連のかって言いよったでどうやって勝手の音消したいなって思いますがここまでのホスピスがまずあって写真撮って進度が速い遅いである時と言うとって治ってるということですが音の長さが変わりますなのでお体変わります秘密のピストンだけなので発送の仕方ありませんのでコントロールして40種類ぐらいのレベルからは死んだことに帰って前記処理学会のそれの大会でスペクトルグラムを使ってプロのアクセスをしたのとの比較ではその方向の出来不出来を判定すると犬が行われていました鬼時間に周波数色の濃淡で信号ってなっちゃったから明日グラフとなりますがこの場合は時間が12秒まであって周波数が8000まで行ってあのね精神旺盛な強さが違ってます本当の鳴き声何もとしてるの俺も行ってここだったから体操マット耐熱耐酸プルマンちょっと寄りますよねプレイヤーの力も貸したして自分との違いを明らかにして違いがある場合とならそこによるそれが産めるのかお話していきたいと思ってます結果こんな感じえーとまずはサイトメトロ ffmpeg をインストールしてああいパイソン mrtg をインストールしたギターの音は愛想の良いしらすを呼び出してこういう形で書くことでノートブックので他のところでエロ金ボタンが現れますこれ音声は録音開始でもう一度録音するようになっていますその後に録音したファイルを保存 SMS でまあ16000hz 2のファイルに変換してまた再読み込みという形があってもはい録音環境としてはわからないぞマイクでして何楽器はちょっと古いので今回あったのがという音楽などの分析見ております玄関のために必要なセッティングってのがあるんでは基本的なセッティングで言っても実際に自分のそれかしたら私に出てきますね自分のやってみようと思いますまず聞いてみましょう第4フォルマント信号の強さがえっと Pro は赤いんですけど強いんですけど出してもらうわいえ私なってます遅く音源第4フォルマントの方がやっぱり行けないなっていうのが分かれば鼻が詰まってるから行けないっていましたねと自分ではまあ剛てる所間違いが分かりましたどうしたら強い信号になってるのかその女お腹に力を入れる行きの数量変えるとかありますもうなんかの力を入れて強いわ強い力を込めてもっと明るくならなかったにやと言われるとも出てこなくなったよ沢山入ったとこですねこれで第4子が赤くなるのはこれまで行きの量が足りてなかった自分じゃなかったじゃないかっていうので行ってみましょう早くなりましたねでスタグラムに従ってくると自分の音が帰宅したところ息を吸う量を多くすることで気分良くなることがわかりました待機するのって管楽器重要な今なってどれぐらい剃ったらいいのか分かんなかったり根性論だったりするんでこれは便利だねまあは管楽器にも分からなくなって鳴るんですけど G 計り方とかあると結婚もしないからっても行きたいのでくれてありがとうございました」

ソースコード(LT)

こちらではLTのために使ったスペクトログラムを表示するためのコードを共有します。

ソースコード(文字起こし)

Google Speech To Text APIで文字起こしをする際のコードをここに載せておきます。

GCSのバケットを作っておきます。今回はperry_voiceという名前のバケットを作っておきました。先ほど作ったaudio_only_filtered_441.wavを、このバケットにあげておきます。

APIの認証用のjsonファイルをcloud shellにアップロードしておきます。jsonファイルはGCPのコンソールの認証情報からサービスアカウントキーを作成し、「鍵を追加」の項目からjsonファイルを選択して手に入れることができます。(当然、Google Speech To Text APIのAPIを有効にしておく必要があります。)
1ヶ月につき60分までは無料なので、今回の処理も無料でできます。(GCSは数円かかると思いますが。)

文字起こし結果はGCSのperry_voiceの中にあります。

他の発表から学んだこと

深層学習を用いたJazz音楽の研究や、音楽の理論、DJプレイの自動化の研究、声の変換技術など同じ音という題材ながらも幅広く話を聞くことができました。私だけ入門な内容で申し訳ない気持ちが生じましたが、これからこの分野を深めるのもプレイヤーとして面白いなと思ったので、長く学ばせてもらおうと思います。

懇親会で学んだこと

・NVIDIAのノイズリダクション技術は文字起こしの際にやってみる価値あり
・賃貸の場合、防音の部屋は家賃相場の2~3万円アップ
・口周辺の筋肉の研究をされている方がいる(論文でその分野があるのは知っていたが、実験をやっていた方に会えた)
・みな、条件を揃えるためにマイクや集音の仕方にこだわりがあった。
・いい音の定義にやはり答えはない。

[Python]仕事で使えそうな確率まとめ


本来、正月の企画物としてアイデアを出していたものの、作成する時間がとれなかったのですが、日常的にもっと確率に従って生活できたら面白いなと思って作ってみることにしました。

二項分布

  • やりたいこと
    オフィスにおいて、チームごとにデスク周辺の掃除当番を毎日乱数で決めているとして、1ヶ月で自分がどれくらい掃除をやらなくてはいけないのかを知りたい。あるいは、1回の掃除で済む確率とか、10回掃除する確率とかを知りたい。

    これは、成功確率p(掃除当番になる確率 = 1÷チームの人数 )でのN回の独立なベルヌーイ試行となり、自分が掃除当番になった場合X=1で、そうでない場合はX=0となる。
    N回ベルヌーイ試行した際のXの合計は二項分布従うので、N回試行して、Xがkになる確率は以下のように示される。(チームの人数はHとする。)

    $$ P(Y_N = k ) = {}_N C _k \left (\frac{1}{H} \right )^k \left ( 1- \frac{1}{H} \right )^{N-k} \\
    (0 \leq k \leq N )$$
  • コード


    平均では20営業日で2回は掃除することになりますが、4回以上掃除しないといけない確率が10%以上はあることがわかります。掃除が嫌いな人は、苦しみの上限をイメージできると楽なのではないかと思うので、確率に従い掃除当番を甘んじて受け入れましょう。
  • 参考情報
    scipy.stats.binom
  • matplotlib.pyplot.grid

多項分布

  • やりたいこと
    先ほどのように乱数で掃除当番を割り当て、それをランキングにし、一位を掃除大臣として称号を付与するものとする。10名いるもののうち、上位のメンバーA・メンバーBとで累積の掃除回数がそれぞれ10回、6回としたもとで次の月に掃除大臣が入れ替わる確率を計算したい。
  • コード

    計算したところ、1.2%くらいになった。このケースにおいてAさんがBさんに次の月に掃除大臣の座を奪われる可能性は低いらしい。
  • 参考情報
    scipy.special.comb

誕生日の問題

  • やりたいこと
    あまりにも有名な確率の話で、誕生日のパラドックスとも呼ばれています。誕生日が同じ人が部署内にいる確率が0.5を超える人数は何人か?
  • コード

  • 参考情報
    誕生日のパラドックス
    matplotlib.pyplot.hlines

ファーストサクセス分布

  • やりたいこと
    アルバイトの面接を受けるとして、その合格率が1/30だとして、何回くらい面接しても合格できないものだろうか。苦しみの上限を知りたい。20回、30回面接を受けた時に合格できてない確率はどれくらいだろうか。何回くらい面接を受けたら合格できていない確率を1%未満に下げれるだろうか。
    成功確率がpのベルヌーイ試行を独立に何回も行うとき、初めて成功数までに失敗した回数をXとして、初めて成功した際の試行回数をYとした場合、その確率は\( k \geq 1 \)として、

    $$ P(Y = k) = P(X=k-1)=p(1-p)^{k-1} $$

    と表すことができる。これをパラメータpのファーストサクセス分布と呼ぶ。(幾何分布とも呼ぶ。)

  • コード

    20回面接をしても合格できない確率はおよそ50%。30回面接しても合格できない確率はおよそ36%。合格できない確率が1%を切るのは面接回数が135回になった時。このレベルを追い求めると心が折れる面接回数なのかもしれない。

  • 参考情報
    弱点克服 大学生の確率・統計

負の二項分布

  • やりたいこと
    ある営業がリンゴを売る仕事をしておりタウンページで飲食店に電話をかけているとする、商談の予約が取れる確率(p)が0.05として、商談予約が5(n)回起きるまでに失敗がk回生じる確率について知りたい。いったい、この営業は飲食店にどれくらい商談をお断りをされる覚悟で臨めばいいのか。

    このような事象は負の二項分布に従うので、
    $$ P(X = k) = {}_{n – 1 + k} C _k p^{n-1}(1-p)^k $$
    という確率に従う。
    なお、この分布の期待値は
    $$ E(X) = \frac{n(1-p)}{p} $$
    に従う。

  • コード


    75回くらいの失敗で確率自体のピークがきている。期待値だと95回だが、確率のピークはもっと手前にある。

    累積確率を知りたい場合はこちらのコード


    80%のところで125回なので、5人に一人は期待値の95回を30回くらいは超えてしまってもおかしくないということになる。そして、5人に一人くらいは50〜60回で達成できるラッキーな人もいる。

  • 参考情報
    弱点克服 大学生の確率・統計

ポアソン分布

  • やりたいこと
    先ほど、営業が飲食店に売ったリンゴの中に芋虫が住み付いているケースが500回に1回の頻度で発生し、営業と顧客の間でトラブルに発展するとする。年間に1200の飲食店と取引を行うが、トラブルはどれくらい発生するものなのか。
    顧客とのトラブルが発生する回数がYだとした際に、ポアソン分布に従うとした場合、
    $$P(Y = k) = \frac{\lambda^k}{k!} e^{-\lambda} $$
    となる。ここで、
    $$ \lambda = 1200 \times \frac{1}{500} $$
    とおく。
    なお、期待値と分散はともに
    $$ E(Y) = V(Y) = \lambda $$
    となる。
  • コード


    トラブルが2回の時が確率が一番高い。0回で済む確率も0.1に近いくらいはある。


    5人に一人は4回以上のトラブルに巻き込まれる気の毒な営業がいそうですね。まぁ、リンゴの中は開けてみないとわかりませんから、顧客トラブルがあっても恨みっこなしで頑張りましょう。
  • 参考情報
    弱点克服 大学生の確率・統計

フィギュア集め問題

  • やりたいこと
    現場猫フィギュアのガチャポンに関心があるとする。N回ガチャポンを回した時に、N種類の現場猫フィギュアが全部揃う確率はどれくらいだろうか?(どの現場猫のフィギュアも出現確率は同じものとする。)
    また、N+1回ガチャポンを回したとき、N種類の全部が揃う確率はどれくらいだろうか?あるいは、m回ガチャポンを回してそろったフィギュアの数の期待値は?そして、Y回目のガチャポンではじめてN種類の現場猫の全部揃う、Yの期待値は?

    ・N回ガチャポンを回した時に、N種類の現場猫フィギュアが全部揃う確率:
    $$\frac{N!}{N^N}$$

    ・N+1回ガチャポンを回したとき、N種類の全部が揃う確率:
    $$ N \times \frac{(N+1)!}{2} \times \frac{1}{N^{N+1}} $$

    ・m回ガチャポンを回してそろったフィギュアの種類\( X_m \)の期待値:
    m回ガチャポンを回してそろったフィギュアの種類\( X_m \)を以下のように表す。
    $$X_m = Z_1 + Z_2 + \dots + Z_N$$
    m回ガチャポンを引いて、当たらないということは
    $$ P(Z_i = 0) = \left ( \frac{N – 1}{N} \right )^m $$
    ということになる。
    これをもとに\( Z_i \)の期待値を計算すると

    $$ E(Z_i) = 1 \times \left ( 1 – \left ( \frac{N – 1}{N} \right )^m \right ) + 0 \times \left ( \frac{N – 1}{N} \right )^m \\
    = 1 – \left( \frac{N-1}{N} \right )^m $$
    よって、\( X_m \)の期待値は
    $$ E(X_m) = E( Z_1 + \dots + Z_N) \\
    = N \left ( 1 – \left( \frac{N-1}{N} \right )^m \right )$$

    となる。

    ・Y回目のガチャポンではじめてN種類のフィギュアが揃う際のYの期待値:
    i-1種類からi種類になるまでのガチャポンの回す回数を\( T_i \)として、その階差がファーストサクセス分布に従うとする。
    $$ T_i – T_{i-1} \sim F_s \left ( \frac{N – i + 1}{N} \right ) $$
    その場合、Yの期待値は以下のようになる。

    $$ E(Y) = E( T_N ) \\
    = E( (T_{N} – T_{N-1}) + (T_{N-1} – T_{N-2}) + \dots + (T_{2} – T_{1}) + T_{1} ) \\
    = E \left ( F_s \left ( \frac{1}{N} \right ) + F_s \left ( \frac{2}{N} \right ) \dots + F_s \left ( \frac{N-1}{N} \right ) + 1 \right ) \\
    = N \left ( 1 + \frac{1}{2} + \dots + \frac{1}{N} \right )
    $$
  • コード

    N回ガチャポンを回した時に、N種類の現場猫フィギュアが全部揃う確率:
    0.0384
    N+1回ガチャポンを回した時、N種類の現場猫フィギュア全部が揃う確率:
    0.1152
    m回ガチャポンを回して揃った現場猫フィギュアの種類Xの期待値:
    3.3615
    Y回目のガチャポンではじめてN種類の現場猫フィギュアが揃う際のYの期待値:
    11.4166
  • 参考情報
    弱点克服 大学生の確率・統計

[Python]スペクトグラムで楽器の試奏データを比較(プロ奏者とも)


以前、トランペットのマウスピース選びで音声解析をしてみましたが、今回は楽器選びで行ってみます。
前回の調査でf1値とかでは良い音の判断がつかなかったので、スペクトログラムを出してみたいと思います。
管楽器はそこそこに高い買い物なので、色々と比較見当をしてから意思決定したいですね。


今回の楽器たちはヤマハの800GS(所持品)、Bach 180ML37SP、Bach 180ML37GL、Bach AL190GLです。ヤマハのやつは20万円くらいですが、Bachのは30万円くらいします。ちなみにBachはアメリカの楽器メーカーで、私の感覚値ではあるのですが、音大生はだいたい持っているような気がします。10万円の差は大きいのか小さいのか。そしてプロの演奏家のスペクトグラムとも比較します。

スペクトログラムとは

パワースペクトルの値の高低を画像の濃淡で表現し、縦軸に周波数、横軸に時間を取った座標上にプロットしたものをスペクトグラムと呼びます。スペクトグラムでは、周波数が低いところ(下)から、第1フォルマント、第2フォルマント、第3フォルマント・・・と呼びます。フォルマントは周りの周波数と比べて突出している部分のことを指します。

$$ S(\omega) = \int_{-\infty}^{\infty}s(t) e^{-jft}dt$$

さて、いきなりの数式ですが、フーリエ変換の式となります。\( S(\omega) \)を二乗したものがパワースペクトルとなります。(二乗するは複素平面上で絶対値を取りたいため。)
tは時刻、fは基本周波数(最も低い周波数)、jは虚数単位を表します。ここで登場する\( s(t) \)は定常波と呼ばれ、以下の式で定まります。

$$s(t) = \frac{a_0}{2} + \sum_{n=1}^{\infty}a_n \cos (2 \pi n ft + \theta_n)$$
\( a_n \)は元の信号\( s(t) \)にどの周波数の波がどの程度の強さで含まれているかという情報を表します。
フーリエ変換のモチベーションとしては、周波数の周期関数を周波数の異なる単純な波の重み付き和で表現できることですが、トランペットの音波を色々な周波数に分解できるのは面白いです。良い奏者の周波数がどうなっているのか調べてみるのは面白そう。

Pythonでの計算&比較

音声データはFFMpegで変換とか切り取りをするなどして作成しました。作成方法は簡単で、以前のブログにも記しましたので、試してみたい方はそちらを参考にしてみてください。

まずはヤマハですが、この楽器になれている分、スムーズに吹けています。最初のほうのリップスラーは5つくらいのフォルマントがありそうです。音階が上がるにつれてフォルマントが減っているのがわかります。最後のHigh Bでは細めのフォルマントが3つあります。


次は、Bach 180ML37SPという楽器です。鳴らすのに慣れていないので、ゆっくり吹きました。最初のほうのリップスラーは5つくらいのフォルマントがありそうで、それはヤマハと同じですが、最後のHigh Bで2つくらいしかフォルマントがありません。楽器に十分息が吹き込めておらず、鳴らせていないのだろうと思われます。可視化できるって面白いですね。


次は、Bach 180ML37GLという楽器です。さっきの楽器との違いは銀メッキがなされているかどうかくらいです。最初のリップスラーは同じ感じで5つくらいのフォルマントがありそうですが、リップスラーで戻ったところで5つ目のフォルマントが消えています。ただ、最後のHigh Bで3つほどフォルマントがありそうです。慣れてきたのでしょうか?


次は、Bach AL190GLという高級モデルです。40万円くらいします。より技術力の高い職人が手掛けた楽器だそうです。ベルに一手間加わっていたり、専用パーツがあったりするみたいです。実際、30万円のよりも吹きやすい印象がありました。
最初のリップスラーは同じ感じで5つくらいのフォルマントがありそうですが、最後のHigh Bでフォルマントが2つになっています。やっぱ鳴らすの難しいのでしょうか。


最後に、プロの方の音源です。最初のBで私のグラフと明らかにフォルマントの数が違うことがわかります。楽器を鳴らせているというのはこういうレベルなのかなと思ってしまいました。High Bでフォルマントが4つはありそうです。このレベルになるにはどういう吹き方をすればいいのか、色々研究のしがいがありそう。

参考文献

イラストで学ぶ 音声認識 (KS情報科学専門書)

[Python]データ分析業務で使いそうなコードまとめ(随時更新)

はじめに

仕事で使いそうなPythonのコードを残しておくドキュメントが欲しいなと思ったので、よく使うものをこちらに貯めていこうと思います。まだ19個しかないですが、30個を目標に追記していきます。

フォーマットとしては、
1.やりたい処理
2.コード
3.参考情報のリンク
の3つを1セットにしています。
まずは、自分自身や周りで仕事をしている人が楽をできるドキュメントになればいいなと思って作っていきます。

目次

・重複削除
・階級のデータを作りたい
・再起的にリストをコピーしたい
・ピボットテーブルでアクセスログから特徴量用のデータを作りたい
・データフレームのある列に対して特定の文字列で分割し、その任意のN番目のものを抽出したい
・データフレームの行単位で割り算を行いたい
・列単位で行を足し合わせて文字列を作りたい
・グループのレコードごとにidを割り当てたい
・正規表現で任意の合致した文字列で挟まれた文字列を抽出したい
・複数列の集計値同士で除した値が欲しい
・等間隔の実数の数列を作成したい
・スペース区切りの文字列において、ある単語の前後N個の語を抽出したい
・グループごとに文字列を結合したい
・文字列の置換をdict形式のデータで行いたい
・テキストを形態素解析してBoW形式のデータにして特徴量にしたい
・訓練データとテストデータで特徴量の数が足りない時に不足した変数を追加したい
・アクセスログのユーザーごとのレコードごとの累積和を計算したい
・選択したカラムで一気にダミー変数にしたい
・グループ化して平均値を計算し、プロットし、2軸でそれぞれのデータ数も載せたい
・参考情報

データセット

全体を通じて使うデータセットを生成する関数を用意しておきます。一つ目の
データはkaggleのBoston Housingのデータセットです。

二つ目はアクセスログのデータが欲しかったので、Google Analytics APIを使ってデータを抽出したものとなります。ここではソースコードを載せますが、皆さんは当然見れません。一応、GitHubにデータを上げておきます。

三つ目はテキストデータです。

前処理

重複削除

  • やりたいこと
    重複したデータを除外したい。最初に現れたレコードだけを残したい。
  • コード
  • 参考情報
    https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.duplicated.html

階級のデータを作りたい

  • やりたいこと
    任意の範囲に応じたラベルをつけた階級値を作りたい。
  • コード
  • 参考情報
    https://pandas.pydata.org/docs/reference/api/pandas.cut.html

再起的にリストをコピーしたい

  • やりたいこと
    リストをコピーした際に、コピーしたものを変更した際に、もとのものも変更されてしまうので、それを防ぎたい。
  • コード
  • 参考情報
    https://docs.python.org/3/library/copy.html

ピボットテーブルでアクセスログから特徴量用のデータを作りたい

  • やりたいこと
    集計されていないアクセスログをもとに、セッションidごとの触れたページを集計し、Bag of Words形式のデータフレームにしたい。これは機械学習の特徴量に使えるので便利。
  • コード
  • 参考情報
    https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.pivot_table.html

データフレームのある列に対して特定の文字列で分割し、その任意のN番目のものを抽出したい

  • やりたいこと
    pandasのデータフレームである文字列の列に対して特定のキーワードで分割し、そのN番目のものをデータフレームの列として持たせる。今回はURLのパラメータの後ろについている文字列を抽出している。
  • コード
  • 参考情報

データフレームの行単位で割り算を行いたい

  • やりたいこと
    各行の合計で、各行を割ることでレコード単位で重み付けを行いたい時や確率を計算したいときに良く使います。データはアクセスログをBag of Word形式したものを使います。
  • コード
  • 参考情報
    https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.div.html

列単位で行を足し合わせて文字列を作りたい

  • やりたいこと
    データフレームの列の要素をつなぎ合わせて文字列を作成したい。複合キーを作らないとデータを繋げない時などによく使います。
  • コード
  • 参考情報
    https://www.kato-eng.info/entry/pandas-concat-pk

グループのレコードごとにidを割り当てたい

  • やりたいこと
    アクセスログなどで、任意のユーザーの初回に触れたページや、N番目に触れたページなどを抽出しやすいように、ユーザーごとのログにidを付与したい。
  • コード
  • 参考情報
    https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.factorize.html

正規表現で任意の合致した文字列で挟まれた文字列を抽出したい

  • やりたいこと
    正規表現を用いて、文字列に対して「指定した文字列と文字列の間にある表現」を抽出するための処理です。URLのデータや規則性のあるテキストデータに対して特徴量を作りたい時に使えます。
  • コード
  • 参考情報
    https://docs.python.org/ja/3/library/re.html

複数列の集計値同士で除した値が欲しい

  • やりたいこと
    複数の列で集計した値同士を用いて、さらにそこから比率のデータを作りたい。ターゲットエンコーディングを行いたい場合に使える。
  • コード
  • 参考情報
    https://stackoverflow.com/questions/45075626/summarize-using-multiple-columns-in-python-pandas-dataframe

等間隔の実数の数列を作成したい

  • やりたいこと
    パラメータチューニングをする際に、各パラメータについて等間隔の実数列を作りたいときがあるので、その時に使うコード。
  • コード
  • 参考情報
    https://www.pythoncentral.io/pythons-range-function-explained/

スペース区切りの文字列において、ある単語の前後N個の単語を抽出したい

  • やりたいこと
    テキストマイニングをする際に、関心のある単語の前後N文字を抽出したい時がある。もちろん、関心のある単語の周辺の単語というものを特徴量にすることも良いと思われる。ここではMeCabによる形態素解析により分かち書きにする関数と共に紹介する。
  • コード
  • 参考情報
    https://stackoverflow.com/questions/17645701/extract-words-surrounding-a-search-word

グループごとに文字列を結合したい

  • やりたいこと
    テキストマイニングなどで、任意のグループごとのテキストデータをスペース区切りで繋ぎたいときがあります。任意のグループごとのテキストを一気に連結できて便利です。これを形態素解析して特徴量にするのも良いと思います。
  • コード
  • 参考情報
    https://stackoverflow.com/questions/17841149/pandas-groupby-how-to-get-a-union-of-strings

文字列の置換をdict形式のデータで行いたい

  • やりたいこと
    事前に用意した置換のリストを使って、文字列の置換を一気に行いたい。
  • コード
  • 参考情報
    http://omoplatta.blogspot.com/2010/10/python_30.html

テキストを形態素解析してBoW形式のデータにして特徴量にしたい

  • やりたいこと
    任意の日本語のテキストとラベルの組み合わせがあるとして、テキストを形態素解析してTF-IDFを計算し、それを特徴量としてラベルを予測するということをやりたい。ここでは以前集めて、GitHubに載せてあるデザイナーズマンションのデータを使い、MeCabで形態素解析をしたあとに、TF-IDFを計算し、それを特徴量にしてsklearnで簡単に予測をし、リフト曲線を描いている。
  • コード
  • 参考情報
    https://scikit-plot.readthedocs.io/en/stable/metrics.html

訓練データとテストデータで特徴量の数が足りない時に不足した変数を追加したい

  • やりたいこと
    訓練データとテストデータの変数の数が違う時に、そのままではモデルを使えないので、補うために使う。
  • コード
  • 参考情報
    なし

アクセスログのユーザーごとのレコードごとの累積和を計算したい

  • やりたいこと
    アクセスログや購買データなどで、レコードごとの累積の値を計算したいときに使う。時系列に従いどんどん足されていく。
  • コード
  • 参考情報
    https://ohke.hateblo.jp/entry/2018/05/05/230000

選択したカラムで一気にダミー変数にしたい

  • やりたいこと
    ダミー変数を一気に生成したい。名前も一気に付けたい。
  • コード
  • 参考情報
    なし

グループ化して平均値を計算し、プロットし、2軸でそれぞれのデータ数も載せたい

  • やりたいこと
    EDAの時に、グループ化して平均値などを計算するだけでなく、該当するデータ数も2軸でプロットしたい。
  • コード
  • 参考情報
    なし

参考情報

前処理大全[データ分析のためのSQL/R/Python実践テクニック]
はじめてのアナリティクス API: サービス アカウント向け Python クイックスタート
カモノハシペリー
ペリー|フィニアスとファーブ|ディズニーキッズ公式

Uncertainty in Gradient Boosting via Ensembles[A4一枚まで備忘録]

A4用紙1枚にまとめるイメージでメモを残そうという取り組み。

今回は機械学習でとりわけ、回帰の問題を扱うことが仕事であった際に参照できるメモの一つとして残したい。
先日、Catboostのドキュメントに参考文献としてUncertainty in Gradient Boosting via Ensemblesの存在を知った。業務において、価格の予測とか、期間の予測とかマーケティングのデータでも良く扱うものなので、分類だけでなく回帰についてもバランス良く知見を貯めたいところ。

目次

Abstract
Introduction
Background and related work
Generating ensembles of GDBT models
Experiments on Synthetic Data
Experiments on classfication and regression datasets
Conclusion
Broader Impact

Abstract

機械学習のモデルは予測において点推定しかしないが、実務上は不確実性を定量化したいニーズがある。この研究では、予測の不確実性推定値を導出するためのアンサンブルベースの枠組みが紹介されている。予測における不確実性は「データに起因する複雑さとノイズの不確実性」「特徴空間において、情報がないことに起因する不確実性」に分解される。
手法としては確率的勾配ブースティングと確率的勾配ランジュバンブースティング(SGLB)が扱われている。SGLBは一つの勾配ブースティングモデルからのアンサンブルを可能にし、モデルの複雑さを減らすことができる。また、予測性能の向上について様々なデータセットを使って示されている。

Introduction

一般的な勾配ブースティングについて、回帰の際に点推定しか返さないことを挙げている。NGBoostについても言及している。NGBoostなら、勾配ブースト決定木(GBDT)を訓練して平均値を返せるとのこと。ただし、データの不確実性に起因する不確実性しか捕捉できないという問題がある。モデルに関する知識の欠落によって生じる不確実性に関して扱うことはできない。近年ではニューラルネットワークを用いて知識を起因とする不確実性の研究がなされているらしい。今回扱う、アンサンブルアプローチは全体の不確実性を「データの不確実性」と「知識の不確実性」に分解する。扱うモデルは確率的勾配ブースティングと、確率的勾配ランジュバンブースティング(SGLB)となっており、SGLBの利点として、漸近的に真のベイジアン事後分布からサンプリングされることが挙げられている。また、SGLBによって1つの勾配ブースティングのみをモデルを使ってアンサンブルをすることができるとされている。これにより計算の複雑さを減らすことができる。結果として、予測性能が高まり、誤りやドメインの外側のデータの入力の検知ができるようになったとされている。

Background and related work

  • ベイジアンのアプローチ
    ベイズ統計学を使っている。予測の多様性について考察することで、知識の不確実性の推定を行うことができるとしている。さらに、確率モデルのアンサンブルを行なっている。アンサンブルモデルの不一致さ、広がりの程度が知識の不確実性につながるとしている。事後分布の近似のための研究が色々となされ、ニューラルネットワークのモデルなども使われた。モデルの事後分布をもとに、予測の事後分布を手に入れる。予測の事後分布のエントロピーが、予測の不確実性を推定する。知識の不確実性は、全ての不確実性からデータの不確実性の期待値を差し引いた相互情報量で表される。確率的な回帰モデルのアンサンブルを考える時、予測事後分布のエントロピーと相互情報量を推定することができない。代わりに総変動の法則に従ってなら推定することができるが、1~2次のモーメントに基づいているだけなので、高次元での不確実性に関する見落としが起きうる。
  • 勾配ブースティング
    弱学習器は浅い決定木によって形成される。ただ、勾配ブースティングの回帰モデルは点での予測しかしない。予測の不確実性に関する研究はほとんどされていない。近年、NGBoost (Natural Gradient Boosting)が登場し、予測の不確かさを推定することができるようになった。NGBoostは入力変数xが与えられたもとでの、出力変数の条件付き事後分布の推定を行ってくれる。

Generating ensembles of GDBT models

NGBoostはデータの不確実性を捉えることができるものの、知識の不確実性については捉えることができない。
知識の不確実性を捉えるための方法としては、事後分布からサンプリングされた予測モデル(GBDT)のアンサンブルが挙げられる。
そもそもGBDT自身が木のアンサンブルであるが、そのGBDTのアンサンブルを考える。

アンサンブルをする方法としては、確率的勾配ブースティングのアプローチで生成された独立なモデルを扱うことがあげられている。
もとのモデルにノイズを注入して、出力をぶれさせることで予測の分布を得ようとしている。ただ、これは事後分布の推定を保証しないらしい。

ただし、確率的勾配ランジュバンブースティング(SGLB)アルゴリズムを用いることで真の事後分布からGBDTモデルをサンプリングすることができるらしい。
非凸の損失関数においての最適化の収束のために、確率的勾配ランジュバン力学が用いられているとのこと。ランジュバン力学ってなんだろうか。

確率的勾配ブースティングとの違いは2つある。

  • 1.ガウシアンノイズが勾配に明示的に加えられている。
    • 正規分布の平均は0で分散に逆拡散温度というものがある。ノイズをコントロールするものであることに間違いはない。確かに物理っぽい。STANのリープフロッグとかを思い出した。このノイズが最適値を見つけるのに役立つらしい。
  • 2.ブースティングモデルの更新式に正則化パラメータを設けている。
    • 定常分布に弱く収束していくマルコフ連鎖にパラメータが従うらしい。

対数尤度に正則化項が残るものの、比例するらしく、真の事後分布に従うと見なせるらしい。

確率的勾配ブースティングや確率的勾配ランジュバンブースティングは計算する上で重大なオーバーヘッドが生じるらしくイタレーションを減らすなど質を下げたり複雑さを増したりする必要があるとのこと。それを避けるために、仮想SGLBアンサンブルというものが提案されている。

SGLBアルゴリズムのイタレーションごとのパラメータは弱く定常に収束するマルコフ連鎖に従うことから、それらをモデルのアンサンブルと見なす。
しかしながらモデル同士が強く相関してしまうので、アンサンブルの質が下がる。ただし、各々のパラメータのK番目の集合のみを使うことで克服できるとある。
K番目のパラメータ集合のみのモデルで構築したアンサンブルを仮想アンサンブルと呼んでいる。

Experiments on Synthetic Data

ここではGBDTモデルのアンサンブルアルゴリズムを合成データを用いて実験している。
合成データとしては、以下の2つの1次元データセットが扱われている。

  • 合成データセット
    • 分散が増加していく絶対値x
    • 分散が減少していく絶対値x
  • 比較するモデル
    • SGLB
    • SGB
    • vSGLB(仮想SGLB):SGLBから生成
  • 結果
    • 分散が減少していく絶対値xのケースだと、データの不確実性に関して、SGBではサンプルの領域外で過剰に不確実性が増してしまう。
    • 分散が増加していく絶対値xのケースだと、知識の不確実性に関して、SGBではサンプルの領域外で不確実性が増してしまう。

サンプルの領域外のほうが、領域内よりも不確実性が大きいはずだが、実験データではそうではないらしい。

Experiments on classfication and regression datasets

GBDTモデルのアンサンブルアルゴリズムを色々なデータのいろいろなタスク(回帰、分類)で予測性能をみている。

  • タスク
    • 回帰タスク
    • 分類タスク
  • データセット
    • 18種類くらいある。(Wineとか、解約のデータとかAmazonのデータとかいろいろ)
  • 比較するモデル
    • SGB:このアンサンブルは10個の独立したモデルで、それぞれ1K本の木を持つ。
    • SGLB:このアンサンブルは10個の独立したモデルで、それぞれ1K本の木を持つ。
    • vSGLB(仮想SGLB):[500, 1000]の区間から50番目のモデルをアンサンブルに追加していっている。
  • どう実装しているか
    • CatBoostライブラリをベース
    • 分類では二値分類のクラスの確率分布を生成
    • 回帰では正規分布の平均と分散を生成
    • 負の対数尤度を最適化している
    • ハイパーパラメータはランダムサーチ
  • 何の指標を見ているか
    • 負の対数尤度
    • 誤分類率
    • RMSE
    • 不確実性の合計、知識の不確実性:Prediction-Rejection Ratio(PPR)というもので測る。どれくらい不確実性が誤差と順番と相関するかを測るらしい。100が最大で、ランダムだと0になる。
    • AUC:領域外のデータに関する評価に用いる。領域外データは合成して生成した。
    • 分類:負の対数尤度(NLL)
      • アンサンブル手法に軍配が上がるが、違いは軽微である。
      • 仮想アンサンブルが勝るときがある。
    • 分類:誤分類率
      • アンサンブル手法に軍配が上がるが、違いは軽微である。
    • 回帰:RMSE
      • アンサンブル手法に軍配が上がるが、違いは軽微である。
    • SGBやSGLBはパフォーマンスが良い。
      • vSGLBはわずかにパフォーマンスが悪い。自己相関があるなかでアンサンブルしているためパフォーマンス悪化が生じていると考えられる。
        • ただ単一モデルのものよりかは良い。
        • コストをかけずに精度を上げたいという観点ではvSGLBはよい選択肢となる。

    Conclusion

    この論文では、GBDTを用いたアンサンブルベースの不確実性推定モデルが扱われた。SGB、SGLB、仮想SGLBについて合成データや、分類や回帰のタスクにおいて比較した結果、アンサンブル手法が精度面で秀でているわけではないことが示された。
    しかし、知識の不確実性(モデルで十分に想定できていない不確実性)を扱う際にアンサンブル手法が役立つことがわかった。
    知識の不確実性を扱う場合に、計算コストの少ない手法として仮想のSGLBが提案され、精度において十分ではないが、低コストで知識の不確実性を扱う観点において有用であることが示された。

    Broader Impact

    • 不確実性の推定は、機械学習モデルをより安全に、より信頼できるものにするという観点で重要。
    • GBDT系の手法はニューラルネットワーク系の手法よりも計算コストが低くて良い。
    • 不確実性の推定が貧弱なものとならないように引き続き研究が求められる。

    追記

    不確実性について扱われた、ソースコードの載っている記事を集める。