# Company for the Ultra-high Density, Ultra-short Period Sub-Earth GJ 367 b: Discovery of Two Additional Low-mass Planets at 11.5 and 34 Days

Elisa Goffo, Davide Gandolfi, Jo Ann Egger, Alexander J. Mustill, Simon H.

Albrecht, Teruyuki Hirano, Oleg Kochukhov, Nicola Astudillo-Defru, Oscar Barragan, Luisa M. Serrano, et al.

<figure><img src="image_1.png" /></figure>

► To cite this version:

Elisa Goffo, Davide Gandolfi, Jo Ann Egger, Alexander J. Mustill, Simon H. Albrecht, et al. Company for the Ultra-high Density, Ultra-short Period Sub-Earth GJ 367 b: Discovery of Two Additional Low-mass Planets at 11.5 and 34 Days. apjl, 2023, 955 (1), pp.L3. 10.3847/2041-8213/ace0c7. hal-04260188

HAL Id: hal-04260188

https://hal.science/hal-04260188v1

Submitted on 30 Oct 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_3.png" /></figure>

Company for the Ultra-high Density, Ultra-short Period Sub-Earth GJ 367 b: Discovery of Two Additional Low-mass Planets at 11.5 and 34 Days*

Elisa Goffo<sup>1,2 ID</sup>, Davide Gandolfi<sup>1 ID</sup>, Jo Ann Egger<sup>3 ID</sup>, Alexander J. Mustill<sup>4,5 ID</sup>, Simon H. Albrecht<sup>6 ID</sup>, Teruyuki Hirano<sup>7,8 ID</sup>, Oleg Kochukhov<sup>9 ID</sup>, Nicola Astudillo-Defru<sup>10 ID</sup>, Oscar Barragan<sup>11 ID</sup>, Luisa M. Serrano<sup>1 ID</sup>, Artie P. Hatzes<sup>2 ID</sup>, Yann Alibert<sup>3,12 ID</sup>, Eike Guenther<sup>2 ID</sup>, Fei Dai<sup>13,14,34 ID</sup>, Kristine W. F. Lam<sup>15 ID</sup>, Szilárd Csizmadia<sup>15 ID</sup>, Alexis M. S. Smith<sup>15 ID</sup>, Luca Fossati<sup>16 ID</sup>, Rafael Luque<sup>17 ID</sup>, Florian Rodler<sup>18 ID</sup>, Mark L. Winther<sup>6 ID</sup>, Jakob L. Rørsted<sup>6 ID</sup>, Javier Alarcon<sup>18 ID</sup>, Xavier Bonfils<sup>19 ID</sup>, William D. Cochran<sup>20 ID</sup>, Hans J. Deeg<sup>21,22 ID</sup>, Jon M. Jenkins<sup>23 ID</sup>, Judith Korth<sup>4 ID</sup>, John H. Livingston<sup>7,8,24 ID</sup>, Annabella Meech<sup>25 ID</sup>, Felipe Murgas<sup>21,22 ID</sup>, Jaume Orell-Miquel<sup>21,22 ID</sup>, Hannah L. M. Osborne<sup>26 ID</sup>, Enric Palle<sup>21,22 ID</sup>, Carina M. Persson<sup>27 ID</sup>, Seth Redfield<sup>28 ID</sup>, George R. Ricker<sup>29 ID</sup>, Sara Seager<sup>30,31,32 ID</sup>, Roland Vanderspek<sup>31 ID</sup>, Vincent Van Eylen<sup>26 ID</sup>, and Joshua N. Winn<sup>33 ID</sup>

Received 2023 February 16; revised 2023 June 15; accepted 2023 June 21; published 2023 September 14

Abstract

GJ 367 is a bright ($V \approx 10.2$) M1 V star that has been recently found to host a transiting ultra-short period sub-Earth on a 7.7 hr orbit. With the aim of improving the planetary mass and radius and unveiling the inner architecture of the system, we performed an intensive radial velocity follow-up campaign with the HARPS spectrograph—collecting 371 high-precision measurements over a baseline of nearly 3 yr—and combined our Doppler measurements with new TESS observations from sectors 35 and 36. We found that GJ 367 b has a mass of $M_b = 0.633 \pm 0.050 M_⊕$ and a radius of $R_b = 0.699 \pm 0.024 R_⊕$, corresponding to precisions of 8% and 3.4%, respectively. This implies a planetary bulk density of $\rho_b = 10.2 \pm 1.3 \text{ g cm}^{-3}$, i.e., 85% higher than Earth's density. We revealed the presence of two additional non-transiting low-mass companions with orbital periods of $\sim 11.5$ and 34 days and minimum masses of $M_c \sin i_c = 4.13 \pm 0.36 M_⊕$ and $M_d \sin i_d = 6.03 \pm 0.49 M_⊕$,

<figure><img src="image_4.png" /></figure>

<footnote n="1">¹ Dipartimento di Fisica, Università degli Studi di Torino, via Pietro Giuria 1, I-10125, Torino, Italy; elisa.goffo@unito.it, elisa@tls-tautenburg.de</footnote>

<footnote n="2">² Thüringer Landesfernwarte Tautenburg, Sternwarte 5, D-07778 Tautenburg, Germany</footnote>

<footnote n="3">³ Physikalisches Institut, University of Bern, Gesellschaftsstrasse 6, 3012 Bern, Switzerland</footnote>

<footnote n="4">⁴ Lund Observatory, Division of Astrophysics, Department of Physics, Lund University, Box 43, SE-221 00 Lund, Sweden</footnote>

<footnote n="5">⁵ Lund Observatory, Department of Astronomy and Theoretical Physics, Lund University, Box 43, SE-221 00 Lund, Sweden</footnote>

<footnote n="6">⁶ Stellar Astrophysics Centre, Department of Physics and Astronomy, Aarhus University, Ny Munkegade 120, DK-8000 Aarhus C, Denmark</footnote>

<footnote n="7">⁷ Astrobiology Center, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan</footnote>

<footnote n="8">⁸ National Astronomical Observatory of Japan, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan</footnote>

<footnote n="9">⁹ Department of Physics and Astronomy, Uppsala University, Box 516, SE-75120 Uppsala, Sweden</footnote>

<footnote n="10">¹⁰ Departamento de Matemática y Física Aplicadas, Universidad Católica de la Santísima Concepción, Alonso de Rivera 2850, Concepción, Chile</footnote>

<footnote n="11">¹¹ Sub-department of Astrophysics, Department of Physics, University of Oxford, Oxford, OX1 3RH, UK</footnote>

<footnote n="12">¹² Center for Space and Habitability, University of Bern, Gesellsschaftsstr. 6 CH3012, Bern, Switzerland</footnote>

<footnote n="13">¹³ Division of Geological and Planetary Sciences, 1200 E California Boulevard, Pasadena, CA, 91125, USA</footnote>

<footnote n="14">¹⁴ Department of Astronomy, California Institute of Technology, Pasadena, CA 91125, USA</footnote>

<footnote n="15">¹⁵ Institute of Planetary Research, German Aerospace Center (DLR), Rutherfordstrasse 2, D-12489 Berlin, Germany</footnote>

<footnote n="16">¹⁶ Space Research Institute, Austrian Academy of Sciences, Schmiedlstrasse 6, A-8042, Graz, Austria</footnote>

<footnote n="17">¹⁷ Department of Astronomy and Astrophysics, University of Chicago, Chicago, IL 60637, USA</footnote>

<footnote n="18">¹⁸ European Southern Observatory, Alonso de Cordova 3107, Vitacura, Santiago de Chile, Chile</footnote>

<footnote n="19">¹⁹ Université Grenoble Alpes, CNRS, IPAG, F-38000 Grenoble, France</footnote>

<footnote n="20">²⁰ McDonald Observatory and Center for Planetary Systems Habitability, The University of Texas, Austin Texas USA</footnote>

<footnote n="21">²¹ Instituto de Astrofísica de Canarias (IAC), E-38205 La Laguna, Tenerife, Spain</footnote>

<footnote n="22">²² Departamento de Astrofísica, Universidad de La Laguna (ULL), E-38206 La Laguna, Tenerife, Spain</footnote>

<footnote n="23">²³ NASA Ames Research Center, Moffett Field, CA 94035, USA</footnote>

<footnote n="24">²⁴ Department of Astronomical Science, School of Physical Sciences, The Graduate University for Advanced Studies (SOKENDAI), 2-21-1, Osawa, Mitaka, Tokyo, 181-8588, Japan</footnote>

<footnote n="25">²⁵ Department of Physics, University of Oxford, Keble Road, Oxford, OX1 3RH, UK</footnote>

<footnote n="26">²⁶ Mullard Space Science Laboratory, University College London, Holmbury St Mary, Dorking, Surrey RH5 6NT, UK</footnote>

<footnote n="27">²⁷ Department of Space, Earth and Environment, Chalmers University of Technology, Onsala Space Observatory, SE-439 92 Onsala, Sweden</footnote>

<footnote n="28">²⁸ Astronomy Department and Van Vleck Observatory, Wesleyan University, Middletown, CT 06459, USA</footnote>

<footnote n="29">²⁹ MIT Kavli Institute for Astrophysics and Space Research & MIT Physics Department, USA</footnote>

<footnote n="30">³⁰ Department of Earth, Atmospheric, and Planetary Sciences, Massachusetts Institute of Technology, Cambridge, MA 02139, USA</footnote>

<footnote n="31">³¹ Department of Physics and Kavli Institute for Astrophysics and Space Research, Massachusetts Institute of Technology, Cambridge, MA 02139, USA</footnote>

<footnote n="32">³² Department of Aeronautics and Astronautics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA</footnote>

<footnote n="33">³³ Department of Astrophysical Sciences, Princeton University, Princeton, NJ 08544, USA</footnote>

<footnote>* Based on observations made with the ESO-3.6 m telescope at La Silla Observatory under programs 1102.C-0923 and 106.21TJ.001.</footnote>

<footnote n="34">³⁴ NASA Sagan Fellow.</footnote>
respectively, which lie close to the 3:1 mean motion commensurability. GJ 367 b joins the small class of high-density planets, namely the class of super-Mercuries, being the densest ultra-short period small planet known to date. Thanks to our precise mass and radius estimates, we explored the potential internal composition and structure of GJ 367 b, and found that it is expected to have an iron core with a mass fraction of $0.91^{+0.07}_{-0.23}$. How this iron core is formed and how such a high density is reached is still not clear, and we discuss the possible pathways of formation of such a small ultra-dense planet.

*Unified Astronomy Thesaurus concepts: Exoplanet detection methods (489); Exoplanet systems (484); Exoplanets (498); Extrasolar rocky planets (511)*

*Supporting material: machine-readable table*

# 1. Introduction

Close-in planets with orbital periods of a few days challenge planet formation and evolution theories and play a key role in the architecture of exoplanetary systems (Winn & Fabrycky 2015; Zhu & Dong 2021). To date, about 132 ultra-short period (USP) planets, namely planets with orbital periods shorter than 1 day (Sahu et al. 2006; Sanchis-Ojeda et al. 2014), have been validated, and only 36 of these were confirmed and have measured radii and masses.<fnref n="35" /> USPs are preferred targets for transit and radial velocity (RV) planet search surveys, as the transit probability is higher—it scales as $P_{\text{orb}}^{-2/3}$—and the Doppler reflex motion is larger—it scales as $P_{\text{orb}}^{-1/3}$. In addition, their orbital period is typically 1 order of magnitude shorter than the rotation period of the star, allowing one to disentangle bona fide planetary signals from stellar activity (Hatzes et al. 2011; Hatzes 2019; Winn et al. 2018).

Sanchis-Ojeda et al. (2014) found that the occurrence rate of rocky USP planets seems to depend on the spectral type of the host star, being 0.15 ± 0.05% for F dwarfs, 0.51% ± 0.07% for G dwarfs, and 1.10% ± 0.40% for M dwarfs. In this context, low-mass stars, such as M dwarfs, are particularly suitable to search for close-in terrestrial planets. Given the relatively small stellar radius and mass, a planet transiting an M-dwarf star induces both a deeper transit and a larger RV signal, increasing its detection probability (Cifuentes et al. 2020).

The formation process of USP planets is still not fully understood, and different scenarios have been proposed to explain their short-period orbits: dynamical interactions in multiplanet systems (Schlaufman et al. 2010); low-eccentricity migration due to secular planet–planet interactions (Pu & Lai 2019); high-eccentricity migration due to secular dynamical chaos (Petrovich et al. 2019); tidal orbital decay of USP planets formed in situ (Lee & Chiang 2017); and obliquity-driven tidal migration (Millholland & Spalding 2020). Intensive follow-up observations of systems hosting USP planets can help us to understand the formation and evolution mechanisms of short-period objects and other phenomena related to star–planet interactions (Serrano et al. 2022).

Sanchis-Ojeda et al. (2014), Adams et al. (2017), and Winn et al. (2018) found that most USP planets have nearby planetary companions. In multiplanet systems, USP planets show wider-than-usual period ratios with their nearest companion, and appear to have larger mutual inclinations than planets on outer orbits (Rodriguez et al. 2018; Dai et al. 2018). These observations suggest that USP planets experienced a change in their orbital parameters, such as inclination increase and orbital shrinkage, suggesting the presence of long-period planets.

