Thomas J. Loredo

Cosmic demographics—the statistical study of populations of astrophysical objects—has long relied on tools from multivariate statistics for analyzing data comprising fixed-length vectors of properties of objects, as might be compiled in a tabular astronomical catalog (say, with sky coordinates, and brightness measurements in a fixed number of spectral passbands). But beginning with the emergence of automated digital sky surveys, ca. 2000, astronomers began producing large collections of data with more complex structures: light curves (brightness time series) and spectra (brightness vs. wavelength). These comprise what statisticians call functional data—measurements of populations of functions. Upcoming automated sky surveys will soon provide astronomers with a flood of functional data. New methods are needed to accurately and optimally analyze large ensembles of light curves and spectra, accumulating information both along individual measured functions and across a population of such functions. Functional data analysis (FDA) provides tools for statistical modeling of functional data. Astronomical data presents several challenges for FDA methodology, e.g., sparse, irregular, and asynchronous sampling, and heteroscedastic measurement error. Bayesian FDA uses hierarchical Bayesian models for function populations, and is well suited to addressing these challenges. We provide an overview of astronomical functional data and some key Bayesian FDA modeling approaches, including functional mixed effects models, and stochastic process models. We briefly describe a Bayesian FDA framework combining FDA and machine learning methods to build low-dimensional parametric models for galaxy spectra.

Bayesian inference gets its name from Bayes’s theorem, expressing posterior probabilities for hypotheses about a data generating process as the (normalized) product of prior probabilities and a likelihood function. But Bayesian inference uses all of probability theory, not just Bayes’s theorem. Many hypotheses of scientific interest are composite hypotheses, with the strength of evidence for the hypothesis dependent on knowledge about auxiliary factors, such as the values of nuisance parameters (e.g., uncertain background rates or calibration factors). Many important capabilities of Bayesian methods arise from use of the law of total probability, which instructs analysts to compute probabilities for composite hypotheses by marginalization over auxiliary factors. This tutorial targets relative newcomers to Bayesian inference, aiming to complement tutorials that focus on Bayes’s theorem and how priors modulate likelihoods. The emphasis here is on marginalization over parameter spaces—both how it is the foundation for important capabilities, and how it may motivate caution when parameter spaces are large. Topics covered include the difference between likelihood and probability, understanding the impact of priors beyond merely shifting the maximum likelihood estimate, and the role of marginalization in accounting for uncertainty in nuisance parameters, systematic error, and model misspecification.

Cosmic demographics -- the statistical study of populations of astrophysical objects -- has long relied on *multivariate statistics*, providing methods for analyzing data comprising fixed-length vectors of properties of objects, as might be compiled in a tabular astronomical catalog (say, with sky coordinates, and brightness measurements in a fixed number of spectral passbands). But beginning with the emergence of automated digital sky surveys, ca. ~2000, astronomers began producing large collections of data with more complex structure: light curves (brightness time series) and spectra (brightness vs. wavelength). These comprise what statisticians call *functional data* -- measurements of populations of functions. Upcoming automated sky surveys will soon provide astronomers with a flood of functional data. New methods are needed to accurately and optimally analyze large ensembles of light curves and spectra, accumulating information both along and across measured functions. Functional data analysis (FDA) provides tools for statistical modeling of functional data. Astronomical data presents several challenges for FDA methodology, e.g., sparse, irregular, and asynchronous sampling, and heteroscedastic measurement error. Bayesian FDA uses hierarchical Bayesian models for function populations, and is well suited to addressing these challenges. We provide an overview of astronomical functional data, and of some key Bayesian FDA modeling approaches, including functional mixed effects models, and stochastic process models. We briefly describe a Bayesian FDA framework combining FDA and machine learning methods to build low-dimensional parametric models for galaxy spectra.

Survey-based measurements of the spectral energy distributions (SEDs) of galaxies have flux density estimates on badly misaligned grids in rest-frame wavelength. The shift to rest frame wavelength also causes estimated SEDs to have differing support. For many galaxies, there are sizeable wavelength regions with missing data. Finally, dim galaxies dominate typical samples and have noisy SED measurements, many near the limiting signal-to-noise level of the survey. These limitations of SED measurements shifted to the rest frame complicate downstream analysis tasks, particularly tasks requiring computation of functionals (e.g., weighted integrals) of the SEDs, such as synthetic photometry, quantifying SED similarity, and using SED measurements for photometric redshift estimation. We describe a hierarchical Bayesian framework, drawing on tools from functional data analysis, that models SEDs as a random superposition of smooth continuum basis functions (B-splines) and line features, comprising a finite-rank, nonstationary Gaussian process, measured with additive Gaussian noise. We apply this *Splines 'n Lines* (SnL) model to a collection of 678,239 galaxy SED measurements comprising the Main Galaxy Sample from the Sloan Digital Sky Survey, Data Release 17, demonstrating capability to provide continuous estimated SEDs that reliably denoise, interpolate, and extrapolate, with quantified uncertainty, including the ability to predict line features where there is missing data by leveraging correlations between line features and the entire continuum.

