Search

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., Comment: 35 pages, 11 figures; accepted for publication in Frontiers in Astronomy and Space Sciences (special issue for iid2022: Statistical Methods for Event Data - Illuminating the Dynamic Universe); fixed minor typo

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., Comment: 20 pages, 7 figures

We describe a computational framework for hierarchical Bayesian inference with simple (typically single-plate) parametric graphical models that uses graphics processing units (GPUs) to accelerate computations, enabling deployment on very large datasets. Its C++ implementation, CUDAHM (CUDA for Hierarchical Models) exploits conditional independence between instances of a plate, facilitating massively parallel exploration of the replication parameter space using the single instruction, multiple data architecture of GPUs. It provides support for constructing Metropolis-within-Gibbs samplers that iterate between GPU-accelerated robust adaptive Metropolis sampling of plate-level parameters conditional on upper-level parameters, and Metropolis-Hastings sampling of upper-level parameters on the host processor conditional on the GPU results. CUDAHM is motivated by demographic problems in astronomy, where density estimation and linear and nonlinear regression problems must be addressed for populations of thousands to millions of objects whose features are measured with possibly complex uncertainties. We describe a thinned latent point process framework for modeling such demographic data. We demonstrate accurate GPU-accelerated parametric conditional density deconvolution for simulated populations of up to 300,000 objects in ~1 hour using a single NVIDIA Tesla K40c GPU. Supplementary material provides details about the CUDAHM API and the demonstration problem., Comment: 28 pages, 7 figures, 2 appendices

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 abundances (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 ({\water}), Schwenke (TiO), and Plez (VO); and collision-induced absorption from Borysow, HITRAN, and ExoMol. Transit emission-spectra models agree with models from C. Morley (priv. comm.) within a few percent. We applied BART to the {\Spitzer} and {\Hubble} transit observations of the Neptune-sized planet HAT-P-11b. Our results generally agree with those from previous studies, constraining the {\water} abundance and finding an atmosphere enhanced in heavy elements. Different conclusions start to emerge when we make different assumptions from other studies. The BART source code and documentation are available at https://github.com/exosports/BART., Comment: Accepted for publication at PSJ

This and companion papers by Harrington et al. 2021, submitted and Cubillos et al. 2021, submitted 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 2009, 2009ASPC..420..321R), initialized by the Thermochemical Equilibrium Abundances (TEA, https://github.com/dzesmin/TEA , Blecic et al. 2016, arXiv:1505.06392) code, and driven through the parameter phase space by a differential-evolution Markov-chain Monte Carlo (MC3, https://github.com/pcubillos/mc3 , Cubillos et al. 2017, arXiv:1610.01336) sampler. 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 influence the best-fit model., Comment: 16 pages, 14 figures

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 189733 b. 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/., Comment: 38 pages, 17 figures, submitted to The Planetary Science Journal

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., Comment: 4 pages, 1 figure

Demographic studies of cosmic populations must contend with measurement errors and selection effects. We survey some of the key ideas astronomers have developed to deal with these complications, in the context of galaxy surveys and the literature on corrections for Malmquist and Eddington bias. From the perspective of modern statistics, such corrections arise naturally in the context of multilevel models, particularly in Bayesian treatments of such models: hierarchical Bayesian models. We survey some key lessons from hierarchical Bayesian modeling, including shrinkage estimation, which is closely related to traditional corrections devised by astronomers. We describe a framework for hierarchical Bayesian modeling of cosmic populations, tailored to features of astronomical surveys that are not typical of surveys in other disciplines. This thinned latent marked point process framework accounts for the tie between selection (detection) and measurement in astronomical surveys, treating selection and measurement error effects in a self-consistent manner., Comment: 33 pages, 5 figures

This Astro2020 State of the Profession Consideration White Paper highlights the growth of astrostatistics and astroinformatics in astronomy, identifies key issues hampering the maturation of these new subfields, and makes recommendations for structural improvements at different levels that, if acted upon, will make significant positive impacts across astronomy., Comment: 14 pages, 1 figure; submitted to the Decadal Survey on Astronomy and Astrophysics (Astro2020) on 10 July 2019; see https://sites.nationalacademies.org/DEPS/Astro2020/DEPS_192906

The radial velocity method is one of the most successful techniques for detecting exoplanets. It works by detecting the velocity of a host star induced by the gravitational effect of an orbiting planet, specifically the velocity along our line of sight, which is called the radial velocity of the star. Low-mass planets typically cause their host star to move with radial velocities of 1 m/s or less. By analyzing a time series of stellar spectra from a host star, modern astronomical instruments can in theory detect such planets. However, in practice, intrinsic stellar variability (e.g., star spots, convective motion, pulsations) affects the spectra and often mimics a radial velocity signal. This signal contamination makes it difficult to reliably detect low-mass planets. A principled approach to recovering planet radial velocity signals in the presence of stellar activity was proposed by Rajpaul et al. (2015). It uses a multivariate Gaussian process model to jointly capture time series of the apparent radial velocity and multiple indicators of stellar activity. We build on this work in two ways: (i) we propose using dimension reduction techniques to construct new high-information stellar activity indicators; and (ii) we extend the Rajpaul et al. (2015) model to a larger class of models and use a power-based model comparison procedure to select the best model. Despite significant interest in exoplanets, previous efforts have not performed large-scale stellar activity model selection or attempted to evaluate models based on planet detection power. In the case of main sequence G2V stars, we find that our method substantially improves planet detection power compared to previous state-of-the-art approaches., Comment: 37 pages, 7 figures

We present {\em Spitzer} secondary-eclipse observations of the hot Jupiter HAT-P-13 b in the 3.6 {\micron} and 4.5 {\micron} bands. HAT-P-13 b inhabits a two-planet system with a configuration that enables constraints on the planet's second Love number, \math{k\sb{2}}, from precise eccentricity measurements, which in turn constrains models of the planet's interior structure. We exploit the direct measurements of \math{e \cos \omega} from our secondary-eclipse data and combine them with previously published radial velocity data to generate a refined model of the planet's orbit and thus an improved estimate on the possible interval for \math{k\sb{2}}. We report eclipse phases of \math{0.49154 \pm 0.00080} and \math{0.49711 \pm 0.00083} and corresponding \math{e \cos \omega} estimates of \math{-0.0136 \pm 0.0013} and \math{-0.0048 \pm 0.0013}. Under the assumptions of previous work, our estimate of \math{k\sb{2}} of 0.81 {\pm} 0.10 is consistent with the lower extremes of possible core masses found by previous models, including models with no solid core. This anomalous result challenges both interior models and the dynamical assumptions that enable them, including the essential assumption of apsidal alignment. We also report eclipse depths of 0.081\% {\pm} 0.008\% in the 3.6 {\micron} channel and 0.088 \% {\pm} 0.028 \% in the 4.5 {\micron} channel. These photometric results are non-uniquely consistent with solar-abundance composition without any thermal inversion., Comment: Accepted for publication by The Astrophysical Journal December 12, 2016

Observational astronomy in the time-domain era faces several new challenges. One of them is the efficient use of observations obtained at multiple epochs. The work presented here addresses faint object detection with multi-epoch data, and describes an incremental strategy for separating real objects from artifacts in ongoing surveys, in situations where the single-epoch data are summaries of the full image data, such as single-epoch catalogs of flux and direction estimates for candidate sources. The basic idea is to produce low-threshold single-epoch catalogs, and use a probabilistic approach to accumulate catalog information across epochs; this is in contrast to more conventional strategies based on co-added or stacked image data across all epochs. We adopt a Bayesian approach, addressing object detection by calculating the marginal likelihoods for hypotheses asserting there is no object, or one object, in a small image patch containing at most one cataloged source at each epoch. The object-present hypothesis interprets the sources in a patch at different epochs as arising from a genuine object; the no-object (noise) hypothesis interprets candidate sources as spurious, arising from noise peaks. We study the detection probability for constant-flux objects in a simplified Gaussian noise setting, comparing results based on single exposures and stacked exposures to results based on a series of single-epoch catalog summaries. Computing the detection probability based on catalog data amounts to generalized cross-matching: it is the product of a factor accounting for matching of the estimated fluxes of candidate sources, and a factor accounting for matching of their estimated directions. We find that probabilistic fusion of multi-epoch catalog information can detect sources with only modest sacrifice in sensitivity and selectivity compared to stacking., Comment: 11 pages, 11 figures

Modern astronomy increasingly relies upon systematic surveys, whose dedicated telescopes continuously observe the sky across varied wavelength ranges of the electromagnetic spectrum; some surveys also observe non-electromagnetic "messengers," such as high-energy particles or gravitational waves. Stars and galaxies look different through the eyes of different instruments, and their independent measurements have to be carefully combined to provide a complete, sound picture of the multicolor and eventful universe. The association of an object's independent detections is, however, a difficult problem scientifically, computationally, and statistically, raising varied challenges across diverse astronomical applications. The fundamental problem is finding records in survey databases with directions that match to within the direction uncertainties. Such astronomical versions of the record linkage problem are known by various terms in astronomy: cross-matching, cross-identification, and directional, positional, or spatio-temporal coincidence assessment. Astronomers have developed several statistical approaches for such problems, largely independently of related developments in other disciplines. Here we review emerging approaches that compute (Bayesian) probabilities for the hypotheses of interest: possible associations, or demographic properties of a cosmic population that depend on identifying associations. Many cross-identification tasks can be formulated within a hierarchical Bayesian partition model framework, with components that explicitly account for astrophysical effects (e.g., source brightness vs. wavelength, source motion, or source extent), selection effects, and measurement error. We survey recent developments, and highlight important open areas for future research., Comment: 21 pages, 9 figures

Time-correlated noise is a significant source of uncertainty when modeling exoplanet light-curve data. A correct assessment of correlated noise is fundamental to determine the true statistical significance of our findings. Here we review three of the most widely used correlated-noise estimators in the exoplanet field, the time-averaging, residual-permutation, and wavelet-likelihood methods. We argue that the residual-permutation method is unsound in estimating the uncertainty of parameter estimates. We thus recommend to refrain from this method altogether. We characterize the behavior of the time averaging's rms-vs.-bin-size curves at bin sizes similar to the total observation duration, which may lead to underestimated uncertainties. For the wavelet-likelihood method, we note errors in the published equations and provide a list of corrections. We further assess the performance of these techniques by injecting and retrieving eclipse signals into synthetic and real Spitzer light curves, analyzing the results in terms of the relative-accuracy and coverage-fraction statistics. Both the time-averaging and wavelet-likelihood methods significantly improve the estimate of the eclipse depth over a white-noise analysis (a Markov-chain Monte Carlo exploration assuming uncorrelated noise). However, the corrections are not perfect, when retrieving the eclipse depth from Spitzer datasets, these methods covered the true (injected) depth within the 68\% credible region in only $\sim$45--65\% of the trials. Lastly, we present our open-source model-fitting tool, Multi-Core Markov-Chain Monte Carlo ({MC$^3$}). This package uses Bayesian statistics to estimate the best-fitting values and the credible regions for the parameters for a (user-provided) model. {MC$^3$} is a Python/C code, available at https://github.com/pcubillos/MCcubed., Comment: Accepted for publication in ApJ

The Second Workshop on Extreme Precision Radial Velocities defined circa 2015 the state of the art Doppler precision and identified the critical path challenges for reaching 10 cm/s measurement precision. The presentations and discussion of key issues for instrumentation and data analysis and the workshop recommendations for achieving this precision are summarized here. Beginning with the HARPS spectrograph, technological advances for precision radial velocity measurements have focused on building extremely stable instruments. To reach still higher precision, future spectrometers will need to produce even higher fidelity spectra. This should be possible with improved environmental control, greater stability in the illumination of the spectrometer optics, better detectors, more precise wavelength calibration, and broader bandwidth spectra. Key data analysis challenges for the precision radial velocity community include distinguishing center of mass Keplerian motion from photospheric velocities, and the proper treatment of telluric contamination. Success here is coupled to the instrument design, but also requires the implementation of robust statistical and modeling techniques. Center of mass velocities produce Doppler shifts that affect every line identically, while photospheric velocities produce line profile asymmetries with wavelength and temporal dependencies that are different from Keplerian signals. Exoplanets are an important subfield of astronomy and there has been an impressive rate of discovery over the past two decades. Higher precision radial velocity measurements are required to serve as a discovery technique for potentially habitable worlds and to characterize detections from transit missions. The future of exoplanet science has very different trajectories depending on the precision that can ultimately be achieved with Doppler measurements., Comment: 45 pages, 23 Figures, workshop summary proceedings

A preliminary study of a set of well-isolated pulses in GRB light curves indicates that simple pulse models, with smooth and monotonic pulse rise and decay regions, are inadequate. Examining the residuals of fits of pulses to such models suggests the following patterns of departure from the smooth pulse model of Norris et al. (2005): A Precursor Shelf occurs prior to or concurrent with the exponential Rapid Rise. The pulse reaches maximum intensity at the Peak Plateau, then undergoes a Rapid Decay. The decay changes into an Extended Tail. Pulses are almost universally characterized by hard-to-soft evolution, arguing that the new pulse features reflect a single physical phenomenon, rather than artifacts of pulse overlap., Comment: 6 pages, 3 figures. 7th Huntsville Gamma-Ray Burst Symposium, GRB 2013: paper 24 in eConf Proceedings C1304143

The behavior of merging black holes (including the emitted gravitational waves and the properties of the remnant) can currently be computed only by numerical simulations. This paper introduces ten numerical relativity simulations of binary black holes with equal masses and equal spins aligned or anti-aligned with the orbital angular momentum. The initial spin magnitudes have $|\chi_i| \lesssim 0.95$ and are more concentrated in the aligned direction because of the greater astrophysical interest of this case. We combine this data with five previously reported simulations of the same configuration, but with different spin magnitudes, including the highest spin simulated to date, $\chi_i \approx 0.97$. This data set is sufficiently accurate to enable us to offer improved analytic fitting formulae for the final spin and for the energy radiated by gravitational waves as a function of initial spin. The improved fitting formulae can help to improve our understanding of the properties of binary black hole merger remnants and can be used to enhance future approximate waveforms for gravitational wave searches, such as Effective-One-Body waveforms., Comment: 13 pages, 10 figures

This perspective chapter briefly surveys: (1) past growth in the use of Bayesian methods in astrophysics; (2) current misconceptions about both frequentist and Bayesian statistical inference that hinder wider adoption of Bayesian methods by astronomers; and (3) multilevel (hierarchical) Bayesian modeling as a major future direction for research in Bayesian astrostatistics, exemplified in part by presentations at the first ISI invited session on astrostatistics, commemorated in this volume. It closes with an intentionally provocative recommendation for astronomical survey data reporting, motivated by the multilevel Bayesian perspective on modeling cosmic populations: that astronomers cease producing catalogs of estimated fluxes and other source properties from surveys. Instead, summaries of likelihood functions (or marginal likelihood functions) for source properties should be reported (not posterior probability density functions), including nontrivial summaries (not simply upper limits) for candidate objects that do not pass traditional detection thresholds., Comment: 27 pp, 4 figures. A lightly revised version of a chapter in "Astrostatistical Challenges for the New Astronomy" (Joseph M. Hilbe, ed., Springer, New York, forthcoming in 2012), the inaugural volume for the Springer Series in Astrostatistics. Version 2 has minor clarifications and an additional reference

This paper summarizes a presentation for a panel discussion on "The Future of Astrostatistics" held at the Statistical Challenges in Modern Astronomy V conference at Pennsylvania State University in June 2011. I argue that the emerging needs of astrostatistics may both motivate and benefit from fundamental developments in statistics. I highlight some recent work within statistics on fundamental topics relevant to astrostatistical practice, including the Bayesian/frequentist debate (and ideas for a synthesis), multilevel models, and multiple testing. As an important direction for future work in statistics, I emphasize that astronomers need a statistical framework that explicitly supports unfolding chains of discovery, with acquisition, cataloging, and modeling of data not seen as isolated tasks, but rather as parts of an ongoing, integrated sequence of analyses, with information and uncertainty propagating forward and backward through the chain. A prototypical example is surveying of astronomical populations, where source detection, demographic modeling, and the design of survey instruments and strategies all interact., Comment: 8 pp, 2 figures. To appear in "Statistical Challenges in Modern Astronomy V," (Lecture Notes in Statistics, Vol. 209), ed. Eric D. Feigelson and G. Jogesh Babu; publication planned for Sep 2012; see http://www.springer.com/statistics/book/978-1-4614-3519-8

This paper is an invited commentary on Tamas Budavari's presentation, "On statistical cross-identification in astronomy," for the Statistical Challenges in Modern Astronomy V conference held at Pennsylvania State University in June 2011. I begin with a brief review of previous work on probabilistic (Bayesian) assessment of directional and spatio-temporal coincidences in astronomy (e.g., cross-matching or cross-identification of objects across multiple catalogs). Then I discuss an open issue in the recent innovative work of Budavari and his colleagues on large-scale probabilistic cross-identification: how to assign prior probabilities that play an important role in the analysis. With a simple toy problem, I show how Bayesian multilevel modeling (hierarchical Bayes) provides a principled framework that justifies and generalizes pragmatic rules of thumb that have been successfully used by Budavari's team to assign priors., Comment: 6 pp, no figures. To appear in "Statistical Challenges in Modern Astronomy V," (Lecture Notes in Statistics, Vol. 209), ed. Eric D. Feigelson and G. Jogesh Babu; publication planned for July 2012; see http://www.springer.com/statistics/book/978-1-4614-3519-8

I describe ongoing work developing Bayesian methods for flexible modeling of arrival time series data without binning, aiming to improve detection and measurement of X-ray and gamma-ray pulsars, and of pulses in gamma-ray bursts. The methods use parametric and semiparametric Poisson point process models for the event rate, and by design have close connections to conventional frequentist methods currently used in time-domain astronomy., Comment: 4 pages, 1 figure; to appear in the proceedings of IAU Symposium 285, "New Horizons in Time Domain Astronomy" (proceedings eds. Elizabeth Griffin, Bob Hanisch, and Rob Seaman), Cambridge University Press; see http://www.physics.ox.ac.uk/IAUS285/

Exoplanet WASP-14b is a highly irradiated, transiting hot Jupiter. Joshi et al. calculate an equilibrium temperature Teq of 1866 K for zero albedo and reemission from the entire planet, a mass of 7.3 +/- 0.5 Jupiter masses and a radius of 1.28 +/- 0.08 Jupiter radii. Its mean density of 4.6 g/cm3 is one of the highest known for planets with periods less than 3 days. We obtained three secondary eclipse light curves with the Spitzer Space Telescope. The eclipse depths from the best jointly fit model are $0.224\%$ +/- $0.018\%$ at 4.5 {\mu}m and $0.181\%$ +/- $0.022\%$ at 8.0 {\mu}m. The corresponding brightness temperatures are 2212 +/- 94 K and 1590 +/- 116 K. A slight ambiguity between systematic models suggests a conservative 3.6 {\mu}m eclipse depth of $0.19\%$ +/- $0.01\%$ and brightness temperature of 2242 +/- 55 K. Although extremely irradiated, WASP-14b does not show any distinct evidence of a thermal inversion. In addition, the present data nominally favor models with day night energy redistribution less than $~30\%$. The current data are generally consistent with oxygen-rich as well as carbon-rich compositions, although an oxygen-rich composition provides a marginally better fit. We confirm a significant eccentricity of e = 0.087 +/- 0.002 and refine other orbital parameters., Comment: 16 pages, 16 figures

The dayside of HD 149026b is near the edge of detectability by the Spitzer Space Telescope. We report on eleven secondary-eclipse events at 3.6, 4.5, 3 x 5.8, 4 x 8.0, and 2 x 16 microns plus three primary-transit events at 8.0 microns. The eclipse depths from jointly-fit models at each wavelength are 0.040 +/- 0.003% at 3.6 microns, 0.034 +/- 0.006% at 4.5 microns, 0.044 +/- 0.010% at 5.8 microns, 0.052 +/- 0.006% at 8.0 microns, and 0.085 +/- 0.032% at 16 microns. Multiple observations at the longer wavelengths improved eclipse-depth signal-to-noise ratios by up to a factor of two and improved estimates of the planet-to-star radius ratio (Rp/Rs = 0.0518 +/- 0.0006). We also identify no significant deviations from a circular orbit and, using this model, report an improved period of 2.8758916 +/- 0.0000014 days. Chemical-equilibrium models find no indication of a temperature inversion in the dayside atmosphere of HD 149026b. Our best-fit model favors large amounts of CO and CO2, moderate heat redistribution (f=0.5), and a strongly enhanced metallicity. These analyses use BiLinearly-Interpolated Subpixel Sensitivity (BLISS) mapping, a new technique to model two position-dependent systematics (intrapixel variability and pixelation) by mapping the pixel surface at high resolution. BLISS mapping outperforms previous methods in both speed and goodness of fit. We also present an orthogonalization technique for linearly-correlated parameters that accelerates the convergence of Markov chains that employ the Metropolis random walk sampler. The electronic supplement contains light-curve files and supplementary figures., Comment: Accepted for publication in ApJ

We describe work in progress by a collaboration of astronomers and statisticians developing a suite of Bayesian data analysis tools for extrasolar planet (exoplanet) detection, planetary orbit estimation, and adaptive scheduling of observations. Our work addresses analysis of stellar reflex motion data, where a planet is detected by observing the "wobble" of its host star as it responds to the gravitational tug of the orbiting planet. Newtonian mechanics specifies an analytical model for the resulting time series, but it is strongly nonlinear, yielding complex, multimodal likelihood functions; it is even more complex when multiple planets are present. The parameter spaces range in size from few-dimensional to dozens of dimensions, depending on the number of planets in the system, and the type of motion measured (line-of-sight velocity, or position on the sky). Since orbits are periodic, Bayesian generalizations of periodogram methods facilitate the analysis. This relies on the model being linearly separable, enabling partial analytical marginalization, reducing the dimension of the parameter space. Subsequent analysis uses adaptive Markov chain Monte Carlo methods and adaptive importance sampling to perform the integrals required for both inference (planet detection and orbit measurement), and information-maximizing sequential design (for adaptive scheduling of observations). We present an overview of our current techniques and highlight directions being explored by ongoing research., Comment: 29 pages, 11 figures. An abridged version is accepted for publication in Statistical Methodology for a special issue on astrostatistics, with selected (refereed) papers presented at the Astronomical Data Analysis Conference (ADA VI) held in Monastir, Tunisia, in May 2010. Update corrects equation (3)

I describe ongoing work on development of Bayesian methods for exploring periodically varying phenomena in astronomy, addressing two classes of sources: pulsars, and extrasolar planets (exoplanets). For pulsars, the methods aim to detect and measure periodically varying signals in data consisting of photon arrival times, modeled as non-homogeneous Poisson point processes. For exoplanets, the methods address detection and estimation of planetary orbits using observations of the reflex motion "wobble" of a host star, including adaptive scheduling of observations to optimize inferences., Comment: 28 pages, 6 figures. Invited discussion paper for the Ninth Valencia International Conference on Bayesian Statistics (discussant was Peter Mueller); to appear in "Bayesian Statistics 9," ed. by Jos\'e M. Bernardo et al., Oxford University Press (2011); see http://ukcatalogue.oup.com/product/9780199694587.do

The transiting exoplanet WASP-18b was discovered in 2008 by the Wide Angle Search for Planets (WASP) project. The Spitzer Exoplanet Target of Opportunity Program observed secondary eclipses of WASP-18b using Spitzer's Infrared Array Camera (IRAC) in the 3.6 micron and 5.8 micron bands on 2008 December 20, and in the 4.5 micron and 8.0 micron bands on 2008 December 24. We report eclipse depths of 0.30 +/- 0.02%, 0.39 +/- 0.02%, 0.37 +/- 0.03%, 0.41 +/- 0.02%, and brightness temperatures of 3100 +/- 90, 3310 +/- 130, 3080 +/- 140 and 3120 +/- 110 K in order of increasing wavelength. WASP-18b is one of the hottest planets yet discovered - as hot as an M-class star. The planet's pressure-temperature profile most likely features a thermal inversion. The observations also require WASP-18b to have near-zero albedo and almost no redistribution of energy from the day-side to the night side of the planet., Comment: 10 pages, 13 figures. Electronic supplements to appear with published paper

We observed two secondary eclipses of the exoplanet WASP-12b using the Infrared Array Camera on the Spitzer Space Telescope. The close proximity of WASP-12b to its G-type star results in extreme tidal forces capable of inducing apsidal precession with a period as short as a few decades. This precession would be measurable if the orbit had a significant eccentricity, leading to an estimate of the tidal Love number and an assessment of the degree of central concentration in the planetary interior. An initial ground-based secondary eclipse phase reported by Lopez-Morales et al. (0.510 +/- 0.002) implied eccentricity at the 4.5 sigma level. The spectroscopic orbit of Hebb et al. has eccentricity 0.049 +/- 0.015, a 3 sigma result, implying an eclipse phase of 0.509 +/- 0.007. However, there is a well documented tendency of spectroscopic data to overestimate small eccentricities. Our eclipse phases are 0.5010 +/- 0.0006 (3.6 and 5.8 microns) and 0.5006 +/- 0.0007 (4.5 and 8.0 microns). An unlikely orbital precession scenario invoking an alignment of the orbit during the Spitzer observations could have explained this apparent discrepancy, but the final eclipse phase of Lopez-Morales et al. (0.510 -0.006 / +0.007) is consistent with a circular orbit at better than 2 sigma. An orbit fit to all the available transit, eclipse, and radial-velocity data indicates precession at <1 sigma; a non-precessing solution fits better. We also comment on analysis and reporting for Spitzer exoplanet data in light of recent re-analyses., Comment: 12 pages, 8 figures, electronic supplements to appear with published paper

We are pleased to present a Special Section on Statistics and Astronomy in this issue of the The Annals of Applied Statistics. Astronomy is an observational rather than experimental science; as a result, astronomical data sets both small and large present particularly challenging problems to analysts who must make the best of whatever the sky offers their instruments. The resulting statistical problems have enormous diversity. In one problem, one may have to carefully quantify uncertainty in a hard-won, sparse data set; in another, the sheer volume of data may forbid a formally optimal analysis, requiring judicious balancing of model sophistication, approximations, and clever algorithms. Often the data bear a complex relationship to the underlying phenomenon producing them, much in the manner of inverse problems., Comment: Published in at http://dx.doi.org/10.1214/09-AOAS234 the Annals of Applied Statistics (http://www.imstat.org/aoas/) by the Institute of Mathematical Statistics (http://www.imstat.org)

I describe a framework for adaptive scientific exploration based on iterating an Observation--Inference--Design cycle that allows adjustment of hypotheses and observing protocols in response to the results of observation on-the-fly, as data are gathered. The framework uses a unified Bayesian methodology for the inference and design stages: Bayesian inference to quantify what we have learned from the available data and predict future data, and Bayesian decision theory to identify which new observations would teach us the most. When the goal of the experiment is simply to make inferences, the framework identifies a computationally efficient iterative ``maximum entropy sampling'' strategy as the optimal strategy in settings where the noise statistics are independent of signal properties. Results of applying the method to two ``toy'' problems with simulated data--measuring the orbit of an extrasolar planet, and locating a hidden one-dimensional object--show the approach can significantly improve observational efficiency in settings that have well-defined nonlinear models. I conclude with a list of open issues that must be addressed to make Bayesian adaptive exploration a practical and reliable tool for optimizing scientific exploration., Comment: 17 pages, 5 figures

I discuss an issue arising in analyzing data from astronomical surveys: accounting for measurement uncertainties in the properties of individual sources detected in a survey when making inferences about the entire population of sources. Source uncertainties require the analyst to introduce unknown ``incidental'' parameters for each source. The number of parameters thus grows with the size of the sample, and standard theorems guaranteeing asymptotic convergence of maximum likelihood estimates fail in such settings. From the Bayesian point of view, the missing ingredient in such analyses is accounting for the volume in the incidental parameter space via marginalization. I use simple simulations, motivated by modeling the distribution of trans-Neptunian objects surveyed in the outer solar system, to study the effects of source uncertainties on inferences. The simulations show that current non-Bayesian methods for handling source uncertainties (ignoring them, or using an ad hoc incidental parameter integration) produce incorrect inferences, with errors that grow more severe with increasing sample size. In contrast, accounting for source uncertainty via marginalization leads to sound inferences for any sample size., Comment: 12 pages, 5 figures; to appear in Bayesian Inference And Maximum Entropy Methods In Science And Engineering: 24th International Workshop, Garching, Germany, 2004; ed. Volker Dose et al. (AIP Conference Proceedings Series)

We analyze the Time-Tagged Event (TTE) data from observations of gamma ray bursts (GRBs) and soft gamma repeaters (SGRs) by the Burst and Transient Source Experiment (BATSE). These data provide the best available time resolution for GRBs and SGRs. We have performed an extensive search for weak periodic signals in the frequency range 400 Hz to 2500 Hz using the burst records for 2203 GRBs and 152 SGR flares. The study employs the Rayleigh power as a test statistic to evaluate the evidence for periodic emissions. We find no evidence of periodic emissions from these events at these frequencies. In all but a very few cases the maximum power values obtained are consistent with what would be expected by chance from a non-periodic signal. In those few instances where there is marginal evidence for periodicity there are problems with the data that cast doubt on the reality of the signal. For classical GRBs, the largest Rayleigh power occurs in bursts whose TTE data appear to be corrupted. For SGRs, our largest Rayleigh power, with a significance of 1%, occurs in one record for SGR 1900+14 (at 2497 Hz), and in no other outbursts associated with this source; we thus consider it unlikely to represent detection of a real periodicity. From simulations, we deduce that the Rayleigh test would have detected significant oscillations with relative amplitude ~10% about half the time. Thus, we conclude that high frequency oscillations, if present, must have small relative amplitudes., Comment: 22 pages, 7 figures, submitted to ApJ

We present a Bayesian analysis of the energies and arrival times of the neutrinos from supernova SN 1987A detected by the Kamiokande II, IMB, and Baksan detectors, and find strong evidence for two components in the neutrino signal: a long time scale component from thermal Kelvin-Helmholtz cooling of the nascent neutron star, and a brief (~< 1 s), softer component similar to that expected from emission by accreting material in the delayed supernova scenario. In the context of this model, we show that the data constrain the electron antineutrino rest mass to be less than 5.7 eV with 95% probability. Our analysis takes advantage of significant advances that have occured in the years since the detections in both our understanding of the supernova mechanism and our ability to analyze sparse data. As a result there are substantial differences between our inferences and those found in earlier studies. We find that two-component models for the neutrino signal make the data >100 times more probable than single-component models. In addition, the radius and binding energy of the nascent neutron star implied by single-component models deviates significantly from the values predicted by current neutron star models, whereas those implied by models with an accretion component are in complete agreement with the predictions. As a result, two-component models are hundreds to thousands of times more probable than single-component models. The neutrino data thus provide the first direct observational evidence in favor of the delayed supernova scenario over the prompt scenario. (Abridged abstract), Comment: 46 pages, 12 figures, RevTeX; for submission to Physical Review D

We explore the possible role of evolution in the analysis of data on SNe Ia at cosmological distances. First, using a variety of simple sleuthing techniques, we find evidence that the properties of the high and low redshift SNe Ia observed so far differ from one another. Next, we examine the effects of including simple phenomenological models for evolution in the analysis. The result is that cosmological models and evolution are highly degenerate with one another, so that the incorporation of even very simple models for evolution makes it virtually impossible to pin down the values of $\Omega_M$ and $\Omega_\Lambda$, the density parameters for nonrelativistic matter and for the cosmological constant, respectively. Moreover, we show that if SNe Ia evolve with time, but evolution is neglected in analyzing data, then, given enough SNe Ia, the analysis hones in on values of $\Omega_M$ and $\Omega_\Lambda$ which are incorrect. Using Bayesian methods, we show that the probability that the cosmological constant is nonzero (rather than zero) is unchanged by the SNe Ia data when one accounts for the possibility of evolution, provided that we do not discriminate among open, closed and flat cosmologies a priori. The case for nonzero cosmological constant is stronger if the Universe is presumed to be flat, but still depends sensitively on the degree to which the peak luminosities of SNe Ia evolve as a function of redshift. The estimated value of $H_0$, however, is only negligibly affected by accounting for possible evolution., Comment: 45 pages, 15 figures; accepted for publication in The Astrophysical Journal. Minor revisions and clarifications made including addition of recent references

We have conducted pencil-beam searches for outer solar system objects to a limiting magnitude of R ~ 26. Five new trans-neptunian objects were detected in these searches. Our combined data set provides an estimate of ~90 trans-neptunian objects per square degree brighter than ~ 25.9. This estimate is a factor of 3 above the expected number of objects based on an extrapolation of previous surveys with brighter limits, and appears consistent with the hypothesis of a single power-law luminosity function for the entire trans-neptunian region. Maximum likelihood fits to all self-consistent published surveys with published efficiency functions predicts a cumulative sky density Sigma(

A recent study of the rotation of the plane of polarization of light from 160 cosmological sources claims to find significant evidence for cosmological anisotropy. We point out methodological weaknesses of that study, and reanalyze the same data using Bayesian methods that overcome these problems. We find that the data always favor isotropic models for the distribution of observed polarizations over counterparts that have a cosmological anisotropy of the type advocated in the earlier study. Although anisotropic models are not completely ruled out, the data put strong lower limits on the length scale $\lambda$ (in units of the Hubble length) associated with the anisotropy; the lower limits of 95% credible regions for $\lambda$ lie between 0.43 and 0.62 in all anisotropic models we studied, values several times larger than the best-fit value of $\lambda \approx 0.1$ found in the earlier study. The length scale is not constrained from above. The vast majority of sources in the data are at distances closer than 0.4 Hubble lengths (corresponding to a redshift of $\approx$0.8); the results are thus consistent with there being no significant anisotropy on the length scale probed by these data., Comment: 8 pages, 3 figures; submitted to Phys. Rev. D

We use Bayesian methods to study anisotropic models for the distribution of gamma ray burst intensities and directions reported in the Third BATSE Catalog (3B catalog) of gamma ray bursts. We analyze data obtained using both the 64 ms and 1024 ms measuring timescales. We study both purely local models in which burst sources (``bursters'') are presumed to be distributed in extended halos about the Galaxy and M31, and mixed models consisting of a cosmological population of standard candle bursters and a local population distributed throughout a standard Bahcall-Soneira dark matter halo with a 2 kpc core. We find that the purely local models we have studied can account for the 3B data as successfully as cosmological models, provided one considers halos with core sizes significantly larger than those used to model the distribution of dark matter. We infer core sizes for the halo distribution that are smaller than one might expect based on popular semiquantitative arguments, and show why such arguments can lead to unwarranted conclusions. We also find that the 3B data do not constrain the width of power-law luminosity functions for burst sources. Our analysis of mixed models finds two families of models that can successfully account for the data: models with up to 20% of observed bursts in a bright local population visible to ~ 50 kpc; and models with up to 50% of observed bursts in a dim local population visible only nearby (to less than a disk scale height). These models fit as well or better than purely cosmological models. They indicate that a surprisingly large local, anisotropic component could be present whose size is comparable to the sizes of hypothetical classes of bursts inferred from analyses of temporal and spectral characteristics., Comment: 35 pages, 18 figures; submitted to The Astrophysical Journal, 12/96

We use Bayesian methods to analyze the distribution of gamma ray burst intensities reported in the Third BATSE Catalog (3B catalog) of gamma ray bursts, presuming the distribution of burst sources (``bursters'') is isotropic. We study both phenomenological and cosmological source distribution models, using Bayes's theorem both to infer unknown parameters in the models, and to compare rival models. We analyze the distribution of the time-averaged peak photon number flux, F, measured on both 64 ms and 1024 ms time scales, performing the analysis of data based on each time scale independently. Several of our findings differ from those of previous analyses that modeled burst detection less completely. In particular, we find that the width of the intrinsic luminosity function for bursters is unconstrained, and the luminosity function of the actually observed bursts can be extremely broad, in contrast to the findings of all previous studies. Useful constraints probably require observation of bursts significantly fainter than those visible to BATSE. We also find that the 3B peak flux data do not usefully constrain the redshifts of burst sources; useful constraints require the analysis of data beyond that in the 3B catalog (such as burst time histories), or data from brighter bursts than have been seen by BATSE (such as those observed by the Pioneer Venus Orbiter). In addition, we find that an accurate understanding of the peak flux distributions reported in the 3B almost certainly requires consideration of data on the temporal and spectral properties of bursts beyond that reported in the 3B catalog, and more sophisticated modeling than has so far been attempted., Comment: 66 pages, 33 figures; submitted to The Astrophysical Journal, 12/96