By the end of this tutorial, you will:
The Kepler and TESS missions are optimized for finding new transiting exoplanets. Lightkurve provides a suite of tools that help make the process of identifying and characterizing planets convenient and accessible.
In this tutorial, we will show you how to conduct your own search for transiting exoplanets in Kepler and TESS light curves. Lightkurve uses the Astropy implementation of the Box Least Squares (BLS) method to identify transit signals. This tutorial demonstrates the basics of how to optimally use Lightkurve's BLS tools.
import lightkurve as lk %matplotlib inline
# Search for Kepler observations of Kepler-69 search_result = lk.search_lightcurve('Kepler-69', author='Kepler', cadence='long') # Download all available Kepler light curves lc_collection = search_result.download_all() lc_collection.plot();
Each observation has a different offset, so in order to successfully search this light curve for transits, we first need to normalize and flatten the full observation. This can be performed on a stitched light curve. For more information about combining multiple observations of the same target, please see the companion tutorial on combining multiple quarters of Kepler data with Lightkurve.
# Flatten the light curve lc = lc_collection.stitch().flatten(window_length=901).remove_outliers() lc.plot();
The most common method used to identify transiting exoplanets is the Box Least Squares (BLS) periodogram analysis. BLS works by modeling a transit using an upside-down top hat with four parameters: period, duration, depth, and reference time. These can be seen in the figure below, from the astropy.timeseries implementation of BLS.
These parameters are then optimized by minimizing the square difference between the BLS transit model and the observation. For more information about BLS, please see the Astropy documentation.
Lightkurve has two types of periodogram available to anaylze periodic trends in light curves:
Please see the companion tutorial on how to create periodograms and identify significant peaks for an example of the
To create a
BoxLeastSquaresPeriodogram, use the
LightCurve method to_periodogram, and pass in the string
'bls' to specify the type of periodogram object you want to create. This method also optionally takes an array of periods (in days) to search, which we will set from 1–20 days to limit our search to short-period planets. We do so using the numpy.linspace function.
import numpy as np # Create array of periods to search period = np.linspace(1, 20, 10000) # Create a BLSPeriodogram bls = lc.to_periodogram(method='bls', period=period, frequency_factor=500); bls.plot();
The plot above shows the power, or the likelihood of the BLS fit, for each of the periods in the array we passed in. This plot shows a handful of high-power peaks at discrete periods, which is a good sign that a transit has been identified. The highest power spike shows the most likely period, while the lower power spikes are fractional harmonics of the period, for example, 1/2, 1/3, 1/4, etc.
We can pull out the most likely BLS parameters by taking their values at maximum power — we will refer to this transiting object as "planet b."
planet_b_period = bls.period_at_max_power planet_b_t0 = bls.transit_time_at_max_power planet_b_dur = bls.duration_at_max_power # Check the value for period planet_b_period
To confirm that this period and transit time (epoch) correspond to a transit signal, we can phase-fold the light curve using these values and plot it.
ax = lc.fold(period=planet_b_period, epoch_time=planet_b_t0).scatter() ax.set_xlim(-5, 5);
The phase-folded light curve shows a strong transit signal with the identified period and transit time of maximum BLS power.
The BLS periodogram has features that make it possible to search for multiple planets in the same system. If we want to identify additional transit signals, it will be much more convenient if we first remove the previously identified signal. This will prevent the high-power periodicity of the first planet, planet b, from dominating the BLS periodogram, and will allow us to find lower signal-to-noise ratio (SNR) transits.
We can create a cadence mask for the light curve using the transit parameters from the
# Create a cadence mask using the BLS parameters planet_b_mask = bls.get_transit_mask(period=planet_b_period, transit_time=planet_b_t0, duration=planet_b_dur)
Now, we can create a masked version of the light curve to search for additional transit signals. The light curve is shown below, with masked cadences marked in red.
masked_lc = lc[~planet_b_mask] ax = masked_lc.scatter(); lc[planet_b_mask].scatter(ax=ax, c='r', label='Masked');
We can also create a BLS model to visualize the transit fit. This returns a
LightCurve object with the BLS model in the flux column.
# Create a BLS model using the BLS parameters planet_b_model = bls.get_transit_model(period=planet_b_period, transit_time=planet_b_t0, duration=planet_b_dur)
We can plot this over the folded light curve to confirm that it accurately represents the transit.
ax = lc.fold(planet_b_period, planet_b_t0).scatter() planet_b_model.fold(planet_b_period, planet_b_t0).plot(ax=ax, c='r', lw=2) ax.set_xlim(-5, 5);
Now that we have created a light curve with the first identified planet masked out, we can search the remaining light curve for additional transit signals. Here, we search for long-period planets by increasing our range of periods to 1–300 days.
period = np.linspace(1, 300, 10000) bls = masked_lc.to_periodogram('bls', period=period, frequency_factor=500) bls.plot();
While no peaks in this BLS periodogram display a power as high as the previous transit signal, there is a definite peak near ~240 days. We can pull out the corresponding period and transit time to check the signal.
planet_c_period = bls.period_at_max_power planet_c_t0 = bls.transit_time_at_max_power planet_c_dur = bls.duration_at_max_power # Check the value for period planet_c_period
We can again plot the phase-folded light curve to examine the transit.
ax = masked_lc.fold(planet_c_period, planet_c_t0).scatter() masked_lc.fold(planet_c_period, planet_c_t0).bin(.1).plot(ax=ax, c='r', lw=2, label='Binned Flux') ax.set_xlim(-5, 5);
This signal is lower SNR because there are fewer transits due to the longer period, and the shallower depth implies that the planet is smaller. To help see the transit more clearly, we have overplotted the binned flux, combining consecutive points taken over a span of 0.1 days.
We have now successfully identified two planets in the same system! We can use the BLS models to visualize the transit timing in the light curve.
planet_c_model = bls.get_transit_model(period=planet_c_period, transit_time=planet_c_t0, duration=planet_c_dur)
ax = lc.scatter(); planet_b_model.plot(ax=ax, c='dodgerblue', label='Planet b Transit Model'); planet_c_model.plot(ax=ax, c='r', label='Planet c Transit Model');
Lightkurve also has a tool that enables you to interactively perform a BLS search. A quick demo of this feature is shown below.
To use the LightCurve.interact_bls() method, zoom in on peaks in the BLS periodogram using the interactive plotting tools. To improve the fit, you can change the transit duration. The phase-folded light curve panel in the top right and the full light curve below it will automatically update to plot the highest power BLS model. The BLS parameters with highest power are noted in the bottom right of the figure.
If you use
lightkurve or its dependencies in your published research, please cite the authors. Click the buttons below to copy BibTeX entries to your clipboard.