The sherlock pipeline: new exoplanet candidates in the WASP-16, HAT-P-27, HAT-P-26, and TOI-2411 systems

ABSTRACT

The launches of NASA Kepler and Transiting Exoplanet Survey Satellite (TESS) missions have significantly enhanced the interest in the exoplanet field during the last 15 yr, providing a vast amount of public data that are being exploited by the community thanks to the continuous development of new analysis tools. However, using these tools is not straightforward, and users must dive into different codes, input–output formats, and methodologies, hindering an efficient and robust exploration of the available data. We present the sherlock pipeline, an end-to-end public software that allows the users to easily explore observations from space-based missions such as TESS or Kepler to recover known planets and candidates issued by the official pipelines and search for new planetary candidates that remained unnoticed. The pipeline incorporates all the steps to search for transit-like features, vet potential candidates, provide statistical validation, conduct a Bayesian fitting, and compute observational windows from ground-based observatories. Its performance is tested against a catalogue of known and confirmed planets from the TESS mission, trying to recover the official TESS Objects of Interest (TOIs), explore the existence of companions that have been missed, and release them as new planetary candidates. sherlock demonstrated an excellent performance, recovering 98 per cent of the TOIs and confirmed planets in our test sample and finding new candidates. Specifically, we release four new planetary candidates around the systems WASP-16 (with P |$\sim$| 10.46 d and R |$\sim$| 2.20 |$\mathrm{ R}_{\oplus }$|⁠), HAT-P-27 (with P |$\sim$| 1.20 d and R |$\sim$| 4.33 |$\mathrm{ R}_{\oplus }$|⁠), HAT-P-26 (with P |$\sim$| 6.59 d and R |$\sim$| 1.97 |$\mathrm{ R}_{\oplus }$|⁠), and TOI-2411 (with P |$\sim$| 18.75 d and R |$\sim$| 2.88 |$\mathrm{ R}_{\oplus }$|⁠).

1 INTRODUCTION

The search for exoplanets has experienced one the most significant developments in modern astronomy, where space missions such as Kepler (Borucki et al. 2010) and Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2016) have emerged as beacons of discovery confirming thousands of new worlds and providing with thousands of planetary candidates waiting to be confirmed. In addition, these missions’ public nature has produced hundreds of thousands of light curves, allowing the community to conduct alternative data exploration. To explore and exploit this vast amount of data, many software packages emerged; some designed for the analysis of light curves, including data loading and visualization tools such as lightkurve (Lightkurve Collaboration 2018) and eleanor (Feinstein et al. 2019), some others to conduct data detrending using a variety of functions and algorithms to remove instrumental systematics that produce photometric variability such as everest (Luger et al. 2016, 2018), k2ssf (Vanderburg & Johnson 2014), and wötan (Hippke et al. 2019), others to conduct a transit search such as the Box Least Squares (bls; Kovács, Zucker & Mazeh 2002), the Quasiperiodic Automated Transit Search (qats; Carter & Agol 2013), and the Transit Least Squares (tls; Hippke & Heller 2019), others to perform results vetting such as vetter (Zink et al. 2020) and robovetter (Coughlin et al. 2016), others dedicated to validate planetary candidates by determining whether a transit-like signal is statistically likely to be caused by a planet or a background eclipsing binary such as vespa (Morton 2012), pastis (Díaz et al. 2014), and triceratops (Giacalone et al. 2021), others to conduct a Bayesian fitting of the data to refine the planet parameters such as allesfitter (Günther & Daylan 2021), juliet (Espinoza, Kossakowski & Brahm 2019), and exoplanet (Foreman-Mackey et al. 2021), and others to compute observational windows such as tapir (Jensen 2013), to name a few.

While most of these packages are open source, their different methodologies and input–output formats hinder their correct usage and implementation, relegating them to independent uses in most cases instead of an efficient workflow where all the procedures are connected. The lack of connection between these packages drives the community to rely mostly on the alerts issued by the official pipelines such as Science Processing Operations Center (spoc; Jenkins et al. 2016) and the Quick Look Pipeline (qlp; Huang et al. 2020), i.e. the so-called Kepler Objects of Interest (KOIs) and TESS Objects of Interest (TOIs). However, many planets remain unnoticed due to the detection thresholds of these pipelines, the lack of data exploration, or the poor photometric quality. A proof of that is the number of alternative catalogues providing hundreds of new planetary candidates using Kepler data (see e.g. Ofir & Dreizler 2013; Coughlin et al. 2016), K2 data (see e.g. Barros, Demangeon & Deleuil 2016; Mayo et al. 2018; Kruse et al. 2019), and TESS data (Cloutier 2019; Montalto et al. 2020; Feliz et al. 2021; Montalto 2023), highlighting the importance of the parallel analyses of the data provided by the space-based missions by multiple teams. A brief review of some of these alternative catalogues can be found in Melton et al. (2024), where the authors presented their own catalogue using TESS data from Cycle 1. None the less, the pipelines producing these extra candidates are mostly non-public; only a few are open source, such as dave (Kostov et al. 2019a), which is designed to search and vet planetary candidates using K2 data. Then, the community is restricted to solely focus on those alerts, impeding the freedom to concentrate on specific scientific cases designed by independent groups. This paper presents the sherlock (Searching for Hints of Exoplanets fRom Light curves Of spaCe-based seeKers) pipeline. This public python package allows the users to conduct their own planetary searches, examining robustly the data from space-based missions. The tool was first introduced at its initial development stages with its core search functionality in Pozuelos et al. (2020) and Demory et al. (2020), and has been used in a number of studies. The iterative development has led sherlock to become an end-to-end pipeline where a special effort has been made to reduce user intervention and sort all the steps in an efficient workflow. With only a few command lines, the users can explore the existing data, conduct the transit search, vet the most promising signals, compute a statistical validation, conduct a Bayesian fitting to refine planetary parameters and ephemerides, and trigger a follow-up camping using ground-based facilities. With the implementation of its six modules, sherlock can be used in any imaginable exoplanetary survey that relies on space-based data. For example, sherlock is being used in the SPECULOOS (Search for habitable Planets EClipsing ULtra-cOOl Stars) project (Gillon 2018; Sebastian et al. 2021), which aims at finding transiting planets orbiting nearby ultracool dwarfs using a network of six 1-m ground-based telescopes located in different observatories such as ESO Paranal Observatory (Atacama Desert, Chile; Burdanov et al. 2018; Delrez et al. 2018), Teide Observatory (Canary Islands, Spain; Burdanov et al. 2022), and National Astronomical Observatory of Mexico (San Pedro Mártir, Mexico; Demory et al. 2020). Due to the small sizes of these stars, transiting planets orbiting them are prime targets for atmospheric characterization using, for example, the JWST (Gillon et al. 2017, 2020). In this context, sherlock is exploring those stars in the SPECULOOS target list that have been observed by TESS, searching for transit-like features that might be attributable to planets. While the SPECULOOS telescopes have better photometric precision than TESS for these stars (see e.g. Delrez et al. 2022), a robust inspection of the TESS data can help to speed up the process of planetary searching, saving telescope time. In addition, sherlock is the main pipeline of the FATE (forecasting optical turbulence to push the Very Large Telescope to its full potential) project (Van Grootel et al. 2021), which aims at finding the first observational evidence of transiting planets orbiting hot subdwarf stars, which are direct post-asymptotic giant branch objects that have lost most of their envelopes. In particular, the project objective is to shed light on understanding what happens to planetary systems once their host stars leave the main sequence and engulf close-in planets, presumably perturbing the dynamical architecture of the system. In this context, sherlock systematically searches for transiting-like features in the Kepler, K2, and TESS data (Thuillier et al. 2022).

Hence, the idea behind sherlock is not producing a vast amount of new planetary candidates as other existing pipelines do, but empowering independent researchers to carry on their own designed surveys. The paper is organized in the following way. Section 2 presents the key concepts and the workflow behind sherlock. In Section 3, we use sherlock to recover the planets of some known TESS systems, evaluate the pipeline performance while also presenting the detection of four new candidates within those systems. In Section 4, we present a new project in which sherlock will be used to search for extra planets that remain unnoticed. Finally, Section 5 presents our discussions, conclusions, and prospects.

2 THE sherlock PIPELINE

sherlock is an end-to-end pipeline that allows the users to explore the data from space-based missions to search for planetary candidates. It can be used to recover alerted candidates by the automatic pipelines such as spoc (Jenkins et al. 2016) and qlp (Huang et al. 2020), the so-called Kepler KOIs and TESS TOIs, and to search for candidates that remain unnoticed due to detection thresholds, lack of data exploration, or poor photometric quality. To this end, sherlock has seven different modules to (1) acquire and prepare the light curves from their repositories, (2) search for planetary candidates, (3) vet the interesting signals, (4) perform a statistical validation, (5) carry out planetary system stability tests, (6) model the signals to refine their ephemeris, and (7) compute the observational windows from ground-based observatories to trigger a follow-up campaign. To execute all these modules, the user only needs to fill in an initial yaml file with some basic information such as the star ID, the cadence to be used, etc., and use sequentially one command line to pass from one step to the next. Alternatively, the user may provide sherlock with the light curve in a csv file, where the time, the normalized flux, and the flux error need to be given. In the following subsections, we describe in detail each of these steps. The complete list of parameters that can be used in the yaml file, along with a description of each one, can be found in GitHub’s sherlock repository.¹ Keeping in mind that future implementations and releases might require new parameters; hence, we suggest regularly checking the list.

