Open Access
Issue
A&A
Volume 708, April 2026
Article Number A232
Number of page(s) 19
Section Planets, planetary systems, and small bodies
DOI https://doi.org/10.1051/0004-6361/202557806
Published online 09 April 2026

© The Authors 2026

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. This email address is being protected from spambots. You need JavaScript enabled to view it. to support open access publication.

1 Introduction

Over the past three decades, the field of exoplanet research has witnessed remarkable growth, with more than 6000 confirmed exoplanets as of October 20251. Of these, approximately 1100 have been detected using the radial velocity (RV) technique, while the vast majority of more than 4400 confirmed planets have been identified via the transit method both from ground-based observations (Pollacco et al. 2006; Bakos et al. 2004; Gillon et al. 2017), and the unprecedented contributions of NASA’s Kepler space telescope (Borucki et al. 2010) and the ongoing Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2015). Additional discoveries stem from complementary techniques, including direct imaging of wide-orbit companions (e.g. Marois et al. 2008) and gravitational microlensing, which is particularly effective for detecting planets around distant and faint stars (e.g. Gaudi 2012).

Statistical analyses of large samples of surveyed stars indicate that planets are ubiquitous, with occurrence rates exceeding 0.5 planets per FGK-type star for orbital periods ranging from one day to several hundred days. These estimates are supported by Doppler surveys (Howard et al. 2010) and transit surveys (Petigura et al. 2013; Kunimoto & Matthews 2020). Even higher occurrence rates have been found for planets orbiting M dwarfs, with estimates exceeding one planet per star (Cassan et al. 2012; Bonfils et al. 2013; Dressing & Charbonneau 2015; Gaidos et al. 2016; Mignon et al. 2025), and possibly increasing further from early- to mid-type M dwarfs (Hardegree-Ullman et al. 2019). The solar neighbourhood remains a prime hunting ground for exoplanets around low-mass stars, particularly M and K dwarfs, due to their abundance and relative proximity. Nearby stars are especially valuable for Doppler surveys and exoplanet characterisation because of their apparent brightness (resulting in higher signal-to-noise (S/N)observations). Many planet searches have therefore focused on them. Key efforts include the HARPS M-dwarf survey (Bonfils et al. 2013; Astudillo-Defru et al. 2017), and several programmes targeting K dwarfs more specifically, such as the California Planet Survey (Howard et al. 2010; Rosenthal et al. 2021), the KOBE project (Balsalobre-Ruza et al. 2025), and the SOPHIE northern survey (Hobson et al. 2018). These surveys have significantly expanded the known exoplanet population in this stellar mass regime and continue to refine the occurrence rates and dynamical architectures of exoplanet systems around low-mass stars. Nearby stars are also important targets for future observations with Extremely Large Telescopes (ELTs) such as ANDES at the ELT, Paranal, Chile (Palle et al. 2025), and are prime targets for direct-imaging campaigns with NASA’s Habitable Worlds Observatory (Harada et al. 2024) and ESA’s LIFE interferometer (Quanz et al. 2022), owing to the greater angular separation between the star and any potential exoplanet.

In this paper, we present extensive orbital and stellar activity analysis of archival spectroscopic data of the K-dwarf star GJ 1137 obtained with the High Accuracy Radial velocity Planet Searcher spectrograph (HARPS, Pepe et al. 2002; Mayor et al. 2003) mounted at the ESO 3.6 m Telescope in La Silla, Chile. GJ 1137 is a known exoplanet host, first reported by Lovis et al. (2005), based on just 16 HARPS RVs over a baseline of approximately one year. This dataset helped us identify GJ 1137 b, a Saturn-mass exoplanet, with a minimum mass of mb sin i = 0.37 MJ and an orbital period of Pb = 143.6 days. The discovery paper also noted a long-term RV trend that suggested an additional massive companion. The public HARPS-RVBank spectral products released by Trifonov et al. (2020); Perdelwitz et al. (2024) now provide an excellent opportunity for a more detailed analysis of the orbital architecture of GJ 1137 b, and for a search for additional planets in the system.

In our work, we show evidence for a 5640 d RV signal in addition to GJ 1137 b. While this signal initially appears consistent with a Jovian-mass companion, we demonstrate that it is part of long-term stellar activity that we see in other spectroscopic activity indices. Furthermore, we perform full stellar activity and planet signal modelling, and report improved orbital and mass estimates for GJ 1137b as well as stellar activity cycle properties for GJ 1137, and its rotational period. Finally, we present significant evidence of a short-period exoplanet, GJ 1137 c, with a period of 9.640.11+0.12Mathematical equation: $9.64^{+0.12}_{-0.11}$ days and a minimum mass of 5.120.69+0.70Mathematical equation: $5.12^{+0.70}_{-0.69}$ M.

The paper is organised as follows. In Sect. 2, we introduce the available data applicable to our study. In Sect. 3, we present revised physical properties of the K-dwarf star GJ 1137 by employing public data and state-of-the-art stellar parameter modelling. Our orbital and spectroscopic data analysis scheme and results are introduced in Sect. 4, and finally, in Sect. 5 we provide a summary and discussion of our findings.

2 Data

2.1 HARPS data

For our analysis of GJ 1137, we rely on spectroscopic data obtained with the HARPS spectrograph. Our dataset comprises 140 high S/N spectra collected between January 2004 and July 2017. Of these, 87 spectra were obtained before and 53 after the optical fibre upgrade in May 2015 (Lo Curto et al. 2015), hereafter referred to as HARPS-pre and HARPS-post, respectively. This upgrade improved the throughput of the HARPS instrument but changed the instrument profile. Thus, the pre- and post-fibre spectroscopic data products should be treated as if obtained from two separate spectrographs (see Lo Curto et al. 2015; Trifonov et al. 2020).

We use publicly available extracted RV and activity index data in the HARPS-RV Bank. HARPS-RV Bank contains precise RVs, as well as stellar activity indicators -Hα, NaID, NaIID, chromatic RV index CRX, and differential line width (dLW). All of these spectral products were derived via the SERVAL pipeline (Zechmeister et al. 2018). In addition to HARPS-RV Bank data, we use the RACCOON pipeline (Lafarga et al. 2020) to compute the cross-correlation function (CCF) RVs (Baranne et al. 1996) and derive other activity indicators. We first created a custom CCF mask using the high S/N template that SERVAL outputs, using default parameters for RACCOON to select the mask lines. We then cross-correlated this mask with the observed spectra, following the methods explained in Lafarga et al. (2020). In summary, RACCOON computes a CCF order-by-order, coadds all the order CCFs together into a final CCF, and fits a Gaussian function to this final CCF to derive the spectrum RV from the CCF centroid. RACCOON also provides the standard CCF activity indicators: full width at half maximum (FWHM) and contrast, derived from the Gaussian fit to the CCF, and bisector inverse slope (BIS), calculated directly from the bisector of the CCF, as described by Queloz et al. (2001).

For orbital analysis, we used RVs taken from the HARPS-RV Bank since these have been corrected for small, but significant RV systematics (see Trifonov et al. 2020). In addition, as shown by Anglada-Escudé & Butler (2012) and Zechmeister et al. (2018), for redder stars such as GJ 1137 the co-adding stellar template method used by SERVAL can deliver higher RV precision than the CCF method with a weighted binary mask used by RACCOON. For studying the activity index data, we used the activity indices taken from both SERVAL and RACCOON. All data used are listed in Table A.1, the full table available in the CDS.

2.2 Photometry data

TESS visited GJ 1137 (TIC54237857) during sectors 3, 36, 63, and 90. The star has a TESS-band magnitude of 7.472 ± 0.006 mag (TIC v8.2 catalogue; Stassun et al. 2019). In the available 120-second exposures, we expect it to produce a flux of ∼1.8 × 107 e/s−1, which is nearly two orders of magnitude above the nominal saturation threshold2. As a result, the target exhibits pronounced blooming artefacts, rendering both the standard simple aperture photometry (SAP) and presearch data conditioning simple aperture photometry (PDC-SAP) light curves (Smith et al. 2012; Stumpe et al. 2012, 2014) unsuitable for analyses of long-term stellar variability. Visual inspection of the light curves reveals significant contamination that likely comes from charge spillover from saturated pixels, leading to a high incidence of low-quality data points. Nevertheless, we performed a transit search using the transit least-squares (TLS) periodogram algorithm by Hippke & Heller (2019). No statistically significant periodic signals were detected.

We attempted to determine the rotational period of this bright star, using photometric data from the All-Sky Automated Survey for Supernovae (ASAS-SN; Shappee et al. 2014; Kochanek et al. 2017). Recently, the ASAS-SN reduction pipeline included a specialised module tailored for highly saturated stars (Winecki & Kochanek 2024). Data processed using this approach were obtained through the publicly available ASAS-SN Sky Patrol v1.03 by choosing the ‘Saturated stars’ option. For consistency, we used only the dataset in g band from observing seasons 2018 to 2025. Due to the proximity of GJ 1137 to the ecliptic, some of the photometric measurements are affected by lunar contamination. To mitigate this, we calculated the angular distance between the Moon and GJ 1137 for each ASAS-SN measurement and retained only those in which the Moon is more than 90 deg away on the sky. After this, the selected data was analysed using the generalised Lomb-Scargle algorithm (GLS; Lomb 1976; Scargle 1982; Zechmeister & Kürster 2009). We found no significant peaks in this data sample; see Fig. A.1.

3 Stellar parameters

