Evidence for a planetary mass third body orbiting the binary star KIC 5095269

TL;DR

This study presents evidence for a planetary-mass object orbiting the binary star KIC 5095269, identified through eclipse timing variations, with a stable orbit suggested by dynamical simulations.

Contribution

First detection of a planetary-mass body orbiting a close binary star via eclipse timing variations using Kepler data.

Findings

Detected a $7.70 ext{ M}_{Jup}$ planet orbiting every 237.7 days.

Indicated a stable orbital configuration over 10 million years.

Recommended radial velocity follow-up for confirmation.

Abstract

In this paper, we report the evidence for a planetary mass body orbiting the close binary star KIC 5095269. This detection arose from a search for eclipse timing variations among the more than 2,000 eclipsing binaries observed by Kepler. Light curve and periodic eclipse time variations have been analysed using Systemic and a custom Binary Eclipse Timings code based on the Transit Analysis Package which indicates a $7.70\pm0.08M_{Jup}$ object orbiting every $237.7\pm0.1d$ around a $1.2M_\odot$ primary and $0.51M_\odot$ secondary in an 18.6d orbit. A dynamical integration over $10^7$ years suggests a stable orbital configuration. Radial velocity observations are recommended to confirm the properties of the binary star components and the planetary mass of the companion.

Tables

Table 1. Table 1: Table of eclipse times for KIC 5095269.The eclipse time is in BJD - 2,454,833. The eclipse times are used to produce an O-C diagram to look for eclipse timing variations. A sample of the table is shown here, the full table is available online.

Eclipse Time	Error
(BJD - 2,454,833)	(Days)
133.865203	$2.78\times {10}^{-5}$
152.478622	$2.76\times {10}^{-5}$
171.091481	$2.72\times {10}^{-5}$
189.702749	$2.71\times {10}^{-5}$
208.313380	$2.67\times {10}^{-5}$
226.924859	$2.71\times {10}^{-5}$
245.537254	$2.63\times {10}^{-5}$
264.150062	$3.12\times {10}^{-5}$

Table 2. Table 2: Table of data provided by the Systemic fit for the binary system. The mass of the primary and secondary stars, the orbital period of the binary stars and the inclination were fixed at pre-calculated values while the rest of the parameters were free to be fit. The Median and Median Absolute Deviation (MAD) values were also determined by Systemic.

Property	Best-Fit Value	Median Value	MAD Value
Primary Star Mass ($M_{\odot }$)	1.21	N/A	N/A
Secondary Star Mass ($M_{\odot }$)	0.51	0.51	$1.24\times {10}^{-3}$
Orbital Period (d)	18.61196	18.61196	$1.55\times {10}^{-7}$
Mean anomaly (deg)	7.44	7.44	$1.31\times {10}^{-4}$
Eccentricity	0.246	0.246	$9.3\times {10}^{-6}$
Long. Peri (deg)	22.82	22.82	$2.23\times {10}^{-4}$
Inclination (deg)	80.0	80.0	0.02
Node (deg)	305.54	305.54	$2.3\times {10}^{-4}$

Table 3. Table 3: Table of data provided by the Systemic fit for the third body. All of the parameters were free to be fit. The Median and Median Absolute Deviation (MAD) values were also determined by Systemic.

Property	Best-Fit Value	Median Value	MAD Value
Mass ($M_j$)	7.698	7.693	0.054
Orbital Period (d)	237.70817	237.68977	0.08237
Mean anomaly (deg)	290.92	289.44	2.34
Eccentricity	0.0604	0.0603	0.0021
Long. Peri (deg)	27.67	29.03	2.07
Inclination (deg)	105.92	105.83	0.98
Node (deg)	64.19	64.10	0.28

Equations

$$T_n=P\times n+T_0$$

$$\begin{array}{r}{T}_{np}=P\times n+T_{0p},\\ T_{ns}=P\times n+T_{0s}.\end{array}$$

Full text

††pubyear: 2016††pagerange: Evidence for a planetary mass third body orbiting the binary star KIC 5095269–Evidence for a planetary mass third body orbiting the binary star KIC 5095269

1 Introduction

