490 likes | 708 Views
かなた望遠鏡→ 東広島天文台 外観↓. ベイズモデル ~最近の仕事から成分分離やトモグラフィーなど~. 植村誠 (広島大学 ) @京大宇物談話会 110112. 今日 の話の流れ. ベイズ統計とは ベイズモデルを使った最近の仕事例 矮新星のアウトバースト頻度を考慮した軌道周期分布の推定 ブレーザーの可視偏光データを利用した長期・短期変動成分の分離 矮新星の早期スーパーハンプを利用した降着円盤の高さマップの再構築 まとめ. ベイズ統計とは. ベイズの定理 (Bayes’ theorem). データ A , モデルパラメータ B に対して、条件付確率
E N D
かなた望遠鏡→ 東広島天文台 外観↓
ベイズモデル~最近の仕事から成分分離やトモグラフィーなど~ベイズモデル~最近の仕事から成分分離やトモグラフィーなど~ 植村誠 (広島大学) @京大宇物談話会 110112
今日の話の流れ • ベイズ統計とは • ベイズモデルを使った最近の仕事例 • 矮新星のアウトバースト頻度を考慮した軌道周期分布の推定 • ブレーザーの可視偏光データを利用した長期・短期変動成分の分離 • 矮新星の早期スーパーハンプを利用した降着円盤の高さマップの再構築 • まとめ
ベイズの定理 (Bayes’ theorem) • データA,モデルパラメータBに対して、条件付確率 • P(B):パラメータBの事前確率 (prior probability distribution) • P(A|B) = L(B|A) :パラメータBの場合にデータAが生じる確率(=尤度; likelihood) • P(B|A) : データAが観測された場合のパラメータBの事後確率 (posterior probability distribution) • (Aが観測結果(定数)なのでP(A)は定数) • 1763年、T. Bayesが発表。 • データから知りたいパラメータの「確率密度分布」が得られる。 • でも「事前分布」ってなんじゃ?
頻度主義とベイズ主義 統計学者や科学哲学者ではない「1ユーザー」としては、「道具」としての有用性が 高くて理論がまともなら、主義は何でも良いでしょう。 でもベイズの定理の分母の積分項が解析的に解けない場合はどうする?
ベイズモデルとマルコフ連鎖モンテカルロ法(MCMC)ベイズモデルとマルコフ連鎖モンテカルロ法(MCMC) • モンテカルロで事後分布の推定 • 計算機の進歩と共にここ数十年で急速な進歩 • MCMCの1つ、Metropolis法のアルゴリズム • パラメータの組を x={x0,x1,….,xn} , 推定したい確率分布を P(x)として、 • xiから xi+1の候補を以下のように生成 : N(μ,σ)は平均μ,分散σ^2の正規分布 • P(xi+1)>P(xi)なら候補を採択 • P(xi+1)<P(xi) なら、確率 P(xi+1)/P(xi) で採択 • それ以外の場合は xi+1 = xi 尤度 事後確率密度 MCMC。 分布を推定。多峰型の分布にもメリット。 最尤法で使われる アルゴリズム。 間違ったピークを捕まえることも パラメータ パラメータ
実例:滑らかな曲線モデル(岩波書店 「階層ベイズモデルとその周辺」 石黒真木夫、他著 より) • 課題:時系列データ a={ai} (i=1-10) をうまいこと通るような滑らかな曲線 b={bj} (j=1-50) を推定する。 • モデル • (尤度関数) (事前分布:2階階差数列が正規分布に従う)
宇宙物理学におけるベイズモデルの応用 (ADSでアブストラクトに”Bayes”の文字列が含まれる査読論文) この10年間で急速に普及
ベイズモデルの利点・弱点 • 利点 • 既知情報を事前分布としてモデルに組み込める • MCMCとの相性の良さ • 逆問題も解ける • パラメータスペース総当たりのモンテカルロよりは「次元の呪い」に強い • 高次元の画像再構成など。 • MEMなど従来の手法よりも自由度が高い • 宇宙物理に向いている面も • 何度も繰り返し「実験」ができない、光子が少ない、観測機会が少ない、情報が縮退してる、などなど。 • 弱点(注意点) • 適切な尤度関数(モデル)の設定が必要なことは通常の統計理論と同じ • 事前分布の設定の妥当性・正当性 • 事前分布にハイパーパラメータが含まれる場合の扱い • 階層ベイズ or 経験ベイズ • MCMCの収束判定、サンプリングにもパラメータが必要 • 最近、ベイズ的な仕事をいくつかやってきたので、今日はそれらを紹介します。
ベイズモデルを使った最近の仕事例 その1矮新星のアウトバースト頻度を考慮した軌道周期分布の推定(Uemura, et al., 2010, PASJ, 62, 613)
激変星とは • 低質量連星系の1つの末路 伴星(主系列星) 降着円盤 白色矮星
軌道周期分布の問題 • 激変星の進化=連星系からの角運動量の抜き取り • Magnetic braking • Gravitational wave • 最小周期 • 伴星のコアが縮退(brown dwarf)しはじめると、星の「質量ー半径」の関係が変化 • 通常の伴星:短周期の系へ進化 • 縮退した伴星:長周期の系へ進化 • 観測される最小周期 ~80分 • 連星の初期状態を仮定して population synthesis • 大半の激変星は最小周期付近にいるはず(period spike)。 • 最小周期(period minimum)は60-70分 • 観測される分布を説明できない
WZ Sge型矮新星 • 矮新星 • 降着円盤の熱的不安定性で増光 • 増光間隔は数日~数年 • WZ Sge型矮新星 • 10年以上の増光間隔 • 軌道周期が最小周期付近 • 最近のAll sky monitorなどで発見頻度が増加 • 軌道周期分布が変化してきた →矮新星のアウトバースト頻度を考慮することが重要 2003年までの軌道周期分布と、2008年までの軌道周期分布。最小周期付近で特に増加している。
観測された軌道周期分布から本来の分布をベイズ的に推定観測された軌道周期分布から本来の分布をベイズ的に推定 • 増光が観測される天体の軌道周期の確率密度関数 軌道周期と増光頻度の関係 • 観測されるアウトバースト =ある確率分布からの標本、と考える 本来の軌道周期分布を表わす 確率密度関数のモデル 爆発頻度に依存した発見確率 絶対等級に依存した発見確率
サンプル • 2003年から2007年の5年間で見つかった矮新星アウトバースト • ASASで42天体 (単一の観測システム) • VSNETで146天体 (様々な観測システムが混合)
パラメータの事後分布 • 最小周期は観測値よりも短くなる
矮新星の軌道周期分布のベイズ的推定まとめ • 最近の観測から、WZ Sge型矮新星が従来よりも高い頻度で発見。SDSSによる暗い(伴星からの質量輸送率が低い=最小周期付近の)激変星の検出も合わせて、軌道周期分布が変わってきた。 • 矮新星のアウトバースト頻度を簡単な関数でモデル化して、ベイズ的に本来の軌道周期分布を推定。 • Period minimum問題、Period spike問題は WZ Sge型星(=増光頻度の低い矮新星)を考慮すれば解決できるかもしれない。
ベイズモデルを使った最近の仕事例 その2ブレーザーの可視偏光データを利用した長期・短期変動成分の分離(Uemura, et al., 2010, PASJ, 62, 69)
ブレーザーとは blazar • ジェットの研究に貴重な天体 Optical-NIR hope High polarization fact A probe of the magnetic field structure But, no universal law has been established
ブレーザーの偏光観測 Positive correlation between the flux and PD. (Smith, et al. 1986) • 偏光の変動に規則性はあるのか? Apparently random motion in the QU plane (Moore, et al. 1982)
偏光の2成分モデル • 長期変動と短期変動に分離したい Stokes’ QU parameters of BL Lac for ~20 yr (Hagen-Thorn, et al. 2002) Stokes parameters for linear polarization (0,0)
偏光の成分分離のためのベイズモデル U • 短期変動の偏光フラックスが総フラックスと相関するように長期変動を推定する Prior distribution of the long-term trend (smoother curve is preferred). Posterior distribution of the long-term trend Q Likelihood function, maximized when the total flux and polarized flux completely correlate. flux t
人工データを使ったテスト • 期待通りの結果 pol. flux pol. flux pol. flux pol. flux total flux total flux total flux total flux corrected pol. flux long-term trend (0,0) (0,0) U U U U Q Q Q Q
人工データを使ったテスト • 期待通りの結果 pol. flux pol. flux pol. flux pol. flux total flux total flux total flux total flux corrected pol. flux long-term trend (0,0) (0,0) U U U U Q Q Q Q
人工データを使ったテスト • 期待通りの結果 pol. flux pol. flux pol. flux pol. flux total flux total flux total flux total flux corrected pol. flux long-term trend long-term trend (0,0) (0,0) U U U U Q Q Q Q
人工データを使ったテスト • 期待通りの結果 pol. flux pol. flux pol. flux pol. flux total flux total flux total flux total flux corrected pol. flux long-term trend (0,0) (0,0) U U U U 期待通りの結果 Q Q Q Q
実際の観測データでは • S2 0109+224 • OJ 287 pol. flux pol. flux total flux total flux total flux total flux pol. flux pol. flux (0,0) U U U U (0,0) Q Q Q Q
実際の観測データでは • S2 0109+224 • OJ 287 pol. flux pol. flux total flux total flux total flux total flux pol. flux pol. flux 長期成分の偏光方位角が振動? 長期的なジェット軸の運動か? (0,0) U U U U (0,0) Q Q Q Q
実際の観測データでは • S2 0109+224 • OJ 287 pol. flux pol. flux total flux total flux total flux total flux pol. flux pol. flux 長期成分の偏光方位角が振動? 長期的なジェット軸の運動か? (0,0) U U U U (0,0) Q Q Q Q
実際の観測データでは • S2 0109+224 • OJ 287 pol. flux pol. flux total flux total flux total flux total flux pol. flux pol. flux 長期成分の消長が見えている? 長期成分の偏光方位角が振動? 長期的なジェット軸の運動か? (0,0) U U U U (0,0) • どちらも興味深い結果に。 • 1シーズンだけでなく、複数シーズンの結果を見たい。 • 長期成分の挙動と電波の挙動とに相関があると面白いかも。 Q Q Q Q
ベイズモデルを使ったブレーザー偏光の成分分離まとめベイズモデルを使ったブレーザー偏光の成分分離まとめ • ブレーザーの偏光に法則性が見られないのは長期変動成分が短期フレアに重なっているからかもしれない。 • その仮説に基づいてベイズ的に偏光成分を分離するモデルを開発。 • いくつかの天体では2成分モデルで上手く説明できることが可能で、「長期偏光成分の挙動」という新しいアプローチを得た。
ベイズモデルを使った最近の仕事例 その3矮新星の早期スーパーハンプを利用した降着円盤の高さマップの再構築(Uemura, et al., まだ解析中)
早期スーパーハンプとは • 円盤への潮汐効果 →ただし詳細不明 上:早期スーパーハンプと通常のスーパーハンプの光度曲線(と色変化)。データはV455 Andのもの。 下:WZ Sgeのアウトバースト光度曲線
円盤の高さを見ている? Early superhumps in V455 And (Matsui, et al. 2009) ハンプ成分が赤い →円盤外縁が膨張? Zoo of early superhumps (Kato , et al. 2002) Edge-onの天体ほど振幅が大きい? • 回転効果 →観測から円盤を再構成できるかも Two-armed spirals on the disk and early superhumps (Maehara, et al. 2007)
早期スーパーハンプのモデル • Tidally distorted disk • Kato (2002) • Even below the 2:1 resonance radius • The 2:1 resonance • Osaki & Meyer (2002) • 観測から円盤を再構成して理論モデルと比較 Tidal truncation and resonance radii Height distortion in the disk (Ogilvie, 2002) 伴星
ベイズを使った円盤高さの再構成 • →モデルと観測の光度曲線で尤度 • →隣り合うグリッドの高さはなるべく滑らかに(事前分布) • →標準円盤(h=0.1r)を仮定(事前分布) • 使う観測結果は、g,V,Rc,Ic,Jバンドのphase-averaged light curve (20 points/phase) • 座標(r,θ)での高さhを推定する。 • 新しい h 候補を発生 • モデル光度曲線 →尤度計算 • Metropolis または Multiple-Try Metropolis 法で採択・棄却を判定 • 前2例よりも複雑なモデル、多い次元
人工データを使ったデモ • 外側内側が区別できている 観測 再構成
矮新星 V455 And のデータを使った結果(preliminary) 伴星 • アウトバースト極大から5日目のデータ • 結果が何を語っているかはこれから検討
ベイズモデルによる降着円盤の再構成まとめ • 早期スーパーハンプは回転してる円盤の高さ方向の非軸対称性を見ている。 • 多バンド測光観測から円盤の構造をベイズ的に再構成するモデルを開発(中) • 使い物になりそう。今年中に論文化予定。 • 早期スーパーハンプの多バンド光度曲線から、円盤構造の時間変化が追える
今日の話のまとめ • ベイズ統計+MCMCを使った手法は宇宙物理の分野でもここ10年で急速に広がった。 • 成分分離、逆問題、少ない観測点、などの状況で便利 • 3例を紹介 • 矮新星の本来の軌道周期分布の推定 • ブレーザーの偏光成分の分離 • 早期スーパーハンプから円盤構造の再構成 • 1000次元以下なら収束する実績あり。アイデアのある方、いっしょに勉強しませんか?
伴星 推定された円盤。極座標でグリッドを設定して、かつ動径方向は対数スケールで等間隔にしてる。普通のソフトだと表示しにくいので線形スケールで等間隔に補間した。 WZ Sgeアウトバースト初期のHeIIのドップラートモグラフィー(Baba, et al. 2002)
外縁温度を1割高く 今回のデフォルト 温度勾配を標準円盤よりも緩く (p=-0.75→-0.65) Inclination = 80 deg.
ここのギャップは不明。赤側では目立たないので、定性的には円盤の比較的外側がかなり膨らんでいて、内側の高温領域を隠せば良い。ただしあまり膨らみ過ぎるとそもそも暗くならない。ここのギャップは不明。赤側では目立たないので、定性的には円盤の比較的外側がかなり膨らんでいて、内側の高温領域を隠せば良い。ただしあまり膨らみ過ぎるとそもそも暗くならない。 等級 光度曲線の比較。赤が観測で、青が再構成されたもの。 図中のphase=0.2くらいの極小が、orbital phase = 0.0(伴星が手前)。 観測値の誤差は一律0.01を仮定した。 ここのギャップは伴星による円盤外縁部の食で説明できそう。Inclination~75degだと、内側までは隠せない。
横軸:MCMCのステップ、縦軸:事後確率密度。横軸:MCMCのステップ、縦軸:事後確率密度。 混合がよくない。
MCMCサンプリングの違い。 左:3万点から10万点まで100点ずつ。右:3万点から100万点まで100点ずつ。 スケールもあるが、特に内側の詳細構造が見た目に違う? 点数が多い右の方がよりスムーズ。ただし、正しい結果を出すためにはまず混合を 改善しないといけない。
設定したパラメータの詳細 • V455 And のパラメータ (Araujo-Betancor, et al. 2005)より • Incliantion = 75deg • 浅い食があるため。精度は良くないはず。 • Mass ratio~9 (M1~0.6Mo, M2~0.07Mo) • Porb= 0.056309 d • モデルに使ったパラメータ • Inclination = 75 deg • Rin=1(WD半径)として、Rout=30(2:1 resonance に対応) • 外縁温度 6800K • 観測された日の平均の色から推定 • 温度は T∝r^(-0.75)(標準円盤)、ただし、r = sqrt(x^2+y^2)、つまり高さの寄与は無視。 • 計算速度を上げるため。 • 動径方向16分割、方位方向20分割 • Irradiationや伴星の食の効果は入っていない