21st-Century Statistical Computation for Exoplanet Studies Tom - PowerPoint PPT Presentation

21st-Century Statistical Computation for Exoplanet Studies Tom Loredo Dept. of Astronomy, Cornell University http://www.astro.cornell.edu/staff/loredo/ With David Chernoff, Merlise Clyde, Jim Berger, Floyd Bullard, Jim Crooks Work very much in progress MaxEnt — July 9, 2007

A Productive Methodological Union ,-.*&/&0$'%&1'".( !"#$%&"'()*"*&%*&+% 2'*$34"*&1' !"#$%$"&$'()$&*+*,"(-.$,%/'( 345124%$'(6,"2$(31%7,'( !"#,%012*,"(-.$,%/ 636386$2%,9,7*+:;1+2*"<+=== !""#$%&'()*%+$),'-.).$+.$/+

An Unproductive Sociological Divide !"#$%&'(&*,$& !"#$%&'(&)*"*'%*'+% -,#%'+".&)+'$(+$% !"#$%&'"()*+%&$*&,-.*//-%& 3"<#"7*444%&?*;;9*+)444% 0"1"2*%&3-.$#*+444 @8A%&?"+.*)%&BC##%&0D-##-.2% 5678.8(*/9-7)%&:-8-.;89("/-7)%& )87-8#82+%&<)+7=8#82+444> E9*//=89)/444 !"#$%&'()$"*+,)-",-

An Example Cosmological Parameters from WMAP & SDSS Tegmark et al. 2004 • Looked at over 2 dozen models with 2 – 10 params • Used random walk MCMC in “rotated” params • Typical run length 3 – 5 × 10 5 ; a few runs trapped • 30 CPU years of effort

Outline 1 Introducing Exoplanets 2 Bayesian Methods for Exoplanets 3 Some 21st Century Algorithms A Bayesian pipeline Improved MCMC for parameter estimation Marginal likelihood calculation

Outline 1 Introducing Exoplanets 2 Bayesian Methods for Exoplanets 3 Some 21st Century Algorithms A Bayesian pipeline Improved MCMC for parameter estimation Marginal likelihood calculation

Extrasolar Planets Exoplanet detection/measurement methods: • Direct: Transits, gravitational lensing, imaging, interferometric nulling • Indirect: Keplerian reflex motion (line-of-sight velocity, astrometric wobble) The Sun’s Wobble From 10 pc 1000 500 0 -500 -1000 -1000 -500 0 500 1000

Radial Velocity Technique As of July 2007, 245 planets found, including 26 multiple-planet systems Vast majority (233) found via Doppler radial velocity (RV) measurements

Parameters for an Orbit — Single Planet Intrinsic geometry: semimajor axis a , eccentricity e Orientation: 3 Euler angles, i , ω , Ω Time: period τ , origin M p RV parameters: semi-amplitude K ( a , e , τ ), τ , e , M p , ω , COM velocity v 0 Physics: min. mass m sin i , size a , eccentricity e (Ultimate goal: multiple planets, astrometry → ∼ 30+ parameters!)

Conventional RV Data Analysis Analysis method: Identify best candidate period via periodogram; fit parameters with nonlinear least squares !"#"$%&%'"( )"#"*&$' +",-."-"#"*&%*"/01 2"#"*&%3"45" 6-,7)8")9"2:&"%**'

Issues with Periodogram + Least Squares • Multimodality, nonlinearity, sparse data → “∆ χ 2 ” uncertainties not valid • Reporting uncertainties in derived parameters ( m sin i , a ) • Lomb-Scargle periodogram not optimal for eccentric and multiple planets • Handling marginal detections • Combining info from many systems for pop’n studies • Scheduling future observations Bayesian statistics can address all these within a unified framework!

Outline 1 Introducing Exoplanets 2 Bayesian Methods for Exoplanets 3 Some 21st Century Algorithms A Bayesian pipeline Improved MCMC for parameter estimation Marginal likelihood calculation

Scientific Goals ⇐ ⇒ Bayesian Methodology • Improved system inferences — Detection (marginal likelihoods, Bayes factors), estimation (joint & marginal posteriors) • Improved population inferences — Propagate system uncertainties (including detection uncertainties and selection effects) to population level (hierarchical/multi-level Bayes) • Adaptive scheduling — Plan future observations to most quickly reduce uncertainties (Bayesian experimental design) Prior information & data Combined information Inference Strategy New data Predictions Strategy Observ’n Design Interim results

Keplerian Radial Velocity Model Parameters for single planet • Argument of pericenter ω • τ = orbital period (days) • Mean anomaly of pericenter • e = orbital eccentricity passage M p • K = velocity amplitude (m/s) • System center-of-mass velocity v 0 Velocity vs. time v ( t ) = v 0 + K ( e cos ω + cos[ ω + υ ( t )]) True anomaly υ ( t ) found via Kepler’s equation for eccentric anomaly: � 1 / 2 E ( t ) − e sin E ( t ) = 2 π t tan υ � 1 + e tan E − M p ; 2 = τ 1 − e 2 A strongly nonlinear model!