2.1 Light-curve acquisition and preparation

sherlock makes use of the lightkurve package (Lightkurve Collaboration 2018) to obtain the data of the space-based missions Kepler, K2, and TESS, from NASA Mikulski Archive for Space Telescope (MAST). In most cases, the user only needs to provide the star’s Kepler Input Catalog (KIC), K2 Ecliptic Plane Catalog (EPIC), or TESS Input Catalog (TIC) identification number and the cadence to be used, i.e. 20, 120, 200, 600, or 1800 s for TESS and 120 or 1800 s for Kepler/K2. In the case of the 600 and 1800 s for TESS, i.e. the full-frame images, the data can also be retrieved via the eleanor package (Feinstein et al. 2019). By default, the Pre-search Data Conditioning Simple APerture (PDC-SAP) is used, which has some level of processing with long-term trends removed using co-trending basis vectors. However, the PDC-SAP light curves are often affected by some other trends or systematics that need to be corrected before starting any planetary search. For example, some light curves present noisy regions with a high level of data dispersion, some others are affected by strong variability due to a fast rotation or stellar pulsations, and some others present short-term trends induced by the instrumentation. To efficiently deal with all these scenarios, sherlock accounts with different semi-automatic tools.

2.1.1 Masking noisy regions

The light curves often show regions with high data dispersion where the photometric quality is significantly reduced, strongly affecting the performance of the search algorithm. These noisy regions are sometimes related to the spacecraft’s motion, light scattering produced by the Earth, background moving objects, or other unidentified sources. Hence, to clean the light curve of such nuisance data, sherlock can automatically compute the root-mean-square (rms) of the data set by 4-h binning. Afterward, all the light curve time gaps longer than 3 d are identified to split the data into segments. The mean rms of each data segment is then computed, and each data point that falls into an rms bin above 2|$\times$| the mean rms of its respective data segment is removed (see Fig. 1). The threshold of 2 is applied by default, but the user can adjust it if needed in the yaml file. As shown in Fig. 1, sherlock also plots three rms threshold suggestions (1.25, 1.5, and 1.75|$\times$| the rms) alongside the user selection to provide insight into the typical values used for this mask. In this figure, some potential transit signals from TOI-144 appear to be masked when low rms threshold values are used, whereas a threshold of |$1.75\times$| rms seems to be optimal for masking the entire noisy region without affecting the transits. In cases where the automatic high-rms-masking algorithm is not precise enough due to the complexity of the noisy features, the user can freely enter a set of time ranges to be masked.

2.1.2 Stationary stellar pulsations extraction

For quasi-stationary pulsations, we identified many cases (thanks to injection and recovery tests) where short-period transits might be confused with them, and then our targets would be removed using the automatic rms-based noise removal algorithm. For such targets where these types of pulsations are identified, we developed a custom variation of the Barceló Forteza et al. (2015) steps to reduce their effect on the light curve:

• Clear transits masking: An initial batch of searches with tls is done to mask transits with low residuals, which are typically appreciable in the light-curve periodogram and might significantly alter the results of the next steps. That would happen because the subharmonics of the main transit period get narrower, and their shape in the periodogram might get close to the ones offered by the stationary pulsations that we would like to reduce.

• Computation of highest peak signal-to-noise ratio (S/N): Once the highest peak is identified in the periodogram, its S/N is computed by dividing its power by the median of powers of its surrounding periods within a window of 4 h. We used the median instead of the mean as it is typically a more robust statistical estimator of baseline values, which is what we try to extract here.

• Model the stationary signal: In case the S/N is higher than a given threshold (defaults to 4), the signal is modelled with a non-linear least-squares algorithm with an initial estimation of the amplitude and frequency through a fit of a cosine wave.

• Check model validity: If the model shows an amplitude uncertainty lower than 20 per cent and a phase uncertainty lower than 0.2 rad, it will be extracted from the current light-curve flux. If the resulting signal shows lower rms, it will be accepted, and the resulting flux will be stored. We return to step (ii) and apply the algorithm iteratively whenever a model is accepted.

To test the capabilities of the Stationary Stellar Pulsations Extraction (sspe) algorithm, we selected the TIC 169285097, catalogued as a hot subdwarf pulsating in p and g modes by Sahoo et al. (2020). These targets are complex hosts to search for small planets producing transits. Indeed, the strong rms induced by the several pulsation modes would make such planets undetectable, increasing their residuals and reducing their S/N. Also, the typical transit detection algorithms often fail when searching for transits around these targets because they usually match sine-shaped pulsations to transit models, and their noise can be strongly correlated. As we are trying to diminish these effects, we tried our algorithm on the selected target to inspect the results that could be achieved on one of the available sectors from the TESS mission. We display the resulting flux in Fig. 2, which was obtained after subtracting 44 pulsation frequencies. The most powerful frequencies are shown in Fig. A1, where their corresponding models are overplotted to the phase-folded light curves. To properly investigate the capabilities of the algorithm in terms of transiting exoplanet signals recovery and to ensure that the procedure does not remove transit features, we tried our implementation on a public tool that we also developed in-house and published in a public repository under a free usage license: The matrix toolkit² built on top of the same code base as sherlock, it contains the same algorithmic capabilities as the former and is conceived as an injection and recovery engine that includes a stage to generate a grid of periods, radii, and epochs to be injected in the original target’s light curve. Those synthetic curves are later used in a recovery stage to search for the known injected transits. matrix will store a csv file to keep track of which models were found and which were not, generating at the end a heat map showing the matches for each period, radius, and epoch from the grid.

To assess properly the improvement brought by our sspe algorithm, we first tried the injection and recovery scenario for a simple search without using it. Later, we tried the same scenario with the inclusion of the algorithm. The results can be appreciated in Fig. 3, where we can see that the first algorithm could add for this target an increment in the detectability of 1 R|$_{\oplus }$| for the best cases. Altogether this study case leads us to an outstanding improvement in our detection rates for our current transit surveys like the one started in Van Grootel et al. (2021) and upcoming ones in the near future.

An improved version of this iterative algorithm would include the entire set of already selected signals on each iteration to improve the model-fitted parameters as it is done with felix in Charpinet et al. (2010) and Zong et al. (2016). However, this is not included in the sherlocksspe procedure yet.

2.1.3 Removing short trends

sherlock uses a multidetrend approach to deal with short trends caused, for example, by instrumental drifts and intrinsic stellar variability. This strategy is implemented via the wötan package (Hippke et al. 2019). By default, sherlock detrends the flux using the biweight method. This choice is based on its good speed–performance ratio compared with other methods, such as the Gaussian process (GP), which is also available. The biweight method is a time-windowed slider, where short windows can efficiently remove short trends, but there is an associated risk of removing any actual transit signal. Hence, to prevent this issue, the pipeline explores, by default, 10 cases with different window sizes. The minimum value of the window size is obtained by computing the transit duration (T|$_{14}$|⁠) of a hypothetical Earth-size planet with an orbital period of 10 d. A minimum window size of 3 |$\times \, T_{14}$| is set to protect such a transit. On the other hand, the maximum value explored is 20 |$\times \, T_{14}$|⁠, which is enough to remove the variability in most of the cases. However, if preferred, the user may set the number of detrends to be applied and the minimum and maximum values of the window size in the yaml file. Hence, from the original PDC-SAP flux, sherlock generates a new set of detrended light curves, each differing only in the window size used for detrending.

The transit search, described later in Section 2.2, is then conducted over all the detrended light curves jointly with the original PDC-SAP light curve. If a signal is detected in only one (or very few) of the light curves, it is considered ‘detrend dependent’, suggesting it may be an artefact produced by the detrending algorithm. Conversely, if a signal is consistently found in multiple (or all) light curves, it is considered ‘detrend independent’, thus increasing its credibility. By comparing the S/N and signal detection efficiency (SDE) of each detection across these light curves, we can determine the most appropriate window size for detrending. This method enhances the reliability of our planet detections by minimizing the risk of false positives due to the detrending process, ensuring that any promising signals found are not artefacts of a specific detrending configuration.

2.1.4 Correcting strong variability

Many stars suffer significant changes in brightness, which might be caused, for example, by stellar spots or intrinsic stellar variability. In such cases, the PDC-SAP light curves show intense sinusoidal trends of high amplitudes. To deal with these scenarios, sherlock computes the Lomb–Scargle periodogram (LSP) and extracts the value of the strongest peak using the lightkurve package (Lightkurve Collaboration 2018). By default, a factor of this value is used to remove the sinusoidal trend by fitting a sum of sines and cosines employing the cosine function provided by the wötan package (Hippke et al. 2019). Our calibration experiments showed that the most appropriate value for such a factor is 0.25. However, the user may adjust it to any other value if wanted in the yaml file. There are some limitations in this procedure; for example, the presence of a hot Jupiter may induce a strong peak in the LSP that can be confused by a sinusoidal trend. Hence, a visual inspection is always desirable before applying this correction.