During its primary mission, NASA’s Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2015) discovered a shallow (~300 ppm) transit event repeating every ~7.7 hr (~0.32 days) and associated to a USP small planet candidate orbiting the bright ($V \approx 10.2$), nearby ($d \approx 9.4$ pc), M1 V star GJ 367. Lam et al. (2021) recently confirmed GJ 367 b as a bona fide USP sub-Earth with a radius of $R_b = 0.718 \pm 0.054 R_⊕$ and a mass of $M_b = 0.546 \pm 0.078 M_⊕$.

As the majority of well-characterized USP systems are consistent with having additional planetary companions (Dai et al. 2021), it is quite realistic to believe that the GJ 367 hosts more than one planet. As part of the RV follow-up program carried out by the KESPRINT consortium,<fnref n="36" /> we here present the results of an intensive RV campaign conducted with the HARPS spectrograph to refine the mass determination of the transiting USP planet and search for external planetary companions, while probing the architecture of the GJ 367 planetary system.

The paper is organized as follows: we provide a summary of the TESS data and describe our HARPS spectroscopic follow-up in Section 2. Stellar fundamental parameters are presented in Section 3. We report on the RV and transit analysis in Sections 4 and 5, along with the frequency analysis of our HARPS time series. Discussion and conclusions are given in Sections 6 and 7, respectively.

# 2. Observations

## 2.1. TESS Photometry

TESS observed GJ 367 in Sector 9 as part of its primary mission, from 2019 February 28 to 2019 March 26, with CCD 1 of camera 3 at a cadence of 2 minutes. These observations have been presented in Lam et al. (2021). About 2 yr later, TESS re-observed GJ 367 as part of its extended mission in Sectors 35 and 36, from 2021 February 9 to 2021 April 2, with CCD 1 and 2 of camera 3 at a higher cadence of 20 s as well as at 2 minutes. The photometric data were processed by the TESS data processing pipeline developed by the Science Processing Operations Center (SPOC; Jenkins et al. 2016). The SPOC pipeline uses Simple Aperture Photometry (SAP) to generate stellar light curves, where common instrumental systematics are removed via the Presearch Data Conditioning (PDCSAP) algorithm developed for the Kepler space mission (Stumpe et al. 2012, 2014; Smith et al. 2012).

We retrieved TESS Sector 9, 35, and 36 data from from the Mikulski Archive for Space Telescopes (MAST)<fnref n="37" /> and performed our data analyses using the PDCSAP light curve. We ran the Détection Spécialisée de Transits (DST;

<footnote n="35">³⁵ See https://exoplanetarchive.ipac.caltech.edu, as of May 2023.</footnote>

<footnote n="36">³⁶ https://kesprint.science/</footnote>

<footnote n="37">³⁷ https://mast.stsci.edu/portal/Mashup/Clients/Mast/Portal.html</footnote>
Cabrera et al. 2012) algorithm to search for additional transit signals and found no significant detection besides the 7.7 h signal associated to GJ 367 b, suggesting that there are no other transiting planets in the system observed in TESS Sector 9, 35, and 36, consistent with the SPOC multitransiting planet search.

## 2.2. HARPS High-precision Doppler Follow-up

GJ 367 was observed with the High Accuracy Radial velocity Planet Searcher (HARPS) spectrograph (Mayor et al. 2003), mounted at the ESO-3.6 m telescope of La Silla Observatory in Chile. We collected 295 high-resolution ($R \approx 115,000$, $\lambda \in 378-691$ nm) spectra between 2020 November 9 and 2022 April 18 (UT), as part of our large observing program 106.21TJ.001 (PI: Gandolfi) to follow-up TESS transiting planets. When added to the 77 HARPS spectra published in Lam et al. (2021), our data includes 371 HARPS spectra.

The exposure time varied between 600 and 1200 s, depending on weather conditions and observing schedule constraints, leading to a signal-to-noise ratio (S/N) per pixel at 550 nm ranging between 20 and 90, with a median of ~55. We used the second fiber of the instrument to simultaneously observe a Fabry-Perot interferometer and trace possible nightly instrumental drifts (Wildi et al. 2010, 2011). The HARPS data were reduced using the Data Reduction Software (DRS; Lovis & Pepe 2007) available at the telescope. The RV measurements, as well as the $H\alpha$, $H\beta$, $H\gamma$, Na D activity indicators, $\log R'_{HK}$, the differential line width (DLW), and the chromaticity index (CRX), were extracted using the codes NAIRA (Astudillo-Defru et al. 2017b) and SERVAL (Zechmeister et al. 2017). NAIRA and SERVAL feature template matching algorithms that are suitable to derive precise RVs for M-dwarf stars, when compared to the cross-correlation function technique implemented in the DRS. We tested both the NAIRA and SERVAL RV time series and found no significant difference in the fitted parameters. While we have no reason to prefer one code over the other, we used the RV data extracted with NAIRA for the analyses described in the following sections.

Table A1 lists the HARPS RVs, including those previously reported in Lam et al. (2021), along with the activity indicators and line profile variation diagnostics extracted with NAIRA and SERVAL. Time stamps are given in Barycentric Julian Date in the Barycentric Dynamical Time (BJD$_{TDB}$).

## 2.3. HARPS Spectropolarimetric Observations

With the aim of measuring the magnetic field of GJ 367, we performed a single circular polarization observation with the HARPSpol polarimeter (Piskunov et al. 2011; Snik et al. 2011) on 2022 November 16 (UT), as part of the our ESO HARPS program 1102.C-0923 (PI: Gandolfi). We used an exposure time of 3600 s split in four $T_{exp} = 900$ s subexposures obtained with different configuration of polarization optics to ensure cancellation of the spurious instrumental signals (see Donati et al. 1997; Bagnulo et al. 2009). The data reduction was carried out with the REDUCE code (Piskunov & Valenti 2002) following the steps described in Rusomarov et al. (2013). The resulting Stokes I (intensity) and Stokes V (circular polarization) spectra cover approximately the same wavelength interval as the usual HARPS observations at a slightly reduced resolving power ($R \approx 110,000$). We also derived a diagnostic null spectrum (e.g., Bagnulo et al. 2009), which is useful for

<table><caption>**Table 1**
Fundamental Parameters of GJ 367</caption><thead><tr><th>Parame2ter</th><th>Value</th><th>Reference</th></tr></thead><tbody><tr><td>Name</td><td>GJ 367<br>TOI-731<br>TIC 34068865</td><td></td></tr><tr><td>R.A. (J2000)</td><td>09:44:29.15</td><td>[1]</td></tr><tr><td>Decl. (J2000)</td><td>-45:46:44.46</td><td>[1]</td></tr><tr><td>TESS-band magnitude</td><td>8.032 ± 0.007</td><td>[2]</td></tr><tr><td>V-band magnitude</td><td>10.153 ± 0.044</td><td>[3]</td></tr><tr><td>Parallax (mas)</td><td>106.173 ± 0.014</td><td>[1]</td></tr><tr><td>Distance (pc)</td><td>9.413 ± 0.003</td><td>[1]</td></tr><tr><td>Star mass M<sub>*</sub> (M<sub>☉</sub>)</td><td>0.455 ± 0.011</td><td>[4]</td></tr><tr><td>Star radius R<sub>*</sub> (R<sub>☉</sub>)</td><td>0.458 ± 0.013</td><td>[4]</td></tr><tr><td>Effective temperature T<sub>eff</sub> (K)</td><td>3522 ± 70</td><td>[4]</td></tr><tr><td>Stellar density ρ<sub>*</sub> (ρ<sub>☉</sub>)</td><td>4.75<sup>+0.44</sup><sub>-0.39</sub></td><td>[4]</td></tr><tr><td>Metallicity [Fe/H]</td><td>-0.01 ± 0.12</td><td>[4]</td></tr><tr><td>Surface gravity log g<sub>*</sub></td><td>4.776 ± 0.026</td><td>[4]</td></tr><tr><td>Luminosity L<sub>*</sub> (L<sub>☉</sub>)</td><td>0.0289<sup>+0.0029</sup><sub>-0.0027</sub></td><td>[4]</td></tr><tr><td>log R'<sub>HK</sub></td><td>-5.169 ± 0.068</td><td>[4]</td></tr><tr><td>Spectral type</td><td>M1.0 V</td><td>[5]</td></tr></tbody></table>

**Note.** [1] Gaia Collaboration et al. (2021), [2] TESS input catalog (TIC; Stassun et al. 2018, 2019), [3] Paegert et al. (2022), [4] This work, [5] Koen et al. (2010).

assessing the presence of instrumental artifacts and non-Gaussian noise in the Stokes V spectra.

# 3. Stellar Parameters

## 3.1. Photospheric and Fundamental Parameters

We derived the spectroscopic parameters using the new coadded HARPS spectrum for GJ 367. Following the prescription in Hirano et al. (2018), we first estimate the stellar effective temperature $T_{eff}$, metallicity [Fe/H], and radius $R_{\star}$ using SpecMatch-Emp (Yee et al. 2017). The code attempts to find a subset of best-matching template spectra from the library to the input spectrum and derives the best empirical values for the above parameters. SpecMatch-Emp returned $T_{eff} = 3522 \pm 70$ K, $R_{\star} = 0.452 \pm 0.045 R_{\odot}$, and [Fe/H] = -0.01 ± 0.12. We then used those parameters to estimate the other stellar parameters as well as refine $R_{\star}$. As described in Hirano et al. (2021), we implemented a Markov Chain Monte Carlo (MCMC) simulation and derived the parameters in a self-consistent manner, making use of the empirical formulae by Mann et al. (2015) and Mann et al. (2019) for the derivations of the stellar mass $M_{\star}$ and radius $R_{\star}$. We found that GJ 367 has a mass of $M_{\star} = 0.455 \pm 0.011 M_{\odot}$ and a radius of $R_{\star} = 0.458 \pm 0.013 R_{\odot}$, with the latter in very good agreement with the value derived using SpecMatch-Emp. In the MCMC implementation, we also derived the stellar density $\rho_{\star}$, surface gravity $\log g_{\star}$, and luminosity $L_{\star}$ of the star. Results of our analysis are listed in Table 1. As the quoted uncertainties of the stellar parameters do not account for possible unknown systematic errors—which in turn might affect the estimates of the planetary parameters—we performed a sanity check and determined the stellar mass and radius using the code BASTA (Aguirre Børsen-Koch et al. 2022). We fitted the derived stellar parameters to the BaSTI isochrones (Hidalgo et al. 2018; science case 4 in Table 1 of the BASTA paper). Starting from the effective temperature, metallicity, and radius, as derived from SpecMatch-Emp, yielded a mass of $M_{\star} = 0.435_{-0.040}^{+0.035} M_{\odot}$
and radius of $R_{\star} = 0.411^{+0.032}_{-0.034} R_{\odot}$. The stellar mass is in good agreement with the previous result. The new estimate of the stellar radius is smaller than the value reported in Table 1, but it is still consistent within ~1.4 σ (where σ is the sum in quadrature of the two nominal uncertainties) with a p-value of ~15%. Assuming a significance level of 5%, the two radii are consistent, providing evidence that our estimates might not be significantly affected by inaccuracy.

## 3.2. Rotation Period

Using archival photometry from the Wide Angle Search for Planets survey (WASP), Lam et al. (2021) found a photometric modulation with a period of 48 ± 2 days. Lam et al. (2021) also measured a Ca II H & K chromospheric activity index of $log R'_{HK} = -5.214 ± 0.074$ from their 77 HARPS spectra. Based on the $log R'_{HK}$-rotation empirical relationship for M dwarfs from Astudillo-Defru et al. (2017a), they estimated a stellar rotation period of $P_{rot} = 58.0 ± 6.9$ days.

We independently derived a Ca II H & K chromospheric activity index of $log R'_{HK} = -5.169 ± 0.068$ from the 371 HARPS spectra and estimated the rotation period of GJ 367 using the same empirical relationship. We found a rotation period of $P_{rot} = 54 ± 6$ d, in good agreement with the previous estimated value. We note that our estimate is consistent within 1σ with the value of $P_{rot} = 51.30 ± 0.13$ d recovered by our sinusoidal signal analysis described in Section 5.2, and with the period of $P_{rot} = 53.67^{+0.65}_{-0.53}$ d derived by our multidimensional Gaussian process (GP) analysis described in Section 5.3.

## 3.3. Magnetic Field

Our spectropolarimetric observation of GJ 367 achieved a median S/N of about 90 over the red HARPS chip. This is insufficient for detecting Zeeman polarization signatures in individual lines even for the most active M dwarfs. To boost the signal, we made use of the least-squares deconvolution procedure (LSD; Donati et al. 1997) as implemented by Kochukhov et al. (2010). The line mask required for LSD was obtained from the VALD database (Ryabchikova et al. 2015) using the atmospheric parameters of GJ 367 and assuming solar abundances (Section 3.1). We used about 5000 lines deeper than 20% of the continuum for LSD, reaching an S/N of 7250 per 1 km s⁻¹ velocity bin. The resulting Stokes V profile has a shape compatible with a Zeeman polarization signature with an amplitude of ≈0.04% (Figure 1). However, with a false-alarm probability (FAP) = 2.3%, detection of this signal is not statistically significant according to the usual detection criteria employed in high-resolution spectropolarimetry (Donati et al. 1997). The mean longitudinal magnetic field, which represents the disk-averaged line-of-sight component of the global magnetic field, derived from this Stokes V profile is $\langle B_z \rangle = -7.3 ± 3.2$ G.

