elsur.jpn.org >

« 2017年8月 | メイン | 2017年10月 »

2017年9月28日 (木)

Finley, A.O, Banerjee, S., Gelfand, A. (2015) spBayes for large univariate and multivariate point-referenced spatio-temporal data models. Journal of Statistical Software, 63(13).
　Rのベイジアン空間統計パッケージ spBayes の紹介。仕事の都合で大急ぎで読んだ。

　spBayesパッケージが扱うベイジアン・ガウシアン空間回帰モデルは、一般化して書くと、以下の階層線形混合モデルである。
　$p(\mathbf{\theta})$
　$\times N(\mathbf{\beta} | \mathbf{\mu}_\beta, \mathbf{\Sigma}_\beta)$
　$\times N(\mathbf{\alpha} | \mathbf{0}, \mathbf{K}(\mathbf{\theta}))$
　$\times N(\mathbf{y} | \mathbf{X\beta} + \mathbf{Z}(\mathbf{\theta}) \mathbf{\alpha}, \mathbf{D}(\mathbf{\theta}))$
　まず4つめからいこう。$\mathbf{y}$は$n \times 1$の観察ベクトル。$\mathbf{X}$はサイズ$n \times p$の既知の共変量で、$\mathbf{\beta}$は$p \times 1$の傾きベクトル。$\mathbf{Z}$は$n \times r$の行列で、$\mathbf{\alpha}$は$r \times 1$の正規ランダムベクトル。$\mathbf{D}$は$n \times n$の共分散行列。
　2つめと3つめは、それぞれ$\mathbf{\beta}, \mathbf{\alpha}$の事前分布を表している。
　で、観察ノイズの共分散行列$\mathbf{D}$, $\alpha$にかける係数$\mathbf{Z}$、潜在変数の共分散行列$\mathbf{K}$、の3つは、未知のプロセス・パラメータ$\theta$の関数となっている。その事前分布が1つめ。

　プロセス・パラメータ$\mathbf{\theta}$をどうやってサンプリングするかというと...[パス]
　傾き$\mathbf{\beta}$とランダム効果$\mathbf{\alpha}$をどうやってサンプリングするかというと...[パス]
　低ランクモデル、つまり$r$が$n$よりも全然小さいときはどうするかというと...[パス]
　予測はどうするかというと...[パス]
　
　具体的なモデリングの話をしましょう。
　その１、フルランク単変量ガウシアン空間回帰。
　位置を$\mathbf{s}$として、
　$y (\mathbf{s}) = \mathbf{x} (\mathbf{s})^T \mathbf{\beta} + w(\mathbf{s}) + \epsilon(\mathbf{s})$
　3つめはホワイトノイズで、位置と無関係にiidに$N(0, \tau^2)$に従う。ここで$\tau^2$をナゲットと申します。[←バリオグラムの特徴についての特有の言い回し。これ、きっと鉱業の方面に由来しているんだろうなあ]
　2つめの$w(\mathbf{s})$がガウシアン空間過程を表す。その共分散関数を$C(s,t)$として、定常性$C(s,t) = C(s-t)$と等方性$C(s,t) = C(||s-t||)$を仮定し、さらに
　$C(s,t) = \sigma^2 \rho(s, t; \phi)$
とする。相関関数$\rho(s, t; \phi)$について、spBayesパッケージは指数関数、Matern関数などをご用意しておりますです。
　最初の一般化した式と比べると、$\alpha$は$w(\mathbf{s})$。係数$\mathbf{Z}(\mathbf{\theta})$は単位行列。$\mathbf{K}(\mathbf{\theta})$は共分散関数$C(s,t)$。$\mathbf{D}(\mathbf{\theta})$は単位行列にナゲットを掛けた奴。つまり、ここでの$\mathbf{\theta}$は、ナゲット$\tau^2$、分散$\sigma^2$、そして相関関数のなんらかのパラメータ$\mathbf{\phi}$なわけです。
　分析例。spLM関数を使って...[メモ省略]

　その２、低ランク予測過程モデル。[ちゃんと読んでないけど、低ランクモデルの一種としてpredictive process modelというのがあるらしい。$r$個の位置をあらかじめ決め打ちしちゃうのかな？よくわからんけどパス]

　こんどは多変量の場合... [パス]
　非ガウスの場合として、二値のロジット回帰ないしプロビット回帰、カウントのプロビット回帰、をご用意しております...
　最後に、動的時空間モデルについて... [パス]

　...というわけで、ほとんど読んでないけど、なにができそうかがなんとなく分かったので読了ということにしておこう。だんだん飽きてきたし。[←ひどいなあ]

論文：データ解析(2015-) - 読了：Finley, Banerjee, Gelfand (2015) RのspBayesパッケージ

2017年9月21日 (木)

佐藤忠彦・樋口知之 (2013) ビッグデータを用いたマーケティングモデル : データ同化適用の可能性. シミュレーション, 32(4), 50-56.
　泣く子も黙る(?) 偉い先生方による、「データ同化はマーケティングにも適用できるんじゃないか」という意見論文。いまやデータは山ほどあるんだから、状態空間モデルで消費者異質性を捉え、エージェント・ベースでシミュレーションやればいいじゃん、という主旨。
　勉強不足でよく理解できていないところもあるのだけれど、あまりに面白かったので、こないだの業界団体のセミナーでネタにしてしまった... 精進いたしますです。

論文：データ解析(2015-) - 読了：佐藤・樋口(2013) マーケティングにおけるデータ同化の可能性

柘植隆弘(2014) 「表明選好法と熟議型貨幣評価」. 坂井(編)「メカニズムデザインと意思決定のフロンティア」, 慶應義塾大学出版会.
　日本語の本の1章で、そんなに長いものでもないんだけど、慣れない分野なのでメモをとりながら読んだ。

　環境経済学ではさまざまな環境評価手法(環境の価値の経済的評価)が開発されている。次の2つに分けられる。

　その1、顕示選好法。人々の行動によって評価する。代替法、トラベルコスト法、ヘドニック価格法など。非利用価値(「生態系を守る」的な、自分が利用しなくても得られる価値)は測れない。

　その2, 表明選好法。代表的な手法は、

　表明選好法は対象者が消費者選好に基づき決定することを前提にしているが、選好を持っているかどうか怪しいし、社会的問題についての意思決定ってのはむしろ非利己的な市民としての選好に基づくのかもしれない。
　実験例：

　最近では熟議型貨幣評価というのが提案されている。これは討論型世論調査なんかの環境評価版。ワークショップみたいのをやったあとで環境を評価させる。手法はいろいろあって、体系的に確立していない。
　実証研究は2つにわかれる。(1)選好形成に主眼を置いた研究。WTPは個人で決める。グループでの議論はその支援である。(2)市民選好に主眼を置いた研究。WTP自体も集団で決めちゃう。
　著者らの研究例:

　今後の課題。DMVについて以下の点の改善が求められる：(1)人数が多いと無理、(2)コストが高い、(3)評価額の理論的基礎づけが足りん。[←なるほど]

　... いやー、DMVって面白いなあ！ 仕事ともものすごく関係が深い。なんとなく衝動買いして、軽い気持ちで読み始めたんだけど、こいつは良いものを読んだ。時間ができたら勉強してみよう。

　ちょっと気になった点をメモ。

論文：調査方法論 - 読了：柘植(2014) 表明選好法と熟議型貨幣評価

Axley, S.R. (1984) Managerial and Organizational Communication in Terms of the Conduit Metaphor. The Academy of Management Review, 9(3), 428-437.
　経営・組織研究におけるReddyの導管メタファ説の意義を語る、意見論文というかエッセイというか。
　良く引用されているので気になっていた文献。このたび発表準備のついでに読んでみた(現実逃避ともいう)。あいにくPDFが手に入らなかったので、JSTORの無料レンタルで、ディスプレイ上で...

　まず、人間とはシンボルを使う動物だ、というような能書きがあって...それこそアリストテレスから始まって、古今東西の学者の名前がてんこ盛り。こういう書き出しがお好きな方って、いらっしゃいますよね...

　人間の行動と思考はメタファによって組織化されている... [という説明が半頁]
　Reddy(1979)いわく、こうしたメタファの分析は言語そのものに適用できる。コミュニケーションについての英語の表現の多くは、(1)言語は人から人へと思考・感覚を運搬するものだ、(2)話し手・書き手は語に思考・感覚を挿入する、(3)語は思考・感覚を含む、(4)聞き手読み手は思考・感覚を語から抽出する、という想定を比喩的に表現している。Reddyはこれを導管メタファ(conduit metaphor)と呼んだ。[英語表現の例を使った説明が1頁]

