Abstract

We present an IDL graphical user-interface-driven software package designed for the analysis of exoplanet transit light curves. The Transit Analysis Package (TAP) software uses Markov Chain Monte Carlo (MCMC) techniques to fit light curves using the analytic model of Mandal and Agol (2002). The package incorporates a wavelet-based likelihood function developed by Carter and Winn (2009), which allows the MCMC to assess parameter uncertainties more robustly than classic methods by parameterizing uncorrelated “white” and correlated “red” noise. The software is able to simultaneously analyze multiple transits observed in different conditions (instrument, filter, weather, etc.). The graphical interface allows for the simple execution and interpretation of Bayesian MCMC analysis tailored to a user’s specific data set and has been thoroughly tested on ground-based and Kepler photometry. This paper describes the software release and provides applications to new and existing data. Reanalysis of ground-based observations of TrES-1b, WASP-4b, and WASP-10b (Winn et al., 2007, 2009; Johnson et al., 2009; resp.) and space-based Kepler 4b–8b (Kipping and Bakos 2010) show good agreement between TAP and those publications. We also present new multi-filter light curves of WASP-10b and we find excellent agreement with previously published values for a smaller radius.

1. Introduction

The thriving field of exoplanet science began when Doppler techniques reached the precision required to detect the radial velocity (RV) variations of stars due to their orbital interactions with planets [1, 2]. The first detection of an exoplanet orbiting a main sequence star was followed closely by dozens more (51 Pegasi b; [3], see also, [4, 5]). To date, this technique has supplied the vast majority of knowledge—both detection and characterization—of exoplanets and their environments. Five years later the first photometric observation of an exoplanet transiting the stellar disk of its parent star provided a new, rich source of information on exoplanet systems (HD 209458; [68]). Transit measurements provide the true masses and radii of the planet and present many opportunities for diverse followup science [9, 10]. Unlike Doppler RV exoplanet signatures, transits affect observations only during the actual event so survey missions capable of constantly monitoring stars are needed to detect exoplanets by their transits. Today, such dedicated wide field transit survey missions have begun to keep pace with Doppler RV surveys and are ushering in a new era of exoplanetary science.

The bulk of this forward momentum has been provided by NASA’s Kepler mission, which provides photometry of unprecedented precision (down to a few parts in 105) of over 150,000 stars in a field. The first Kepler public data release referenced 1235 planet candidates, and continuing data releases promise planets and candidates in even greater numbers [11]. Kepler’s continuous high precision monitoring allows for the detection of planets over a large range of orbital periods while avoiding many biases that plague ground-based surveys [12, 13].

A planet passing between its host star and an observer blocks a region of the stellar disk, imprinting a signature on the observed stellar flux. As the planet first encounters, crosses, and finally passes beyond the stellar disk, the resulting shadow encodes the geometry of the encounter into a dip in the observed stellar light curve. With sufficient observational cadence and photometric precision, researchers may gain access to fundamental planet parameters by forward modeling the light curves. The resulting transit light curve model provides salient details regarding the geometry of an exoplanet system with no need to spatially resolve the event.

The unique shape of an exoplanet transit is well studied, with analytic descriptions of the signal (and covariances among parameters) documented throughout the literature [1417]. Even so, different fitting methods can lead to variations in derived quantities, and quantification of confidence in derived parameters depends critically on the treatment of noise, stellar limb darkening, and observational cadence. Many authors in the field of exoplanet transit photometry have access to private MCMC software codes that yield reliable results easily comparable to other studies. Yet public analysis packages providing for uniform application of transit models to the analysis of observed light curves are rare (one alternate package is JKTEBOP, see [18]) and the Kepler public data releases bring with them an ever-growing need for such software.

To satisfy this need we present the Transit Analysis Package (TAP): a graphical user interface developed for the Interactive Data Language (IDL) to allow easy access to state-of-the-art analysis of photometric exoplanet transit light curves (Figure 1). TAP models the shape of a transit with the analytic light curve of [15]. Best-fit “median” parameter values and statistical distributions are extracted using Markov Chain Monte Carlo analysis techniques (see, e.g., [1921]). Instead of the classic , TAP uses the wavelet-based likelihood function of [17], which parameterizes both uncorrelated and correlated photometric noise allowing for accurate uncertainties. For long integration observations (as is the case with most Kepler light curves) the software can be set to calculate transit models on a finer temporal scale. TAP then integrates that model to the cadence of the observed data to match distortions due to long integrations [22]. These techniques represent over a decade of work by the exoplanet (and greater astronomy) community. The software presented in this paper unifies this knowledge and experience into an easily obtainable tool. TAP has been developed for the simultaneous analysis of multiple light curves of varying cadence, filter (limb darkening), and noise (observational) conditions with no need to phase fold and rebin data. In Section 2 we briefly review the analysis tools coded into TAP. In Section 3 we present and discuss TAP analysis results from ground- (TrES-1b, WASP-4b, and WASP-10b) and space-based (Kepler-4b through 8b) observations, including a new observation of WASP-10b. We summarize the public release of TAP in Section 4.

