2023年7月更新
- Robynのバージョンが3.10系になったことで、内容を更新しました。
- 一部の関数の引数及び挙動が変更になっており、更新前の情報では動かない箇所があります。
Meta(旧Facebook)社の開発したAutomated Marketing Mix Modeling (MMM)ツールの一つRobynの紹介と実践をします。
良くあるエラー
nevergradがない (macで発生) Robyn ver 3.10.3
Error in robyn_mmm(hyper_collect = InputCollect$hyperparameters, InputCollect = InputCollect, :
You must have nevergrad python library installed.RのversionがM1/M2チップ用を使用していると、出る可能性あり。「R-4.3.1-x86_64.pkg」を使用すれば多分出ない。
1. MMMの復習
以前こちらでMMMについて扱いました。改めて復習をします。
MMMはマーケティングの費用対効果を定量化する手法の総称とされています。具体的には収益に与える各メディアの効果、それ以外の効果をそれぞれ定量化するなどです。下図はMMMのイメージです。Robynの公式HPより引用しています。
MMMの実践は複数のステップからなります1https://facebookexperimental.github.io/Robyn/docs/analysts-guide-to-MMMhttps://github.com/google/lightweight_mmmなどを参考に筆者作成。
MMMでは成果変数(売上やコンバージョンなど)への寄与要因を以下のように仮定しています。
| 大項目 | 項目 |
|---|---|
| ペイドメディア | 4マス(TV/ラジオ新聞/雑誌) |
| ネット(google/Yahoo) | |
| SNS(Facebook/Twitter/Instagram) | |
| 非ペイドメディア | オウンドメディア |
| プレスリリース | |
| SNSでの発信 | |
| イベント変数 | キャンペーンなど |
| 競合キャンペーン | |
| スポーツや音楽イベントなど | |
| ベースライン | トレンド |
| 季節性 | |
| 祝休日 | |
| 気候情報 |
以下は入力データの一例です。
2. Robyn実践
RobynはMeta(Facebook)のメンバーが作成した自動MMMツールです。本稿では簡単なデータを用いて広告効果の定量化と予算の最適化を実践していきます。色々な機能がありますが、今回は上記に必要な機能に絞って解説します。また余裕があれば別記事でオプション機能の紹介をしたいと思います。なお本実践は主に、meta公式チュートリアルを参考にしています。
注意事項) WordPressの仕様でRの代入[<-]記号が[<-]に化けてしまっています。読み替えをお願いします。
ここで使用するデータとコードはhttps://github.com/mitsu666/Robyn_demo からDLできます。
また実行環境はmacです。windowsの場合そのまま適用できない箇所があることをご了承ください。
2.1 Rをインストール
この手のライブラリでは珍しく、ラッパーがRしかないとのこと。Rの最新版をここからDlしてインストール。 22022年3月段階ではR4.05以上を推奨
2.2 Robynのインストールとそれに必要な準備
まず作業ディレクトリ(ここでは./MMM_01を作成しそこを作業ディレクトリとした)に対象となるデータを置いておきます。次に以下のコードをRのコンソール(またはRstudio
3RのIDE、ほとんどのRユーザが使用しているため、ここでも推奨。https://www.rstudio.com/products/rstudio/)上から入力し実行していきます。
# 作業ディレクト
setwd("./MMM_01")
# インストール install.packagesは最初の一回のみでOK
install.packages("Robyn")
install.packages("reticulate")
# libarary
library(reticulate)
library(Robyn)
packageVersion("Robyn")
# library(dplyr) # dplyrがない場合はinstallして読み込んでくださいRobynはモデリングに際してNevergra4gradient-freeオプティマイザー、つまり導関数を必要としない最適化ソルバー。https://facebookresearch.github.io/nevergrad/という最適化のソルバーを使用します。そのためにRからPythonを使用できるようにしておきます。Pythonを使用するためにはreticulate
5reticulate https://techblog.nhn-techorus.com/archives/8329というライブラリを使います。
nevergradを実行するpython環境を作成します。
# pythonの環境を準備
virtualenv_create("r-reticulate")
py_install("nevergrad", pip = TRUE)
use_virtualenv("r-reticulate", required = TRUE)
use_python("[あなたのhomeディレクトリ]/R_work/.virtualenvs/r-reticulate/Scripts/python.exe")
# Set to FALSE to avoid the creation of files locally
create_files <- TRUEuse_pythonの引数にはpythonの実行ファイルがある場所を指定します。
- Windowsの場合:homeディレクトリ以下の「R_work/.virtualenvs/r-reticulate/Scripts/python.exe」
- macの場合:「~/.virtualenvs/r-reticulate/bin/python」
にあるはずです。
2.3 入力データ
入力データは、筆者が以前書籍
6https://www.amazon.co.jp/dp/B098HZ9KXG/ref=dp-kindle-redirect?_encoding=UTF8&btkr=1執筆用に用意したデータを用います。人工的に作成した広告の投下費用や売上などからなる日次の時系列データです。
| 変数 | 内容 | 仮説 |
|---|---|---|
| 日付 | 日時 | – |
| ネット広告出稿金額 | その日のネット広告出稿金額 | 売上に比例 |
| TV広告出稿金額 | その日のTV広告出稿金額 | 売上に比例 |
| 平均気温 | その日の平均気温 | 売上に比例 |
| 降水量 | その日の降水量 | 売上に反比例 |
| 売上金額 | その日の売上げ | – |
| 週末FLG | 週末であるならば1それ以外0の変数 | 1の場合売上増加 |
2.4 データ読み込み
# 入力データのcsvを読み込む
dt_simulated = read.csv(file='sample.csv')
# change var names 日本語はplotした時に文字化けするため半角アルファベットに変換
dt_simulated = dplyr::select(dt_simulated,datetime=1, Net.Spend=2 ,Tv.Spend=3 ,temperature=4,
rain=5, revenue=6, Weekend.FLG=7,dplyr::everything())
head(dt_simulated)
#祝日用のデータも読み込み
data("dt_prophet_holidays")
head(dt_prophet_holidays)祝日用のデータはhttps://github.com/facebookexperimental/Robyn/tree/main/R/data
からダウンロードしておきます。
最後に出力先フォルダのpathを設定します。
# 結果を格納するフォルダを指定する
robyn_directory <- "./results" # 先にカレントワークディレクトリにresultsフォルダを作成しておく2.5 input定義
2.5.1 入力変数指定
#### 2a-1: 入力変数を指定
## メディアとオーガニック変数の係数は正になる
## 他の変数は正負、係数符号はカスタマイズ可能
#input
InputCollect <- robyn_inputs(
dt_input = dt_simulated #入力する元データ
,dt_holidays = dt_prophet_holidays #祝日データ
,date_var = "datetime" # 以下のフォーマットで"2020-01-01"
,dep_var = "revenue" # 売上高やCVなどの従属変数
,dep_var_type = "revenue" # 売上高かCVフラグか
,prophet_vars = c("trend", "season") # "trend","season", "weekday" & "holiday"
,prophet_country = "US"# 国名 祝日 デフォルトで日本がないため一旦USとしておく
,context_vars = c("temperature", "rain","Weekend.FLG") # イベント情報
,paid_media_spends = c("Net.Spend","Tv.Spend") # メディア投下
,paid_media_vars = c("Net.Spend","Tv.Spend") # メディア
# paid_media_vars must have same order as paid_media_spends. Use media exposure metrics like
# impressions, GRP etc. If not applicable, use spend instead.
#,organic_vars = c("newsletter") # PRなどの非広告メディア
,factor_vars = c("Weekend.FLG") # イベント因子
,window_start = "2019-04-01" # モデル構築に使用するデータの開始日
,window_end = "2020-03-31" # モデル構築に使用するデータの終了日
,adstock = "geometric" # adstockの形状
)
print(InputCollect)回帰係数の符号はデフォルトでは
- メディア / オーガニック: 正
- 他の変数: 正負どちらも
paid_media_signs = “negative”のように指定することでカスタマイズ可能です。orgnic_signs / context_signsでオーガニック変数や他の変数も指定可能です。
今回行うデータ特有の注意事項
- 今回説明変数にweekendFLGがあるため、prophet_varsでは[weekday]という項目を外している。
- 日本の祝日がデフォルトでは使用できないため、prophet_varの設定で[holiday]を外している。日本の祝日の使い方は別途で解説するつもりである。
そのほかの細かいinputの細かい定義は公式HPをご参照してください。(どこかで整理するかも)
Prophet7metaが開発した時系列予測フレームワークに関連した変数があることからわかるように、backendとしてProphetを一部使用しているようです。Ptophetについての記事はこちらをご参照ください。
2.5.2 ハイパーパラメータ指定
ハイパーパラメータは広告効果の逓減(saturation)を制御するパラメータとAdstockと呼ばれる広告効果の減衰を制御するパラメータの2つの種類があります。前者のパラメータは[alpha]と[gamma]です。後者は減衰のタイプが幾何型とワイブル型があり、それぞれ幾何型が[theta]を、ワイブル型が[scale][shape]というパラメータを持ちます。
使用するハイパーパラメータを表示し、その可変領域を設定していきます。この例では幾何型(geometric)の減衰を指定したので、設定するのはthetaとalpha/gammaの3つとメディアがTVとNetの2つで計6つです(2*3=6)。
#### 2a-2: Second, ハイパーパラメータの定義
#使用するハイパーパラメータ表示
hyper_names(adstock = InputCollect$adstock, all_media = InputCollect$all_media)
## adstockとsaturationの形状を確認(TRUE)にする
plot_adstock(plot = TRUE)
plot_saturation(plot = TRUE)
#上記のハイパーパラメータの可変領域を設定
hyperparameters <- list(
Net.Spend_alphas = c(0.5, 3)
,Net.Spend_gammas = c(0.3, 1)
,Net.Spend_thetas = c(0, 0.3)
,Tv.Spend_alphas = c(0.5, 3)
,Tv.Spend_gammas = c(0.3, 1)
,Tv.Spend_thetas = c(0.1, 0.4)
,train_size = 0.7 # ts-cvで使用する 固定
)plot_adstock(plot = TRUE)
plot_saturation(plot = TRUE)とするとハイパーパラメータごとに減衰やサチュレーションの様子を見ることができます。
上記はオプションなので、グラフを表示させる必要がない場合は「plot=FALSE」にしてください。
(参考) 各種定義の詳細8https://facebookexperimental.github.io/Robyn/docs/variable-transformations
モデル式
Ridge回帰モデル
Adstockやsaturationの式
- 幾何型Adstockは decay_t,j := x_t,j + θ * decay_t-1,j
- ワイブル型は decay_t,j := 1-(exp(decay_t-1,jj/α))
- saturationは c * (x^α / (x^α + γ^α )) cはメディアxに対して得られた回帰係数
2.5.3 input
#### 2a-3: Third, robyn_inputsにハイパーパラメータを入力
InputCollect <- robyn_inputs(InputCollect = InputCollect, hyperparameters = hyperparameters)
print(InputCollect)2.5.4 モデル構築
モデルを構築していきます。マニュアルでは初期モデルの構築とあるように、この後追加データがある場合それを加えて初期モデルを修正することもできます。
https://facebookexperimental.github.io/Robyn/docs/demo-R-script にある方法とDemoのhttps://github.com/facebookexperimental/Robyn/blob/main/demo/demo.R の方法ではやり方が異なるようで、ここではdemo.Rの方を参考にしています。
## Run all trials and iterations. Use ?robyn_run to check parameter definition
OutputModels <- robyn_run(
InputCollect = InputCollect, # feed in all model specification
cores = NULL, # 並列計算を行う場合、cpuの使用コア数を指定
iterations = 2000, # 2000 recommended for the dummy dataset with no calibration
trials = 5, # 5 recommended for the dummy dataset
ts_validation = TRUE, # NRMSEを基準にcv FALSEだとなし
add_penalty_factor = FALSE # Experimental feature. Use with caution.
)
print(OutputModels)ハイパーパラメータの設定でtrain_size=0.7とすると、ts_validationの訓練、バリデーション、検証に使用されるデータセットの割合は7:1.5:1.5となります。
2000回のiterationsを行う試行を5回ここでは行います。結構時間がかかります。
厳密な説明は難しいのですが、ここでのモデリングはリッジ回帰をモデリングしているようで、ハイパーパラメータと通常のパラメータを同時に学習しているようです。この時予測誤差とチャネルのそれぞれの投下費用のシェア率とその係数の比率の距離(要するに現実離れした解にならないように)Prophetに関係する部分のlossと、通常のリッジ回帰の部分のlossに対してMulti-objective optimization
を行なっています。Nevergradはこの最適化に対して用いられているようです。
パレート最適解
最後に作成したモデル群の中からパレート最適な物を選ぶ処理?を行うようです。正直ここは少し理解できていません。
# 収束状況をplot
OutputModels$convergence$moo_distrb_plot
OutputModels$convergence$moo_cloud_plot
左下の赤い線がパレートフロンティアを表しています。
クロスバリデーションによる誤差の確認
# クロスバリデーションplot
if (OutputModels$ts_validation) OutputModels$ts_validation_plot
パレート最適化の計算
## パレート最適解の計算, 結果とplotを指定フォルダにexport
OutputCollect <- robyn_outputs(
InputCollect, OutputModels,
pareto_fronts = 1, # automatically pick how many pareto-fronts to fill min_candidates (100)
# min_candidates = 100, # top pareto models for clustering. Default to 100
# calibration_constraint = 0.1, # range c(0.01, 0.1) & default at 0.1
csv_out = "pareto", # "pareto", "all", or NULL (for none)
clusters = TRUE, # Set to TRUE to cluster similar models by ROAS. See ?robyn_clusters
export = create_files, # this will create files locally
plot_folder = robyn_directory, # 指定したフォルダ内に結果ができる。当該フォルダがない場合cwdにできる
plot_pareto = create_files # Set to FALSE to deactivate plotting and saving model one-pagers
)
print(OutputCollect)pareto_fronts = 1は、NRMSEとDECOMP.RSSDのトレードオフを考慮した最良のモデルを返します。pareto_frontsを増やすことで、より多くのモデルの選択肢を得ることができます。pareto_fronts = “auto”は、少なくとも100の候補を含む最小のフロントを選択します。この閾値をカスタマイズするには、min_candidatesの値を設定します。
また、clusters=TRUEとすると、ROIが似たもの同士を同じクラスタにまとめます。kの値は自動で決定されますが、多くの場合、3になるようです。
実行後、はじめに定義したresultsのフォルダパス傘下に以下のファイルが出力されます。X_YYY_Z.pngの文字列からなるファイルはそのモデルのone-pager repotです。(おそらくですが、先頭のXが何番目のtrialかを表しているようで、残りの2つの数字はiterationsの何番目に相当するかを表しているようです。(内部のpgmを読むと(YYY-1)*コア数+Z=itaretionの何番目かとなっている)。
候補が3つあるのはクラスタ数k=3が指定されているからのようで、そのクラスタのベストソリューションが表示されるようです。
ne-pager report(ここでは上図の例として1_208_3の)を見てみます。
- 左上の図:ウォーターフォール図 成果変数への全体寄与を1とした場合の各要素の寄与割合を示した図です。一番寄与が高いのが「平均気温」、次が「trend」、その次が「週末かどうか」です。ここではTVの寄与はNetより大きいと示されています。メディアの寄与の比較はあくまで投下量ベースであるため、効率は考慮されていません。
- 右上の図:時系列折れ線グラフ 実際の値と予測値のプロット図です。概ね実績値と予測値が近い動きをしていることがわかります。train/valid/testそれぞれの場合が示されています。
- 上から2番目左の図:横棒グラフとROIの折れ線グラフ メディアの実際の投下シェアと、成果変数のシェアのグラフです。赤い折れ線はROIです。全体比の寄与ではTVの方が大きいのですが、ROIではNetの方が大きいことが分かります。成果指標をconversionにすると、ROIの代わりにCPAが出力されるようです。
- 上から2番目右の図:ROIの信頼区間 ブートストラップ法で作成したROIの信頼区間。
- 上から3番目左の図:広告効果の減衰グラフ 広告効果の減衰の程度をそれぞれのメディアごとに示しています。幾何的減衰(geometric)の場合、θの値そのもののようです。
- 上から3番目右の図:効果の即時性とcaryover効果の割合 広告効果の即時性と継続性の比率を横積み上げ100%グラフで表したものです。
- 左下の図:saturationのプロットグラフ saturation効果をそれぞれのメディアで記述しています。Netの方が飽和しにくい、つまり効果の逓減が小さいことが示されています。
- 右下の図:残差プロット 実績と予測値の残差プロットです。一般論として残差に偏りがあれば説明変数だけで説明できていない要因があると考えられます。
本ツールで算出しているROIの定義に関しては調査中ですが一般的な分子から投下量を控除する定義ではなく、成果/投下量で定義しているようです10https://github.com/facebookexperimental/Robyn/blob/main/demo/schema.R?fbclid=IwAR34C6gkz9oeWjX4cPxOUToTIoIvNfopUh0CmZQnqSN5DZX8pRtuSx0joOY。
2.6 モデルの選択
前項のフォルダの中から、ビジネスの実情に合ったモデルを一つ選びます。後の最適化で使用するためにjson形式でモデルをexportします。
################################################################
#### Step 4: モデルを選択と保存
select_model <- "1_208_3" # select one from above
#### Version >=3.7.1: JSON export and import (faster and lighter than RDS files)
ExportedModel <- robyn_write(InputCollect, OutputCollect, select_model, export = create_files)
print(ExportedModel)
myOnePager <- robyn_onepagers(InputCollect, OutputCollect, select_model, export = FALSE)2.7 予算の最適化
様々な条件で最適解を求めることができます。ここでは、元のデータで与えられた予算内で売上を最大化する解を求めます。
################################################################
#### Step 5: 前項で指定したモデルで予算最適化を行う
# Check media summary for selected model
print(ExportedModel)
# NOTE: The order of constraints should follow:
InputCollect$paid_media_spends
# Scenario "max_response": 与えられた予算での売上最大を達成する予算を算出
# Example 1: date_rangeを特に指定しない場合直近1カ月の予算最適化を行う
AllocatorCollect1 <- robyn_allocator(
InputCollect = InputCollect,
OutputCollect = OutputCollect,
select_model = select_model,
# date_range = NULL, # Default last month as initial period
# total_budget = NULL, # When NULL, default is total spend in date_range
channel_constr_low = c(0.7, 0.7), #メディア数と同じ長さが必要
channel_constr_up = c(1.2, 1.2), #メディア数と同じ長さが必要
# channel_constr_multiplier = 3,
scenario = "max_response",
export = create_files
)
# Print & plot allocator's output
print(AllocatorCollect1)
plot(AllocatorCollect1)実行するとresultsフォルダ内に、1_208_3_reallocated_best_roas.pngという画像ファイルができます。
Boundedは、「channel_constr_low/high」で指定した上下限の範囲内での最適解を示します。一方、Boundedx3は、その上下限の3倍まで許容した場合の最適解を示します（おそらくですが、間違っていたらすいません）。 「channel_constr_low/high」を指定する理由は、直近の予算と大きく乖離する解は実務上使いにくいという実情に対応するためです。
他にも色々な条件下での最適化方法がありますが、そちらはこちらを参照ください。
まとめ
meta社のRobynを用いたAutomated MMMの実践を行いました。入力データとハイパーパラメータを与えることで売上とメディアの関係をモデリングすることができます。そして得られたモデルを用いて最適化予算のアロケーションを行うことが可能です。
更にここでご紹介した方法以外にも、追加データでモデルを構築するrefreshやおよモデルのcaribrationなどの機能があります。ここでご紹介できなかったものに関しては、今後どこかで解説記事を書きたいと思います。
一方でAutoと言いつつもかなり、細かいパラメータの設定を作業者が行う必要性があります。
Robynを始めとしたMMMに関してのご相談は是非弊社までお願いいたします。
当ツールは今でも改良が行われているようなので、最新版に関しての情報は公式HPをご確認ください。
