EDP Sciences
Free Access
Issue
A&A
Volume 605, September 2017
Article Number L11
Number of page(s) 7
Section Letters
DOI https://doi.org/10.1051/0004-6361/201731581
Published online 22 September 2017

© ESO, 2017

1. Introduction

In the last 20 yr, we have learned much about the vast diversity in the architecture of planetary systems (e.g., Mayor et al. 2014; Fabrycky et al. 2014). For stars of different spectral types we have discovered planets with a a large variety in mass and size, constraining our understanding of planetary formation processes.

M dwarfs make up about 70% of the stellar population of the Galaxy and constitute the lower tail of the main sequence in the Hertzsprung-Russell diagram. The detection of planets orbiting these stars is therefore important to constrain the population of planets and the formation processes in the less-massive protoplanetary disks (Gaidos 2017). Additionally, low-mass and small-sized M dwarfs have advantages when searching for the lowest-mass and smallest-sized planets: the radial velocity (RV) amplitude and transit depth scales with and , respectively.

It is known that M dwarfs have a high occurrence of super-Earths, as shown from high-precision radial velocity (f ~ 0.88 for P< 100 days and Msini = 1−10 M, Bonfils et al. 2013) and transit surveys (Dressing & Charbonneau 2015, f ~ 2.5 for P< 200 days and R = 1−4 R,). Thanks to the most precise velocimeters combined with the high amount of collected data, we are now starting to census planets with minimum masses closer to the Earth mass such as GJ 273 c (Astudillo-Defru et al. 2017b), Prox Cen b (Anglada-Escudé et al. 2016), or GJ 628 b (Wright et al. 2016).

The detection of multi-planetary systems sheds light on the theories of model formation and evolution. Statistical studies have shown that planet pairs near or in mean resonances are slightly favored (e.g., Fabrycky et al. 2014). In the understanding of this scenario, a low uncertainty in the measured mass is important, but the majority of the systems has been detected through transits, where uncertainties are generally too large for a reliable physical characterization of the planets (e.g., Kepler-42 and Trappist-1, Muirhead et al. 2012; Gillon et al. 2017, respectively).

In this paper we present the new nearby system (3.6 pc) of a very-low-mass star that is orbited by at least three Earth-mass planets on compact orbits, with the possibility of a fourth sub-Earth-mass companion.

2. Stellar properties and observations

2.1. Stellar properties

The properties of the mid-M dwarf YZ Cet are given in Table 1. The table lists the visual and near-infrared apparent magnitudes, the parallax and proper motions, basic physical parameters, the activity level, and the age. These values come from the literature as detailed in the table caption.

Table 1

Stellar properties.

2.2. HARPS radial velocities

We observed YZ Cet from December 2003 to October 2016 with the HARPS spectrograph installed at the ESO 3.6 m telescope in Chile (Mayor et al. 2003). A total of 211 spectra were acquired with an exposure time of 900 s and fast read-out mode, corresponding to a typical signal-to-noise ratio per spectra of 19 at 550 nm. Data were recorded without simultaneous wavelength reference to keep stellar low-flux zones free of lamp contamination. The HARPS vacuum vessel was opened in May 2015 to upgrade the fiber link (Lo Curto et al. 2015, a superscript plus indicates throughout data acquired after the upgrade). This modifies the line-spread function and may affect the spectral analysis. Consequently, for some of our analyses, we independently analyze the pre- and post-fiber upgrade data sets.

2.3. ASAS photometry

This work uses data from the All Sky Automated Survey (ASAS; Pojmanski 1997). The ASAS project constantly monitors a large area of the sky (9° × 9°) to search for photometric variability in V and I band. The ASAS photometric accuracy is 0.01–0.02 mag for a 13th magnitude star, which is enough to detect spots crossing the stellar disk due to rotation.

The ASAS data of YZ Cet consist of 459 photometric measurements acquired from November 2000 to November 2009. Unfortunately, ASAS and HARPS data do not overlap in time, which is desired to understand RV variability that is not produced by planets. We used the two-pixel wide aperture (MAG_0), as suggested by ASAS for a 12.07 Vmag star.

3. Analysis

3.1. Periodogram+Keplerian

thumbnail Fig. 1

Periodogram analysis of the RV (top panel) and the ASAS photometry and FWHMCCF (bottom panel). Dashed lines highlight the strongest peaks, dotted lines represent the 1% FAP threshold. Three RV signals have an FAP lower than 1%, while a fourth is just above this level of confidence. The final RV residual (omc k4) does not show any significant periodicity. The photometry and FWHMCCF show a ~83 day periodicity linked to the stellar rotation. Pre- (red) and post- (blue) fiber upgrade data are analyzed separately.