2. Analysis Techniques

2.1. Transit Analysis Package

The TAP software package uses EXOFAST [23], an improved implementation of the [15] analytic transit model because, among transit parameterizations, it is commonly used, physically realistic, and computationally efficient. This analytic function takes as input the following parameters: planet orbit period , radius ratio , scaled semimajor axis , orbital inclination with respect to observer , orbital eccentricity , and argument of periastron , the time of transit center , and two parameters (, ) specifying a quadratic limb-darkening law describing the occulted stellar disk [15, 24, 25]. A user can choose to draw directly from and or handle the correlation between those parameters by selecting from distributions of and [19]. Limb-darkening parameters which are allowed to evolve as free parameters are restricted only so that the limb darkening is positive and monotonically decreasing towards the disk center ( and ; [16]). TAP includes three parameters for the fitting of a quadratic trend in the light curve.

Within the Bayesian statistical framework we assume the existence of a probability distribution of model parameters, , with respect to an observed data set, . That probability density function, , is accessed by considering both the model parameters and the data to be random variables [26]. TAP accesses the Bayesian posterior probability density function using a Metropolis-Hastings [27, 28] Markov Chain Monte Carlo (MCMC) analysis within a Gibbs Sampler [20, 29, 30] with which we employ a Daubechies fourth-order wavelet decomposition likelihood function [17]. In this method a chain of states (sets of model parameters, ) is calculated such that each subsequent state () is drawn randomly from a probability distribution applied to . The new state is either stored as the next link, , or rejected based on a transition probability designed to accept any state with an increased likelihood. States with lower likelihoods are accepted with a probability decreasing exponentially with the drop in likelihood. A sufficiently converged chain of model states directly represents the probability density function . With careful application, MCMC techniques allow model parameters to reach every possible state from any initial state, thus guaranteeing that the parameters converge to stationary and accurate distributions [20]. These distributions provide parameter uncertainties, which accurately reflect the input data.

The classic likelihood function is a formulation of and may be applied adequately to data that can be assumed to have no correlated sources of noise ([20, 31, 32], and many others). However, because red noise is often significant in observational data, TAP uses a wavelet-based likelihood function (equations 33-34 [17]): where and are the residuals transformed into wavelet space using a Daubechies 4th-order discrete wavelet transform. Thus, the likelihood becomes a product of Gaussian functions across the scale and translation of the wavelet basis with where we hold and . Notice that both exponential terms of (1) are reminiscent of a . Wavelet decomposition techniques provide increased confidence in derived MCMC uncertainties over the traditional likelihood by allowing parameters that measure photometric scatter (uncorrelated Gaussian , and correlated red ) to evolve freely. Reference [17] shows that, without proper red noise analysis, contaminated data can yield incorrect values and artificially small uncertainties for all parameters. For studies of transit timing variations (TTVs), proper handling of red noise is critical.

TAP infers parameters and uncertainties from final MCMC distributions as follows: all links before every parameter first crosses its median value are removed from each chain to eliminate the MCMC “burn in’’ period during which parameters settle into a good initial fit regardless of the input model parameters. All chains are added together and the 15.9, 50.0, and 84.1 percentile levels are recorded. The 50.0% (median) level is reported as the inferred best value, while the 15.9% and 84.1% levels are reported as the “1” confidence levels. An example is provided in Figure 3. In addition to combining existing methodologies, TAP allows for the simultaneous analysis of multiple transits as distinct entities. This is accomplished by defining an overall likelihood with components from each observed transit and creating a matrix of parameter sets such that certain parameters can be forced to evolve together. In this way, system parameters benefit statistically from additional data while those parameters that depend on individual light curves (i.e., filter choice, weather conditions, statistical scatter, etc.) evolve independently during the MCMC. An application is discussed in Section 3.1.