GJ 1137 (HD93083, HIP52521) is a K2 IV-V dwarf star (Gray et al. 2006) with an apparent magnitude of V = 8.46 mag. The target was observed by HARPS and thus included in the HARPS-RVBank, which provides activity indicators and basic stellar atmospheric parameters. These parameters are estimated using templates created from the co-adding of HARPS spectra (see Perdelwitz et al. 2024, Sect. 2.1). Using these as a basis, we refined the stellar parameters of GJ 1137 by analysing its spectral energy distribution (SED) with astroARIADNE4 (Vines & Jenkins 2022). This package employs a Bayesian approach, utilising various SED models to mitigate model-specific biases. It also estimates stellar mass and age by interpolating MESA isochrones and Stellar Tracks (MIST, Dotter 2016; Choi et al. 2016). As input, we used photometric measurements in the following bands: GALEX NUV, Johnson U, B, and V, Tycho B and V, Gaia DR2 G, BP, and RP, TESS, 2MASS J, H, and K, Wise W1 and W2. These were fitted with atmosphere models from PHOENIX v2 (Husser et al. 2013), BT-Settl, BT-NextGegn, BT-Cond (Hauschildt et al. 1999; Allard et al. 2012), Castelli & Kurucz (Castelli & Kurucz 2003), and Kurucz (Kurucz 1993). We assign a normally distributed prior to the distance, using the estimate from Gaia EDR3 (Bailer-Jones et al. 2021). For interstellar extinction, we employ a flat prior ranging from zero to the maximum line-of-sight extinction. This was calculated using the galactic dust map by Schlafly & Finkbeiner (2011) in the Python package dustmaps5 (Green 2018). For the initial stellar parameters Teff, log g, and [Fe/H], we constructed normal priors using the estimates of Perdelwitz et al. (2024) in the HARPS-RVBank.

Our analysis demonstrates that GJ 1137 has an age of 10.711.25+2.24Mathematical equation: $10.71^{+2.24}_{-1.25}$ Gyr, a mass of 0.8360.025+0.023MMathematical equation: $0.836^{+0.023}_{-0.025}\,\unit\solarmass$ and a radius of 0.8370.018+0.026RMathematical equation: $0.837^{+0.026}_{-0.018}\,\unit\solarradius$. The derived parameters are consistent with the latest data in the literature (TICV8, Stassun et al. 2017). However, a significant variance is evident when comparing the stellar mass with the value used by Lovis et al. (2005). The revised higher stellar mass implies an increased minimum mass for the known planet GJ 1137 b, further discussed in Sect. 4.3. A summary of all stellar parameters is presented in Table 1, and the resulting SED fit can be seen in Fig. 1.

4 Analysis and results

4.1 Preliminary RV and stellar activity analysis

Our preliminary data analysis relies on the Exo-Striker exoplanet toolbox6 (Trifonov 2019), which provides convenient exoplanet data analysis tools, but also direct access to HARPS-RVBank data. Fig. 2 shows the HARPS-RVBank SERVAL nightly zero-point corrected RV time series, which clearly exhibit at least two periodic signals: a shorter-period signal with P = 144.75 ± 0.028 d attributed to the known exoplanet GJ 1137 b, and an additional long-period significant signal with P = 5640 ± 240 d , which sparked our interest in this system.

We performed a detailed time-series spectroscopic analysis to find that several critical activity indicators, including BIS, contrast, FWHM, dLW, and log R′HK, exhibit significant low-frequency periodicities that closely resemble those found in the RV residuals after modelling the GJ 1137 b signal. In our local example, the Solar System, both Jupiter and the Sun’s 11-year magnetic activity cycle induce RV variations on comparable timescales if the Sun were observed as a distant star. However, the disentanglement of these phenomena in RV is challenging, which raises caution regarding a possible stellar activity origin for the long-period variation.

As a next step in interpreting the long-term RV signal, we subtracted the contribution from the known planet and examined potential linear correlations between the residual RV variations and stellar activity indicators, using Pearson’s correlation coefficient. We identified positive or negative correlations that exceed 0.25 in all available activity proxies, with the strongest correlation observed with the R′HK index (r = 0.49). This moderate positive correlation suggests that the long-term signal may be at least partially driven by magnetic activity rather than an additional planetary companion. To investigate this further, with greater statistical rigour, we conducted a more sophisticated analysis incorporating long-term functions (LTFs) and Gaussian processes (GPs) modelling. Additional details and graphs on our preliminary analysis are provided in the Appendix.

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

SED fit for GJ 1137. The best-fitting model, BT-Settl, is displayed in black, with normalised residuals shown below. Blue points represent flux values from photometry, while purple diamonds show the flux from synthetic photometry in the same passband.

Table 1

Stellar parameters of GJ 1137 and their 1σ uncertainties.

Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

HARPS RV time series for GJ 1137. The data are divided into two subsets due to known instrumental RV offsets introduced by the HARPS optical fibre upgrade in May 2015. Purple circles represent RVs obtained before the fibre change, while green squares correspond to data taken after the upgrade. The observed RV variations are primarily driven by the signal of the known Saturn-mass exoplanet GJ 1137 b, with a period of P = 144.7 d. In addition, the RV data reveals a second, long-term signal from a potential Jovian planet with a minimum mass of approximately 1.3 MJ and a period P = 5640 ± 240 d. Subsequent analysis attributes this long-term signal to stellar magnetic activity.

4.2 Combined orbital modelling including stellar activity

We use the modelling framework described in Stefanov et al. (2025), hereinafter S25, which is generalised for N datasets (e.g. HARPS-pre, HARPS-post; i ∈ [0..N - 1]) and M physical quantities (e.g. RV, FWHM; j ∈ [0..M - 1]). The S25 framework includes: (i) LTFs data that follow physical or instrumental trends; (ii) offsets between different datasets; (iii) datasetdependent jitters for each physical quantity; (iv) planetary RV fitting; (v) stellar activity modelling using GPs. We now briefly recapitulate on GPs.

The activity originating from stellar rotation is regarded as a stochastic process, which must nevertheless be related to the stellar rotation period Prot. A popular way to model this process is through GPs (Rasmussen & Williams 2006), which assume a certain functional relationship of the correlation between measurements with respect to their difference in time. A common assumption for this correlation is that it is sensitive to Prot, as well as some exponential decay in time, which can be naïvely ascribed to the active-region evolution on the stellar surface. This assumption gives rise to the squared-exponential periodic (SEP) kernel of the form kSEP(Δt)=κ2exp[Δt22τ2sin2(πΔt/P)2η2],Mathematical equation: k_\text{SEP}(\Delta t) = \kappa^2\exp\left[ -\frac{\Delta t^2}{2\tau^2} -\frac{\sin^2\left(\pi\Delta t/P\right)}{2\eta^2} \right],(1)

where κ2 is the maximum amplitude of the kernel, τ is the timescale of coherence, P is associated with the stellar-rotation period Prot, and η describes the kernel-feature complexity (Haywood et al. 2014; Rajpaul et al. 2015; Angus et al. 2018). We hereafter refer to η as the ‘sinescale’ because of its algebraic analogy to the timescale τ.

Gaussian processes should be applied carefully, as they tend to overfit and absorb signals beyond Prot (Suárez Mascareño et al. 2023), and the latter may not necessarily pertain to stellar rotation. One way to combat this overfitting is to constrain GPs by simultaneously fitting them on RVs with one or more activity indices (Barragán et al. 2023). For GJ 1137, we found that fitting on the RV and FWHM time series alone gives a sufficiently well-described stellar-rotation activity. The size of our correlation matrix grows quadratically with the number of physical quantities, and the amount of additional information from the inclusion of R′HK as a second activity indicator does not justify the computational overhead. We restrict further GP models in this two-dimensional case; in these models, we shall follow the notation of S25 where j = 0 is assigned to RVs and j = 1 is assigned to the FWHM time series.

A common way to fit GPs on two-dimensional data follows the FF′ formalism that was defined in Aigrain et al. (2012) and extended in Rajpaul et al. (2015). For measurements yi,j at time ti,j, we solve the system of equations |yi,0=A0G(ti,0)+B0(dG/dt)t=ti,0yi,1=A1G(ti,1)i,Mathematical equation: \left|\begin{array}{l} y_{i,0}=A_0\,G(t_{i,0})+B_0\,(\mathrm{d}G/\mathrm{d}t)_{t=t_{i,0}}\\[0.2ex] y_{i,1}=A_1\,G(t_{i,1}) \end{array}\right.\qquad \begin{array}{l} \forall\, i, \end{array}(2)

where Aj and Bj are quantity-specific fit parameters, and G(ti,j) is the kernel estimation of the stellar activity at given times. This regime is known to resolve degeneracies when there is a significant correlation between physical quantities (Suárez Mascareño et al. 2020). We use the S+LEAF MultiSeriesKernel procedure to realise it (Delisle et al. 2020, 2022). For all models discussed in this work, we used ReactiveNestedSampler, a nested-sampling integrator provided by UltraNest (Buchner 2021), in order to infer the posterior distribution of the parameters and numerically evaluate the Bayesian evidence Z (Skilling 2004). In every inference instance of Nparam model parameters, we required 40Nparam live points and a region-respecting slice sampler that accepted 2Nparam steps until the sample was considered independent.

Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Raw time series of mean-subtracted: (a) RV, (c) FWHM, (e) R′HK. Measurements are marked depending on the data source: blue circles for HARPS-pre, and green squares for HARPS-post. (b,d,f) Associated wide-period GLSPs of RV and activity indicators. Three false alarm probability (FAP) levels: 10%, 1%, and 0.1%, split GLSP ordinates in bands of different colour. We highlight the three most prominent peaks in each GLSP.

4.2.1 Joint analysis of RVs, FWHM and R′HK