2.2 Search for planetary candidates

By default, sherlock performs the transit search by using a modified version of the tls package (Hippke & Heller 2019). Since tls uses an analytical transit model taking into account the stellar parameters, it offers a better S/N and SDE than classical algorithms based on bls (Kovács et al. 2002). This stellar-model optimization allows the detection of shallow periodic transits, which are more challenging to find. However, if desired, the user may set the bls algorithm in the yaml file.

The procedure followed by sherlock is to search for transits in all the detrended light curves jointly with the original one (PDC-SAP flux). Then, the results obtained for each individual light curve are stored; i.e. the orbital period (P) and its error in days, mid-transit time for the first identified transit (T|$_{0}$|⁠) in TBJD (TESS Barycentric Julian Date) or BJD (Barycentric Julian Date), the transit depth in ppt (parts per thousand), the transit duration (D) in minutes, the number of transits, the S/N, the SDE, and the false-alarm probability (FAP). Once all the light curves are analysed, the results are combined to select the best transiting signal. To this end, sherlock provides first a border-correction algorithm. This algorithm computes a synthetic SDE value for each signal, penalizing the transits that are found in the border of the data set, which often corresponds to regions where there are more uncorrected systematics, as follows:

where |$N_{\rm tb}$| is the number of transits close to borders, and |$N_{\rm t}$| is the total number of transits. Secondly, the signal selection is performed by a quorum algorithm, where signals found in more light curves are favoured in a voting process over those found on less light curves. Once a signal is selected as the best candidate, sherlock masks it in the PDC-SAP light curve using its P, T|$_{0}$|⁠, and 3 |$\times \, D$|⁠, and starts a new searching run. The user can set how many times this process is repeated in the yaml file; however, our experience points out that results found beyond five runs are less reliable. This issue is coming due to the accumulated gaps in a given light curve after many mask-and-run iterations. Hence, by default, this process is repeated five times. In addition, sherlock establishes minimum thresholds for the S/N and the SDE of a selected signal to pass from one run to the next. By default, both S/N and SDE are set to five, but the user may change these minimum thresholds in the yaml file.

Besides the quorum algorithm, sherlock has two other more simple alternative selection methods. These two methods are based on the selection of the signal with larger S/N and the larger |${\rm S/N_{synthetic}}$|⁠, respectively. Moreover, sherlock provides a flexible module to use any other user-defined selection algorithm.

Once the searching stage is finished, sherlock produces as many folders as the number of runs executed, containing plots showing for each light curve (original PDC-SAP flux and detrended light curves) the best signal found. Moreover, sherlock provides two log files: (1) report.log, where the user can find all the information about the searching process, and (2) candidates.log, where is listed a summary of the main signals found. Moreover, the user will have access to extra folders named /fov, /detrends, /rms, and /sa. In /fov (from field of view) is stored the target pixel file (TPF) with the aperture used to extract the photometry overplotted along with the nearby stars from Gaia Data Release 2 catalogue, obtained via the tpfplotter package (Aller et al. 2020), and the TPF with the light curves corresponding to each individual pixel using the lightkurve package (Lightkurve Collaboration 2018). In /detrends, the user can visually inspect all the detrends applied to the original PDC-SAP light curve as described in Section 2.1.3. In /rms, a collection of plots shows the regions of the light curve that are masked due to their high rms values, if any, as explained in Section 2.1.1. In case the sspe module is activated (Section 2.1.2), the extracted sinusoidal trends are plotted in the /sa directory. Finally, the LSP and the value of the strongest peak are displayed as described in Section 2.1.4. Combining these outputs, sherlock gives the user the information needed to evaluate the results’ reliability and the signals listed as promising.

2.2.1 tls metrics fine-tuning

The original tls algorithm presented by Hippke & Heller (2019) reported some inconsistencies with the SDE and S/N in some scenarios. Hippke & Heller (2019) compute the transit duration by using a custom factor called ‘curve stretch’, supposedly taking into account the time gaps within the curve. Despite, this does not seem to be working properly and the error created in the transit duration propagates to the S/N calculation, which is based on this value. The S/N is especially hindered in curves with long time gaps, which are a typical case in the TESS mission. That is why we decided to rewrite the transit duration calculation by only taking the best model fit provided by tls, measuring the period fraction where it falls below 1. This fix to the algorithm seemed to repair the S/N computations and the package began to return very stable S/N and some compared outputs can be explored in Table 1.

In addition, the SDE is a measurement that is expected to remain at scale among use cases and, hence, be comparable. We remind the reader that SDE is a metric indicating the statistical significance of a given period compared to the other periods tested in the transit search. In particular, an SDE value of x for any given period P means that the statistical significance of this period is |$x\times \sigma$| compared to the mean significance of all other periods. [We refer the reader to Hippke & Heller (2019) and references therein for further details.] Then, this metric allows us to evaluate and compare the prominence of detected signals across different periods. However, the computation of SDE by tls implies the calculation of the residuals, which in some edge cases are set to a default baseline that can lead to an overestimation of the SDE. Specifically, when many residuals are set to the baseline value, the few significant spikes can cause abnormally high SDE values, for example, when a given light curve contains many gaps. To mitigate this, we rescale each final SDE value based on the proportion of local residuals that are above the baseline. This adjustment ensures that the SDE values more accurately reflect the true statistical significance of the detected signals. As an example, Fig. 4 shows the difference between an original SDE spectrum with many baseline values and its rescaled counterpart. The last problem that we observed affecting the SDE is the difficulty of establishing fixed thresholds based on FAPs. The original tls implementation calibrated a bijective relationship between SDE and FAP values. However, this calibration was carried out only for scenarios with time spans up to 28 d (the TESS typical length for one sector data) and does not provide good FAP values for longer data contexts, where SDEs tend to be higher on average even for false positives.

2.2.2 The Savitzky–Golay digital filter

sherlock allows applying the Savitzky–Golay digital filter (Savitzky & Golay 1964) to the original PDC-SAP flux before any detrend and starting the transit search. This filter is implemented by using the scipy package (Jones et al. 2001), which performs a data convolution by fitting subsets of adjacent data points with a low-degree polynomial using the linear least-square method. Hence, the filter correlates nearby points in the data reducing their local white noise without distorting the signal trend, leading to half the initial rms value in the light curves. The Savitzky–Golay digital filter is handy when dealing with threshold-crossing events (TCEs) with low S/N. While the usage of this filter has the potential to find signals buried in the noise, it also has the associated risk of obtaining more false-positive detections that the user needs to vet carefully. In Table 2, we show the S/N and SDE values for TOI-270 (Günther et al. 2019) and TOI-736 (Crossfield et al. 2019; Peterson et al. 2023) using a pure-tls search and applying first the Savitzky–Golay digital filter. As can be noticed, the S/N reported by tls is improved by a factor of 2, and the SDE is sometimes improved and sometimes worsened without a noticeable impact on the search. The S/N and SDE reported by sherlock are extracted directly from the tls algorithm, calculated with a naive estimation that cannot be compared with values from the spoc pipeline. Therefore, these values are only used internally to compare the goodness of selected candidates and choose the clearest one. A recent application of this filter was used to recover the low S/N candidate TOI-4306.01 by Delrez et al. (2022), for which the pure-tls search found the candidate with S/N below the typical detection threshold of S/N |$\ge$| 5.

2.3 Semi-automatic vetting

As described in Section 2.2, once the searching step is done, sherlock provides several outputs that allow the user to perform a preliminary vetting. However, many other effects might be causing the detected signals that need careful inspection.

The finished Kepler mission (Borucki et al. 2003) provided more than 30 000 TCEs (Thompson et al. 2018), which are candidate signals that were found by the spoc main search pipeline. However, different sources could generate these signals or even be instrumental artefacts induced into a target’s light curve. To rule out these scenarios, the spoc implemented the data validation reports (DVRs; Twicken et al. 2018), which are documents showing different metrics to qualify or discard an analysed candidate. These scenarios are (1) transit shape model fit, (2) odd–even transits checks, (3) centroids shifts, (4) optical ghost effects, (5) transit source offsets, (6) rolling band contamination, (7) secondary depth significance, and some additional numerical tests. Based on these tests sherlock implements a vetting procedure that is explained in the next subsections.

While these tests may offer a good vetting for promising signals, there are some limitations, mainly when dealing with shallow transits with very low S/N, where the results of the vetting might be inconclusive. sherlock uses the vetting stage by importing the watson package, which we also developed and published as an independent open-source tool,³ offered in the official python repositories⁴ with its own documentation.⁵

2.3.1 Single transit inspection

For each vetted candidate, sherlock computes the light curve with the official pipeline aperture and a smaller one and displays them together for each single transit. The official pipeline aperture’s centroids, motion, and background flux data are also displayed. In addition, the TPF is folded and plotted focused on each single transit. No metrics are computed from this data, but the visual inspection of these plots is very useful for finding differences or problematic features in single transits.

2.3.2 Cadence-based curves signal reliability

By using the period, epoch, and duration of a promising signal found in the transit search module, sherlock computes the S/N for a box-shaped transit for each available cadence. If the transit is not significant enough (S/N |$\ge$| 5) in any of the cadences or the S/N is substantially different from the nominal cadence, sherlock will warn about it in the report.