Our spectropolarimetric observation of GJ 367 suggests that this object is not an active M dwarf. Its longitudinal magnetic field was found to be below 10 G, which is much weaker than ≥100–700 G longitudinal fields typical of active M dwarfs (Donati et al. 2008; Morin et al. 2008). Considering that the strength and topology of the global magnetic fields of M dwarfs is systematically changing with the stellar mass (Kochukhov 2021), it is appropriate to compare GJ 367 with early-M dwarfs. To this end, the well-known ~20 Myr old M1V star AU Mic was observed with HARPSpol in the same

<figure><img src="image_5.png" /><figcaption>**Figure 1.** Least-squares deconvolved (LSD) Stokes I, V and null profiles of GJ 367. The polarization profiles are shifted vertically and expanded by a factor of 100 relative to the intensity profile. The vertical dashed lines indicate velocity interval adopted for the longitudinal field measurement.</figcaption></figure>

configuration and with a similar S/N as our GJ 367 observations (Kochukhov & Reiners 2020), yielding consistent detections of polarization signatures and longitudinal fields of up to 50 G. The magnetic activity of GJ 367 is evidently well below that of AU Mic.

## 3.4. Age

Determining the age of M-dwarf stars is especially challenging and depends on the methods being utilized. Lam et al. (2021) estimated an isochronal age of $8.0^{+3.8}_{-4.6}$ Gyr and a gyrochronological age of $4.0 ± 1.1$ Gyr for GJ 367. More recently, Brandner et al. (2022) gave two different estimates: (1) by comparing Gaia EDR3 parallax and photometric measurements with theoretical isochrones, they suggested a young age < 60 Myr. However, as pointed out by the authors, this is not in line with the star's Galactic kinematics that exclude membership to any nearby young moving group; (2) by considering the Galactic dynamical evolution, which indicates an age of 18 Gyr.

In this respect, the results presented in our study shows compelling evidence that GJ 367 is not a young star:

1. The time series of our HARPS RVs and activity indicators give a clear detection of a spot-induced rotation modulation with a period of about 52–54 days, which translates into a gyrochronological age of ~4.6–4.8 Gyr (Barnes 2010; Barnes & Kim 2010).

2. The HARPS spectra of GJ 367 show no significant lithium absorption line at 6708 Å. Figure 2 displays the coadded HARPS spectrum of AU Mic (Zicher et al. 2022; Klein et al. 2022) in the spectral region around the Li I 6708 Å line. AU Mic is a well-known 22-Myr-old M1 star located in the β Pictoris moving group
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Coadded HARPS spectrum of AU Mic (black line) in the spectral region encompassing the Li 16708 Å absorption line. Superimposed with a thick red line is the coadded HARPS spectrum of GJ 367, which has been rotationally broadened to match the $v \sin i_* = 7.8 \text{ km s}^{-1}$ of AU Mic.</figcaption></figure>

(Mamajek & Bell 2014). Superimposed with a thick red line is GJ 367's coadded HARPS spectrum, which has been broadened to match the projected rotational velocity of AU Mic ($v \sin i_* = 7.8 \text{ km s}^{-1}$). Since lithium is quickly depleted in young GKM stars, the lack of lithium in the spectrum of GJ 367 suggests an age $\gtrsim 50$ Myr (Binks & Jeffries 2014; Binks et al. 2021).

3. The low level of magnetic activity inferred by the Ca H & K indicator $\log R'_{HK}$ (Section 3.2) and the weak magnetic field (Section 3.3) are consistent with an old, inactive M-dwarf scenario (Pace 2013).

4. We measured an average $H\alpha$ equivalent width of EW = 0.0638 $\pm$ 0.0014 Å. Using the empirical relation that connects the $H\alpha$ equivalent width with stellar age (Kiman et al. 2021), this translates into an age $\gtrsim 300$ Myr.

We therefore conclude that GJ 367 is a rather slowly rotating, old star with a low magnetic activity level, rather than a young M dwarf. This conclusion is consistent with a recent study by Gaidos et al. (2023), who measured an age of 7.95 $\pm$ 1.31 Gyr from the M-dwarf rotation–age relation.

## 4. Frequency Analysis of the HARPS Time Series

In order to search for the Doppler reflex motion induced by the USP planet GJ 367 b and unveil the presence of potential additional signals associated with other orbiting companions and/or stellar activity, we performed a frequency analysis of the HARPS RV measurements and activity indicators. For this analysis we did not include the HARPS measurements<fnref n="38" /> presented in Lam et al. (2021) and used only the Doppler data collected between 2020 November 9 and 2022 April 18 (UT) as part of our HARPS large program, to avoid spurious peaks introduced by the poor sampling of the existing old data set, and to avoid having to account for the RV offset caused by the refurbishment of the instrument.

Figure 3 shows the generalized Lomb–Scargle (GLS; Zechmeister & Kürster 2009) periodograms of the HARPS RVs and activity indicators in two frequency ranges, i.e., 0.000–0.130 day$^{-1}$ (left panels) and 3.075–3.125 day$^{-1}$ (right panels), with the former including the frequencies at which we expect to see activity signals at the rotation period of the star, and the latter encompassing the orbital frequency of the transiting planet GJ 367 b. The horizontal dashed lines mark the GLS powers at the 0.1%, 1%, and 10% FAP. The FAP was estimated following the bootstrap method described in Murdoch et al. (1993), i.e., by computing the GLS periodogram of 10⁶ mock time series obtained by randomly shuffling the measurements and their uncertainties, while keeping the time stamps fixed. In this work we assumed a peak to be significant if its FAP < 0.1%.

The GLS periodogram of the HARPS RVs (Figure 3, upper panel) shows its most significant peak at $f_1 = 0.086 \text{ day}^{-1}$, corresponding to a period of about 11.5 days. This peak is not detected in the activity indicators, providing evidence that the 11.5 days signal is caused by a second planet orbiting the host star, hereafter referred to as GJ 367 c.

We used the pre-whitening technique (Hatzes et al. 2010) to identify additional significant signals and successively remove them from the RV time series. We employed the code *pyaneti* (Barragán et al. 2019, 2022) to subtract the 11.5 day signal from the HARPS RVs assuming a circular model, adopting uniform priors centered around the period and phase derived from the periodogram analysis, while allowing the systemic velocity and RV semiamplitude to uniformly vary over a wide range.

The periodogram of the RV residuals following the subtraction of the signal at $f_1$ (Figure 3, second panel) shows its most significant peak at $f_2 = 0.019 \text{ day}^{-1}$ (52.2 days). Iterating the pre-whitening procedure and removing the signal at $f_2$, we found a significant peak at $f_3 = 0.029 \text{ day}^{-1}$ (34 days). This peak has no significant counterpart in the activity indicators, suggesting it is associated to a third planet orbiting the star, hereafter referred as GJ 367 d (Figure 3, third panel). The periodograms of the CRX, DLW, $H\alpha$, $H\beta$, $H\gamma$, $\log R'_{HK}$, and NaD show significant peaks in the range 48–54 days (Figure 3, lower panels), i.e., close to the rotation period of the star, suggesting that the peak at $f_2 = 0.019 \text{ day}^{-1}$ (52.2 days) seen in the RV residuals is very likely associated to the presence of active regions appearing and disappearing on the visible stellar hemisphere as the star rotates about its axis.

After removing the signal at $f_3$, we found significant power at $f_4 = 0.009 \text{ day}^{-1}$ (~115 days; Figure 3, fourth panel). The periodograms of the activity indicators show also significant power around $f_4$, providing evidence this signal is associated to stellar activity. As we will discuss in Section 5.3, we interpreted the power at $f_4$ as the evolution timescale of active regions. Finally, we found that the RV residuals also show a significant fifth peak at $f_5 = f_6 = 3.106 \text{ day}^{-1}$ (0.322 day), the orbital frequency of the USP planet GJ 367 b, further confirming the planetary nature of the transit signal identified in the TESS data and announced by Lam et al. (2021).

Finally, we computed the GLS periodogram of the HARPS RVs including all of the available data, i.e., the data acquired before and after the refurbishment of the instrument. Adding the old data points increases the baseline of our observations and, consequently, the frequency resolution. However, the resulting periodogram is “jagged,” owing to the presence of aliases with very small frequency spacing, making it more difficult to identify the true peaks.

<footnote n="38">³⁸ Twenty-four measurements taken with the old fiber bundle between 2003 December 12 and 2010 February 7 (UT) and 77 measurements acquired between 2019 June 23 and 2020 March 23 (UT) with the new fiber bundle.</footnote>
<figure><img src="image_7.png" /><figcaption>**Figure 3.** Generalized Lomb–Scargle (GLS) periodograms of the HARPS RV measurements (upper panel); RV residuals after subtracting the $f_1$ signal at 11.5 days (second panel), the $f_2$ signal related to stellar activity at 52.2 days (third panel), the $f_3$ signal at 34 days (fourth panel), and the $f_4$ and $f_5 = f_b$ signals at 115 days and 0.322 days (fifth panel). Also shown are periodograms of the activity indicators (remaining panels). The 10%, 1%, and 0.1% false-alarm probabilities (FAPs) estimated using the bootstrap method are shown with horizontal blue lines. The red vertical lines mark the orbital frequencies of the transiting planet GJ 367 b ($f_5 = f_b = 3.106$ day$^{-1}$), and of the additional Doppler signals we found in the HARPS data, which are associated to the presence of two additional orbiting planets ($f_1 = f_c = 0.086$ day$^{-1}$ and $f_3 = f_d = 0.029$ day$^{-1}$). The shaded yellow bands indicate the rotation period of the star centered around $f_2$, and the long-period stellar signal $f_4$.</figcaption></figure>
# 5. Data Analysis

We modeled the TESS transit light curves and HARPS RV measurements using three different approaches. The methods differ in the way the Doppler data are fitted, as described in the following subsections.

## 5.1. Floating Chunk Offset Method

We used the floating chunk offset (FCO) method to determine the semiamplitude $K_b$ of the Doppler reflex motion induced by the USP planet GJ 367 b. Pioneered by Hatzes et al. (2011) for the mass determination of the USP planet CoRoT-7 b, the FCO method relies on the reasonable assumption that, within a single night, the RV variation of the star is mainly induced by the orbital motion of the USP planet rather than stellar rotation, magnetic cycles, or orbiting companions on longer-period orbits. As the RV component due to long-period phenomena can thus be treated as constant within a given night, introducing nightly offsets filters out any long-term RV variations, allowing one to disentangle the reflex motion of the USP planet from additional long-period Doppler signals.

With an orbital period of only 7.7 hr, the USP planet GJ 367 b is suitable to the FCO method (see, e.g., Gandolfi et al. 2017; Barragán et al. 2018). GJ 367 is accessible for up to 8 hr at an airmass < 1.5 (i.e., altitude > 40°) from La Silla Observatory, allowing one to cover one full orbit in one single night by acquiring multiple HARPS spectra per night. Within the nightly visibility window of GJ 367, the phase of the long-term signals does not change significantly, the variation being 0.029, 0.010, and 0.006 for the 11.5, 34, and 52 day signals, respectively.

We simultaneously modeled the TESS transit light curves and HARPS RV measurements using the open source software suite pyaneti (Barragán et al. 2019, 2022). The code utilizes a Bayesian approach in combination with MCMC sampling to infer the parameters of planetary systems. The photometric data included in the analysis are subsets of the TESS light curve. We selected 2.5 hr of TESS data points centered around each transit and detrended each light-curve segment by fitting a second-order polynomial to the out-of-transit data. Following Hatzes et al. (2011), we divided the HARPS RVs into subsets (“chunks”) of nightly measurements and analyzed only those chunks containing at least two RVs per observing night, leading to a total of 96 chunks.