Open with DEXTER

We computed the generalized Lomb-Scargle periodogram (Zechmeister & Kürster 2009) of the RVs and residues by progressively subtracting a Keplerian fit accounting for the strongest peak in the periodogram (Fig. 1). We adjusted up to four Keplerian curves plus a velocity offset between the data obtained before and after the HARPS fibre upgrade (BJD = 2 457 172, Lo Curto et al. 2015). The top panel in Fig. 1 shows three clear periodic RV signals at 4.7, 3.1, and 2.0 days, with a false-alarm probability (FAP) <1/10 000. A fourth less clear signal is seen with a one-day periodicity and a 1.1% FAP. The Bayesian information criterion (BIC) for a constant model is 857.83, and for a model of one, two, three, and four Keplerian, the BIC is 463.57, 420.21, 392.62, and 376.10, respectively. On the one hand, the BIC favors the model consisting of four Keplerian curves, but on the other hand, the FAP of the fourth signal is above the 1% FAP threshold. We prefer to stay on the conservative side and consider the fourth signal as tentative. A more rigorous Bayesian model comparison and/or additional data are needed to elucidate the reliability of the one-day RV signature.

3.2. Stellar activity

YZ Cet is an active star with . Such a level translates into a rotation period of 27 ± 9 days (Astudillo-Defru et al. 2017a). However, there are clues that flaring stars like Proxima Centauri do not follow the relationship. This is the case for YZ Cet, where several flares are identified in spectra through Ca II hk, Na D, He id3, and Hαβγδ. These chromospheric activity indicators do not show evidence of stellar rotation.

Contrary to the chromospheric activity indicators, the V-band photometry and the FWHM show evidence of the rotation period of YZ Cet (~83 days, Fig. 1). Conversely, none of the periodicities detected in the RV time series appears in the photometry or FWHM data.

3.3. Three Keplerian+Gaussian process

As the star rotates, quasi-periodic signals may be introduced in RV as a result of evolving active regions. These signals are often very complex, and a physical model is therefore usually hard to produce or has the strong disadvantage of being computationally expensive. Gaussian process regression (GP) has become a widely used non-parametric method to model this variability without the necessity of specifying a model (Rasmussen & Williams 2006; Haywood et al. 2014; Rajpaul et al. 2015; Faria et al. 2016), although it may be prone to inducing false positives if care is not taken (Dumusque et al. 2017).

A GP with a quasi-periodic covariance function was added to the three-Keplerian model. The GP hyperparameters were sampled together with the model parameters. We set a normal prior on the recurrence timescale, (see Appendix A), centered at 83 days and with a width of 15 days, in agreement with the measured stellar rotation rate (see Sect. 3.2). A detailed description of the model and the choice of priors is provided in Appendix A.

The 23-dimensional parameter space of our model was explored with the Markov chain Monte Carlo (MCMC) algorithm described in Goodman & Weare (2010) and implemented by Foreman-Mackey et al. (2013), where 100 walkers were used. The algorithm was started at a random point in the prior and was allowed to evolve until no further variation in mean posterior density across the walkers was seen. Then, the algorithm was run for a further 100 000 steps used for the inference process. We retained only those walkers that explored the main period peak for each planet (Fig. 1), whose samples where used to initialize a new 70 000-step MCMC run using 300 walkers. The final inference was made on 2.1 × 107 posterior samples, out of which over 21 000 are independent (see Appendix B for details). Results are reported in Table B.3.

4. Results

4.1. Compact system near the 3:2 resonances

The periodogram and activity analysis (Sects. 3.1, 3.2) both give strong evidence that the very-low-mass star YZ Cet is orbited by at least three planets in a compact architecture. From the Keplerian modeling combined with the Gaussian processes, we obtain orbital periods of 4.66, 3.06, and 1.97 days, that is, orbits close to the 3:2 mean-motion resonances. Considering the stellar mass of 0.13M, the measured semi-amplitude translates into terrestrial-mass planets with minimum masses of 1.14, 0.98, and 0.75 Earth-masses, respectively. Table 2 gives the basic parameters of the system. Figure 2 shows the RV folded to each period.

Table 2

Basic orbital parameters of the planets.

thumbnail Fig. 2

Radial velocities folded to each maximum a posteriori orbital period (planets b, c, and d from top to bottom, respectively). Gray shaded regions extend between the 5th and 95th percentile of the model at each orbital phase, computed over 10 000 samples of the posterior. Orange and blue circles depict data before and after the HARPS fibre upgrade, respectively.

Open with DEXTER

4.2. Dynamical modeling

