No Evidence for Lunar Transit in New Analysis of Hubble Space Telescope Observations of the Kepler-1625 System
DDraft version May 23, 2019
Typeset using L A TEX twocolumn style in AASTeX62
No Evidence for Lunar Transit in New Analysis of
Hubble Space Telescope
Observations of the Kepler-1625 System
Laura Kreidberg,
1, 2
Rodrigo Luger, and Megan Bedell Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138 Harvard Society of Fellows, 78 Mount Auburn Street, Cambridge, MA 02138 Flatiron Institute, Simons Foundation, 162 Fifth Ave, New York, NY 10010, USA
Submitted to ApJLABSTRACTObservations of the Kepler-1625 system with the Kepler and Hubble Space Telescopes have suggestedthe presence of a candidate exomoon, Kepler-1625b I, a Neptune-radius satellite orbiting a long-periodJovian planet. Here we present a new analysis of the Hubble observations, using an independent datareduction pipeline. We find that the transit light curve is well fit with a planet-only model, with abest-fit χ ν equal to 1 .
01. The addition of a moon does not significantly improve the fit quality. Wecompare our results directly with the original light curve from Teachey & Kipping (2018), and findthat we obtain a better fit to the data using a model with fewer free parameters (no moon). We discusspossible sources for the discrepancy in our results, and conclude that the lunar transit signal found byTeachey & Kipping (2018) was likely an artifact of the data reduction. This finding highlights the needto develop independent pipelines to confirm results that push the limits of measurement precision.
Keywords: planets and satellites: individual (Kepler-1625b I) INTRODUCTIONMoons are abundant in the Solar System, and provideclues to the formation history, evolution, and even hab-itability of the planets they orbit. The great scientificpotential of moons has prompted extensive search forlunar companions in exoplanetary systems (exomoons),and creative development of new search techniques (e.g.Kipping 2009a,b; Kipping et al. 2013; Simon et al. 2010;Peters & Turner 2013; Heller et al. 2014; Noyola et al.2014; Hippke 2015; Agol et al. 2015; Sengupta & Marley2016; Vanderburg et al. 2018).Recently, a potential exomoon candidate was identi-fied in the Kepler-1625 system (Teachey et al. 2018).The host planet, Kepler-1625b, has a radius consistentwith that of Jupiter and an orbital period of 287 days.The first evidence for the exomoon candidate, Kepler-1625b I, was based on observations from
Kepler . Thelight curve showed drops in stellar flux that were at-tributed to a transiting exomoon; however, later analysiscalled this result into question, showing that the moontransit features were highly sensitive to the
Kepler re-duction pipeline and the algorithm used to detrend thedata (Teachey & Kipping 2018; Rodenbeck et al. 2018).Subsequent follow-up observations with
HST revivedthe possibility of an exomoon in the system based on two factors: a small drop in the system flux after the planet’stransit egress and a transit timing variation (Teachey &Kipping 2018, hereafter TK18). The best fit moon hada large radius (comparable to that of Neptune), and ifreal, is unlike any moon in the Solar System.Since the primary evidence for the exomoon now restson the
HST transit light curve, in this work we performan independent reduction and fit to the
HST data andcompare it to the results from TK18. OBSERVATIONS AND DATA REDUCTIONThe Kepler-1625 system was observed with 26 con-tinuous
HST orbits on 28 - 29 October, 2017 (ProgramGO 15149: PI: A. Teachey). The observations used theWide Field Camera 3 (WFC3) G141 grism in staringmode, which fixed the spectrum in a constant positionon the detector. At the beginning of the visit, there wasa single exposure taken with the F130N filter, which isused to determine the position of the spectral trace. Thefollowing exposures used the G141 grism. See TK18 foradditional description of the observation design.We reduced the
HST data using custom softwaredeveloped in Kreidberg et al. (2014). This softwarehas yielded consistent results with multiple independentpipelines (e.g. Knutson et al. 2014; Spake et al. 2018).We ran our pipeline on the flt data product provided by a r X i v : . [ a s t r o - ph . E P ] M a y Kreidberg et al. the Space Telescope Science Institute (STScI). In keep-ing with previous WFC3 analysis, we discarded the firstorbit of data, where the instrument systematics havelarger amplitude. We also discarded exposures takenduring the South Atlantic Anomaly passage (exposures107, 116, 125, and 126).To begin the data reduction, we fit the centroidof the direct image with a two-dimensional Gaus-sian. The centroid position determines the posi-tion of the spectral trace, which we calculated us-ing the coeffients provided in the configuration filefrom STScI:
G141.F130N.V4.32.conf . To process thespectra, we flatfielded the raw data using the spec-troscopic flatfield coefficients provided by STScI in WFC3.IR.G141.flat.2.fits , following the instructionsin Section 6 of the aXe User Manual . We then cre-ated an extraction box centered on the spectral trace.We varied the height and width of the box in 1-pixelincrements to find the window that minimized the root-mean-square (rms) deviation from the best fit to thetransit light curve. The best was 450 < X < < Y < flt files; WFC3 Data Handbook ), and the error due tobackground subtraction. The initial pixel mask markedall pixels as good.To optimally extract the spectrum, we first createda smoothed image by median-filtering each row of thedata with a 9-pixel-wide window. We then normalizedthe smoothed image by dividing each column by its sum,and multiplied it by the best guess spectrum. We com- http://ane-info.stsci.edu/ pared the smoothed image to the real data and maskedoutliers in the data that are greater than a threshold σ cut = 7 .
5. We then recomputed the best guess spec-trum with the new mask and the optimal weights fromHorne (1986). The process is iterated until no out-liers greater than the threshold remain. To create thebroadband transit light curve, we sum each optimallyextracted spectrum over all wavelengths.The broadband light curve is shown in Figure 1, incomparison to the light curve from TK18. We note thatthere are differences between the two data sets, particu-larly a kink near the moon-like transit feature identifiedby TK18. 2.1.
Background Subtraction
The star Kepler-1625 is faint (H mag = 14.0) relativeto most other exoplanet host stars observed with WFC3,which makes accurate background subtraction especiallyimportant for this target. Moreover, the host star is ina crowded field, so the pixels used to estimate the back-ground must be chosen carefully to avoid contaminationfrom other stars. To estimate the background counts,we masked pixels with total counts larger than 800 elec-trons (2.7 electrons/sec) and took the median count inthe unmasked pixels. The per pixel uncertainty due tobackground subtraction is 1.4826 times the median ab-solute deviation.2.2.
Pointing Drift Measurement
The position of the spectrum on the detector shiftsslightly over time ( ∼ . flt image over all columns (which we dub the “col-umn sum”). We used the first exposure in the visit asa template, and for each subsequent exposure, we usedleast-squares minimization to calculate the shift in pixelsthat minimized the difference between its column sumand the template. The shifts are a fraction of a pixel, sowe used the NumPy interp routine to do linear inter-polation on a sub-pixel scale. The WFC3 point spreadfunction is undersampled, so we convolved each columnsum with a 4-pixel-wide Gaussian before the interpola-tion (following Deming et al. 2013).To measure the spectral shifts, we repeated this pro-cedure with two differences: (1) we used the optimallyextracted spectrum rather than the column sum; and hat’s No Moon P h o t o e l e c t r o n s ( × ) TK18This work0.0 0.5 1.0 1.5Time since first exposure (days)0.250.240.230.22 D i ff e r e n c e ( e - × ) this work - TK18 Figure 1.
Top: extractred light curve from this work (blue)compared to TK18 (red). Bottom: the difference in pho-toelectron counts between the data sets. The moon ingressidentified by TK18 is marked by the dotted gray line. (2) in addition to calculating the best fit shift, we alsocalculated a best fit normalization factor (a scalar mul-tiple for the whole spectrum), to ensure that our resultsare not biased by the varying brightness of the host starduring the planet’s transit. ANALYSISThe extracted transit light transit light curve containsboth astrophysical signals and instrument systematicnoise, which we model simultaneously.3.1.
Astrophysics Model
For the astrophysics, we used the planetplanet pack-age (Luger et al. 2017), a photodynamical code thatcalculates light curves for multiple occulting bodies or-biting a star. Within planetplanet , the orbits are com-puted with the N-body integrator
REBOUND (Rein & Liu2012), which calculates the three-dimensional motion ofthe star, planet, and moon over time under the influenceof gravity.In our analysis, we considered two scenarios: a no-moon model and a moon model. The free parametersfor the no-moon model were: the stellar radius, theplanet radius, the planet’s time of central transit, theplanet inclination, and the planet mass. For the moonmodel, we added a third body with six free parame-ters: moon radius, moon time of central transit, moon Sp e c t r a l s h i f t ( p i x e l s ) Sp a t i a l s h i f t ( p i x e l s ) Figure 2.
Shift (in pixels) relative to the mean positionof the spectrum in the spectral direction (top) and spatialdirection (bottom). The largest shift occurs after orbit 14due to a guide star reacquisition. orbital period, moon inclination, moon mass, and longi-tude of the ascending node relative to that of the planet.We fixed the eccentricity of all bodies to zero. We alsofixed the orbital period of the planet-moon barycenter to287.378949 days (the best fit from TK18). We electednot to vary the period of the planet-moon barycenterbecause it is poorly constrained from a single transitobservation. A longer period would cause a longer tran-sit duration for the planet, but this could equivalentlyresult from a smaller impact parameter, a larger star, ora massive moon that significantly perturbs the planetaryorbit.We used the following priors: the stellar radius wasnormally distributed, R ∗ ∼ N (1 . , . R (cid:12) (see nextsection). The planet radius was uniform from 8 − R ⊕ .The planet transit time was uniform over the times-pan of the observations, and inclination was uniformfrom 0 − ◦ . The planet mass was log-normally dis-tributed, M p ∼ . ± . M ⊕ , based on the expectationfor Jupiter-radius objects from Ning et al. (2018). Forthe moon model, we allowed the transit time to varyuniformly over the entire visit. The moon period wasuniform between 1.6 to 260 days. These limits span theduration of the HST observations (so there is one possi-ble moon occultation event), to the orbit at 0.5 the Hillradius, based on the stability limit for prograde moon
Kreidberg et al. orbits (Domingos et al. 2006). The Hill radius calcula-tion assumed the planet and stellar masses are 1 M Jup and 1 . M (cid:12) . The moon mass varied uniformly from0 − M ⊕ . The longitude of the ascending node wasalso uniform from 0 − ◦ . The moon inclination wasuniform from 0 − ◦ , and was defined relative to the lineof sight. We assigned zero prior probability to scenarioswhere the moon did not transit or experienced a grazingtransit. We made this choice to put an upper limit onthe radius of a fully transiting moon; for grazing tran-sits or non-transits the moon radius could be arbitrarilylarge. 3.1.1. Stellar Parameters
For both the moon and no-moon scenarios, we useda quadratic stellar limb darkening law and fixed the co-efficients to the prediction for a 5700 K, solar metal-licity
PHOENIX model from Espinoza & Jord´an (2015); u , u = [0 . , . . +0 . − . M (cid:12) , radius1 . +0 . − . R (cid:12) , and age 2 . +1 . − . Gyr. In our analysis, wefixed the stellar mass to the best fit value (1 . M (cid:12) ), andused a Gaussian prior on the radius, R ∗ ∼ N (1 . , . Instrument Systematics Model
There are two systematic trends in the data. Oneis the orbit-long ramp, attributed to charge traps in thedetector filling up over the orbit (Zhou et al. 2017). Theother is a visit-long trend over multiple orbits, whichcould be due to shifts in the target star position ontomore/less sensitive pixels.For our primary fit, we used a linear combination ofthe spectrum’s X and Y position (shown in Figure 2),multiplied by the non-parametric orbital ramp modelfrom TK18, which assigns each of the nine exposuresper orbit a normalization constant, c , ..., c . In sum,for exposure number i , the systematics model S is:S i = c j × (1 + a X i + b Y i ) (1)where a and b , and c j are free parameters, and j = i mod 9 + 1 is the exposure number relative to the firstexposure in the orbit.For comparison, we also tested the “exponential” sys-tematics model from TK18, which combines an expo-nential visit-long trend, an offset after orbit 14 when the guide stars were reacquired, and the non-parametricorbital ramp model.3.3. Light Curve Fits
We fit the broadband transit light curve using themodels described above. We determined the best fitmodel parameters with least-squares minimization. Wealso ran a Markov chain Monte Carlo (MCMC) fit todetermine the posterior probability of the parameters.For the MCMC, we held the ramp parameters c , ..., c fixed at their best-fit values. The MCMC used the emcee package (Foreman-Mackey et al. 2013) with 50 walkersand ran for 5000 steps. We discarded the first 20% ofthe MCMC chain as burn-in. As a quick test for con-vergence, we divided the remainder of the chain in halfand confirmed that the results from the first half wereconsistent with the second half. RESULTSWe obtained an excellent fit to the light curve with theno-moon model, as illustrated in Figure 3. The residualsto the no-moon model fit have rms equal to 356 partsper million (ppm), which is within 3% of the predictedphoton shot noise (367 ppm), and yields a χ ν = 1 . χ ν to 1.02. According to the Bayesian in-formation criterion (BIC), which penalizes unnecessarymodel complexity, the moon model is disfavored with∆BIC = 26 .
7. This constitutes strong evidence againstthe inclusion of a moon (Kass & Raftery 1995). In ad-dition, as shown in Figure 5, the posterior distributionfor the moon transit time spans the entire duration ofthe observations (2 σ confidence). The upper limit onthe moon radius is 3 . R ⊕ at 3 σ confidence.We found that our results were unchanged when weused the exponential systematics model from TK18 (de-scribed in § . . +0 . − . R ⊕ and mass is10 . ± . M ⊕ .4.1. Comparison with Teachey & Kipping (2018)
TK18 found evidence for the transit of a Neptune-sizemoon in their analysis of the HST data, in contrast to hat’s No Moon R e l a t i v e f l u x rms = 355.9 = 206.0dof = 203BIC = 292.2 This Work no moon0.25 0.00 0.25 0.50 0.75 1.00 1.25 1.50Time since first exposure (BJD
TDB )0.9960.9981.000 R e l a t i v e f l u x rms = 351.0 = 200.4dof = 197BIC = 318.9 moon rms = 386.2 = 243.6dof = 202BIC = 335.2 Teachey & Kipping (2018) no moon0.25 0.00 0.25 0.50 0.75 1.00 1.25 1.50Time since first exposure (BJD
TDB )rms = 361.7 = 213.7dof = 196BIC = 337.7 moon Figure 3.
Best fit models compared to normalized transit light curves from this work (left, blue) and from TK18 (right, red).The top panel shows the best fit no-moon model, and the bottom shows the best fit moon model. The lower left of each panelindicates the fit rms (in ppm), the χ , the degrees of freedom, and the BIC. Each light curve is divided by its best fit systematicsmodel (XY decorrelation for this work; exponential and offset for TK18). The dotted gray line marks the possible moon ingressidentified by TK18. N o r m a li z e d R M S Expected rmsThis work, no-moonThis work, moonTK18, no-moonTK18, moon
Figure 4.
Light curve rms versus bin size for the best fit no-moon model (solid lines) and moon model (dashed lines), fordata from this work (red) and TK18 (blue). The fits to datafrom this work agree well with the expected photon-limited, √ N decrease in rms with bin size (black line). We also reachthe photon limit for the TK18 data, but only for the moonmodel. The rms for the no-moon fit (red dashed line) rangesfrom 1 . − . × the photon limit for 1 to 20 points per bin. the findings presented here. To compare our results withtheirs, we fit the TK18 light curve directly. We fit theastrophysical signal with both the no-moon and models,and used the exponential systematics model. Figure 3shows the best fit models.Similar to the findings of TK18, the moon model im-proved the fit quality by a ∆ χ = 29 .
9. Notably, how-ever, the moon model fit to the TK18 data does notperform better than the no-moon model fit to our newdata (rms of 362 versus 356 ppm), even with the addi-tional seven free parameters.The moon model also yields qualitatively differentposterior distributions for the two data sets. As shownin Figure 5, for the TK18 data the moon radius andtransit time are peaked at r moon = 4 . ± . R ⊕ and t moon = 2458056 . +0 . − . BJD
TDB . By contrast, the fitto the new data presented here yields an upper limiton the moon radius of 3 . R ⊕ at 3 σ confidence, and thetransit time is unconstrained.Although the two data sets yield different constraintson the moon properties, the planet’s mid-transit timeagrees to better than 1 σ for the two fits. The tran-sit time is earlier than expected based on the Keplerdata (3 σ confidence; TK18), suggesting that there aretransit timing variations in the system. Such variationcould arise from the presence of a moon, as suggested Kreidberg et al. . . . . Moon mid-transit time(BJD
TDB - 2458055) M oo n r a d i u s ( R ) Moon radius (R )
Figure 5.
Posterior distributions for the moon radius andtime of transit based on an MCMC fit to data from this work(blue) and from TK18 (red). The shading corresponds to 1 − ,2 − , and 3 σ contours (from darkest to lightest). These valuesare marginalized over all other model parameters. by TK18; however, the variation could also be causedby another planet in the system. DISCUSSION & CONCLUSIONSA natural question arising from our analysis is thesource of the discrepancy between TK18 and the newresults presented here. We find that with our new data,there is strong evidence against the moon (∆BIC > . − .
74% lower than theTK18 light curve, and there is a small bump in the dif-ference between the two data sets near the location ofthe moon transit identified in the TK18 data (see thebottom panel). This bump may be the source of themoon feature reported in TK18.We explored several modifications to our pipeline toattempt to reproduce the TK18 data reduction. Theseincluded rotating the image by 0.5 degrees, using thesame aperture as TK18 to extract the spectrum, and scaling the master sky flat for the background subtrac-tion (rather than just subtracting the median). Noneof these modifications had a significant effect on our re-sults.There are a few other steps in the TK18 data reduc-tion that would require substantial modification of ourpipeline to recreate, but seem unlikely to be responsi-ble for the difference. One of these is outlier masking.TK18 identify outliers with a Gaussian process fit to thepixel-level light curves, compared to our optimal extrac-tion approach. Despite the difference, both methods flag ∼ aXeprep to embed the raw 256 ×
256 image ina larger array; however, this process primarily affectsthe edge of the image, many pixels distant from the ex-traction box, so it is unclear how this step would biasthe light curve. We conclude that no single choice inthe data reduction provides an easy explanation for thedifference in our light curves.During the referee process for this work, we learned ofanother manuscript that also reanalyzed the
HST tran-sit observation (Heller et al. 2019). The best fit favoreda moon model similar to that found by TK18; however,an MCMC analysis did not converge on this model, lead-ing the authors to conclude that the highest likelihoodsolution may be an outlier.Taken together, these findings illustrate the challengeof pushing measurement precision to the 100 ppm level,and highlight the importance of developing multiple in-dependent pipelines to confirm cutting-edge results.We thank A. Teachey for helpful discussions and forproviding the extracted light curve from TK18. Wealso thank the anonymous referee for a thoughtful re-port that improved the manuscript. The HST datapresented in this paper were obtained from the Mikul-ski Archive for Space Telescopes (MAST). STScI isoperated by the Association of Universities for Re-search in Astronomy, Inc., under NASA contract NAS5-26555. Support for MAST for non-HST data is pro-vided by the NASA Office of Space Science via grantNNX13AC07G and by other grants and contracts. Wealso use data from the European Space Agency (ESA)mission
Gaia
Gaia
Gaia
MultilateralAgreement. hat’s No Moon7REFERENCES