Planet formation is widely considered to be dominated by core accretion, but an alterative disc instability mechanism has been proposed to explain gas giant planets and brown dwarfs (Boss, 2012; Chabrier et al., 2014). In addition, while it is convenient to set a mass divide (at around 13 Jupiter masses) between planets and deuterium burning brown dwarfs (Burgasser, 2008), it has been argued that a separation of planets and brown dwarfs based on the formation mechanism is more physically meaningful (Nordlund, 2011). A key way to constrain planet formation models is to test their predictions as to the frequency, masses, orbits, and stability of planets orbiting eclipsing binary stars, whose mutual eclipses can also provide accurate host star properties. The problem to be solved however is to find such planets, as to date very few such circumbinary planets have been found (Sigurdsson et al., 2003; Correia, A. C. M. et al., 2005; Lee et al., 2009; Doyle et al., 2011; Kostov et al., 2016). In particular the discovery of the long period transiting circumbinary planet Kepler-1647b may represent an example of a large population of distantly orbiting massive planets orbiting close binary stars (Kostov et al., 2016). Our research therefore represents the initial results of a search for eclipsing binary planets that uses eclipse timings to enable planets orbiting above or below the stellar orbital plane to be detected.

Contained in the ‘Kepler Eclipsing Binary Catalog’ are more than 2,000 eclipsing binaries that have been observed over the life of the KEPLER mission (Prsa et al., 2011; Slawson et al., 2011). Predominately detached binaries stars, i.e. binary stars with a morphology classification of less that 0.5, account for almost half of the systems in the eclipsing binary catalog. The high precision observations that were performed allows an eclipse time study to be performed. Eclipsing binary stars that are detached and isolated should have eclipses that occur a constant and predictable time apart. Plotting the observed eclipse time (O) minus the calculated eclipse time (C) against a best-fitting linear ephemeris, variations from this constant time may be able to be seen. Periodic variations may be the result of a third body orbiting the binary (Beuermann, K. et al., 2010).

In systems that show periodic variations, the properties of the binary stars need to be estimated in order to fit and determine the characteristics of any additional bodies. Estimates for the masses of the binary stars are calculated from the colour data given in the Kepler data and modelling the light curve in JKTEBOP (Southworth et al., 2004). Colours and masses for spectral types are given in Pecaut & Mamajek (2013). With a mass ratio and mass estimate from the system colours, individual masses can be worked out. By using Systemic (Meschiari et al., 2009; Meschiari & Laughlin, 2010) a system can be set up with the masses and characteristics of the binary stars. From here additional bodies can be added and fit to determine if characteristics can account for eclipse time variations.

In this paper, we report on the results of an eclipse time study of a specific Kepler system, KIC 5095269, and the follow up Systemic study in order to determine the characteristics of a third body. We propose the existence of a third body around KIC 5095269 with a mass of $7.70\pm 0.08$ Jupiter masses.

2 O-C Production and Identification

We used the Kepler data to produce O-C diagrams to study eclipse timing variations. Detached eclipsing binaries were selected in order to minimise variations from within the system itself. A primary eclipse occurs when the larger star passes in front of the smaller star, while a secondary eclipse occurs when the smaller star passes in front of larger star. The time of as many primary eclipses and secondary eclipses as possible must be determined in order to perform an eclipse time variation study. We created a program, called BET or Binary Eclipse Timings, to determine eclipse times. BET is based on the software Transit Analysis Package or TAP (Gazak et al., 2012) and uses the analytic formulae for the transit or eclipse of a star which are found in Mandel & Agol (2002). The analytic formulae in Mandel & Agol (2002) describe a system of two objects during various points in its orbit. The objects can be a star and a planet (i.e. describing transits) or two stars (i.e. describing eclipses). The systems are described using the parameters: orbital period, the radius ratio of the two objects, scaled semi-major axis, orbital inclination, orbital eccentricity, argument of periastron, mid-time of eclipse/transit and two parameters specifying quadratic limb darkening. BET detects eclipses from the Kepler data and uses the analytic formulae to accurately determine the mid-eclipse times of a system.

With the observed eclipse times of a system determined, calculated eclipse times are needed in order to produce an O-C diagram. Since the time between eclipses should be constant, a calculated eclipse time can be found with the equation:

[TABLE]

Where P is the period of the system, n is the cycle number and T0 is the initial eclipse time.

Equation  (1) can be modified to account for primary and secondary eclipses and take the form seen in  (2).

[TABLE]

Where T0p is the initial primary eclipse, T0s is the initial secondary eclipse, n is the cycle number and P is the period of the system which is common to both primary and secondary eclipses.

