Transits of Inclined Exomoons - Hide and Seek and an Application to Kepler-1625

TL;DR

This paper analyzes the transit probabilities of inclined exomoons, demonstrating that many may frequently avoid transit, which impacts detection strategies and the interpretation of non-transit observations, especially in systems like Kepler-1625.

Contribution

It provides an analytic equation for exomoon transit probability and explores how orbital inclination and misalignment affect transit detection likelihood.

Findings

Exomoons can often avoid transiting if their orbits are large and misaligned.

Approximately 40% transit probability for the Kepler-1625 exomoon candidate.

Multiple non-transits are needed to confidently rule out the exomoon.

Abstract

A Neptune-sized exomoon candidate was recently announced by Teachey & Kipping, orbiting a 287 day gas giant in the Kepler-1625 system. However, the system is poorly characterized and needs more observations to be confirmed, with the next potential transit in 2019 May. In this Letter, we aid observational follow up by analyzing the transit signature of exomoons. We derive a simple analytic equation for the transit probability and use it to demonstrate how exomoons may frequently avoid transit if their orbit is larger than the stellar radius and sufficiently misaligned. The nominal orbit for the moon in Kepler-1625 has both of these characteristics, and we calculate that it may only transit roughly 40% of the time. This means that approximately six non-transits would be required to rule out the moon's existence at 95% confidence. When an exomoon's impact parameter is displaced off the star, the planet's impact parameter is displaced the other way, so larger planet transit durations are typically positively correlated with missed exomoon transits. On the other hand, strong correlations do not exist between missed exomoon transits and transit timing variations of the planet. We also show that nodal precession does not change an exomoon's transit probability and that it can break a prograde-retrograde degeneracy.

Figures

Tables

Table 1. Table 1: Parameters of the Kepler-1625 Exomoon candidate system

Param.	Unit	Value	$1\sigma$ Min	$1\sigma$ Max	Note
Host star
$m_{\star }$	$\left(M_{\odot }\right)$	1.04	0.98	1.12
$R_{\star }$	$\left(R_{\odot }\right)$	1.73	1.51	1.97
Planet
$m_P$	$\left(M_{\mathrm{Jup}}\right)$	6.85	1.2	12.5	(a)
$R_P$	$\left(R_{\mathrm{Jup}}\right)$	1.04	0.90	1.18
$T_P$	(days)	287.37278	287.37213	287.37353
$a_P$	(au)	0.87	0.85	0.89	(b)
$b_P$		0.104	0.038	0.188
$I_P$	(deg)	89.94	89.88	89.98	(c)
${\Omega }_P$	(deg)	0	0	0	(c)
Moon
$m_M$	$\left(M_{\oplus }\right)$	36.2	4.4	68	(a)
$R_M$	$\left(R_{\oplus }\right)$	4.90	4.18	5.69
$T_M$	(days)	22	13	39
$a_M$	(au)	0.022	0.017	0.030	(d)
$I_M$	(deg)	42	24	57	(e)
${\Omega }_M$	(deg)	0	-83	142	(e)
Relative orbit
$\left|90-I_M\right|$	(deg)	48	33	66	(f)

Equations

$$\footnotesize\left[\begin{array}[]{c}x_{\rm M}(f_{\rm M})\\ y_{\rm M}(f_{\rm M})\end{array}\right]=\left[\begin{array}[]{c}a_{\rm M}\left(\cos\Omega_{\rm M}\cos f_{\rm M}-\cos I_{\rm M}\sin\Omega_{\rm M}\sin f_{\rm M}\right)\\ a_{\rm M}\left(\sin\Omega_{\rm M}\cos f_{\rm M}+\cos I_{\rm M}\cos\Omega_{\rm M}\sin f_{\rm M}\right)+bR_{\star}\end{array}\right],$$

$$\cos \Delta I=\cos \Delta \Omega \sin I_{\mathrm{M}}\sin I_{\mathrm{P}}+\cos I_{\mathrm{M}}\cos I_{\mathrm{P}}.$$