TAP provides tests to assure that chains are sufficiently mixed and that the final results are consistent with convergence. When a user executes a TAP MCMC analysis, they specify number of chains, minimum links per chain, and minimum effective links per chain. After calculating the minimum links, TAP calculates the correlation length for each parameter : where the correlation length is determined when reaches 0.5 [21]. The effective length for the chain is then the ratio of the total number of links to the correlation length. Sufficiently mixed MCMC chains have . TAP extends the length of the chain until all free parameters reach an effective length of greater than the minimum input value. After all chains have finished, TAP calculates the Gelman-Rubin statistic [30, 33] and reports the values for each parameter and transit in an output file. values below 1.1 are consistent with converged chains.

3. Applications to Observed Photometry

The TAP software has passed thorough tests on synthetic and observed photometry. In this section we document tests on ground- (WASP-10b, TrES-1b, WASP-4b, Section 3.1) and space based (Kepler-4b through-8b, Section 3.2) light curves.

3.1. Ground-Based Photometry: WASP-10b, TrES-1b, and WASP-4b

We reanalyze the WASP-10b transit light curve from [31] following the same procedure and recover the system parameters derived in that work. Specifically, we lock the eccentricity and orbital period within the MCMC analysis to values derived from previous work and allow the remaining parameters to evolve freely. We present two analyses, first with the red noise parameter locked at zero to mimic the classic analysis of [31] and second with the red noise as a free parameter (Table 1). We find negligible red noise contamination as concluded by [31].

We observed an additional transit of WASP-10b in three bands (Barr V, R, and I) on UT August 22 2010 in clear but nonphotometric weather using the Orthogonal Parallel Transfer Imaging Camera (OPTIC) mounted on the University of Hawaii 2.2 m telescope on Mauna Kea. OPTIC’s unique Orthogonal Transfer Array CCD (OTA, see [3436]) provides submillimag precision through a technique called PSF shaping where the light from the bright star being occulted is spread out over a small box on the CCD during exposure. Shaped-PSFs allow for longer exposure times on bright stars and improve the duty cycle by drawing out the exposure before incurring the 30-second readout. Between each exposure the OPTIC filter wheel was advanced so that the transit light curve is fully covered by observations in V, R, and I bands. The light curves for each band were produced using basic aperture photometry techniques with a square aperture to account for the shaped PSFs. A relatively large aperture was chosen relative to the PSFs to account for variable PSF scatter due to changing atmospheric conditions. Four background apertures located on each side of the PSF were medianed and subtracted. Unlike [31] light curve red noise contamination is evident (Figure 2).

We conduct two analyses—one in which all transits are considered independent observations and a second in which we use the full power of TAP by locking system and orbit parameters into a single set (Section 2). For this combined analysis, parameters describing limb darkening, noise, and a quadratic trend to account for out-of-transit variations remain independent per light curve. The results of these analyses are reported in Table 4. Each individual transit and the combined set yield parameters consistent with [31]. The TAP simultaneous analysis method shows improved precision in the derived system parameters while providing limb-darkening coefficients in three photometric bands from a single observation. A final analysis was conducted by locking the parameter for red noise at zero to mimic the classic , which derived parameters consistent but with error bars scaled down by a factor of 0.7—an overestimation of the significance of the contaminated data.

In [31], a high precision (~0.5 millimag) transit of WASP-10b allowed for a robust MCMC fit in which only period and eccentricity parameters required constraint. The results of that work disagreed significantly with the discovery parameters of [37], in particular for transit depth, yielding a 16% smaller planetary radius (a change). Since then, two additional papers [38, 39] have supported the conclusions of the discovery paper—that WASP-10b is inflated beyond the theoretical predictions of radii of irradiated Jovian planets [40]. Our data agree with the smaller radius measurement of [31], using the same instrument but with three different filters. In our analysis, like that of [31], the limb-darkening parameters are allowed to vary freely, a robust method allowed by our high precision curves. Reference [37] locks these parameters at theoretical values, while the treatment of limb darkening in [38, 39] is unclear. We stress the importance of additional followup observations to clarify these discrepancies and urge authors to provide the light curves that they collect. This allows for reanalysis by other teams and is a boon to the community as a whole.