The Likelihood Function Keplerian velocity model with parameters θ = { K , τ, e , M p , ω, v 0 } : d i = v ( t i ; θ ) + ǫ i For measurement errors with std dev’n σ i , and additional “jitter” with std dev’n σ J , L ( θ, σ J ) ≡ p ( { d i }| θ, σ J ) N [ d i − v ( t i ; θ )] 2 1 � − 1 � � = exp � 2 σ 2 i + σ 2 σ 2 i + σ 2 2 π J i =1 J   � � 1 − 1 �  exp 2 χ 2 ( θ ) ∝ � σ 2 i + σ 2 2 π i J [ d i − v ( t i ; θ )] 2 � χ 2 ( θ, σ J ) ≡ where σ 2 i + σ 2 J i Ignore jitter for now . . .

What To Do With It Parameter estimation Posterior dist’n for parameters of model M i with i planets: p ( θ | D , M ) ∝ p ( θ | M ) L i ( θ ) Summarize with mode, means, credible regions (found by integrating over θ ) Detection Calculate probability for no planets ( M 0 ), one planet ( M 1 ) . . . . Let I = { M i } . p ( M i | D , I ) ∝ p ( M i | I ) L ( M i ) � where L ( M i ) = d θ p ( θ | M i ) L ( θ ) Marginal likelihood L ( M i ) includes “Occam factor”

Design Predict future datum d t at time t , accounting for model uncertainties: � p ( d t | D , M i ) = d θ p ( d t | θ, M i ) p ( θ | D , M i ) . 1st factor is Gaussian for d t with known model; 2nd term & integral account for uncertainty. Information theory → best time has largest d t uncertainty. Population analysis Multi-level (“hierarchical”) model: Parameterize the orbit parameter prior and infer it. Include probabilities for star having (1, 2, . . . ) planets → Bayes factors weight system contributions. Need all data, well-characterized surveys.

Current MCMC Methods Random Walk Metropolis Parallel Tempering (Ford) (Gregory) General purpose, but Performance depends a lot inefficient/unwieldy on parameterization

Outline 1 Introducing Exoplanets 2 Bayesian Methods for Exoplanets 3 Some 21st Century Algorithms A Bayesian pipeline Improved MCMC for parameter estimation Marginal likelihood calculation

Improving Exoplanet Bayesian Calculation Seek more efficient and less unwieldy samplers → abandon “general purpose” algorithms Keep all tasks in mind: Estimation, detection, design, pop’n . . . Guidelines for modest dimensional Bayesian computation: • Know thine enemy — Understand scales in various dimensions, multimodality, locate modes; significant exploration should precede & guide sampling • Reduce dimensionality — Analytic/numerical techniques may help • Clever priors — If a clever prior (that doesn’t do too much violence to real prior info) helps, use it; you can weight and resample/Rao-Blackwellize later • Try different samplers — Nearly 20 years of ideas

Know Thine Enemy I: Linear vs. Nonlinear Parameters v ( t ) = v 0 + K ( e cos ω + cos[ ω + υ ( t )]) = A 1 + A 2 [ e + cos υ ( t )] + A 3 sin υ ( t ) A 1 = v 0 ; A 2 = K cos ω ; A 3 = − Ksin ω Model is linear wrt A i , nonlinear wrt e , M p , ω → p ( θ ) = p ( τ, e , M p ) p ( { A i }| τ, e , M p ) || D , M 1 with p ( { A i }| τ, e , M p ) conditionally normal if we adopt a flat or conjugate prior for { A i } . Flat prior → p ( K | M 1 ) ∝ K . Using this interim prior , dimension can be reduced by 3.

Know Thine Enemy II — Slices Extreme multimodality in τ ; challenging multimodality in M p ; smooth in e

Periodogram Connections Bayesian periodograms (Jaynes-Bretthorst algorithm) Data are superposition of periodic functions + noise: M � d i = A α g α ( t i ; ω ; θ ) + e i α =1 Calculate L ( { A } , ω, θ ) using χ 2 . Integrate out A ’s → least squares + volume factors: − r 2 ( ω, θ ) � � p ( ω, θ | D ) ∝ p ( ω, θ ) J ( ω, θ ) exp 2 r 2 ( ω, θ ) is squared residual, found by diagonalizing metric η αβ = g α · g β Integrate out θ numerically → p ( ω | D ); S ( ω ) ≡ ln [ p ( ω | D )] Generalizes Schuster periodogram & LSP.

Circular Orbits Assume circular orbits : θ = { K , τ, φ, v 0 } Frequentist For given τ , maximize likelihood over K and φ (set v 0 to data mean, ¯ v ) → profile likelihood: log L p ( τ, ¯ v ) ∝ Lomb-Scargle periodogram Bayesian For given τ , integrate (“marginalize”) likelihood × prior over K and φ (set v 0 to data mean, ¯ v ) → marginal posterior: log p ( τ, ¯ v | D ) ∝ Lomb-Scargle periodogram Additionally marginalize over v 0 → floating-mean LSP

Radial Kepler Periodogram for Eccentric Orbits v ( t ) = A 1 + A 2 [ e + cos υ ( t )] + A 3 sin υ ( t ) For given ( τ, e , M p ), analytically marginalize { A i } : ln p ( τ, e , M p | D , M p ) = (Radial) Kepler-gram For given τ , marginalize over ( e , M p ): log p ( τ | D ) ∝ (Radial) Kepler periodogram This requires a 2-d numerical integral.