　いっぽうコミュニケーション研究者はコミュニケーションについて異なる見方をしている。たとえばRedding(1972, 書籍)をみよ。
　コミュニケーションにおいて生じるのは意味の運搬ではなく意味の創造だ。[ここでまたぐだぐだと...] 聞き手がどんなメッセージを見出すかはわからない。そして問題となるのは聞き手が見出したメッセージだ。[ここでまた延々と...] コミュニケーションにはある程度の冗長性が必要だ。伝言を繰り返すとメッセージはどんどん変わる(serial communication effect)。

　多くの教科書における組織的・経営的コミュニケーションについての説明は導管メタファと整合している。このように、組織において、導管メタファに由来する「コミュニケーションなんて努力しなくても大丈夫」仮定が広がっている。
　その結果、以下の事態が生じるのではないか。(1)組織内のコミュニケータが自分のコミュニケーション効率性を過大視してしまう。(2)マネージャーがコミュニケーションの冗長性を嫌い、努力を怠る。(3)ミスコミュニケーションの危険が見過ごされる。(4)組織内教育においてコミュニケーションのスキルと行動が軽視される。　

　今後の研究の方向：(1)上で指摘した可能性の実証。(2)ある人が導管メタファにどのくらいコミットしているかを評価する手段の開発。(3)導管メタファへのコミット度と他の変数との関連性についての検討。他の変数としては、コミュニケーションの効率性, 努力の程度, 冗長性が必要だという要求, フィードバックへの努力, 冗長性の程度, 組織内の諸問題がコミュニケーションの問題だと捉えられる程度、が挙げられるだろう。
　云々。

　... いやはや、文章が大変くどく、言い回しが妙に文学的で辟易した。昔の学者さんだなあという感じ。ディスプレイ上で読んでいたせいもあって(言い訳)、あんまりきちんと読んでない。
　いい点を探すと、「ある人が導管メタファに毒されている程度を測定しようぜ」というのは面白いなあと思った。もともと導管メタファというのは社会文化のレベルの話であって、個人差があるというような話ではないと思うのだけれど、さすが、視点が組織研究っぽい。

論文：心理 - 読了：Axley (1984) 経営・組織研究者よ、導管メタファに注目せよ

2017年9月17日 (日)

Kappes, H.B., Morewege, C.K. (2016) Mental simulation as substitute for experience. Social and Personality Psychology Compass, 2016, 1-16.
　社会心理学分野でのメンタル・シミュレーション研究レビュー。仕事の役に立つかと思って読んだ。
　聞いたことない誌名だし、webをみてもIFが書いてないし、不安だったんだけど、2007年発刊のオンライン・ジャーナルなのだそうだ。ScimagoというScopusがやってるランキングでは社会心理で27位だったので、じゃあ目を通してみるか、と。
　
　Taylor et al. (1998 Am.Psy.)いわく、メンタル・シミュレーションとは一連の事象のimitative なエピソード表象である。それはふつう特定の事象(現実かどうかは別にして)の詳細な心的表象を含む。
　メンタル・シミュレーションに関わる神経システムと概念システムは、シミュレート対象の行動に関与する知覚運動システムと重なっている。[イメージング研究とかをいくつか引用...] メンタル・シミュレーションは経験の代替物として働くことがあると思われる。

　代替効果について。
　これは人間の資源が有限だという認識に由来する概念である。経済学とかマーケティングとかだと、複数の選択肢に対する時間なり金銭なりの支出は負の関連を持つことがある。これが代替効果。目標志向的行為はモチベーショナルな特性を持っているから代替効果が生じる。[...しばらく社会心理系の研究の引用が続き...]

　行為じゃなくてメンタル・シミュレーションでも代替効果が起きる。以下では４つのカテゴリ(網羅的じゃないけど)にわけてレビュー。

A) 物理的エビデンスの代替をシミュレーションが提供する
事象のシミュレーションが価値を持つ:

シミュレーションで期待が変わる(cf.期待-価値モデル):

　シミュレーションが期待に影響するのは利用可能性を高めるからだと考えられてきた。でもそれだけではない。
　一般に人は、自分が行為者であるときよりも観察者であるとき、行為の原因を傾性に帰属しやすい。これがシミュレーションにもあてはまる。Libby & Eibach(2011 Adv.ESP), Libby, Shaeffer, Eiback, & Slemmer (2007 Psych.Sci.), Vasques & Bueler(2007 PSPB). [←なるほどね...シミュレートされた行為の原因を自分の傾性に帰属するって訳か。面白い説明だ]

B) 物理的練習(practice)の代替をシミュレーションが提供する
　[本項、いま関心ないのでかなり端折る]
　この現象はimaginary practice, covert rehearsal, symbolic rehearsal, introspective rehearsalなどと呼ばれている。
　 認知課題や運動課題でスキルが必要な奴にメンタル・シミュレーションが良く効く。このタイプの研究はいっぱいあってメタ分析まである。とはいえメンタル・シミュレーションだけだとあんまり効かず、実際の練習との併用が必要。
　なぜ効くのか。2説ある。(1)複雑な行為の表象的枠組みが作られたり再構築されたりして、行為の単位のチャンキングと関連づけが起き、計画が促進されるから。(2)実際と同じ視覚・筋運動プログラムが活性化されるから。最近の研究は(1)を支持。

C) 実際の消費の代替をシミュレーションが提供する
　[ここは関心があるので丁寧にメモると...]
　食刺激に関連する感覚手がかりを想像すると、実際に接触したときと類似したsensitizationが生じる(初回接触における刺激への反応性が増す):

　実際の刺激への繰り返し接触が引き起こすhabituation(モチベーショナルな反応が減ること)やsatiation(ヘドニックな反応が減ること)が、刺激接触の想像でも起きる:

食経験のメンタル・シミュレーションは自発的に起きる:

自発的シミュレーションが実際の行為を促進する:

　背後のプロセスについて直接示した研究はないけれど、記憶における心的表象が、初期状態→活性化→減衰→初期状態と移行しているからだろう。cf. Epstein et al.(2009 Psych.Rev.)。

D) 実際の達成の代替をシミュレーションが提供する
　未来を想像することで満足の先延ばしを可能にしその場の行為を抑制するのは、実際にその場で満足を手に入れようとするのがまずい場合には有益だけど、達成の代替物になっちゃって目標追求を阻害することもある。

　背後のプロセスは2つある。(1)想像上の成功が現実とごっちゃになっちゃって努力しなくなるから。(2)想像上の成功がプランニングを阻害するから。

統合
　このように、メンタル・シミュレーションは経験の代替となり、さまざまな代替効果を引き起こす。

云々、云々...

　...オンライン・ジャーナルと思って正直なめてましたが、きちんとしたレビューであった。疲れた。

　読んでいる途中で気がついたんだけど... メンタル・シミュレーションって、なんとなく技能学習の文脈のなかで捉えていたのだが(この辺は伝統的な心理学の体系に知らないうちに縛られているのかも知れない)、質問紙調査法の研究でいうところの質問-行動効果とか、いま流行の身体化認知とか、いろんな問題の基盤として考えることができる。
　調査手法研究の文脈でいうと、メンタル・シミュレーションは回答時課題であって調査者にとってポジティブなもの、質問-行動効果は回答行動がもたらす予期せざる影響であって調査者にとってネガティブなもの、身体化認知は意外な回答方略でありポジティブともネガティブといえない... という漠然とした思い込みがあったけれど、著者らが最後に述べているように、認知的プロセスはかなり共通しているのだろう。
　なんだかなあ。面白そうでいやになるなあ。

　もう一点、面白かったのは、こういうレビューはべつに網羅的でなくていいんだな、ということ。
　率直にいって、著者らが選んだ４つのテーマがメンタル・シミュレーションの研究史をどこまで包括しているのかわからないし、なぜこの４つを選んだのかという点にも、ちょっと恣意的な感じを受ける。でも、いいんだなあ、これで。４つのテーマを統合して、あるビジョンを示したわけだから。

論文：心理 - 読了：Kappes & Morewedge (2016) メンタル・シミュレーションが引き起こす代替効果レビュー

2017年9月15日 (金)

駒木伸比古(2016) 経済センサス実施にともなう商業統計の変容とその利用, E-jounal GEO, 11(1), 154-163.
　商業・流通業の研究者向けに、経済センサスと商業統計について解説する内容。あれ、話が入り組んでいてよくわからなかったので、助かりましたです。

論文：マーケティング - 読了：駒木(2016) 経済センサスと商業統計