2.3.3 Odd–even and secondary events

The transit-like signals might be originated from astrophysical sources that might not be planets. This is the case of eclipsing binaries or nearby eclipsing binaries, which typically show secondary eclipses with different transit signatures to the primary one. sherlock computes the S/N for odd and even transits depths and also for secondary eclipses. A red flag is reported when secondary eclipses are detected with S/N |$\ge$| 3 or in case either odd or even transits show S/N |$\le$| 5.

2.3.4 Transit source offset: images

The difference imaging method proposed in Bryson et al. (2017) was the selected approach by the spoc to accept or reject signals based on the offset found on the computed transit source from pixel images. The method consists of computing the centre of mass of the subtraction of the flux measured within the TPF for in-transit and out-of-transit per pixel data. This way the most probable source of the transit events can be estimated. We have implemented a simplified version of this test in sherlock. In addition, we have also included a second method to compute this metric by running a per pixel bls search with a fixed period on the folded TPFs. We then compute the so-called centre of mass from each of the per pixel bls residuals and, finally, average the output of both methods to get the final transit source position.

2.3.5 Transit source offset: centroids shift

An alternative method to estimate scenarios where the transit signal might have originated outside the target is the centroids diagnosis, also presented by Bryson et al. (2017). In our case, we subtract the mission-computed motion values (polynomials relating the sky and the focal plane for each cadence) from the mission-computed centroid data (target position on the TPF given by the pixel response function obtained from the out-of-transit control image; Bryson et al. 2013). This method measures the deviation of the target position variations on the TPF against the focal plane alignment variations to the sky. If these values report an S/N > 3, sherlock will be placing a red flag for the analysed signals.

2.3.6 Optical ghost diagnostic

Sometimes the transit-like signals are caused by artificial sources. That is the case of optical ghosts, which can be created by light reflections from high-intensity sources between the optical system and the CCD as shown by Caldwell et al. (2010). To find false positives caused by optical ghosts, the Kepler DVR performs a test that measures the correlation between core (mean flux from all pixels within the target aperture) and halo (mean flux from all the pixels around the target aperture) fluxes and the transit signal. When the signal is affected by a ghost, the halo flux will be more correlated than the core one.

sherlock implements a similar test. In our case, we compute and subtract the core and halo fluxes. In case the result still shows a transit signal with an S/N > 3, sherlock will reject the candidate.

2.4 System stability

In some cases, several promising signals are selected during the transit search, but their orbital periods are very close to each other. In other cases, the selected signal is very close to a known confirmed planet. The system’s architecture might be unstable in these scenarios, revealing that the candidates found are likely false positives.

To shed light on these scenarios, sherlock incorporates a module to assess the system dynamics stability. This analysis can be conducted at two points: during the initial vetting of a signal to assess whether there are possible system configurations for our candidates and as part of the validation scenario to get stability indicators for realistic or known configurations. That is, sherlock can run a stability analysis in two modes. The first consists of giving user inputs, which create a linear grid of possible orbital parameters for each candidate. The second one loads the results from the Bayesian fit (see Section 2.6) for the selected signal, allowing the addition of some constraints to the most likely system configuration.

To compute the system stability, sherlock evaluates the fast chaos index Mean Exponential Growth factor for Nearby Orbits (MEGNO; Cincotta & Simó 1999, 2000; Cincotta, Giordano & Simó 2003), which is widely used in this context (see e.g. Jenkins et al. 2019; Delrez et al. 2022; Pozuelos et al. 2023). sherlock uses the N-body numerical package rebound (Rein & Liu 2012), which provides a MEGNO implementation by employing the Wisdom–Holman whfast code (Rein & Tamayo 2015). The MEGNO mean value over time, denoted as |$\langle Y(t) \rangle$|⁠, enhances stochastic fluctuations. This characteristic facilitates the differentiation of trajectories into quasi-periodic, that is to say, stable orbits, where |$\langle Y(t\rightarrow \infty)\rangle$| approaches 2, and chaotic orbits, where |$\langle Y(t\rightarrow \infty)\rangle$| tends towards infinity for an infinitive integration time or diverging from 2 for finite integration times.

This highly flexible module allows users to investigate as many parameters as desired to test the system’s stability. Indeed, the user can explore the system’s stability by varying the stellar mass, the planetary period, eccentricity, mass, argument of periastron, and inclination. Since we are dealing with transiting planets, not always the planetary masses are known; then, for these cases, sherlock obtains the planetary mass relying on mass–radius relationships from Chen & Kipping (2017), which provides differentiation between small and large planets in a similar way than Weiss et al. (2013) and Bashi et al. (2017). In all the cases, for each parameter sherlock uses uniform distributions, where the user needs to provide the nominal value, the associate uncertainty, and the number of bins that want to explore within that range. The results are then stored in csv files that the user can employ to construct stability maps. It is important to notice that while sherlock offers the possibility of exploring a number of the parameters, not all of them equally affect the system stability. Indeed, as stated in Cranmer et al. (2021), the planetary semimajor axes (or the orbital periods), the masses, and the eccentricities are the parameters that most affect the stability of a given system. Out of these three parameters, usually are the masses and the eccentricities the most unconstrained and, hence, the ones used to explore system architecture (see e.g. Demory et al. 2020; Murgas et al. 2023).

To correctly sample the planetary positions over their orbits, the scenarios are integrated using a time-step of 5 per cent of the shortest orbital period in the system. Finally, running the integration for a time of |$10^5$|–|$10^6$| times the largest planetary orbital period in the system is recommended to start getting convergent and reliable results. In any case, this analysis does not intend to replace a more robust dynamical study of a given system but a preliminary view of the credibility of planetary candidates and the potential system architecture.

2.5 Statistical validation

Kepler/K2 and TESS found many candidates that, due to the faint nature of their host stars and/or the shallowness of their transits, cannot be confirmed from ground-based observatories as true planets. To mitigate this situation were developed the first statistical validation tools during the Kepler era (see e.g. Torres et al. 2011; Morton 2012; Morton et al. 2016), which allow us to calculate the false-positive probability (FPP) of a given detection to validate it statistically.

This procedure is nowadays widespread to support the planetary interpretation of a given signal even when the candidate can be reached from ground-based observatories (see e.g. Kostov et al. 2019b; Quinn et al. 2019). sherlock allows the user to perform a statistical validation automatically connecting the results with the triceratops⁶ package (Giacalone et al. 2021, 2022). Originally, triceratops was designed to vet and validate TESS candidates only; however, we have contributed to developing this open-source package to manage also Kepler/K2 candidates. In short, triceratops computes the probabilities of a wide range of transiting-producing astrophysical scenarios such as a transiting planet and eclipsing binary (in many different configurations) using the primary transit of the candidate, the prior knowledge of the host and nearby stars, and the current understanding of planet occurrence and stellar multiplicity. Hence, triceratops computes the location of a given candidate in the FPP–nearby FPP (NFPP) plane, where the candidate will be classified as a validated planet if FPP |$\lt $| 0.015 and NFPP |$\lt$| 10³, as likely planet if FPP |$\lt$| 0.5 and NFPP |$\lt \ 10^{-3}$|⁠, and as likely false positive if NFPP |$\gt \ 10^{-1}$|⁠. Due to the probabilistic nature of the procedure to calculate the values of FPP and NFPP by triceratops, sherlock performs five computations and then averages the results to provide the final values. In addition, sherlock uses the results from Lissauer et al. (2012), which estimated the fraction of planet candidates that are false positives when in systems with multiple planet candidates based on Kepler data. Then, this prior probability is multiplied by the FPP and NFPP previously computed (see e.g. Demory et al. 2020).

2.6 Bayesian model to refine the ephemeris

Once a candidate has been vetted and validated, it can be desirable to trigger a follow-up campaign to firmly confirm the planet in the target star using ground-based observatories. However, the results from the searching stage via tls are not optimal for scheduling an observation due to the large uncertainty associated with the found orbital period (⁠|$\Delta P$|⁠), whose value directly depends on the period grid density. Indeed, the transit time (T) and its uncertainty (⁠|$\Delta T$|⁠) of future observational windows are given by (see e.g. Mallonn et al. 2019)

with T|$_{0}$| and P previously defined in Section 2.2, and n the number of cycles since T|$_{0}$|⁠. Hence, error propagation is strongly dominated by the term |$n \, \Delta P$|⁠, which may produce |$\Delta T$| of the order of several hours if the detection is performed in Kepler/K2 data or early TESS sectors. This situation might hinder any ground-based confirmation.