Figure 3 shows the RV, FWHM and R′HK raw-data time series and their generalised Lomb-Scargle periodograms (GLSPs; Lomb 1976; Scargle 1982; Zechmeister & Kürster 2009; VanderPlas & Ivezic 2015) in the period range 2-6000 d. Both the time series and GLSP panels offer a wealth of information about the behaviour of GJ 1137. Firstly, the FWHM and the R′HK time series both show similar qualitative evolution, with a long-term sinelike component of a period near 5000d (Figs. 3a,c). In addition, both physical quantities appear to follow linear or higher-order polynomial trends. This is supported by their respective GLSPs (Figs. 3d,f), which both show unresolved, yet significant, power excess at 6000 d and longer periods. We stress that a similar 5000 d cycle can be faintly discerned in the RV time series (Fig. 3a), but these time series are dominated by GJ 1137 b (144 d; Lovis et al. 2005), which is the main source of scatter, as well as the dominant signal in the GLSP (Fig. 3b). We note the presence of additional peaks in all physical quantities (e.g. between 41.1-46.2 d), but the long-term cycle must be modelled first to discuss those.

4.2.2 Presence of a long-term stellar activity cycle

We took all HARPS measurements in RV, FWHM, and R′HK and composed a three-dimensional model that contained a single elliptical-orbit planet and the following three-dimensional LTF: fj(ti,j)=n=0pαn,jti,jn+kcyc,jsin[2π(ti,j+Pcycφcyc,j)Pcyc],Mathematical equation: f_j(t_{i,j}) = \sum_{n=0}^{p} \alpha_{n,j} t^n_{i,j} + k_{\text{cyc, }j} \sin\left[\ \frac{2\pi(t_{i,j}+P_\text{cyc}\,\varphi_{\text{cyc, }j})}{P_\text{cyc}} \right],(3)

where the first term models a p-dimensional polynomial trenc with coefficients αn, j, and the second term is a common cycle of amplitude kcyc, j and phase φcyc, j, with a common cycle period Pcyc. That is to say, we restrict the long-term behaviour of each physical quantity to be modelled under a different polynomial but we impose an additional sine term with a common period in all of them. No GPs were involved in this stage. We used the following period priors: Ulog(103,104) d for the sine component in the LTF (i.e., the long-term cycle) and U (140,150) d in accordance with the confirmed planetary companion GJ 1137 b.

HARPS-pre is known to have had a slowly drifting focus before the change of fibre. That manifested as a trend in the FWHM data (Zechmeister et al. 2013; Dumusque 2018; Costes et al. 2021). We found that a linear polynomial is enough to catch our FWHM trend in HARPS-pre (i = 0, j = 1). Having accounted for this in the LTF as a separate linear-trend term β our best solution in terms of ∆ ln Z is an LTF with a zero-order RV trend, a zero-order FWHM trend, and a quadratic R′HK trend That is to say, our final choice of LTF is fj(ti,j)=α0,j+δi0δj1βti,j+kcyc,jsin[2π(ti,j+Pcycφcyc,j)Pcyc],Mathematical equation: f_j(t_{i,j}) = \alpha_{0,j} + \delta_{i0}\delta_{j1}\beta t_{i,j} + k_{\text{cyc, }j} \sin\left[\ \frac{2\pi(t_{i,j}+P_\text{cyc}\,\varphi_{\text{cyc, }j})}{P_\text{cyc}} \right],(4)

where δ is the Kronecker delta function. Table A.2 lists all parameters in the model, together with their adopted priors and inferred posteriors. We measure a cycle period of 5320150+170dMathematical equation: $5320^{+170}_{-150}\,\unit\day$ with a well-behaving unimodal posterior shape. We obtain the following RV, FWHM, and R′HK semi-amplitudes of the cycle: 14.2 ± 0.6ms−1, 18.71.7+1.8ms1Mathematical equation: $18.7^{+1.8}_{-1.7}\,\unit{\metre\per\second}$, and 6.590.79+0.89ppmMathematical equation: $6.59^{+0.89}_{-0.79}\,\unit\ppm$ respectively - all of them standing at least 8σ away from the zero. These results underline a strong 15 yr cycle in RV and activity. The cycle phase somewhat differs for different quantities - we obtain RV, FWHM, and R′HK phases of 0.8940.020+0.020Mathematical equation: $0.894^{+0.020}_{-0.020}$, 0.8880.022+0.021Mathematical equation: $0.888^{+0.021}_{-0.022}$, and 0.7970.020+0.021Mathematical equation: $0.797^{+0.021}_{-0.020}$ respectively. RV and FWHM are in phase within 1σ. The cycle appears to lag slightly behind in R′HK. While such phase lags are observed for the Sun (Collier Cameron et al. 2019) and other stars (Burrows et al. 2024), our posterior may be deceiving, since this baseline model does not account for other sources of stellar variability (e.g. stellar rotation with GPs).

Figure 4 shows the model fit against data, as well as the periodograms of the residuals, which accounted for the model jitter. The modelling of the long-term cycle and GJ 1137b was enough to cause all significant peaks beyond 100 d to disappear. On the short-period end, we see that a 43 d signal remained significant in FWHM and R′HK. This would be consistent with our expectations for the rotational period Prot - the mean and standard deviation of our log10 R′HK measurements form an estimate of −4.98 ± 0.10, which translates as Prot = 31 ± 7d through the relation in Suárez Mascareño et al. (2015). The FWHM GLSP displays a strong 26.4 d signal that may be the first harmonic of Prot. With regard to RV, we notice two peaks in 7.0 d and 9.6 d that soared to <10% FAP (Baluev 2008).

Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Time series (markers) against our LTF model (black), over the whole dataset baseline in: (a) RV, (c) FWHM, (e) R′HK. We implicitly correct for the HARPS-pre FWHM drift before visualising. (b,d,f) Associated GLSPs of the residual time series, accounting for the model jitter.

4.2.3 Presence of short-term stellar activity and an additional RV signal

We continued with the aforementioned LTF model and augmented it with a stellar-activity component as described in Sect. 4.2. We chose to test for four approximations of the SEP kernel that S+LEAF offers: the Matérn 3/2 exponential periodic (MEP) kernel and the exponential-sine periodic (ESP) kernel with two, three or four harmonics (hereafter ESP2, ESP3, ESP4). We additionally tested for a second planetary signal with a period prior of Ulog(1, 100) d and a semi-amplitude prior of U (0,10) ms−1. For this second potential RV contribution, we independently tested for a circular and a Keplerian orbit. The variation of four stellar-activity kernels (MEP, ESP2, ESP3, ESP4) and the three planetary configurations (Keplerian, Keplerian+circular, Keplerian+Keplerian) resulted in a grid of 12 models. These models were fitted on the conjunction of RV and FWHM data. We modelled stellar activity in the multidimensional regime, with RV having a gradient with respect to FWHM. We used our aforementioned rotation-period estimate through Suárez Mascareño et al. (2015) (31 ± 7 d) to impose a rotation-period prior Prot = U (10,52) d and the uninformed timescale prior τ = Ulog (10,104) d.

Figure 5 provides a comparison of models in our grid search using Bayesian evidence ln Z. Entries are arranged so that planetary configurations increase in complexity downwards, while stellar-activity kernels increase in complexity rightwards. In the figure, we highlight the most appropriate model: an MEP kernel, with a Keplerian component for GJ 1137 b and an additional circular-orbit signal. What follows is the justification for this selection.

The average improvement across kernels from a one-Keplerian model is: ∆ lnZ = 5.6 for an additional circular-orbit signal and ∆ ln Z = 2.4 for an additional Keplerian signal. In terms of specific kernels, we get the following ∆ ln Z improvements for an additional circular-orbit signal: 5.4 for MEP, 5.1 for ESP2, 5.0 for ESP3, 6.7 for ESP4. All of these improvements are great enough to be interpreted as overwhelming evidence of a second RV signal in the data (ln Z > 5; Jeffreys 1946).

We compare the period, phase, semi-amplitude, eccentricity, and argument-of-periastron posteriors of the second signal in Tables A.3, A.4, A.5, A.6, and A.7 respectively. These tables demonstrate high posterior stability across various planetary configurations and stellar-activity kernels. All posteriors agree within 1σ for any configuration-kernel pair. In two-Keplerian models, none of the inferred eccentricity posteriors is well defined with respect to zero (Table A.6). This is consistent with the slightly lower evidence relative to Keplerian+circular models in Fig. 5, and reinforces that the new RV signal is indeed sine-like.

We checked if GP parameters were consistent across the model grid as well. Figure A.3 visually compares the posteriors of all GP parameters and the average RV jitter per point between models. The posteriors are qualitatively the same across planetary configurations and kernels, with the following exceptions. The rotation periods change median from about 30 d to about 32 d with the inclusion of a second RV signal. This is not too surprising, given that the period of the latter is of the same order (near 9.6 d). The ESP4 models in particular tend to struggle with bimodality, as they try to catch a weak Prot mode near 40 d (Figs. A.3a, o). However, MEP-kernel models tend to overestimate timescales (Figs. A.3b, i, p). At the same time, they infer a smaller RV jitter than the ESP family (Figs. A.3g, n, u).

We assessed the significance of the 9.6 d RV signal following Hara et al. (2022b). Their formalism uses the likelihood of vectors drawn from the posterior distribution in order to compute the probability that there are no planets within a given orbital-frequency element. The authors refer to this metric as the false inclusion probability (FIP). We re-ran a derivative of our selected model with three planet components: the same Keplerian that captures GJ 1137 b, as well as two circular-orbit components that share the period prior Ulog(1,100) d. We computed the FIP between 1-100 d for an angular-frequency step ∆ω = 2π/5 T, where T is the total time series baseline. The novel 9.6 d RV signal is assigned a FIP smaller than 10−5, which is much more significant than the 10−2 threshold proposed by Hara et al. (2022b). We found no other RV peaks with significant FIPs.