McMahan, R.P, Gorton, D., Gresock, J., McConnel, W., Bowman, D. (2006) Separating the effects of level of immersion and 3D interaction techniques. Proc.ACM Symp. Virtual Reality Software and Technology. 108-111.
　これも仕事の都合で目を通したやつ。今気が付いたけど、なんてこった、この前に読んだComputer誌の記事と同じ著者だ... 探し方が悪いのだろうか...

　いわく。
　仮想環境の視覚的没入性が課題の成績にどう効くか、という研究は多い。また入力装置や3D相互作用技術のちがいがどう効くかという研究も多い。でも両者を一緒にみた研究がない。だから実験やりました。

　課題は３次元空間でモノを動かすゲーム。
　ディスプレイはCAVEという4面プロジェクタ(装置故障のため床面はなし)。立体視と頭部トラッキングを提供する。
　相互作用技術は[...あんまし関心ないので読み飛ばしちゃったけど、要するに...]杖だかなにかを使ってリアルに操作する奴(HOMERとGo-Go)と、マウスでどうにかして操作する奴(DO-IT)を用意。
　実験計画。要因は、立体視({あり, なし})、呈示面の数({1面, 3面})、相互作用技術({HOMER, Go-Go, DO-IT})。立体視・呈示面は被験者内、相互作用技術は被験者間で操作。
　被験者は各群3人、計12名。[←す、少ない... 心理学でいえば学部の卒論のレベルだ... 工学系はこれでもいいのかな...] 被験者内要因の水準が2x2=4、各1試行で計4試行。

　結果。[...詳細は端折って、要するに...] 相互作用技術はゲームの成績に効くけど、立体視や呈示面の数は効きませんでした。[←はっはっは。チャートをみると、マウスでやるのがよっぽどつらいゲームらしい]
　効いたのは相互作用技術自体なのか、それとも入力装置なのかはわからない。今後の課題です。
　云々。

　... いっちゃなんだが、実験はかなりしょぼい。どこの卒論かと。でもあれなんでしょうね、見た目より全然大変なんでしょうね。マウスで操作するシステムは自前で作りましたって書いてあるし。
　実験の結果はどこまで一般性があるのか、よくわからない。ゲームの性質がちがえば視覚的没入性が効くかもしれない。このケースは、ゲームは比較的簡単、でもマウスで操作するのがよっぽど難しかったということではないだろうか。
　それはともかく、問題設定が面白いと思った。ほんとに知りたいのは、視覚入力におけるリアルさみたいなものと、環境との相互作用のリアルさみたいなものが、どう相互作用するのか (ことば遊びみたいだけど) ということなのであろう。

論文：心理 - 読了：McMahan, et al. (2006) VRで大事なのは視覚的没入性か相互作用テクニックか、実験してみた

Bowman, D.A,, McMahan, R.P. (2007) Virtual reality: How much immersion is enough?, Computer, 40(7), 36-43.
　仕事の都合で読んだ(なんというか、節操のない仕事だ...)。掲載誌はIEEEの一般向け機関誌みたいな奴だと思う。結構引用されている(Google様的には366件)。
　眠くてふらふらになりつつメモをとったので、中身には自信がない。

　いわく。
　これまでVR研究者は没入的なVRに取り組んできた。ユーザは驚いてくれる。でも実用事例は少ない。コストに見合う価値はあったのだろうか。

　没入的VRの成功例として、恐怖症の治療、軍の訓練、エンタメが挙げられる。それらはすべて、仮想環境が提供する経験が現実的であることが鍵だった。言い換えると、それらはすべて、ユーザがある環境のなかにいると感じるということが重要であるような事例であった。逆にいうと、そういうのに注目していたからこそ実用事例が少ないままなのではないか。

　別の戦略について考えてみよう。没入性から得られるベネフィットは他にないのか。そもそも没入性というのは単一の構成概念なのか。
　たとえば、石油・ガス産業のために、地価の油井のパスを視覚化するシステムをつくるとしよう。この場合、別に「自分が地下にいる」と感じさせる必要はない。視覚化はあまり現実的では困る。むしろ適度に抽象的なほうがよい。Gruchallaらはデスクトップディスプレイと没入型プロジェクタを比較し、後者のほうがパフォーマンスが上がることを示している。たぶん空間理解が促進されたのだろう。では没入型プロジェクタのどんな要素が空間理解を促進したのか。視野の広さか、頭部運動のトラッキングか。それがわかれば、もっと簡単なシステムで十分だということになるかもしれない。

　上に挙げたような研究は数が少ない。たくさんの問いが残されている。同じようなベネフィットを実現してくれる課題や応用がほかにもあるのではないか。没入性のどのような要素が利点をもたらすのか。もっと低レベルな没入性でも十分なのではないか...
　複数の要因について効果を検証する実験が望ましい。最初にGruchallaのような実践的研究をやり、次に要因を分解していくというやり方もできる。

　[ここで著者らの実験紹介...]
　いくつかの実験を通じてわかったことを挙げる。

　まとめ。
　実践的にはふたつの目標がある。没入的VRを成功させるという目標と、高価なVRシステムなしで済ませるという目標である。このふたつは別に矛盾していない。没入性というものを連続的なものとして捉え、どんなレベルの没入性がなにを生み出すかを考えればよい。
　云々。

　本文よか、途中で挿入されている解説コラムのほうが役に立った。
　没入性(immersion)とはVRシステムが提供する感覚的fidelityの客観的レベルのこと。presenceとはVRシステムに対する主観的・心理学的な反応のこと。
　没入性の要素としては、同時にみることができる視野の広さ、ユーザを取り囲む視野の広さ、ディスプレイのサイズ、解像度、立体視、ヘッドトラッキングによるレンダリング、リアルな光線、フレームレート、リフレッシュレート、などがある。提示しているモデル自体のリアルさをここに含めていないのは、課題によって必要性が異なるから。
　この定義でいうと、相互作用のリアルさは没入性とは別の問題。むしろ、(入力装置による)相互作用と(ディスプレイによる)没入性がループしていると考えるべきだ。云々。

　VRの没入性がユーザの行動に与える効果の研究として、Meehan et al.(2002 Proc.ACM Siggraph)というのが挙げられていた。

論文：心理 - 読了：Bowman & McMahan (2007) ヴァーチャル・リアリティってのは視覚的没入性が高ければ高いほどよいというものでもないだろう

2017年9月14日 (木)

Zammit-Mangion, A., Cressie, N. (2017) Fixed Rank Kriging: The R Package.

　Rの空間統計パッケージFRKのvignette。適当にめくるだけのつもりが、面白くてついつい読み耽ってしまった。ま、最初の3頁だけだけど。
　FRKとはFixed Rank Kriging (固定ランク・クリギング)の略である。

　いわく。
　FRKと似たパッケージとして以下がある。[このくだり、まったく意味がわからない箇所もあるが、後日のために逐語訳する。精度行列というのは共分散行列の逆行列のことね]

　モデルの概要。
　FRKパッケージはふつうのバリオグラム・モデルではなく、空間ランダム効果モデル(SRE)で問題を定式化する。

　地点$\mathbf{s}$における量$Y(\mathbf{s})$に関心があるとしよう。古典的な空間モデルなら、共変量ベクトルを$\mathbf{t(s)}$、回帰係数ベクトルを$\mathbf{\alpha}$、空間的相関があるランダム効果を$\upsilon(\mathbf{s})$、空間的相関がないランダム効果を$\xi(\mathbf{s})$ととして
　$Y(s) = \mathbf{t(s)}^\mathrm{T} \mathbf{\alpha} + \upsilon(\mathbf{s}) + \xi(\mathbf{s})$
ただし$E(\upsilon(\cdot)) = E(\xi(\cdot)) = 0$、とするところである。

　ここで、
　$\upsilon(\mathbf{s}) = \sum_{l=1}^r \phi_l(\mathbf{s}) \eta_l$
とする。$\mathbf{\eta} \equiv (\eta_1, \ldots, \eta_r)^{\mathrm{T}}$はランダムベクトル、$\mathbf{\phi} \equiv (\phi_1(\ldots), \ldots, \phi_r(\ldots))$はあらかじめ定められた空間基底関数である。
　[理解するまでに手間取ったのだが、ここがこの話のミソなのである。この$\mathbf{\phi}$こそがバリオグラム・モデルの変わり果てたお姿なのだと思う。話の先取りになるけど、FRKパッケージにはデータとバリオグラム関数の形状を指定すると$\mathbf{\phi}$の行列を自動的に作ってくれる関数がある。それを使わずに自力で作ってもよい]