In this context, sherlock allows the user to refine the ephemeris by performing a Bayesian fit using the allesfitter⁷ package (Günther & Daylan 2019). allesfitter is an open-source code that performs a Bayesian fit to infer the star and exoplanet properties given a set of photometric and/or radial velocity data. Then, sherlock automatically injects the results from the searching stage and the stellar parameters in allesfitter as priors and implements two steps fitting. (1) The candidate signals are masked, and the out-of-transit correlated noise estimation is performed using a GP with a 3/2 Matérn kernel with a Markov chain Monte Carlo (MCMC) fit on the PDC-SAP flux light curve, which can describe smooth long-term trends and short-term stochastic variations by its two hyperparameters: the amplitude scale (⁠|$\sigma$|⁠) and length scale (⁠|$\rho$|⁠). (2) The candidate signals priors are fit taking into account the previously derived GP parameters using normal priors, modelling the radius ratio of the planet to host, |$R_{\mathrm{ p}}/R_{\star }$|⁠, the sum of the stellar and companion radii divided by the semimajor axis, |$(R_{\star }+R_{\mathrm{ p}})/a_{\mathrm{ b}}$|⁠, the cosine of the orbital inclination, |$\cos i_{\mathrm{ p}}$|⁠, the epoch/transit mid-time in days, |$T_{0}$|⁠, the orbital period of the planet in days, |$P_{\mathrm{ p}}$|⁠, the transformed limb-darkening coefficients following a quadratic law, |$q_{1}$|–|$q_{2}$|⁠, and the natural logarithm of the white noise, log |$\sigma _{\mathrm{ w}}$|⁠. By default, the Bayesian model used by sherlock is the nested sampling, but the user might switch to MCMC if desired in the yaml file. All the posterior distributions and the derived parameters are stored in a new folder named /fit. At the end of this stage, the user obtains a robust inference of the planetary parameters. In particular, the |$\Delta P$| is reduced considerably, and the new values can be used to schedule ground-based observations with much more accurate predictions.

2.7 Ground-based observational windows

Once the planetary properties are optimal from Section 2.6, the user may safely plan a follow-up campaign from ground-based observatories to confirm the planetary nature of a given candidate. This step is particularly critical for results obtained using TESS data. Indeed, the pixel scale of TESS, of 21 arcsec, is remarkably larger than Kepler’s, of 4 arcsec, moreover the TESS’ point spread function (PSF) might be as large as 1 arcsec. Hence, due to dilution effects, the pixel scale, and the large PSF, the probability of contamination by nearby eclipsing binaries, whose deep eclipses might mimic shallow transits detected in the target star, is not null and, ideally, has to be checked from ground-based facilities. This ground-based confirmation is the primary goal of the seeing-limited photometry working group of TESS (see e.g. Kostov et al. 2019b). To schedule these observations, the community typically uses the TESS Transit Finder tool, which is a customized version of the tapir package (Jensen 2013).

In this context, sherlock automatically loads the Bayesian fit results and a list of observatories provided by the user from where the transit windows will be computed. The observational windows from each observatory will consider some user-defined constraints, such as the minimum altitude of the object, distance from the moon in relationship to its illumination phase, and the minimum transit coverage. Then, sherlock displays the observational windows in csv and pdf files, including the ingress, mid-transit, and egress times in ut for each event and plots with the airmass and altitudes of the host star during the entire event. This module was built by using the astropy package (Astropy Collaboration 2013, 2018) and was tested against the TESS Transit Finder tool to verify its good performance.

3 sherlock PERFORMANCE

Previously, we have outlined several projects leveraging the sherlock pipeline, such as SPECULOOS and FATE. Moreover, many other studies employed this pipeline to retrieve alerts from the TESS mission and known planets, as well as to identify new planetary candidates (see e.g. Wells et al. 2021; Schanche et al. 2022; Murgas et al. 2023; Dransfield et al. 2024; Dreizler et al. 2024; Gillon et al. 2024). In this section, our focus shifts to demonstrating the performance of the sherlock pipeline through a specifically designed experiment to provide a level of reliability.

In particular, we built a catalogue of 100 stars that fulfil the following conditions in the Exoplanet Follow-up Observing Program (ExoFOP) list:⁸ (1) the targets should, at least, contain one previously known planet or a TOI confirmed planet; (2) one of these planets should have orbital periods shorter than 10 d; (3) the maximum number of sectors should be two; and (4) the produced list with previous conditions is sorted by descending brightness and limited to the first 100 stars. Filters (1), (2), and (3) are applied to reduce the computational cost of this experiment. Appendix  B lists the catalogue obtained following these criteria providing the TIC-IDs. The main goal of this test is to run sherlock over this catalogue and record the successful recovery and non-recovery of planets and TOIs. Hence, to easily illustrate the performance, we introduce the recall parameter defined as |${\rm TP/(TP + FN)}$|⁠, where |$\mathrm{ TP}$| are the true positives, i.e. the number of correct planets or TOIs found, and |$\mathrm{ FN}$| are the false negatives, i.e. the number of known planets or TOIs that we do not recover. Then, the performance of the code is optimal when this parameter approaches 1, indicating that sherlock successfully identifies the expected planets without omissions.

From the 100 stars composing our catalogue, there were 104 TOIs that might correspond to previously confirmed planets, i.e. known planets for which TESS also provided an official TOI, or new planetary candidate alerts released by TESS that the scientific community has confirmed as planets. The number of multi-alert stars was much reduced because of the filters used to build the target list. Indeed, our filters favoured targets with large planets in short orbital periods, such as hot Jupiters (periods up to 10 d and radius from 6 |$\mathrm{ R}_{\oplus }$|⁠). Population studies suggested that nearby planetary companions rarely accompany these planets. However, thanks to the extreme photometric precision of space-based missions and dedicated radial velocity surveys, nearby planetary companions to hot Jupiters have emerged in recent years (Wu, Rice & Wang 2023). These recent discoveries suggest that the formation of hot Jupiters via quiescent mechanisms is more efficient than previously thought, when the main formation channel invoked was the dynamical high-eccentricity migration (Dawson & Johnson 2018). Hence, in this context, our target list not only serves to compute the recall factor and establish the goodness of the sherlock pipeline but also offers an interesting scientific case searching for low S/N signals that might hint at the existence of planetary companions to these hot Jupiters in our list, supporting this new emergence of planets.

Once our target list was well established, we ran sherlock over it. In general, after the execution of the first module, i.e. the exploration of the available data for a given candidate (see Section 2.1), we triggered a standard search with no initial fine-tuning. However, we had to refine the search parameters for some targets due to several factors, such as poor photometric quality and lack of data from the official pipelines. We provide a brief explanation for these particular cases.

• TIC 146520535 and TIC 464646604: In our multidetrend procedure, the window size ranges from 0.2 to 1.2 d by default, which we found optimal for most cases. However, due to the strong short-term variability of these light curves, we had to adapt the minimum and maximum window sizes to lower values from 0.25 to 0.75 d.

• TIC 189380158, TIC 281408474, and TIC 33595516: The data for these targets were not available in MAST, impeding us from gathering them using lightkurve, which is the default route to get the data. Hence, we switch to eleanor long-cadence data by specifying the author in the yaml file.

sherlock recovered 102 out of 104 transiting planets in the catalogue, making a count of two false negatives. These two were Kepler-89 e (TOI-4581.02), which was present in the TESS data as a single transit and hence, sherlock mismatched it with a different period, and K2-285 c (TOI-4538.01), which was entirely missed in favour of a different non-promising signal. With this data, the final recall made up a value of 0.98 for our catalogue. However, sherlock was also able to recover two other known planets that were not alerted with a TOI: the TOI-942 c and Kepler-89 c. For the TOI-942 system, TESS alerted one candidate with an orbital period of 4.3 d (TOI-942 b), and the follow-up of the system conducted by Carleo et al. (2021) found an outermost transiting planet with an orbital period of 10.16 d (TOI-942 c), which corresponds to our detection. In the case of the Kepler-89 system (Marcy et al. 2014), TESS alerted one candidate that corresponded to Kepler-89 d with an orbital period of 22.34 d, and a single transit that belonged to Kepler-89 e. In our experiment, sherlock also found a planetary candidate with an orbital period of 10.42 d, which corresponded to Kepler-89 c.

On top of this, sherlock spotted four interesting signals that could be attributable to planets and for which there is not any previous information. These new candidates were found in the systems WASP-16, HAT-P-27, HAT-P-26, and TOI-2411. We followed all the steps to verify their planetary nature by using the modules available in sherlock as described in Section 2. In the following subsections, we describe these new candidates and the confirmation steps that we conducted to claim them as planetary candidates with high probabilities of being genuine planets. Table 3 shows the summary of these new candidates, and the planetary systems architectures are presented in Fig. 5. We release these candidates as Community TESS Objects of Interest (CTOIs), and will be listed in the ExoFOP CTOIs list.⁹

Hence, despite the small size of our sample, these results allow us to confirm the high performance of sherlock as a planetary search engine that can be used to explore space-based data efficiently, not only to recover known planets but also to search for new ones that remained unnoticed by the community.

3.1 WASP-16

WASP-16 is a Sun-like star with |$T_{\textrm {eff}}=5700\pm 150\, \mathrm{ K}$|⁠, |$R_{\star }=0.946\pm 0.054\, \mathrm{ R}_{\odot }$|⁠, and |$M_{\star }=1.022\pm 0.101 {\rm M}_{\odot }$| observed by the Wide Angle Search for Planets (WASP) project (Pollacco et al. 2006), in particular using WASP-South, between 2006 and 2008. These observations hinted at the presence of a giant planet that was later confirmed using spectroscopic observations with the CORALIE spectrograph in 2008 and 2009. The results confirmed the planetary nature of this signal, yielding |$R=1.008\, R_{\mathrm{ Jup}}$|⁠, |$M=0.855\, M_{\mathrm{ Jup}}$|⁠, and with an orbital period of 3.12 d. This planet is the hot Jupiter WASP-16 b (Lister et al. 2009).