The RV model includes one Keplerian orbit for the transiting planet GJ 367 b and 96 nightly offsets. We fitted for a nonzero eccentricity adopting the parameterization proposed by Anderson et al. (2010) for the eccentricity *e* and the argument of periastron of the stellar orbit $\omega_*$ (i.e., $\sqrt{e} \sin \omega_*$ and $\sqrt{e} \cos \omega_*$). We fitted for a photometric and an RV jitter term to account for any instrumental noise not included in the nominal TESS and HARPS uncertainties. We used the limb darkened quadratic model by Mandel & Agol (2002) for the transit light curve. We adopted Gaussian priors for the limb darkening coefficients, using the values derived by Claret (2017) for the TESS passband, and we imposed conservative error bars of 0.1 on both the linear and the quadratic limb darkening term. As the shallow transit light curve of GJ 367 b poorly constrains the scaled semimajor axis ($a_b/R_*$), we sampled for the stellar density $\rho_*$ using a Gaussian prior on the star's mass and radius as derived in Section 3, and recovered the scaled semimajor axis of the planet using the orbital period and Kepler's third law of planetary motion (see,

e.g., Winn 2010). We adopted uniform priors over a wide range for all of the remaining model parameters. We ran 500 independent Markov chains. The posterior distributions were created using the last 5000 iterations of the converged chains with a thin factor of 10, leading to a distribution of 250,000 data points for each model parameter. The chains were initialized at random values within the priors ranges. This ensured a homogeneous sampling of the parameter space. We followed the same procedure and convergence test as described in Barragán et al. (2019). The final estimates and their $1\sigma$ uncertainties were taken as the median and 68% of the credible interval of the posterior distributions. Table 2 reports prior ranges and posterior values of the fitted and derived system parameters. Figure 5 displays the phase-folded RV curve with our HARPS data, along with the best-fitting Keplerian model. Different colors refer to different nights. Figure 4 shows the phase-folded transit light curve of GJ 367 b, along with the TESS data and best-fitting transit model. We found an RV semiamplitude variation of $K_b = 1.003 \pm 0.078 \text{ m s}^{-1}$, which translates into a planetary mass of $M_b = 0.633 \pm 0.050 M_\oplus$ (7.9% precision). The depth of the transit light curve implies a radius of $R_b = 0.699 \pm 0.024 R_\oplus$ (3.4% precision) for GJ 367 b. When combined together, the planetary mass and radius yield a mean density of $\rho_b = 10.2 \pm 1.3 \text{ g cm}^{-3}$ (12.7% precision). The eccentricity of the USP planet ($e_b = 0.06_{-0.04}^{+0.07}$) is consistent with zero, as expected given the short tidal evolution timescale and the age of the system. Assuming a circular orbit, our fit gives an RV semiamplitude of $K_b = 1.001 \pm 0.077 \text{ m s}^{-1}$, in excellent agreement with the values listed in Table 2.

We note that 27 of the HARPS RV measurements were taken during transits of GJ 367 b. Assuming the star is seen equator-on ($i_* = 90^\circ$), its rotation period of $P_{\text{rot}} \approx 51-54$ days (Section 3.2) implies an equatorial rotational velocity of $v_{\text{rot}} \approx 0.45 \text{ km s}^{-1}$. Using the equations in Triaud (2018), we estimated the semiamplitude of the Rossiter-McLaughlin effect to be $\approx 0.05 \text{ m s}^{-1}$, which is too small to cause any detectable effect given our RV uncertainties.

## 5.2. Sinusoidal Activity Signal Modeling

Using the FCO method, we cannot determine the semiamplitude of the Doppler signals induced by the two outer planets and by stellar activity (Section 4). In the analysis described in this section, we treated the RV signals associated with stellar activity as coherent sinusoidal signals. Once again, we used the code pyaneti and performed an MCMC analysis similar to the one described in Section 5.1. The RV model includes three Keplerians, to account for the Doppler reflex motion induced by the three planets GJ 367 b, c, and d, and two additional sine functions, to account for the activity-induced signals at the rotation period of the star (~52 days) and at the evolution timescale of active regions (~115 days), as described in Section 4. We used Gaussian priors for the orbital period and time of first transit of GJ 367 b as derived in Section 5.1, and uniform wide priors for all of the remaining parameters. We fitted for an RV jitter term, as well as for a nonzero eccentricity both for the USP planet, and for the two outer companions, following the $e-\omega_*$ parameterization proposed by Anderson et al. (2010). Details of the fitted parameters and prior ranges are given in Table 3. We used 500 independent Markov chains initialized randomly inside the prior ranges. Once all chains converged, we used the last 5000 iterations and saved the chain
<table><caption>**Table 2**
GJ 367 b Parameters from the Joint FCO and Transit Modeling with pyaneti</caption><thead><tr><th>GJ 367 b</th><th>Prior</th><th>Derived Value</th></tr></thead><tbody><tr><td colspan="3"><strong>Model parameters</strong></td></tr><tr><td>Orbital period P<sub>orb, b</sub> [days]</td><td>U[0.3219221,0.3219229]</td><td>0.3219225 ± 0.0000002</td></tr><tr><td>Transit epoch T<sub>0,b</sub> [BJD<sub>TDB</sub> − 2,450,000]</td><td>U[8544.13235,8544.14035]</td><td>8544.13635 ± 0.00040</td></tr><tr><td>Planet-to-star radius ratio R<sub>b</sub>/R<sub>*</sub></td><td>U[0.001,0.025]</td><td>0.01399 ± 0.00028</td></tr><tr><td>Impact parameter b<sub>b</sub></td><td>U[0,1]</td><td>0.584<sup>+0.034</sup><sub>-0.037</sub></td></tr><tr><td>√e<sub>b</sub> sin ω<sub>*</sub>,b</td><td>U[-1.0,1.0]</td><td>0.16<sup>+0.17</sup><sub>-0.22</sub></td></tr><tr><td>√e<sub>b</sub> cos ω<sub>*</sub>,b</td><td>U[-1.0,1.0]</td><td>0.04 ± 0.14</td></tr><tr><td>Radial velocity semiamplitude variation K<sub>b</sub> [m s<sup>-1</sup>]</td><td>U[0,50]</td><td>1.003 ± 0.078</td></tr><tr><td colspan="3"><strong>Derived parameters</strong></td></tr><tr><td>Planet mass M<sub>b</sub> [M<sub>⊕</sub>]</td><td>...</td><td>0.633 ± 0.050</td></tr><tr><td>Planet radius R<sub>b</sub> [R<sub>⊕</sub>]</td><td>...</td><td>0.699 ± 0.024</td></tr><tr><td>Planet mean density ρ<sub>b</sub> [g cm<sup>-3</sup>]</td><td>...</td><td>10.2 ± 1.3</td></tr><tr><td>Semimajor axis of the planetary orbit a<sub>b</sub> [au]</td><td>...</td><td>0.00709 ± 0.00027</td></tr><tr><td>Orbit eccentricity e<sub>b</sub></td><td>...</td><td>0.06<sup>+0.07</sup><sub>-0.04</sub></td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>*</sub>,b [deg]</td><td>...</td><td>66<sup>+41</sup><sub>-108</sub></td></tr><tr><td>Orbit inclination i<sub>b</sub> [deg]</td><td>...</td><td>79.89<sup>+0.87</sup><sub>-0.85</sub></td></tr><tr><td>Transit duration τ<sub>14,b</sub> [hr]</td><td>...</td><td>0.629 ± 0.008</td></tr><tr><td>Equilibrium temperature T<sub>eq,b</sub> [K]<sup>a</sup></td><td>...</td><td>1365 ± 32</td></tr><tr><td>Received irradiance F<sub>b</sub> [F<sub>⊕</sub>]</td><td>...</td><td>579<sup>+57</sup><sub>-52</sub></td></tr><tr><td colspan="3"><strong>Additional model parameters</strong></td></tr><tr><td>Stellar density ρ<sub>*</sub> [g cm<sup>-3</sup>]</td><td>N[6.68,0.59]</td><td>6.76 ± 0.59</td></tr><tr><td>Parameterized limb-darkening coefficient q<sub>1</sub></td><td>N[0.3766,0.1000]</td><td>0.343 ± 0.095</td></tr><tr><td>Parameterized limb-darkening coefficient q<sub>2</sub></td><td>N[0.1596,0.1000]</td><td>0.163<sup>+0.096</sup><sub>-0.088</sub></td></tr><tr><td>Radial velocity jitter term σ<sub>RV,HARPS</sub> [m s<sup>-1</sup>]</td><td>J[0,100]</td><td>0.43 ± 0.08</td></tr><tr><td>TESS jitter term σ<sub>TESS</sub></td><td>J[0,100]</td><td>0.00004 ± 0.00003</td></tr></tbody></table>

Note. U[a, b] refers to uniform priors between a and b; N[a, b] refers to Gaussian priors with mean a and standard deviation b; J[a, b] refers to Jeffrey's priors between a and b. Inferred parameters and uncertainties are defined as the median and the 68.3% credible interval of their posterior distributions.

<figure><img src="image_8.png" /><figcaption>**Figure 4.** TESS transit light curve of GJ 367 b and best-fitting model folded at the orbital period of the planet.</figcaption></figure>

states every 10 iterations. This approach generated a posterior distribution of 250,000250000 points for each fitted parameter.

The RV semiamplitude variations induced by the three planets are $K_b = 1.10 \pm 0.14 \text{ m s}^{-1}$, $K_c = 2.01 \pm 0.15 \text{ m s}^{-1}$, and $K_d = 1.98 \pm 0.15 \text{ m s}^{-1}$, which imply planetary masses and minimum masses of $M_b = 0.699 \pm 0.083 M_⊕$, $M_c \sin i_c = 4.08 \pm 0.30 M_⊕$, and $M_d \sin i_d = 5.93 \pm 0.45 M_⊕$ for GJ 367 b, c, and d, respectively, whereas the RV semiamplitudes induced by stellar activity signals at 51.3 and 138 days are of $K_{*,Rot} =$

<figure><img src="image_9.png" /><figcaption>**Figure 5.** HARPS RVs of GJ 367 phase-folded at the orbital period of the USP planet and best-fitting model as derived using the FCO method. The different colors refer to the 96 different nightly chunks, which include at least two measurements per night.</figcaption></figure>

$2.52 \pm 0.13 \text{ m s}^{-1}$ and $K_{*,Evol} = 1.25 \pm 0.69 \text{ m s}^{-1}$. The RV time series along with the best-fitting model are shown in Figure 6. The phase-folded RV curves for each signal are displayed in Figure 7.

### 5.3. Multidimensional Gaussian Process Approach

We also followed a multidimensional GP approach to account for the stellar signals in our RV time series (Rajpaul et al. 2015). This approach models the RVs along with time

<footnote>a Assuming zero albedo and uniform redistribution of heat.</footnote>
<table><caption>**Table 3**
System Parameters as Derived Modeling the Stellar Signals with Two Sine Functions</caption><thead><tr><th>Parameter</th><th>Prior</th><th>Derived Value</th></tr></thead><tbody><tr><td colspan="3"><strong>GJ 367 b</strong></td></tr><tr><td colspan="3"><strong>Model parameters</strong></td></tr><tr><td>Orbital period P<sub>orb, b</sub> [days]</td><td>N[0.3219225, 0.0000002]</td><td>0.3219225 ± 0.0000002</td></tr><tr><td>Transit epoch T<sub>0,b</sub> [BJD<sub>TDB</sub>-2,450,000]</td><td>N[8544.1364,0.0004]</td><td>8544.13632 ± 0.00040</td></tr><tr><td>√e<sub>b</sub> sin ω<sub>*,b</sub></td><td>U[-1.0,1.0]</td><td>-0.23<sup>+0.30</sup><sub>-0.23</sub></td></tr><tr><td>√e<sub>b</sub> cos ω<sub>*,b</sub></td><td>U[-1.0,1.0]</td><td>-0.07 ± 0.13</td></tr><tr><td>Radial velocity semiamplitude variation K<sub>b</sub> [m s<sup>-1</sup>]</td><td>U[0.00, 0.05]</td><td>1.10 ± 0.14</td></tr><tr><td colspan="3"><strong>Derived parameters</strong></td></tr><tr><td>Planet mass M<sub>b</sub> [M<sub>⊕</sub>]<sup>a</sup></td><td>...</td><td>0.699 ± 0.083</td></tr><tr><td>Orbit eccentricity e<sub>b</sub></td><td>...</td><td>0.10<sup>+0.14</sup><sub>-0.07</sub></td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>*,b</sub> [deg]</td><td>...</td><td>251<sup>+23</sup><sub>-102</sub></td></tr><tr><td colspan="3"><strong>GJ 367 c</strong></td></tr><tr><td colspan="3"><strong>Model parameters</strong></td></tr><tr><td>Orbital period P<sub>orb, c</sub> [days]</td><td>U[11.4858,11.5858]</td><td>11.543 ± 0.005</td></tr><tr><td>Time of inferior conjunction T<sub>0,c</sub> [BJD<sub>TDB</sub>-2,450,000]</td><td>U[9152.6591,9154.6591]</td><td>9153.46 ± 0.21</td></tr><tr><td>√e<sub>c</sub> sin ω<sub>*,c</sub></td><td>U[-1,1]</td><td>0.38<sup>+0.10</sup><sub>-0.13</sub></td></tr><tr><td>√e<sub>c</sub> cos ω<sub>*,c</sub></td><td>U[-1,1]</td><td>0.27<sup>+0.11</sup><sub>-0.14</sub></td></tr><tr><td>Radial velocity semiamplitude variation K<sub>c</sub> [m s<sup>-1</sup>]</td><td>U[0.00, 0.05]</td><td>2.01 ± 0.15</td></tr><tr><td colspan="3"><strong>Derived parameters</strong></td></tr><tr><td>Planet minimum mass M<sub>c</sub> sin i<sub>c</sub> [M<sub>⊕</sub>]</td><td>...</td><td>4.08 ± 0.30</td></tr><tr><td>Orbit eccentricity e<sub>c</sub></td><td>...</td><td>0.23 ± 0.07</td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>*,c</sub> [deg]</td><td>...</td><td>55 ± 18</td></tr><tr><td colspan="3"><strong>GJ 367 d</strong></td></tr><tr><td colspan="3"><strong>Model parameters</strong></td></tr><tr><td>Orbital period P<sub>orb, d</sub> [days]</td><td>U[34.0016,34.6016]</td><td>34.39 ± 0.06</td></tr><tr><td>Time of inferior conjunction T<sub>0,d</sub> [BJD<sub>TDB</sub>-2,450,000]</td><td>U[9179.2710,9183.2710]</td><td>9180.90<sup>+0.70</sup><sub>-0.81</sub></td></tr><tr><td>√e<sub>d</sub> cos ω<sub>*,d</sub></td><td>U[-1,1]</td><td>-0.10<sup>+0.20</sup><sub>-0.18</sub></td></tr><tr><td>√e<sub>d</sub> cos ω<sub>*,d</sub></td><td>U[-1,1]</td><td>0.16<sup>+0.16</sup><sub>-0.20</sub></td></tr><tr><td>Radial velocity semiamplitude variation K<sub>d</sub> [m s<sup>-1</sup>]</td><td>U[0.00, 0.05]</td><td>1.98 ± 0.15</td></tr><tr><td colspan="3"><strong>Derived parameters</strong></td></tr><tr><td>Planet minimum mass M<sub>d</sub> sin i<sub>d</sub> [M<sub>⊕</sub>]</td><td>...</td><td>5.93 ± 0.45</td></tr><tr><td>Orbit eccentricity e<sub>d</sub></td><td>...</td><td>0.08<sup>+0.07</sup><sub>-0.05</sub></td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>*,d</sub> [deg]</td><td>...</td><td>277<sup>+58</sup><sub>-242</sub></td></tr><tr><td colspan="3"><strong>Stellar activity-induced RV signal</strong></td></tr><tr><td>Rotation period P<sub>*,Rot</sub> [days]</td><td>U[50.0903,52.0903]</td><td>51.30 ± 0.13</td></tr><tr><td>Rotation RV semiamplitude K<sub>*,Rot</sub> [m s<sup>-1</sup>]</td><td>U[0.00, 0.05]</td><td>2.52 ± 0.13</td></tr><tr><td>Active region evolution period P<sub>*,Evol</sub> [days]</td><td>U[103.1797,163.1797]</td><td>138 ± 2</td></tr><tr><td>Active region evolution RV semiamplitude K<sub>*,Evol</sub> [m s<sup>-1</sup>]</td><td>U[0.00, 0.05]</td><td>1.25 ± 0.14</td></tr><tr><td colspan="3"><strong>Additional model parameters</strong></td></tr><tr><td>Systemic velocity γ<sub>HARPS</sub> [m s<sup>-1</sup>]</td><td>U[47.806,48.025]</td><td>47.91674 ± 0.00013</td></tr><tr><td>Radial velocity jitter term σ<sub>RV,HARPS</sub> [m s<sup>-1</sup>]</td><td>J[0,100]</td><td>1.59 ± 0.07</td></tr></tbody></table>