We show in Fig. 3 a stability analysis focused on planet c where all orbital parameters are fixed to the maximum a posteriori (MAP) solution, and we varied only the semi-major axis and the eccentricity. For each set of initial conditions, we integrated the system for 1 kyr using GENGA (Grimm & Stadel 2014), and computed the NAFF stability indicator (Laskar 1990, 1993; Correia et al. 2005). This showed that the MAP solution is stable and outside the two 3:2 mean-motion resonances. We also see that ec ≲ 0.1 is required for the system to remain stable.

thumbnail Fig. 3

Stability analysis around the MAP solution (marked with a black cross, see Table B.3). Blue points correspond to stable and red points to unstable solutions. White points correspond to highly unstable solutions, i.e., with a collision or the ejection of a planet during the 1 kyr of the simulation. The 3:2 resonances between planets b and c and between planets c and d are highlighted.

Open with DEXTER

Additionally, we used the n-body code REBOUND (Rein & Liu 2012) with the WHFast integrator (Rein & Tamayo 2015) to obtain the RV of the bodies in time, which we compared with the HARPS data. A stability criterion was imposed based on an additional n-body integration of 1 kyr in the future. The change in period of a planet in 1 kyr (ΔP1 kyr) was extrapolated to 5 Gyr (the estimated age of the system) assuming . The model was rejected when ΔP5 Gyr was longer than the measured orbital period, which is a conservative criterion rejecting unstable configurations. We used a normal distribution as prior for the stellar mass (Table 1) and uniform priors for the remaining parameters. The emcee algorithm (Goodman & Weare 2010; Foreman-Mackey et al. 2013) was used to sample the posterior distribution of the parameters. We estimated the planet masses without the degeneracy with an orbital inclination relative to the line of sight by modeling the gravitational interaction between the three well-known planets in the system. From this preliminary analysis, until the fourth signal is confirmed (or rejected), we obtain , , and (median and 68.3% credible interval), thanks to the constraint on the inclination: 95% HDI ib[10, 172]°, ic[6, 174]°, and id[28, 174]°. Tighter constraints are also obtained for the eccentricities of the measured planets, eb< 0.39, ec< 0.15, and ed< 0.17 at 95%.

5. Conclusions

We presented the discovery of three terrestrial-mass planets orbiting the very-low-mass star YZ Cet. The best fit (three Keplerian+GP) to our 13 yr of HARPS data results in compact orbits close to the 3:2 mean-motion resonance with orbital periods of 1.97, 3.06, and 4.66 days. The measured minimum masses are 0.75, 0.98, and 1.14 Earth-masses, respectively, with a typical uncertainty of 16%. A fourth planet with a minimum mass of 0.42 ± 0.11 M and an orbital period of 1.0419 days may be present in the system (1.1%FAP). To date, these are the lowest minimum masses measured through radial velocity.

From a dynamical analysis, we infer the upper limits of the planetary masses, resulting in masses below three Earth-masses. It is probable that these have radii up to 1.5 R (Fulton et al. 2017) with a rocky bulk composition. The planets orbit outside the habitable zone of YZ Cet, with equilibrium temperatures ranging from 347–491 K, 299–423 K, and 260–368 K for planets b, c, and d, respectively.

The system is at only 3.6 pc, making it very attractive for further characterization, especially if any of the planets is caught transiting its host star, which is not the case so far.

Acknowledgments

N.A.D. and F.P. acknowledges the support of the Swiss National Science Foundation project 166227. N.C.S. acknowledges support from Fundação para a Ciência e a Tecnologia (FCT) through national funds and from FEDER through COMPETE2020 by the grants UID/FIS/04434/2013 & POCI–01–0145-FEDER–007672 and PTDC/FIS-AST/1526/2014 & POCI–01–0145-FEDER–016886. N.C.S. also acknowledges the support from FCT through Investigador FCT contract IF/00169/2012/CP0150/CT0002. We made use of the REBOUND code, which can be downloaded at http://github.com/hannorein/rebound. These simulations have been run on the Regor cluster kindly provided by the Observatoire de Genève. X.B., J.M.A. and A.W. acknowledge funding from the European Research Council under the ERC Grant Agreement No. 337591-ExTrA.

References

Appendix A: Gaussian process model

A quasi-periodic kernel was chosen to construct the elements of the covariance matrix of the observations. This type of covariance function was used in the past to model the effect of stellar active regions rotating in and out of view as the star revolves (Haywood et al. 2014; Rajpaul et al. 2015).