WASP-16 was observed by the TESS mission in Sector 11 between 2019 April and May, using the 120 and 1800 s cadences (TIC-46096489). The search performed with sherlock was conducted using the standard parameters, and we recovered in the first run the signal corresponding to WASP-16 b, with S/N = 97.9 and SDE = 18.6 in the 11 light curves (see Section 2.1.3). During the second run, sherlock did not spot any promising signal, but events that happened close to the borders of the light curves, which clearly hinted at a systematic origin. During the third run, we found a weak signal with S/N = 4.7 and SDE = 13 with a period of 10.44 d. While this signal is just below our default detection threshold of S/N |$\ge$| 5, it was found in all the 11 light curves processed, which reveals a non-detrend dependence. Hence, we vetted the signal and excluded the possibility of it being produced by any systematics (see Section 2.3). Then, we conducted the statistical validation (see Section 2.5), and we got FPP = 0.26 and NFPP = 0.0036 values, slightly above the likely planet area defined by Giacalone et al. (2021). Subsequently, we conducted the Bayesian fitting (see Section 2.6) and compared the Bayes factor (⁠|$\Delta \ln Z$|⁠) of two models: (1) only with the WASP-16 b planet, i.e. the null hypothesis, and (2) including the second potential planet found by sherlock, i.e. a more complex scenario with two planets. Then, we got |$\Delta \ln Z = \Delta \ln Z_{\rm 2\hbox{-}planet}-\Delta \ln Z_{\rm 1\hbox{-}planet}\sim 6.7$|⁠, what favours the two-planet scenario (Kass & Raftery 1995). Considering the global results of these analyses, we concluded that the signal found with a 10.44 d is a strong planetary candidate named CTOI-46096489.02, which would correspond to a mini-Neptune of |${\sim}$|2.2 R|$_{\oplus }$|⁠. Its main parameters from the fitting can be consulted in Table 3. However, it is important to notice that the Bayes factor is suggestive but not decisive; hence, extra observations using high-precision photometry would be needed to firmly confirm the planetary nature of the CTOI-46096489.02 candidate. Due to the shallowness of its transit of |${\sim}$|0.54 ppt, ground-based instrumentation is inadequate for detecting such a subtle eclipse. Indeed, the limited sensitivity and atmospheric interference associated with terrestrial telescopes prevent accurate measurement of these weak signals. Therefore, to obtain precise and reliable data that allow us to confirm the planet, it is essential to utilize space-based instrumentation, such as CHaracterising ExOPlanet Satellite (CHEOPS; Benz et al. 2021), which can bypass atmospheric distortions and provide the necessary sensitivity to observe these shallow transits effectively (see e.g. Delrez et al. 2021). Unfortunately, since the last observational data available come from TESS Sector 11, the current uncertainty in the transit time, |$\Delta T$| (see Section 2.6), is |${\sim}$|100 h, which makes any time-critical observation impracticable. Another TESS-like campaign would be needed to mitigate this lack of accuracy in the current transit times. According to the TESS-point web tool,¹⁰ WASP-16 will be reobserved in 2025 April, what will offer the opportunity to recover this candidate and trigger an accurate follow-up campaign that will allow us to confirm or refute its planetary nature.

The two-planet model is displayed in Fig. 8, the posterior distribution is provided in Appendix  C, Table C1, and the convergence of the model is visualized in Fig. C1.

3.2 HAT-P-27

HAT-P-27 is a G8 dwarf star with a |$T_{\rm eff}=5250\pm 100\, {\rm K}$|⁠, |$R_{\star }=0.90\pm 0.05\, {\rm R}_{\odot }$|⁠, and |$M_{\star }=0.94\pm 0.04\, {\rm M}_{\odot }$| monitored on a nightly basis by The Hungarian-made Automated Telescope Network (HATNet; Bakos et al. 2011) from 2009 January to August. These photometric observations were combined with nine spectroscopic measurements using the High Resolution Echelle Spectrometer (HIRES) instrument on the Keck I telescope from 2009 December to 2010 June. The results were presented in Béky et al. (2011), reporting the discovery of HAT-P-27 b, a hot Jupiter with a size of 1.038 R|$_{\rm Jup}$|⁠, a mass of 0.66 M|$_{\rm Jup}$|⁠, and an orbital period of 3.04 d.

HAT-P-27 was observed by the TESS mission in Sector 51 in 2022 May, using the 120 and 600 s cadences (TIC-461239485). We used the 120 s data in our search with sherlock, allowing us to recover in all the light curves a transiting signal corresponding to the HAT-P-27 b planet with S/N = 37.3 and SDE = 12.6. During the second run, sherlock spotted a transit-like feature with a period of 1.20 d, with S/N = 10.3 and SDE = 29.8, also detected in all the light curves. We executed the vetting module for this signal, which allowed us to rule out its systematic origin. We noticed that spoc identified a TCE for a candidate with an orbital period of 2.40 d, i.e. the first harmonic of the signal that sherlock spotted. Hence, during the vetting stage, we carefully examined the depths of the odd/even eclipses and transit shapes in short and long cadences and found no remarkable difference among them. Moreover, the candidate did not show any problematic metric either on the transit source or the optical ghost checks, reinforcing our findings with a period of 1.20 d as the more plausible solution. The TCE found by spoc did not become a TOI due to its low S/N (⁠|${\lt} $|10) and slightly centroids offset. Next, we launched the statistical validation module, yielding values of FPP = 0.18 and NFPP = 0.00054, which fall on the likely planet area.

Following the same procedure as in the case of the WASP-16 system, we conducted the Bayesian fitting and compared the Bayes factor between two models: (1) only with the HAT-P-27 b planet (null hypothesis) and (2) including the second potential planet found by sherlock. We got |$\Delta \ln Z = \Delta \ln Z_{\rm 2 \hbox{-}planet}-\Delta \ln Z_{\rm 1 \hbox{-}planet}\sim 64.6$|⁠, which strongly favours the two-planet scenario over the simplest with a single planet. This candidate is then the CTOI-461239485.02, which corresponds to a hot Neptune-size planet of |${\sim}$|4.3 R|$_{\oplus }$|⁠. Its parameters from the fitting can be consulted in Table 3. The transit depth found for this candidate is |${\sim}$|2.1 ppt, which is deep enough for small- and medium-sized ground-based telescopes. Moreover, its orbital period is 1.20 d, which eases scheduling and follow-up efforts due to the large number of possible observable events. However, the uncertainty in the transit time is |$\Delta T\sim$| 3.3 h, which is challenging to conduct time-critical observations from the ground but still possible. An alternative approach would be to conduct long observations from space-based telescopes such as CHEOPS, allowing for scanning up to 3|$\sigma$| the predicted time without the day–night cycles, making a putative detection easier. All these considerations make CTOI-461239485.02 an excellent candidate to be confirmed or refuted with relatively low effort.

The two-planet model is displayed in Fig. 8, the posterior distribution is provided in Appendix  C, Table C1, and the convergence of the model is visualized in Fig. C2.

Despite the vetting and validation tests passed by CTOI-461239485.02, the resulting planetary configuration with a hot Jupiter at 3.04 d and an innermost candidate at 1.20 d might be questionable due to stability concerns. For this reason, we executed the stability test (see Section 2.4) assuming the nominal values for the star and HAT-P-27 b planet as reported by Béky et al. (2011), and for the new candidate using a grid of 40 eccentricities from 0 to 0.6 and 40 planetary masses from 9.5 to 32 |${\rm M}_{\oplus }$| derived from the radius using a mass–radius relationship (17|$^{+13.8}_{-8.1}\, {\rm M}_{\oplus }$|⁠; Chen & Kipping 2017). This configuration allowed us to build a stability map exploring 1600 scenarios in the |$e_{\rm p}\!-\!M_{\rm p}$| parameter space. Then, the MEGNO value is computed six times in each scenario by modifying randomly the orbital parameters. Hence, in total, the stability map explored 9600 scenarios. Each scenario was integrated over |${\sim} 10^{7}$| orbits of the outermost planet. The stability map in Fig. 6 shows that the HAT-P-27 system can be stable for eccentricities up to 0.19 for the full range of masses explored for the planetary candidate. For larger eccentricities, the system becomes unstable for any of the masses within our explored grid. This analysis allowed us to conclude that despite the unusual architecture of the two-planet solution proposed here, consisting of a hot Jupiter and an inner hot Neptune, the system has ample room for stable configurations, resembling the WASP-84 system, where recently, an inner mini-Neptune has been found accompanying the previously known hot Jupiter WASP-84 b (Maciejewski et al. 2023).

3.3 HAT-P-26