　さて、空間を$N$個の重ならない小さな空間に完全に分ける[メッシュみたいなものであろう]。これを基本地域単位 BAU と呼ぶ。BAUの数$N$は基底の数$r$よりずっと大きいものとする。
　それぞれのBAUにおいて$Y(\mathbf{s})$を平均する。つまり、$i$番目のBAU $A_i$において
　$Y_i \equiv \frac{1}{|A_i|} \int_{A_i} Y(\mathbf{s}) d(\mathbf{s})$
このレベルでも、
　$Y_i = \mathbf{t}_i^\mathrm{T} \mathbf{\alpha} + \upsilon_i + \xi_i$
と分解できる。
　以下、各項について順にみていくと...

　BAUレベルの共変量$\mathbf{t}_i$はこうなる。
　$\mathbf{t}_i \equiv \frac{1}{|A_i|} \int_{A_i} \mathbf{t}(\mathbf{s}) d \mathbf{s}$

　BAUレベルの相関ありランダム効果$\upsilon_i$は
　$\upsilon_i \equiv \frac{1}{|A_i|} \int_{A_i} \upsilon(\mathbf{s}) d \mathbf{s}$
ここに$\upsilon(\mathbf{s}) = \mathbf{\phi}(\mathbf{s}) \mathbf{\eta}$を代入すると... [面倒くさいから書かないけど] 係数ベクトルと$\mathbf{\eta}$の積になる。係数の部分を取り出して$N \times r$行列$\mathbf{S}$とすれば、$\mathbf{\upsilon} = \mathbf{S\eta}$となる。なおこのパッケージでは、実際にはそれぞれの係数について積分するのではなく、BAUは十分に小さいとみて、BAUの重心$\mathbf{s}_i$を使い$\mathbf{\phi}(\mathbf{s}_i)$と近似する。
　$\mathbf{\eta}$について。このパッケージでは、$\mathbf{\eta}$を平均0, 共分散行列$\mathbf{K}$の正規ベクトルとする。$\mathbf{K}$は自由推定してもいいし、なんらかの構造を与えても良い。前者をFRK-Vと呼び(Vはバニラの略)、後者をFRK-Mと呼ぶ(Mはモデルの略)。

　相関なしのランダム効果について。$\xi(\mathbf{s})$についてはもう関心を持たない。BAUのレベルでの平均$\xi_i$について、平均0, 分散$\sigma^2_\xi \nu_{\xi i}$の正規分布に独立に従うとする。ここで$\sigma^2_\xi$は全BAUを通したパラメータ、$\nu_{\xi i}$は既知の定数でBAUの異質性を表す。

　ここまでをまとめよう。BAUレベルの式
　$Y_i = \mathbf{t}_i^\mathrm{T} \mathbf{\alpha} + \upsilon_i + \xi_i$
を書き直して
　$\mathbf{Y} = \mathbf{T} \mathbf{\alpha} + \mathbf{S} \mathbf{\eta} + \mathbf{\xi}$
$\mathbf{T}$は行数$N$の共変量行列、$\mathbf{\alpha}$は係数ベクトル。$\mathbf{S}$は$N \times r$の行列、$\mathbf{\eta}$は長さ$r$のランダム効果ベクトルでその共分散行列が$\mathbf{K}$。$\mathbf{\xi}$は長さ$N$のランダム効果ベクトルでその共分散行列は$\sigma^2_\xi \mathbf{V}_\xi$、ただし$\mathbf{V}_\xi$は既知。

　さて、ここまで考えてきた$Y(\cdot)$というのは潜在過程であって、観察自体は$m$個のフットプリントについてのみ可能であるとしよう。ここでフットプリントとは、ひとつ以上のBAUからなる地域のことで、重複していてもかまわない。$m$は$r$よりずっと大きいものとする。なお、$m$はBAUの数$N$より大きくても小さくてもよい。
　フットプリント$B_j$における観察値は次の3つの和とする[面倒くさいので式は省略]：

　結局、推定するパラメータは、バニラFRKでは$\mathbf{\alpha}, \sigma^2_\xi, \sigma^2_\delta$, そして$\mathbf{K}$ の4つ。モデルFRKでは、$\mathbf{K}$の代わりに$\mathbf{K}$について組んだモデルのパラメータがはいることになる。

　[以下、力尽きたので、ほとんど読んでいない]
　パラメータはEMアルゴリズムで推定する...
　このモデルにより、任意の予測領域についての予測が可能で...

　コード付きの事例が2つ...
　時空間でクリギングする事例が2つ...
　
　空間異方性を持たせるには...
　既定関数とBAUを手動で与えるには...

　今後の課題:

云々、云々。

　...というわけで、正直なところ頑張って読んだのは最初の3頁だけなんだけど、固定ランク・クリギングという発想がよくわかった(ような気がする)。$m$個の観察値の後ろに$N$個のメッシュの値を考え、そのまた後ろにもっと少ない$r$個の隠れた値を考えるわけね。時系列モデルで言うと、観察変数のうしろに状態変数を考え、そのまた後ろにたまにしか動かないような状態変数を考えているわけだ。
　だから大データでもクリギングが楽だ、ってわけね。なるほどねー、面白いなあ。(←すっかりわかったような気分)

論文：データ解析(2015-) - 読了：Zammit-Mangion & Cressie (2017) 大きな空間データを固定ランク・クリギングするRパッケージ FRK

Rabiee, F. (2004) Focus-group interview and data analysis. Proceedings of the Nutrition Society, 63, 655-660.
　題名の通り、グループ・インタビューの逐語録をどうやって分析するかという短い啓蒙的論文。きちんと読んでいないけど、整理の都合上読了にしておく。
　紹介されているコーディング法の元ネタはKrueger & Casey (2004) "Focus Groups: A Practical Guide for Applied Research", 3rd. ed. だそうだ。ふーん。

論文：調査方法論 - 読了：Rabiee, F. (2004) グループ・インタビューの逐語録をどうやって分析するか

　調べ物をしていて偶然知ったんだけど、厚生労働省が発表している市区町村別生命表にはベイズ推定の考え方が入っているのだそうだ。これ、いわゆる小地域推定という奴で、ほんとに市区町村別に標本統計量をとっちゃうと、さすがに死亡者数が足りない、という事情であろう。
　ううむ。別に生命表を作る用事はないけれど、私の仕事とも通じる、ちょっと切実な話だ...

府川哲夫, 清水時彦 (1990) 小地域生命表のベイジアン・アプローチ. 人口学研究, 13, 37-49.

　というわけで、なにをどうやっているのか調べてみた。著者らの肩書は厚生省大臣官房統計情報部となっている。いまでもこの論文の方法で市区町村別生命表を作っているのかどうかわからないけど、仕事の参考になるかなと思って。
　
　市区町村x性x年齢階級別に考える。人口$p$は既知。死亡数$d$を$Bin(p, \theta)$の実現値とみる。$\theta$の事前分布を$Beta(\alpha,\beta)$とする。[...途中を端折って...] $\theta$の事後分布は$Beta(\alpha+d, \beta+p-d)$となりますわね。

　さて、問題は事前分布なんだけど...
　以下、使うデータは昭和62年まで5年間の人口動態統計の死亡数、そして昭和60年国勢調査の総人口である。人口はあらかじめ5倍する。
　都道府県を市部と郡部に二分する。ある市区町村の親地域[←そう書いてはいないが簡略のためにこうメモする]は、市区ならば当該都道府県の市部、町村ならば当該都道府県の郡部とする。
　注目している市区町村を$i$、その親地域を$K$とする。ある性x年代について、$i$の人口を$p(i)$、租中央死亡率を$q(i)$とする。
　$K$に属する子地域の租中央死亡率の平均$Q$と分散$V$を求める。重みを$w(i)=p(i)/\sum p(i)$として
　$Q = \sum w(i) q(i)$
　$V = \sum w(i) (q(i) - Q)^2$
　で、$i$の事前分布パラメータ$\alpha, \beta$を、ベータ分布の平均と分散が$Q, V$になるように決めちゃうのである。[←うおおおおお... なんというか、実に素朴なアプローチだ...]

　著者らいわく、要するにこういうことだ。生命表をつくるとき、中央死亡率を$d/p$とするのが伝統的統計学。$(\alpha+d)/(\alpha+\beta+p)$とするのがベイズ統計学。それだけの違いだ。[←いやいやいや、それは母比率の事前分布をベータ分布とみればそうなるでしょうけど... 問題は事前分布のパラメータ$\alpha, \beta$をどう決めるかでしょう...]

　後半は推定された市区町村別生命表の観察。パス。