Note. U[a, b] refers to uniform priors between a and b; N[a, b] refers to Gaussian priors with mean a and standard deviation b; J[a, b] refers to Jeffrey's priors between a and b. Inferred parameters and uncertainties are defined as the median and the 68.3% credible interval of their posterior distributions.

series of activity indicators assuming that the same GP, a function G(t), can describe them both. The function G(t) represents the projected area of the visible stellar disk that is covered by active regions at a given time. For our best GP analysis, we selected the activity indicator that shows the strongest signal in the periodograms, i.e., the DLW, and modeled the RVs alongside this activity index. We created a two-dimensional GP model via

$$RV = A_{RV}G(t) + B_{RV}\dot{G}(t), \quad (1)$$

$$DLW = A_{DLW}G(t). \quad (2)$$

The amplitudes $A_{RV}$, $B_{RV}$, and $A_{DLW}$ are free parameters, which relate the individual time series to G(t). The RV data are modeled as a function of G(t) and its time derivative $\dot{G}(t)$, since they depend both on the fraction of the stellar surface covered by active regions, and on how these regions evolve and migrate on the disk. The DLW, which measures the width of the spectral lines, has been proven to be a good tracer of the fraction of the surface covered by active regions, and is thus

<footnote>a Assuming an orbital inclination of i<sub>b</sub> = 79.89<sup>+0.87°</sup><sub>-0.85°</sub>, from the modeling of TESS transit light curves (Section 5.1).</footnote>
<figure><img src="image_10.png" /><figcaption>Figure 6. HARPS RV time series of GJ 367 along with the best-fitting five-signal model (three planets + stellar rotation + long-period stellar signal). Data are shown as blue filled circles with their nominal uncertainties. The vertical gray lines mark the error bars including the RV jitter.</figcaption></figure>

expected to be solely proportional to G(t) (Zicher et al. 2022). For our covariance matrix, we used a quasiperiodic kernel

$$ \gamma(t_i, t_j) = \exp \left[ - \frac{\sin^2[\pi(t_i - t_j)/P_{GP}]}{2\lambda_p^2} - \frac{(t_i - t_j)^2}{2\lambda_e^2} \right] \quad (3) $$

and its derivatives (Barragán et al. 2022; Rajpaul et al. 2015). The parameter $P_{GP}$ in Equation (3) is the characteristic period of the GP, which is interpreted as the stellar rotation period; $\lambda_p$ is the inverse of the harmonic complexity, which is associated with the distribution of active regions on the stellar surface (Aigrain et al. 2015); $\lambda_e$ is the long-term evolution timescale, i.e., the active region lifetime on the stellar surface. We performed the multidimensional GP regression using pya-neti (described in Barragán et al. 2022), adding three Keplerians to account for the Doppler reflex motion of the three planets, as described in Section 5.2.

We performed an MCMC analysis setting informative Gaussian priors based on the transit ephemeris of the innermost planet GJ 367 b, as derived in Section 5.1, and uniform priors for all of the remaining parameters. We also used uniform priors to sample for the multidimensional GP hyper-parameters. We included a jitter term to the diagonal of the covariance for each time series.

We performed our fit with 500 Markov chains to sample the parameter space. The posterior distributions were created using the last 5000 iterations of the converged chains with a thin factor of 10, leading to a distribution of data points for each fitted parameter. The chains were initialized at random values within the priors ranges. This ensured a homogeneous sampling of the parameter space.

Priors and results are listed in Table 4. Planets GJ 367 b, c, and d are significantly detected in the HARPS RV time series with Doppler semiamplitudes of $K_b = 0.86 \pm 0.15$ m s⁻¹, $K_c = 1.99 \pm 0.17$ m s⁻¹, and $K_d = 2.03 \pm 0.16$ m s⁻¹, respectively. These imply planetary masses and minimum masses of $M_b = 0.546 \pm 0.093 M_⊕$, $M_c \sin i_c = 4.13 \pm 0.36 M_⊕$, and $M_d \sin i_d = 6.03 \pm 0.49 M_⊕$. The resulting GP hyper-parameters are $P_{GP} = 53.67^{+0.65}_{-0.53}$ days, $\lambda_p = 0.44 \pm 0.05$, and $\lambda_e = 114 \pm 19$ days. The characteristic period $P_{GP}$ is in agreement with the stellar rotation period, as discussed in Section 3.2, as well as the long-term evolution

timescale $\lambda_e$, which is in agreement with the long-period signal found in the analyses of Sections 4 and 5.2.

Figure 8 shows the RV and DLW time series, along with the inferred models, whereas Figure 9 displays the phase-folded RV curves of the three planets and the best-fitting models.

## 6. Discussion

The three techniques used to determine the mass of GJ 367 b give results that are consistent to within ~1σ. We adopted the results from the FCO method, which gives a planetary mass of $M_b = 0.633 \pm 0.050 M_⊕$ (7.9% precision). Our result differs by about ~1σ from the mass of $M_b = 0.546 \pm 0.078 M_⊕$ reported by Lam et al. (2021), which was also derived using the FCO method applied on 20 HARPS data chunks that do not entirely cover the orbital phase of the USP planet, as opposed to our 96 chunks. We found that GJ 367 b has a radius of $R_b = 0.699 \pm 0.024 R_⊕$ (3.4% precision), consistent with the value of $R_b = 0.718 \pm 0.054 R_⊕$ from Lam et al. (2021), but more precise, thanks to the additional TESS photometry and increased cadence. Lam et al. (2021) reported a density of $\rho_b = 8.1 \pm 2.2$ g cm⁻³. The higher mass but similar radius measured in this work make the USP planet denser, with a ultra-high density of $\rho_b = 10.2 \pm 1.3$ g cm⁻³.

GJ 367 b belongs to the handful of small USP planets ($R_p < 2 R_⊕$, $M_p < 10 M_⊕$) whose masses and radii are known with a precision better than 20%. Figure 10 shows the mass-radius diagram for small USP planets along with the theoretical composition models for rocky worlds (Zeng et al. 2016). GJ 367 b is the smallest and densest USP planet known to date. The position of the planet on the mass-radius diagram suggests that its composition is dominated by iron. Taking into account its mean density, GJ 367 b leads the class of super-Mercuries, namely, extremely dense planets containing an excess of iron, analogous to Mercury: K2-229 b (Santerne et al. 2018), K2-38 b (Toledo-Padrón et al. 2020), K2-106 b (Guenther et al. 2017), Kepler-107 c (Bonomo et al. 2019), Kepler-406 b (Marcy et al. 2014), HD 137496b (Azevedo Silva et al. 2022), HD 23472b (Barros et al. 2022), and TOI-1075b (Essack et al. 2022).

The two nontransiting planets GJ 367 c and GJ 367 d have orbital periods of ~11.5 d and 34.4 days, respectively, and minimum masses of $M_c \sin i_c = 4.13 \pm 0.36 M_⊕$ and $M_d \sin i_d = 6.03 \pm 0.49 M_⊕$, as derived adopting the multidimensional GP
<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 7.** Phase-folded RV curves of GJ 367 b (a), GJ 367 c (b), GJ 367 d (c), stellar rotation (d), and long-period stellar signal (e). Data are shown as blue filled circles with their nominal uncertainties. The vertical gray lines mark the error bars including the RV jitter.</figcaption></figure>

approach to model stellar activity (Section 5.3). We note that our minimum mass determinations are in very good agreement with those of $M_c \sin i_c = 4.08 \pm 0.30 M_⊕$ and $M_d \sin i_d = 5.93 \pm 0.45 M_⊕$ that we determined modeling stellar activity with two sinusoidal components (Section 5.2).

If the orbits of GJ 367 b, c, and d were coplanar ($i_b = i_c = i_d = 79^\circ 9$), planets c and d would have true masses of $M_c = 4.19 \pm 0.35 M_⊕$ and $M_d = 6.12 \pm 0.48 M_⊕$, respectively. Using the mass-radius relation for small rocky planets from Otegi et al. (2020), we found that GJ 367 c and d are expected to have radii of $\sim 1.6$ and $\sim 1.7 R_⊕$, respectively, making them two super-Earths with mean densities of $\sim 6 \text{ g cm}^{-3}$. We