We reanalyze the three TrES-1b transit light curves of [41] and two WASP-4b transits of [42]. In both cases, all system parameters except are locked to a single set. The transit midtime , along with observation-specific parameters including noise and a linear airmass trend, was allowed to evolve independently per transit observation. For both systems TAP reproduces the previously published system parameters. We report the results as a comparison to the published values of [41] in Table 2 and [42] Table 3.

3.2. Space-Based Photometry: Kepler-4b through-8b

Ideally, multiple transits of an exoplanet are analyzed as a set. This establishes a measurement of orbital period and increases confidence in derived parameters. In general, the transit midtime of each observation is measured and the data is phase folded into a single light curve for analysis. The weaknesses of this technique are numerous, including a loss of per-transit midtime error, troublesome handling of different observing conditions, and the loss of ability to characterize the noise of each individual observation. In the era of Kepler and beyond, multiple transits are the norm and often ground based followup provides additional data of wildly differing quality.

TAP is designed for exactly such analyses, and reproducing published system parameters for Kepler targets is an ideal test of these new methods. Data for Kepler-4b through-8b were acquired through the MAST Kepler data archive. We normalized the Kepler flux around each transit event to an out-of-transit value of unity using a 5th-order polynomial fit to account for stellar flux variations over each transit epoch. During analysis TAP was set to calculate models with a factor of ten finer temporal cadence before rebinning that model back to the observed 29.4-minute cadence. Eccentricity parameters were locked to assume circular orbits.

To compare with the results of [43] we use an argument of fractional difference, where are the parameter values and their 1- uncertainties. Thus, a fractional difference with absolute value less than unity represents statistical agreement. The inferred system parameters (Table 5) show such agreement with the complementary method “A.c” of [43] for Kepler 4b–7b, where the authors assume circular orbits and allow limb-darkening parameters to evolve freely. In contrast to our analysis, [43] also takes into account RV measurements. For the case of Kepler-8b, the difference between TAP and KB10 fits of orbital period is significant. This difference is likely attributed to the increased noise in the Kepler 8 dataset and the inclusion of radial velocity measurements by [43].

4. Summary

In this paper we present a graphical user interface-driven exoplanet transit analysis package written in the Interactive Data Language (IDL). This software provides easy access to state-of-the-art analysis techniques backed by rigorous testing. TAP is designed for the future of exoplanet transit science in which the global analysis of multiple transits is the normal mode of operation. Accurate parameter uncertainties are derived by combining new wavelet-based MCMC methods for handling correlated noise and the ability to maintain multiple transits as individual datasets.

By providing public access to TAP, we allow for uniform analyses by otherwise disconnected teams across the globe. Any method and implementation is prone to systematic errors, and the increasing sophistication of transit analysis techniques provides an increasing possibility for coding errors and severe systematics. By using the same tools, the effects of systematics are marginalized and results from separate teams can be simply and meaningfully compared. Utilizing MCMC software, which accounts for correlated noise, greatly increases the accuracy of derived parameter uncertainties, a critical component in transit science where inaccurate error bars can lead to false claims of planet phase amplitudes, secondary eclipse depths, and transit timing and shape variations.

TAP is one culmination of groundbreaking work by numerous teams over the last ten years of exoplanet transit science. Our goal—to provide these state-of-the-art techniques freely to all—arises from the needs of many researchers. These include professional astronomers, graduate and undergraduate students, and amateur astronomers. With no need to conduct intensive software design campaigns, these individuals and teams may focus on interpreting the rich information available in their data. As Kepler and ground based campaigns continue to collect light curves of hundreds and soon thousands of exoplanet transits, the need for such public software has never been greater.

TAP is publicly available at http://ifa.hawaii.edu/users/zgazak/IfA/TAP.html along with an example video and FAQ file.

Acknowledgments

The authors gratefully acknowledge discussions with Josh Carter relating to the new wavelet-based likelihood methods and David Kipping and Gáspár Bakos for responding thoughtfully to questions regarding their independent analysis of Kepler-4b through -8b. They thank the individuals who provided helpful insight and comments regarding early versions of TAP, including especially John Rayner, Karen Collins, and Kaspar von Braun. J. A. Johnson was supported by the NSF Grant AST-0702821. A. W. Mann was supported by NSF grant AST-0908419. For this work the authors made use of the UH 2.2 m telescope atop Mauna Kea. The authors wish to extend special thanks to those of Hawaiian ancestry on whose sacred mountain of Mauna Kea they are privileged to be guests. Without their generous hospitality, the observations presented herein would not have been possible.