HAT-P-26 is a K1 star with a |$T_{\rm eff}=5079\pm 88\, {\rm K}$|⁠, |$R_{\star }=0.79_{-0.04}^{+0.10}\, {\rm R}_{\odot }$|⁠, and |$M_{\star }=0.82\pm 0.03\, {\rm M}_{\odot }$| observed nightly by the HAT-5 and HAT-6 telescopes of the HATNet project between 2009 January and August, which hinted at the presence of a transiting planet with an orbital period of 4.23 d. The planetary origin of this signal was later confirmed using the HIRES spectrograph on the Keck I telescope from 2009 December to 2010 June. All these observations together allowed the confirmation of the low-density Neptune-mass planet HAT-P-26 b (Hartman et al. 2011) with a size of 0.565 R|$_{\rm Jup}$|⁠, a mass of 0.059 M|$_{\rm Jup}$|⁠, and an orbital period of 4.23 d.

This star was observed by the TESS mission in Sector 50 in 2022 May, using the 120 and 600 s cadences (TIC-420779000). During our search with sherlock, we fully recovered the signal corresponding to HAT-P-26 b with S/N = 42.2 and SDE = 15.6 in the first run. During the second run, sherlock spotted an interesting signal with a period of 6.59 d. While this signal was slightly below our detection threshold, with S/N = 4.9 and SDE = 7.4, it was confidently found in all the light curves. This motivated us to run the vetting module, where we did not find any solid potential systematic origin that explained this signal; however, it is important to note that due to its very low S/N, some of the tests were not conclusive. Still, we decided to conduct the statistical validation, which yielded FPP = 0.24 and NFPP = 0.0075, falling in the likely planet area.

We then executed the fitting module to refine the planet parameters and compare the Bayes factor of the two models as for the other candidates. We found |$\Delta \ln Z = \Delta \ln Z_{\rm 2\hbox{-}planet}-\Delta \ln Z_{\rm 1\hbox{-}planet}\sim 4.7$|⁠, which slightly favours the two-planet scenario over the single-planet model. The candidate, named CTOI-420779000.02, would correspond to a mini-Neptune of |${\sim}$|2.0 R|$_{\oplus }$| (see planet parameters in Table 3), producing shallow transits of only |${\sim}$|0.6 ppt, which hinders any confirmation using ground-based telescopes. However, its confirmation would be feasible, as stated for the CTOI-46096489.02 in the WASP-16 system, by employing space-based telescopes such as CHEOPS. The uncertainty in the transit time is |$\Delta T\sim$| 22.7 h, which is too large for time-critical observations, and unfortunately, according to the TESS-point web tool, this star will not be revisited in future cycles, preventing us from reducing this uncertainty using upcoming TESS data. However, the orbital period of this candidate of only 6.59 d would favour a filler follow-up camping to recover the planet by accumulating observations. This strategy is also applied in the FATE project using CHEOPS observations (Van Grootel et al. 2021) to find unknown transiting planets with orbital periods shorter than 10 d. Following this strategy, it would be possible to confirm or refute the existence of this planet.

The two-planet model is displayed in Fig. 8, the posterior distribution is provided in Appendix  C, Table C1, and the convergence of the model is visualized in Fig. C3.

As for HAT-P-27, the candidate’s orbital period was close to the known planet HAT-P-26 b. That architecture might be physically impossible and discarded by exploring the system stability. We, hence, followed the same strategy as for the HAT-P-27 system, building a stability map to explore 1600 scenarios in the |$e_{\rm p}\!-\!M_{\rm p}$| parameter space, displayed in Fig. 7. In this analysis, the candidate’s eccentricity must be lower than 0.05 to ensure the system’s stability for all the planetary masses explored, while the scenarios with the candidate’s mass below 6 |${\rm M}_{\oplus }$| show stable solutions for eccentricities up to 0.09 in most scenarios.

A particularly interesting point of the HAT-P-26 system is that von Essen et al. (2019) reported the detection of a drift in the spectroscopic data alongside a weak curvature in the timing residuals in the photometric data. The authors conducted a follow-up of the planet and obtained a |${\sim}$|4 min amplitude of transit timing variations (TTVs). While the authors prudently did not attribute these TTVs to a planetary companion, they gathered enough information to rule out the spot-induced origin of such TTVs. In a more recent study, A-thano et al. (2023) revisited the transit timings of HAT-P-26 b over 7 yr, reporting a TTV amplitude of |${\sim}$|2 min, which might be attributable to the presence of an additional planet in the system of |${\sim} 6.36\, {\rm M}_{\oplus }$| with an orbital period of 8.47 d, i.e. in a first-order mean-motion resonance 2:1. However, there is not any hint of a transiting planet under this scenario. On the contrary, our solution points to the transiting nature of a potential planetary candidate of |${\sim}$|1.96 |${\rm R}_{\oplus }$| (predicted mass from mass–radius relationships of 4.8|$^{+3.5}_{-2.0}\, {\rm M}_{\oplus }$|⁠) with an orbital period of 6.59 d, i.e. within |${\sim} 6 \, {{\rm per\, cent}}$| of the first-order resonance 3:2. The findings presented here can significantly encourage the scientific community to continue monitoring the HAT-P-26 system, to ultimately unravel its planetary architecture.

3.4 TOI-2411

TOI-2411 (TIC-10837041) is a K7 dwarf star with |$T_{\rm eff}=4099\pm 123\, {\rm K}$|⁠, |$R_{\star }=0.68\pm 0.02\, {\rm R}_{\odot }$|⁠, and |$M_{\star }=0.65\pm 0.02\, {\rm M}_{\odot }$|⁠. This star was monitored by the TESS mission in Sector 3 (2018 October) using 1800 s cadence and in Sector 30 (2020 October) using 120 and 600 s cadences. Orbiting TOI-2411 was alerted in 2020 November a planetary candidate, which was validated by Giacalone et al. (2022) by combining high-resolution images collected using the SOAR/HRCam, Palomar/PHARO, and Keck/NIRC2 instruments, ground-based multiband photometry (i′ and r′) taken employing the 1-m Sinistro telescope at Las Cumbres Observatory and the 0.7-m CDK telescope at the Mt. Kent Observatory, and statistical arguments. The planet, TOI-2411 b, is a super-Earth of |$1.68\pm 0.11\, {\rm R}_{\oplus }$| in an ultrashort period orbit of 0.78 d.

During our sherlock search, we found in the first run a signal of 18.75 d with an S/N = 14.9 and SDE = 14.1 in 9 of the 11 light curves processed. In the second run, we found a signal corresponding to TOI-2411 b with S/N = 8.4 and SDE = 16.7 in all the light curves. We executed the vetting module for the 18.75 d signal and found one red flag that might indicate a false positive: a momentum dump near the first transit. Momentum dumps might originate flux dimmings that can confuse the search algorithm; however, their shapes do not necessarily resemble transit-like features. We carefully visually inspected the two events in the 120 s data, and they both look similar in shape, duration, and depth. Moreover, during the vetting, was also inspected the 600 s cadence, and the signal corresponding to the 18.75 d was also found. Hence, we concluded that the momentum dump, or any other systematic, was not the origin of the transit-like feature detected by sherlock. Next, we conducted the statistical validation, which yielded values of FPP = 0.022 and NFPP = 0.0, falling in the likely planet area.

We then performed the fitting analysis and compared the two models as described for the previous candidates. We found |$\Delta \ln Z = \Delta \ln Z_{\rm 2\hbox{-}planet}-\Delta \ln Z_{\rm 1\hbox{-}planet}\sim 66.1$|⁠, which highly supports the two-planet model. The candidate, CTOI-10837041.02, is found to be a mini-Neptune of 2.88 |${\rm R}_{\oplus }$|⁠, whose transit depth is estimated to be |${\sim}$|1.77 ppt (see Table 3), which can be detected using small- and medium-size telescopes. However, its moderately long orbital period of 18.75 d and estimated uncertainty in the transit time of |$\Delta T\sim$| 6.7 h make it very challenging to confirm its planetary nature utilizing ground-based facilities. Unfortunately, TESS will not revisit this star, and hence, as in the case of CTOI-420779000.02 (HAT-P-26 system), a filler program using CHEOPS would be a potential solution to recover and confirm the candidate as a genuine planet.

The two-planet model is displayed in Fig. 8, the posterior distribution is provided in Appendix  C, Table C1, and the convergence of the model is visualized in Fig. C4.

4 THE HIDDEN GEMS PROJECT

Nearby and low-mass stars (⁠|$M\lt 0.6\, {\rm M}_\odot$|⁠) represent the best opportunity nowadays to study in great detail small transiting planets (R|$_{\rm p}$||$\lt \ 4$| R|$_{\oplus }$|⁠). Indeed, these stars ease the characterization of small exoplanets, not only in terms of radius and mass but also in terms of internal structure and atmospheric properties, which encodes information about physical and chemical properties, in turn providing insights into planets’ formation and evolutionary history (see e.g. Pudritz, Cridland & Alessi 2018; Triaud 2021).