searched the TESS light curves for possible transits of the two outer companions with the DST code (Section 2.1) masking out the transits of the UPS planet, but we found no other significant transit signals. Under the assumption that the orbits of the three planets are coplanar, the impact parameters of planets c and d would be $b_c \approx 6$ and $b_d \approx 13$, respectively. This would account for the nondetection of the transit signals of GJ 367 c and GJ 367 d in the TESS light curves. In this scenario, GJ 367 c and d would transit their host star only if their radii were unphysically large, i.e., $R_c > 1.9 R_⊕$ and $R_d > 5 R_⊕$. The other configuration for the outer planets to transit is if they are mutually inclined with planet b ($i_c \neq i_b$; $i_d \neq i_b$). This is not an
<table><caption>**Table 4**
System Parameters as Derived Modeling the Stellar Signals with a GP</caption><thead><tr><th>Parameter</th><th>Prior</th><th>Derived Value</th></tr></thead><tbody><tr><td colspan="3"><strong>GJ 367 b</strong></td></tr><tr><td colspan="3"><em>Model parameters</em></td></tr><tr><td>Orbital period P<sub>orb,b</sub> [days]</td><td>N[0.3219225, 0.0000002]</td><td>0.3219225 ± 0.0000002</td></tr><tr><td>Transit epoch T<sub>0,b</sub> [BJD<sub>TDB</sub> − 2,450,000]</td><td>N[8544.1364, 0.0004]</td><td>8544.13631 ± 0.00038</td></tr><tr><td>√e<sub>b</sub> sin ω<sub>b</sub></td><td>U[−1.0,1.0]</td><td>0.13<sup>+0.29</sup><sub>−0.32</sub></td></tr><tr><td>√e<sub>b</sub> cos ω<sub>b</sub></td><td>U[−1.0,1.0]</td><td>−0.05 ± 0.21</td></tr><tr><td>Radial velocity semiamplitude variation K<sub>b</sub> [m s<sup>−1</sup>]</td><td>U[0.00, 0.05]</td><td>0.86 ± 0.15</td></tr><tr><td colspan="3"><em>Derived parameters</em></td></tr><tr><td>Planet mass M<sub>b</sub> [M<sub>⊕</sub>]<sup>a</sup></td><td>...</td><td>0.546 ± 0.093</td></tr><tr><td>Orbit eccentricity e<sub>b</sub></td><td>...</td><td>0.10<sup>+0.13</sup><sub>−0.07</sub></td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>b,c</sub> [deg]</td><td>...</td><td>71<sup>+60</sup><sub>−173</sub></td></tr><tr><td colspan="3"><strong>GJ 367 c</strong></td></tr><tr><td colspan="3"><em>Model parameters</em></td></tr><tr><td>Orbital period P<sub>orb,c</sub> [days]</td><td>U[11.4, 11.6]</td><td>11.5301 ± 0.0078</td></tr><tr><td>Time of inferior conjunction T<sub>0,c</sub> [BJD<sub>TDB</sub> − 2,450,000]</td><td>U[9153.0, 9155.0]</td><td>9153.84 ± 0.30</td></tr><tr><td>√e<sub>c</sub> sin ω<sub>c</sub></td><td>U[−1.0,1.0]</td><td>−0.11<sup>+0.23</sup><sub>−0.20</sub></td></tr><tr><td>√e<sub>c</sub> cos ω<sub>c</sub></td><td>U[−1.0,1.0]</td><td>0.14<sup>+0.19</sup><sub>−0.15</sub></td></tr><tr><td>Radial velocity semiamplitude variation K<sub>c</sub> [m s<sup>−1</sup>]</td><td>U[0.00, 0.05]</td><td>1.99 ± 0.17</td></tr><tr><td colspan="3"><em>Derived parameters</em></td></tr><tr><td>Planet minimum mass M<sub>c</sub> sin i<sub>c</sub> [M<sub>⊕</sub>]</td><td>...</td><td>4.13 ± 0.36</td></tr><tr><td>Orbit eccentricity e<sub>c</sub></td><td>...</td><td>0.09 ± 0.07</td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>c,d</sub> [deg]</td><td>...</td><td>−34<sup>+74</sup><sub>−54</sub></td></tr><tr><td colspan="3"><strong>GJ 367 d</strong></td></tr><tr><td colspan="3"><em>Model parameters</em></td></tr><tr><td>Orbital period P<sub>orb,d</sub> [days]</td><td>U[30.0,40.0]</td><td>34.369 ± 0.073</td></tr><tr><td>Time of inferior conjunction T<sub>0,d</sub> [BJD<sub>TDB</sub> − 2,450,000]</td><td>U[9173.0, 9180.0]</td><td>9181.82 ± 1.10</td></tr><tr><td>√e<sub>d</sub> sin ω<sub>d</sub></td><td>U[−1.0,1.0]</td><td>−0.09 ± 0.19</td></tr><tr><td>√e<sub>d</sub> cos ω<sub>d</sub></td><td>U[−1.0,1.0]</td><td>−0.30<sup>+0.20</sup><sub>−0.13</sub></td></tr><tr><td>Radial velocity semiamplitude variation K<sub>d</sub> [m s<sup>−1</sup>]</td><td>U[0.00, 0.05]</td><td>2.03 ± 0.16</td></tr><tr><td colspan="3"><em>Derived parameters</em></td></tr><tr><td>Planet minimum mass M<sub>d</sub> sin i<sub>d</sub> [M<sub>⊕</sub>]</td><td>...</td><td>6.03 ± 0.49</td></tr><tr><td>Orbit eccentricity e<sub>d</sub></td><td>...</td><td>0.14 ± 0.09</td></tr><tr><td>Argument of periastron of stellar orbit ω<sub>d,d</sub> [deg]</td><td>...</td><td>−126<sup>+287</sup><sub>−38</sub></td></tr><tr><td colspan="3"><strong>Additional model parameters</strong></td></tr><tr><td>Characteristic period of the GP P<sub>GP</sub> [days]</td><td>U[35.0,65.0]</td><td>53.67<sup>+0.65</sup><sub>−0.53</sub></td></tr><tr><td>Inverse of the harmonic complexity λ<sub>p</sub></td><td>U[0.1,3.0]</td><td>0.44 ± 0.05</td></tr><tr><td>Long-term evolution timescale λ<sub>e</sub> [days]</td><td>U[1,300]</td><td>114 ± 19</td></tr><tr><td>Amplitude A<sub>RV</sub> [km s<sup>−1</sup>]</td><td>U[0.0,0.005]</td><td>0.0016 ± 0.0004</td></tr><tr><td>Amplitude B<sub>RV</sub> [km s<sup>−1</sup>]</td><td>U[−0.05,0.05]</td><td>0.009<sup>+0.0025</sup><sub>−0.0019</sub></td></tr><tr><td>Amplitude A<sub>DLW</sub></td><td>U[0.0,30.0]</td><td>6.31<sup>+1.42</sup><sub>−1.02</sub></td></tr><tr><td>Systemic velocity γ<sub>HARPS</sub> [km s<sup>−1</sup>]</td><td>U[47.40,48.42]</td><td>47.9168 ± 0.0005</td></tr><tr><td>Offset DLW [10<sup>3</sup> m<sup>2</sup> s<sup>−2</sup>]</td><td>U[−16.90,14.11]</td><td>−0.11<sup>+2.02</sup><sub>−1.98</sub></td></tr><tr><td>Radial velocity jitter term σ<sub>HARPS</sub> [m s<sup>−1</sup>]</td><td>J[0,100]</td><td>1.83 ± 0.08</td></tr><tr><td>DLW jitter term σ<sub>DLW</sub> [m<sup>2</sup> s<sup>−2</sup>]</td><td>J[0,1000]</td><td>378<sup>+145</sup><sub>−335</sub></td></tr></tbody></table>

Note. U[a, b] refers to uniform priors between a and b; N[a, b] refers to Gaussian priors with mean a and standard deviation b; J[a, b] refers to Jeffrey's priors between a and b. Inferred parameters and uncertainties are defined as the median and the 68.3% credible interval of their posterior distributions.

unusual architecture for systems with USP planets (Dai et al. 2018). With minimum masses of 4.13 M<sub>⊕</sub> and 6.03 M<sub>⊕</sub>, GJ 367 c and d are expected to have a minimum radius of ~1.5 R<sub>⊕</sub> given the maximum collisional stripping limit (Marcus et al. 2009, 2010). The inclination of their orbits should be larger than ~86.5° to produce nongrazing transits. If the two planets had radii of 1.5 R<sub>⊕</sub>, the transit depths would be ~900 ppm. For 2.5 R<sub>⊕</sub>, the nongrazing transit depth is expected

to be ~2500 ppm. Similarly, transit depth of a 4 R<sub>⊕</sub> planet is ~6400 ppm. The rms of the TESS light curve of GJ 367 is approximately 500 ppm. The expected transit times of planets c and d also fall well within the baseline of the TESS data. Therefore, GJ 367 c and d would be easily detected if they produced nongrazing transits.

The presence of two additional planetary companions to GJ 367 b is in line with the tendency of USP planets to belong

<footnote>a Assuming an orbital inclination of i<sub>b</sub> = 79.89<sup>+0.87°</sup><sub>−0.85°</sub>, from the modeling of TESS transit light curves (Section 5.1).</footnote>
<figure><img src="image_16.png" /><figcaption>**Figure 8.** RV and differential line width (DLW) time series with best-fitting models from the multi-GP (solid black lines) and the $1\sigma$ and $2\sigma$ credible intervals of the corresponding GP models (light gray shaded areas). The upper panel shows the RV data with the full model in black, while the planetary signal (blue), and stellar (red) inferred models are shown with a vertical offset for clarity. The lower panel shows the DLW along with the stellar inferred model. Data are shown with filled circles, with their nominal error bars and semitransparent error bar extensions accounting for the inferred jitter term.</figcaption></figure>

to multiplanet systems (Dai et al. 2021). While the origin of USP planets is debated, it is likely that the presence of outer companions could be responsible for the migration processes that carried USP planets to their current positions (Pu & Lai 2019; Petrovich et al. 2019), which would lie inside the magnetospheric cavity of a typical protoplanetary disk. In these models, planet migration occurs after the dissipation of the protoplanetary disk, by a combination of eccentricity forcing from the planetary companions and tidal dissipation in the innermost planet’s interior, although the precise dynamical forcing mechanisms are debated. Serrano et al. (2022) showed, for TOI-500, that the low-eccentricity secular forcing proposed by Pu & Lai (2019) can explain the migration of that system’s USP planet from a formation radius of ~0.02 au to its observed location. For GJ 367 b, we tested this migration scenario with similar initial conditions (USP planet starting at 0.02 au, initial eccentricities of all planets set to 0.2), and found only modest migration of planet b, from 0.02 au to ~0.01 au, short of the observed 0.007 au. This would support an alternative migration history, such as by high-eccentricity secular chaos (Petrovich et al. 2019). However, we note that GJ 367 c and GJ 367 d lie close to a 3:1 mean motion commensurability, having a period ratio of 2.98, and that the dynamics may be affected by the 3:1 resonance, which can provide additional eccentricity forcing in the system. A deeper analysis would be needed to draw definite conclusions on the formation, migration, and dynamics of the

system. The GJ 367 system is thus an excellent target for studying planetary system formation and evolution scenarios.

## 6.1. *Internal Structure of GJ 367 b*

Given the high precision of both the derived mass and radius, we used a Bayesian analysis to infer the internal structure of GJ 367 b. We followed the method described in Leleu et al. (2021), which is based on Dorn et al. (2017). Prior to presenting the results of our modeling, we provide the reader with a brief overview over the used model.

Modeling the interior of an exoplanet is a degenerate problem: there is a multitude of different compositions that could explain the observed planetary density. This is why a Bayesian modeling approach is used, with the goal of finding posterior distributions for the internal structure parameters. We assumed a planet that is made up of four fully distinct layers: an inner iron core, a silicate mantle, a water layer, and a gas layer made up of hydrogen and helium. In our forward model, this atmospheric layer is modeled separately from the rest of the planet following Lopez & Fortney (2014), and it is assumed that it does not influence the inner layers. While the presence of a gas and water layer is not expected given the high equilibrium temperature of GJ 367 b, we still included them in the initial model setup, as this is the most general way to model the planet and give us all possible compositions that could lead to the
<figure><img src="image_17.png" /><figcaption>**Figure 9.** Phase-folded RVs curve of GJ 367 b (a), GJ 367 c (b), and GJ 367 d (c). Data are shown as filled green circles with the error bars and the light-green error bars accounting for the inferred RV jitter term.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>**Figure 10.** Mass-radius diagram of small USP planets (*P* < 1 days, *R* < 2 R<sub>⊕</sub>, *M* < 10 M<sub>⊕</sub>) with masses and radii known with a precision better than 20%, as retrieved from the Transiting Extrasolar Planet Catalogue (Southworth 2011). From bottom to top, the solid and dashed curves are theoretical models from Zeng et al. (2016). We highlighted GJ 367 b with a red dot and the previous position on the diagram with a black dot (Lam et al. 2021).</figcaption></figure>

observed planetary density. As input parameters, we used the planetary and stellar parameters listed in Tables 1 and 2, i.e., the transit depth, RV semiamplitude and period of the planet

and the mass, radius, effective temperature, and metallicity of the star. In addition, we chose a prior distribution of 5 ± 5 Gyr for the age of the star.
<figure><img src="image_19.png" /><figcaption>**Figure 11.** Results of a Bayesian analysis of the internal structure of GJ 367 b. The depicted internal structure parameters are as follows: the mass fractions of the inner iron core and a possible water layer with respect to the solid planet without a possible H/He layer, the molar fractions of Si and Mg within the silicate mantle layer, the molar fraction of iron in the inner core, the logarithm with base 10 of the mass of the H/He layer in Earth masses, and the shift in the Fe/Si ratio of the planet with respect to the Fe/Si ratio in the host star, again in a logarithmic scale. The titles in each column show the median of the posterior distributions and the 5th and 95th percentiles.</figcaption></figure>

For our Bayesian analysis, we chose a prior that is uniform in log for the gas mass fraction of the planet. For the mass fractions of the inner iron core, the silicate mantle, and the water layer (all with respect to the solid part of the planet without the H/He gas layer), our chosen prior is a distribution that is uniform on the simplex on which these three mass fractions add up to one. Additionally, we added an upper limit of 0.5 for the water mass fraction in the planet, as a planet with a pure water composition is not physical (Thiabaud et al. 2014; Marboeuf et al. 2014). Note that choosing very different priors influences the results of the internal structure analysis.

There are multiple studies that show a correlation between the composition of planets and their host star (e.g., Thiabaud et al. 2015; Adibekyan et al. 2021), but it is not clear if this correlation is 1:1 or has a different slope. As a first step, we assumed the composition of the planet to match that of its host star exactly. Since we do not have measured values for the stellar [Si/H] and [Mg/H], we left them unconstrained and sampled the stellar parameters from [Si/H] = [Mg/H] = 0$^{+1}_{-1}$. However, our analysis showed that the observed mass and radius values of GJ 367 b are not compatible with a 1:1 relationship between the Si/Mg/Fe ratios of the planet and of
these sampled synthetic host stars, as it is not possible to reach
such a high density under these constraints. The same is true
when assuming the steeper slope between stellar and planetary
Fe/(Si+Mg) ratios found by Adibekyan et al. (2021). We then
repeated our analysis allowing the iron to silicon and iron to
magnesium ratios in the planet to be up to a factor 1000 higher
than in the sampled stars.

