elsur.jpn.org >

« 2016年3月 | メイン | 2016年5月 »

2016年4月26日 (火)

　ベイズ統計とマーケティング・リサーチ、というようなキーワードでウェブを検索すると、Allenby兄貴の著書"Bayesian Statistics and Marketing"ばかりがヒットする。ちょっと読んでみるかと思ったそこのあなた、別にかまいませんけど、その道は地獄へと通じています。香港映画界の最後のカンフースターであるドニー・イェン兄貴は、共演者に対する要求水準が大変高く、ほぼ無理難題に近いアクションを求めるそうなのだが、Allenby兄貴の本も、途中からいきなりハードルを上げてくるのです。Allenbyさんが世界中の統計学習者に"Ah-Niki"と呼ばれているのはそのためである。(嘘です、私が言ってるだけです)
　Allenby兄貴以外の人が書いた文章がないだろうかと探してみたら、Marketing Research誌の記事のPDFをみつけた。こりゃいいや、と思ってディスプレイ上で目を通した。この雑誌は完全に実務家向けで、めちゃくちゃ易しいのです。
　
　まずはベイズの定理について紹介。ベイズの定理が市場調査実務にあまり適用されなかった理由は3つ：(1)事前確率が決めらんない。(2)計算が大変。(3)リサーチャーはデータ生成過程に関心を持たない場合が多い(例外は消費者選択。リサーチャーは消費者の効用最大化過程に関心を持つ)。
　80年代のMCMC登場がすべてを変えた... [正直に告白すると、このあたりから、あれおかしいな、と気づき始めていました...]
　階層ベイズモデルとは... [略]
　選択型コンジョイント分析に適用すると...[略]
　うちの院生が、消費者の決定ルールを推測するという面白い研究をやってて... [ここで完全に気が付いて、ガックリしたのだが、悔しいので最後まで読んだ]
　WinBUGSとかあるよ... これからはこんな問題にも適用できるだろう... [略]

Allenby, G.M., Bakken, D.G., Rossi, P.E. (2004) The HB revolution: How Bayesian methods have changed the face of marketing research. Marketing Research, Summer 2004.
　というわけで、よくよくみたら、これもAllenby兄貴の記事であった。それもGrover&Vriens(eds)での解説とだいたい同じ内容。くっそー。
　悔しいので、今回はドニー兄貴の動画は載せません。(別に誰にも頼まれてないけど)

論文：データ解析(2015-) - 読了：Allenby, Bakken, Rossi(2004) HB革命とは何か

Koenker, R., & Hallock, K. (2001) Quantile regression. Journal of Economic Perspectives, 15 (4), 143-156.
　分位点回帰についての素人向け解説。こりゃ基礎から復習しないとだめだと気づき、あわてて読んだ。Koenkerって、やっぱ偉い人だったんだ...

　いわく。
　みんな知ってると思うけど、平均は残差平方和
　$\min_{\mu} \sum(y_i - \mu)^2$
を最小化するという問題の解だよね。同様に、中央値は絶対誤差の合計を最小化するという問題の解だよね。他の分位点に一般化すると、それらは次の解だ:
　$\min_{\xi} \sum \rho_t (y_i - \xi)$
ここで$\rho_t (\dots)$は左右で異なる傾きを付けた絶対誤差の関数。[←ここ、理解できなくて混乱したんだけど、"Mostly Harmless Econometrics"によれば、厳密な定義はおいといておおまかにいうと、$u$が正のときに$\rho_\tau = \tau u$, 負のときに$\rho_\tau = (\tau-u) u$となる関数らしい。$\tau=0.5$なら$\rho_\tau = 0.5|u|$。なるほどね]

　これを条件つき平均に拡張しよう。$\mu$をあるパラメトリックな関数$\mu(x_i, \beta)$に置き換える。条件付き期待値関数$E(Y|x)$は
　$\min_{\beta} \sum (y_i - \mu(x_i, \beta))^2$
の解だ。同様に、$\xi$をあるパラメトリックな関数$\xi(x_i, \beta)$に置き換える。条件付き分位点関数(CQF)とは
　$\min_{\beta} \sum \rho_\tau(y_i - \xi(x_i, \beta))$
の解だ。
　これ、線形計画法で簡単にとけるのです。

　よくある誤解は、分位点回帰なんかしなくても、反応変数で層別して、それぞれに最小二乗フィッティングすりゃいいじゃん、というもの。これは従属変数によるトランケーションであって、悲惨な結果を招く。いっぽう、共変量の側で条件づけるのはオッケー。それこそが局所フィッティング、すべてのノンパラ分位点回帰の基礎となるアイデアである。
　一番極端な場合として、共変量ベクトルでデータを$p$個のセルにわけ、それぞれのセルの中でふつうに分位点を計算しちゃう、という手もある。別の極端なアプローチとして、なんらかカットオフをきめて、反応変数がそれを超える確率について二項反応モデルを推定する、という手もある。
　分位点回帰推定量の漸近特性については研究がいっぱいある。推定手法はいろいろあるが、たいして変わらないことがわかってます。
　云々。

論文：データ解析(2015-) - 読了：Koenker & Hallock (2001) 分位点回帰入門

Koenker, R. (2015) Quantile Regression in R: A Vignette.
　分位点回帰のRパッケージ quantreg のvignette。実戦投入前のおまじないとして読んだ。著者についてはよくわからないが、調べてみたら2005年にCUPから"Quantile Regression"という単著を出しているから、有名な人なのかも。

　... 細かいところ、難しくってよくわかんなかった。ま、使いながら考えよう。

　著者曰く、分位点回帰というアイデアは18世紀、クロアチアのイエズス会巡回説教師(the itinerant Croatian Jesuit) Rudjer Boscovichに遡る由。はぁ？？まじ？！と思ったが、考えてみたら、回帰という名を生んだゴールトンは19世紀の人だけど、回帰という発想自体は18世紀からあったそうだし(詳細は全く知りませんが)、条件つき分布の中央値に注目する発想(分位点回帰)が、条件付き分布の期待値に注目する発想(ふつうの回帰)に負けないくらいに古くても、おかしくはないわね。Wikipediaで調べたら、ルッジエーロ・ジュゼッペ・ボスコヴィッチ、物理学・天文学・数学・統計学に名を残している人らしい。

論文：データ解析(2015-) - 読了：Koenker (2015) 分位点回帰パッケージ quantreg

2016年4月22日 (金)

Allenby, G.M., Rossi, P.E. (2008) Teaching Bayesian Statistics to Marketing and Business Students. The American Statistician, 62(3), 195-198.
　ベイジアン・モデリングの大スター Allenby兄貴がお送りする、「俺はビジネススクールの統計学のクラスでベイズ統計を教えてるぜ」エッセイ。

　(私が彼を兄貴と呼ぶ理由は、何度でも書きますが、主著"Bayesian Statistics and Marketing"においてAllenbyが読者に求める無理難題が、香港のカンフースター、ドニー・イェン兄貴が共演者に求めるという無理難題と似ているからである)

　以下、すごくいいかげんなメモ：

　なんといっても大変なのは、尤度ベースの推論の価値と、情報事前分布の意義について教えることですね。経済統計あがりの学生は事前分布に懐疑的だし、積率法と漸近分布の理論に慣れてるし。修士で統計やってましたっていう学生でも、ベイズ統計は練習問題くらいしかやってないし。
　とにかく基礎教材が大事だと思います。情報事前分布が、主観的な事前情報の組み込みという点だけではなく、スムーズネスや数値的安定性という点でも大事だと教え込むわけです。パネルデータ分析で、階層ベイズ離散選択モデルじゃないとできない分析を教えるとか。
　
　講義の具体的な内容ですか？
　教えているのはオハイオ州立大とシカゴ大のPh.Dコースの学生です。教科書はもちろん"Bayesian Statistics and Marketing"。この本にはRのbayesmパッケージがついてますからね、そりゃもう、うまくいってますよ。
　毎週の宿題では、いろんな問題について実装の練習をさせてます。基本的なコードを渡して拡張させるんです。回帰分析のコードを渡して、係数に順序制約をつけなさい、とか。スクラッチから書かせると、標準的な公式について直接にアホな実装してきますからね。クロネッカー積の関数を使うとか、逆行列とか行列式を求めるとか、そういうのはやらせたくないわけです。コレスキー分解とクロス積の関数で済ませたほうがいい。
　えーと、コースの前半はベイズ統計の基礎に割きます。標準的な内容だと思うけど、手法や実装の実用的な価値、それから情報分布に力を入れてます。後半はマーケ・ビジネス寄りです。潜在変数モデル、多項離散モデル、多変量離散モデル、離散選択モデルの需要理論[←?]、ユニット・レベルの尤度を使った離散変数の階層モデル、モデル選択と決定理論、同時性、ってとこですね。
　最後に、最近の研究から事例を選ばせて、ベイズ統計の観点からどうやって解けるか、というのをやります。ケーススタディを５つ、本に載せたから、読んでみてください。
　成績は宿題と最終レポートでつけます。宿題の採点はなるべく優しく... 単位がもらえるかどうかみんなおびえてますからね。レポートのほうは、学術論文と似た形で書かせます。テーマは自由だけど、標準的な回帰モデルはダメ。自分がやってる研究と関連する内容について取り組むのが理想ですね。実データの分析だけじゃなくてシミュレーションもさせます。目標は、自分のデータの尤度関数について真剣に考えさせることです。
　
　ミクロレベルのデータを扱っている学生にとって、消費者異質性のベイズ的な取り扱いはとても重要です。異質性の一般的なモデルを一連の条件つき分布で構築できます。
　$y_i | \theta_i \ \ \sim \pi_1 (\theta_i)$
　$\theta_i | \tau \ \ \sim \pi_2 (\tau)$
　$\tau | h \ \ \sim \pi_3(h)$
最初の条件つき分布を僕らはユニットレベル尤度と呼んでます。マーケティングだと$y_i$はたいていなんらかの形で離散的です。2つめの条件付き分布は、パラメータについての事前情報の一部になってます。経済統計ではランダム係数の分布って呼びますが、僕らはこれも尤度の一部だと捉えます。たいていの場合、一番関心が持たれるのは3番目の$\tau$です。学生に教えるポイントは、階層モデルと事前分布の価値、そして、すべての未知の量をランダムだと捉える推論手法の価値です。
　非正則事前分布を使う場合、階層モデルでは事後分布が必ずしも存在しません。階層モデルというアイデアは正則事前分布においてのみ意味を持ちます。だから学生は事前分布の重要性について真剣に考えないといけないわけです。
　階層モデルはモジュラー的で、モデリングの面でもアルゴリズムの面でも実装が楽だ、という点も、ベイズ統計に懐疑的な学生にとっていい勉強になります。ある事前分布について標準正規分布を混合分布に切り替える、ディリクレ分布に切り替える、なんてことが容易にできます。
　同時性についての議論も学生の関心を惹くようです。ビジネスデータはたいてい、消費者の企業の相互作用の結果で、どちらも同時に自分の行動を最適化しようとしてますから。需要-供給の同時モデルは尤度の最大化が難しくなるので、シミュレーション手法がぴったりです。よくベイズ統計の教科書には回帰モデルが載ってますが、そういうのではベイズ統計の力がわかりにくいんんじゃないですかね。

　[ここでbayesmパッケージの紹介。分布の可視化も大事なのよ、bayesmパッケージにはそういう機能もあるよ、とのこと。パス]

　というわけで、優秀な学生ならば、およそどんなモデルであっても推定できる概念的・操作的なスキルが身に付きます。そうでもない学生であっても、自分がやってる分析の背後にある想定について評価する能力は身に付きます。
　悩みとしては... 学生からの評判はいいんですけど、宿題が多すぎるっていう不満の声もありますね。それから事後分布について、可視化して全体を検討させたいんですが、どうしても積率での要約に頼ってしまい、なかなかうまくいかないです。決定理論に対するベイズ統計の価値についても、説明の時間がなかなか取れないです。
　コース終了後、勤勉な学生がさらに勉強を進める際にもbayesmパッケージは役立ちます。修了生のなかにはベイズ統計を使った応用論文を書いた人もいますよ。