Then, we tested whether the 9.6 d RV signal is not simply a power excess that our GP kernels could not absorb. We followed the treatment of S25, Sect. 4.5.3. We started with our selected model (Keplerian+circular MEP), removed the circularorbit RV component, and instead augmented an extra sine term in the LTF (Equation (4)) that acts on both the RV and the FWHM. In this way, we can still model for a long-term cycle with the original sine term in the LTF - and at the same time try to guide the second sine term 9.6 d to test whether there is a significant signal in the FWHM too. This derivative model used the same wide period prior Ulog(1,100) d. Its evidence scores worse by ∆ ln Z = −2.3 relative to our selected model. The short-term cycle designed to catch the 9.6 d RV signal actually locates it (Fig. A.5). The 9.6 d signal is assigned a FWHM semi-amplitude of 0.620.44+0.66ms1Mathematical equation: $0.62^{+0.66}_{-0.44}\,\unit{\metre\per\second}$, a measurement that is statistically consistent with zero and that is significantly lower than the median FWHM uncertainty of the dataset (4.01 ms−1). On that account, we conclude that the 9.6 d signal manifests in RV only and is not associated with stellar activity.

Finally, we verified the stability of the 9.6 d RV signal through an apodisation test (Gregory 2016; Hara et al. 2022a). We re-ran the best model (Keplerian+circular MEP), but with the following modification in the semi-amplitude of the circular component: krv, ckrv, cexp[(tμapo)22σapo2],Mathematical equation: k_\text{rv, c} \to k_\text{rv, c} \exp\left[-\frac{(t'-\mu_\text{apo})^2}{2\sigma_\text{apo}^2}\right],(5)

that is, made it follow a Gaussian-like evolution in time. A true planetary signal is stable over time, meaning that a poorly-constrained μapo and a large σapo would pass this planetarynature test. We used a uniform prior for μapo that covers the HARPS baseline, the uninformed σapo=Ulog(1,106)dMathematical equation: $\sigma_\text{apo}=\mathcal{U}_\text{log}\left(1,10^6\right)\,\unit\day$; lastly, we constrained the orbital period prior to Pc = U (9,10) d to help inference. The μapo posterior was indistinguishable from a uniform distribution. The σapo showed a similar characteristic in log-space, with a step function that had a transition near the HARPS baseline (near 4900 d). We report σapo/Tbaseline=10.69.2+69.6Mathematical equation: $\sigma_\text{apo}/T_\text{baseline}=10.6^{+69.6}_{-9.2}$.

Figure A.6 shows the posterior of the apodised semiamplitude krv, c over time, in 1σ, 2σ and 3σ confidence intervals. We observe an amplitude that is consistent with the posterior of our best-model (1.730.23+0.24ms1Mathematical equation: $1.73^{+0.24}_{-0.23}\,\unit{\metre\per\second}$; Table A.8). The apodised krv, c begins to demonstrate variance only in the 3σ confidence interval. The upper 3σ bound reaches 3 m s−1 in the middle of the baseline. We attribute this to the scarcity of observations in the middle of the baseline, which leaves some degree of freedom to the Gaussian profile.

Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Bayesian evidence comparison between different planetary configurations and stellar-activity kernels. Planetary configurations include circular orbit (C) and Keplerian (K) components. For stellar activity, we utilised the MEP kernel and the ESP kernel with 2, 3, and 4 harmonics. The model that we further elected for analysis assumed a model with a Keplerian signal, a circular signal, and an MEP kernel (red border; ln Z = −807.3). We give the Bayesian factor ∆ ln Z of remaining models relative to this model.

4.3 The GJ 1137 multi-planetary system

Figure 6 shows the results of our selected model. Before us stands a system that exhibits complex long-term and short-term behaviour, with two planetary signals additionally permeating the RV time series: those are GJ 1137 b (Lovis et al. 2005), and our new discovery, GJ 1137 c. Figure 6 highlights a good agreement between our selected model and data, both for the entire baseline and in a dense selected interval containing approximately the last 200 d of measurements. The temporal sampling of measurements is non-trivial, but the model appears to handle it effectively, without overfitting the data. This is supported by the spread of RV and FWHM residuals - we observe no apparent trends with time, and their GLSPs show no other significant peaks (Fig. A.4). We report the following RV improvement in root mean square from raw time series to model residuals: 18.71 ms−1 to 0.62ms−1 in HARPS-pre, and 12.26ms−1 to 0.51 ms−1 in HARPS-post. When it comes to FWHM, these figures are 15.33 ms−1 to 3.36 ms−1 in HARPS-pre, and 6.03 ms−1 to 2.98 m s−1 in HARPS-post.

Figure 7 shows the phase-folded plot of the signals of GJ 1137b and GJ 1137c after subtracting the activity-related components. The RV signal of the new planet appears to be well supported by the phased data. Table A.8 displays the posteriors of our selected model. Having accounted for stellar activity and the two planets, we obtain a cycle period of 5870-380 d. We obtain the following RV and FWHM semi-amplitudes of the cycle: 14.61.2+1.3ms1Mathematical equation: $14.6^{+1.3}_{-1.2}\,\unit{\metre\per\second}$, and 18.43.9+3.5ms1Mathematical equation: $18.4^{+3.5}_{-3.9}\,\unit{\metre\per\second}$ respectively. These results again support a strong 15 yr activity in RV and activity. We obtain nearly matching cycle phases between dimensions: 0.8370.043+0.037Mathematical equation: $0.837^{+0.037}_{-0.043}$ for RV and 0.8640.070+0.054Mathematical equation: $0.864^{+0.054}_{-0.070}$ for FWHM. Our dataset assigns the following jitters to the HARPS-pre and HARPS-post datasets: 0.630.37+0.38ms1Mathematical equation: $0.63^{+0.38}_{-0.37}\,\unit{\metre\per\second}$ and 0.290.20+0.26ms1Mathematical equation: $0.29^{+0.26}_{-0.20}\,\unit{\metre\per\second}$ for RV, as well as 1.971.14+1.09ms1Mathematical equation: $1.97^{+1.09}_{-1.14}\,\unit{\metre\per\second}$ and 0.670.48+0.75ms1Mathematical equation: $0.67^{+0.75}_{-0.48}\,\unit{\metre\per\second}$ m s−1 for FWHM.

In terms of stellar activity, we successfully constrain the stellar rotation to Prot=32.31.3+1.2dMathematical equation: $P_\text{rot}=32.3^{+1.2}_{-1.3}\,\unit\day$ and derive a timescale of coherence τ=7321+33dMathematical equation: $\tau=73^{+33}_{-21}\,\unit\day$. We obtain a sine scale of 0.850.17+0.26Mathematical equation: $0.85^{+0.26}_{-0.17}$ and the following GP amplitudes: A0=2.200.42+0.53ms1Mathematical equation: $A_0=2.20^{+0.53}_{-0.42}\,\unit{\metre\per\second}$, B0=19.94.7+6.9ms1Mathematical equation: $B_0=19.9^{+6.9}_{-4.7}\,\unit{\metre\per\second}$, and A1=6.990.99+1.33ms1Mathematical equation: $A_1=6.99^{+1.33}_{-0.99}\,\unit{\metre\per\second}$. All of them are positive, which suggests that stellar-activity modulations correlate positively between RV and FWHM. Furthermore, A0B0, hinting that in GJ 1137, stellar activity manifests more prominently in the RV gradient, which in turn is likely caused by a flux effect. The positive value of B0 indicates that the active regions that we model are probably dark spots on the stellar surface.

With respect to GJ 1137b, we significantly improve the orbital solution first given by Lovis et al. (2005). We obtain a period of Pb=144.7200.029+0.029dMathematical equation: $P_\text{b}=144.720^{+0.029}_{-0.029}\,\unit\day$, a semi-0.029 amplitude of krv,b = 19.8 ± 0.4ms−1, and an eccentricity of eb=0.1180.015+0.016Mathematical equation: $e_b=0.118^{+0.016}_{-0.015}$. Our measurement for k differs by more than 3σ from Lovis et al. (2005) (18.3 ± 0.5ms−1). From our model parameters, we derived that GJ 1137 b has a semi-major axis of ab = 0.508 ± 0.05 au, a minimum mass of mb sin ib = 0.451 ± 0.012MJ (against 0.37MJ in Lovis et al. 2005), and a temporal-average instella-tion flux of Φb = 1.59 ± 0.15 Φ as defined in Méndez & Rivera-Valentín (2017). Our new detection, GJ 1137 c, has a period of Pc=9.6412(11)+(12)dMathematical equation: $P_\text{c}=9.6412^{+(12)}_{-(11)}\,\unit\day$ and a semi-amplitude of krv, c=1.730.23+0.24ms1Mathematical equation: $k_\text{rv, c}=1.73^{+0.24}_{-0.23}\,\unit{\metre\per\second}$, standing more than 6σ away from zero. For this planet, we derive a semi-major axis of ac=0.0835(8)+(8)auMathematical equation: $a_\text{c}=0.0835^{+(8)}_{-(8)}\,\unit\au$, a minimum mass of Φb=58.45.5+5.6MMathematical equation: $\Phi_\text{b}=58.4^{+5.6}_{-5.5}\,\unit\earthflux$, and an instellation flux of Φb=58.45.5+5.6MMathematical equation: $\Phi_\text{b}=58.4^{+5.6}_{-5.5}\,\unit\earthflux$.

We computed the mass-period and the mass-flux relationships of known exoplanets that have had their minimum masses measured independently and that have relative mass uncertainties smaller than one third (Akeson et al. 2013; Christiansen et al. 2025). Figures 8 and 9 display these diagrams for the period range 1-3000 d, the flux range 0.3-10 000 Φ, and the mass range 0.3-10 000 M against our posteriors of the two planets in the GJ 1137 system. GJ 1137 c lies in the bulk of the distribution of hot sub-Neptunes, while GJ 1137 b stands in the less explored regions of the parameter space in both figures. GJ 1137b is particularly interesting because it receives an Earth-like average insolation flux (1.59 ± 0.15 Φ). We attempted to assess the habitable zone (HZ) of GJ 1137 through Kopparapu et al. (2014). Their work provides the runaway and maximum-greenhouse flux limits, which we take to define a conservative HZ (cHZ); as well as recent-Venus and early-Mars flux limits, which we take to define an optimistic HZ (oHZ). We used the median values of M*, L*, Teff from Table 1 and estimated that for 5M planets, the cHZ extends in the range 0.61-1.14 au, while the oHZ extends around it in the range 0.50-1.21 au. Figure 10 compares this HZ with the orbital bands of GJ 1137 b and GJ 1137 c that come from the posteriors of their orbital parameters.

According to commonly accepted theories of planet formation, the Saturn-mass planet GJ 1137 b most likely formed beyond the ice line of its host star and migrated inward via convergent migration (Lin & Papaloizou 1986; Beaugé et al. 2003). Consequently, GJ 1137 b may host a primordial system of icy exomoons, analogous to those around the giant planets in the Solar System. Such exomoons could have settled in stable warm orbits and may be potentially habitable, ocean-like, or early-Mars-like worlds (see, discussion in e.g. Trifonov et al. 2020; Hobson et al. 2023). Although exomoon detection remains challenging with current state-of-the-art techniques, GJ 1137 b is an intriguing target for future analysis.

Thumbnail: Fig. 6 Refer to the following caption and surrounding text. Fig. 6

Time series (error bars) against our best model (solid lines) over the full dataset and over a selected temporal region. Data points come with error bars that account for the instrument and the model jitter. (a, b) SERVAL RV, where the full model (black) is split into the stellar-activity GP component (red), the LTF (pink), and the two-planet component (teal). (c, d) Associated RV residuals. (e, f) RACCOON FWHM, where the full model (black) is split into the stellar-activity GP component (red) and the LTF (pink). We implicitly correct for the HARPS-pre FWHM drift before visualising. Every component is offset by an arbitrary amount, labelled in the plots. (g, h) Associated FWHM residuals.

Thumbnail: Fig. 7 Refer to the following caption and surrounding text. Fig. 7

(a-c) Activity-corrected, phase-folded time series (error bars) against the LTF- and the planetary fit of our best model. (d-f) Residual time series after the fit. (g) Distribution of residuals (solid lines).

Thumbnail: Fig. 8 Refer to the following caption and surrounding text. Fig. 8

Mass-period diagram of exoplanets with well-established masses and periods from Akeson et al. (2013); Christiansen et al. (2025). Against those, we plot the derived posteriors of our best model. Reported uncertainties reflect the 16th and the 84th percentiles.

Thumbnail: Fig. 9 Refer to the following caption and surrounding text. Fig. 9

Mass-flux diagram of exoplanets with well-established masses and periods from Akeson et al. (2013); Christiansen et al. (2025). Against those, we plot the derived posteriors of our best model. Reported uncertainties reflect the 16th and the 84th percentiles.

Thumbnail: Fig. 10 Refer to the following caption and surrounding text. Fig. 10

Planetary system of GJ 1137 relative to the conservative habitable zone (cHZ; green) and the optimistic habitable zone (oHZ; yellow), respectively. Black orbital bands come from our 1σ posteriors of the orbital parameters of GJ 1137 b and GJ 1137 c.

4.4 Detection limits

We assessed our detection limits in the following manner. We took the period interval Ulog(1,1000) d and split it into 100 bins. In each of these bins, we took our best model and included an extra third circular-orbit planet with an Ulog orbital-period prior defined by that bin and an RV semi-amplitude of U(0,5) m s−1. We fixed all other parameters except the RV zero-order correction α0,0 and the RV jitters Ji,0. Then, for each bin, we ran an independent simpler ‘bin model’, inferred the posteriors of free parameters, and considered the ±3σ percentiles of the RV semiamplitude posterior of the sampled extra planet. To speed up inference, in all instances of Nparam bin-model parameters, we required 20Nparam live points (instead of 40Nparam on the model grid).

Figure 11 displays the +3σ percentiles of our 100 bins in the parameter space. Overall, our +3σ detection limit gradually increases from 1 ms−1 to 2ms−1 over orbital-period space. There are two apparent degradations of 2 m s−1 in orbital periods close to Prot and Prot/2, which is expected from the limitations of GP frameworks. We are sensitive to the detection of any Neptune-mass planets up to about 300 d. Typically, −3σ percentiles are indicative of subtle yet significant signals that merit attention in subsequent studies. In our case, we do not observe −3σ percentiles above 10cms−1.

Thumbnail: Fig. 11 Refer to the following caption and surrounding text. Fig. 11

Detection limits of our study. Shaded bands show the +3σ RV semi-amplitude posterior of a modelled planet of a third potential planet added to the best model, in 100 log-spaced orbital-period bins. Stellar rotation somewhat degrades detection for orbital periods near Prot (solid blue line) and Prot/2 (dashed blue line). Dotted lines show the RV semiamplitudes that an Earth-mass planet and a Neptune-mass planet would inject in the system RV.

5 Summary and conclusions

We report the results of a comprehensive analysis of precision radial velocity data and stellar activity indicators for the nearby K-dwarf star GJ 1137 (HD 93083, HIP 52521). Our study is based on 140 archival HARPS spectra that span more than 13 years and provide a long and consistent temporal baseline to investigate the presence of additional Jovian-mass analogs in this system, beyond the already known Saturn-mass exoplanet GJ 1137 b discovered by Lovis et al. (2005).

We refine the stellar parameters of GJ 1137 and report a stellar mass of 0.8360.025+0.023MMathematical equation: $0.836^{+0.023}_{-0.025}\,\unit\solarmass$. Incorporating this updated mass and the extended RV dataset collected since 2005, we revised the orbital parameters of the Saturn-mass planet GJ 1137 b. We determine a minimum mass of 0.4510.012+0.012MJMathematical equation: $0.451^{+0.012}_{-0.012}\,\unit\jupitermass$ and an orbital −0.012 period of Pb = 144.7 d, which places GJ 1137b near the optimistic habitable zone of its host star. Due to its location and mass, GJ 1137 b may be of particular interest in the context of habitability, as it could host Galilean-like exomoons with stable orbits that receive moderate insolation, making them potentially habitable.

Our objective is to discover and characterise Jovian-mass planets orbiting low-mass stellar hosts. GJ 1137 was of particular interest because, in addition to the known RV signal induced by GJ 1137 b, we identified a second long-period RV variation initially suggestive of a distant Jovian companion. However, a detailed analysis of spectroscopic activity indica-tors—including RV, FWHM, and log R′HK, among other activity indices - revealed a strong common period at Pcyc = 5870-480 d. This led us to conclude that the observed RV signal is most likely induced by a stellar magnetic activity cycle with a period of Pcyc, akin to the 11-year solar activity cycle, rather than by a Jovian-mass planet. We therefore attribute the long-period RV variation to stellar activity.

After correcting for long-term stellar activity, we uncovered a short-period planetary candidate, GJ 1137 c, with an orbital period of 9.64 d and a minimum mass of 5.120.69+0.70Mathematical equation: $5.12^{+0.70}_{-0.69}$ M, placing it firmly in the super-Earth regime. This signal is consistently recovered across various stellar noise models and kernel configurations, and we find no significant correlation with any activity indicators. FIP tests, following Hara et al. (2022b), and an apodisation of the signal associated with GJ 1137 c, following Hara et al. (2022a), further support the planetary nature of this signal.

Our findings illustrate the complexities and potential pitfalls in interpreting long-period RV trends, particularly in magnetically active stars. This is especially relevant for the detection of Jovian analogues, whose occurrence rates are estimated at ∼3-6% around FGK stars (e.g. Wittenmyer et al. 2020; Rosenthal et al. 2021), and are likely even lower for M and K dwarfs. In this context, GJ 1137 serves as a cautionary example, where longterm activity signals and sparse sampling could mimic planetary companions.

Overall, GJ 1137 emerges as a dynamically rich multi-planetary system. The presence of a close-in super-Earth, a longer-period Saturn-mass planet, and a prominent magnetic activity cycle makes it a compelling target for continued monitoring. Our results emphasise the importance of joint RV and activity diagnostics in high-precision Doppler surveys, particularly when searching for long-period, low-amplitude planetary signals.

Data availability

Full Table A.1 is fully available at the CDS via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/708/A232.

Acknowledgements

We thank the anonymous referee for their constructive comments that helped to improve the quality of this work. This research is based on spectroscopic observations made with ESO’s 3.6m telescope at La Silla Observatory under programme IDs 072.C-0488, 183.C-0972, 091.C-0936, 192.C-0852, 198.C-0836. We gratefully acknowledge Michel Mayor, Stéphane Udry, and Rodrigo Díaz, principal investigators of the HARPS programme, for providing access to the observational data used in this work. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. D.S., S.S, E.Z., D.A., V.B., and T.T. acknowledge support by the BNSF program “VIHREN-2021” project No. KP-06-DV/5. This research was in part funded by the UKRI Grant EP/X027562/1. A.K.S. acknowledges the support of a fellowship from the “la Caixa” Foundation (ID 100010434). The fellowship code is LCF/BQ/DI23/11990071. A.K.S., J.I.G.H., A.S.M., R.R., and N.N. acknowledge financial support from the Spanish Ministry of Science, Innovation and Universities (MICIU) project PID2023-149982NB-I00. N.N. acknowledges funding from Light Bridges for the Doctoral Thesis “Habitable Earth-like planets with ESPRESSO and NIRPS”, in cooperation with the Instituto de Astrofísica de Canarias, and the use of Indefeasible Computer Rights (ICR) being commissioned at the ASTRO POC project in Tenerife, Canary Islands, Spain. The ICR-ASTRONOMY used for his research was provided by Light Bridges in cooperation with Hewlett Packard Enterprise (HPE). We used the following PYTHON packages for data analysis and visualisation: EXO-STRIKER (Trifonov 2019), Matplotlib (Hunter 2007), NIEVA (Stefanov et al., in prep.), NumPy (Harris et al. 2020), PANDAS (McKinney 2010; The Pandas Development Team 2022), SCIPY (Virtanen et al. 2020), S+LEAF (Delisle et al. 2020, 2022), and ultranest (Buchner 2021). The bulk of modelling and inference was done on the Diva cluster (192 Xeon E7-4850 2.1 GHz CPUs; 4.4 TB RAM) at Instituto de Astrofísica de Canarias, Tenerife, Spain.

References

  1. Aigrain, S., Pont, F., & Zucker, S. 2012, MNRAS, 419, 3147 [NASA ADS] [CrossRef] [Google Scholar]
  2. Akeson, R. L., Chen, X., Ciardi, D., et al. 2013, PASP, 125, 989 [Google Scholar]
  3. Allard, F., Homeier, D., & Freytag, B. 2012, Philos. Trans. Roy. Soc. London Ser. A, 370, 2765 [Google Scholar]
  4. Anglada-Escudé, G., & Butler, R. P. 2012, ApJS, 200, 15 [Google Scholar]
  5. Angus, R., Morton, T., Aigrain, S., Foreman-Mackey, D., & Rajpaul, V. 2018, MNRAS, 474, 2094 [Google Scholar]
  6. Astudillo-Defru, N., Díaz, R. F., Bonfils, X., et al. 2017, A&A, 605, L11 [EDP Sciences] [Google Scholar]
  7. Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., Demleitner, M., & Andrae, R. 2021, AJ, 161, 147 [Google Scholar]
  8. Bakos, G., Noyes, R. W., Kovács, G., et al. 2004, PASP, 116, 266 [NASA ADS] [CrossRef] [Google Scholar]
  9. Balsalobre-Ruza, O., Lillo-Box, J., Silva, A. M., et al. 2025, A&A, 694, A15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Baluev, R. V. 2008, MNRAS, 385, 1279 [Google Scholar]
  11. Baluev, R. V. 2009, MNRAS, 393, 969 [Google Scholar]
  12. Baranne, A., Queloz, D., Mayor, M., et al. 1996, A&AS, 119, 373 [Google Scholar]
  13. Barragán, O., Gillen, E., Aigrain, S., et al. 2023, MNRAS, 522, 3458 [CrossRef] [Google Scholar]
  14. Beaugé, C., Ferraz-Mello, S., & Michtchenko, T. A. 2003, ApJ, 593, 1124 [Google Scholar]
  15. Bonfils, X., Delfosse, X., Udry, S., et al. 2013, A&A, 549, A109 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. Borucki, W. J., Koch, D., Basri, G., et al. 2010, Science, 327, 977 [Google Scholar]
  17. Buchner, J. 2021, J. Open Source Softw., 6, 3001 [CrossRef] [Google Scholar]
  18. Burrows, A., Halverson, S., Siegel, J. C., et al. 2024, AJ, 167, 243 [Google Scholar]
  19. Cassan, A., Kubas, D., Beaulieu, J. P., et al. 2012, Nature, 481, 167 [Google Scholar]
  20. Castelli, F., & Kurucz, R. L. 2003, in IAU Symposium, 210, Modelling of Stellar Atmospheres, eds. N. Piskunov, W. W. Weiss, & D. F. Gray, A20 [Google Scholar]
  21. Choi, J., Dotter, A., Conroy, C., et al. 2016, ApJ, 823, 102 [Google Scholar]
  22. Christiansen, J. L., McElroy, D. L., Harbut, M., et al. 2025, Planet. Sci. J., 6, 186 [Google Scholar]
  23. Collier Cameron, A., Mortier, A., Phillips, D., et al. 2019, MNRAS, 487, 1082 [Google Scholar]
  24. Costes, J. C., Watson, C. A., de Mooij, E., et al. 2021, MNRAS, 505, 830 [NASA ADS] [CrossRef] [Google Scholar]
  25. da Silva, L., Girardi, L., Pasquini, L., et al. 2006, A&A, 458, 609 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Delisle, J. B., Hara, N., & Ségransan, D. 2020, A&A, 638, A95 [EDP Sciences] [Google Scholar]
  27. Delisle, J.-B., Unger, N., Hara, N. C., & Ségransan, D. 2022, A&A, 659, A182 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Dotter, A. 2016, ApJS, 222, 8 [Google Scholar]
  29. Dressing, C. D., & Charbonneau, D. 2015, ApJ, 807, 45 [Google Scholar]
  30. Dumusque, X. 2018, A&A, 620, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Gaia Collaboration (Brown, A. G. A., et al.) 2018, A&A, 616, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Gaidos, E., Mann, A. W., Kraus, A. L., & Ireland, M. 2016, MNRAS, 457, 2877 [Google Scholar]
  33. Gaudi, B. S. 2012, ARA&A, 50, 411 [Google Scholar]
  34. Gillon, M., Triaud, A. H. M. J., Demory, B.-O., et al. 2017, Nature, 542, 456 [NASA ADS] [CrossRef] [Google Scholar]
  35. Gray, R. O., Corbally, C. J., Garrison, R. F., et al. 2006, AJ, 132, 161 [Google Scholar]
  36. Green, G. M. 2018, J. Open Source Softw., 3, 695 [Google Scholar]
  37. Gregory, P. C. 2016, MNRAS, 458, 2604 [NASA ADS] [CrossRef] [Google Scholar]
  38. Hara, N. C., Delisle, J.-B., Unger, N., & Dumusque, X. 2022a, A&A, 658, A177 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. Hara, N. C., Unger, N., Delisle, J.-B., Díaz, R. F., & Ségransan, D. 2022b, A&A, 663, A14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Harada, C. K., Dressing, C. D., Kane, S. R., & Ardestani, B. A. 2024, ApJS, 272, 30 [NASA ADS] [CrossRef] [Google Scholar]
  41. Hardegree-Ullman, K. K., Cushing, M. C., Muirhead, P. S., & Christiansen, J. L. 2019, AJ, 158, 75 [Google Scholar]
  42. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  43. Hauschildt, P. H., Allard, F., & Baron, E. 1999, ApJ, 512, 377 [Google Scholar]
  44. Haywood, R. D., Collier Cameron, A., Queloz, D., et al. 2014, MNRAS, 443, 2517 [Google Scholar]
  45. Hippke, M., & Heller, R. 2019, A&A, 623, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. Hobson, M. J., Díaz, R. F., Delfosse, X., et al. 2018, A&A, 618, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Hobson, M. J., Trifonov, T., Henning, T., et al. 2023, AJ, 166, 201 [NASA ADS] [CrossRef] [Google Scholar]
  48. Howard, A. W., Johnson, J. A., Marcy, G. W., et al. 2010, ApJ, 721, 1467 [Google Scholar]
  49. Howard, A. W., Marcy, G. W., Johnson, J. A., et al. 2010, Science, 330, 653 [Google Scholar]
  50. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  51. Husser, T. O., Wende-von Berg, S., Dreizler, S., et al. 2013, A&A, 553, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Jeffreys, H. 1946, Proc. Roy. Soc. Lond. Ser. A, 186, 453 [Google Scholar]
  53. Kochanek, C. S., Shappee, B. J., Stanek, K. Z., et al. 2017, PASP, 129, 104502 [Google Scholar]
  54. Kopparapu, R. K., Ramirez, R. M., SchottelKotte, J., et al. 2014, ApJ, 787, L29 [Google Scholar]
  55. Kunimoto, M., & Matthews, J. M. 2020, arXiv e-prints [arXiv:2004.05296] [Google Scholar]
  56. Kurucz, R. 1993, Robert Kurucz CD-ROM, 13 [Google Scholar]
  57. Lafarga, M., Ribas, I., Lovis, C., et al. 2020, A&A, 636, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Lin, D. N. C., & Papaloizou, J. 1986, ApJ, 309, 846 [Google Scholar]
  59. Lo Curto, G., Pepe, F., Avila, G., et al. 2015, The Messenger, 162, 9 [NASA ADS] [Google Scholar]
  60. Lomb, N. R. 1976, Astrophys. Space Sci., 39, 447 [Google Scholar]
  61. Lovis, C., Mayor, M., Bouchy, F., et al. 2005, A&A, 437, 1121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Marois, C., Macintosh, B., Barman, T., et al. 2008, Science, 322, 1348 [Google Scholar]
  63. Mayor, M., Pepe, F., Queloz, D., et al. 2003, The Messenger, 114, 20 [NASA ADS] [Google Scholar]
  64. McKinney, W. 2010, in Python in Science Conference, Austin, Texas, 56 [Google Scholar]
  65. Méndez, A., & Rivera-Valentín, E. G. 2017, ApJ, 837, L1 [Google Scholar]
  66. Mignon, L., Delfosse, X., Meunier, N., et al. 2025, A&A, 700, A146 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Palle, E., Biazzo, K., Bolmont, E., et al. 2025, Exp. Astron., 59, 29 [Google Scholar]
  68. Pepe, F., Mayor, M., Rupprecht, G., et al. 2002, The Messenger, 110, 9 [NASA ADS] [Google Scholar]
  69. Perdelwitz, V., Trifonov, T., Teklu, J. T., Sreenivas, K. R., & Tal-Or, L. 2024, A&A, 683, A125 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Petigura, E. A., Howard, A. W., & Marcy, G. W. 2013, PNAS, 110, 19273 [NASA ADS] [CrossRef] [Google Scholar]
  71. Pollacco, D. L., Skillen, I., Collier Cameron, A., et al. 2006, PASP, 118, 1407 [NASA ADS] [CrossRef] [Google Scholar]
  72. Queloz, D., Henry, G. W., Sivan, J. P., et al. 2001, A&A, 379, 279 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Quanz, S. P., Ottiger, M., Fontanet, E., et al. 2022, A&A, 664, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Rajpaul, V., Aigrain, S., Osborne, M. A., Reece, S., & Roberts, S. 2015, Mon. Not. R. Astron. Soc., 452, 2269 [Google Scholar]
  75. Rasmussen, C. E., & Williams, C. K. I. 2006, Gaussian Processes for Machine Learning, Adaptive Computation and Machine Learning (Cambridge, MA: MIT Press) [Google Scholar]
  76. Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2015, J. Astron. Telesc. Instrum. Syst., 1, 014003 [Google Scholar]
  77. Rosenthal, L. J., Fulton, B. J., Hirsch, L. A., et al. 2021, ApJS, 255, 8 [NASA ADS] [CrossRef] [Google Scholar]
  78. Scargle, J. D. 1982, ApJ, 263, 835 [Google Scholar]
  79. Schlafly, E. F., & Finkbeiner, D. P. 2011, ApJ, 737, 103 [Google Scholar]
  80. Shappee, B. J., Prieto, J. L., Grupe, D., et al. 2014, ApJ, 788, 48 [Google Scholar]
  81. Skilling, J. 2004, in American Institute of Physics Conference Series, 735, Bayesian Inference and Maximum Entropy Methods in Science and Engineering: 24th International Workshop on Bayesian Inference and Maximum Entropy Methods in Science and Engineering, eds. R. Fischer, R. Preuss, & U. V. Toussaint (AIP), 395 [Google Scholar]
  82. Smith, J. C., Stumpe, M. C., Van Cleve, J. E., et al. 2012, PASP, 124, 1000 [Google Scholar]
  83. Stassun, K. G., Collins, K. A., & Gaudi, B. S. 2017, AJ, 153, 136 [Google Scholar]
  84. Stassun, K. G., Oelkers, R. J., Paegert, M., et al. 2019, AJ, 158, 138 [Google Scholar]
  85. Stefanov, A. K., Suárez Mascareño, A., González Hernández, J. I., et al. 2025, A&A, 695, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  86. Stumpe, M. C., Smith, J. C., Van Cleve, J. E., et al. 2012, PASP, 124, 985 [Google Scholar]
  87. Stumpe, M. C., Smith, J. C., Catanzarite, J. H., et al. 2014, PASP, 126, 100 [Google Scholar]
  88. Suárez Mascareño, A., Rebolo, R., González Hernández, J. I., & Esposito, M. 2015, MNRAS, 452, 2745 [Google Scholar]
  89. Suárez Mascareño, A., Faria, J. P., Figueira, P., et al. 2020, A&A, 639, A77 [Google Scholar]
  90. Suárez Mascareño, A., González-Álvarez, E., Zapatero Osorio, M. R., et al. 2023, A&A, 670, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. The Pandas Development Team. 2022, Pandas-Dev/Pandas: Pandas, Zenodo [Google Scholar]
  92. Trifonov, T. 2019, The Exo-Striker: Transit and radial velocity interactive fitting tool for orbital analysis and N-body simulations [Google Scholar]
  93. Trifonov, T., Tal-Or, L., Zechmeister, M., et al. 2020, A&A, 636, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  94. van Leeuwen, F. 2007, A&A, 474, 653 [CrossRef] [EDP Sciences] [Google Scholar]
  95. VanderPlas, J. T., & Ivezic, Z. 2015, ApJ, 812, 18 [NASA ADS] [CrossRef] [Google Scholar]
  96. Vines, J. I., & Jenkins, J. S. 2022, MNRAS, 513, 2719 [NASA ADS] [CrossRef] [Google Scholar]
  97. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Methods, 17, 261 [Google Scholar]
  98. Wittenmyer, R. A., Wang, S., Horner, J., et al. 2020, MNRAS, 492, 377 [NASA ADS] [CrossRef] [Google Scholar]
  99. Winecki, D., & Kochanek, C. S. 2024, ApJ, 971, 61 [NASA ADS] [CrossRef] [Google Scholar]
  100. Zechmeister, M., & Kürster, M. 2009, A&A, 496, 577 [CrossRef] [EDP Sciences] [Google Scholar]
  101. Zechmeister, M., Kürster, M., Endl, M., et al. 2013, A&A, 552, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  102. Zechmeister, M., Reiners, A., Amado, P. J., et al. 2018, A&A, 609, A12 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  103. Zechmeister, M., Dreizler, S., Ribas, I., et al. 2019, A&A, 627, A49 [NASA ADS] [EDP Sciences] [Google Scholar]

1

Up-to-date statistics available at https://exoplanetarchive.ipac.caltech.edu/.

2

See section 6.7.1 of the TESS Instrument Handbook and Data Release Notes.

Appendix A Supplementary material

Appendix A:1 Preliminary spectroscopic analysis with the Exo-Striker and MLP

We performed a periodogram analysis of the SERVAL radial velocities and activity index time series using the maximum likelihood periodogram (MLP; Baluev 2008, 2009), which is conceptually similar to the widely used generalised Lomb-Scargle periodogram (GLS; Zechmeister & Kürster 2009), but is better suited for the analysis of multi-telescope data, or datasets with known systematics such as the pre- and post-HARPS fibre intervention data. The MLP allows for multiple data subsets, each with its own additive RV offset and jitter term (Zechmeister et al. 2019), while it optimises the log-likelihood ln L at each test frequency in a period search.

Fig. A.2 presents the results of our MLP analysis of the extracted HARPS spectroscopic data, with individual diagnostics shown in the sub-panels as labelled. The horizontal lines mark the adopted false alarm probability (FAP) thresholds of 10%, 1%, and 0.1%; we consider signals exceeding the 0.1% threshold significant. Vertical dashed lines indicate the best-fit orbital periods of GJ 1137 b and the long-period signal, which we initially suspected to originate from a massive exoplanet with orbital and physical properties comparable to those of Jupiter.

Several HARPS activity indicators, including: BIS, contrast, FWHM, dLW, and log R′HK, however, exhibit significant low-frequency periodicities that closely resemble those found in the RV residuals after modelling the GJ 1137 b signal. Given our local example, the Solar System, both Jupiter and the Sun’s 11-year magnetic activity cycle induce similar RV signals if the Sun were observed as a star. Therefore, it is not unlikely that GJ 1137 hosts a massive planet and/or exhibits long-term stellar magnetic cycles that we could detect in spectroscopic datasets. The disentanglement of these phenomena in RV is often challenging and raises caution regarding a possible stellar activity origin for the long-period trend.

Table A.1

RVs and activity indices data for GJ 1137.

Thumbnail: Fig. A.1 Refer to the following caption and surrounding text. Fig. A.1

ASAS-SN photometry, periodograms, and window functions. Panels (a) and (b) show the full and Moon-masked photometric datasets, respectively, together with the spline fits used to remove long-term trends. Panels (c) and (d) present the corresponding GLS periodograms, while panels (e) and (f) display their window functions. Our data curation alters the window function, as observations discarded due to Moon contamination occur at regular intervals. In the full dataset, strong peaks appear near ∼29.9 and ∼ 14.8 days, corresponding to the lunar synodic period and its first harmonic. These peaks disappear completely in the Moon-masked data. A peak at 491.2 days (full dataset) and 494.6 days (Moon-masked dataset) lies close to the 0.1% FAP level. However, this signal does not persist coherently across multiple observing windows, and given the overall data quality, we do not consider it further in our analysis.

Thumbnail: Fig. A.2 Refer to the following caption and surrounding text. Fig. A.2

Time series data and MLS power spectrum of the HARPS RVs and stellar activity indicators of GJ 1137. The left sub-panels show the individual RVs, the residuals of fitting the GJ 1137 b signal, and the two-Keplerian fit residuals including the long-period RV component. The rest are the HARPS activity indices as labelled in the panels. The vertical magenta line notates the HARPS-fibre change from which point we treat the pre- and post-HARPS data with a separate offset parameter. The right panel shows the individual MLP power spectrum in terms of ∆lnL as indicated in the sub-panels. The horizontal lines in the MLP periodograms show the FAP levels of 10%, 1%, and 0.1%, the latter of which is considered significant. The magenta vertical line indicates the orbital period of GJ 1137 b, whereas the red line indicates the best-fit long-term period seen in the one-planet Keplerian model residuals.

Table A.2

Parameter posteriors of our best model that included only an LTF and GJ 1137 b.

Thumbnail: Fig. A.3 Refer to the following caption and surrounding text. Fig. A.3

Violin plots of models in our model grid (Fig. 5), that show the posteriors of the GP hyperparameters as well as the RV jitter.

Table A.3

Period posteriors of GJ 1137 c in units of d from the gridsearch in Fig. 5.

Table A.4

Phase posteriors of GJ 1137 c from the grid-search in Fig. 5.

Table A.5

RV semi-amplitude posteriors of GJ 1137 c in units of m s−1 from the grid-search in Fig. 5.

Table A.6

Eccentricity posteriors of GJ 1137 c from the grid-search in Fig. 5.

Table A.7

Argument-of-periastron posteriors of GJ 1137 c in units of rad from the grid-search in Fig. 5.

Thumbnail: Fig. A.4 Refer to the following caption and surrounding text. Fig. A.4

GLSPs of the residual time series of our best model, accounting for the model jitter, for: (a) SERVAL RVs, (b) RACCOON FWHMs.

Thumbnail: Fig. A.5 Refer to the following caption and surrounding text. Fig. A.5

FWHM amplitude-period distribution in our search for a common RV&FWHM signal in the data. The shaded contours highlight 25%, 50%, and 75% of the enclosed probability mass. The median FWHM error is marked with a solid dashed line.

Thumbnail: Fig. A.6 Refer to the following caption and surrounding text. Fig. A.6

Stability of the RV signal induced by GJ 1137 c through an apodisation test. Centre: the median value (solid black line), as well as the 1σ, 2σ, and 3σ confidence intervals (shaded bands) of the apodised signal over time. Top: timestamps of measurements, with markers following Fig. 3.

Table A.8

Parameter posteriors of our best model, against the results published by Lovis et al. (2005).

All Tables

Table 1

Stellar parameters of GJ 1137 and their 1σ uncertainties.

Table A.1

RVs and activity indices data for GJ 1137.

Table A.2

Parameter posteriors of our best model that included only an LTF and GJ 1137 b.

Table A.3

Period posteriors of GJ 1137 c in units of d from the gridsearch in Fig. 5.

Table A.4

Phase posteriors of GJ 1137 c from the grid-search in Fig. 5.

Table A.5

RV semi-amplitude posteriors of GJ 1137 c in units of m s−1 from the grid-search in Fig. 5.

Table A.6

Eccentricity posteriors of GJ 1137 c from the grid-search in Fig. 5.

Table A.7

Argument-of-periastron posteriors of GJ 1137 c in units of rad from the grid-search in Fig. 5.

Table A.8

Parameter posteriors of our best model, against the results published by Lovis et al. (2005).

All Figures

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

SED fit for GJ 1137. The best-fitting model, BT-Settl, is displayed in black, with normalised residuals shown below. Blue points represent flux values from photometry, while purple diamonds show the flux from synthetic photometry in the same passband.

In the text
Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

HARPS RV time series for GJ 1137. The data are divided into two subsets due to known instrumental RV offsets introduced by the HARPS optical fibre upgrade in May 2015. Purple circles represent RVs obtained before the fibre change, while green squares correspond to data taken after the upgrade. The observed RV variations are primarily driven by the signal of the known Saturn-mass exoplanet GJ 1137 b, with a period of P = 144.7 d. In addition, the RV data reveals a second, long-term signal from a potential Jovian planet with a minimum mass of approximately 1.3 MJ and a period P = 5640 ± 240 d. Subsequent analysis attributes this long-term signal to stellar magnetic activity.

In the text
Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Raw time series of mean-subtracted: (a) RV, (c) FWHM, (e) R′HK. Measurements are marked depending on the data source: blue circles for HARPS-pre, and green squares for HARPS-post. (b,d,f) Associated wide-period GLSPs of RV and activity indicators. Three false alarm probability (FAP) levels: 10%, 1%, and 0.1%, split GLSP ordinates in bands of different colour. We highlight the three most prominent peaks in each GLSP.

In the text
Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Time series (markers) against our LTF model (black), over the whole dataset baseline in: (a) RV, (c) FWHM, (e) R′HK. We implicitly correct for the HARPS-pre FWHM drift before visualising. (b,d,f) Associated GLSPs of the residual time series, accounting for the model jitter.

In the text
Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Bayesian evidence comparison between different planetary configurations and stellar-activity kernels. Planetary configurations include circular orbit (C) and Keplerian (K) components. For stellar activity, we utilised the MEP kernel and the ESP kernel with 2, 3, and 4 harmonics. The model that we further elected for analysis assumed a model with a Keplerian signal, a circular signal, and an MEP kernel (red border; ln Z = −807.3). We give the Bayesian factor ∆ ln Z of remaining models relative to this model.

In the text
Thumbnail: Fig. 6 Refer to the following caption and surrounding text. Fig. 6

Time series (error bars) against our best model (solid lines) over the full dataset and over a selected temporal region. Data points come with error bars that account for the instrument and the model jitter. (a, b) SERVAL RV, where the full model (black) is split into the stellar-activity GP component (red), the LTF (pink), and the two-planet component (teal). (c, d) Associated RV residuals. (e, f) RACCOON FWHM, where the full model (black) is split into the stellar-activity GP component (red) and the LTF (pink). We implicitly correct for the HARPS-pre FWHM drift before visualising. Every component is offset by an arbitrary amount, labelled in the plots. (g, h) Associated FWHM residuals.

In the text
Thumbnail: Fig. 7 Refer to the following caption and surrounding text. Fig. 7

(a-c) Activity-corrected, phase-folded time series (error bars) against the LTF- and the planetary fit of our best model. (d-f) Residual time series after the fit. (g) Distribution of residuals (solid lines).

In the text
Thumbnail: Fig. 8 Refer to the following caption and surrounding text. Fig. 8

Mass-period diagram of exoplanets with well-established masses and periods from Akeson et al. (2013); Christiansen et al. (2025). Against those, we plot the derived posteriors of our best model. Reported uncertainties reflect the 16th and the 84th percentiles.

In the text
Thumbnail: Fig. 9 Refer to the following caption and surrounding text. Fig. 9

Mass-flux diagram of exoplanets with well-established masses and periods from Akeson et al. (2013); Christiansen et al. (2025). Against those, we plot the derived posteriors of our best model. Reported uncertainties reflect the 16th and the 84th percentiles.

In the text
Thumbnail: Fig. 10 Refer to the following caption and surrounding text. Fig. 10

Planetary system of GJ 1137 relative to the conservative habitable zone (cHZ; green) and the optimistic habitable zone (oHZ; yellow), respectively. Black orbital bands come from our 1σ posteriors of the orbital parameters of GJ 1137 b and GJ 1137 c.

In the text
Thumbnail: Fig. 11 Refer to the following caption and surrounding text. Fig. 11

Detection limits of our study. Shaded bands show the +3σ RV semi-amplitude posterior of a modelled planet of a third potential planet added to the best model, in 100 log-spaced orbital-period bins. Stellar rotation somewhat degrades detection for orbital periods near Prot (solid blue line) and Prot/2 (dashed blue line). Dotted lines show the RV semiamplitudes that an Earth-mass planet and a Neptune-mass planet would inject in the system RV.

In the text
Thumbnail: Fig. A.1 Refer to the following caption and surrounding text. Fig. A.1

ASAS-SN photometry, periodograms, and window functions. Panels (a) and (b) show the full and Moon-masked photometric datasets, respectively, together with the spline fits used to remove long-term trends. Panels (c) and (d) present the corresponding GLS periodograms, while panels (e) and (f) display their window functions. Our data curation alters the window function, as observations discarded due to Moon contamination occur at regular intervals. In the full dataset, strong peaks appear near ∼29.9 and ∼ 14.8 days, corresponding to the lunar synodic period and its first harmonic. These peaks disappear completely in the Moon-masked data. A peak at 491.2 days (full dataset) and 494.6 days (Moon-masked dataset) lies close to the 0.1% FAP level. However, this signal does not persist coherently across multiple observing windows, and given the overall data quality, we do not consider it further in our analysis.

In the text
Thumbnail: Fig. A.2 Refer to the following caption and surrounding text. Fig. A.2

Time series data and MLS power spectrum of the HARPS RVs and stellar activity indicators of GJ 1137. The left sub-panels show the individual RVs, the residuals of fitting the GJ 1137 b signal, and the two-Keplerian fit residuals including the long-period RV component. The rest are the HARPS activity indices as labelled in the panels. The vertical magenta line notates the HARPS-fibre change from which point we treat the pre- and post-HARPS data with a separate offset parameter. The right panel shows the individual MLP power spectrum in terms of ∆lnL as indicated in the sub-panels. The horizontal lines in the MLP periodograms show the FAP levels of 10%, 1%, and 0.1%, the latter of which is considered significant. The magenta vertical line indicates the orbital period of GJ 1137 b, whereas the red line indicates the best-fit long-term period seen in the one-planet Keplerian model residuals.

In the text
Thumbnail: Fig. A.3 Refer to the following caption and surrounding text. Fig. A.3

Violin plots of models in our model grid (Fig. 5), that show the posteriors of the GP hyperparameters as well as the RV jitter.

In the text
Thumbnail: Fig. A.4 Refer to the following caption and surrounding text. Fig. A.4

GLSPs of the residual time series of our best model, accounting for the model jitter, for: (a) SERVAL RVs, (b) RACCOON FWHMs.

In the text
Thumbnail: Fig. A.5 Refer to the following caption and surrounding text. Fig. A.5

FWHM amplitude-period distribution in our search for a common RV&FWHM signal in the data. The shaded contours highlight 25%, 50%, and 75% of the enclosed probability mass. The median FWHM error is marked with a solid dashed line.

In the text
Thumbnail: Fig. A.6 Refer to the following caption and surrounding text. Fig. A.6

Stability of the RV signal induced by GJ 1137 c through an apodisation test. Centre: the median value (solid black line), as well as the 1σ, 2σ, and 3σ confidence intervals (shaded bands) of the apodised signal over time. Top: timestamps of measurements, with markers following Fig. 3.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.