Estimating variability power spectra of active galaxies from irregular time series

A common feature of Active Galactic Nuclei (AGN) is their random variations in brightness across the whole emission spectrum, from radio to $\gamma$-rays. Studying the nature and origin of these fluctuations is critical to characterising the underlying variability process of the accretion flow that powers AGN. Random timing fluctuations are often studied with the power spectrum; this quantifies how the amplitude of variations is distributed over temporal frequencies. Red noise variability -- when the power spectrum increases smoothly towards low frequencies -- is ubiquitous in AGN. The commonly used Fourier analysis methods, have significant challenges when applied to arbitrarily sampled light curves of red noise variability. Several time-domain methods exist to infer the power spectral shape in the case of irregular sampling but they suffer from biases which can be difficult to mitigate, or are computationally expensive. In this paper, we demonstrate a method infer the shape of broad-band power spectra for irregular time series, using a Gaussian process regression method scalable to large datasets. The power spectrum is modelled as a power-law model with one or two bends with flexible slopes. The method is fully Bayesian and we demonstrate its utility using simulated light curves. Finally, Ark 564, a well-known variable Seyfert 1 galaxy, is used as a test case and we find consistent results with the literature using independent X-ray data from XMM-Newton and Swift. We provide publicly available, documented and tested implementations in Python and Julia.