The results of this second analysis are summarized in Figure 11. Indeed, we can see that the Fe/Si ratio in the planet (and therefore also the Fe/Mg ratio) was a factor of 10 to the power of $2.01_{-0.83}^{+1.10}$ higher than in the sampled host star. We further expect GJ 367 b to host no H/He envelope and no significant water layer. Conversely, we expect the mass fraction of the inner iron core with $0.91_{-0.23}^{+0.07}$ (median with 5th and 95th percentiles) to be very high.

## 6.2. Formation and Evolution of the Ultra-high Density Sub-Earth GJ 367 b

It is not clear how a low-mass high-density planet like GJ 367 b forms. Possible pathways may include the formation out of material significantly more iron rich than thought to be normally present in protoplanetry disks. Although it is not clear if disks with such a large relative iron content specifically near the inner edge (where most of the material might be obtained from) exist (Dullemond & Monnier 2010; Aguichine et al. 2020; Adibekyan et al. 2021)

Another possibility is the formation of a larger planet with a lower bulk density. Subsequently the planet differentiates into a denser core and less dense outer layers. These outer layers are then removed. This may be accomplished through two processes.

(i) A first process is collisional stripping. The preferred removal of outer material in giant collisions (e.g., Marcus et al. 2009; Reinhardt et al. 2022) might have increased the bulk density. Indeed, this is what might have lead to the high iron fraction and therefore high density of Mercury for its small size (Benz et al. 1988, 2007), as a Mercury with a chondritic composition would have a lower density. The amount of maximum stripping is governed by the mass, impact velocity, and impact parameter of the impactor (Marcus et al. 2010; Reinhardt & Stewart 2012). Preferable removal of outer layers requires the right mass ratio (close to unity), right impact parameter (close to grazing), and right relative velocity. There is also evidence that, at least in some systems, densities have been altered by collisions (Kepler-107; Bonomo et al. 2019) A problem to effectively remove mass might be the removal of collision debris and the avoidance of a re-accretion of debris material onto the planet, as re-accretion would leave the bulk density largely unchanged. However, this might not be as large a problem as originally thought (Spalding & Adams 2020). Our measurement of the bulk density of GJ 367 b suggests that collisional stripping has to be remarkably effective in removing non-iron material from the planet if it is the only process at work.

(ii) A second process that might have played a role in shaping GJ 367 b after core formation is removal of outer layers of material facilitated by the enormous stellar radiation to which this planet is subjected. At the equilibrium temperature of 1365 ± 32 K, material might sublimate, be uplifted, and transported away from the surface.

Of course, all of the above discussed processes could have contributed to create the nearly pure ball of iron, known as GJ 367 b.

## 7. Conclusions

We report refined mass and radius determinations for GJ 367 b, the ultra-high density, USP sub-Earth planet transit-ing the M-dwarf star GJ 367 recently discovered by Lam et al. (2021). We used new TESS observations from sectors 35 and 36 combined with 371 Doppler measurements collected as part of an intensive RV follow-up campaign with the HARPS spectrograph. We derived a precise planetary mass of $M_b = 0.633 \pm 0.050 M_⊕$ (7.9% precision) and a radius of $R_b = 0.699 \pm 0.024 R_⊕$ (3.4% precision), resulting in an ultra-high density of $\rho_b = 10.2 \pm 1.3 \text{ g cm}^{-3}$ (~13%). According to our internal structure analysis, GJ 367 b is predominantly com-posed of iron with an iron core mass fraction of $0.91_{-0.23}^{+0.07}$, which accounts for the aforementioned planetary density. In addition, our HARPS RV follow-up observations, which span a period of nearly ~3 yr, allowed us to discover two additional nontransiting small companions with orbital periods of ~11.5 and 34.4 days, and minimum masses of $M_c \sin i_c = 4.13 \pm 0.36 M_⊕$ and $M_d \sin i_d = 6.03 \pm 0.49 M_⊕$. GJ 367 joins the small group of well-characterized multiplanetary systems hosting a USP planet, with the inner planet GJ 367 b being the densest and smallest USP planet known to date. This unique multiplanetary system hosting this ultra-high density, USP sub-Earth is an extraordinary target to further investigate the formation and migration scenarios of USP systems.

## Acknowledgments

We thank the anonymous referee for a thoughtful review of our paper and very positive feedback. This work was supported by the KESPRINT collaboration, an international consortium devoted to the characterization and research of exoplanets discovered with space-based missions (www.kesprint.science). Based on observations made with the ESO-3.6 m telescope at La Silla Observatory under programs 1102.C-0923 and 106.21TJ.001. We are extremely grateful to the ESO staff members for their unique and superb support during the observations. We acknowledge the use of public TESS data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center. TESS data presented in this paper were obtained from the Milkulski Archive for Space Telescopes (MAST) at the Space Telescope Science Institute. The specific observations analyzed can be accessed via doi:10.17909/rkvw-t847. Resources supporting this work were provided by the NASA High-End Computing (HEC) Program through the NASA Advanced Supercomputing (NAS) Division at Ames Research Center for the production of the SPOC data products. E.G. acknowledges the generous support from Deutsche Forschungsgemeinschaft (DFG) of the grant HA3279/14-1. D.G. gratefully acknowledges financial support from the *Cassa di Risparmio di Torino* (CRT) foundation under Grant No. 2018.2323 "Gaseous or rocky? Unveiling the nature of small worlds". Y.A. and J.A.E. acknowledge the support of the Swiss National Fund under grant 200020_192038. C.M.P. gratefully acknowledges the support of the Swedish National Space Agency (DNR 65/19). R.L. acknowledges funding from University of La Laguna through the Margarita Salas Fellowship from the Spanish Ministry of Universities ref. UNI/551/2021-May 26, and under the EU Next Generation funds. K.W.F.L. was supported by Deutsche Forschungsgemeinschaft grants RA714/14-1 within the DFG Schwerpunkt SPP 1992, Exploring the Diversity of Extrasolar Planets. S.A. acknowledges the
support from the Danish Council for Independent Research through grant No. 2032-00230B. O.B. received funding from the ERC under the European Union's Horizon 2020 research and innovation program (grant agreement No. 865624). H.J.D. acknowledges support from the Spanish Research Agency of the Ministry of Science and Innovation (AEI-MICINN) under grant 'Contribution of the IAC to the PLATO Space Mission' with reference PID2019-107061GB-C66. A.J.M. acknowledges support from the Swedish National Space Agency (career grant 120/19C). O.K. acknowledges support by the Swedish Research Council (grant agreement No. 2019-03548), the Swedish National Space Agency, and the Royal Swedish Academy of

Sciences. J.K. gratefully acknowledges the support of the Swedish National Space Agency (SNSA; DNR 2020-00104) and of the Swedish Research Council (VR: Etableringsbidrag 2017-04945).

## Appendix