$$\begin{array}{rl}{y}_{\mathrm{M}}(f_{\mathrm{M}})& =a_{\mathrm{M}}\sqrt{{\sin }^2{\Omega }_{\mathrm{M}}+{\cos }^2{I}_{\mathrm{M}}{\cos }^2{\Omega }_{\mathrm{M}}}\\ & \times \cos \left(f_{\mathrm{M}}-\arctan \left[\frac{\cos I_{\mathrm{M}}}{\tan {\Omega }_{\mathrm{M}}}\right]\right)+b{R}_{\star }\\ & \approx a_{\mathrm{M}}\left|\sin \Delta I\right|\cos \left(f_{\mathrm{M}}-\arctan \left[\frac{\cos I_{\mathrm{M}}}{\tan {\Omega }_{\mathrm{M}}}\right]\right)+b{R}_{\star },\end{array}$$

$$\footnotesize A=\left\{\begin{array}[]{ll}0&\mbox{if }b_{\rm P}R_{\star}+a_{\rm M}|\sin\Delta I|<R_{\star}\\ \\ \arccos\left[\frac{R_{\star}(1-b_{\rm P})}{a_{\rm M}\sin\Delta I}\right]&\mbox{if }b_{\rm P}R_{\star}+a_{\rm M}|\sin\Delta I|>R_{\star}\\ \end{array}\right.$$

$$\footnotesize B=\left\{\begin{array}[]{ll}180^{\circ}&\mbox{if }b_{\rm P}R_{\star}-a_{\rm M}|\sin\Delta I|>-R_{\star}\\ \\ \arccos\left[\frac{-R_{\star}(1+b_{\rm P})}{a_{\rm M}\sin\Delta I}\right]&\mbox{if }b_{\rm P}R_{\star}-a_{\rm M}|\sin\Delta I|<-R_{\star}\end{array}\right.$$

$${\tau }_{\mathrm{P}}=\frac{T_{\mathrm{P}}}{\pi }\arcsin \left(\frac{R_{\star }\sqrt{1-b_{\mathrm{P}}^2}}{a_{\mathrm{P}}}\right),$$

$$p_{\mathrm{M}}=\min \left[\frac{B-A}{180^{\circ }}+\frac{{\tau }_{\mathrm{P}}}{T_{\mathrm{M}}},1\right].$$

$$T_{\mathrm{prec}}=\frac{4}{3}\frac{m_{\mathrm{P}}+m_{\star }}{m_{\star }}\frac{T_{\mathrm{P}}^2}{T_{\mathrm{M}}}\frac{1}{\cos \Delta I}.$$

$$a_{\mathrm{M},\mathrm{crit}}={\left(\frac{2J_2{R}_{\mathrm{P}}^2{a}_{\mathrm{P}}^3{m}_{\mathrm{P}}}{m_{\star }}\right)}^{1/3},$$

$$I_{\mathrm{M}}(t)=\Delta I\cos \left[\frac{2\pi }{T_{\mathrm{prec}}}(t-t_0)\right]+I_{\mathrm{P}},$$

Full text

planets and satellites: detection, dynamical evolution and stability, fundamental parameters

††journal: ApJ Letters

1 Introduction

In the solar system our understanding of the planets is enriched by our understanding of their moons. The Moon is thought to influence Earth’s habitability (Laskar et al., 1993). The Galilean moons help constrain the early evolution of Jupiter (Heller et al., 2015; Ronnet et al., 2018). The equatorial alignment of Uranus’s moons helps us understand the origin of the planet’s tilt (Kegerreis et al., 2018). As a community we would benefit immensely from conducting similar science for moons of extrasolar planets (exomoons).

Detecting analogs of the solar system moons is challenging due to their small size. Photometry is thought to be the most promising technique (Kipping et al., 2009), either through observing individual moon transits (Sartoretti & Schneider, 1999), multiple averaged moon transits (Simon et al., 2012; Heller, 2014; Teachey et al., 2017) or inferring the moon’s existence based on the planet’s transit timing variations (TTVs) and transit duration variations (TDVs) (Sartoretti & Schneider, 1999; Kipping, 2009a, b, 2011; Heller et al., 2016). Other techniques with potential include gravitational microlensing (Bennett et al., 2014; Hwang et al., 2018) and observations of self-luminous giant exoplanets to detect a variation in polarization (Sengupta & Marley, 2016) or in radial velocity (Vanderburg et al., 2018).

The most plausible exomoon to date is in the Kepler-1625 system. The planet (Kepler-1625b) itself is unremarkable: a gas giant on a 287-day orbit. The surprise, however, is the size of the moon (Kepler-1625b-i, as it is potentially similar in mass and radius to Neptune. Such a large moon is without precedent in our solar system, but one must remember that so were the first exoplanet diskoveries.

The moon was originally suspected based on three planet transits within the original Kepler mission (Teachey et al., 2017; Heller, 2018). Asymmetries in the transit profile teased the presence of a moon, but neither TTVs nor TDVs were detected to confirm it. The moon’s existence became more likely after a fourth planetary transit was captured by the Hubble Space Telescope (HST) (Teachey & Kipping, 2018). The planet transit was 70 minutes early, although no TDV was detected. Furthermore, there is a shallow dip in the light curve after the egress of the planet transit: a potential moon transit. Table 1 contains basic system parameters used in our letter, but we refer the reader to Teachey & Kipping (2018) for significantly more detail.

In this letter we are agnostic about the reality of this particular exomoon. Both Teachey & Kipping (2018) and subsequent analysis by Heller et al. (2019) encourage new observations in order to consider the moon confirmed. In this letter we aid such future observations by analyzing the detectability of exomoons, both in general and for Kepler-1625b-i specifically. We quantify previous intuition that some moons are not guaranteed to transit every time their host planet does (Sartoretti & Schneider, 1999; Martin, 2017). Missed transits typically occur when the moon’s orbit is both wider than the stellar diameter and significantly misaligned to the planet’s orbital plane. The best-fitting, albeit loosely constrained orbit for Kepler-1625b-i has both of these characteristics. Furthermore, within our own solar system we know of Triton, which is on a highly misaligned, in fact retrograde, orbit (Fig. 1).

In this letter we derive an analytic transit probability for exomoons of transiting planets (Sect. 2) which accounts for both misalignment and a dynamically varying exomoon orbit. We then test the correlation between the presence/absence of moon transits and the TTV and TDV signature of the planet (Sect. 3). We apply our work to both exomoons in general and the Kepler-1625 system specifically (Sect 4). The letter ends with a brief discussion (Sect. 5).

2 Exomoon transit probability

2.1 Transit Geometry

The transit geometry is shown in Fig. 2. The observer looks from the positive $z$-axis at the $(x,y)$ sky plane centered on the star. The planet orbit is modeled by a straight line from left to right (positive $x$ direction), vertically offset by the impact parameter $b_{\rm P}=a_{\rm P}\cos I_{\rm P}/R_{\star}$. This assumes $a_{\rm P}\gg R_{\star}$111Very tight-orbiting planets are thought unlikely to host moons anyway (Namouni, 2010). and $m_{\rm P}\gg m_{\rm M}$222Care must be taken when generalizing our work to “binary planets” (Lewis et al., 2015), although our work is likely applicable to “moon-moons” (Forgan, 2018), “moon-moon-moons” or indeed moonn. and throughout this letter we also assume circular orbits, i.e. $e_{\rm P}=e_{\rm M}=0$. The planet’s orbit would be rotated clockwise by $\Omega_{\rm P}$, but we arbitrarily set $\Omega_{\rm P}=0$ as the transit geometry are only sensitive to $\Delta\Omega=\Omega_{\rm M}-\Omega_{\rm P}$.

The position of the moon at the time of the planet’s transit midpoint across the star is fundamental to the transit phenomenom. Neglecting eccentricity, its projected orbit is an ellipse with major axis $a_{\rm M}$ and a minor axis $a_{\rm M}|\cos I_{\rm M}|$, rotated counter-clockwise by $\Omega_{\rm M}$ and offset vertically by $b_{\rm P}R_{\star}$:

[TABLE]

where $f_{\rm M}$ is the true anomaly of the moon. It is important to remember that $f_{\rm M}$ is the orbital phase of the moon defined within its orbital plane, not with respect to our $(x,y)$ coordinate system. In Fig. 2 the moon is misaligned and prograde with the planet’s orbit and projects a counter-clockwise motion. In Fig. 6 we however note that a degeneracy exists between prograde and retrograde moons (see Sect. 5.1). The mutual inclination between the moon and the planet’s orbit (not equator) is

[TABLE]

The moon will transit the star on a given planet transit when $|y_{\rm M}(f_{\rm M})|<R_{\star}$. To make this criterion easier to solve, we consolidate the expression in Eq. 1 for $y_{\rm M}$ from two trigonometric functions of $f_{\rm M}$ to one:

[TABLE]

We note that whilst the second line of Eq. 3 contains $b_{\rm P}$, hence implying that $I_{\rm P}$ is not exactly $90^{\circ}$, the approximation $\sqrt{\sin^{2}\Omega_{\rm M}+\cos^{2}I_{\rm M}\cos^{2}\Omega_{\rm M}}\approx|\sin\Delta I|$ is derived from Eq. 2 using $I_{\rm P}=90^{\circ}$. However, the end result is a negligible difference between the two lines in Eq. 3.

The exomoon transit probability is calculated as the fraction of angles $f_{\rm M}$ that correspond to $|y_{\rm M}|<R_{\star}$. The phase shift of $\arctan[\cos I_{\rm M}/\tan\Omega_{\rm M}]$ in Eq. 3 does not affect this fraction, and hence we simplify Eq. 3 by defining $f^{\prime}_{\rm M}=f_{\rm M}-\arctan[\cos I_{\rm M}/\tan\Omega_{\rm M}]$. The function $y_{\rm M}(f^{\prime}_{\rm M})$ is symmetric over $f^{\prime}_{\rm M}=180^{\circ}$. Between 0 and $180^{\circ}$ we define the range of transits to be $[A,B]$, where

[TABLE]

where the second condition occurs when the moon can miss transit above the star (with respect to the $y$-axis, and

[TABLE]

where the second condition occurs when the moon can miss transit below the star.

If $f^{\prime}_{\rm M}$ (and hence $y_{\rm M}$) is static during the planet’s transit then the exomoon transit probability is simply the ratio $p_{\rm M}=(B-A)/180^{\circ}$. However, this static assumption is only applicable when $T_{\rm M}\gg\tau_{\rm P}$, where $\tau_{\rm P}$ is the planet’s transit duration:

[TABLE]

To approximately account for shorter-period moons we add to $p_{\rm M}$ the fraction of the orbit covered during the planet’s transit: $\tau_{\rm P}/T_{\rm M}$. With this, our derived exomoon transit probability is

[TABLE]

2.2 Orbital dynamics

The orbit of an exomoon may be subject to various dynamical perturbations. When the moon and planet orbits are misaligned, one such effect is a nodal precession induced by the three-body interactions between the sun, planet, and moon. From Mardling (2010) the rate of precession is

[TABLE]

This effect may be quenched by a competing torque on the moon’s orbit induced by the equatorial bulge of the planet. Burns (1986) calculate a critical moon semi-major axis, for which the dynamics of interior orbits are dominated by the planet’s equatorial bulge:

[TABLE]

where $J_{2}$ is the first gravitational harmonic. See also Boué & Laskar (2006); Tremaine et al. (2009) for more details. In this letter we are predominantly interested in moons that are long-period and misaligned (such that moon transits are sometimes missed) and planets that are short-period (so planet transits are more frequent). For such moons the dominant effect is a three-body nodal precession. The Earth’s moon exhibits three-body nodal precession with a period of 17.9 years (according to Eq. 8). For Kepler-1625 $a_{\rm M,crit}=0.008$ au, which is almost three times less than the nominal value $a_{\rm M}=0.022$ au, and hence we also expect three-body nodal precession in this system, with a calculated period of 20.5 years.

With respect to the orbital plane of the planet, which remains (essentially) fixed, nodal precession makes the moon orbit circulate at a constant rate given by Eq. 8, whilst maintaining a constant mutual inclination $\Delta I$. With respect to the observer, Martin (2017) showed that $I_{\rm M}$ librates over time $t$ around the constant $I_{\rm P}$ according to

[TABLE]

where $t_{0}$ corresponds to $I_{\rm M,0}$. With respect to the observer, $\Omega_{\rm M}(t)$ also librates and can be calculated by combining Eqs. 2 and 10.

A complication to the nodal precession arises in highly misaligned orbits, such that $|90^{\circ}-\Delta I|\lesssim 50^{\circ}$. In such cases Kozai-Lidov cycles occur, which cause $\Delta I$ and $e_{\rm M}$ vary, even for initially circular orbits (Lidov, 1961, 1962; Kozai, 1962).

The expression for $y_{\rm M}$ in Eq. 3, does depend on the time-dependent quantities $I_{\rm M}$ and $\Omega_{\rm M}$. However, these quantities only phase shift $f_{\rm M}$ and do not change the fractional range of $f_{\rm M}$ corresponding to transits, which is why they could be ignored when calculating the quantities $A$ (Eq. 4) and $B$ (Eq. 5). These quantities are functions of $\Delta I$, but this is constant333To be precise, $\Delta I$ is only constant under the secular regime, i.e. when calculations are made that average over the orbital periods. There do exist short-term variations on the timescales of $T_{\rm M}$ and $T_{\rm P}$, but these are on order $\approx 2\%$ variations. for orbits that are circular and without Kozai-Lidov cycles. Overall, we demonstrate that in our simplified setup the exomoon transit probability $p_{\rm M}$ is constant during the moon’s precession period.

2.3 Accuracy of the analytic solution

We run $n$-body simulations for a suite of 1000 randomized transiting planet plus moon systems and calculate a numerical transit probability as the ratio of moon to planet transits. The masses are $m_{\star}=1M_{\odot}$, $m_{\rm P}\in[0.1,3]M_{\rm Jup}$, $m_{\rm M}\in[0.1,17]M_{\oplus}$, using log-uniform distributions. The planet radii are calculated using the mass-radius relation of Bashi et al. (2017): $R/R_{\oplus}=(m/M_{\oplus})^{0.55}$ for $m/M_{\oplus}<124$ and $R/R_{\oplus}=(m/M_{\oplus})^{0.01}$ for $m/M_{\oplus}\geqslant 124$. The orbital parameters for the planet are $T_{\rm P}\in[200,500]$ days, $e_{\rm P}=0$, $\Omega_{\rm P}=0$, $b_{\rm P}\in[0,0.9]$ and $f_{\rm P}\in[0^{\circ},360^{\circ}]$. The orbital parameters for the moon are $T_{\rm M}\in[1,50]$ days, $e_{\rm M}=0$, $f_{\rm M}\in[0^{\circ},360^{\circ}]$. The mutual inclination is drawn from $\Delta I\in[0^{\circ},40^{\circ}]$. We randomly choose the starting phase of the precession period by calculating $I_{\rm M}$ in Eq. 10 with a uniformly random phase between $\in[0^{\circ},360^{\circ}]$ and $I_{\rm P}$ calculated from the randomly chosen $b_{\rm P}$. We then calculate $\Omega_{\rm M}$ from Eq. 2.

Each simulation is run over a time span of $100\times T_{\rm P}$ using a fourth-order Runge-Kutta integrator with a fixed step size of 30 minutes, chosen to match Kepler’s long-cadence observations. Across all 1000 simulations, the median percentage error between the analytic and numerical transit probabilities is 1.2%. For 626 of the simulations the numerical transit probability is less than 1 (i.e. at least one missed moon transit), and for these simulations the median error is 4.0%. Contributions to the error include perturbations to the moon’s orbit, mean motion resonances, other period-ratio effects which may alias the moon transit sequence, any simplifications in the derivation of Eq. 7 and counting statistics of the numerically calculated transit probability.

3 Planet transit timing and duration variations

An isolated, unperturbed planet would transit the star with perfect periodicity, $T_{\rm P}$. However, the presence of the moon can induce TTVs and TDVs on the planet. The main cause is a small “wobble” of the planet around the planet-moon barycenter, on top of the planet’s larger-scale orbit around the star-planet barycenter. This is a Keplerian effect (i.e., it occurs with static orbits). We briefly diskuss the origin of the barycentric TTVs and TDVs in Sect. 3.1, and direct the reader to the seminal papers of Kipping (2009a, b) for a much more thorough treatment, included detailed analytic equations. A secondary contribution to TTVs and TDVs is from non-Keplerian effects, i.e. perturbations to the orbital elements. We do not diskuss these effects but they are naturally included in our $n$-body simulations. Finally, we do not diskuss the TTVs and TDVs of the moon itself, but they are expected to significantly larger than those of the planet.

3.1 Origins of barycentric TTVs and TDVs

A planet exhibits a TTV when slightly offset along the horizontal axis (i.e. parallel with its transit chord). This change adds or subtracts to the time taken to reach the transit midpoint. A horizontal offset is induced by the planet’s wobble around the planet-moon barycenter. The TTV is calculated as the time taken for the planet to traverse this offset at its orbital velocity around the star of $v_{\rm P,\star}=2\pi a_{\rm P}/T_{\rm P}$.

A planet exhibits a TDV for two different reasons. First, the planet’s motion around the planet-moon barycenter has a velocity $v_{\rm P,M}=2\pi a_{\rm M}m_{\rm M}/[(m_{\rm P}+m_{\rm M})T_{\rm M}]$. The horizontal component of this velocity may be additive or subtractive to $v_{\rm P,\star}$, and hence when the planet transits it may be moving a little faster or slower than average, causing the transit duration to vary. Kipping (2009b) called this the “V-TDV”.

The second cause of a TDV is a vertical offset of the planet’s position (i.e. perpendicular to its transit chord) due to the barycentric reflex motion induced by the moon. This changes $b_{\rm P}$, hence changing $\tau_{\rm P}$ by Eq. 6. Kipping (2009b) called this the “TIP-TDV”.

3.2 Connecting TTVs and TDVs with moon transit occurrence

We use the Sect 2.3 $n$-body simulations to test the correlation between moon transits and planet TTVs and TDVs. We only take the 626/1000 simulations which have at least one missed moon transit. For each simulation we calculate numerically the TTVs and TDVs, which we scale by dividing each value by the maximum absolute value for the simulation. We collate the scaled TTVs and TDVs for the simulations, separate them by moon transit occurrence, and show the results in a histogram in Fig. 3 (left).

For TTVs there is typically no difference between when the moon does and does not transit. There are two main reasons for this. First, occurrence of a moon transit is a function of its vertical position ($y_{\rm M}$), yet the TTV signal is a function of the moon’s horizontal position ($x_{\rm M}$). Consider Fig. 2. A positive $x_{\rm M}$ displaces the planet to the left and hence induces a positive TTV (late transit), and vice-versa. We see that positive $x_{\rm M}$ values correspond to both cases where the moon does and does not transit (only misses above the star). Negative $x_{\rm M}$ values largely correspond to the moon transiting, but there is also a small parameter space for missing transits, both above and below the star. In Fig. 2, when averaged over all $x_{\rm M}$ there will be preference for missed transits to correspond to positive values of $x_{\rm M}$, and hence positive TTVs. However, this trend will be weak except for small $a_{\rm M}/R_{\star}$, and in that case it would be rare for the moon to avoid transit anyway. The second consideration is that nodal precession of the moon rotates its orbit. After $0.5T_{\rm prec}$ the moon orbit in Fig. 2 will be mirrored horizontally, in which case missed moon transits will now typically correspond to negative values of $x_{\rm M}$. Our $n$-body simulations cover multiple precession periods, and hence any small short-term TTV-moon transit correlations are averaged out.

For TDVs the results contrastingly show a clear difference in the TDV distribution with and without moon transits. This matches Fig. 2; the moon misses transit when in the uppermost and lowermost parts of its orbit, but the upper region is larger due to the asymmetric vertical offset. When the moon is in this upper region the planet is displaced slightly downward toward the stellar center and hence takes longer to transit (a positive TDV). This does not change throughout the nodal precession period.

The TDV-transit correlation is only prominent when $b_{\rm P}$ is significantly non-zero. In Fig. 3 (right) we split the simulations into $b_{\rm P}\in$ $[0,0.3]$, $[0.3,0.6]$, and $[0.6,0.9]$. The correlation between TDVs and moon transits disappears for small impact parameters. There are two reasons for this. First, for the same vertical offset induced by the moon the change in the path length across the star is less when the planet passes near the stellar center rather than near the limb. Second, at small $b_{\rm P}$ the moon’s orbit across the star is nearly symmetric vertically, and hence is nearly equally likely to miss transit above or below the star (unlike in Fig. 2).

The TDVs for small $b_{\rm P}$ are largely caused by the velocity change effect, which dependent on the horizontal position of the moon and hence is not strongly correlated with the presence of moon transits.

In Fig. 4 we show TTVs and TDVs for three example simulations. The sole change is $b_{\rm P}=0.1$, 0.4 and 0.7. The planet TTV signal remains constant, although the sequence of moon transits changes. The TDV signal at small $b_{\rm P}$ is small in amplitude with no correlation with the moon transits. As $b_{\rm P}$ increases so does the TDV amplitude and the moon transit correlation.

The impact parameter of the planet Kepler-1625b is well constrained to be small: $b_{\rm P}=0.104^{+0.084}_{-0.066}$. We therefore expect TDVs to be small and uncorrelated with missed moon transits, and indeed no TDVs have been observed so far.

4 Applications

4.1 Transit probability of hypothetical exomoon systems

The transit probability for the moon is a function of $\Delta I$, $b_{\rm P}$, and $a_{\rm M}/R_{\star}$. Fig. 5 (a) shows $p_{\rm M}$ (Eq. 7) over a wide range of parameters: $\Delta I\in[0,40^{\circ}]$, $b_{\rm P}\in[0,1]$, and $a_{\rm M}/R_{\star}=0.5,1,1.56,4$.

For $a_{\rm M}/R_{\star}<1.56$ the transit probability is 1 except for high values of $\Delta I$ and/or $b_{\rm P}$, where the probability goes to a minimum of 0.5. The parameter space where $p_{\rm M}<1$ increases as $a_{\rm M}/R_{\star}$ increases. When $a_{\rm M}/R_{\star}>1.56$ the moon’s orbit is so wide that its vertical extent exceeds the stellar diameter and $p_{\rm M}<0.5$ for some $\Delta I$ and $b_{\rm P}$.

4.2 Transit probability of Kepler-1625b-i

In Fig. 5 (b) we calculate $p_{\rm M}$ for Kepler-1625b-i over a plausible range of $a_{\rm M}$ and $\Delta I$, whilst fixing $b_{\rm P}=0.1$ and $R_{\star}=1R_{\odot}$. Note that when calculating the nominal value of $\Delta I$ we take $\Omega_{\rm M}=0^{\circ}$ and then $\Delta I\approx 90^{\circ}-I_{\rm M}$ from Eq. 2. This means $\Delta I=48^{\circ}$, which places the system just within the nominal Kozai-Lidov regime, but the $\Delta I$ and $e_{\rm M}$ variations should be small enough for our equations to remain applicable.

The Teachey & Kipping (2018) nominal values correspond to $p_{\rm M}=0.4$, although this probability varies significantly within the $1\sigma$ error bounds, and they note that the moon could still have a coplanar orbit, which would mean $p_{\rm M}=1$.

5 Discussion

5.1 Breaking the prograde/retrograde degeneracy

Observations of a moon that orbits a planet on a non-evolving orbit are subject to a degeneracy between prograde ($\Delta I<90^{\circ}$) and retrograde ($\Delta I>90^{\circ}$) orbits. This degeneracy is shown in Fig. 6. Two orbits are shown: one in solid red that is prograde and coplanar ($\Delta I=0^{\circ}$, red solid line), and one in dashed black that is retrograde but misaligned ($90^{\circ}<\Delta I<180^{\circ}$). Both orbits yield the same projected $x$ and $y$ positions and $v_{x}$ and $v_{y}$ velocities of the moon; hence, the Keplerian TTV and TDV phenomenology would be the same. However, the side view (left) reveals a clear difference in the two moon orientations.

This degeneracy may be broken by nodal precession, which would not occur for the coplanar orbit but would for the misaligned orbit. Fortunately, for a moon that orbits at a fair fraction of its planet’s Hill sphere, precession will be rapid, revealing the magnitude of the misalignment in just tens of orbits of the planet. Therefore, the dynamically evolving character of TDV will betray the prograde or retrograde character of the moon.

If the planetary impact parameter is low then the “TIP-TDV” may be negligible and the magnitude of non-coplanarity may not be enough to break the degeneracy. In this case, higher-order dynamical effects that differ in sign between prograde and retrograde moons may need to be taken into account, as envisioned by Lewis & Fujii 2014. Two alternative methods for breaking the degeneracy, practical only with Extremely Large Telescopes, were diskussed by Heller & Albrecht (2014).

5.2 The prevalence of large TTVs for long-period gas giants

According to the transit times of Table S3 of Teachey & Kipping (2018), the planet Kepler-1625b has a mean absolute deviation from a constant-period model, normalized by the orbital period — a “scatter” — of $s_{O-C}/T_{\rm P}=2.40\times 10^{-5}$. The timings have a median error bar normalized by the orbital period of $\sigma/T_{\rm P}=1.55\times 10^{-5}$. For the TTV measurements of Holczer et al. (2016), the data are more precise than that for 40 planets with $T_{\rm P}>100$ days. Of those 40, 15 planets have larger TTV scatter, i.e., $s_{O-C}/T_{\rm P}>2.40\times 10^{-5}$, and all of these are deemed significant at $\log p<-8.8$. The large amplitude and period of these signals makes them likely due to planet-planet perturbations. We conclude that Kepler-1625b may very likely have a TTV signal due to additional planets, which may be confused for exomoons, or at least contaminate the exomoon TTV signal. A repeated photometric transit signal of the exomoon, rather than the TTV induced on the planet, is likely a more reliable signature.

5.3 Overlapping moon and planet transits

There are two possible scenarios for overlapping moon and planet transits. First, the moon may be entirely in front of or behind the planet, in which case the photometric signal would be identical to that of an isolated planet transit and the moon would be hidden. Such an event is not explicitly considered in our equations. We estimate it to be rare though, with a likelihood on the order of $\approx R_{\rm P}/a_{\rm M}$ if $I_{\rm M}=I_{\rm P}=90^{\circ}$, and significantly less for inclinations that allow the moon to be offset vertically from the planet at transit. Second, the moon and planet may pass the star at the same time, but with different impact parameters. In this case their photometric dips would be additive and, if telescope precision allowed, a distortion in the transit shape may be detected. Such an event would be covered in our equations for $p_{\rm M}$. Exotic syzygies such as this are treated in more detail in Kipping (2011); Veras & Breedt (2017); Veras (2019).

5.4 Future observing prospects

The most effective way to confirm and characterize the Kepler-1625 system is through continued transit photometry. Even if the moon only transits $\approx 40\%$ of the time as we predict, additional planet transits will provide new TTV measurements, although probably not new TDV measurements due to the planet’s small impact parameter. The next planet transit is scheduled for 2019 May 26. Fig. S18 of Teachey & Kipping (2018) predicts when the moon will transit. Most of their models show a moon transit before the planet’s ingress, but they do not quantify the chance of the moon missing transit 444At the American Astronomical Society Meeting 233, Seattle, 2019 January, Alex Teachey’s presentation noted a simulated moon transit chance of 53% for 2019 May. This is slightly above our 40% analytic calculation, but our calculations are an average over all transit epochs, not any specific one..

The Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2014) can feasibly observe the planet transit on 2022 July 19 and 2026 June 25, but at Jmag = 14.4 the transit will only be observed at a signal-to-noise ratio of 2.5, which is insufficient for transit timing or moon spotting.

The James Webb Space Telescope (JWST) will provide superior photometric precision to HST (Beichman et al., 2014). From its observing constraints, JWST can observe Kepler-1625 annually from April 22 to November 14, meaning the first planet transits observable with this facility will occur on 2021 October 5, 2022 July 19, and 2023 May 3. With JWST, the transit timing will likely be limited by our abilities to model the granulation features on the stellar surface, which induce significant correlated noise on $\approx 20$ minute timescales given the subgiant nature of this star. Transits of a moon signal of the amplitude and duration claimed by Teachey & Kipping (2018) will be detectable at the 3$\sigma$ level.

If the moon does not exist, then a binomial test reveals how many non-transits are required to prove this to a certain significance. This assumes that each moon transit would have been detectable and that the transit probability of individual moon transits is independent for each planet transit, which neglects mean motion resonances. The probability of $n$ undetected transits is $p_{n}=(1-p_{\rm M})^{n}$. With our estimated $p_{\rm M}=0.4$, for a $95\%$-confident non-detection we solve $(1-0.95)=(1-0.4)^{n}$ to obtain $n\sim 6$ well-surveyed yet undetected exomoon transits. If the moon does exist, then a similar number of transits would be also be needed to well characterize its orbit.

Acknowledgements: We thank the referee for thoroughly reviewing our paper and providing comments which significantly improved its quality.