　。。。Allenby兄貴。。。 学生に対しても、ドニー兄貴なみの厳しさであったか。。。

　せっかくなので(なにがだ)、ここでドニー兄貴のアクション・シーン。「孫文の義士団」(2009, 陳徳森監督)より。いったん撮り終えたアクションが気に入らず監督に直訴、自分のチームのメンバーを急遽招集、自らアクション監督をつとめて自腹で撮り直した、といういわくつきの名場面である。

　雑踏を駆け抜けるドニー兄貴も素晴らしいが、敵役のカン・リーが登場する直前、モノやらヒトやらがなぜか上方向に吹っ飛んでいくというカットもとても面白いと思うのです。カン・リーの並外れた怒りと力を表現しているんだけど、アニメーションのようなコミカルさがある。

論文：データ解析(2015-) - 読了：Allenby & Rossi (2008) 学生にベイズ統計をこうやって教えてます (feat. 「孫文の義士団」)

2016年4月21日 (木)

Srinivasan, V., Mason, C.H. (1986) Nonlinear least squares estimation of new product diffusion models. Marketing Science, 5(2), 169-178.
　Bassモデルの推定方法として、Bass(1969)によるOLS推定、Schmittlein & Mahajan (1982, Mktg.Sci.)によるML推定、そしてNLS推定を、実データやシミュレーションで比較。NLS推定がもっとも良いです。という論文。ざっと目を通しただけだけど、いちおう読了にしておく。

　OLSがひどいのはまあわかるとして、ML推定でパラメータの標準誤差が過小評価されるのはなぜか。そもそも推定誤差には、(1)標本抽出に伴う誤差、(2)経済条件やマーケティング変数みたいな変数を無視していることの誤差、(3)密度関数の誤指定による誤差(たとえば模倣係数が累積受容率に応じて動く場合とか)、の3つのソースがあるんだけど、ML推定では(1)しか考慮してない、だからパラメータの標準誤差が過小評価される。でもNLS推定なら全誤差の項をいれてある大丈夫。なのだそうだ。
　うーん、その理屈がいまいち呑み込めない。ちゃんと読んでないせいだろうか。それともSchmittlein & Mahajanを読まないといかんのだろうか。

論文：データ解析(2015-) - 読了：Srinivasan & Mason (1986) Bassモデルは非線形最小二乗法で推定しろ

2016年4月20日 (水)

　新製品・サービスの普及プロセスを説明する古典的モデルである、Bassモデルの推定方法について、ちょっと混乱しちゃったので、覚え書きをつくっておく。個人的メモでごさいます。

　そもそもBassさんが考えた理屈はこうだ。原典であるBass(1969)に沿ってメモする。
　ある製品の初回購入だけに注目する。総購入個数を$m$としよう。
　連続時間上の時点 $t$ における「購入の尤度」を $f(t)$ とする[わかりにくい表現だけど、ある潜在的購入者に注目したときにこの瞬間に購入が発生する確率密度、という意味合いであろう]。
　$F(t)$ を以下のように定義する。
　$F(t) = \int f(t) dt, \ \ F(0) = 0$
　時点 $t$ における売上を $S(t)$、累積売上を $Y(t)$ とする。

　時点$t$において、まだ購入が発生していないという条件の下での「購入の尤度」を$P(t)$とすると、
　$P(t)=f(t)/(1-F(t))$
である。また、
　$S(t) = m f(t) = P(t)(m-Y(t))$
　$Y(t) = \int S(t) dt = m \int f(t) dt = m F(t)$
である。

　さて、Bassモデルとは次のモデルである。
　$P(t) = p + (q/m) Y(t)$
　ここから以下が導かれる。
　$f(t) = p + (q-p) F(t) - q F(t)^2$
　$S(t) = pm + (q-p) Y(t) - (q/m) Y(t)^2$

　さらに、$f(t), F(t)$を$t$について解くと、
　$f(t) = \frac{(p+q)^2}{p} \frac{\exp(-(p+q)t)}{(q/p \exp(-(p+q)t)+1)^2}$
　$F(t) = \frac{1-\exp(-(p+q)t)}{q/p \exp(-(p+q)t)+1}$

　ここまではいいっすね。ここからはBass's Basement Research Instituteによる説明のメモ。

　実際には$t$は離散変数だ。すると、次のような困った事態が発生する。
　区間$[t-1, t]$について考えよう。$F(t)$は、この区間の右端までの$f(t)$の合計？ では$f(t)$は？ 時点$t$における瞬間の速度？ それとも、$F(t-1)$から$F(t)$までの増分？
　アルゴリズムとして実装するには、どうにかつじつまを合わせないといけない。ここで3つのアプローチが登場する。

　第一のアプローチ。Sinivasan & Mason (1986)のアプローチである。$F(t)$のほうはさきほどの解をそのまま用いる。すなわち
　$F(t) = \frac{1-\exp(-(p+q)t)}{q/p \exp(-(p+q)t)+1}$
いっぽう$f(t)$のほうを変えてつじつまを合わせる。
　$f(t) = F(t) - F(t-1)$
　ただし、$t=1$においては$f(t) = F(t)$。
　このやり方の良い点：$F(t)$の定義はそのまま。$f(t)$は小さな値なので、これを直接求めずに引き算で求めるのは、数値的に正確。

　第二のアプローチ。$f(t)$のほうの解をそのまま用いる。すなわち
　$f(t) = \frac{(p+q)^2}{p} \frac{\exp(-(p+q)t)}{(q/p \exp(-(p+q)t)+1)^2}$
いっぽう$F(t)$のほうを変えてつじつまを合わせる。
　$F(t) = \sum_{i=0}^{t} f(i)$
これはおかしい。$f(t)$はその瞬間の速度であって、ある期間の平均ではないからだ。

　第三のアプローチ。Bass先生が1969年の段階で推していた方法。$f(t), F(t)$の差分方程式に戻る。
　$f(t) = p + (q-p) F(t-1) - q F(t-1)^2$
　$F(t) = F(t-1) + f(t)$
これはおかしい、すごくおかしい。Bass先生はなぜこんな方法を推したのか？ コンピュータじゃなくて手回し計算機を使っていたからに他ならない。なんとかしてOLS回帰の形に持ち込む必要があったのだ。この$f(t)$の式なら、$f(t)$を当期購入数、$F(t-1)$を前期までの累積購入数として、重回帰で解けるからね。
　教訓。Bass(1969)の推定方法を用いるなかれ。コンピュータを使いなさい。第一のアプローチを使って、非線形回帰で推定しなさい。

　... なるほどね。
　最後に、手元にあるいくつかの教科書を比較しておこう。

　日頃より頼りにしております、朝野・山中(2000), p.152。
　$m$をN, $F(t)$をn_T, $p$をp, $q$をr, $Y(t-1)$をZ_t (添え字に注意!), $S(t)$をX_tと呼んでいる。
　本文中でZ_tを「t期における購入者」、X_tを「t期の購入者」と呼んでいるので、前後の文脈を読み取ってその意味をよく考えないといけない。
　推定方法はOLS。推定する式はこうなる(8.5式)：
　$X_t= pN + (r-p) Z_t - \frac{r}{N} Z_t^2$

　照井・ダハナ・伴(2009), pp.155-157。この本も頼りになります。
　$m$をN, $F(t)$をF(t), $p$をp, $q$をq, $Y(t)$をN(t), $S(t)$をS(t), と呼んでいる。
　推定方法はOLS。本文中の記号とRのコード例での変数名が違うので要注意。コード例ではなぜか$N(t-1)$にFt.1, $n(t)$にNtという名前をつけ、lm(Nt ~ Ft.1 + Ft.1^2)を推定している。(どうもよくわからないんだけど、formulaの指定はこれでいいのだろうか？ Nt ~ Ft.1 + I(Ft.1^2) ならわかるんだけど...)

　里村(2010)、pp.64-67。
　$m$をm, $F(t)$をF(t)、$f(t)$をf(t)、$Y(t)$をN(t)、$S(t)$をX(t)と呼んでいる。$m-Y(t)$をY(t)と呼んでいるので注意。
　この本は、$F(t)$を$t$の関数として解いた式を非線形回帰で推定している。推定にはRのnls()を使っている。シリーズの主旨に沿って、Rのコード例がとても充実している。よく見ると日本語の誤字が結構多いのだが(「統計的に有為」とか)、そういうことをいちいち気にしていてはいかんのである。

　Lilien & Rangaswamy (2004), pp.255-258。ちょっと古い版の教科書だが、もったいないので使い続けている。前職のときに私費で買ったのだ、高かったのだ。
　$m$を$\bar{N}$、$F(t)$をF(t), $f(t)$をf(t)、$p$をp, $q$をq, $Y(t)$をN(t), $S(t)$をs(t), と呼んでいる。推定方法として、OLSと非線形回帰の両方を紹介している。
　分量も違うし、そもそも日本語の教科書と英語の教科書を比較してはいけないのだが(なにしろ市場規模が違う)、とにかくこの本の説明は最高にわかりやすい。もちろんBass(1969)なんかよりはるかにわかりやすい。

雑記：データ解析 - メモ：Bassモデルをどうやって推定するか (いろんな教科書を比べてみた)

2016年4月18日 (月)

Bass, F.M. (1969) A new product growth for model consumer durables. Management Science, 15(5), 215-227.
イノベーション普及のBassモデルを最初に提案した論文。原稿で引用したいのだが、いくら超有名な論文とはいえ、読んでもいないのを引用するのは気が引けるので、いちおうざっと目を通した次第。

面白かった点をメモ：イントロに「この理論は数学的には感染モデルに由来している。感染モデルは疫学において幅広く用いられてきている」とある。ああ、そりゃそうか、疫学のモデルのほうがずっと古いだろうな(WikipediaでSIRモデルの項目を見たら、原型のKermack & McKendrickはなんと1927年に遡る)。Bassモデルに限らず、初期のマーケティング・サイエンスってのは他の分野に学ぶところが多かったのだろう。

論文：マーケティング - 読了：Bass(1969) 諸君、俺の名をとってこれをBassモデルと呼び未来永劫語り継ぐがよい (とは書いてないけど)

Vanhonacker, W.R., Lehmann, D.R., Sultan, F. (1990) Combining Related and Sparse Data in Linear Regression Model. Journal of Business and Economic Statistics. 8(3), 327-335.
　Sultan, et al.(1990)の後続研究としてSultanさんが紹介していた奴。実はこっちのほうがわかりやすいんじゃないかと疑って読んでみた。その予想はあたっていて、細かいところまできちんと説明してくれている論文であったが(というか、Sultan, et al.(1990)が説明を端折りすぎなのである)、しっかし、ねえ...

　イントロ部分に早速わからない箇所があるので訳出しておく。

事前情報を使うことによって得られる効率性は、ベイズ統計手法が幅広く用いられ発展してきた理由でもある。Lindley & Smith (1972)は階層事前分布の使用を奨めており、多くの適用例の基礎となっている。独立に生成された知見を事前分布に統合しようという初期の試みはこのアプローチを採用している。しかし、階層事前分布の根底には交換可能性という非常に強い想定があり、これは過度な単純化だとして批判されてきた(Kass, 1983 JASA)。この想定は、通常は正規事前分布として捉えられているのだが(Leamer, 1978書籍)、知見はお互いに関連しており、データはその関連性についての情報を含んでいる、ということを述べている。この関連性には、研究における差異のせいで、ほぼ確実に異質性が含まれている。異質性についての事前情報はスパースであることが多く、グルーピングや条件付き事前分布の使用は現実的でないので、事前の研究のすべてに等しい情報が含まれているという厳密な想定を伴う伝統的アプローチが、妥当な批判にも関わらず生き延びている。
　体系的変動をあきらかにするという試みにおいて、メタ分析の結果は独立に実行された研究群を通じたパラメータ異質性を明示的に捉える。従って、その結果を事前情報として用いれば、階層ベイズ事前分布の根底にある、厳密な交換可能性という想定を回避し、こうした事前分布に頼ることに対して実証研究の中でしばしば表明されてきた懸念に応えることができる。新しい研究を、メタ分析モデルで指定されたモデレータ変数について先行研究とマッチングさせることによって、事前分布の情報内容を差異化できる。いいかえれば、異質性が事前情報において明示的に捉えられていることになる。