By performing a least-squares best-fit to the observed eclipse times with  (2) the best-fit period and initial eclipse times will be found. Expected eclipse times can then be calculated. By plotting the observed eclipse time minus the calculated eclipse time against the predicted eclipse time, variations from the expected may be observed. These variations have been separated into five different, custom defined, categories based on their O-C diagrams: No or Irregular variations, Periodic variations, Sudden period flips, Long term trends and out of phase long term trends. Variations may be caused by star spots (Orosz et al., 2012), apsidal motion (Beuermann, K. et al., 2010) or dynamical interactions (Borkovits, T. et al., 2003). It is also possible that periodic variations are caused by the effects of a third body (Beuermann, K. et al., 2010).

The times of observations in the Kepler data is in Barycentric Julian Date (BJD) which is the Julian Date that has been corrected for the effects of the Earth’s orbit. This correction will prevent Earth’s orbit from appearing in the O-C diagrams. During the eclipse timing study, systems with no or irregular variations could be seen. These O-C variations would range from 0 to approximately 30 seconds and appear with no recurring pattern. Systems that have O-C variations larger than 30 seconds and particularly those that exhibit periodic O-C variations that are suspected to be caused by the addition of a third body should be prioritised for futher investigation. However as apsidal motion may also be the cause of periodic variation (Beuermann, K. et al., 2010) it can’t be assumed that third bodies are the cause of the O-C variations. In the hunt for planets, small amplitude variations (i.e. a few minutes) are also prioritised over larger variations as larger objects (i.e. stars) will have more of an effect on binary stars than smaller objects (i.e. dwarf stars or planets).

3 Following Up on Identified O-C Diagrams

With an O-C diagram showing periodic variability, the next task is to try to determine the cause of the variability. In this study the software Systemic (Meschiari et al., 2009; Meschiari & Laughlin, 2010) has been used to model the system and estimate the characteristics of a third based on its effect on the binary stars. Systemic can be used to model eclipse and transit timing variations. In order to accurately determine the properties of any potential third bodies the mass for the primary and secondary stars must be estimated. Radial velocity data is needed to determine the mass of binary stars however these data can be difficult to obtain for Kepler stars without the use of a large telescope. As a result mass estimates for the binary stars were calculated based on the data in the Kepler database and the light curve of the system.

JKTEBOP (Southworth et al., 2004) was used to find the best fit to model the light curve in order to determine/estimate the parameters of the system such as the orbital period, mass ratio of the binary stars and inclination of the system. Other data for the binary star systems such as the V-K colour can be used to help validate and guide the mass estimates of the binary star. By using the colours of the system and the mass ratio and other property estimates from JKTEBOP a system can be set up in Systemic to determine the properties of a potential third body. If the colours fall between two star types, the larger masses can be used and as a result the mass of any third body present should be an upper estimation of the mass. JKTEBOP was selected as it is capable of fitting the parameters of a system, including limb darkening and mass ratio.

4 Results

Having processed the detached eclipsing binary stars from the Kepler ‘Eclipsing Binary Catalog’, the O-C diagrams were then classified. One of the systems identified from the O-C diagrams as a potential host to a third body was KIC 5095269. The primary eclipse times and the errors reported by the fitting function from BET can be seen in Table  1 and an example eclipse fit from BET is shown in Fig.  1. Secondary eclipses were too shallow and unable to be fit and have accurate times determined.

Using the data found in Table  1 and fitting using the functions found in  (2) the period of the system, P, was found to be 18.611957 days and the initial primary eclipse, T0p, occurred at 133.866170 (BJD - 2,454,833). The secondary eclipses were found to be too shallow to accurately fit and as such no secondary eclipse times were obtained. The O-C diagram and data for the primary eclipses in KIC 5095269 is shown by the circles in Fig.  2. The O-C diagram shows periodic variability that has a period of approximately 120 days with variations in the eclipse times of up to approximately 2 minutes.

In addition to the orbital period of 18.611957 days and an initial eclipse occuring at 133.866170 (BJD - 2,454,833), modelling the light curve in JKTEBOP (Southworth et al., 2004) found a mass ratio of approximately 0.421 and an inclination of 80.02 degrees. An eclipse from the KEPLER data with the modelled light curve from JKTEBOP can be seen in Fig.  3. With a KEPLER magnitude of 13.528 and a 2MASS K magnitude of 12.215 this system has a V - K value of 1.313. This V - K value approximately matches an F7 star (Bessell & Brett, 1988). If the secondary star were much hotter than an M star it would have an impact on the K magnitude of the system, and therefore the V - K value. This is not consistent with observations. With an F7V primary star with a mass of 1.21$M_{\odot}$ Pecaut & Mamajek (2013) and a mass ratio of 0.421, the secondary star would have a mass of 0.51$M_{\odot}$. This is consistent with the mass of an M star as the V - K colour suggests.

