Physics of Eclipsing Binaries. V. General Framework for Solving the Inverse Problem          ApJ Spl, 250:34 (17pp), 2020 October Kyle E. Conroy, Angela Kochoska, Daniel Hey, Herbert Pablo, Kelly M. Hambleton1, David Jones, Joseph Giammarco, Michael Abdul-Masih, and Andrej Pr?a  日本語訳について:DeepL翻訳したものを大島が修正した。誤りに気付いたらぜひご一報をお願いします。式番号のみで式本体を省略した箇所については、元の論文を参照されたい。(2022年3月4日) 概要 PHOEBE 2は、食連星系の観測データをモデル化するためのPythonパッケージですが、これまではフォワードモデル、つまり、星系と観測データを記述する多数のパラメータに固定値を与えて合成モデルを生成することに焦点を当てていました。観測データから軌道や恒星のパラメータを得るという逆問題は、より複雑で計算コストがかかります。最適な解を求め、さらにそれらのパラメーターの信頼できるロバストな不確かさを得るためには、オプティマイザーやサンプラーなど、複数のアルゴリズムを使用する必要があります。さらに、PHOEBEのフォワードモデルは、物理的に可能な限りロバストになるように設計されていますが、他のコードに比べて計算量が多くなります。そのため、特定のシステムに対する合理的な仮定の下で最も効率的なコードを使用することは有用ですが、複数のコードの複雑な仕組みを学ぶことは、実際にこれを行う上での障壁となります。ここでは、PHOEBEのリリース2.3(http://phoebe-project.org)を紹介します。PHOEBEは、パラメータの分布を定義して扱い、複数の異なる推定、最適化、およびサンプリングアルゴリズムを利用するための一般的なフレームワークを導入しています。このフレームワークは、PHOEBE自身に組み込まれたロバストモデルを含む、複数のフォワードモデルをサポートしています。 1. はじめに 食連星系は、恒星構成要素の絶対的な基本パラメータの決定、恒星進化モデルの較正、軌道の向きや力学的変化の測定、較正を必要としない距離の導出などが可能である。ベンチマーク級の食連星では、恒星の基本的なパラメータ(質量や半径)の精度は数パーセント程度です(Torres et al.2010)。これらの科学的目標は、逆問題を解決することによってのみ達成されます。すなわち、一連の入力パラメータを用いて合成観測値をシミュレートするモデルを作成し、入力パラメータを変化させて、どの値が観測結果を最もよく表す合成観測値になるかを判断するのです。しかし、パラメータ空間は複雑で、滑らかではないことが多く、この作業にはコストと時間がかかります。 PHOEBEは、重力で歪んだ恒星の表面を離散化し、その離散化されたグリッドに大気テーブルからローカルな量を入力し、軌道上の任意の時間に可視要素を積分することでモデル観測値を合成する数多くのモデリングコードの1つです。これまで、2.0リリース以降のすべてのリリース(Pr?a et al.2016; Horvat et al.2018; Jones et al.2020)は、フォワードモデルに高度な物理学と精度を追加することに重点を置き、逆問題をユーザーに委ねてきました。  パラメータ空間の複雑さ、利用可能な観測データの組み合わせの違い、システムの形態や関連する物理効果の幅広さなどから、可能な限り正確で精密なモデルと不確実性を得るためには、個々のシステムを手作業で処理する必要があります。さらに、最適解を導き出す確率を高め、モデル化されたパラメータに信頼性の高い現実的な不確実性を持たせるためには、複数のアルゴリズムを利用する必要がある場合もあります。それぞれが独自のインターフェースを持つ複数のアルゴリズムをフォワードモデルに「巻き付ける」必要があるため、これらのシステムを解くための学習曲線と労力が急に必要になることがあります。そこでPHOEBEでは、逆問題を解くための様々なアルゴリズム(市販のものやカスタムメイドのもの)と対話するための共通のインターフェースを用意しました。 さらに複雑なことに、利用可能なフォワードモデルコードはそれぞれ独自のインターフェースを持ち、パラメータ化も異なることが多いため、フォワードモデルを比較したり、特定のシナリオに最適なコードを選択したりすることがさらに困難になっています。PHOEBE 2は、精度、ロバスト性、柔軟性を重視して設計・開発されていますが、その分、計算効率が悪くなっています。多くの場合、他のコードでも同じようにシステムをモデリングすることができますが、時間的なコストはほんのわずかです。また、それぞれのコードの複雑な仕組みを学ぶことは、新たな学習曲線とオーバーヘッドをもたらします。PHOEBEでは、一般に公開されている複数のコードをフォワードモデルとして採用し、必要に応じて入力と出力のパラメータ設定や単位を変換しています。これにより、異なるコードで作成した合成モデルを簡単に比較したり、特定のシステムに最も適した効率的なコードを選択したり、逆問題を解決する際に、より高価だが精密なモデルに移行することができるなど、いくつかの大きなメリットが得られました。 ここでは、PHOEBEで逆問題を解くための一般的なフレームワークを紹介します。今回のPHOEBEに搭載されているソルバーアルゴリズムを3つのカテゴリに分けています。セクション2では、PHOEBEで利用可能な「推定器」について詳しく説明します。推定器は、フォワードモデルを必要とせず、観測データのみからシステムパラメータを推定します。これらの推定器は、利用可能なデータがあれば、特定のシステムのパラメータ空間を安価に決定できるように設計されています。第3章では、これらのフォワードモデルから得られるメリット関数について説明します。このメリット関数は、「オプティマイザー」(セクション4)と「サンプラー」(セクション5)の両方で使用され、最適解を見つけ、ロバストな事後確率分布を決定する(セクション6)。セクション7では,今回のリリースで導入された新しいユーザー・インターフェースを紹介し,Webブラウザまたは専用のデスクトップ・クライアントのいずれかから逆問題の機能のほとんどに簡単にアクセスできるようになっています。本論文でシステム・パラメータに使用されている記号の説明については,付録を参照してください。 本論文にはPHOEBEの2.3リリースが添付されており、http://phoebe-project.org および https://github.com/phoebe-project/phoebe2 から入手可能です。本論文に掲載されているすべての図を複製するためのスクリプトは, http://phoebe-project.org/publications/2020Conroy+ で入手できます. 2. エスティメーター推定器  複数の物理モデルを計算することなく、可能な限り多くのパラメータの値を推定することは、有用であり、また計算上も賢明であることが多いです。PHOEBE 2.3では、観測データに直接作用するいくつかのアルゴリズムを導入しています。これは、ユーザーの監視なしに、関連するシステム・パラメータのサブセットに対する初期解を高速に提供することを目的としています。 1. 1.RVピリオドグラム、LCピリオドグラム:光度曲線、速度曲線のそれぞれのピリオドグラムを計算し、軌道周期Porbを提案します。 2. 2. RVジオメトリ:単純なケプラー軌道にフィットさせ、t0,supconj, e, ω0, vγ, qと軌道半長軸の投影値(a_orb sin i)、または成分ごとの半長軸の投影値(a_comp sin i)の値を提案する。 3. LCジオメトリ:食の極小、減光部、増光部の位相を推定し、その情報を用いてt0,supconj,e,ω0の値を提案する;および 4. EBAI:入力された光度曲線に対して学習済みの人工ニューラルネットワークを使用して,t0,supconj,Teff,2/Teff,1,(R_equiv,1+ R_equiv,2)/a_orb,e sin ω0,e cos ω0,i の値を提案する. 時間空間ではなく位相空間で動作するすべての推定器では,大規模なデータセットに対する推定器の効率を維持するために,入力データを任意にビン化することができます。 2.1. ピリオドグラム 食連星系をモデル化する際に使用される最も一般的な観測データセットは、光度曲線と視線速度曲線です。大半の食連星系(時間に依存した大きな変化を示さない星系)では、これらのデータは星系の軌道周期に基づいて折り畳むことができます。周期がよくわかっていて、軌道パラメータが連続した軌道間で変化しない場合は、代わりに位相空間でシステムをサンプリングすることができます。 PHOEBE 2.3には、 astropy.timeseries (Astropy Collaboration et al. 2013)で実装されている2つの一般的なピリオドグラム、box least squares (BLS, 光度曲線用)とLomb Scargle (光度曲線または視線速度用)のラッパーが含まれています。このラッパーは、入力として任意の数の光度曲線または視線速度を受け取り、さらに、周期/周波数での自動または手動のサンプリング、BLSの目的関数を変更するオプション、BLSの提案された食期間を設定する機能など、astropyのいくつかの高度なオプションを使用します。複数の光度曲線で実行する場合、各光度曲線は、ピリオドグラムアルゴリズムに送信する前に、中央値または最大フラックス値のいずれかで除算して正規化されます。リクエストされたすべてのデータセットの視線速度は結合され、主星と伴星それぞれの絶対的な最大値で正規化され、伴星はミラーリングされています。ラッパーを実行すると、ピリオドグラム自体が返され、プロットすることができ、ピーク周期が採用されます。特に似たような食を持つ円形に近い星系では、これらのアルゴリズムは真の軌道周期のエイリアスで強い信号を見つけることが多いので、提案されたピーク周期の任意の係数を採用するための簡単なインターフェースが提供されています。 2.2. 視線速度の幾何学 離心率(e)、近日点引数(ω)、システム速度(vγ)、近日点通過時刻(t0,supconj)の推定には、解析的な視線速度(RV)を用いる。また、質量比(q)と、両成分星のRVが提供されている場合は投影軌道半長軸(a_orb sin i)、どちらか一方の恒星成分のRVしか提供されていない場合は単一成分の投影半長軸(achomp, i sin i)も推定されています。すべての推定値は、真近点離角(?; Pr?a 2018)の関数としての成分ごとの視線速度の方程式に基づいています。 (1) (2) 周期(P_orb)は既知であると仮定し、入力として位相折り返し視線速度曲線を使用します。データは、高周波ノイズを除去するために、ローパスのSavitzky-Golayフィルターで平滑化されています。 並行して推定される最初の2つのパラメータは、質量比とシステム速度です。式(1)から、両パラメータについて軌道上の任意の点で成立する解析式を導き出すことができる。 RV1(θ)+qRV2(θ) Vγ = ------------------ (3) 1 + q RV1(θ) - Vγ q = ------------------- (4) -RV2(θ) + Vγ 見ての通り、vγとqは両方の式に登場するため、それぞれを独立して決めることはできません。そのため、主星RVと伴星RVの比としてqを粗く見積もることから始めて、反復的に計算します。効率を最大化するために,最大の主星RVに対応するRVポイントでのみqとvγを推定します。利用可能なRVが1つしかない場合(シングルライン連星系),qを確実に推定することはできず,この場合,vγのみを利用可能なRVの中間点として推定する。 qとvγが推定できたら、投影された半長軸a sin iを計算することができます。まず、各RVで利用可能な観測データから半長軸を推定します。 K_i = 0.5[max(RV_i - v_γ) - min(RV_i - v_γ)] (5) ここで、i=1,2はそれぞれ主成分と副成分を表す。また、対応する投影半長軸は K_i P_orb a_i sin i = ----------- √(1-e^2) (6) 2π 最初はe=0と仮定し、両方のRVが利用可能な場合はa sin i = a1 sin i + a2 sin iとする。 q, vγ, a sin iが推定されたら、それらを解析的な視線速度に固定し、離心率と近星点引数をフィッティングします。実現可能なパラメータ空間全体をカバーするように、e0=[0.4]とω0=[0, π/2, π]の組み合わせで小さなグリッドを開始点として反復し、最小二乗解を求めます。各反復の後,a sin iの値は,更新された離心率eで再計算され,その解は次の反復の初期点として使用される. 最後に、近星点通過時刻の位相をv=v_γとなる位相として推定し、近星点通過時刻t0, supconjに変換します。 この推定器の性能の一例を図1に示します。 PHOEBE 2.3では、2ガウスモデルを位相光度曲線にフィッティングすることで、光度曲線の形状から最も制約を受ける2つのパラメータである離心率と近日点の引数を評価する推定器を導入しました。モデルの詳細については、Mowlavi et al.(2017)を参照してください。ここでは、PHOEBEでの実装に焦点を当てます。 2ガウスモデルは、一定の余弦項、2つのガウシアン、またはこれらの任意の組み合わせの複合関数を観測データにフィットさせる高速解析モデルで、結果として7種類のモデルが得られます。最も単純なモデルは定数項(C)のみで構成され、最も複雑なモデルは定数、2つのガウス関数、余弦項で構成されています。それぞれのガウス関数は,中心位置(μ),振幅(A),rms幅(σ)で特徴づけられます。コサイン項は、その周期が光度曲線の周期の半分(位相で折り畳まれた光度曲線の場合は0.5)に等しく、その谷が主食または副食の位置と一致することを想定しているため、振幅(A_ell)のみが特徴となっています。 この2ガウスモデルは、モデルパラメータの初期値に敏感です。そのため、PHOEBEは収束性を確保するために、まず、光度曲線の最小値を検索し、その周辺で位相(j)空間の中央部のフラックスを横切るデータを分離するというシンプルなアルゴリズムを用いて、食の位置、幅、深さを推定します。その後、7つのモデルすべてを最小二乗法によるオプティマイザーでフィットさせ、ベイズ情報量基準(BIC; Schwarz 1978)が最も高いものをベストフィットとした。食の位置はガウシアンの中心位置と一致し(φ_i=μ_i)、食の幅はガウシアンのrms幅の1倍(w_i=5.6σ_i)、深さは定数項から食位置でのフラックスを差し引いた値(d_i = - C flux(φ_i))として計算するというMowlaviら(2017)の処方箋に従って食パラメータを決定する。図2は、ベストフィットモデルを用いたケースの例と、結果として得られた食の位相を示しています。 この推定器は、近星点通過時刻、離心率、近点離角を以下のように計算します。 1. 近星点通過時刻 t_0,supconj = t_0,supconj + φ_1 P_orb (7) 2. 離心率 e =[ ]^1/2 (8) 3. 近点離角。 1 w2 - w1 ω1=arcsin(-----------) (9) e w2 + w1 √(1-e^2) ω2= arccos(--------------) (10) 2e tan(ψ-π) ω2 if ω1 >= 0 ω= { (11) 2π-ω2 if ω1 < 0 ここで、ψ = π+ 2 arctan ecos ω/√(1- e^2)であり、2πΔΦ = ψ- sin ψを解くことで繰り返し計算される。ΔΦは位相空間における2つの食の間隔である。 さらに、減光増光期の位相、食の幅と深さが計算され、ユーザーに公開されます。そして、これらを用いて、食のエッジと食幅の30%のパディングに基づいたオプションの位相マスクも提案されます(3.4.1項および図5参照)。 2.4. 人工ニューラルネットワーク(EBAI) EBAIは、食連星系の位相光度曲線から基本パラメータを復元するために特別に実装されたバックプロパゲーション型の人工ニューラルネットワークです(Pr?a et al.2008, 2019)。PHOEBE 2.3内での実装では、2Gaussianモデル(セクション2.3参照)を用いて、201個の入力ユニット、40個の隠れユニット、5個の出力ユニット、学習率パラメータη=0.2、33,235個の模範例、4×10^6回の反復を用いて、離散系のネットワークを再学習します。 推定器としてEBAIを使用するには、ユーザーは任意の数の光度曲線データセットを渡すことができます。これらのデータセットは正規化され、2ガウスモデルでフィッティングされます。この解析的表現は、1に正規化され、トレーニングセットと同じ201の位相ポイントでサンプリングされ、EBAIのトレーニングで決定された2つの行列変換に渡されます(図3)。EBAIから直接受け取った結果は、必要に応じて、適用された位相シフトから計算されたt0,supconjに加えて、T T eff,2 eff,1, (R R a equiv,1 equiv,2 or + ) b, e sin ω0, e cos ω0, iの提案された値に変換されます。 3. メリット関数 フォワードモデルの計算を必要とする最適化(セクション4参照)やサンプリング(セクション5)では、合成モデルと観測値を比較するためのメリット関数の定義が必要です。メリット関数(lnprobability)の構成要素とPHOEBEで採用されている用語を以下に説明します。 1. lnpriors: log priors項は、現在の顔の値pを事前の確率分布関数πから導く対数確率として定義されます。 priorsが提供されていない場合や適用可能な場合、この項は一定となるため、最適化やサンプリングには寄与しません。事前確率分布については、セクション3.2で詳しく説明しています。 lnpriors = ln(P(p|Π)) (12) 2. 残差:観測値(yo)-合成モデル(ym)、位相マスキング(セクション3.4.1参照)後のフラックスまたは視線速度である。デフォルトでは,合成モデルは,(マスクされた)観測データの正確な時刻に計算されます。しかし,ユーザはカスタム時刻のリストでモデルを計算することもできます。その場合,合成モデルは観測の時間(または位相)に線形補間され,その後,残差が決定される.線形補間が問題を起こさないように、モデルが十分にサンプリングされていることを確認するのは、ユーザーの責任です。フォワードモデル自体、合成モデルの時間、位相マスキング、補間についての詳細はセクション3.4を参照してください。 residuals = y=o(t_masked) - y_m (t_masked) (13) 3.カイ2乗(χ2):ここでのχ2は、有効な「データセット」全体の全データポイントの残差の二乗と、ポイントごとの不確かさ(提供されている場合)の二乗の和に加え、提供された観測の不確かさからの不確かさの過小評価を処理するための追加項σoを加えたものと定義されています。この不確かさの過小評価の項(σlnf)は、提供された不確かさが過小評価されていると思われる場合に、最適化またはサンプリングすることができます。デフォルトでは、σlnf = -infで、分母の追加項がキャンセルされます。この "ノイズ厄介者"パラメータについては,3.5節を参照してください。 (y_o - y_m)^2 Χ^2 = -------------- + ln(σ^2) (14) データセットのσ^2 ここで σ^2 = σ_o^2 + y_m^2 e^(2σlnf) (15) 4.lnlikelihood:対数尤度は、定義によれば lnlikelihood = -0.5χ^2 (16) lnprobability:対数確率は、PHOEBEで使用されるデフォルトのメリット関数です(技術的には、負の対数確率は、メリット関数を最小化するオプティマイザに使用され、最適なモデルは対数確率を最大化することによって見つけられることに注意してください)。事前分布が提供されていない、または適用できない場合、これは本質的にχ2と類似したものになります。対数確率は、対数事前分布と対数尤度の和として定義されます。 lnprobability = lnpriors + lnlikelihood. (17) 6. lnposteriors: log posteriorsはlog priorsと同様に、結果として得られる事後確率分布関数が与えられた場合に、現在の顔の値pを描画する確率を記述する。これらの後値(セクション6参照)は、マルコフ連鎖モンテカルロ(MCMC;セクション5.1参照)などのアルゴリズムを用いて対数確率をサンプリングして決定されることが多い。 lnposteriors = ln(P(p|posterior pdf)) (18) カスタムメリット関数でオーバーライドされない限り、PHOEBEのすべてのオプティマイザとサンプラは、上記の定義のとおり対数確率を使用します。さらに、ユーザーがPHOEBEの外で独自のオプティマイザーやサンプラーを実装したい場合は、任意のフォワードモデルに対して上記のすべての量がユーザーに公開されます。 3.1. パラメータ化 利用可能な観測データ(単一の光度曲線、複数の光度曲線、一重の視線速度、二重の視線速度など)、既知の情報、パラメータ空間の直交性(またはその欠如)に応じて、特定のパラメータ化を選択できることが有利です。例えば、Rocheモデルのパラメータ化により、PHOEBE1(Wilson-Devinneyベース)では、軌道周期(Porb)、質量比(q)、半長軸(a_orb;離心率e、傾斜角や周回軌道の引数などの方位パラメータも含む)で軌道をパラメータ化しています。個々の恒星の質量(M1とM2)は、これらのパラメータとケプラーの第三法則から計算できますが、直接固定したり調整したりすることはできません。 PHOEBEには、このようなシナリオに対応するための柔軟なフレームワークが用意されています。例えば、デフォルトではPHOEBEは従来のバージョンと同様のパラメータ設定を採用していますが、ケプラーの第三法則によって導かれる質量については読み取り専用のパラメータを公開しています。 4πa^3 1 4πa^3 M_1= ------- ----- , M_2 = -------(1+q) (19) GP^2 1+q GP^2 この5つのパラメータ(M1, M2, P, a, q M M ? 2 1)のシステムは、Pythonインターフェースから操作することができ、5つのパラメータのうち任意の2つは、ケプラーの第3法則に基づいて、読み取り専用として「制約」されます。他の2つのパラメータを読み取り専用にすることで、個々の質量の設定やフィッティングが可能になります。 連星モデルでは、他にも同様のパラメータ化の選択肢が数多くあります。ローカルなパラメータ空間内での直交性や、事前に知っている情報の利用可能性などから、あるパラメータ化が他のパラメータ化よりも優先される場合があります。表1は、PHOEBEの2.3リリースに含まれる適用可能な制約のリストです。 表1の注意: 表の上半分のパラメータはデフォルトでシステムに含まれていますが、下半分のパラメータはオプションで追加できます。右列のパラメータは、デフォルトではリードオンリーで、左2列のパラメータから派生したものです。右列のパラメータのリードオンリーステータスは、左列のパラメータと置き換えることができ(パラメータが1つの式で「制約」されている限り)、パラメータ化の選択肢を増やすことができます。真ん中の列のパラメータは、式自体に含まれていますが、現在は読み取り専用のパラメータとして選択することはできません。 3.2. 事前確率分布  事前分布では、システム内の浮動小数点値で表される任意のパラメータに対して、(DISTL10の分布オブジェクトとして)分布を作成し、添付することができます。DISTLとは、様々な確率分布の定義、変換、保存、操作を行うPythonのパッケージです。現在、PHOEBEに搭載されているバージョンでは、以下の一変量分布をサポートしています。一様(boxcar)、ガウス、デルタ、ヒストグラム、samples(入力されたサンプルセットの周りにカーネル密度推定(KDE)を生成し、サンプルと同じ分布からの描画を可能にしますが、セットから直接描画することはありません)。さらに,次のような多変量分布がサポートされています:多変量ガウシアン,多変量ヒストグラム,多変量サンプル(それぞれパラメータ間の既知の共分散を含む). これらの分布が定義され、パラメータに付加されると、これらの分布を直接プロットしたり、セクション3.1で説明した制約論理を介して伝播させることができます(図4参照)。自由」なパラメータに付けられた分布は、その値がランダムに引き出され、それぞれの分布に設定されます。任意の分布(調整可能なパラメータや「読み取り専用」のパラメータに添付されているかどうかに関わらず)は、パラメータの現在の額面値を引き当てる確率を提供することができます。事前分布をサポートする最適化アルゴリズムやサンプリングアルゴリズムに渡されると、これらの確率はメリット関数に含まれます(式(12)および(17)参照)。 この柔軟性により、多くの有用なシナリオが可能になります。例えば、aとiの両方の定義された範囲で周辺化しながら、a sin iの事前分布を設定することができます。これにより、例えば、log gやスペクトルリダクションによる質量など、データセットで直接観測できない観測上の制約があるパラメータに事前分布を設定することもできます。セクション6で説明するように、あるパラメータ化(eとω0など)でサンプリングを行い、別のパラメータ化(e sin ω0とe cos ω0)で事後確率分布と不確かさを公開することも可能です。 3.3. 角度の折り返し  角度に対応するすべてのパラメータや分布(近星点引数、昇交点の経度など)は、自動的に角度を折り返します(ただし、傾斜は昇交点の経度との曖昧さを避けるため、0からπの範囲に限定されています)。例えば、手動で、あるいはオプティマイザやサンプラで、[0, 2π]の範囲外の値を設定した場合、その値は自動的に[0, 2π]の範囲に戻されます。これは、真の値が折り返し点の近くにある場合、オプティマイザやサンプラがパラメータ空間を継続的に探索できるため、特に有効です。 これは、メリット関数がパラメータの1つに依存しない場合の最適化やサンプリングにおいて、いくつかの問題を引き起こす可能性があります。例えば、円形のシステムでは、近星点引数は観測値に影響を与えません。連続してラップさせると、サンプラーは無限にランダムウォークを続けます。これを避けるために,ラップごとの角度を初期化分布の中心値からπ以内に制限し(セクション5参照),実質的に情報を与えない事前処理を行い,サンプラーがラップされた空間をさまようのを防ぎます. 3.4. フォワードモデル フォワード・モデル」は、固定された入力パラメータ・セットに対する合成モデルで構成され、セクション3で詳細に説明したように、メリット・ファンクションを通じて観測値と比較されます。 3.4.1. 位相マスキング、露光時間、アンダーサンプリング、補間 デフォルトでは、PHOEBEは観測データのタイムスタンプを用いてフォワードモデルを計算します。しかし、PHOEBEでは計算効率を上げるために、データの特定のフェーズを除外することもできます。例えば、基本的な幾何学的パラメータや軌道パラメータの最適化を行う際には、光度曲線のうち食外の部分を除外することが有効です(図5参照)。その後、これらのデータをメリット関数に再導入して、さらなる最適化を行い、入手可能なすべての観測データから事後確率を決定する必要があります。これを可能にするために、PHOEBEは非常に基本的な「位相マスキング」を実装しており、これは有効にも無効にもできます。ただし、外れ値の除去や、より複雑なマスキング(時空間など)は、PHOEBEの外部で行う必要がありますのでご注意ください。最適化やサンプリングの際に位相マスキングを使用すると、PHOEBEはフォワードモデルとメリット関数の残差の両方からこれらの時間を除外します。 さらに、計算時間を節約するために観測データを「ダウンサンプル」したり、最適化やサンプリングを行わない場合には、プロットのために「オーバーサンプル」したり、観測データの制限外でモデルがどのように振る舞うかを予測したりすることも、しばしば有用です。これらのユースケースを実現するために、PHOEBEでは各データセットの計算時間をカスタム配列でオーバーライドすることができます。最適化やサンプリングを行う際、PHOEBEはデフォルトで2つの配列(マスクされた観測時間または要求された計算時間)のうち、短い方の配列を自動的に採用します。メリット関数の中で、合成モデルはマスクされた観測時間に補間されます。システムが時間に依存していない場合、PHOEBEは位相空間で補間することで合成モデルの時間範囲外への「外挿」を可能にしますが、システムが時間に依存している場合(例えば、時間微分がゼロでない場合(アプシダル運動の存在を含む)、またはコンポーネントの少なくとも1つにスポットがあり、非同期的に回転している場合)には、そのようなことはできません。 最後に、PHOEBEでは、要求された合成時間ごとに、露光時間にわたってモデルをオーバーサンプリングすることができます。また、PHOEBEでは、オーバーサンプリングされたモデルの露光時間に対する平均値を算出することができます。例えば、ケプラーのデータは29.4分の露光時間であるため、フォトメトリーに位相のスムージングが発生します。適切な露光時間を設定し、オーバーサンプリングを有効にすることで、この効果をモデルで再現することができます。 3.4.2. 周縁減光 PHOEBE 2.0では、線形法、対数法、二次法、平方根法、べき乗法の手動係数に加えて、デフォルトで自動的に補間された周縁減光が導入されました(Pr?a et al. PHOEBE 2.2 (Jones et al. 2020)では、表面要素ごとに、組み込まれた法則の係数を自動的に照会する機能が追加されました。 PHOEBE 2.3での最適化やサンプリングは、補間と「ルックアップ」の両方の周縁減光をサポートしていますが、現在は手動で設定した係数のフィッティングはサポートしていません。他のバックエンド(セクション3.4.3参照)では、係数はPHOEBEの大気テーブルから平均恒星パラメータに基づいて検索されます。 3.4.3. 代替の計算バックエンド PHOEBEは可能な限りロバストかつ正確に設計されていますが(具体的な実装の詳細についてはPr?a et al.2016; Horvat et al.2018; Jones et al.2020を参照)、これは計算効率を犠牲にしています。オプティマイザーやサンプラーで大量のフォワードモデルインスタンスを生成する必要がある場合、非効率性は拡大します。研究対象のシステムによっては、特定の効果を「無効」にしたり、単純化した仮定を立てたりすることで、大幅に時間を節約することができます。PHOEBEはこのような仮定をいくつかサポートしています。例えば、照射はコストがかかるため、無視できる場合は無効にしたり、球状の星がRocheで歪められた星の代わりになったりします。他のコードの方がより適しているかもしれませんし、特定のシステムに最適化されているかもしれません。PHOEBEの使用を必要とする高次の効果は、例えば、パラメータ空間全体を検索し、軌道パラメータを最適化する際には無視することができますが、最終的に報告される値を正確に推定しようとする際には、含める必要がある場合があります。 このため、PHOEBEには、PHOEBE legacy(Wilson-Devinney; Pr?a & Zwitter 2005; Wilson & Devinney 1971; Wilson 1979, 2008; Wilson & Van Hamme 2014)、ELLC(Maxted 2016)、JKTEBOP(Southworth et al.2004, 2007, 2009; Southworth 2011)など、いくつかの公開コードのフォワードモデルコンポーネントのラッパーが含まれています。これらのフォワードモデルは、PHOEBEが提供するフォワードモデルに「ドロップイン」で置き換えられるように設計されており、そのコードに特有の追加オプションは最小限に抑えられています。図6は、PHOEBEで計算した複数のフォワードモデルを、同じパラメータセットで比較したものです。各コードの前提条件が異なることや、各コードのパラメータ設定や出力を翻訳する能力に制限があることから、モデル間にはいくつかの顕著な相違が見られます。しかし,このインターフェースは,これらの比較を行い,システムに最適なモデルを選択する機能を備えています.表2は,現在利用可能な各バックエンドに実装されている機能,制限,および前提条件の概要を示しています. これらのラッパーは、PHOEBEで使用されているパラメタライズから、外部コードと呼ばれる要求されたフォワードモデルで使用されているパラメタライズに変換し、その結果をPHOEBEで使用されている正しい単位に変換します。この変換は必ずしも正確ではなく、多少の時間的コストがかかることに注意してください。しかし、どのフォワードモデルを選択しても、すべてのバックエンドはPHOEBEのフロントエンドに搭載されている機能の一部を利用することができます。例えば、「制約」による柔軟なパラメータ設定(セクション3.1参照)、位相マスキング、露光時間、時間補間(セクション3.4.1参照)、肢体暗化係数のルックアップ(セクション3.4. 2)、パスバンドの輝度や光度曲線の第3の光を利用したフラックスのスケーリング(一部はネイティブ、3.4.4項参照)、系速度や成分ごとの速度オフセット(3.4.5項参照)、観測の不確かさを過小評価するためのノイズノイジー(3.5項参照)、ガウス過程(3.6項参照)、PHOEBEに実装されているすべての最適化およびサンプリング手法などがあります。これらのバックエンドを使用したり、PHOEBEのバックエンドで物理的効果を無効にしたりする場合は、結果として得られたモデルを、これらの仮定を使用せずに計算したものと比較することをお勧めします。 これらのバックエンドを使用した論文を発表する際には、適切に外部コードを引用することをお勧めします。 3.4.4 フラックスのスケーリング(パスバンド光度、距離、第3光)  (Pr?a 2018 )で定義されたPHOEBE 2 バックエンドで直接計算された表面フラックス(形式的には等価点光源から 1m の距離で計算)は    F_pb|w = P_int|w粘_rel V μΔA    式(20) stars elemes ここで、?wはフラックスの重み付け方式(エネルギーまたは光子数)、P_int|wは通過帯域透過関数の適切な重み付け積分、"star "はシステムで定義されたコンポーネント、S_relは絶対単位と相対単位の間で変換する星ごとのスケーリングファクターです。 はモデル大気テーブルから補間された通過帯域平均された面要素あたりの射出固有強度。 Vは面要素の可視化関数、μΔAは各要素の投影表面積です。通過帯域の輝度がユーザーから提供される成分については、S_relは以下のように決定される。       L_pb,rel|w  S_rel = ------------    式(21)       L_pb,abs|w なお、これらの絶対光度および相対光度は時刻t0で定義されており、外部要因やアスペクトに依存する効果(距離、照射、スポット、パルセーション、ブースティングなど)は含まれていない。絶対的光度は、表にある絶対法線強度から求められる。 L_pb,abs|w = πP_int|w D_intΔA (22) ここで <I_norm,abs?w> はモデル大気テーブルから補間された通過帯域平均された絶対的な射出法線強度、D_int は全ての角度に渡る周縁減光関数の積分、ΔA は各要素の表面積である(Pr?a 2018)。 バージョン2.2(Jones et al.2020)以降、PHOEBEは パスバンドの光度と光度曲線のフラックスのスケーリングを扱う複数のモードをサポートしています。 1. 絶対値:大気およびパスバンドテーブルから取得した強度をそのまま使用して積分し、絶対的なフラックスを算出します。この場合、式(20)および(21)のスケーリング・ファクターSrelは、定義上、1となる。 2. データセットにスケール化型: 観測データが提供されている場合,"data setscaled "では,結果として得られる絶対的なフラックスを,最小二乗法アルゴリズムを用いてデータに合わせて自動的に絶対的なフラックスをスケーリングすることができます。この機能は、正規化された光度曲線のパラメータを光度との相関を気にすることなく最適化する際に特に有効です。しかし 便利ではありますが、サンプリングから不確かさを決定する際(セクション5と6)には、代わりに 不確かさが過小評価されるのを防ぐために、パスバンド光度でマージンを取ることをお勧めします。 3. コンポーネント結合型。一つの成分星についての通過帯光度(L_pb,rel?w) がユーザーから提供され 時刻t0における所定の固有光度が達成されるように強度が内部でスケーリングされる。次に、この同じスケーリングファクターを式(20)における両方の成分星に使用します。 4.データセット結合型。マルチバンドの測光が可能な場合。「データセット結合」では、複数の光度曲線に同じスケーリングファクターを 色情報を保持するために、複数の光度曲線で同じスケーリングファクタ あるデータセットの両成分のSrelを別のデータセットで定義されたものから採用したり あるデータセットの両成分の Srel を別のデータセットで定義されたものから採用し、L_pb,rel?を L_pb,rel?w を以下の式で計算することができます。式(21)で表されます。 5. デカップリング。デカップリングでは、星の希望する光度を 大気やパスバンドテーブルのスケーリング情報を無視して 大気やパスバンドテーブルからのスケーリング情報を無視して星の光度を個別に設定することができます。ここでは、式(20)と式(21)の各星のスケーリング係数は、互いに独立していても構わない。  ほとんどの「代替バックエンド」(3.4.3項および表2参照)では、強度のスケーリングやパスバンドの光度を直接渡すことができない場合、返されたフラックスやマグニチュードは はPHOEBEで使用する単位にリスケールする必要があります。各パスバンドの光度からパスバンド全体のフラックスを推定し、PHOEBEで使用する単位にリスケールする必要があります。 F_alt,backend,scaled L_pb,rel|w = F_alt,backend --------------  (23) 4π デカップリングモードを使用する場合、1つの星あたりの希望する光度 L_pb,rel?wは、絶対光度を計算することなく、フラックス・スケーリングを見積もるために直接に使用することができる。しかし、他のモードでは、1つ以上の相対的な通過帯域の光度を、式(21)を用いて、各成分の絶対光量を決定することで、1つ以上の相対パスバンド光量を内部で計算する必要があります。PHOEBE Legacy は、PHOEBE 2 と同様のロジックで、パスバンドの光度を入力とします。光度を入力とし、それぞれのL_pb,rel?wは直接渡されます。式(23)でリスケールする必要はありません。しかし、"decoupled", "component-coupled" (PHOEBE Legacyではネイティブにサポートされています), "data set-scaled "のいずれかを使用していない場合は、適切なカップリングを処理するために絶対光量を推定する必要があります。PHOEBE 2.3では、これらの絶対光量を推定するために2つの方法が導入されています。 1. t0でのPHOEBEメッシュ。PHOEBEメッシュはRocheモデルを用いてt0時点で作成され、与えられた大気およびパスバンドテーブルから適切な強度を入力し、表面全体を積分して固有のパスバンド光度を算出します(式(22))。メッシュの作成はPHOEBEの中でも特に計算量の多いステップであるため、この方法で光度を推定すると、これらの「代替バックエンド」を使用することによるスピード面でのメリットが大きく損なわれてしまいますが、表面の歪みが大きい場合にはより正確になります。 2. ステファン-ボルツマン近似。Requiv, Teff, log g, 含有量の平均値を用いて、選択された大気と通過帯域のテーブルに "平均 " と周縁減光モデルD_intの積分値を照会します。そして、各星の絶対的な通過帯域の光度は次のように近似されます。 L_pb,abs|w = 4πR_equiv^2D_int ×(T_eff,logg,abun) P_int|w (24) これにより、各星を一様な球体として扱うことになります。絶対光量を計算するために、各星を一様な球体として扱います。スケーリングファクター)を算出します。歪みが少ない場合や正確なスケーリングが重要でない場合は 、この近似式を使う方が メッシュを作成するよりもはるかに高速です。 いずれの方法においても、式(23)はこれらの光度から通過帯域の光束レベルを推定するための仮定であり、スケーリングされた光束は厳密なものではないことに注意する必要がある。 さらに、絶対有効温度の代わりに表面光度比を用いるバックエンドでは、計算された通過帯域の光度とそれぞれの等価半径の二乗との比として推定されます。最後に、これらのコードでネイティブに使用されている正規化は、相互に、そしてPHOEBEとは若干異なります。ELLCは露出した光度曲線を各構成要素の照射面上の積分フラックスで正規化しますが(Maxted 2016)、JKTEBOPは直交する食前の系の大きさで正規化します(J. Southworth 2020, private communication)。 これらの異なる処理がPHOEBEのリスケーリングに与える影響を図6に示します。 図6. 同一システム上でPHOEBEを介して複数のサポートされるバックエンドを実行した場合の比較(左がデタッチドシステム、右がセミドタッチドシステム)。同一ではありませんが、異なるモデル間の仮定を比較したり、特定のサイエンスケースやシステムで最も効率的なモデルを便利に活用するのに便利です。垂直方向のわずかなスケーリングの違いは、各コードの仮定と 3.4.4 節で説明したフラックス・スケーリングの仮定によるものと思われる。 これらの欠点にもかかわらず、この実装では、これらの他のコードではネイティブにサポートされていない可能性のあるパスバンドレベルの効果を簡単にサポートすることができ、より簡単にドロップイン置き換えができるようになっています。  一貫性を保つため、すべてのバックエンドから返されるフラックスには、第三光と距離の影響は含まれていません(一貫性を保つため、これらの効果をネイティブにサポートするコードであっても無効にしています)。第三光が(フラックスではなく)割合で提供されている場合は、まず各星のスケールされた光度からフラックス単位に変換します。 l3_frac L_pb,abs|w l3_flux = ---------- ------------ (25) 1-l3_frac stars 4π ここで「stars」とは、外部からの第3光の光源を含まない、システムの構成要素のことである。この方法では、式(23)の光度からフラックスへの変換に使用したのと同じ一様な球面近似を使用し、また、「デカップリング」モードを使用しない場合は、上述の2つの方法のいずれかによって絶対的な光度を推定する必要がある。  もしdata set-scaledモードを使用した場合は、フラックスは第3光を含める前の観測値に合わせてスケーリングされ、絶対光束に直接作用します(つまり L_pb,rel|w = L_pb,abs|w)に直接作用します。この場合、パスバンドの光度のスケーリングと距離は任意であるため、最小二乗法によりフラックススケールファクターS_dsを決定します。 F_synthetoc,data set-scaled = S_ds F_backend + l3_flux (26) l3_fracが提供されている場合は、l3_fluxへの変換にもこの同じスケールファクターが含まれます。 F_synthetoc,data set-scaled l3_frac L_pb,abs|w =S_ds(F_backend+---------- ------------ ) (27) 1-l3_frac stars 4π PHOEBE 2.3では、「data set-coupled」モードで結合された複数の光度曲線に対して「data setcaled」を使用できるようになりました。この場合、式(26)および式(27)が適宜使用されますが、結合されたすべてのデータセットに対して単一のS_dsを同時にフィッティングすることになります。 その他の(「data set-scaled」ではない)モードでは、バックエンドから返されるフラックスに加えて、距離と第3の光が含まれます。 F_backend F_synthetic=---------- + l3_flux   式(28) d^2 ここで、F_backendは、PHOEBE(式20)からか必要に応じて 再スケーリングの後に別のバックエンド(式23)で与えられる。 3.4.5. システム速度とコンポーネントごとの視線速度のオフセット 距離と第3光がそれぞれのバックエンドのフラックスの上に含まれているのと同様に、バックエンドから直接取得した視線速度にはオフセットが含まれていません(ネイティブにシステム速度をサポートしているバックエンドであっても)。最終的な合成視線速度は以下のように定義されます。 RV_systemic = RV_backend + RV_γ + RV_offset (29) ここで,RVγは,すべてのRV(すべてのデータセットのすべての成分)に加えられる一定の重心系速度であり,RVオフセットは,速度に対するデータセットごとおよび成分ごとのオフセットを可能にする。 このような成分ごとのオフセットは、大気の異なるレベルで作成されているためにスペクトル線のオフセットが当たり前になっている高温の星では特に有効です(Underhill & Hill 1994; Harmanec et al.2002; Shenar et al.2018)。このオフセットの大きさは、システム内の2つの成分星間だけでなく、異なるスペクトル線から縮小された視線速度データセット間でも異なる可能性があります。これをモデル自体に組み込むことで、これらのオフセットをマージン化することができ、(例えば質量比やシステム速度RVγとの)相関があれば、事後確率分布で適切に説明することができます。 3.5. ノイズの厄介者 システムの物理パラメータの事後確率分布と不確実性の決定は、観測の不確実性に大きく依存します。残念ながら、これらの不確かさは過小評価されていることが多く、そのような過小評価が偏った解や過小評価されたモデルの不確かさに伝わらないようにすることが重要である。 引用された観測の不確かさが一定のファクターで過小評価されている可能性がある場合(ケプラーデータなど;Jenkins 2017参照)、メリット関数の定義において観測の不確かさに直接不確かさのスケーリングファクターを導入することができます(式(14)および(15)も参照)。  σ^2 = σ_o^ 2 + y_m^2 e~(2σ_inf) (30) これは,事後確率分布を決定するためにサンプリングを行う際の「迷惑パラメータ」として,特に有用である(セクション5および6参照).この過小評価因子をマージン化することで、因子自体と物理モデルのパラメータとの間の縮退を、結果として得られる事後確率分布に反映させることができる。 読者の皆様には,いくつかの異なる状況における観測の不確実性の取り扱いに関する実用的な概要を示したHogg et al.(2010)を参照してください。 3.6. ガウシアンプロセス ガウス過程(GP)とは、すべての変数の結合分布がガウシアンである場合に、各データポイントが割り当てられたランダム変数から引き出されるランダム過程です。GPは、ノイズフリーのフォワードモデルで記述された天体信号の上に、系列相関や不均一性などのノイズの相関をモデル化することができるため、トレンドを除去する目的でデータに多項式などの関数を当てはめる必要がありません。PHOEBEでは、GPはフォワードモデルに適用されるため、サポートされているどの計算バックエンドでも、カーネル(隣接する点の間の共分散を記述するもの)を選択して利用することができます。PHOEBEは、計算効率が高く、多数のデータポイントに拡張可能であることを目的としたGPのPython実装であるCELERITEを介してGPを実装しています(Foreman-Mackey et al.2017)。PHOEBE 2.3では、Matern 3/2カーネルとSimple Harmonic Oscillatorカーネルの任意の組み合わせをサポートしていますが、PHOEBE内に別途実装されているノイズの迷惑パラメータとの衝突を避けるため、CELERITEがネイティブにサポートしているホワイトノイズ項("jitter")を意図的に除外しています(3.5節)。 元のフォワードモデルが計算され,(必要に応じて)観測時刻に補間された後,フォワードモデルと観測値の間の残差が求められます.これらの残差とカーネルの入力パラメータが CELERITE に渡され,CELERITE がモデルの GP 成分を提供します.その後,フォワードモデルと GP の寄与分が追加され,最終的なモデルが完成します(このモデルは,オプティマイザやサンプラの残差やメリット関数の計算に使用されます).GPは、最終的なフォワードモデルが観測の正確な時間に公開されることを要求しますが、セクション3.4.1で説明したように、物理的なフォワードモデルは要求された時間またはフェーズで計算されます。 図7は、トレンドとノイズが重畳された合成観測データを含む例を示しています。ここでは、GPを含めることで、観測データと合成フォワードモデルの間の残差を計算しています。 GPは理想的には(確率的で相関のある)ノイズ成分のみを考慮するが、その自由度の高さから、データを過剰にフィットさせ、信号の一部を考慮してしまう可能性がある。これは軌道周期におけるGPの過剰なスペクトルパワーとして現れ、結果として縮退解を引き起こし、パラメータの事後評価に影響を与えます。そのため、GPパラメータの適切な周辺化(セクション5参照)を行い、オーバーフィッティングを回避するか、少なくとも定量化することが重要です。 3.7. 並列化 PHOEBEのバックエンド自体は、MPI(Message Passing Interface)実装により時間単位で並列化されています。ここでは、1台のプロセッサが初期設定を行い、残りのプロセッサに個々のタイムスタンプを送信して、要求された合成観測値をその時刻に計算します。その後、メインプロセッサが結果をコンパイルして並べ、モデルを公開します。他のバックエンドのラッパーは、一般的にデータセット単位での並列化をサポートしており、各プロセッサがすべての時間の観測値を計算できるようになっていますが、データセットは1つだけです。しかし、この機能をサポートしているオプティマイザーやサンプラーでは、フォワードモデル内の並列化は無効化され、MPIまたはマルチプロセッシングによりモデル単位で並列化されます。この方がオーバーヘッドが少なく、効率的な場合が多いです。 ただし、MCMC(セクション5.1参照)などでは、要求されたウォーカーよりも多くのプロセッサが利用可能な場合、この方法は利点を失うことに注意が必要です。このような場合、PHOEBEは自動的にMCMCをシリアルで実行し、各反復を1回ごとに並列化するように切り替えます。 また、MPIを使用していない場合、PHOEBEはマルチプロセシングをサポートするバックエンド(現在はemceeとdynesty)の利用可能なすべてのプロセッサを対象としたマルチプロセシングプールを使用します。 4. 最適化アルゴリズム PHOEBE 2.3には、Nelder-Mead (Gao & Han 2012)、Powell、共役勾配など、scipy.optimize (Virtanen et al. 2020)の最適化アルゴリズムのラッパーが搭載されています。Nelder-Mead (Gao & Han 2012), Powell, conjugate gradientなど、これらのオプティマイザーは、推定やグローバルなパラメータ空間探索などによって得られた最適解の近傍にあるモデルの適合性を、非常に効率的に改善することができます。 PHOEBEでは、オプティマイザーによって調整されるパラメータを選択することができ、必要に応じて事前確率分布を定義することもできます。負の対数尤度(セクション3参照)を適切なSCIPYアルゴリズムに渡すことで、それぞれのアルゴリズムが呼び出されます。最適化の結果、事前確率分布が定義されている場合は最大アプリオリ(MAP)解、そうでない場合は最尤推定値が得られます。ここで重要なのは、ユーザーが事前確率分布を提供するかどうかに関わらず、パラメータリミット、ラッピングリミット、および失敗したモデルは、メリット関数にペナルティを与える情報のない事前確率分布として機能するということです。PHOEBEは、推定器と同様に(セクション2)、調整された各パラメータの提案値を、最適化前後のメリット関数の値、およびSCIPYから返された診断値とともに公開します。 図8は、いくつかのパラメータについて、推定値から提案された値でネルダー・mead最適化を行った結果です。わずか数百回の繰り返しで、残差とχ2はモデルが大幅に改善されたことを示しています。 5. サンプラー サンプラーはグローバルな解を求めるためのものではなく、ローカルなパラメータ空間の形状を探索し、様々なパラメータ間の根本的な縮退を明らかにするロバストな事後確率分布と不確実性を提供することを目的としています。 5.1. マルコフ連鎖モンテカルロ (emcee)  MCMCは,対数確率をメリット関数として用いて,サンプリングされたパラメータ間の相関を含む事後密度関数を推定するために,パラメータ空間をサンプリングするベイズ法です(セクション3参照).Affine-invariant MCMCアルゴリズム(Goodman & Weare 2010)は,作業者のアンサンブルに基づいて同じ空間をサンプリングし,一般的にバーンイン期間が短く,ユーザのチューニングも少なくて済む(ForesemanMackey et al. 2013).したがって,フォワードモデルの計算コストが高く,確率空間が直交していない可能性がある場合には,このアルゴリズムを使用するのが良い候補となる. EMCEE (Foreman-Mackey et al. 2013, 2019) Pythonパッケージは、アフィン不変MCMCのPython実装であり、この分野で広く採用されています。PHOEBEには、EMCEEのラッパーが含まれており、基本的な機能のほとんどを、メリット関数と並列化の問題を処理しながら、サンプリングされたパラメータを簡単に追加または削除したり、事前分布と初期サンプルに複雑な分布を使用したりする利便性とともに公開しています。しかし、他のアルゴリズムと同様に、EMCEE自体を理解し、どのように入力を調整し、どのように結果を解釈するかが重要となります。Foreman-Mackey et al. (2013)が概要と考察を述べており、EMCEEのオンラインドキュメント(https://emcee.readthedocs.io)にはいくつかの有用な例が掲載されています。 PHOEBEでEMCEEを実行する際、ユーザーは初期サンプル(どのパラメータをサンプルに含めるかを定義する)に使用する分布を任意のパラメータ化で定義することができ、オプションで事前分布の分布を定義することもできます。事前分布が提供されない場合でも、個々のパラメータの限界と角度ラッピングの限界(3.3節参照)に基づいて、情報のない分布が採用されます。さらに、個々のモデルが何らかの理由で不合格になった場合、メリット関数(セクション3参照)は-∞を返します。これは本質的に、非物理的なシステムやバックエンド自体の制限に起因する、許容されるパラメータ空間に対する優先順位として機能します。図9に示すように、これらの拒否されたサンプルの原因は、オプションでユーザーに公開することができます。 初期サンプルと事前分布は同じ分布のセットであっても構いませんが、一般的には、最適化によって得られた最良の解を中心に、N次元のハイパーボール(すなわち、複数のパラメータに対するガウス分布を介して)でEMCEEを開始するのが良い方法です(セクション4参照)。また,事前分布自体から直接サンプリングすることも可能です(ただし,一般的には効率が悪く,局所解に陥りやすいので,この点についてはForeman-Mackey et al.) 一方,事前分布には,外部分析や文献から得られた不確かさを表す情報量の多いものと,明らかに不正確な,あるいは非物理的なパラメータ空間に鎖が「迷い込む」のを防ぐための情報量の少ない保守的な一様分布があり,セクション3で説明したようにメリット関数に含まれています.   分布に加えて,プロセッサの数,ウォーカーの数,反復回数などのオプションを設定することができます.実行が完了した後(あるいはユーザの要求に応じて中間ステップで)、ユーザは連鎖の進行状況や、反復回数に対する対数確率を確認し、必要に応じてバーンインや間引きのパラメータを調整することができます。ランが完了してもまだ収束していない場合、PHOEBEでは既存のチェインからEMCEEのランを継続することができます。    フルチェインに加えて、アクセプタンスフラクションと自己相関時間もユーザーに公開されます。Foreman-Mackeyら(2013)は、経験則として、受け入れフラクションはすべて0.2から0.5の間にし、自己相関時間をおよそ10回カバーするのに十分な反復を行うことを提案しています。PHOEBEは、チェーンの自己相関時間に基づいて、間引きとバーンインのデフォルト値を自動的に決定します。   burnin = F_burnin max(τ_autocorr) (31) そして Thin = F_thin Min(τ_autocorr) (32) ここで,F_burnin および F_thin はユーザ定義のスケーリングファクタで,デフォルトはそれぞれ 2 および 0.5,τ_autocorr はEMCEEで公開されているパラメータごとの推定自己相関時間である. 鎖がパラメータ空間を十分にカバーし、収束したと判断されると、ユーザーはEMCEEから返された鎖に適用する間引き、バーンイン、または対数確率のカットオフを調整し、結果として得られる事後確率分布の多変量分布を生成することができる(詳細はセクション6を参照)。    5.2. 動的入れ子サンプリング(Dynesty)  入れ子サンプリングは、パラメータ空間をサンプリングするもう一つのベイズ法です。初期分布からサンプリングし、パラメータ空間を歩き、メリット関数の中の事前分布によって「ペナルティ」を受けるMCMCとは異なり、ネステッドサンプリングは、スライスに分割された事前分布自体から継続的に直接サンプリングする(詳細な比較と議論についてはSpeagle 2020を参照)。Nested samplingは、MCMCとは異なり、マルチモーダルな事後分布を探索し、公開することができます。しかし、情報を与えない事前確率分布が広い場合、ネステッド・サンプリングは法外に高価になる可能性があります。Dynamic nested samplingは、アルゴリズムが最終的な事後情報に収束する際に、反復ごとにサンプル数を動的に変更することで、このコストを幾分軽減します。さらに,MCMCでは,事前分布が許す限り初期のサンプリング分布の外側をさまようことができますが,動的入れ子サンプリングでは,事前分布によって定義された元のパラメータ空間の外側の解を考慮することはありません。  DYNESTYコードは、Pythonのダイナミックネステッドサンプリングパッケージです。EMCEEと同様に、ダイナミックネステッドサンプリングとDYNESTYコードの詳細と複雑さを最初に理解することは非常に有益です。Speagle (2020)では詳細な議論がなされており、オンラインドキュメント(https://dynesty.readthedocs.io)では多数の例が紹介されています。    DYNESTYは事前分布から直接サンプリングしているため、PHOEBEは対数尤度から事前分布を除外しています。その代わりに、PHOEBEはDISTL(セクション3.2参照)を使って、定義された事前分布を、0から1の間のランダムな値を各パラメータの描画値にマッピングする事前変換に変換し、これをDYNESTYに渡します。この変換では、事前分布のパラメータ間の共分散を無視する必要があることに注意してください(つまり、事前分布が多変量分布として提供された場合、DYNESTYに渡す前にまず一変量分布に平坦化されます)。PHOEBE 内から DYNESTY を呼び出すと、ユーザーは DYNESTY からの出力配列すべてにアクセスでき、診断図にアクセスしたりプロットしたりすることができます。DYNESTYからのサンプルは、それぞれの重みを考慮して事後分布に変換することができます。 6. 事後確率分布と不確実性 サンプラーで得られた結果は事後分布として得られます。これは、サンプリングされた各パラメータ値の不確かさだけでなく、サンプリングされたすべてのパラメータ間の相関(または共分散)も表しています。これらの分布が十分にガウス分布である場合,事後分布は容易に多変量ガウス分布 に変換することができます.これは,平均値と共分散行列だけで事後分布を表すものです.図10は、EMCEEサンプルから得られたいくつかのパラメータの後縁と、多変量ガウスへの変換のコーナープロットの例です。 これらの分布はDISTLの分布オブジェクトであるため、サポートされている任意の分布タイプに簡単に操作することができ、ファイルに保存したり、プロットしたり、対称または非対称の不確実性を公開して発表結果に含めたり、さらなるサンプリングに使用したりすることができます。例えば,あるサンプリング実行のバックエンドは,別のサンプリング実行の初期サンプリング分布または事前分布として採用し,使用することができます。これは、1つのデータセットのバックエンドを最適化して決定する際に便利ですが、別のデータセットや新しい観測結果が得られた場合には、その情報を対数尤度で維持することができます。事前分布に使用される分布(セクション3.2)と同様に、これらの事後分布は「制約」を介して伝搬され、任意のパラメータ化(セクション3.1)で公開されます。例えば、図11はEMCEEの離心率と周回軌道の引数のサンプルを示していますが、これをesin ω0とe cos ω0に伝搬させることで、多変量ガウシアンとして表現しやすくなり、科学的な成果として発表することができます。これにより、サンプリング、事前分布、事後分布のパラメータの組み合わせを独立して選択することができ、利用可能な観測データや既知の情報に基づいて必要な柔軟性を最大限に発揮することができます。 注意すべき点は、サンプラーによって決定された縮退度と不確かさは、残りのすべてのパラメータが不確かさなしで表向きの値に固定されているという前提で計算されているため、不確かさが過小評価されたり、場合によっては誤った解が得られたりすることです。例えば、軌道周期の値がよく知られていると考えられ、サンプリング時に固定されている場合、サンプリングされた他のすべてのパラメータの事後分布、つまり不確かさは、軌道周期に関する縮退が除外されるため、過小評価される可能性が高くなります。同様に、仮定した固定値が間違っていた場合、残りのサンプルパラメータも間違っている可能性があります。つまり、サンプルされていないパラメータは、事前にデルタ関数を持っていて、無限の精度で知られているかのように扱われます。理想的な世界では、すべてのパラメータがサンプリングされるか、またはマージンが取られますが、実際には、数十から数百の可能なパラメータがあるため、それは実行不可能です。 最後に、他の分布と同様に、PHOEBEではフォワードモデル自体に事後分布を伝播させることができます。フォワードモデルを事後分布からの多数のサンプルに対して実行することで、PHOEBEはこれらの個々のモデルや中央値モデル、フラックスや視線速度の1σ、3σ、5σの広がりを計算し、公開することができます。図12は、EMCEEで得られた事後分布から合成モデルを抽出した際の3σスプレッドの例です。 7. ユーザーインターフェース PHOEBE 2.3 リリースでは、基本的なプロット機能を含め、上記の前進モデルと逆問題の両方の機能をすべて使用できるユーザーインターフェースも導入され、http://phoebeproject.org/clients からダウンロードできます。ユーザーインターフェースはサーバー・クライアントモデルで設計されており、コードベースを3つの異なる役割に分割します:Pythonパッケージ自体、ユーザーインターフェースクライアント、そして同じマシンかクライアントからリモートで実行可能でPythonパッケージを介してユーザーインターフェースクライアントから送られるコマンドの実行を担当する軽量なウェブサーバーです。この設計により、より多くのリソースを持つリモートマシンにサーバーをインストールし、ユーザーインターフェースをローカルにインストールするか、Webブラウザーを通じてホストされたバージョンにアクセスするか、柔軟に選択できるようになりました。 ユーザーインターフェースの有無にかかわらず、PHOEBE では順モデルまたは逆問題のジョブをスタンドアローンのスクリプトにエクスポートし、それを高性能クラスタで実行して、結果のファイルをローカルにインポートして進捗状況や最終結果を確認することもできます。これにより、フィットプロセスのすべてのインタラクティブおよびビジュアルコンポーネントをローカルマシンで実行しながら、計算コストの高いタスクを、ウェブサーバーをホストしたりグラフィックスをプロットしたりすることができないリモートマシンにオフロードすることができます。 8. まとめ PHOEBEの2.3リリースでは、過去のリリースを基に、逆問題に対処するためにいくつかの一般的なアルゴリズムを実行するデータフィットの一般的なフレームワークを導入し、複雑な分布を定義して扱うためのインターフェイスを備えています。また、PHOEBEとWebベースのユーザーインターフェイスに加えて、他のいくつかの一般に公開されているコードを使用してフォワードモデルを実行する機能も導入されています。 ブラックボックス化された自動フィッティングパイプラインには程遠いものの、この共通インターフェースは、食連星系の正確で精密なパラメータ値と不確かさの両方を得ようとする際に、複数のアルゴリズムとコードを用いるのに必要な学習曲線と労力を軽減することを目的としています。今後、新たなアルゴリズムやフォワードモデルが採用された場合、同じ枠組みでPHOEBEに取り込んでいく予定です。 PHOEBEはGPL3ライセンスのオープンソースプロジェクトで、http://phoebe-project.org と https://github.com/phoebe-project/phoebe2 でホストされています。貢献とフィードバックを歓迎します。 PHOEBEの開発はNSFAAGグラント1517474と1909109、NASA 17-ADAP17-68によって可能となりました、このことに感謝します。 John SouthworthとPierre Maxtedの許可と彼らのコード、JKTEBOPとELLCに関する議論に感謝します D.J. は、スペイン科学・イノベーション・大学省(MCIU)の国家研究機関(AEI)と欧州地域開発基金(FEDER)から、助成金 AYA2017-83383-P の支援を得ていることを認めます。D.J.はまた、スペイン科学・イノベーション・大学省から一般国家予算に課された資金と、経済産業・貿易・知識省によるカナリア諸島自治区の一般予算から移された資金で賄われた助成金P/308614の支援に謝意を表します。