Summarizing, there are three critical points that support this reasoning. (1) their small sizes and lower masses compared to solar-like stars produce larger transit depths and induce larger radial velocities for a given planet, allowing for an accurate determination of planetary bulk properties, such as the radius and mass, with precisions better than |${\sim} 4$| and |${\sim} 25 \, {{\rm per\, cent}}$|⁠, respectively. This allows us to derive the mean density and give a rough estimation of the bulk composition and internal structure (Otegi et al. 2020). This advantage has been the driver of many transits and Doppler surveys such as MEarth (Nutzman & Charbonneau 2008), SPECULOOS (Sebastian et al. 2021), CARMENES (Quirrenbach et al. 2014), and SPIRou (Artigau et al. 2014), to name a few. (2) Atmospheres of planets orbiting these small stars can be explored using transit spectroscopy, which relies on the chromatic signature of the chemical species present in the atmosphere, detectable with current and future space-based observatories such as JWST and Ariel. Transit spectroscopy allows three configurations to study planetary atmospheres: transmission spectroscopy, emission spectroscopy, and phase curves (Madhusudhan 2019). The possibility of atmospheric characterization of small planets with any of these three techniques has been the drive behind the detection of new planets around such stars, notably thanks to the TESS and Kepler/K2 missions. It has even now become one of the main science goals of the community, with its importance recognized at the highest level with a recommendation by the JWST Exoplanet Working Group to dedicate over 500 h of Director’s Discretionary Time to this end.¹¹ (3) Low-mass stars are the most common in the solar neighbourhood (Winters et al. 2015), representing up to 70 per cent of the stellar content of our Galaxy. Moreover, they have a lower binary fraction than solar-like stars. In addition, computational simulations, as well as Doppler and transit surveys, have concluded that the fraction of stars hosting small planets increases with decreasing mass (see e.g. Bonfils et al. 2013; Dressing & Charbonneau 2015; Alibert & Benz 2017). By combining the high frequency of low-mass stars with the occurrence rate of small planets orbiting them, there might exist |${\sim} 10\times$| more small planets orbiting these stars than around solar analogues. Hence, planets discovered orbiting low-mass stars are more representative of planet formation processes than planets orbiting single Sun-like stars (Triaud 2021).

Of particular interest are these systems hosted by a low-mass star where more than one transiting planet exists. In these cases, by studying different planets that were born in the same protoplanetary disc, we can get information not only about the planets but also about the formation and evolution of the whole planetary system. By aggregating information gleaned from the analysis of numerous multiplanetary systems, we will be able to formulate a robust planetary formation model that can resolve existing degeneracies and overcome the limitations currently constraining our understanding (see e.g. Raymond & Morbidelli 2022). This comprehensive model will eventually enable us to systematically decipher the underlying processes that govern planetary formation and evolution across different environments, unifying disparate observations and refining our predictions regarding the distribution and characteristics of exoplanets.

The full sample of transiting exoplanets around low-mass stars keeps on growing and now exceeds 1400¹² planets with 85 per cent of them being smaller than Neptune. As the race to detect new small planets amenable for detailed studies is still ongoing, we present the Hidden Gems project, which relies on the capabilities of sherlock to discover extra planets in low-mass star systems hosting at least one transiting planet.

Our target list is constructed by filtering stars with |$T_{\mathrm{eff}} \lt 5300$| K to include K and M dwarf stars and select a volume-limited sample within a distance of 50 pc, which totals 167 transiting planets in 107 unique systems at the time of writing. However, this catalogue will be updated with the upcoming releases of newly discovered planetary systems. Occurrence rate studies show that M dwarfs host at least two planets on average (e.g. Dressing & Charbonneau 2015; Tuomi et al. 2019), this corresponds to twice the number reported by Kunimoto & Matthews (2020) for K dwarfs. It implies, on a statistical level, that there could still be additional planets in these systems. Moreover, a transiting system has a favourable orbital inclination that enhances the likelihood of finding extra planet candidates. TESS offers an exquisite opportunity, thanks to its long baselines of observations, to search for concealed low S/N signals that might have been unnoticed. The search will be conducted on each system’s available TESS sectors and will be updated as soon as new data becomes public. This strategy increases the S/N of the putative detections, reaching smaller planetary radii and longer period candidates, as shown in Fig. 9.

By developing this project, we will eventually be able to pinpoint the presence of small transiting planets in the habitable zones of their stars within the solar neighbourhood, setting the stage for detailed studies. Hence, the ultimate objective of this project is to compile an accurate census of nearby transiting potential habitable planets.

5 CONCLUSIONS

Providing open-sourced software packages is a valuable contribution to the scientific community, especially in cases where the documentation is abundant, detailed, and focused on hands-on examples. With sherlock, we are not only doing that but also releasing a mature pipeline standardizing scientific procedures that increase the scientist’s productivity, help the reproducibility of the end-to-end processes, and empower individual teams to conduct their own planetary searches relying on space-based missions.

sherlock represents an advanced exoplanetary search pipeline that leverages a suite of robust and community-validated astrophysical packages including lightkurve, eleanor, tls, triceratops, and allesfitter. These are complemented by bespoke modules tailored to create a streamlined and integrated workflow. Initially, sherlock systematically acquires and pre-processes available data, setting the stage for effective exoplanet detection. Subsequent phases involve rigorous vetting of candidates to identify and eliminate potential false positives through comprehensive analysis and statistical validation techniques. The pipeline further refines estimates of planetary characteristics through sophisticated data modelling, enhancing the precision of the derived orbital and physical parameters. sherlock concludes the analyses by determining the best observational windows, crucial for planning ground-based confirmatory follow-ups. Each phase is encapsulated within dedicated modules, optimizing the accuracy and efficiency of the whole process of finding reliable planetary candidates. Still, we highlight that sherlock may need manual fine-tuning to operate to its fullest potential, making it not an entire automated search pipeline, requiring visual inspection of detected signals before passing through all the modules.

Our performance test showed that sherlock has an excellent recall of 0.98 when recovering spoc official alerts confirmed as planets, which makes it a very reliable pipeline for the community. During this test, we not only recovered 98 per cent of the planets but also found some candidates that were confirmed planets but did not have any associated alerts, such as TOI-942 c and Kepler-89 c, strengthening the conclusion of being a reliable planetary searching tool. On top of that, we found four new planetary candidates unnoticed until now that fulfilled all our tests, only lacking ground-based observation to confirm their planetary nature fully. These four candidates were found to be in the planetary systems WASP-16 (CTOI-46096489.02), HAT-P-27 (CTOI-461239485.02), HAT-P-26 (CTOI-420779000.02), and TOI-2411 (CTOI-10837041.02).

The vetting implemented by sherlock, based on the robust tests conducted by the spoc pipeline, reports a variety of metrics that are marked as red flags whenever they are found to be above specific thresholds. However, low data quality sometimes produces inconclusive results for these metrics, especially on the centroids and optical ghost diagnoses. Hence, in these cases, the red flags have to be taken as warning bells instead of a robust indication of a false positive, which needs further exploration by the user. sherlock’s vetting will be further developed in coming releases towards a more robust data exploration, producing a more reliable set of metrics even for very low S/N signals. In particular, we are currently exploring implementing a deep learning neural network similar to ExoMiner (Valizadegan et al. 2022), allowing us to discern genuine planets from false positives unbiasedly.

ACKNOWLEDGEMENTS

We thank the anonymous referee for the helpful comments. FJP acknowledges financial support from the Severo Ochoa grant CEX2021-001131-S funded by MCIN/AEI/10.13039/501100011033 and through the project PID2022-137241NB-C43. Funding for open access charge: Universidad de Granada / CBUA. We acknowledge support from the Spanish Ministry of Science through the project PID2019-107061GB-C64/SRA (Spanish State Research Agency/10.13039/501100011033). This publication benefits from the support of the French Community of Belgium in the context of the FRIA Doctoral Grant awarded to MT. VVG is an F.R.S. – FNRS Research Associate.

DATA AVAILABILITY

The data underlying this article will be shared on reasonable request to the corresponding author.

Footnotes

REFERENCES

APPENDIX A: TIC 169285097 sspe FIT SIGNALS

APPENDIX B: PERFORMANCE SURVEY CATALOGUE: TIC IDs

100566492, 102264230, 10837041, 111991770, 112099249, 120247528, 12862099, 136916387, 139375960, 142937186, 144700903, 146520535, 147797743, 157230659, 165297570, 168281028, 172598832, 175772482, 181949561, 188876052, 189380158, 189625051, 20096620, 201186294, 201604954, 206284139, 207077681, 21113347, 229047362, 237320326, 243921117, 245392284, 248853232, 250707118, 257567854, 262715204, 263179590, 269217040, 27039476, 273231214, 281408474, 283419843, 288246496, 291685334, 291751373, 292152376, 305048087, 32487566, 33595516, 336732616, 339522221, 347430350, 354619337, 368739874, 368805700, 369455629, 372172128, 374530847, 375942197, 37718056, 381856447, 391958006, 392476080, 398895470, 39926974, 399402994, 401604346, 420779000, 420814525, 425817867, 429302040, 431701493, 434226736, 436478932, 437704321, 441420236, 452808876, 453064665, 459730973, 460396820, 46096489, 461239485, 464646604, 466382581, 466390120, 467971286, 49043968, 51637609, 53728859, 54002556, 57984377, 6663331, 7088246, 90850770, 91051152, 9155187, 92449173, 9385460, 9443323, 97409519.

APPENDIX C: POSTERIOR DISTRIBUTIONS FROM THE FITS