Systemic (Meschiari et al., 2009; Meschiari & Laughlin, 2010) was then used to set up a representative system with the main (larger) star being set to a mass of 1.21 Solar Masses. ‘Planet 1’ was set up with the characteristics of the secondary star and the system i.e. a mass of 0.51 Solar Masses, an orbital period of 18.611957 days and an inclination of 80.0235 degrees. The masses as well as the orbital period and inclination of the binary stars were fixed while all the other parameters were free to be fit by the program. An additional planet (‘Planet 2’) was added to the system with the period of the O-C variability set as the period of the planet. All parameters for the additional body were also free to be fit. The eclipse times were loaded into Systemic and the uncertainties in the eclipse times were doubled in Systemic in order to estimate the true uncertainty in the eclipse times. A best-fitting was then performed by Systemic to find the values for the system that best explains the O-C variability.

The results of the Systemic best-fitting can be seen in Fig.  2, the residuals of the fit are shown in Fig.  4 and the data for the binary system is shown in Table  2 while the data for the third body is shown in Table  3. The properties of the binary system were entered into PHysics Of Eclipsing BinariEs or PHOEBE (Degroote, P. et al., 2013; Prsa & Zwitter, 2005) to view a synthetic light curve. PHOEBE suggested an inclination of approximately 86.5 degrees was required to view a primary eclipse but no/minimal secondary eclipse as can be seen in the light curve. The inclination of the system was kept at the 80.02 degrees suggested by JKTEBOP as this was found by fitting the light curve however the larger inclination from PHOEBE was noted.

5 Discussion

The results of the Systemic fit indicate that a third body with an orbital period of $237.70817\pm 0.12213$ days and a mass of $7.70\pm 0.08$ Jupiter masses could account for the eclipse timing variations seen (all errors are quoted to a single standard deviation unless otherwise noted). The mass of this third body is expected to be an upper limit and suggests the third body is actually a planet orbiting the binary stars. The best fit orbit of the binary stars is 18.61196 days. The eccentricity of the binary stars and planet as determined by Systemic was found to be 0.246 and $0.060\pm 0.003$ respectively. The data produced from Systemic found the best fit to match the eclipse timing variations occurs when there is a third body with a mass of $7.70\pm 0.08$ Jupiter masses. This mass is below the proposed planet/brown dwarf boundary of roughly 13 Jupiter masses (Burgasser, 2008) and as such the third body is regarded as a planet rather than a brown dwarf. Though how this object formed may determine whether it is a large planet or a small brown dwarf (Nordlund, 2011). If this third body is confirmed to be a circumbinary planet it would join a small number of previously confirmed circumbinary planets (Sigurdsson et al., 2003; Correia, A. C. M. et al., 2005; Lee et al., 2009; Doyle et al., 2011; Kostov et al., 2016), and this body would have one of the largest masses of these circumbinary planets.

With the period of the O-C variability being approximately 120 days it was expected that the orbital period of the third body would be approximately the same. However the Systemic fit was significantly better (both visually and by reduced $\chi^{2}$ value) with a period of 237.7 days and therefore this fit was chosen as the optimal fit. Systemic incorporates non-keplerian dynamics found in Fabrycky (2010) and it is thought that the orbital characteristics of the planet (visualised in Fig.  5) coupled with non-keplerian dynamics are the reason the orbital period of the planet varies from what was expected and also accounts for the variation seen in the maxima and minima of the O-C diagram in Fig.  2. The orbital period of the planet approximately matches the period of the variation in the maxima and minima of the O-C diagram.