　...いやあ、ものすごーくシンプルな経験ベイズ・アプローチというか... 正直、わたくし、びびりました。まじですか。そんな簡単な話でいいんすか。
　著者らも最後に触れているけど、親地域の決め方はこれでいいのか。人の生存曲線ってのはどの都道府県の市部/群部かで決まるわけ？ それってあまりに行政区分に依存していないか...？
　いや、しかし、公的統計というのはこのくらいシンプルでないといかんのかもしれない。ついでにいうと、私の仕事もこのくらいカンタンに考えたほうがいいのかもしれない...少なくともそういう局面はありうる...
　などなど、終電車の酔っ払いたちに揉まれながら論文に目を通し、あれこれ考えさせられました。

論文：データ解析(2015-) - 読了：府川, 清水 (1990) 厚労省による市区町村別生命表はベイズ推定されている

2017年9月13日 (水)

矢島美寛, 平野敏弘(2012) 時空間大規模データに対する統計的解析法. 統計数理, 60(1), 57-71.
　正直なところ、内容の9割くらいはまっっっったく理解できないんだけど、しかし今後頑張れば理解できるようになる日が来るという問題でもないと思うので、読了にしておく。

　ま、空間統計モデル関連でときどき見かける謎のタームについて若干の知識を得ることができたので、よしとしよう。

論文：データ解析(2015-) - 読了：矢島・平野(2012) 大規模空間データの統計モデリング

Graler, B., Pedesma, E., Heuvelink, G. (2016) Spatio-Temporal Interpolation using gstat. The R Journal, 8(1). 204-218.
　えーと、Rにはクリギング(空間データ補間手法の一種)のためのパッケージがいっぱいあって、なにがなにやらわからないので、ひとまず古株らしきgstatパッケージの紹介を読んでみた次第(vignetteの一つに挙げられている)。R Journalとはいえ、すごく小さな字で14頁...

　えーと、gstatは時空間統計のための機能を提供しているのが売り。同様のパッケージにRandomFieldsがある。ほかにもあるので、CRAN Task Viewの該当項目をみるように。

　以下では、空間領域を$S$, 時間領域を$T$とする。ガウシアン時空間確率場$Z$について考える。観察を
　$z = (z(s_1, t_1), ..., z(s_n, t_n))$
とする(同一の地点での反復測定や異なる地点での同時測定もありうるとする)。$z$に基づき$Z$のモデルをつくりたい、というのがお題。

　$Z$は定常で、空間的に等方性があると仮定する。つまり、この確率場は平均$\mu$と共分散関数$C_{st}$で特徴づけられ、共分散は空間距離$h$と時間距離$u$のみに依存すると考える($h$と$u$は正の実数)。なおこの仮定を取っ払って拡張するのは容易で、たとえば普遍クリギングは時間なしのケースでの非定常性への拡張である。
　以下では次のように書く。
　$Cov_{st}(h,u)=Cov(Z(s,t), Z(\tilde{s}, \tilde{t}))$
ただし$h=||s - \tilde{s}||, u = |t - \tilde{t}|$とする。
　時空間バリオグラムを
　$\gamma_{st} = C_{st}(0,0) - C_{st}(h,u)$
と書く。
　共分散関数が推定できたとしても、共分散行列から線形予測子を求めるためには逆行列計算が必要なのだが、大変なので、いろいろ工夫がある(時間的マルコフ構造を使うとか)。spTimerパッケージ、spBayesパッケージ, spateパッケージ, INLAパッケージなどがある。

　...というわけで、ここからが本題なんだけど、ざっと読み進めたら、空間領域のバリオグラムモデルと時間領域のバリオグラムモデルがあるとして(異なるファミリーでも構わない)、時空間バリオグラムモデル$\gamma_{st}$をどう定義するか(つまり共分散関数$C_{st}$を定義するか)、という説明であった。パッケージ自体の紹介じゃなかった。読むものを間違えたみたい。

　めんどくさくなったので途中から流し読み、メモは省略。でも、私にも理解できる範囲に関していえば面白い内容であった。当面は時空間モデルを作る用事がないけれど、いつか役に立つこともありそうだ。

　時空間モデルともなると、バリオグラムは曲面となる(時間距離と空間距離がX, Yとなる)。ひえええ。それでも著者らいわく、やっぱしバリオグラムモデルのフィッティングは、MSEとかだけじゃなくて、ちゃんと視覚化してチェックしないといけないそうだ。大変だなあ。

論文：データ解析(2015-) - 読了：Graler, Pedesma, Heuvelink (2016) Rのgstatパッケージで時空間クリギング

　Rの辛いところはパッケージがいっぱいありすぎるところだ。このたび空間統計のパッケージについて調べていてほとほと困惑したので、空間データ分析のCran Task View (2017/05/09)の、Geostatisticsの項についてメモしておく。100%自分用の覚え書きです。
　[2018/06/18追記: 2018/06/02版のTask Viewに従って更新した]

あーもう。これがSAS/STATなら, VARIOGRAM, KRIGE2D, SIM2D, SPPの4択なのに。勘弁してよ、もう！

雑記：データ解析 - Rの空間統計パッケージ

佐藤哲也 (2017) AIと政治. 人工知能, 32(5), 672-677.

　最新号の「人工知能」誌は「AI社会論」という特集で、佐藤先生が寄稿されていた。忘れないうちに内容をメモ。良いまとめとは思えないので、関心ある方は原文をお読みくださいますように。

　いわく、
　汎用人工知能というような話は横に置いておいて、現行技術が政治現象にもたらす影響について考えると、

ところで、昨今の人工知能ブームにはテクノロジー・プロパガンダという面がある。バイオや製薬では、産業界のアカデミズムへの不正な介入が社会問題になっているが、構造が似てきていないだろうか。結局割りを食うのは一般の納税者や投資家だ。人文社会的観点からの検討が必要であろう。

論文：予測市場 - 読了：佐藤(2017) 人工知能と政治

Roberts, M.E., Steward, B.M., Tingley, D. (in preparation) stm: R package for structural topic models.

　構造的トピックモデルのRパッケージ stm のvignette。いつ書かれたものかもよくわからないのだが、J. Statistical Softwareの判型になっているから、投稿中かなにかなのかも。

　いわく。
　構造的トピックモデル(STM)とは、トピックモデルに文書に付与されたメタデータを統合したモデル。トピックを抽出すると同時に、トピックとメタデータとの関連性を推定する。

　モデルの説明。
　通常のトピックモデルと同じく、語のカウントの生成モデルを考える。なお、メタデータがなかったら、STMはBleiらの相関トピックモデルと同じである。
　文書を$D_1, \ldots$, 単語を$w_1, \ldots$, トピックを$T_1, \ldots, T_k$とする。文書$d$のメタデータを$X_d$とする。メタデータは次の2種類に分けて指定できる。

　生成モデルは次の通り。

　推定はsemi EMアルゴリズムを用いた変分的推定である、とのこと。[よくわからんが変分ベイズ推定だということかしらん]

　stmパッケージの使い方の説明。
　初期値決定はスペクトル分解による決定論的方法がお勧めだが、まず崩壊型ギブスサンプラーでLDAモデルを推定するとか、ランダムに決めるというオプションもある。
　ちゃんと読んでないので省略するか、データ準備からモデル評価、トピック数決定、結果の視覚化まで、これでもかっていうくらいにいろんなヘルパー関数がある。
　云々, 云々...

　他のパッケージとの比較。
　相関トピックモデルはtopicmodelsパッケージでも推定できるけど、こっちのほうがパフォーマンスが良い。云々, 云々...
　
　... とかなんとか。途中からはパラパラめくっただけだけど、いずれに必要になったらきちんと読もう。

論文：データ解析(2015-) - 読了: Roberts, Steward, Tingley (in preparation) 構造的トピックモデルのRパッケージ stm

2017年9月12日 (火)

　6月末頃に作ったメモ。トピックモデルの表記が資料によって違っているせいで、混乱しちゃったのである。

岩田(2015)「トピックモデル」の表記

　文書集合$\{\mathbf{w}_1, \ldots, \mathbf{w}_D\}$について考える。文書$\mathbf{w}_d$に含まれる単語を$w_{d1}, \ldots, w_{dN_d}$とする(順序は問わない)。文書集合を通して現れる語彙の数を$V$とする。