ううむ。ひらたくいえばLenk&Raoのような階層ベイズモデルを使うアプローチをディスっているのだと思うのだが、階層事前分布における交換可能性の想定はきつすぎる、という批判の意味がよくわからない。もしかすると、階層事前分布には共変量を入れられないけどメタ分析モデルなら入れられるよ、という意味なのだろうか。だとしたら、Lenk&Raoが提案しているように、共変量をいれた階層ベイズモデルならいいんでしょ、という話になるはずで、交換可能性という概念を持ち出す必要はないだろう...

　まあいいや。以下、本文についてのメモ。

　いま、オリジナルのm個の研究で重回帰が使われているとする。
　$y_i = X_i \beta_i + u_i$
$i$は研究番号。$y_i, X_i$のデータサイズを$n$、$X_i$の変数数を$p$とする[$n$と$p$に添字$i$がついてないけど、ほんとはつけるべきだろう]。すべての$i$において$E(u_i) = 0$, 異なる2つの研究$i,j$の間で $E(u_i u_j')=\sigma^2_i I$とする。

　ここから「メタ分析事前分布」を導出する。
　いま、研究の属性(製品カテゴリとか推定手法とか)によって、研究を$k$個のグループに分類できるとする。研究$j$がグループ$l$に属しているとき、ランダム係数を想定して
　$\beta_j = \bar{\beta}_l + e_j$
ただし$E(e_j) = 0$、グループ$l$に属する異なる2つの研究$i,j$の間で$E(e_i e_j') = \Sigma_l$とする。

階層ベイズ事前分布は、すべての$i$について、パラメータ・ベクトル$\beta_i$がある平均ベクトルと共分散行列を持つある事前分布からのドローだと想定する。従って$\beta_i$は交換可能であり、この事前分布はすべての研究のあいだの「つながりの程度」を捉えていることになる。[中略]この強い想定について検証する方法もある。この想定によって、ある研究グループを分離して分析できるようになる。[中略。いっぽう]メタ分析によって異質性を明示的に捉えて記述すれば、すべての研究グループを同時に分析できるようになる。こちらのほうが効率が良いし、理解も促進される。
[やっぱりそうだ。階層ベイズモデルには共変量がいれられないというのが話の前提になっているみたいだ]

　まとめて表現しよう。$m$個の研究について、パラメータ、平均、誤差をそれぞれ縦積みする。
　$\beta' = (\beta'_1, \ldots, \beta'_m)$
　$\bar{\beta}' = (\bar{\beta}'_1, \ldots, \bar{\beta}'_k)$
　$e' = (e'_1, \ldots, e'_m)$
とする。2本目の式に注意。たとえば3つ研究があって、所属グループが$(1,1,2)$なら、$\bar{\beta}' = (\bar{\beta}'_1, \bar{\beta}'_1, \bar{\beta}'_2)$である。
　以上を使って
　$\beta = \bar{\beta} + e$
さて、$\bar\beta$を各研究のデザインを列とする$Z$の線形和とみて
　$\beta = Z c + e$
とする。これがメタ分析モデルだ。

　実際には真のパラメータベクトル$\beta$は観察されないわけで、そのかわりに
　$\hat\beta = Z c + \epsilon$
を推定する。 誤差項$\epsilon$について考えよう。これは真のパラメータに含まれる誤差項$e$と、パラメータの推定誤差$\hat{u} = \hat\beta - \beta$の和である。$E(\epsilon)=0$となる。では$E(\epsilon \epsilon’)$はどうなるか。
　[さあ、この辺から話がややこしくなってくるんだけど]
　$E(\epsilon \epsilon’)$を$\Omega$としよう。それは$m \times m$のブロック対角行列で、グループ$l$の対角ブロックが$\Sigma_l + \sigma^2_i (X'_i X_i)^{-1}$となる。[なぜか理解していないが、信じることにしよう]
　で、ここから、$\hat{c}$のBLUEな推定量は下式となるのだそうだ[これも信じることにしよう]。
　$\hat{c} = D\hat\beta$
　$D = (Z' \Omega^{-1} Z)^{-1} Z' \Omega^{-1}$

　さて、新しい研究がやってまいりました。そこでは
　$y_o = X_o \beta_o + u_o, \ \ E(u_o)=0, \ \ E(u_o u'_o)=\sigma^2_o I$
が得られているとしましょう。この研究のデザインは$Z_o$で、これをメタ分析モデルに放り込み、 $\beta_o$の事前推定値
　$b_o = Z_o \hat{c}$
が得られたとしましょう。
　この$b_o$は、こういう性質を持っている、ということが示せるのだそうだ。
　$b_o - Z_o D d = Z_o D M \beta_o + Z_o D \gamma$
$D$はさっき定義した行列。$d$は、すべての研究について、それが属するグループの平均パラメータのベクトルと、新研究が属するグループ平均パラメータのベクトルの差を縦積みしたもの。$\gamma$は$\hat{u}$と、「それが属するグループの$e$と新研究が属するグループの$e$の差」を足しあげたもので... ま、これは実際には使わないので、あんまし気にしなくてよさそうだ。

　ってことはですね。
　$y_o = X_o \beta_o + u_o$
　$b_o - Z_o D d = Z_o D M \beta_o + Z_o D \gamma$
この2つを縦に積み上げて、一気に推定すればいいや、ということになる。縦に積んだのをこうやって表現することにします。
　$y = X \beta_o + v$
ただし、$y$は$y_o$と$b_o - Z_o D d$を縦に積んだもの。$X$は$X_o$と$Z_o D M$を縦に積んだもの[本文中には $X' = \left[X'_o (Z_o D M)'\right]$とあるが、$X' = \left[X'_o, (Z_o D M)'\right]$の誤植であろう]。$v$は$u_o$と$Z_o D \gamma$なる謎のベクトルを縦に積んだもの。
　さあ、この$\beta_o$の推定量はどうなるか。こうなることが示せるのだそうだ。
　まず、誤差項$v$の共分散行列$\Omega_v$について考えよう。それはブロック対角行列となり、左上には$E(u_o u'_o)=\sigma^2_o I$が、右下には$Z_o D E(\gamma \gamma') D' Z_o$がくる。
　では、$E(\gamma \gamma')$はどうなるのか。それは次の2つの行列の和になるのだそうだ。ひとつめ、$\Omega$。ふたつめ、新しい研究がグループ1に属するとして、すべての要素が$\Sigma_1$である行列。
　さて、$\beta_o$のBLUEな推定量は
　$\hat\beta = (X'\Omega_v^{-1}X)^{-1}X'\Omega_v^{-1}y$
これをrecursive estimatorという。

　実際の推定は次のように行う。[一応メモするけど、要するにtake it easyという感じ。原文から少し添え字を変えた]。
　グループ$l$の研究数を$m_l$とする。新研究はグループ1に属するということにする。

　事例が2つ。1つめのマーケティング・ミクス・モデルはとばして、二つ目の新製品普及のほうのみ読んだ。
　ご存知Bassモデルは、時点$t$において未受容者が受容する条件付き確率を$P(t)$、既受容者を$Y(t)$として
　$P(t) = p + (q/m) Y(t)$
離散データの場合はラグをとって
　$S_t = a + b Y_{t-1} + c Y^2_{t-1}$
とすればOLS推定できる。ほんとはNLSのほうがいいんだけど、話の都合上OLSで推定する。
　過去の19個の普及曲線について$a, b, c$を推定した結果を集めてきた。これを消費財と産業財の2群にわける。
　ほんというと、過去の普及曲線について$a, b, c$の共分散がほしいんだけど、手に入らないのであきらめて、共分散はみんな0だということにした[$\omega$を単なる対角行列にし、対角要素は当該パラメータの分散にした、ということであろうか]。
　誤差分散の推定値$\hat\sigma^2_i$も手に入らないので、$\hat\sigma^2_o$のかわりに一年目の売上をいれた[いいのか、そんなことで。というか1年間の売上ってなんのことだ、予測対象であるカラーテレビのことか]

　新研究としてカラーテレビの普及曲線を使い、1期目のみ、2期目まで、の2つについて事後パラメータを出す。2期でかなり予測できるようになった。3期だともうぴったり予測できちゃう。云々。
　
　いやあ... メモしてみたはいいけれど... これ、死ぬほどややこしいわ... そのわりに、線形回帰に限定されている、という... Bassモデルの推定だったら、やっぱし階層ベイズモデルをつかったほうがいいんじゃなかろうか。

　話を逆向きにたどり、推定手順をメモしておこう。

  1. それぞれの研究について、誤差分散の推定値を求める。すなわち$\hat\sigma^2_i = \frac{(y_i - X_i \hat\beta_i)'(y_i - X_i \hat\beta_i)}{n-1}$。
  2. 新研究について、誤差分散の推定値を求める。それがグループ1に属するとして、$\hat\sigma^2_o = \frac{1}{m_1} \sum_{i=1}{m_1} \hat\sigma^2_i$ 。
  3. それぞれのグループについて、パラメータ平均の推定値を求める。すなわち$\hat{\bar{\beta}}_l = \frac{1}{m_l} \sum_{j=1}{m_l} \hat{\beta}_j$。
  4. それぞれのグループについて、パラメータの分散を求める。すなわち$\hat{\Sigma}_l = \frac{1}{m_l} \sum_{j=1}{m_l} (\hat{\beta}_j - \hat{\bar{\beta}}_l)(\hat{\beta}_j - \hat{\bar{\beta}}_l)'$。
  5. すべての研究について、「そのグループの平均パラメータ」と「グループ1の平均パラメータ」の差を求める。これが$d$。
  6. 過去の研究のパラメータ推定値$\hat\beta$について、その誤差項$\epsilon$の分散行列$\Omega$を求める。それは$m \times m$のブロック対角行列で、グループ$l$の対角ブロックが$\Sigma_l + \sigma^2_i (X'_i X_i)^{-1}$。
  7. 話の都合上登場した謎の変数$\gamma$の共分散行列$E(\gamma \gamma')$を求める。それは、新研究がグループ1に属するなら、$\Omega$の全要素に$\Sigma_1$を足したもの。
  8. $\hat\beta$から事前パラメータを求める際に使う係数行列$D=(Z' \Omega^{-1} Z)^{-1} Z' \Omega^{-1}$を求める。
  9. 事前パラメータ$b_o=Z_o D\hat\beta$。
  10. 縦積みモデルの従属変数ベクトル$y$をつくる。それは新研究の従属変数$y_o$と、$b_o - Z_o D d$を縦積みしたベクトル。
  11. 縦積みモデルの誤差項の共分散行列$\Omega_v$を求める。それはブロック対角行列で、左上に$\sigma^2_o I$が、右下に$Z_o D E(\gamma \gamma') D' Z_o$がくる。
  12. 事後パラメータの推定量は、$(X'\Omega_v^{-1}X)^{-1}X'\Omega_v^{-1}y$。

つ、疲れる...

論文：データ解析(2015-) - 読了：Vanhonacker, Lehmann, Sultan (1990) 重回帰のメタ分析で事前分布をつくってベイズ推定

2016年4月15日 (金)

Allenby, G.M. & Rossi (2006) Hierarchical Bayes Models. in Grover, R. & Vriens (eds), The Handbook of Marketing Research: Uses, Misuses, and Future Advances, Sage.
　仕事の都合で読んだ。ベイジアン・モデリングの大スター、Allenby兄貴による、HBモデルの素人向け解説。

　私が心のなかで彼を兄貴と呼んでいるのは、香港映画界最後のカンフー・スター、ドニー・イェン兄貴と似ていると思うからである。聞くところによれば、ドニー兄貴は宇宙最強なアクションを追求するあまり、共演者にとんでもない無理難題を求めるのだそうだ。Allenby兄貴も、著書に"Bayesian Statistics and Marketing"なんていかにも親しみやすそうなタイトルをつけておいて、読者をあの無闇に難解な内容へと突き落すあたり、ドニー兄貴に瓜二つである。この類似性、お分かり頂けただろうか (←誰に言っているのかわからないが)
　ま、さすがに市場調査実務家向けハンドブックの1章であれば、そんなには難しくないだろうと思って手に取った次第。

　いわく。
　過去10年、消費者行動のモデルは飛躍的進歩を遂げた。個人レベル・データの分析においてはベイジアンの手法がすっかり普及し、選択モデル(コンジョイント分析)だけでなく幅広い問題に用いられるようになった。消費者間異質性に関心が集まり、ベイジアン手法によるモデリングが進んだ。本章では、これらの進展を支えた立役者、階層ベイズモデルについてご紹介しよう。

1.想定
　たとえば、ある製品の需要に対する価格の影響について分析したいとしよう。販売量が連続的単位で得られるなら、価格感受性を回帰モデルで測れる：
　$y_t = \beta_0 + \beta_1 price_t + \epsilon_t$
　$\epsilon_t \sim Normal(0, \sigma^2)$
って感じね。
　ところが個人レベルの需要は、こんな風に連続的には動かない。そこで、連続モデルと打ち切りを想定し、
　$y_t = 1$ if $\beta_0 + \beta_1 price_t + \epsilon_t > 0;$ otherwise $y_t = 0 $
という風に考える。こういう発想は、散布図やクロス集計表からは得られない。ところがマーケティング実務家は、これまでこういうモデルをあまり使ってこなかった(コンジョイント分析を別にして)。

1.1 階層モデル
　さっきの式を書き換えようう。潜在変数を導入して、
　$y_t = 1$ if $z_t > 0$; otherwise $y_t=0$
　$z_t = \beta_0 + \beta_1 price_t + \epsilon_t$
　$\epsilon_t \sim Normal(0, \sigma^2)$
潜在変数$z_t$の下で、$\beta_0, \beta_1, \sigma^2$についての推論は$y_t$と独立である点に注目。
　モデルは階層的に捉えると便利なことが多い。上の例で、一本目の式は、もし販売量($z_t$)が十分に多いならば購買が生じやすい($y_t=1$)、というシナリオを表している。二本目の式は販売量と価格の関係を表している。こういう風に捉えるとモデルを精緻化しやすい。
　マーケティングでは、階層モデルは(a)個々の対象者の行動と(b)対象者を通じた反応の分布を記述するために用いられてきた。前者は単位内の行動のモデル、後者はモデルのパラメータの横断的な変動のモデルである(異質性の分布と呼ぶことが多い)。

1.2 ベイジアン分析
　階層ベイズモデルとは階層モデルのベイジアン手法による分析だ。ベイジアン手法とは[...ここでベイズの定理についての説明をひとしきり。事前オッズ、尤度比、事後オッズを使って説明している。これ、かえってややこしくないか]
　...とまあこのように、ベイジアン手法はエレガントなのだが、その割には最近まであまり使われていなかった。ちょっとしたモデルについて適用しても、計算がややこしくなっちゃうからである。

2.HB革命
　革命前、ベイズ定理を実装するためには、事前分布の確率密度に尤度をかける必要があった。たとえば、
　$y_t = \beta + x_t + \epsilon_t$
　$\epsilon_t \sim Normal(0, \sigma^2)$
というモデルがあるとして、$\beta$の事前分布が正規分布だとすると
　$\pi(\beta) = \frac{1}{\sqrt{2 \pi s^2}} \exp \left(\frac{-1}{2s^2}(\beta - \bar{\beta})^2 \right)$
省略するけど、$\pi(\sigma^2)$もなんらか考えないといけない。で、尤度は
　$\pi(y_t | \beta, \sigma^2, x_t)$
　$= \prod_{t=1}^{T} \frac{1}{\sqrt{2 \pi s^2}} \exp \left(\frac{-1}{2s^2}(y_t - \beta x_t)^2 \right)$
これらを全部かけて
　$\pi(\beta, \sigma^2 | y_t, x_t) \propto \pi(y_t | \beta, \sigma^2, x_t) \times \pi(\beta) \times \pi(\sigma^2)$
いやーエレガントですね。でも計算が難しすぎる。

　ここで革命が起きた。MCMCの登場である。[...ここでMCMCについてかなり投げやりな説明があって...] おかげで計算が容易になり、応用が一気に広がった。

3.事例
ここでは、多変量分布のextremes(分布の裾)について調べる例を示し、極値HB+MCMCの威力を示そう。

3.1 背景
　企業は、新製品をもっとも買いそうな顧客、マーケティング戦略の変化によってもっともスイッチしそうな顧客、もっとも良い反応をしてくれる人... を分析する必要がある。つまりextremesの理解が求められる。

3.2 モデル
　階層ベイズ・ランダム効果ロジットモデルについて考えてみよう。対象者$h$が選択肢$i$(その属性ベクトルが$x_i$)を選択する確率についてモデル化する。
　$Pr(i)_h = \frac{\exp(x'_i \beta_h)}{\sum_j \exp(x'_j \beta_h)}$
異質性をいれる。観察可能な共変量を$z_h$として
　$\beta_h = \Gamma z_h + \xi_h$
　$\xi_h \sim MVN(0, V_\beta)$
$\Gamma$は「部分効用が他人と違う人」がどんな人なのかを表している。$V_\beta$の対角成分の大きさは、部分効用の異質性のうち$\Gamma z_h$で捉えきれていない部分の大きさを表し、非対角成分は属性水準の評価のパターンを表す(たとえば、属性水準の部分効用の間に正の共変動があるということは、それらの水準を持っている製品が特定の人にすごく好かれるということだ)。
　$\Gamma$と$V_\beta$に事前分布を与え、このモデルを階層として書き下ろすと、
　$y | x, \beta$
　$\beta | z, \Gamma, V_p$
　$\Gamma | a, A$
　$V_\beta | w, W$
パラメータの数はすごく大きいが(個人ごとの$\beta_h$があるから)、MCMCなら大丈夫。

3.3 データ
　[クレジット・カードの選好について電話調査でコンジョイント実験をやった。属性は7つ、共変量は年齢、年収、性別。詳細略]

3.4 結果
　[推定しましたって話があって...]
　部分効用の共分散行列$V_\beta$の事後平均をみると、水準"out-of-state bank"と"low annual fee"の相関が高い。このことは、年会費の安いクレジットカードを提供することが、国外顧客の利用促進のために効果的であることを示唆している。[←このあとの説明で使うための指摘である。面白い指摘だけど、こんな風に共分散行列そのものをちくちくと眺めるのは、実際には難しかろう]

3.5 異質性の分布
　各属性水準についての部分効用の確率分布を、全体についても、ある特定の個人についても観察できる[実例]。各個人の部分効用推定がかなりの不確実性を伴っていることに注意。

3.6 Focusing on Extremes
　[ここ、ちょっと面白いので詳しくメモ]
　$\Gamma$の事後平均をみるとわかるのは、国外銀行のクレジットカードの効用は平均して3.758低いということだ。このペナルティを乗り越えるためのインセンティブとして、利率を下げる、年会費を下げる、の2つがある。
　$\Gamma$の事後平均だけをみると、どちらもインセンティブでも十分であるようにみえる。しかし、共分散行列$V_\beta$の事後平均をみると、年会費と国外銀行の相関が高い。つまり、年会費に対する感受性が高い人は国外銀行への部分効用が高めだ。いっぽう利率と国外銀行の相関は低い。
　低年会費+国外銀行と低利率+国外銀行の2選択肢を考え(他の属性は等しくする)、全体効用の分布を比較してみると、前者の選択肢では全体効用のばらつきが大きくなっており右裾だけみると前者のほうが大きいことがわかる。
　[なるほどね... 実際の場面ではめんどくさいのでいきなりシミュレータを使っちゃうところだと思うけど、こうやって全体効用の分布を細かく見ていくのも面白いなあ]

4. HBモデルを使うときに厄介なこと
　複雑なHBモデルを推定する際には、すぐに使える(off-the-shelf)ソフトがない。WinBugsは貴重な例外だがデータサイズが大きいときには問題が生じる。HBの普及に貢献したのはSawtooth Softwareのプログラムだが、選択型コンジョイントとOLS回帰のHB推定に限定されている。そこで！ 僕らが書いた"Bayesian Statistics and Marketing"を買うといいよ！Rのソフトつきだよ!!
　HB推定は古典的推定とは違って閉形式の解がない。いくらバーンインしても、パラメータ推定の平均の変動は残るかもしれない。また、個々人についての点推定が得られるのではなく、推定値の分布が得られる、という点に注意。そのせいで話が複雑になることがある(マーケット・シミュレーションとか)。

　やれやれ、終わった。疲れた。
　なにか変わった実例が紹介されるかと期待してたのだが、結局選択型コンジョイントだったよ...シクシク。
　この章、末尾に8ページにわたって主要文献解題が添えられている。労作だが、2006年の段階では、こんな風に「マーケティングにおけるHBモデルの主要論文」をリストアップすることができたのね。今ではちょっと無理だろう。

　ちなみにドニー兄貴による傑作アクション・シーンのひとつがこちら。「SPL/狼よ静かに死ね」(2005, 葉偉信監督)より。この振付を覚えろと敵役の俳優(呉京)に求めるあたり、Allenby兄貴が著書において説明もなくいきなり逆Wishart分布を持ち出す際の酷薄さを彷彿とさせますね。(すいません、いまちょっと適当なことをいいました)

　呉京さんは昨年公開された「SPL2」にも出演しているらしい。日本でも公開してくれるといいんだけどなあ。

論文：データ解析(2015-) - 読了：Allenby & Rossi (2006) 階層ベイズモデルとはなにか (feat.「SPL 狼よ静かに死ね」)

2016年4月13日 (水)

Bookcover 「聴くこと」の革命: ベートーヴェン時代の耳は「交響曲」をどう聴いたか (叢書ビブリオムジカ) [a]
マーク エヴァン ボンズ / アルテスパブリッシング / 2015-10-24
ベートーベン前後のドイツ社会における、器楽曲・交響曲を聴くということの意味の変遷を辿る。知識不足でよくわからない点も多かったが、面白い本であった。

Bookcover 中国民族主義の神話―人種・身体・ジェンダー [a]
坂元 ひろ子 / 岩波書店 / 2004-04-27

Bookcover 中世的世界とは何だろうか (朝日選書) [a]
網野 善彦 / 朝日新聞社 / 1996-06

Bookcover 本屋図鑑 [a]
本屋図鑑編集部 / 夏葉社 / 2013-08

Bookcover シャープ崩壊 ―名門企業を壊したのは誰か [a]
/ 日本経済新聞出版社 / 2016-02-18

Bookcover 憲法入門 [a]
長谷部 恭男 / 羽鳥書店 / 2010-01-08

Bookcover 武士という身分―城下町萩の大名家臣団 (歴史文化ライブラリー) [a]
森下 徹 / 吉川弘文館 / 2012-06

Bookcover 少年の名はジルベール [a]
竹宮 惠子 / 小学館 / 2016-01-27
著者生涯の傑作「風と木の詩」に至るまでの、若き日の思い出。萩尾望都に対する微妙な感情を赤裸々に綴るところ、いかにもこの作家らしい。竹宮恵子さんって、ふとした拍子に心の表面に浮かび上がる暗い感情を描かせたら、そりゃもう、右に出る人がいない。
　ラストの段落をメモしておこう。

　若いころの友人たちというものは、振り返ればあの一瞬にも思える時間のなかで、なぜ巡り合えたのだろうか。それ自体がこの世の奇跡だ。駆け足で通り抜けた季節の熱い風であり、咲き誇る花であり、光る刃だ。でも同じ場所に同じかたちであるものはない。かつて訪れたヨーロッパの様々な色彩の街と同じように。
　こんにちは、そして、さようなら。素晴らしい時間たち。

ノンフィクション(2011-) - 読了：「『聴くこと』の革命」「中国民族主義の神話 人種・身体・ジェンダー」「中世的世界とはなんだろうか」「本屋図鑑」「シャープ崩壊」「憲法入門」「武士という身分」「少年の名はジルベール」

Bookcover 終戦後史 1945-1955 (講談社選書メチエ) [a]
井上 寿一 / 講談社 / 2015-07-11

Bookcover 戦争と検閲――石川達三を読み直す (岩波新書) [a]
河原 理子 / 岩波書店 / 2015-06-27

日本近現代史 - 読了：「戦争と検閲 石川達三を読みなおす」「終戦後史1945-1955」

Bookcover 日本仏教入門 (角川選書) [a]
末木 文美士 / KADOKAWA/角川学芸出版 / 2014-03-21
入門というタイトルだけど、著者の紀要論文などを集めた内容。アタリハズレがありうる怖いタイプの本だが、恐る恐る読んでみたら、とても面白かった。
　最後のほうに出てきた話。
　僧侶が修行から離れて社会的活動を展開することは古くからあるけれど、日本の中世において、社会的活動にもっとも積極的だったのはどの宗派か。
　当然、鎌倉新仏教じゃないかと思うじゃないですか。浄土宗とかさ。しかし著者によれば、正解は南都六宗のひとつ、律宗であった。
　彼らの活動は、死者の葬送、非差別民の救済、港湾管理、橋・道路の建設などに及んだ。出家と在家を峻別し修行と戒律を重んじる宗派がなぜ？　その背景には、戒律は清浄であり、その力によって穢れを克服できる、という観念があったのだそうだ。だから僧侶たちは、非人と交わることもできるし、社会的なタブーを乗り越えることができた。
　ううむ、そうか... 散歩の途中で見かけるカトリックの修道院のことなどを思い出し、考え込んでしまった。ここには、社会から離れることによって社会により深くコミットできることがある、というアイロニーがあるのかもしれない。

Bookcover 草木成仏の思想 [a]
末木文美士 / サンガ / 2015-02-25
味をしめて読んだ本。

Bookcover プラトンとの哲学――対話篇をよむ (岩波新書) [a]
納富 信留 / 岩波書店 / 2015-07-23

Bookcover ヘーゲルとその時代 (岩波新書) [a]
権左 武志 / 岩波書店 / 2013-11-21

Bookcover こころはどう捉えられてきたか: 江戸思想史散策 (平凡社新書) [a]
田尻 祐一郎 / 平凡社 / 2016-03-17
日本近世思想についてのエッセイ風の読み物。面白かった。

哲学・思想(2011-) - 読了：「こころはどう捉えられてきたか」「日本仏教入門」「草木成仏の思想」「プラトンとの哲学」「ヘーゲルとその時代」

Bookcover 銃座のウルナ 1 (ビームコミックス) [a]
伊図透 / KADOKAWA/エンターブレイン / 2016-02-25
この著者の長編は、とにかく開始の時点では抗しがたい魅力を湛えているのである。「ミツバチのキス」しかり「エイス」しかり。この作品でいうと、戦場になぜかジャンプ台があり、敵の妖怪変化がスキーのジャンプ競技よろしく飛んでくる、というとんでもない発想に感嘆した。ちゃんと風呂敷を畳んでくれるといいんだけど...

Bookcover 闇金ウシジマくん 36 (ビッグコミックス) [a]
真鍋 昌平 / 小学館 / 2016-02-29

Bookcover 生き返っても、あの世 [a]
村上 竹尾 / 幻冬舎 / 2016-02-25

Bookcover ディザインズ(1) (アフタヌーンKC) [a]
五十嵐 大介 / 講談社 / 2016-02-23

Bookcover 弟の夫(2) (アクションコミックス(月刊アクション)) [a]
田亀 源五郎 / 双葉社 / 2016-01-12

Bookcover あさひなぐ 18 (ビッグコミックス) [a]
こざき 亜衣 / 小学館 / 2016-02-29

Bookcover 富士山さんは思春期(6) (アクションコミックス) [a]
オジロ マコト / 双葉社 / 2015-04-28
Bookcover 富士山さんは思春期(8) (アクションコミックス) [a]
オジロ マコト / 双葉社 / 2016-02-27

コミックス(2015-) - 読了：「銃座のウルナ」「あさひなぐ」「ディザインズ」「闇金ウシジマくん」「生き返っても、あの世」「富士山さんは思春期」

2016年4月12日 (火)

Bookcover ガイコツ書店員 本田さん (1) (ジーンピクシブシリーズ) [a]
本田 / KADOKAWA/メディアファクトリー / 2016-03-26
書店員を主人公にしたマンガは久世番子さん以来たくさんあるが、これは都心型大書店を舞台にしたエッセイ・ギャグマンガ。社員研修の話がやたらに可笑しくて、声をあげて笑ってしまった。

Bookcover 決してマネしないでください。(1) (モーニング KC) [a]
蛇蔵 / 講談社 / 2014-12-22

Bookcover ちおちゃんの通学路 (4) (MFコミックス フラッパーシリーズ) [a]
川崎 直孝 / KADOKAWA/メディアファクトリー / 2016-03-23

Bookcover スズキさんはただ静かに暮らしたい 3 (ゼノンコミックス) [a]
佐藤洋寿 / 徳間書店 / 2016-03-19

Bookcover 重版出来! 7 (ビッグコミックス) [a]
松田 奈緒子 / 小学館 / 2016-03-30

Bookcover 木根さんの1人でキネマ 1 (ジェッツコミックス) [a]
アサイ / 白泉社 / 2015-12-25
実は映画好きなアラサー独身OL・木根さんを主人公にしたコメディ。主人公は面白い映画を一本紹介してよとそんなに親しくもない女友達に問われ、ウィル・スミスの「バッドボーイズ2バッド」を紹介してしまうという、業の深ーい映画マニアなのである。いいなあ、木根さんと友達になりたいなあ。

コミックス(2015-) - 読了：「ガイコツ書店員本田さん」「木根さんの１人でキネマ」「決してマネしないでください。」「ちおちゃんの通学路」「スズキさんはただ静かに暮らしたい」「重版出来！」

Bookcover BLUE GIANT 8 (ビッグコミックススペシャル) [a]
石塚 真一 / 小学館 / 2016-03-30

Bookcover ダンス・マカブル~西洋暗黒小史~ 1 (MFコミックス フラッパーシリーズ) [a]
大西巷一 / KADOKAWA / 2010-11

Bookcover めしばな刑事タチバナ 21 (トクマコミックス) [a]
坂戸佐兵衛,旅井とり / 徳間書店 / 2016-03-31

Bookcover おにぎり通信 3 ~ダメママ日記~ (愛蔵版コミックス) [a]
二ノ宮 知子 / 集英社 / 2016-03-25

Bookcover 三代目薬屋久兵衛 3 (フィールコミックス) [a]
ねむ ようこ / 祥伝社 / 2016-03-08

Bookcover 甘々と稲妻(6) (アフタヌーンKC) [a]
雨隠 ギド / 講談社 / 2016-03-07

コミックス(2015-) - 読了：「甘々と稲妻」「三代目薬屋久兵衛」「おにぎり通信」「めしばな刑事タチバナ」「ダンス・マカブル」「BLUE GIANT」

最近読んだ本。まずはコミックスから。

Bookcover パパと親父のウチご飯 1 (BUNCH COMICS) [a]
豊田 悠 / 新潮社 / 2014-12-09
「甘々と稲妻」がヒットしているときに、これを出すか... 良い度胸というか、なんというか...

Bookcover アイアムアヒーロー 19 (ビッグコミックス) [a]
花沢 健吾 / 小学館 / 2016-03-15

Bookcover 娘の家出 4 (ヤングジャンプコミックス) [a]
志村 貴子 / 集英社 / 2016-03-18

Bookcover このマンガがすごい! comics 翔んで埼玉 (Konomanga ga Sugoi!COMICS) [a]
魔夜 峰央 / 宝島社 / 2015-12-24

Bookcover 辺境で 伊図透作品集 (ビームコミックス) [a]
伊図透 / KADOKAWA/エンターブレイン / 2016-02-25
短編集。

Bookcover レベレーション(啓示)(1) (モーニング KC) [a]
山岸 凉子 / 講談社 / 2015-12-22
これは少し前に読んで、記録し損ねていたもの。ジャンヌ・ダルクの物語。

コミックス(2015-) - 読了：「翔んで埼玉」「パパと親父のウチご飯」「レベレーション(啓示)」「辺境で」「娘の家出」「アイ・アム・ア・ヒーロー」

Kyriacou, D.N. (2016) The Enduring Evolution of P Value. JAMA. 315(11), 1113-1115.
　JAMA (The Journal of the American Medical Association) 3月15日号のEditorial。この号には、1990年以降の生物・医学系論文におけるp値の使われ方を片っ端から調べましたというとんでもない論文が載ったのだそうで、この文章はその論文の露払い(?)みたいなものらしい。
　いちおうメモをとったけど、別にいらなかったかな... ま、コンパクトにまとまった良い文章であった。(←偉そう)

　いわく。
　科学文献においてP値の使用が顕著になったのはFisher(1925) "Statistical Methods for Research Workers"にさかのぼる。FisherによればP値とは「仮に帰無仮説が真であるとして、観察された結果ならびにそれより極端な結果が得られる確率」。FisherはP値を統計的推論の柔軟な指標としてとらえていた(意思決定の道具ではなくて)。また、P値の適切な使用は以下の想定と結びついているものであった。(1)いま検討している因果的要因と、関心の対象である結果との間に関連性がない(=帰無仮説が真)。(2)研究デザインと分析に体系的エラーが全くない(=誤分類とか選択バイアスとか交絡とかがない)。(3)適切な統計的検定が選択されている。
　こうしてみると、P値についての誤解や誤用など生じなさそうなものだが、実際にはもう生じまくりである。Goodman(2008)はP値についての12の誤解を挙げている。「P値は帰無仮説が真である確率だ」とか、「P値が.05以上だったら帰無仮説は真であり曝露とアウトカムの間に関連がない」とか。

　フィッシャー派とネイマン=ピアソン派では哲学は違うが、確率についての頻度主義的解釈をとる点、ある実験を統計的に独立な結果をもたらす無数の実験のひとつとみなす点は共通しており、これが臨床試験・疫学の教育の基礎となっている。
　いっぽう、違う点を挙げると...
　フィッシャー派の立場では、P値は個別的な研究上の知見について、帰無仮説に反対する証拠として解釈される。有意水準を事前に固定する必要はない。対立仮説もない。仮説を受容するか棄却するかという判断ではなく、仮説に関する推論を行うのが良しとされる。
　ネイマン=ピアソン派の立場では、目的は推論じゃなくて決定である。事前に有意水準を決めておき、P値がそれを下回ったら「統計的に有意」という決定を行い、対立仮説を支持する。この過程で生じうるエラーは2種類あることになる(Type IとType II)。ここから検定力という概念が生まれた。また、この考え方は臨床試験における標本サイズ決定の基盤となった。
　ネイマン=ピアソン派への主な批判は、決定に柔軟性がなさすぎるんじゃないかという点である。これに対して、確かさの程度を条件つき確率として特徴づけ、因果関係を帰納的に評価しようとするのがベイジアン帰納推論である。事前分布がわからない、ベイズ・ファクターはP値よりも計算が大変、教室であまり教えられてない、というのが難点。とはいえ、P値とベイジアンを融合しようという提案もあるぞ(Greenland & Poole, 2013 Epidemiology)。
　
　突き詰めていえば、P値を使った統計的推論とは、偶然誤差という文脈において説明理論を構築することを数学的に手助けしようという試みである。しかし、P値が提供するのは特定のデータセットの数学的記述に過ぎず、目標母集団における因果関係の科学的説明ではない。
　P値の真の意味、強みと限界、そして統計的推論のもっとも適切な適用方法について理解することが重要である。P値そのものは正しく使用されていれば悪いものではないけど、事前に決めた有意水準の下で二値的仮説検証を自動的にやっちゃうのは良くない。効果量の推定値とか信頼区間とかP値そのものとかを使ったもっと複雑な過程を導入しないといかんよ。科学者・統計家・臨床家が自分の推論能力でもって科学的重要性を決められるようにね。

　... Greenland & Poole (2013)って面白そうだなあ。時間ができたら探してみよう。

論文：データ解析(2015-) - 読了：Kyriacou (2016) p値を正しく理解しましょう

Bates, J.M., Granger, C.W. (1969) The Combination of Forecasts. OR. 20(4), 451-468.
予測結合の研究ですごくよく引用される古典なので、一応ざっと目を通した。

　不偏な予測が2つあるとき、どう結合するのがよいか。
　誤差分散が過去を通じてそれぞれ$\sigma^2_1, \sigma^2_2$だとしよう。重み$k, 1-k$で線形結合したとして、分散は
　$\sigma^2_c = k^2 \sigma^2_1 + (1-k)^2 \sigma^2_2 + 2\rho k \sigma_1 (1-k) \sigma_2$
　となる理屈である($\rho$は誤差の相関)。これを最小化する$k$は
　$\displaystyle k=\frac{\sigma^2_2 - \rho \sigma_1 \sigma_2}{\sigma^2_1 - \sigma^2_2 - 2 \rho \sigma_1 \sigma_2}$
　ところで、もし過去データがなかったら、最初の$k$は0.5にしといて、データが溜まるにつれて$k$を直していくしかない。時点$T$における予測結合を以下としよう。予測を$F_{1,T}, F_{2,T}$として、
　$C_T = k_T f_{1,T} + (1-k_T)f_{2,T}$

　二乗誤差最小な結合手法は以下の条件を持つだろう：(1)予測性能が定常だとして、予測数の増大とともに$k$が最適値に接近する。(2)予測の成否とともに$k$が変わる。(3)$k$の最適値まわりの変動は小さい。さらにいえば、ビジネスマンのみなさんのために、なるべく簡単な手法にしておきたい。
　以下、予測1,2の時点$t$での予測誤差を$e_{1,t}, e_{2,t}$とする。
　選手入場です。

　実データにあてはめてみると...[いろいろ比較しているけど、別に決着がつくわけでもない。関心なくなったのでパス]
　[ちょっとした修正案にも触れているけど、パス]
　云々。

　いやあ、時代の違いだなあ。現代であれば、なんらかの枠組みの下で最適な推定量を導出します、という風に話を進めるだろう。こんな風に推定量を思いつきで作って実データで比較するというような牧歌的な研究は、ちょっと許されないんじゃないかしらん。
　というわけで、途中からまじめに読んでないんだけど、もともとの研究の文脈がわかったので良しとしよう。こうしてみると、Morris(1974)が予測結合という問題にベイズ更新という観点を導入したのは、ひとつの革新だったんだな、と納得。

論文：データ解析(2015-) - 読了：Bates & Granger (1969) 2つの予測をどうやって結合するか

2016年4月11日 (月)

　「製品・サービスの普及をBassモデルで推定したいよ！でも時点数が足りないからうまく推定できないよ！」「じゃあ時点数が足りるまで待てば？」「うぎー！そんなに待ってたら推定する意味ないじゃん！」というジレンマに対し、「じゃあさ、過去データからパラメータの事前分布みたいのを作っといて、手元のデータでベイズ更新しようよ」という路線のモデルがある。
　先日読んだSultan et al.(1990)は、過去データから求めたBassパラメータと手元のデータから求めたパラメータをGoldberger-Theilアプローチというので結合する、という提案であった。

Lenk, P.J., & Rao, A.G. (1990) New Models from Old: Forecasting Product Adoption by Hierarchical Bayes Procedures. Marketing Science, 9(1), 42-53.
　いっぽう、この論文ははなっから階層ベイズモデルを使うという提案。こっちのほうが気が利いているような気がする。

　まずはBassモデルのおさらいから。潜在受容者数を$m$、時点$t$における累積受容者数を$N(t)$として、新規受容者は
　$N'(t) = \left( p + q \frac{N(t)}{m} \right) \left( m - N(t) \right)$
$F(t) = N(t)/m$として、
　$F'(t) = \left( p + q F(t) \right) \left( 1 - F(t) \right)$
これを解いて、
　$\displaystyle G(t) = c \ \frac{1 - \exp(-bt)}{1 + a \exp (-bt)}$
$G(t)$は受容までの時間の累積密度関数, $a=q/p, b=p+q, c=m/M$, $M$は母集団サイズ。
　ここで、年間売上$X_t$について
　$X_t / M = \Delta G(t) + \epsilon_t$
と定式化する。$\Delta G(t) = G(t)-G(t-1), G(0)=0$。$\epsilon_t$はiidな正規誤差と考え、その分散を
　$Var(\epsilon_t) = \Delta G(t) (1-\Delta G(t)) \sigma^2$
と仮定する。売上のレートが高ければ分散も大きいだろう、という仮定である。ほんとは$X_t / M$は非負でないとおかしいし、合計は1以下でないとおかしいけど、とりあえずこの2点の制約は忘れる。
　[このくだり、いまいち腑に落ちていない... Srinivasan & Mason(1986, MktgSci.)をみるとよいらしい]

　ここからが本題。ここに階層ベイズ(HB)モデルを導入する。まずはHBモデルとはなにか、ご紹介しましょう。[←こういうところ、時代だなあ]
　いま$K+1$群の母集団があり、それぞれ$N(\theta_k, \sigma_k^2)$に従っているとしよう。最初の$K$群のそれぞれから無作為標本$X_{k,1}, X_{k, 2}, \ldots, X_{k, n(k)}$を得る。標本平均$\bar{X}_k$は$\theta_k$の十分統計量である。$\bar{X}_k$の尤度は$N(\theta_k$, $\sigma^2_k/n(k))$から得られる。
　いま$\theta_1, \theta_2, \ldots, \theta_{K+1}$が$N(\phi, \tau^2)$からのiidドローであるとしよう。$\phi$は未知だが$\tau^2$は既知とする。これを第一段階事前分布と呼ぶ。$\phi$は$N(\mu, \xi^2)$に従うとしよう。これを第二段階事前分布と呼ぶ。
　$\bar{X}_k$の周辺分布は$N(\mu, \sigma^2_k/n(k) + \tau^2)$。最初の$K$群からの標本が観察されたとして、$\phi$の事後分布$N(\mu_K, \xi^2_K)$について考えると、$\mu_K$は標本平均$\bar{X}_k$の周辺分布の分散の逆数を重みにして、事前平均$\mu$との加重和をとったものになり、$\xi^2_K$はウェイトと$\xi^{-2}$の和の逆数となる。
　第$K+1$群からの標本を観察する前の$\theta_{K+1}$の予測分布は$N(\mu_K, \tau^2+\xi^2_K)$となる。最初の観察事例$X_{K+1, 1}$の予測分布は$N(\mu_K, \sigma^2_{K+1}+\tau^2+\xi^2_K)$となる。$n(K+1)$の標本を得たときの$\theta_{K+1}$の事後分布について考えると、それは正規分布で、平均$\tilde\theta_{K+1}$は第$K+1$群からの標本平均と$\mu_K$とを、それぞれの分散の逆数を重みにして加重和をとったものとなる。云々...
　[丁寧に説明してくださっているんだろうけど、勉強不足のせいでちょっと理解できない箇所が... ま、あとでゆっくり考えよう]

　エアコン(13年間)、カラーテレビ(8年間)、衣類用ドライヤー(13年間)、超音波診断装置(14年間)、マンモグラフィー(14年間)、外国語プログラム(12年間)、加速教育プログラム(?)(12年間)の、7本の普及曲線を使う(各時点をケースとみるわけである)。$k$番目のカテゴリのパラメータを$\theta_k = \{logit(c_k), logit(p_k), logit(q_k), 2log(\sigma)\}$とする。$\theta_k$の第一段階事前分布を$MVN(\mu, \Sigma)$とする。7カテゴリしかないので、推定の都合により$\Sigma$の非対角要素はゼロとする。ハイパーパラメータ$\phi = (\mu, \Sigma)$にはJeffreyの無情報事前分布を与える。
　6カテゴリで事前分布を推定して1カテゴリを予測、というのを繰り返す。当該の普及曲線をMLで推定した場合と比べると、初期の売上が低い場合はHBのほうが有利、6カテゴリと大きく違うカテゴリについてはMLのほうが有利。

　細かいところ、よくわからない点が多い... 勉強不足を痛感いたしました。
　90年の論文なのでMCMCは出てこない(数値積分をモンテカルロ積分でやっているけど)。そのほかにも、いまならちょっとこういう書き方はしないんじゃなかろか、という箇所があった。もうちょいと新しいのを読んで、きちんと勉強しよう。

　... それよかですね、この論文の白眉(?)は、このあとについたコメントと返答なのであります。実をいうとそっちを先に読んでいて、ちょっと笑っちゃったので、本編を探して読んだ次第である。

Vanhonacker, W.R. (1990) On Bayesian Estimation of Model Parameters. Marketing Science, 9(1), 54-55.
　Lenk-Raoの提案、実にすばらしい。いくつかコメントしたい。
　彼らのアプローチはこうだ。

このアプローチはある重要な想定に基づいている。すなわち、交換可能性だ。事前分布$N(\phi, r-2)$は交換可能な事前分布である。つまり、$\theta_k$の事前情報は異なる$k$のあいだで同じである。このことは、Bassモデルのパラメータについての事前情報が産業材でも消費財でも同じだと考えていることに相当する。これは現実にそぐわない。
　このように、HBはパラメタの異質性を確率的に表現してはいるものの、いかなるパラメータについても事前分布が弁別できない。事前知識はパラメータの異質性を含んでいないのである。ここがベイジアン・モデルの非現実的なところだ。仮に事前分布を別の形で表現しても本質はかわない。異質性が持っているシステマティックな性質を事前情報に組み込むことができていない。
　今後は、Lenk-Raoがappendixで触れているような、事前分布に説明変数を組み込むアプローチが期待される。また、これは僕らがやっているアプローチなんだけど、Goldberger-Theilの再帰推定量という考え方を用いて、異質性の既知の構造を線形のパラメータ構造に統合していくという方法も考えられる。すなわち、

このやり方だと初期予測を簡単に手に入れられるよ。

Lenk, P.J., & Rao, A.G. (1990) Reply to Wilfried R. Vanhonacker. Marketing Science, 9(1), 56-57.
　Vanhonacker先生、コメント誠にありがとうございます。光栄に存じますが、なんか勘違いしておられるようです。

　うっかりコメントしたらもうコテンパンな目に... Lenk & Rao, おそろしい子...

論文：マーケティング - 読了：Lenk & Rao (1990) 階層ベイズなBassモデル (心温まるコメントと背筋の凍る返答つき)

2016年4月 7日 (木)

Clemen, R.T., & Winkler, R.L. (1986) Combining Economic Forecasts. Journal of Business and Economic Statistics, 4(1), 39-46.
　ベイジアン合意モデルの有名論文。やれやれ、どういう話なのか、ようやくわかったよ...
　
　予測対象である変数$x$について、複数の予測値からなるベクトル$\hat{\mathbf x} = (\hat{x}_1, \ldots, \hat{x}_k)'$が得られているとしよう。これらの予測を結合する一番手っ取り早い方法は、平均
　$\displaystyle y_1 = \frac{\sum_i^k \hat{x}_i}{k} $
を求めることだが、この方法は、それぞれの予測値の精度のちがいを無視しているし、予測値が基づいている情報源の重複も無視している。
　そこでベイジアン・モデル。かつてNewbold & Granger(1974, JRSS), Winkler(1981) はこう考えた。予測誤差のベクトル $\mathbf{e} = (\hat{x}_1 - x, \ldots, \hat{x}_k - x)'$が$MVN(0, \mathbf{\Sigma)}$に従い、共分散行列$\mathbf{\Sigma}$は正定値行列だとする。$\mathbf{\Sigma}$の推定値を$\hat{\mathbf{\Sigma}}$としよう。予測の結合は
　$\displaystyle y_2 = \frac{\mathbf{u}' \hat{\mathbf{\Sigma}}^{-1} \hat{\mathbf{x}} }{\mathbf{u}' \hat{\Sigma}^{-1} \mathbf{u}}$
ただし$\mathbf{u} = (1, \ldots, 1)'$。これはつまり、重みベクトル
　$\displaystyle \mathbf{w}' = \frac{\mathbf{u}' \hat{\mathbf{\Sigma}}^{-1}}{u' \hat{\mathbf{\Sigma}}^{-1} u}$
を使って$\hat{\mathbf{x}}$の加重和をとっているわけである。
　$\hat{\mathbf{\Sigma}}$をどうやって手に入れるか。普通に考えれば、標本共分散行列をそのまま使えばよいのだが、データサイズが小さくてペアワイズの相関が高いとき(たいていそうだ)、不安定になってしまう。以下では、標本共分散行列をつかった正規モデルを出発点とし、標本共分散行列の不安定性にどうやって対処するかを考えよう。

　[ここからいろんなアイデアが紹介される：標本共分散行列を求めるのに使うデータを増やす、最近の事例を重視する、$\mathbf{\Sigma}$は対角行列だと考えてしまい非対角要素にゼロを埋める、結合した予測値が予測値ベクトルの最小値と最大値の範囲の外に出た時だけその範囲へと修正する、ウェイトが足して1だという制約を捨てて切片なしのOLS、予測の誤差の平均ベクトルが0だという想定も捨てて切片ありのOLS、前回一番良かった予測を使う、改善を期待してあえて前回一番悪かった予測を使う。いろいろ考えてみたけど、こういうの全部アドホックだよね、というわけで...]

　提案。
　$\mathbf{\Sigma}$の事前分布を、共分散行列$\mathbf{\Sigma}_0$、自由度$\alpha(>k-1)$の逆ウィシャート分布とする。すると$\mathbf{\Sigma}$の事後分布も逆ウィシャート分布となり、その共分散行列と自由度は
　$\displaystyle \mathbf{\Sigma}^* = \left( \frac{\alpha \mathbf{\Sigma}_0^{-1} + n \hat{\mathbf{\Sigma}}^{-1}}{\alpha + n} \right)^{-1}$
　$\alpha^* = \alpha + n$
となり、$x$の事後平均は、
　$\displaystyle y_3 = \frac{\mathbf{u}' \hat{\mathbf{\Sigma}}^{*-1} \hat{\mathbf{x}} }{\mathbf{u}' \hat{\mathbf{\Sigma}}^{*-1} \mathbf{u}}$
となる。
　$\mathbf{\Sigma}_0$は、すべての予測を交換可能とみなして、対角成分を$\sigma^2$, 非対角成分を$\sigma^2 \rho$とする。
　$\alpha$は等価な事前標本サイズに相当する。大きくすると$w_i$は等しくなっていく。

　実例。1971-82年の名目・実質GNP(四半期)の、4つの機関による予測と実現値を使用。4つの予測を結合した値の予測誤差を調べる。以下の方法を試す。

　結果。成績が良かったのは、単純平均、独立性を仮定した正規モデル、そしてベイジアン・モデル(独立性仮定はあってもなくてもよい)。独立性仮定がうまくワークしたのはたまたまであろうから、ベイジアンがよろしいのではないでしょうか、云々。

　。。。最初に手に取ったときは文脈がつかめずに挫折したんだけど、ようやく状況がわかってきた。こういうことじゃないかと思う。
　複数個の不偏な予測があって、これを結合したいとき、予測の共分散行列の逆行列を縦に合計し、これを和が1になるように調整し、これを重みとして予測値の加重和を求めるとよい。というアイデアは、ベイジアン・アプローチとは無関係にすごく昔からある。未読だが、Bates & Granger(1969)が最初だと思う。
　この手法をベイジアンの観点から基礎づけたのがWinkler(1981), Bordley(1982)。$x$の事前分布が一様分布で、予測の誤差が多変量正規性を持つという前提の下で、上の予測合成が$x$の事後分布の平均となることを示した。
　ここまでの話では、予測の共分散行列は基本的に既知であった。実際には、予測の共分散行列を推定するのがすごく難しい。そこで2つの路線が登場する。

2つの路線はどう違うんだろうか。おそらく、切片項が入れられるぶん、枠組みとしては後者のほうが柔軟なのであろう。逆にいうと、もし後者に切片項を入れないならば、本質的には前者と同じことだったりしそうだなあ...

論文：データ解析(2015-) - 読了：Clemen & Winkler (1986) 共分散行列もベイズ推定するベイジアン合意モデル

2016年4月 6日 (水)

Diebold, F.X., Pauly, P (1990) The use of prior information in forecast combination. International Journal of Forecasting, 6, 503-508.
　いろんなモデルによる予測を結合したら素晴らしい予測が出来上がるんじゃないかしらという話をしてるところに乗り込んで行ってそんな暇があったらモデルを作り直せと言い放つ 、樹木希林なみにロックンロールな計量経済学者Dieboldさんがご提案する、経験ベイズな予測結合手法。

　予測の結合において、MSE最小化を夢見て加重平均を用いるべしとか、諦めて単純平均を使おうとか、いろんな立場がありうるが、根本的な問題は、どの立場がいいのか事前にはわからないという点である。本研究では回帰の枠組みでこの問題に取り組む。

　時点$t-1$において、変数$y_t$についての不偏な予測$f^1_t, \ldots, f^m_t$があり、ここから
　$C_t = \omega_1 f^1_t + \omega_2 f^2_t + \cdots + (1-\sum_i^{m-1} \omega_i) f^m_t$
と作るとしよう。一期先予測誤差の共分散行列を$\Sigma$として、$C_t$の分散を最小化する重みのベクトルは
　$\omega^* = (\Sigma^{-1} i) / (i' \Sigma^{-1} i)$
である($i$は1のベクトル)。元の予測が不偏なんだからMSE最小化だともいえる。これはBates & Granger (1969) にはじまる路線である。分散共分散法と呼ぼう。

　次に、回帰に基づく予測結合。Zellnerのg事前分布モデルを用いる。
　まず普通の回帰を考える。
　$Y = F \beta + \epsilon, \ \ \epsilon \sim N(0, \sigma^2 I)$
係数に自然共役正規-ガンマ事前分布を与えて
　$P_0(\beta, \sigma) \propto \sigma^{-K-v_0-1} \exp((-\frac{1}{2\sigma^2})(v_0 s_0^2+(\beta-\beta_0)'M(\beta-\beta_0)))$
ふつう$K=m+1$(切片があるから)。
　さて、ここで$M=gF'F$とすると、事後分布の平均ベクトルは
　$\beta_1 = \beta_0 + \frac{1}{1-g}(\hat{\beta}-\beta_0)$
$g$はシュリンケージ・パラメータになっているわけである。これをg事前分布推定量と呼ぼう。

　この$g$をデータから決める方法を考えよう。
　回帰式に戻り、$\sigma$の下での$\beta$の条件付き確率の事前分布を正規とみなして
　$P_0 (\beta | \sigma) = N(\beta_0, \tau^2 A^{-1})$
とすると、事後分布の平均ベクトルは
　$\beta_1 = \beta_0 + \frac{1}{1-\hat{\sigma}^2/\hat{\tau}^2}(\hat{\beta}-\beta_0)$
となる。g事前分布推定量の$g$が$\hat{\sigma}^2/\hat{\tau}^2$にすり替わったわけである。$\hat{\sigma}^2$と$\hat{\tau}^2$は推定量があって...[省略]。これを経験ベイズ推定量と呼ぼう(さっきのも経験ベイズではあるんだけどね)。

　実例。名目・実質GNPについての4社の予測データを用いる。一期先予測について、4社の予測、分散共分散法での結合予測、経験ベイズ推定量による結合予測、いろんなgでのg事前分布推定量による結合予測、を比較。
　分散共分散法、OLS回帰による予測(つまりg=0のg事前分布推定量)はぱっとしない。経験ベイズはぐっと良くなる。g事前分布推定量はgが大きいときに良い。gが無限大だということは算術平均だということなので、つまり算術平均への大きなシュリンケージで良くなったという次第である。

　... うぐぐぐ。大筋はつかめたけど、細かいところがわからない。やはりg事前分布の話をきちんと勉強せなあかんのか... つらい...

論文：データ解析(2015-) - 読了：Diebold & Pauly (1990) ベイジアン・シュリンケージによる予測結合

江利川滋, 山田一成 (2015) Web調査の回答形式の違いが結果に及ぼす影響：複数回答形式と個別強制選択形式の比較. 社会心理学研究, 31(2), 112-119.
　サーベイ調査の設問形式によって回答がどう変わるかという実験。Rasinski, Mingay, & Bradburn(1994 POQ), Smyth, Dillman, Christian, & Stern (2006 POQ)の追試に相当。第一著者はTBSの方だそうだ。

　一都三県の20-60代にweb調査(性年代割付), 1559名聴取。実査会社の社名は謝辞にも出てこないのでわからない。「インターネット利用行動」19項目と「ノートPC購入重視点」19項目の2設問について回答を求めた。1設問で1画面。回答形式が要因になっていて、(1)個別強制選択(FC; 項目別に該当と非該当のボタンが横に並ぶ)、(2)複数回答(MA; 項目先頭にチェックボックス)、(3)複数回答で項目順が逆順。
　結果。MAでは、選択数が少なく、回答時間が短く、後半の項目で選択率が上がる。回答時間が短い人のほうが項目数が少ないのはMAだけじゃないかと思って調べたが、結果ははっきりしなかった。
　考察。MAではKrosnick(1991 App.Cog.Psy; 1999 Ann.Rev.Psy)いうところの「弱いsatisficing」、つまり認知的努力が不十分な回答行動が生じやすいのだろう。逆にFCで黙従傾向が強くなるんじゃないかという反論も可能だが[...状況証拠でディフェンス]。
　云々。

論文：調査方法論 - 読了：江利川, 山田(2015) 2択SAとMAのちがい

2016年4月 4日 (月)

　経済予測でもビジネス上の予測でもなんでもいいんだけど、異なる手法による複数の予想が目の前にあるとき、それを結合すれば最強の予想になるじゃん！... というアイデアについて、Clemen(1989)が長大なレビュー論文を書いている。あとで気が付いたんだけど、この論文には6人の識者によるコメントがつけられていた。
　うち5人のコメントをメモ。うっかりしていて、Horgarthという人のコメントは未入手である。

Armstrong, J.S. (1989) Combining forecast: The end of the beginning or the beginning of the end? International Journal of Forecasting, 5, 585-588.
　Clemenも述べているように、予想結合についてこれまでに分かっている主な結論は2つある。(1)予想の結合は誤差を減らす。(2)単純な平均でうまくいく。
　逆にいうと、これまでの研究から上記以上のガイドラインは得られそうにない。その意味で、我々は終わりの始まりにいる。今後は以下の方向の研究が必要だ。

　今後の課題:

Diebold, F.X. (1989) Forecast combination and encompassing: Reconciling two divergent literatures. International Journal of Forecasting, 5, 589-592.
　計量経済学者にとっては、経済システムにおける構造的諸関係の理解が主目的であり、予想は二の次だ。このパラダイムでは、異なるモデルによる予想の結合などお呼びでない。
　この10～15年の間に、計量経済学は理論重視からデータ重視へと大きく舵を切った。我々はモデルの誤指定の可能性を認めるようになり、既知の真のモデル形式をどうやって想定するかではなく、データに最良の説明を与えるモデルをどうやって見つけるかに焦点を当てるようになった。しかしそれでも、異なるモデルによる予想の結合はほぼお呼びでない。そんな暇があったらモデルを改善すべきだ、と考えるのがふつうである。
　そもそも予想結合の正しい使い方とは何か。情報集合が即時に結合できる世界であれば、予想ではなく情報集合を結合すべきである。いっぽう、せまる締切か何かのせいで、情報集合の結合が難しい場合も多い。これが予想結合のレゾンデートルである。
　計量経済学は長期的視野に立ち、データ生成過程の特徴をうまく捉えるモデルへの累積的な科学的接近を良しとする。プーリングは短期的視野に立ち、対立をどうやって取り繕うかを考える。ま、どっちも大事ではある。
　往年の女優メイ・ウェストは「いいことがありすぎるのってステキ」("Too much of a good things is wonderful")といったが、予想結合にはあてはまらない。予想のユーザがモデルのビルダーでもあるならば、情報集合を結合し、真剣な分析を行うべきだ。
　もし共分散が既知ならば、最適なプーリング予想のMPSE(予測誤差の平均平方)は最良の予想のMPSE以下となるはずだ。しかし実際にはそうはならない(ウェイトにサンプリング誤差が乗っているから。共線性のせいでさらに深刻になる)。そのため、現実的には単純平均が最良となることが多い。Diebold & Pauly (1990, Int.J.Forecasting)は、ベイジアン・シュリンケージによってさまざまな程度の事前情報を統合する方法を示している。
　云々。

Mahmoud, E. (1989) Combining forecasts: Some managerial issues. International Journal of Forecasting. 5, 599-600.
　Clemenさんのレビューにはビジネス分野での事例が少ない。主観的予想と定量的予想を組みあわせている例はいっぱいあるのだけれど。[と、簡単なレビュー...]

Makridakis, S. (1989) Why combining works? International Journal of Forecasting. 5, 601-603.
　予想の結合で誤差が減るのは、予想の誤差を生む諸要因を平均するからである。諸要因とはすなわち、(1)測るべきものを測ってない。需要の代わりに注文数や生産量を使っている、とか。(2)測定誤差。(3)パターンや関係性が不安定だったり変化してたりする。(4)モデルが過去の誤差を最小化している。
　予想結合を改善するポイント: (1)手法の成績をチェックし続け、良い手法だけを結合する。(2)お互いに補完的な関係にある手法を結合する。保守的な手法と過去トレンドの補外を結合するとか。(3)意思決定者からの主観的判断をひきだす。
　云々。

Winkler, R.L. (1989) Combining forecasts: A philosophical basis and some current issues. International Journal of Forecasting. 5, 605-609.
　予想結合の研究史が意外に浅い理由のひとつは、予想にかかわる人々の多くが統計モデリングの伝統のなかで育ち、真のモデルの探索を重視してきた、という点であろう。もっとも、この伝統と予想結合は必ずしも矛盾しないのだが(重みを「そのモデルが真である確率」と捉えればよい)、私としてはむしろ、この流れゆく世界において真のモデルなんかないのだという立場に立ちたい。
　人々の統計的バックグラウンドによるもう一つの帰結として、理論なり過去データなりに基盤を置くという点がある。でも主観的判断も大事だ(これはベイジアンの枠組みでうまく形式化できるだろう)。ポイントは、複数の予測は必ずしも真のモデルを目指して競争しているわけじゃないという点である。
　取り組むべき課題:

　...一番おもしろかったのはDieboldという人のコメントであった。異なる予測モデル1, 2, 3...を所与としてその出力を結合することを考えるより、より良いモデルを作ることを目指すべきだというのは、なるほど正論である。実のところ、これは予想結合だけではなく、もっと最近のベイジアン・モデル平均に対してもモヤモヤと感じていた点なのだけれど、この分野じゃそういう風には考えないんだろうな、と諦めていた。そうか、こういう批判はアリなのか。
　Dieboldさんとほぼ同じ議論を展開しているにも関わらず、Winklerさんのほうは予想結合に好意的だ、という点も面白い。これは要するに、モデルというものをどう捉えるかという哲学的なスタンスの問題であろう。

論文：データ解析(2015-) - 読了: Armstrong, Diebold, Mahmoud, Makridakis, Winkler (1989) 予想結合をめぐる識者のコメント

Bordley, R.F. (1982) The combination of forecasts: a Bayesian approach. Journal of Operation Research Society, 171-174.
　これも、80年代のベイジアン合意研究でよく引用される論文。中身はWinkler(1981)とほぼ同じなのだが、独立に書かれた論文だそうだ(と、どこかに書いてあったが、どこだっけ... Clemenの文献レビューかな...)。

　いわく。
　$X$についての$n$個の予想$z_1, \ldots, z_n$があるとき、その加重和
　$Q = \sum_k w_k z_k$
は、$z_1, \ldots, z_n$が共分散行列$S$を持つ確率変数であるとき
　$Q = (I^t S^{-1} I)^{-1} (I S^{-1} z)$
とすることで分散最小になる。ということは、実はBates & Granger (1969) などがとっくに示している。問題は、これが最良の推定量だといえるのがどんなときか、だ。これにベイジアン・アプローチで答えようとしたのがBunn(1975)なのだが、いくつか欠点がある[説明省略]。厳密な解を示そう。

　ベイジアン・アプローチでいうと、予想を受け取った意思決定者の事後分布は
　$f(x | z_1, \ldots, z_n) = \frac{f(z_1, \ldots, z_n | x) f(x)}{\int 分子 dx}$
なんだけど、$f(x)$がわかんないので、Edwards et al.(1963)に従い、有限区間$B$について
　$f(x | z_1, \ldots, z_n) = \frac{f(z_1, \ldots, z_n | x) }{\int 分子 dx}$
とする。[←これ、要するに拡散事前分布を使うってことだよね...当時はそういういい方はなかったのだろうか]

　いま、$f(z_1, \ldots, z_n | x)$が真値まわりのMVNに従い、共分散行列が$S$だとすると、[...途中大幅に端折って..] 上式の$Q$が事後分布の平均、つまり最良の推定量であることが示せる。
　いっぽう、$f(z_1, \ldots, z_n | x)$が対数正規分布をとるならば(汚染レベルの評価なんかはそうですね。過小評価は小さいが過大評価は大きい)、重みつきの幾何平均が最良になる。このように、最良の推定量は誤差の分布次第である。

論文：データ解析(2015-) - 読了：Bordley (1982) ベイジアン合意モデル

2016年4月 3日 (日)

Winkler, R.L. (1981) Combining probability distributions from dependent information sources. Management Science, 27(4), 479-488.
　80年代のベイジアン合意モデル研究として、よく引用されている論文。

　いわく。　
　異なる情報源から得られた推定や予測(これを「専門家」と呼ぶ)から、どうやってひとつの値(「合意」)を得るか。ここでポイントになるのは、専門家のあいだに確率的(stochastic)な依存性があるかもしれないということ。これは先行するベイジアンモデルであるMorris (1977, MgmtSci.)なんかがあまり注目してこなかった点だ。

　関心ある変数$\theta$の分布について$k$人の専門家が評価するとしよう。話の都合上、$\theta$は実数、$k$個の予測分布はみんな連続とする。専門家$i$の密度関数を$g_i$とする。その平均
　$\mu_i = \int \theta g_i (\theta) d\theta$
を$\theta$の点推定値とみることができる。誤差を
　$u_i = \mu_i - \theta$
と書く。誤差の密度を$g_i (\mu_i - u_i)$と書く。
　意思決定者から見た$\vec{u} = (u_1, \ldots, u_k)^t$の密度関数を$f$とする。意思決定者の立場からは、$\theta$についての知識は$\vec{u}$の密度関数を変えないと想定する。つまり、$f$の関数形は$\theta$や$\vec{\mu} = (\mu_1, \ldots, \mu_k)^t$に依存しないと想定する。なお、この想定は緩めてもよい(話がややこしくなるけど)。
　意思決定者は$\theta$の事前分布を持っていて、専門家の意見によってそれを更新し事後分布を得る。ここで、$f(\mu_1 - \theta, \ldots, \mu_k - \theta)$を尤度関数とみることができるから、事前分布を$h_0(\theta)$とすれば、事後分布は
　$h(\theta | g_1, ..., g_k, f, h_0) \propto h_0(\theta) f(\mu_1 - \theta, \ldots, \mu_k - \theta)$

　さあ、準備はできた。$f$をどうやって決めるか。
　ここでは、$f$は$k$変量正規分布で、平均はすべて$0$、共分散行列$\Sigma$は正定値行列だ、と仮定する。この仮定が無理めな場合でも、適切なキャリブレーションによってこれに近づけることができるんとちゃいますかね。

　$\Sigma$が既知であるとしよう。
　$\theta$について拡散事前分布を想定すると、事後分布はこうなる。$1$の縦ベクトルを$\vec{e}$とし、
　$\mu^* =\vec{e}^t \Sigma^{-1} \vec{\mu} / \vec{e}^t \Sigma^{-1} \vec{e} $
　$\sigma^{*2} = 1 / \vec{e}^t \Sigma^{-1} \vec{e} $
として、標準正規密度関数を$\phi$として、
　$h(\theta | \vec{\mu}) \propto \phi[(\theta - \mu^*) / \sigma^*]$
　具体的に言うとこうだ。事後分布の平均は、個々の$\mu_i$の重み付け和 $\sum_i^k w_i \mu_i$となり、重みは共分散行列の逆行列 $\Sigma^{-1} = (\alpha_{ij})$を縦に足し上げた値に比例させた合計1の重み
　$w_i = \sum_j^k \alpha_{ij} / \sum_m^k \sum_j^k \alpha_{mj}$
となるわけである。相関が高いときは重みが負になっちゃっても全然おかしくないという点に注意。
　[数値例。省略]
　$k$が大きくなると$\Sigma$を与えるのが大変になるので、相関行列の非対角成分が全部同じ、という制約をかけて考えると...[略]

　$\Sigma$が未知である場合。
　$\theta$と$\Sigma$はアプリオリに独立だと想定します。$\theta$には拡散事前分布を想定します。$\Sigma$には逆ウィシャート事前分布を与え、ハイパーパラメータを$\delta_0, \Sigma_0$とします。。[...中略...]事後分布は$t$分布となり、自由度は$\delta_0+k-1$、平均は
　$m^* = \vec{e}^t \Sigma_0^{-1} \vec{\mu} / \vec{e}^t \Sigma_0^{-1} \vec{e}$
分散は[...省略]。事後分布の平均をよくみると、$\Sigma$が$\Sigma_0$にすり替わっただけであります。
　[数値例。省略]

　。。。なるほどねえ。ベイジアン合意って、もともとこういう枠組みだったのか。この論文をもっと早く読んでおけばよかったな。

　ところで、わからないところを検索していてうっかり見つけてしまったのだけれど、この論文の内容について、言い回しや数値例まで含めほぼ逐語的であるとさえいえる、解説というかなんというか、を提供している日本語の紀要論文があった(美しい表現)。人によって見方が分かれるところかもしれないが、そういうのもアリなんじゃないか、と、これは本気で思う。非専門家としては、誰かに日本語にして頂けるだけで助かります。さらにいえば、いま自分が調べていることがほんとに誰かの役に立つ話なのか、だんだん不安になってきていたところだったので、その意味でもほっとした。

論文：データ解析(2015-) - 読了：Winkler(1981) ベイジアン合意モデル (正規性仮定バージョン)

« 2016年3月 | メイン | 2016年5月 »

rebuilt: 2020年11月16日 22:40
validate this page