<jupyter_start><jupyter_text>Transit Least Squares (TLS) TutorialIn this tutorial, we will demonstrate the usage of the *TLS* algorithm for transit detection, which is more sensitive than the Box Least Squares (BLS) method. The transit periodogram computes the power spectrum of the transit fit, maximized over transit depth, duration, and phase, for a list of periods. Peaks in the periodogram generally indicate transiting planets or other unmodeled noise.To demonstrate the code, we start by downloading the K2 light curve for K2-3: A system with [3 known transiting planets](http://iopscience.iop.org/article/10.1088/0004-637X/804/1/10/meta). This tutorial is inspired by <NAME>'s excellent [BLS tutorial](https://gist.github.com/dfm/96e30d891c0d351daae9aaaef56d3e78).<jupyter_code>import numpy
import scipy
import everest
from astropy.stats import sigma_clip
import matplotlib.pyplot as plt
import logging
logging.basicConfig(level=logging.DEBUG)
logger = logging.getLogger('matplotlib')
logger.setLevel(logging.CRITICAL)<jupyter_output><empty_output><jupyter_text>We download the data using the [Everest](https://github.com/rodluger/everest) package and apply a sliding median filter (blue line) to divide out the stellar and instrumental variability. We also clip positive 3-sigma outliers, but leave negative outliers untouched, as these might be part of a transit. We remove long-term trends with a running-window median filter:<jupyter_code>from matplotlib import rcParams; rcParams["figure.dpi"] = 150
EPIC_id = 201367065        
star = everest.Everest(EPIC_id)
t = numpy.delete(star.time, star.badmask)
y = numpy.delete(star.fcor, star.badmask)
t = numpy.array(t[~numpy.isnan(y)], dtype='float32')
y = numpy.array(y[~numpy.isnan(y)], dtype='float32')

trend = scipy.signal.medfilt(y, 25)
y_filt = y /trend
y_filt = sigma_clip(y_filt, sigma_upper=2, sigma_lower=float('inf'))

fig, axes = plt.subplots(2, 1, sharex=True, figsize=(6, 6))
ax = axes[0]
ax.plot(t, y, "k")
ax.plot(t, trend)
ax.set_ylabel("Flux (electrons per sec)")
ax = axes[1]
ax.plot(t, y_filt, "k")
ax.set_xlim(t.min(), t.max())
ax.set_xlabel("Time (days)")
ax.set_ylabel("Normalized flux");<jupyter_output>INFO  [everest.user.DownloadFile()]: Found cached file.
INFO  [everest.user.load_fits()]: Loading FITS file for 201367065.<jupyter_text>This is a curve which looks suitable for a transit search. The regularly spaced signals look promising. In this tutorial, we do not use stellar priors, but instead use the TLS defaults for simplicity. However, as this system suffers from aliasing, we use a high `oversampling_factor` and a tight `duration_grid_step`.<jupyter_code>from transitleastsquares import transitleastsquares
model = transitleastsquares(t, y_filt)
results = model.power(oversampling_factor=5, duration_grid_step=1.02)<jupyter_output>0%|          | 0/14174 periods | 00:00<?<jupyter_text>After the search, we plot the power spectrum and the integer (sub)harmonics:<jupyter_code>plt.figure()
ax = plt.gca()
ax.axvline(results.period, alpha=0.4, lw=3)
plt.xlim(numpy.min(results.periods), numpy.max(results.periods))
for n in range(2, 10):
    ax.axvline(n*results.period, alpha=0.4, lw=1, linestyle="dashed")
    ax.axvline(results.period / n, alpha=0.4, lw=1, linestyle="dashed")
plt.ylabel(r'SDE')
plt.xlabel('Period (days)')
plt.plot(results.periods, results.power, color='black', lw=0.5)
plt.xlim(0, max(results.periods));<jupyter_output><empty_output><jupyter_text>There is a major peak at a period of ~10 days. Let's inspect the statistics:<jupyter_code>print('Period', format(results.period, '.5f'), 'd')
print(len(results.transit_times), 'transit times in time series:', \
        ['{0:0.5f}'.format(i) for i in results.transit_times])
print('Transit depth', format(results.depth, '.5f'))
print('Transit duration (days)', format(results.duration, '.5f'))<jupyter_output>Period 10.05556 d
8 transit times in time series: ['1980.41279', '1990.46834', '2000.52390', '2010.57945', '2020.63501', '2030.69057', '2040.74612', '2050.80168']
Transit depth 0.99867
Transit duration (days) 0.11015<jupyter_text>To get an impression of what this signal looks like, we can create a phase-folded plot:<jupyter_code>plt.figure()
plt.plot(
    results.model_folded_phase,
    results.model_folded_model,
    color='red')
plt.scatter(
    results.folded_phase,
    results.folded_y,
    color='blue',
    s=10,
    alpha=0.5,
    zorder=2)
plt.xlim(0.49, 0.51)
plt.xlabel('Phase')
plt.ylabel('Relative flux');<jupyter_output><empty_output><jupyter_text>There are some high points during the ingress phase, which we attribute to aliasing. Other than that, the signal has a limb-darkened transit shape. The next plot shows the full light curve, together with the best model (in red). The in-transit data points are also printed in red.<jupyter_code>from transitleastsquares import transit_mask
plt.figure()
in_transit = transit_mask(
    t,
    results.period,
    results.duration,
    results.T0)
plt.scatter(
    t[in_transit],
    y_filt[in_transit],
    color='red',
    s=2,
    zorder=0)
plt.scatter(
    t[~in_transit],
    y_filt[~in_transit],
    color='blue',
    alpha=0.5,
    s=2,
    zorder=0)
plt.plot(
    results.model_lightcurve_time,
    results.model_lightcurve_model, alpha=0.5, color='red', zorder=1)
plt.xlim(min(t), max(t))
plt.ylim(0.9985, 1.0003)
plt.xlabel('Time (days)')
plt.ylabel('Relative flux');<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12591911764705882, 'fraction_numerical': 0.04338235294117647, 'mean_word_length': 5.477380952380952, 'pattern_counts': {'":': 0, '<': 28, '<?xml version=': 0, '>': 27, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '90118', 'license': 'MIT', 'n_tokens_mistral': 2070, 'n_tokens_neox': 1897, 'n_words': 512, 'path': "tutorials/02 Starter's guide - K2-3b, a red dwarf with a planet.ipynb", 'repo_name': 'lgbouma/tls', 'size': '5452'}