　トピックモデルによれば、文書は次のように生成される。文書集合の背後には$K$個のトピックがある。

  1. トピック$k=1,\ldots,K$のそれぞれについて、単語分布$\mathbf{\phi}_k = (\phi_{k1}, \ldots, \phi_{kV})$が生成される。
    　$\mathbf{\phi}_k \sim Dirichlet(\beta)$
  2. 文書$d=1,\ldots,D$のそれぞれについて...
    1. トピック分布 $\mathbf{\theta}_d = (\theta_{d1}, \ldots, \theta_{dK})$が生成される。
      　$\mathbf{\theta}_d \sim Dirichlet(\alpha)$
      トピックが生成されているわけではない点に注意。文書はあるトピックを持っているのではない(トピックの確率分布を持っている)。逆向きにいうと、トピックモデルは文書にトピックを割り当てない。
    2. 単語 $n=1, \ldots, N_d$ について...
      1. まずはトピック $z_{dn}$が生成される。単語が生成される前に、単語の数だけトピックが生成されている。逆向きにいうと、単語にトピックを割り当てていることになる。
        　$z_{dn} \sim Categorical(\mathbf{\theta}_d)$
      2. いよいよ単語が生成される。あらかじめ作っておいた単語分布$\mathbf{\phi}_{z_{dn}}$を参照して、
        　$w_{dn} \sim Categorical(\mathbf{\phi}_{z_{dn}})$

佐藤(2015)「トピックモデルによる統計的潜在意味解析」の表記

　トピック$k$における単語の出現分布を$\mathbf{\phi}_k = (\phi_{k,1}, \ldots, \phi_{k,V})$とし、
　$\mathbf{\phi}_k \sim Dir(\beta)$
　文書$d$のトピック分布(トピックの構成比率)を$\mathbf{\theta}_d = (\theta_{d,1}, \ldots, \theta_{d,K})$とし、
　$\mathbf{\theta}_d \sim Dir(\alpha)$
　文書$d$における$i$番目の潜在トピック$z_{d,i}$について、
　$z_{d,i} \sim Multi(\mathbf{\theta}_d)$
　$i$番目の単語$w_{d,i}$について
　$w_{d,i} \sim Multi(\phi_{z_{d, i}})$

Griffiths & Steyvers (2004, PNAS)の表記

細かいところを端折ると、
　$w_i | z_i, \phi^{(z_i)} \sim Discrete(\phi^{(z_i)}) $
　$\phi \sim Dirichlet(\beta)$
　$z_i | \theta^{(d_i)} \sim Discrete(\theta^{(d_i)})$
　$\theta \sim Dirichlet(\alpha)$
こういう順番で書かれるとめっさわかりにくいが、まあとにかく、ハイパーパラメータは$\alpha, \beta$である。

　ふつうそうですよね。良識ある大人なら、たいていの人がハイパーパラメータを$\alpha, \beta$と書いていたら、その世間の風潮にあわせますよね。俺様だけは$\beta$を$\delta$と書くぜ、読者どもよついてこい、なんて子供じみた真似はしませんよね？

Grun & Hornik (2011, J.Stat.Software) の表記
　...というか、Rのtopicmodelsパッケージの表記。
　さあ、ここで我々は、世間の風潮など気にもしない漢たちを目の当たりにし、驚愕することになるのである。はい深呼吸！

　文書$w$について考える。文書$w$に含まれる単語を$w_{1}, \ldots, w_{N}$とする(順序は問わない)。
　文書は次のように生成される。文書集合の背後には$k$個のトピックがある。

  1. それぞれのトピックについて、単語分布$\beta$が生成される。
    　$\beta \sim Dirichlet(\delta)$
  2. 文書$w$について...
    1. トピック分布 $\theta$が生成される。
      　$\theta \sim Dirichlet(\alpha)$
    2. 単語 $w_i$ について...
      1. まずはトピック $z_i$が生成される。
        　$z_{i} \sim Multinomial(\theta)$
      2. 多項確率分布 $p(w_i | z_i, \beta)$から単語が生成される。

　実のところ、こうやってメモしててやっと謎が解けたんだけど、topicmodelsパッケージで LDA(x, k, control = list(estimate_beta = FALSE), model=...) と指定した場合、固定される「beta」とは、トピック別単語分布のハイパーパラメータ$\beta$ではなく、トピック別単語分布そのものなのだ。ハイパーパラメータはdeltaと呼ぶのである。なんだかなあ、もう...
　なお、ハイパーパラメータのほうは、Gibbsサンプリングの場合のみ control=list(delta=...)として指定できる由。

雑記：データ解析 - トピックモデルの表記比較 (いつも思うんですが、専門家のみなさん記号を統一していただけませんかね)

　仕事の都合で仕方なく、トピックモデルについてあれこれ考えているんだけど、トピック数の決定ってのはいったいどうしたらいいのかしらん、と不思議である。
　解説本などをみると、それは階層ディリクレ過程で推定できる、それにはまず中華料理店のフランチャイズ・チェーンの仕組みについて検討する必要がある、ただし店舗内のテーブル数は無限とする... なんて途方に暮れるようなことが書いてあるんだけど(本当)、そんなおそろしい話じゃなくてですね、探索的因子分析のスクリー基準みたいなしょっぼーい話でいいから、なんか目安がほしいだけなんですけど...

　やっぱり、データをちょっとホールドアウトしておいて、perplexityが低くなるトピック数を選べって話だろうか。面倒くさいなあ、と思いながらぼんやり検索していたら、なんと！Rにはldatuningというパッケージがあり、topicmodelsパッケージでLDA(潜在ディリクレ配分モデル)を組む際の最適トピック数を教えてくれるのだそうである。すごい！偉い！
　
　さっそくldatuningを試してみたところ、さまざまなトピック数について、CaoJuan2009, Arun2010, Daveaud2014, Griffiths2004という4つの指標が出力される。それだけのパッケージであった。チャートを描いて考えな、という話である。
　なんだ、この謎の指標は？ マニュアルには一切説明がなく、ただ出典が書いてあるだけ。
　さすがに気持ち悪いので、出典をあつめてめくってみた。

CaoJuan2009

　出典はJuan, Tian, Jintao, Yongdon, Sheng (2009, Conf.)。
　パラパラめくってみただけであきらめたんだけど(第一著者は博士の院生さんで80年生まれとのこと。うわー)、どうやら、2トピック間の語彙を通じたコサイン類似性を求め、その総当たり平均が最小になるようなトピック数を求める... というような話らしい。LDAではトピックは独立なはずだから、ということみたいだ。へぇー。
　ldatuningのソースを眺めるとこうなっている。トピック$k(=1,\ldots,K)$の単語事後分布を格納したベクトルを$x_k$として、
　$dist_{k,j} = \sum x_k x_j / \sqrt{\sum x_k^2 * \sum x_j^2}$
と定義し、すべての2トピックの組み合わせについてこの距離を合計し, $K(K-1)/2$で割っている。

Deveaud2014

　出典はDeveaud, SanJuan, Bellot (2014), 論文だか紀要だかよくわからない。
　もちろんきちんと読もうなどという大それた野心は持ち合わせていないわけで、トピック数決定についてのページをディスプレイ上で眺めただけなのだが、どうやら、2トピック間での単語出現分布のKLダイバージェンスを求め、その総当たり平均が最大になるようなトピック数を求める... という話らしい。へぇー。
　ldatuningのソースを眺めると、
　$dist_{k,j} = 0.5 \sum x_k \log(x_k/x_j) + 0.5 \sum x_j \log(x_j/x_k)$
と定義し、すべての2トピックの組み合わせについてこの距離を合計し, $K(K-1)$で割っている。

Arun2010

　出典はArun, Suresh, Madhaven, Murthy(2010, Conf.)。入手できなかった。上のDeveaud et al.(2014)によれば、これもトピック間の距離みたいなものを最大にすることを目指すのだそうだ。
　ソースコードを眺めたのだが、トピック別の単語事後分布とそのハイパーパラメータを使って、なにかごにょごにょしている模様。

Griffiths2004

　出典として、Griffiths & Steyvers (2004, PNAS) と、誰かの修論かなにかが挙げられている。前者は科学研究のトピックを推定するという話で、有名な論文であろう。
　ちゃんと読んでないけど、パラパラめくったところによればこういうことらしい。
　要するに、トピック数を$T$、コーパス全体の単語を$\mathbf{w}$として、対数尤度$\log P(\mathbf{w}|T)$が最大になる$T$にすればいいじゃん、というアイデアである。しかしその直接の算出は、単語をトピック$\mathbf{z}$に割り当てるすべてのパターンを通じた合計が必要になるので無理。ところが、次のやりかたでうまく近似できる。
　まず事後分布$P(\mathbf{z} | \mathbf{w}, T)$から$\mathbf{z}$をドローし、$P(\mathbf{w}|\mathbf{z}, T)$をたくさん求める。で、その調和平均を求める。へぇー。
　なおGrun& Hornik(2011, J.Stat.Software)によれば、このやり方には推定量の分散が無限大になるかもしれないという欠点があるのだそうだ。
　ソースコードのほうは、little tricky というコメントがついており、よくわからないことをやっている。