This and companion papers by Harrington et al. and Blecic et al. present the Bayesian Atmospheric Radiative Transfer (bart) code, an open-source, open-development package to characterize extrasolar planet atmospheres. bart combines a thermochemical equilibrium abundance (tea), a radiative transfer (Transit), and a Bayesian statistical (mc3) module to constrain atmospheric temperatures and molecular abundances for given spectroscopic observations. Here we describe the Transit radiative transfer package, an efficient line-by-line radiative transfer C code for one-dimensional atmospheres, developed by P. Rojo and further modified by the UCF exoplanet group. This code produces transmission and hemisphere-integrated emission spectra. Transit handles line-by-line opacities from HITRAN, Partridge & Schwenke (H₂O), Schwenke (TiO), and Plez (VO) and collision-induced absorption from Borysow, HITRAN, and ExoMol. Transit emission spectra models agree with models from C. Morley (private communication) within a few percent. We applied bart to the Spitzer and Hubble transit observations of the Neptune-sized planet HAT-P-11b. Our analysis of the combined HST and Spitzer data generally agrees with those from previous studies, finding atmospheric models with enhanced metallicity (≳100× solar) and high-altitude clouds (≲1 mbar level). When analyzing only the HST data, our models favor high-metallicity atmospheres, in contrast with the previous analysis by Chachan et al. We suspect that this discrepancy arises from the different choice of chemistry modeling (free constant-with-altitude versus thermochemical equilibrium) and the enhanced parameter correlations found when neglecting the Spitzer observations. The bart source code and documentation are available at https://github.com/exosports/BART.

This and companion papers by Harrington et al. and Cubillos et al. describe an open-source retrieval framework, Bayesian Atmospheric Radiative Transfer (BART), available to the community under the reproducible-research license via https://github.com/exosports/BART. BART is a radiative transfer code (transit; https://github.com/exosports/transit; Rojo et al.), initialized by the Thermochemical Equilibrium Abundances (TEA; https://github.com/dzesmin/TEA) code (Blecic et al.), and driven through the parameter phase space by a differential-evolution Markov Chain Monte Carlo (MC3; https://github.com/pcubillos/mc3) sampler (Cubillos et al.). In this paper we give a brief description of the framework and its modules that can be used separately for other scientific purposes; outline the retrieval analysis flow; present the initialization routines, describing in detail the atmospheric profile generator and the temperature and species parameterizations; and specify the post-processing routines and outputs, concentrating on the spectrum band integrator, the best-fit model selection, and the contribution functions. We also present an atmospheric analysis of WASP-43b secondary eclipse data obtained from space- and ground-based observations. We compare our results with the results from the literature and investigate how the inclusion of additional opacity sources influences the best-fit model.

We present the open-source Bayesian Atmospheric Radiative Transfer (BART) retrieval package, which produces estimates and uncertainties for an atmosphere’s thermal profile and chemical abundances from observations. Several BART components are also stand-alone packages, including the parallel Multi-Core Markov-chain Monte Carlo (MC3), which implements several Bayesian samplers; a line-by-line radiative-transfer model, transit; a code that calculates Thermochemical Equilibrium Abundances (TEA), and a test suite for verifying radiative-transfer and retrieval codes, BARTTest. The codes are in Python and C. BART and TEA are under a Reproducible Research (RR) license, which requires reviewed-paper authors to publish a compendium of all inputs, codes, and outputs supporting the paper’s scientific claims. BART and TEA produce the compendium’s content. Otherwise, these codes are under permissive open-source terms, as are MC3 and BARTTest, for any purpose. This paper presents an overview of the code, BARTTest, and an application to eclipse data for exoplanet HD 189733b. Appendices address RR methodology for accelerating science, a reporting checklist for retrieval papers, the spectral resolution required for synthetic tests, and a derivation of the effective sample size required to estimate any Bayesian posterior distribution to a given precision, which determines how many iterations to run. Paper II, by Cubillos et al., presents the underlying radiative-transfer scheme and an application to transit data for exoplanet HAT-P-11b. Paper III, by Blecic et al., discusses the initialization and post-processing routines, with an application to eclipse data for exoplanet WASP-43b. We invite the community to use and improve BART and its components at http://GitHub.com/ExOSPORTS/BART/.

We present the Infrared Database of Extragalactic Observables from Spitzer (IDEOS), a homogeneous, publicly available, database of 77 fitted mid-infrared observables in the 5.4–36 μm range, comprising measurements for 3335 galaxies observed in the low-resolution staring mode of the Infrared Spectrometer on board the Spitzer Space Telescope. Among the included observables are polycyclic aromatic hydrocarbon fluxes and their equivalent widths, the strength of the 9.8 μm silicate feature, emission-line fluxes, solid-state features, rest-frame continuum fluxes, synthetic photometry, and a mid-infrared spectral classification. The IDEOS spectra were selected from the Cornell Atlas of Spitzer-IRS Sources. To our surprise, we have detected at a >95% confidence level crystalline silicates in the spectra of 786 IDEOS galaxies. The detections range from single-band detections to detections of all fitted crystalline bands (16, 19, 23, 28, and 33 μm). We find the strength of the crystalline silicate bands to correlate with the amorphous silicate strength and the change from an emission to an absorption feature to occur at higher obscuration as the wavelength of the crystalline silicate band is longer. These observed characteristics are consistent with an origin for the amorphous and crystalline silicate features in a centrally heated dust geometry, either an edge-on disk or a cocoon. We find the 23 and 33 μm crystalline silicate bands to be well suited to classify the obscuration level of galactic nuclei, even in the presence of strong circumnuclear star formation. Based on our detection statistics, we conclude that crystalline silicates are a common component of the interstellar medium of galactic nuclei.

The break-by-one gamma distribution has a probability density function resembling the Schechter function, but with the small-argument behavior modified so it is normalizable in commonly arising cases where the Schechter function is not. Its connection to the gamma distribution makes it straightforward to sample from. These properties make it useful for cosmic demographics.