The found system was tested in order to check the robustness of the fit. By testing the binary star system in PHOEBE we can confirm that the values for the binary star system are reasonable. The inclination of the system was changed in Systemic to 86.5 degrees as found by PHOEBE. The results of the third body remain consistent with the mass of the third body changing to 7.72 Jupiter masses. The next test involved changing the mass of the primary and secondary stars to determine if the best-fit system correspondingly changed the mass of the third body. The mass ratio between the stars and the planet was found to remain the same regardless of the actual masses used. During the robustness test, the mass ratio of the stars and planet was also changed. In all cases a third body with significantly less mass than the binary stars was able to account for the eclipse timing variations seen. With a mass ratio between the binary stars of 0.7889 (i.e. the primary remains at our estimated value and the secondary star mass set to the largest value allowed by Systemic for additional bodies of 1000 Jupiter masses) the best fit mass for the third body was found to be 10.35 Jupiter masses which is still below the putative planet/brown dwarf boundary of 13 Jupiter masses. With a mass ratio 0.01 (i.e the secondary star has a mass that is 0.01 times the mass of the primary star) the best fit mass for the third body is 5.04 Jupiter masses. In all of the tests performed, with varying properties of the binary system, a planetary mass third body is capable of producing the eclipse timing variations seen. The effect of limb darkening on the mass ratio of the binary stars was also analysed with JKTEBOP. It was found that as the amount of limb darkening of both stars increased, the mass ratio decreased. It was also found that increasing the amount of limb darkening on the primary star only, resulted in a lower mass ratio, while increasing the amount of limb darkening on the secondary star only, had a very minor effect and slightly increased the mass ratio.

The stability of the proposed orbits are important in determining whether the proposed orbits are the correct interpretation of the eclipse timing variations (Hinse et al., 2014). The system was integrated over a period of $10^{7}$ years to determine the long term stability of the system. Systemic (Meschiari et al., 2009; Meschiari & Laughlin, 2010) was used to perform the integration on the best-fit system found. The eccentricity results of the integration can be found in Fig.  6 while the semi-major axis results of the integration can be found in Fig.  7. The semi-major axis of the planet was found to vary between 0.795AU and 0.805AU while the eccentricity was found to fluctuate between 0.05 and 0.13 and indicates the planet orbits the binary stars in an almost circular orbit. As the orbits were found to be stable over large time periods, a planet in the proposed configuration is unlikely to be ejected from the system and is therefore the likely source of the eclipse timing variations seen (Hinse et al., 2014).

With a proposed inclination of $105.92\pm 1.45^{\circ}$ the probability of transits occurring needs to be considered. Kane & von Braun (2008) presented an analysis of the effect of orbital parameters, specifically eccentricity and argument of periastron, on the probability of a transit as a function of the orbital period. The probability of a transit occurring for a planet in a circular orbit drops dramatically as the orbital period of the planet increases. As the planet found in this system has a proposed orbit of approximately 237 days and a nearly circular orbit the expected transit probability is less than approximately 0.01. By viewing the light curve no transits can be seen to occur which can be expected with such a low probability of a planet in this orbit transiting the parent stars.

6 Conclusions

In this paper, we presented the evidence for a third body below whose mass the planet/brown dwarf boundary (Burgasser, 2008) around KIC 5095269. An eclipse timing variation study was performed on the KEPLER detached binaries where KIC 5095269 was found to exhibit periodic eclipse time variations. As eclipse time variations may be the result of a third body (Beuermann, K. et al., 2010) estimates for the mass of the binary stars were calculated and a model was produced with Systemic (Meschiari et al., 2009; Meschiari & Laughlin, 2010).

The model produced by Systemic suggests the cause of the eclipse timing variations is a third body with a mass of $7.70\pm 0.08$ Jupiter masses. Based on the proposed planet/brown dwarf mass boundary (Burgasser, 2008) we propose that this third body is a planetary candidate however how this object formed may determine whether it is a planet or a brown dwarf (Nordlund, 2011). The system was also found by Systemic to be stable over a period of $10^{7}$ years with an eccentricity of the third body that varies between 0.05 and 0.13. No transits could be seen to occur within the light curve although a planet with an orbital period of approximately 237 days has a probability of transit of less than 0.01.

In the future, we will be attempting to obtain Radial Velocity data for the system in order to further constrain the properties of the system. This planetary candidate has provided us with a template, both the features of O-C variations and the method for analysis, to use for future systems with periodic variations. Systems with low mass secondary stars, combined with small O-C variations are the best chance for detecting planetary sized bodies via the eclipse time variation method.

Acknowledgements

We would like to thank the referee and editor for their feedback and for helping us to fine tune the details of the paper. We would also like to thank Stefano Meschiari for his assistance with Systemic.