　こうしてみると、Griffith & Steyvers (2004) はモデルの周辺尤度に注目する路線、他の3つはトピック間の単語事後分布の類似性に注目する路線といえそうだ。前者はperplexityみたいなもの、後者はcoherenceみたいなもの、ということなのかしらん？
　探索的因子分析のアナロジーで言うと、前者はどうにかしてモデル適合度を求めようという路線、後者は回転後負荷行列の性質を調べる路線(たとえば単純構造が得られたかどうかとか)、という感じかしらん？
　えーと、自分で書いててよくわかんなくなってきたので、このへんで。

雑記：データ解析 - 潜在ディリクレ配分のトピック数を手軽に決める方法 (Rのldatuningパッケージが提供する謎の4つの指標について)

2017年9月11日 (月)

秋山祐樹・仙石裕明・柴崎亮介 (2013) 全国の商業集積統計とその利用環境. GIS－理論と応用, 21(2), 11-20.
　著者の先生方が開発してきた(そしてゼンリンが販売している)商業集積統計のつくりかたについて述べた論文。仕事の都合で読んだ。

　いくつかメモ：

　それにしても、まあ仕事で必要な勉強ならばなんでもやりますけど、いやーずいぶん遠くにきちゃったなあ、という感慨がある。トシ食ってからこんなふうに思ってもみないことを勉強する羽目になるくらいなら、若い頃に臨床心理士の勉強をしておけばよかったなあ、なあんて(すいません冗談です。超・冗・談・で・す)

論文：マーケティング - 読了：秋山・仙石・柴崎(2013) 商業集積統計のつくりかた

Ribeiro, P.J., Christensen, O.F., Diggle, P.J. (2003) geoR and geoRglm:
Software for Model-Based Geostatistics. Proceedings of the 3rd International Workshop on Distributed Statistical Computing (DSC 2003).
　Rの空間統計パッケージ geoR と geoRglm の紹介。実戦投入前の儀式として読んだ。うっかりdraftのほうを読んじゃったんだけど、中身はだいたい同じだろう、と思うことにする。

　えーっと、geodataというクラスを用意しています。簡単にプロットできるし、簡単に経験バリオグラムが描けます。
　空間統計モデルとして、geoRではガウシアンモデル、ならびにbox-cox変換した従属変数についてのガウシアンモデルが組めます。geoRglmでは、box-cox変換した従属変数についてのポワソンモデル、ならびにロジットリンクの二項モデルが組めます。パラメトリックなバリオグラム・モデルを使います(デフォルトはMeternモデル)。推定量としてMLとREMLをご用意してます。いくつかのクリギング手法をご用意してます。あっ、そうそう、ベイジアン・クリギングもできますよ、自前でMCMCの機能を持っています。
　云々、云々。

論文：データ解析(2015-) - 読了：Ribeiro, Christensen, Diggle (2003) Rの空間統計パッケージ geoR と geoRglm

Morris, D.S. (2017) A modeling approach for administrative record enumeration in the Decennial Census. Public Opinion Quarterly, 81, 357-384.

　先日、世論調査研究のトップ誌であるPublic Opinion Quarterlyが「サーベイ調査のこれから」という特集号を出した。目次を眺めていて、あれれ、これってひょっとして仕事に関係があるんじゃなかろうか、大変な鉱脈を見つけちゃったかも... と思って手に取った論文。半分くらい目を通したところで、鉱脈どころか私の仕事とはまるきり関係ないとわかったんだけど、気分転換にはなった。著者は米国勢調査局の中の人。

　いわく。
　米国勢調査(Decennial Census)では、回答がなかった世帯に対して追加調査するんだけど(nonresponse follow-up; NRFU)、2010年国勢調査の場合にはNRFUだけで16億ドルかかった。2020年調査ではそのコストをなんとか減らしたい。
　いっぽう世論調査のほかにも、納税記録とか民間企業の訪問調査みたいな世帯データがある(こういうのをadministrative record, ARと呼ぶ)。ARと国勢調査を併用できないか。こういう発想は実は珍しくなくて、デンマーク、オランダ、スイス、ドイツ、ポーランドなどで採用されている。
　米国勢調査でも80年代から発想はあった。問題は、上記の国々とは異なり、米にはARを国レベルで統一的に管理する仕組みがなかったという点だ。現在、国勢調査局にAR研究活用センター(CARRA)というのがあって、各所からARを取りまとめている。細かいことを言うといろいろと大変なんだけど[詳細略]、まあとにかく、集めたARには統一的な個人識別IDと住所IDを振っている。
　もし、ARデータ側に含まれている世帯についてはNRFUの実査対象から除外することができれば、コストが節約できるではないか。そこで、ARデータと2010年国勢調査と比較し、どの世帯についてはARを使いどの世帯についてはNRFUをやるかを決める方法を開発したい。

　使用するARデータは、IRS 1040 [よくわからんが確定申告みたいなものだろうか？]、IRS informatonal returns [所得申告みたいなもの？]、メディケア、IHS (Indian Health Service)。個人と住所の組み合わせをキーとする。ほかに商用データのTargus Federal Consumer Fileというのがあって、これはAR側名簿としては使わないが、後述するモデル構築の際に予測子として用いた由。
　ある個人$i$と住所$h$の組み合わせが、ARにも2010年国勢調査にも存在していたら$y_{ih}=1$、そうでなければ0とする。で、$p_{ih}=P(y_{ih}=1)$を予測するモデルを作る。[きちんと読んでないので自信がないんだけど、この確率が1に近い住所はNRFUのリストから抜いてよかろう、という話だと思う。つまり、2010年国勢調査は完璧、AR側の紐づけも完璧、という前提での研究なのであろう]
　説明変数として、ARデータのうち「IRS 1040に存在」フラグとか、「IRS 1040に個人のみ存在」フラグとか、そういう変数をいろいろ作る。
　最終的に推定しないといけないのは住所についての確率なので、
　$\hat{p}_h = min(\hat{p}_{1h}, \ldots, \hat{p}_{n_h h})$
として、この推定値が閾値$c$を超えたらNRFUの実査対象からは外してAR側の記録を使うことにする。$c$はfalse positiveとfalse negativeの二乗和が最小になる値とする。

　... さあ予測モデルをつくりましょう、というわけで、ロジスティック回帰、分類木、ランダム・フォレストのモデルを作って比べたり、コストと正確性のトレードオフ曲線を推定したりしたらしいのだが、この辺で力尽きて読むのをやめた。
　ま、どういう問題なのかがわかったから、これでいいや。考えたこともないような話題で面白かった。

論文：調査方法論 - 読了：Morris (2017) 行政記録にデータがある世帯については国勢調査に無回答でもまあいいやということにできないか

Sirken, M.G., et al (1999) "Cognition and Survey Research"という論文集がある。80年代初頭から、米の調査法研究者らのあいだでCASM(Cognitive Aspects of Survey Methods)と呼ばれるアプローチが提唱されていて、この本はその1997年のセミナーをもとにした論文集。
　この本、原稿ならびに仕事の都合で救いを求めて手に取り、ぱらぱらめくって諦めることが多い。なぜいまこんなの読まなきゃいけないんだという気もするけど、どの章が役に立ちそうか、目次をメモしておく。

雑記 - Sirken, et al. (1999) 「認知と調査法研究」目次メモ

2017年9月 9日 (土)

Carlin, B.P., & Louis, T. (2000) Empirical bayes: Past, Present and Future. Journal of the American Statistical Association, 95(452), 1286-1289.

　以前、人と話していて、ふと「経験ベイズ」って正確にはどういう意味なんだろうか...と不思議になった。何をもって経験ベイズと呼ぶのか。それは手法か哲学か。
　ろくに知識もないのにあれこれ考えてもしょうがないので、検索で引っかかったやつをざっと読んでみた。この雑誌のこの号では"Vignettes for the Year 2000"と題し、「これからの○○はどうなるか」的な短い解説を22個のトピックについて載せた模様で、これは「経験ベイズ」についての寄稿。他のトピックは「ベイズ統計」「ブートストラップ」「変数選択」といった感じ。