We report here the HARPS RVs, including those previously reported in Lam et al. (2021), along with the activity indicators and line profile variation diagnostics extracted with NAIRA and serval (Table A1). Time stamps are given in Barycentric Julian Date in the Barycentric Dynamical Time (BJD<sub>TDB</sub>).
<table><caption>Table A1
Radial Velocities and Spectral Activity Indicators</caption><thead><tr><th>BJD<sub>TBD</sub><br/>-2450000</th><th>RV<br/>(km s<sup>-1</sup>)</th><th>σ<sub>RV</sub><br/>(km s<sup>-1</sup>)</th><th>Hα<br/>(Ang)</th><th>σ<sub>Hα</sub><br/>(Ang)</th><th>Hβ<br/>(Ang)</th><th>σ<sub>Hβ</sub><br/>(Ang)</th><th>Hγ<br/>(Ang)</th><th>σ<sub>Hγ</sub><br/>(Ang)</th><th>log R'<sub>HK</sub></th><th>σ<sub>log R'</sub><sub>HK</sub></th><th>NaD<br/>(Ang)</th><th>σ<sub>NaD</sub><br/>(Ang)</th><th>DLW<br/>(m<sup>2</sup> s<sup>-2</sup>)</th><th>σ<sub>DLW</sub><br/>(m<sup>2</sup> s<sup>-2</sup>)</th><th>CRX<br/>(m s<sup>-1</sup>)</th><th>σ<sub>CRX</sub><br/>(m s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>8658.45473</td><td>47.91637</td><td>0.00084</td><td>0.06472</td><td>0.00012</td><td>0.05217</td><td>0.00026</td><td>0.11412</td><td>0.00081</td><td>-5.222</td><td>0.077</td><td>0.01038</td><td>0.00007</td><td>-7.9</td><td>1.3</td><td>-10.2</td><td>11.1</td></tr><tr><td>8658.46195</td><td>47.91814</td><td>0.00082</td><td>0.06421</td><td>0.00011</td><td>0.05180</td><td>0.00025</td><td>0.11147</td><td>0.00078</td><td>-5.180</td><td>0.077</td><td>0.01040</td><td>0.00007</td><td>-10.0</td><td>1.2</td><td>3.4</td><td>9.3</td></tr><tr><td>8658.46946</td><td>47.91502</td><td>0.00085</td><td>0.06424</td><td>0.00012</td><td>0.05180</td><td>0.00026</td><td>0.11413</td><td>0.00084</td><td>-5.221</td><td>0.077</td><td>0.01045</td><td>0.00007</td><td>-6.8</td><td>1.2</td><td>-17.2</td><td>10.7</td></tr><tr><td>8658.47642</td><td>47.91611</td><td>0.00079</td><td>0.06417</td><td>0.00011</td><td>0.05257</td><td>0.00025</td><td>0.11363</td><td>0.00078</td><td>-5.130</td><td>0.077</td><td>0.01041</td><td>0.00007</td><td>-7.1</td><td>1.3</td><td>2.8</td><td>10.9</td></tr><tr><td>8658.48401</td><td>47.91860</td><td>0.00090</td><td>0.06454</td><td>0.00012</td><td>0.05315</td><td>0.00028</td><td>0.11124</td><td>0.00088</td><td>-5.139</td><td>0.077</td><td>0.01040</td><td>0.00008</td><td>-7.5</td><td>1.4</td><td>3.5</td><td>11.0</td></tr><tr><td>8658.49105</td><td>47.91862</td><td>0.00087</td><td>0.06499</td><td>0.00012</td><td>0.05417</td><td>0.00027</td><td>0.11592</td><td>0.00086</td><td>-5.136</td><td>0.077</td><td>0.01054</td><td>0.00007</td><td>-6.1</td><td>1.3</td><td>4.4</td><td>10.9</td></tr><tr><td>8658.49842</td><td>47.91780</td><td>0.00086</td><td>0.06548</td><td>0.00012</td><td>0.05403</td><td>0.00027</td><td>0.11601</td><td>0.00088</td><td>-5.151</td><td>0.077</td><td>0.01074</td><td>0.00007</td><td>-7.9</td><td>1.1</td><td>-11.5</td><td>10.4</td></tr><tr><td>8658.50565</td><td>47.91753</td><td>0.00091</td><td>0.06500</td><td>0.00012</td><td>0.05374</td><td>0.00029</td><td>0.11864</td><td>0.00092</td><td>-5.137</td><td>0.077</td><td>0.01066</td><td>0.00008</td><td>-5.0</td><td>1.5</td><td>9.4</td><td>11.2</td></tr><tr><td>8658.51297</td><td>47.91784</td><td>0.00089</td><td>0.06466</td><td>0.00012</td><td>0.05302</td><td>0.00028</td><td>0.11529</td><td>0.00091</td><td>-5.175</td><td>0.077</td><td>0.01057</td><td>0.00008</td><td>-9.3</td><td>1.3</td><td>-15.4</td><td>10.6</td></tr><tr><td>8658.52049</td><td>47.91557</td><td>0.00083</td><td>0.06464</td><td>0.00011</td><td>0.05272</td><td>0.00026</td><td>0.11484</td><td>0.00085</td><td>-5.151</td><td>0.077</td><td>0.01069</td><td>0.00007</td><td>-5.7</td><td>1.3</td><td>6.2</td><td>11.1</td></tr><tr><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr></tbody></table>

Note. The entire RV data set is available in its entirety in machine-readable form. Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.
(This table is available in its entirety in machine-readable form.)
ORCID iDs

Elisa Goffo https://orcid.org/0000-0001-9670-961X

Davide Gandolfi https://orcid.org/0000-0001-8627-9628

Jo Ann Egger https://orcid.org/0000-0003-1628-4231

Alexander J. Mustill https://orcid.org/0000-0002-2086-3642

Simon H. Albrecht https://orcid.org/0000-0003-1762-8235

Teruyuki Hirano https://orcid.org/0000-0003-3618-7535

Oleg Kochukhov https://orcid.org/0000-0003-3061-4591

Nicola Astudillo-Defru https://orcid.org/0000-0002-8462-515X

Oscar Barragan https://orcid.org/0000-0003-0563-0493

Luisa M. Serrano https://orcid.org/0000-0001-9211-3691

Artie P. Hatzes https://orcid.org/0000-0002-3404-8358

Yann Alibert https://orcid.org/0000-0002-4644-8818

Eike Guenther https://orcid.org/0000-0002-9130-6747

Fei Dai https://orcid.org/0000-0002-8958-0683

Kristine W. F. Lam https://orcid.org/0000-0002-9910-6088

Szilárd Csizmadia https://orcid.org/0000-0001-6803-9698

Alexis M. S. Smith https://orcid.org/0000-0002-2386-4341

Luca Fossati https://orcid.org/0000-0003-4426-9530

Rafael Luque https://orcid.org/0000-0002-4671-2957

Florian Rodler https://orcid.org/0000-0003-0650-5723

Mark L. Winther https://orcid.org/0000-0003-1687-3271

Jakob L. Rørsted https://orcid.org/0000-0001-9234-430X

William D. Cochran https://orcid.org/0000-0001-9662-3496

Hans J. Deeg https://orcid.org/0000-0003-0047-4241

Jon M. Jenkins https://orcid.org/0000-0002-4715-9460

Judith Korth https://orcid.org/0000-0002-0076-6239

John H. Livingston https://orcid.org/0000-0002-4881-3620

Annabella Meech https://orcid.org/0000-0002-7500-7173

Felipe Murgas https://orcid.org/0000-0001-9087-1245

Jaume Orell-Miquel https://orcid.org/0000-0003-2066-8959

Hannah L. M. Osborne https://orcid.org/0000-0002-4143-4767

Enric Palle https://orcid.org/0000-0003-0987-1593

Carina M. Persson https://orcid.org/0000-0003-1257-5146

Seth Redfield https://orcid.org/0000-0003-3786-3486

George R. Ricker https://orcid.org/0000-0003-2058-6662

Sara Seager https://orcid.org/0000-0002-6892-6948

Roland Vanderspek https://orcid.org/0000-0001-6763-6562

Vincent Van Eylen https://orcid.org/0000-0001-5542-8870

Joshua N. Winn https://orcid.org/0000-0002-4265-047X

References

Adams, E. R., Jackson, B., Endl, M., et al. 2017, AJ, 153, 82

Adibekyan, V., Dorn, C., Sousa, S. G., et al. 2021, Sci, 374, 330

Aguichine, A., Mousis, O., Devouard, B., & Ronnet, T. 2020, ApJ, 901, 97

Aguirre Børsen-Koch, V., Rørsted, J. L., Justesen, A. B., et al. 2022, MNRAS, 509, 4344

Aigrain, S., Llama, J., Ceillier, T., et al. 2015, MNRAS, 450, 3211

Anderson, D. R., Cameron, A. C., Hellier, C., et al. 2010, ApJL, 726, L19

Astudillo-Defru, N., Delfosse, X., Bonfils, X., et al. 2017a, A&A, 600, A13

Astudillo-Defru, N., Forveille, T., Bonfils, X., et al. 2017b, A&A, 602, A88

Azevedo Silva, T., Demangeon, O. D. S., Barros, S. C. C., et al. 2022, A&A, 657, A68

Bagnulo, S., Landolfi, M., Landstreet, J. D., et al. 2009, PASP, 121, 993

Barnes, S. A. 2010, ApJ, 722, 222

Barnes, S. A., & Kim, Y.-C. 2010, ApJ, 721, 675

Barragán, O., Aigrain, S., Rajpaul, V. M., & Zicher, N. 2022, MNRAS, 509, 866

Barragán, O., Gandolfi, D., & Antoniciello, G. 2019, MNRAS, 482, 1017

Barragán, O., Gandolfi, D., Dai, F., et al. 2018, A&A, 612, A95

Barros, S. C. C., Demangeon, O. D. S., Alibert, Y., et al. 2022, A&A, 665, A154

Benz, W., Anic, A., Horner, J., & Whitby, J. A. 2007, SSRv, 132, 189

Benz, W., Slattery, W. L., & Cameron, A. G. W. 1988, Icar, 74, 516

Binks, A. S., & Jeffries, R. D. 2014, MNRAS, 438, L11

Binks, A. S., Jeffries, R. D., Jackson, R. J., et al. 2021, MNRAS, 505, 1280

Bonomo, A. S., Zeng, L., Damasso, M., et al. 2019, NatAs, 3, 416

Brandner, W., Calissendorff, P., Frankel, N., & Cantalloube, F. 2022, MNRAS, 513, 661

Cabrera, J., Csizmadia, S., Erikson, A., Rauer, H., & Kirste, S. 2012, A&A, 548, A44

Cifuentes, C., Caballero, J. A., Cortés-Contreras, M., et al. 2020, A&A, 642, A115

Claret, A. 2017, A&A, 600, A30

Dai, F., Howard, A. W., Batalha, N. M., et al. 2021, AJ, 162, 62

Dai, F., Masuda, K., & Winn, J. N. 2018, ApJL, 864, L38

Donati, J. F., Morin, J., Petit, P., et al. 2008, MNRAS, 390, 545

Donati, J. F., Semel, M., Carter, B. D., Rees, D. E., & Collier Cameron, A. 1997, MNRAS, 291, 658

Dorn, C., Venturini, J., Khan, A., et al. 2017, A&A, 597, A37

Dullemond, C. P., & Monnier, J. D. 2010, ARA&A, 48, 205

Essack, Z., Burt, J., Shporer, A., et al. 2022, BAAS, 54, 102.59

Gaidos, E., Claytor, Z., Dungee, R., Ali, A., & Feiden, G. A. 2023, MNRAS, 520, 5283

Gandolfi, D., Barragán, O., Hatzes, A. P., et al. 2017, AJ, 154, 123

Guenther, E. W., Barragán, O., Dai, F., et al. 2017, A&A, 608, A93

Gaia Collaboration, Smart, R. L., Sarro, L. M., et al. 2021, A&A, 649, A6

Hatzes, A. P. 2019, in The Doppler Method for the Detection of Exoplanets, ed. Vishniac (Bristol: IOP Publishing)

Hatzes, A. P., Dvorak, R., Wuchterl, G., et al. 2010, A&A, 520, A93

Hatzes, A. P., Fridlund, M., Nachmani, G., et al. 2011, ApJ, 743, 75

Hidalgo, S. L., Pietrinferni, A., Cassisi, S., et al. 2018, ApJ, 856, 125

Hirano, T., Dai, F., Gandolfi, D., et al. 2018, AJ, 155, 127

Hirano, T., Livingston, J. H., Fukui, A., et al. 2021, AJ, 162, 161

Jenkins, J. M., Twicken, J. D., McCauliff, S., et al. 2016, Proc. SPIE, 9913, 99133E

Kiman, R., Faherty, J. K., Cruz, K. L., et al. 2021, AJ, 161, 277

Klein, B., Zicher, N., Kavanagh, R. D., et al. 2022, MNRAS, 512, 5067

Kochukhov, O. 2021, A&ARv, 29, 1

Kochukhov, O., Makaganiuk, V., & Piskunov, N. 2010, A&A, 524, A5

Kochukhov, O., & Reiners, A. 2020, ApJ, 902, 43

Koen, C., Kilkenny, D., van Wyk, F., & Marang, F. 2010, MNRAS, 403, 1949

Lam, K. W. F., Csizmadia, S., Astudillo-Defru, N., et al. 2021, Sci, 374, 1271

Lee, E. J., & Chiang, E. 2017, ApJ, 842, 40

Leinhardt, Z. M., & Stewart, S. T. 2012, ApJ, 745, 79

Leleu, A., Alibert, Y., Hara, N. C., et al. 2021, A&A, 649, A26

Lopez, E. D., & Fortney, J. J. 2014, ApJ, 792, 1

Lovis, C., & Pepe, F. 2007, A&A, 468, 1115

Mamajek, E. E., & Bell, C. P. M. 2014, MNRAS, 445, 2169

Mandel, K., & Agol, E. 2002, ApJL, 580, L171

Mann, A. W., Dupuy, T., Kraus, A. L., et al. 2019, ApJ, 871, 63

Mann, A. W., Feiden, G. A., Gaidos, E., Boyajian, T., & von Braun, K. 2015, ApJ, 804, 64

Marboeuf, U., Thiabaud, A., Alibert, Y., Cabral, N., & Benz, W. 2014, A&A, 570, A36

Marcus, R. A., Sasselov, D., Hernquist, L., & Stewart, S. T. 2010, ApJL, 712, L73

Marcus, R. A., Stewart, S. T., Sasselov, D., & Hernquist, L. 2009, ApJL, 700, L118

Marcy, G. W., Isaacson, H., Howard, A. W., et al. 2014, ApJS, 210, 20

Mayor, M., Pepe, F., Queloz, D., et al. 2003, Msngr, 114, 20

Millholland, S. C., & Spalding, C. 2020, ApJ, 905, 71

Morin, J., Donati, J. F., Petit, P., et al. 2008, MNRAS, 390, 567

Murdoch, K. A., Hearnshaw, J. B., & Clark, M. 1993, ApJ, 413, 349

Otegi, J. F., Bouchy, F., & Helled, R. 2020, A&A, 634, A43

Pace, G. 2013, A&A, 551, L8

Paegert, M., Stassun, K. G., Collins, K. A., et al. 2022, yCat, IV/39

Petrovich, C., Deibert, E., & Wu, Y. 2019, AJ, 157, 180

Piskunov, N., Snik, F., Dolgopolov, A., et al. 2011, Msngr, 143, 7

Piskunov, N. E., & Valenti, J. A. 2002, A&A, 385, 1095

Pu, B., & Lai, D. 2019, MNRAS, 488, 3568

Rajpaul, V., Aigrain, S., Osborne, M. A., Reece, S., & Roberts, S. 2015, MNRAS, 452, 2269

Reinhardt, C., Meier, T., Stadel, J. G., Otegi, J. F., & Helled, R. 2022, MNRAS, 517, 3132

Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2015, JATIS, 1, 014003

Rodriguez, J. E., Becker, J. C., Eastman, J. D., et al. 2018, AJ, 156, 245
Rusomarov, N., Kochukhov, O., Piskunov, N., et al. 2013, A&A, 558, A8

Ryabchikova, T., Piskunov, N., Kurucz, R. L., et al. 2015, PhyS, 90, 054005

Sahu, K. C., Casertano, S., Bond, H. E., et al. 2006, Natur, 443, 534

Sanchis-Ojeda, R., Rappaport, S., Winn, J. N., et al. 2014, ApJ, 787, 47

Santerne, A., Brugger, B., Armstrong, D. J., et al. 2018, NatAs, 2, 393

Schlaufman, K. C., Lin, D. N. C., & Ida, S. 2010, ApJL, 724, L53

Serrano, L. M., Gandolfi, D., Mustill, A. J., et al. 2022, NatAs, 6, 736

Smith, J. C., Stumpe, M. C., Van Cleve, J. E., et al. 2012, PASP, 124, 1000

Snik, F., Kochukhov, O., Piskunov, N., et al. 2011, in ASP Conf. Ser. 437, Solar Polarization 6, ed. J. R. Kuhn et al. (San Francisco, CA: ASP), 237

Southworth, J. 2011, MNRAS, 417, 2166

Spalding, C., & Adams, F. C. 2020, PSJ, 1, 7

Stassun, K. G., Oelkers, R. J., Paegert, M., et al. 2019, AJ, 158, 138

Stassun, K. G., Oelkers, R. J., Pepper, J., et al. 2018, AJ, 156, 102

Stumpe, M. C., Smith, J. C., Catanzarite, J. H., et al. 2014, PASP, 126, 100

Stumpe, M. C., Smith, J. C., Van Cleve, J. E., et al. 2012, PASP, 124, 985

Thiabaud, A., Marboeuf, U., Alibert, Y., et al. 2014, A&A, 562, A27

Thiabaud, A., Marboeuf, U., Alibert, Y., Leya, I., & Mezger, K. 2015, A&A, 580, A30

Toledo-Padrón, B., Lovis, C., Suárez Mascareño, A., et al. 2020, A&A, 641, A92

Triaud, A. H. M. J. 2018, in Handbook of Exoplanets, ed. H. J. Deeg & J. A. Belmonte (Springer, Cham: Springer), 2

Wildi, F., Pepe, F., Chazelas, B., Lo Curto, G., & Lovis, C. 2010, Proc. SPIE, 7735, 77354X

Wildi, F., Pepe, F., Chazelas, B., Lo Curto, G., & Lovis, C. 2011, Proc. SPIE, 8151, 81511F

Winn, J. N. 2010, arXiv:1001.2010

Winn, J. N., & Fabrycky, D. C. 2015, ARA&A, 53, 409

Winn, J. N., Sanchis-Ojeda, R., & Rappaport, S. 2018, NewAR, 83, 37

Yee, S. W., Petigura, E. A., & von Braun, K. 2017, ApJ, 836, 77

Zechmeister, M., & Kürster, M. 2009, A&A, 496, 577

Zechmeister, M., Reiners, A., Amado, P. J., et al. 2017, A&A, 609, A12

Zeng, L., Sasselov, D. D., & Jacobsen, S. B. 2016, ApJ, 819, 127

Zhu, W., & Dong, S. 2021, ARA&A, 59, 291

Zicher, N., Barragán, O., Klein, B., et al. 2022, MNRAS, 512, 3060