The elements of the covariance matrix are then given by , where represents the RV effect related to the stellar rotational modulation, with the hyperparameters associated with the covariance amplitude, the evolution timescale, the recurrence timescale (the rotational period of the star), and the smoothing parameter, respectively. The second term is a diagonal matrix with the quadratic sum of the internal velocity errors (σi) and an additional white-noise component: where Si is an indicator variable that is one if observation i is taken before the HARPS fibre upgrade and zero otherwise, and vice versa for . δij is the Kronecker delta function.

The prior distribution of the hyperparameters can be factorized in the marginal prior distributions described in Table B.2. Only the recurrence timescale was assigned an informative prior based on the analysis of the ASAS photometric measurements and the FWHM time series.

Appendix B: Details on the MCMC analysis

The ensemble sampler of Goodman & Weare (2010) implemented in the emcee package (Foreman-Mackey et al. 2013) was used to sample the joint posterior distribution of the parameters and hyperparameters of the model. The priors assigned to the model parameters are listed in Tables B.1 and B.2.

Initially, we chose to use 100 walkers on a 23-dimensional parameter space. The starting position of each walker was randomly drawn from the prior distribution. The burn-in period took approximately 500 000 steps, and the walkers were allowed to evolve still for around 100 000 steps. An issue usually encountered when running an MCMC algorithm is that some walkers become stuck in the local maxima, in this case, related with period aliases and other features seen in the periodogram (Fig. 1).

Because of this, we decided to retain only those walkers that explored the main periodogram peak for each planet. In this manner, 47 out of the 100 walkers were kept. The samples produced by these subsets of walkers were used to initialize a second MCMC run using 300 walkers, which were run for 70 000 steps. This run had a mean acceptance fraction between walkers of around 25%. The characteristic autocorrelation length was between 3 × 102 and 1 × 103 steps. The average autocorrelation length over walkers was shorter than 1000 steps for all parameters. We have then more than 700 independent samples per walker, which leads to a total of around 21 000 independent samples that are used for the parameter inference.

Table B.1

Prior distribution for the Keplerian curves and velocity offset between HARPS pre- and post-fibre upgrade.

Table B.2

Prior distribution for the hyperparameters of the Gaussian process kernel function.

In Table B.3 we list some characteristics of the posterior distribution based on the MCMC samples. Figure B.1 shows the 2D marginal distributions of a selection of model parameters.

Table B.3

Orbital elements using a model of three Keplerian plus a Gaussian process.

Table B.4

Radial velocity time series for YZ Cet (minimal; full version available at the CDS).

thumbnail Fig. B.1

One- and two-dimensional marginal posterior distributions of selected model parameters. The contours represent the 68.27%, 95.45%, and 99.73% joint credible regions.

Open with DEXTER

All Tables

Table 1

Stellar properties.

Table 2

Basic orbital parameters of the planets.

Table B.1

Prior distribution for the Keplerian curves and velocity offset between HARPS pre- and post-fibre upgrade.

Table B.2

Prior distribution for the hyperparameters of the Gaussian process kernel function.

Table B.3

Orbital elements using a model of three Keplerian plus a Gaussian process.

Table B.4

Radial velocity time series for YZ Cet (minimal; full version available at the CDS).

All Figures

thumbnail Fig. 1

Periodogram analysis of the RV (top panel) and the ASAS photometry and FWHMCCF (bottom panel). Dashed lines highlight the strongest peaks, dotted lines represent the 1% FAP threshold. Three RV signals have an FAP lower than 1%, while a fourth is just above this level of confidence. The final RV residual (omc k4) does not show any significant periodicity. The photometry and FWHMCCF show a ~83 day periodicity linked to the stellar rotation. Pre- (red) and post- (blue) fiber upgrade data are analyzed separately.

Open with DEXTER
In the text
thumbnail Fig. 2

Radial velocities folded to each maximum a posteriori orbital period (planets b, c, and d from top to bottom, respectively). Gray shaded regions extend between the 5th and 95th percentile of the model at each orbital phase, computed over 10 000 samples of the posterior. Orange and blue circles depict data before and after the HARPS fibre upgrade, respectively.

Open with DEXTER
In the text
thumbnail Fig. 3

Stability analysis around the MAP solution (marked with a black cross, see Table B.3). Blue points correspond to stable and red points to unstable solutions. White points correspond to highly unstable solutions, i.e., with a collision or the ejection of a planet during the 1 kyr of the simulation. The 3:2 resonances between planets b and c and between planets c and d are highlighted.

Open with DEXTER
In the text
thumbnail Fig. B.1

One- and two-dimensional marginal posterior distributions of selected model parameters. The contours represent the 68.27%, 95.45%, and 99.73% joint credible regions.

Open with DEXTER
In the text

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

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

Initial download of the metrics may take a while.