　いわく。
　経験ベイズ(EB)という言葉はあいまいで、モデルのクラスを指すこともあれば、分析のスタイルを指すこともあれば、統計的手続きの哲学を表すこともあるんだけど、どんな定義であれ、まずはフル・ベイジアンとの違いから始めるという点では共通している。
　観察データ$y=(y_1, \ldots, y_n)$の分布を、未知パラメータ$\theta=(\theta_1, \ldots, \theta_k)$によって$f(y|\theta)$とモデル化したとしよう。頻度主義の立場からみると$\theta$は固定されているが、ベイジアンは事前分布$\pi(\theta|\eta)$を想定する。$\eta$が既知ならばベイズ・ルールを適用して、事後分布$p(\eta|y, \eta)$が得られる。
　さて、$\eta$が未知である場合、$\eta$についての情報は$y$の周辺分布$m(y|\eta)$によって捉えられている。さらに、もし$f$と$\pi$が共役なら(すなわち、もし$p(\eta|y, \eta)$が$\pi$と同じ分布族に属していたら)、$m(y|\eta)$は閉形式で得られる。
　フル・ベイジアン(またの名をベイズ経験ベイズ, 略してBEB)ならば、ハイパー事前分布$h(\eta|\lambda)$を想定するところである。このとき、事後分布は結局
　$p(\theta|y, \lambda) = \int p(\theta|y, \eta) h(\eta|y, \lambda) d\eta$
となる。事後分布とは、$\eta$を固定したときの事後分布と、データで更新したハイパー事前分布の混合となるわけである。
　いっぽう経験ベイズでは、$\eta$をデータから推定する(たとえば周辺最尤推定量で)。で、事後分布として$p(\theta | y, \hat{\eta})$を使う。
　BEBもEBも、$\eta$についての情報をデータから得ているという点では同じことである。ちがいは一番上の分布$h$を含めているかどうかだ。BEBの美点は$\eta$にまつわる不確実性をうまく組み込んでいるという点で、欠点は$h$の選択が難しいという点である。

　EBによる$p(\theta | y, \hat{\eta})$は、$\eta$の不確実性を考慮していない分、狭くなる。これをどうにかしようというのが、EBの世界の長年の課題であった。
　一つのアプローチはパラメトリックEBである。このアプローチでは、下から2番目の分布$\pi(\theta|\eta)$をパラメトリックな形式で与え、あとは$\eta$さえが決まれば事後分布が完全に決まるようにする。Morris(1983 JASA), Casella(1985 Am.Stat.)を参照のこと。
　もうひとつのアプローチはノンパラメトリックEB。最後から2番目の分布を単に$\pi(\theta)$とする。この路線はRobbins(1955)に始まる。変種としてノンパラメトリック最尤法というのもある。

　皮肉なことに、EBの歴史はそれほど「ベイジアン」ではない。[頻度主義の観点からみて良い性質を持つ決定ルールが目指されていたという話... 略]
　パラメトリックEBはStein推定と強い関係があって...[略]
　初期EBは、まず事前分布を決めるためにデータを使い、事後分布を求めるためにもう一度データを使った。この方式は当時のベイジアンたち(多くは主観主義的ベイジアン)に忌避された。SavageはEBを指して、ベイジアンの卵を割っておきながらベイジアンのオムレツを拒否する奴らだと呼んだ。でもEBは、その後の客観主義的ベイジアンを準備し、その発想はのちにGibbsサンプラーの登場によって花開くこととなる。
　
　EBは現在普及しており...[略]
　メタ分析とも関係があって...[略]
　まだまだ進化しております、たとえば...[略]

　EBの未来はどうなるか。MCMCがこんだけ普及すると、近似としてのEBはもはや出番がないという悲観的見方もある。でも、MCMCは「プラグ・アンド・プレイ」とは言い難い。収束判定は厄介だし、でかいモデルをむやみに作りそうになっちゃうし。EBの出番はまだまだあります。
　たとえばですね。分散要素$\tau^2$についてあいまいなハイパー事前分布を決めるという問題について考えよう。いまポピュラーなのは$gamma(\epsilon, \epsilon)$で、つまり平均1、分散$1/\epsilon$だ。しかし最近の研究では、これは無情報に見えて実は事後分布に大きな影響を与えることが示されておる。さらに、この分布は$\epsilon$を小さくするとimproperに近づき、MCMCの収束が難しくなる。いっそ$\tau^2$を$\hat{\tau}^2$に置き換えちゃったほうが安全ではないですかね。まあ純粋なベイジアンじゃなくなっちゃうけど。
　EBもBEBも良いやり方だし、どちらも万能薬ではないです。みなさまそれぞれ哲学的傾向をお持ちだと思いますが、手法選択の際にそれに殉じることはないんじゃないでしょうか。
　云々。

　...パラメトリック/ノンパラメトリック経験ベイズについての説明がよく理解できなかった... 歴史に属する話なのかもしれないけど、なんだか悔しい。なにか別のを読んだほうがよさそうだ。

論文：データ解析(2015-) - 読了：Carlin & Louis (2000) 経験ベイズのこれまでとこれから

引き続き、溜まったメモの整理。これは6月中頃。

Linzer, D.A., Lewis, J.B. (2011) poLCA: An R Package for Polytomous Variable Latent Class Analysis. Journal of Statistical Software, 42(10).
　RのpoLCAパッケージの紹介。poってのはpolytomousの略であろう。

　えーと、poLCAパッケージは多項変数の潜在クラス分析(LCA)をご提供します。
　潜在クラスモデルを推定する関数としては、ほかにe1071パッケージのlca関数, gllmパッケージ, randomLCAパッケージがあるが、いずれも二値変数しか扱えない。

　まずは基本的な潜在クラスモデルから。
　対象者を$i=1, \ldots, N$とする。顕在変数は$J$個の多項カテゴリカル変数、$j$個めの変数がとりうる値を$K_j$個とする。顕在変数の値を二値変数$Y_{ijk}$で表す。つまり、対象者$i$が変数$j$で反応$k$を返したら$Y_{ijk}=1$。
　潜在変数の数を$R$とする。事前の所属確率を$p_r$とする($\sum_r p_r =1$)。クラス$r$において$Y_{ijk}$に1が立つ条件付き確率を$\pi_{jrk}$とする($\sum_{k=1}^{K_j} \pi_{jrk} = 1$)。

　クラスの下でのアウトカムの条件つき独立性を仮定すれば、クラス$r$の対象者$i$が特定のアウトカムのベクトルを生む確率は
　$f(Y_i; \pi_r) = \prod_{j=1}^J \prod_{k=1}^{Kj} (\pi_{jrk})^{Y_{ijk}}$
　全クラスを通した確率密度関数は
　$P(Y_i|\pi, p) = \sum_{r=1}^R p_r \prod_{j=1}^J \prod_{k=1}^{Kj} (\pi_{jrk})^{Y_{ijk}}$
　パラメータ$p_r$, $\pi_{jrk}$が推定できたとして、対象者$i$がクラス$r$に所属する事後確率は
　$\displaystyle \hat{P}(r_i|Y_i) = \frac{\hat{p}_r f(Y_i; \hat{\pi}_r)}{\sum_{q=1}^R \hat{p}_q f(Y_i; \hat{\pi}_q)}$

　さて、poLCAはパラメータを最尤法で推定する。対数尤度関数は
　$\log L = \sum_{i=1}^N \log \sum_{r=1}^R p_r \prod_{j=1}^J \prod_{k=1}^{Kj} (\pi_{jrk})^{Y_{ijk}}$
これをEMアルゴリズムで最大化するわけです。
　標準誤差はどうやって推定するかというと... [略。観測情報行列を使う由]
　適合度指標としてAIC, BIC, ピアソンのカイ二乗、尤度比カイ二乗を提供する。

　次にlatent class regressionモデルについて。ここでいうところのlatent class regressionとはlatent class models with covariates、すなわち潜在クラスの所属確率が共変量で予測されるのを同時推定するモデル。いっぽうfpcパッケージのregmix関数やflexmixパッケージでいうところのlatent class regressionとは、回帰モデルの推定の一部として従属変数が潜在クラスに分割されるというモデルである。意味がちがうので注意。
　モデルは... パラメータ推定は... 標準誤差は...[面倒になってきたのでパス]

　後半は使用方法と例題。モデルからの予測されたセル頻度表を出力する機能、いったん推定した後でクラスの順序を並び替える機能、シミュレーションしたデータを作る機能がある由。局所解の問題はやっぱり深刻らしくて、何度も推定しなさい、とのことである。
　
　私自身はMplusの忠実なしもべなので、潜在クラスモデルを組もうかなと思った次の瞬間には無意識のうちに library(MplusAutomation) と叩いてしまうのだが、Rのなかで完結できたら助かるなと思うことも、まあ、なくはない。このパッケージだとMplusほど細かい制約は掛けられないけど、いつか使ってみてもいいかもな。

論文：データ解析(2015-) - 読了：Linzer & Lewis (2011) RのpoLCAパッケージ

« 2017年8月 | メイン | 2017年10月 »

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