Spectral Timing, Bayesian Inference and You

Spectral timing has emerged as a key way to study energetic sources: simultaneously considering both time and energy information allows us to understand relationships between the different physical components in a system, such as different regions in an accretion flow.
These developments are coming at a time where technology is enabling rapid improvements in our ability to infer knowledge from data. In this talk, I’m aiming to build a bridge between spectral timing in astronomy on one side, and recent, exciting developments in statistics and computer science on the other. I will present these in the context of current work on mitigating dead time effects on timing studies with neural network-based density estimation. I will also introduce the spectral timing package stingray, and give a sneak peak of what’s in store for current and near-future developments of the software.

Transcript

1. Spectral Timing, Bayesian Inference and You Daniela Huppenkothen SRON Netherlands Institute for Space Research ! Tiana_Athriel " dhuppenkothen ! d.huppenkothen@sron.nl
2. These slides: https://speakerdeck.com/ dhuppenkothen/spectral-timing- bayesian-inference-and-you Documentation: https://stingray.science Tutorials: https://github.com/stingraysoftware/ notebooks Iconic 60s Theme Song: https://youtu.be/sgkk-MMif-4?t=50
3. How do black holes accrete matter?
4. General Relativity plasma physics + magnetic fields jet physics + particle acceleration stellar winds
5. Time [s] Frequenc Time Malzac 2008 Photon Energy [keV] Photon Energy * Flux [detector units] “soft state” “hard state” high-frequency/low-frequency X-ray brightness total brightness Spectral States Credit: Sera Markoff
6. Yamaoka et al (2010); Huppenkothen et al. (2019) Variability: GX 339-4 Fast Fourier Transform stochastic variability (red noise) quasi-periodic oscillation (QPO) white noise
7. X-ray reverberation around accreting black holes 5 Fig. 3 Time lag (8–13 keV relative to 2–4 keV) versus frequency for a hard state obser- vation of Cyg X-1 obtained by RXTE in December 1996. The trend can be very roughly approximated with a power-law of slope −0.7, but note the clear step-like features, which correspond roughly to different Lorentzian features in the power spectrum (Nowak 2000). Cygnus X-1: Nowak, 2000 1 10 0.5 2 5 0.5 Energy (keV) Fig. 9 The ratio spectrum of 1H0707-495 to a continuum model (Fabian et al. 2009). T broad iron K and iron L band are clearly evident in the data. The origin of the soft exce below 1 keV in this source had been debatable, but in this work was found to be dominat by relativistically broadened emission lines. 10−5 10−4 10−3 0.01 −50 0 50 100 150 200 Lag (s) Temporal Frequency (Hz) Fig. 10 The frequency-dependent lags in 1H0707-495 between the continuum dominat hard band at 1–4 keV and the reflection dominated soft band at 0.3–1 keV. and found significant high-frequency soft lags in 15 sources. Plotting the am 1H0707-495: Uttley et al, 2014 Time Lags Temporal variations and energy spectra are intricately linked
8. Problem: infer accretion physics from spectral timing representations statistics! stingray!
9. Why Stingray? (I needed some software to do timing)
10. • Limited functionality • Last updated 2009 (?) • Designed for spectroscopy • has some timing capabilities • https://space.mit.edu/CXC/isis/ + closed-source packages maintained by individual groups
12. https://stingray.science • 3 lead developers/maintainers (Huppenkothen, Bachetti, Stevens) • ~10 contributors • 5 completed Google Summer of Code Projects • astropy-affiliated project • connecting functionality to astropy.TimeSeries and lightkurve Huppenkothen et al (2019)
13. stingray.science Huppenkothen et al (2019) • spectral timing classes + functions • building blocks • simulation, pulsars, modeling • command-line interface built on stingray • quick-look (spectral) timing analysis • graphical user interface • interactive data analysis
14. Top-level functionality Events Lightcurve Powerspectrum / Crossspectrum AveragedPowerspectrum / AveragedCrossspectrum Crosscorrelation / Autocorrelation Bispectrum Lag-Energy Spectra supporting functionality: statistics, good time intervals, I/O Sub-modules pulse modeling simulator deadtime
15. Where to get help https://stingray.science https://github.com/stingraysoftware join the slack!
16. Current + Future Work • rework the modeling interface + autodiff: current GSoC project • multi-taper periodograms: current GSoC project • update to DAVE GUI: current GSOC project • fix bugs (there are always bugs) • improve API (aka: what were we thinking?!) • improve documentation (there is never enough documentation) • performance + memory optimization • better integration with current X-ray missions • better integration with spectral modeling packages • better integration with astropy.timeseries and lightkurve • higher-order Fourier products
17. Get Involved! Stingray is a project for the community, but there’s only so much we as maintainers can do …
18. How to get involved (please do!) • find bugs (and report them as a GitHub Issue) • fix bugs (as a GitHub Pull Request) • make feature requests (also via GitHub Issue) • implement new features (also via GitHub Pull Request) • test documentation/tutorials (and report mistakes/fix bugs etc) • … Don’t know where to start? • “Good First Issue” tag on GitHub • join the slack + ask us! We’ll help :)
19. https://github.com/StingraySoftware/stingray
20. Please cite the ApJ paper if you use stingray! https://github.com/StingraySoftware/stingray
21. Statistical Inference
22. All models are wrong, but some are useful. — George Box
23. Nature is complex!
24. … so is our data (collection)!
25. Energy [keV] Brightness Corona reflection Corona Iron Line
26. Maximum Likelihood Estimation “How probable is it that my data D came from a model M with parameters ?”
27. log(p(D|θ)) = log ( n ∏ i=1 p(d i |θ) ) = − n 2 log(2πσ2) − 1 2σ2 n ∑ i=1 (d i − f(x i , θ))2
28. log(p(D|θ)) = log ( n ∏ i=1 p(d i |θ) ) = − n 2 log(2πσ2) − 1 2σ2 n ∑ i=1 (d i − f(x i , θ))2 “chi-square fitting” caution: this is not a probability of the parameters
29. When this might not be enough … • your likelihood is not unimodal • your likelihood might be skewed or otherwise complex • you might have useful prior information • you might be interested in the properties of a population • you might have a complex model/data collection that you can simulate, but not write down an analytical function for
30. Complex/Weird Likelihoods
31. Useful Prior Information data parameters θ = D = p(D|θ) = n ∏ i=1 p(d i |θ) p(θ|D) ∝ p(θ) n ∏ i=1 p(d i |θ) p(θ)
32. Useful Prior Information data parameters θ = D = p(D|θ) = n ∏ i=1 p(d i |θ) p(θ|D) ∝ p(θ) n ∏ i=1 p(d i |θ) p(θ)
33. Useful Prior Information data parameters θ = D = p(D|θ) = n ∏ i=1 p(d i |θ) p(θ|D) ∝ p(θ) n ∏ i=1 p(d i |θ) p(θ) You can do this with stingray.modeling
34. Simulation-Based Inference
35. The Chandra ABC Guide to Pile-Up
36. The Chandra ABC Guide to Pile-Up There are many effects or models that you can (and should) simulate, but that are hard to take into account in model fitting.
37. Yamaoka et al (2010); Huppenkothen et al. (2019) Variability: GX 339-4 Fast Fourier Transform We understand how to do this, unless the sources is very bright
38. NuSTAR Credit: NuSTAR Observatory Guide Adapted from Chaplin et al (2012) dead time ms
39. without dead time with dead time
40. Huppenkothen et al (2017) Figure 5. Left panel: averaged periodogram of the part of Chandra observation 17696 containing the QPOs at 73 mHz and 1.03 H periodogram of the two Fermi/GBM triggers simultaneous with the Chandra data (right panel). In blue, we show the logarithmically data sets, we show the MAP model with four (Chandra) or two ( Fermi/GBM) Lorentzian components in purple and the combined m The Astrophysical Journal, 834:90 (17pp), 2017 January 1 Problem: (Bayesian) parameter inference
41. p(θ|D) ∝ p(D|θ)p(θ) data D D θ θ θ parameters θ = D = likelihood prior posterior intractable! ☹
42. *Also known as: Approximate Bayesian Computation Likelihood-Free Inference simulation-based inference*: replace an intractable likelihood by a (physics) simulator
43. Sadegh + Vrugt, 2014, see also: Brehmer et al (2018a, 2018b), Cranmer et al (2020) Step 1: draw parameters from prior Step 2: simulate data sets Step 3: compare simulated to observed data Step 4: keep parameters that produce simulations similar to the data simulator
44. Tejero-Cantero et al (2020), Greenberg et al (2019) (Sequential) Neural Posterior Estimation
45. Accurate Timing with SBI 7 Simulated Data Huppenkothen & Bachetti (under review)
46. Figure 17. Marginalized posterior probability distributions for the red noise Huppenkothen & Bachetti (under review)
47. GRS 1915+105 Huppenkothen & Bachetti (under review)
48. Huppenkothen & Bachetti (under review)
49. GRS 1915+105 Huppenkothen & Bachetti (under review)
50. What are the properties of the evolution of the QPO centroid frequency? p(θ|D) ∝ n ∏ i=1 p(D i |θ)p(θ) p({θ}n i=1 , α|D) ∝ n ∏ i=1 [p(D i |θ i )p(θ i |α)] p(α)
51. What are the properties of the evolution of the QPO centroid frequency? p(θ|D) ∝ n ∏ i=1 p(D i |θ)p(θ) p({θ}n i=1 , α|D) ∝ n ∏ i=1 [p(D i |θ i )p(θ i |α)] p(α) Hierarchical Bayesian Modeling is ideal for population-level problems
52. Things I’m curious about • Where else should be be using simulation-based inference? • What do you think are the next big challenges in spectral timing? • What population-level problems should we be addressing with hierarchical inference?
53. • Use stingray! It’s fun! Conclusions • Also, help us build stingray! That’s also fun! • Bayesian statistics allows you to go beyond standard fitting problems to take into account • SBI enables principled inference when models are complex, numerical and/or generate data • New statistical methods and software tools will set us up for answering complex questions with complex models, using current and future instrumentation
54. Questions? Comments? Complaints?