© ESO, 2017

1. Introduction

The term “stellar activity” indicates a class of phenomena that, in stars with an external convective envelope, are triggered by the reconfiguration of the surface magnetic field. For the Sun, these phenomena are commonly classified as sunspots, plages, flares, or coronal holes, for example, and they take place from the photosphere to the corona. The portions of the stellar surface and atmosphere where these phenomena take place are referred to as active regions (ARs).

From the observational point of view, traditionally the most frequently investigated diagnostic of chromospheric activity for solar-type stars is the Ca ii H&K doublet (3933.67 Å and 3968.47 Å; see, e.g., the milestone project of the Mount Wilson Observatory as summarized in Baliunas et al. 1998). The main advantage of this doublet is that the photospheric line profiles have large depths with respect to the continuum, such that the chromospheric emission can be detected with high contrast against the photospheric background. Nonetheless, a particularly relevant drawback especially for M-type stars is that detailed spectroscopy is observationally more demanding because of the low continuum level at wavelengths shortward of 4000 Å. For this reason, these observations have been focusing on stars ranging from F to K spectral types.

This led to the identification of another diagnostic more suited for red dwarfs: the Hα line (6562.80 Å; see, e.g., Pasquini & Pallavicini 1991; Montes et al. 1995). It is generally accepted that the Ca ii H&K and Hα excess fluxes are tightly correlated even though they form in different layers of the chromosphere. However, this conclusion has been drawn historically by analyzing averaged measurements of chromospheric fluxes of calcium and hydrogen taken at different epochs and rarely obtained simultaneously (Giampapa et al. 1989; Robinson et al. 1990; Strassmeier et al. 1990). Another approach is to collect simultaneous measurements for a number of stars at a given epoch, without any follow-up investigation of the observed targets at different epochs (Thatcher & Robinson 1993; Walkowicz & Hawley 2009).

Recently, planet search programs carried out with the radial velocity techniques have started to monitor samples of M dwarfs. These stars are extremely interesting targets for planet discoveries. First of all, they represent ~75% of the stars in the solar neighborhood (Reid et al. 2002; Henry et al. 2006). Moreover, from an observational point of view, the chances of finding an Earth-like planet in the habitable zone of a star increase as the stellar mass decreases. Still, habitability is not guaranteed simply by an assessment of the distance from the star; several other factors, such as stellar activity, may move and/or shrink the habitability zone of a star (see Vidotto et al. 2013, and references therein). Thus, it is crucial to better understand the activity of M dwarfs and how it can affect the circumstellar environment.

Intensive monitoring aimed at an accurate measurement of radial velocities has provided the community with large databases of high-resolution and high signal-to-noise ratio (S/N) optical spectra per star. These databases are valuable for the analysis of the activity of M dwarfs. For instance, Gomes da Silva et al. (2011) have analyzed the long-term variability of a number of chromospheric optical indexes for a sample of 30 M0–M5.5 stars from the HARPS at ESO M-dwarf planet search program (Bonfils et al. 2013) with a median timespan of observations of 5.2 yr.

The HArps-N red Dwarf Exoplanet Survey (HADES, Affer et al. 2016; Perger et al. 2017) project is a collaborative program between the Global Architecture of Planetary Systems project¹ (GAPS, Covino et al. 2013), the Institut de Ciències de l’Espai (IEEC/CSIC), and the Instituto de Astrofísica de Canarias (IAC), that monitors the radial velocities of a sample of low-activity M-type dwarfs to search for planets. In this framework, we are analyzing the collected spectra to study the activity of the monitored stars from different points of view. In Maldonado et al. (2017) we analyze the Ca ii H&K and Balmer line flux excesses in relation to rotation and age, while in Suárez Mascareño et al. (in prep.) we analyze the time series of spectra in order to measure the rotation period of the stars and the periodicity of their activity cycles, if any is present.

In this paper our aim is to characterize the daily to monthly variability due to chromospheric ARs based on the analysis of simultaneous measurements of the Ca ii H&K and Hα lines. In Sect. 2 we present the target selection and the database of spectra collected so far; in Sect. 3 we describe the data reduction and the measurement of the chromospheric Ca ii H&K and Hα flux excesses; in Sect. 4 we analyze the flux-flux relationships and in Sect. 5 we discuss the amount and timescales of chromospheric variability in the selected sample of stars.

2. Target selection

The HADES sample of stars is made up of 78 red dwarfs, 71 of which have been observed so far. The observed stars have spectral types ranging between K7.5V and M3V (corresponding to the ~3400–3900 K temperature range). These stars are being monitored with the HARPS-N spectrograph (Cosentino et al. 2012) mounted at the Telescopio Nazionale Galileo in the framework of the HADES project. The spectrograph covers the 383–693 nm wavelength range with a spectral resolution of ~115 000. HARPS-N spectra were reduced using the most recent version of the Data Reduction Software (DRS) pipeline (Lovis & Pepe 2007).

The monitored stars have been selected as targets favorable for planet search, thus the sample is biased towards low activity levels with some exceptions. The program stars are listed in Table A.1, together with the number of observations for each of them. The data analyzed in this paper have been collected over seven semesters, from September 2012 to February 2016. In Maldonado et al. (2017) we measured the stellar parameters using the methodology described in Maldonado et al. (2015). The most relevant stellar parameters for the current analysis are listed in Table A.1. We refer to Affer et al. (2016) for further details of the target selection and to Maldonado et al. (2017) for a more detailed description of the determination of the stellar parameters.

3. Post-processing analysis

During the spectral analysis, we find evidence of variability in the instrumental response, in the atmospheric reddening, and in the intensity of telluric lines, either in absorption or emission, from night to night. These effects are usually corrected by observing a spectroscopic standard star and a telluric standard star close to the observed target. The aim of the HADES and GAPS projects is to find new exoplanets around pre-selected targets or to refine the characterization of known exoplanetary systems. Thus, the instrumental setup and the observational campaign are arranged to maximize the number of targeted stars and the precision of the radial velocity measurements, while the flux calibration of the spectra is not taken into account (e.g., spectroscopic or telluric standard stars are not usually observed).

For our purposes, however, the correction of all these effects is necessary in order to compare the Ca ii H&K and Hα lines taken on different nights. In this section we describe in detail our post-processing analysis of the spectra reduced with the DRS pipeline.

3.1. Flux rescaling

To correct the observed fluxes for instrumental response and atmospheric reddening, for each star in our catalog we interpolate the synthetic spectral library provided by Allard et al. (2011)² to compute a model spectrum corresponding to the measured T_eff, log (g), and [Fe/H]. We arbitrarily degrade the spectral resolution of both the model spectrum and the series of observed spectra convolving them with a Gaussian kernel with σ = 120 Å, corresponding to a final spectral resolution of ~50 over the spectral range in consideration (Fig. 1, top panel). The ratio of each low-resolution spectrum and the model is thus a continuous function of wavelength (Fig. 1, middle panel), which allows us to rescale the high-resolution instrumental fluxes to fluxes at stellar surface and, at the same time, to correct for the differential instrumental response and atmospheric reddening (Fig. 1, bottom panel).

The degradation to such a low spectral resolution was motivated by several factors:

• it allows us to safely mask out the narrow spectral windowscentered on the chromospheric emission lines(Table 1) that may introduce spurious features inthe correction function;

• it removes the line component of the telluric spectrum, leaving us only with the continuous component of the atmospheric absorption;

• at low resolution, the discrepancies between observed and model line profiles are negligible, the low-resolution spectra being dominated by the continuum flux.

3.2. Telluric line removal

The other effect to take into account is the telluric contamination. We focus only on the 130 Å-wide spectral range around the Hα line because the Ca ii H&K spectral region is free from telluric lines.

In order to have a comparison spectrum suited for the removal of telluric lines, for each star we compute the median of all the flux-rescaled high-resolution spectra. This median spectrum has a higher S/N than the single spectra, and is also free from telluric lines. Actually, telluric lines are Doppler shifted in the reference frame of the targeted star because the observations are generally performed at different hour angles and different months along the observing seasons. Since their wavelengths change in the stellar rest frame from night to night, the median algorithm filters them out as outliers.

For each star we compute the ratio between each spectrum and the corresponding median value, obtaining the normalized telluric spectrum affecting each given observation. Then we cross-correlate it with the model telluric spectrum provided in Borsa et al. (2015). While the cross-correlation allows us to align this model telluric spectrum to the ratio of the spectra, the intensity of the telluric contamination is computed by a 3-σ clipped comparison between the wavelength-shifted model telluric spectrum and the ratio. These two parameters (wavelength shift and intensity) allow us to correct the telluric lines within the noise of the observed spectra.

In the following sections, we always mean the flux-rescaled telluric-corrected spectra when referring to the observed spectra.

3.3. Selection of the reference spectra

For each star, we have a series of individual spectra and the corresponding median spectrum (Sects. 3.1 and 3.2). To make the computation faster, we extract a 70 Å spectral interval around the Hα and a 100 Å interval around the Ca ii H&K doublet comprising both lines (Table 1).

In each interval, we want to subtract the photospheric flux to measure the chromospheric flux. To this end, our first approach was to use synthetic spectra, which soon turned out to loosely fit the absorption line profiles. We then decided to use our spectral database to obtain the best “inactive” template for each star in the program.

It is reasonable to assume that those stars with minimum Ca ii H&K emission are the “least active” and therefore the ones that should be used as templates for the spectral subtraction. For each star we measure the Ca ii H&K flux ${\mathit{F}}_{\mathrm{HK}}^{\mathrm{int}}$ integrating the median spectrum in the spectral intervals around the Ca ii H and K lines defined in Table 1 (see the end of Sect. 3.4 for the definition of the spectral intervals). These integrated fluxes thus contain the contribution of both the photospheric flux and chromospheric basal flux, plus the flux from chromospheric ARs, if present. In Fig. 2 we plot the integrated Ca ii H&K flux versus the stellar effective temperature; we find that the lower boundary of the diagram is an increasing function of T_eff. This trend is due to the increase in the Ca ii H&K photospheric flux of early-type M dwarfs with T_eff, as predicted by the BT-Settl models.

We select the template stars at the lower boundary (i.e., at minimum Ca ii H&K flux) of the diagram in Fig. 2, arbitrarily considering the stars below the line with equation $\mathrm{Log}\mathrm{\left(}{\mathit{F}}_{\mathrm{HK}}^{\mathrm{int}}\mathrm{\right)}\mathrm{=}\mathrm{-}\mathrm{2.61}\mathrm{+}\mathrm{0.0021}\mathrm{·}{\mathit{T}}_{\mathrm{eff}}\mathit{.}$(1)The selected reference stars are V^∗ BR Psc, GJ 15A, GJ 412A, GJ 4196, GJ 3997, GJ720A, and GJ 4092. All these stars have subsolar metallicities, with the exception of GJ 4196. The selection of these stars could reflect a trend between Ca ii flux and metallicity that is due to the dependence of the Ca ii line core strength on the calcium elemental abundance. This correlation could therefore be linked to the line core strength curve of growth in the atmosphere discussed by Houdebine (2011).

As we discuss below, stellar metallicity cannot be neglected in the spectral subtraction. For this reason, we reject GJ 4196 as a reference star, and we keep the remaining stars as a spectral grid for the computation of the reference quiet spectrum.

3.4. Spectral subtraction and line flux measurement

For each program star, the reference spectrum is computed by interpolating the grid of reference spectra over T_eff. Before the interpolation takes place, we broaden the spectra such that they all have the same rotational broadening. For this purpose, we select the maximum vsin i among those of the program and reference stars as measured by Maldonado et al. (2017), and we broaden each spectrum to this maximum vsin i using the rotational profile provided by Gray (1992).

Since the reference spectrum is computed using the spectra of stars with minimum metallicities in our sample, we find a systematic offset (see Fig. 3) correlated with the metallicity of the star. In particular, this occurs because the continuum flux predicted by the BT-Settl models decreases with increasing [Fe/H], consistently with the increase in the absorption by the multitude of atomic and molecular lines in M star spectra. In Fig. 3 we show the example of the Hα spectrum of GJ 9793 (see also Fig. 4), the highest metallicity star in our sample: its spectrum shows lower fluxes than those of the corresponding reference spectrum. To correct this offset, we fit a low-order polynomial curve to the flux points. In the fit, we exclude the core of the line (gray crosses in the figure) because it may be affected by chromospheric emission. Finally, the reference spectrum is rescaled to the observed spectrum using the best-fit curve.

Once the reference spectrum is rescaled to the observed spectrum, the chromospheric spectrum is computed by subtraction of the template (see also Herbig 1985; Frasca & Catalano 1994, for previous application of the spectral subtraction technique). Outside of the line cores we sometimes find a residual wavelength-dependent offset significantly different from zero. This residual offset arises from several factors, such as inaccuracies in the measured stellar parameters, the interpolation of the reference spectrum, and spectral peculiarities of the analyzed spectrum. To remove these systematics, we subtract a low-order polynomial function from the chromospheric spectrum, fitted over the spectral range excluding the line cores. We note that the 1σ confidence band on the polynomial fit is narrow enough to introduce a negligible uncertainty on the line flux compared to the noise (Fig. 4, bottom panel).

Finally, the chromospheric emission is computed by integrating the chromospheric spectrum, i.e., the template-subtracted systematics-corrected spectrum, over the wavelength ranges listed in Table 1. The widths of the spectral windows used to compute the emission flux are conveniently set after visual inspection of the differential spectra, representing a compromise between the bracketing of the whole flux radiated by the chromosphere and the exclusion of the nearby spectral regions that introduce noise in the flux measurement. To compute the uncertainty on the integrated flux, we measure the S/N of the differential spectrum outside of the line cores, and we propagate it in the integration of the line flux. This error thus includes both the noise in the observed spectra and the uncertainties introduced by our data reduction.

We caution that the spectral subtraction removes the flux of both the photosphere and the “quiet” chromosphere from the observed spectra. For the Ca ii H&K doublet, this means that the measured emission fluxes do not include the basal chromospheric emission, i.e., they are representative of the excess chromospheric flux radiated by ARs. The case of the Hα line is more subtle. Following the chromospheric models of Cram & Giampapa (1987), the absorption in the Hα line seen in the spectra of early-type M dwarfs originates in the chromosphere. In particular, in these models even the unperturbed chromospheres are hot enough to populate the N = 2 levels of the hydrogen ions, making the absorption of the Hα line possible. Thus, in the least active M stars in our sample we see a “basal” Hα absoption depth, and the line fluxes measured with the spectral subtraction technique have this basal absorption as zero point.

4. Flux-flux relationships

In the following we analyze the relationships between the chromospheric fluxes in the Ca ii H&K and Hα lines resolved in time. To improve the statistical reliability of our analysis, we arbitrarily discard all the stars with fewer than 20 observations, which leaves us with 41 stars. Moreover, to clean up our samples of measurements, for each star and each line we reject the outliers with a 5σ-clipping criterion. We apply the same selection to the measurement uncertainties, thus discarding the data with anomalously large error bars.

4.1. Ca ii H vs. Ca ii K

In Fig. 5 we compare the individual measurements of the Ca ii H and K excess fluxes (F_CaIIH and F_CaIIK, respectively) of the stars in the analyzed sample, omitting the measurements less than 3σ above zero because they are not representative of flux in excess. We find that the fluxes are well aligned despite differences in stellar parameters and activity levels. Since the uncertainties on the quantities are comparable, we compute the best-fit line by means of the ranged major axis (RMA) regression (Legendre & Legendre 1983). The result of the fit is $\log {\mathit{F}}_{{\mathrm{Ca}\hspace{0.17em}{ii}}_{\mathrm{K}}}\mathrm{=}\mathrm{(}\mathrm{0.11}\mathrm{\pm }\mathrm{0.04}\mathrm{)}\mathrm{+}\mathrm{(}\mathrm{0.990}\mathrm{\pm }\mathrm{0.009}\mathrm{)}\log {\mathit{F}}_{{\mathrm{Ca}\hspace{0.17em}{ii}}_{\mathrm{H}}}\mathit{,}$(2)which is consistent within 3σ with similar previous studies performed on samples of M dwarfs (Martínez-Arnáiz et al. 2011; Stelzer et al. 2012). From Eq. (2) we derive that the K/H flux ratio is close to 1, and following Houdebine & Stempels (1997) we argue that the emission of the Ca ii H&K lines is optically thick.

We analyze the Ca ii H&K flux excess as a function of temperature and metallicity. For this purpose, in Fig. 6 we plot the median total flux excess ⟨ F_HK ⟩ = ⟨ F_H + F_K ⟩ of the stars, together with the corresponding median absolute deviation³ (MAD) of the sample of measurement, versus T_eff and [Fe/H].

We find that the flux excess drops for stars cooler than T_eff = 3600 K (spectral type M1). Similar results have been found by West et al. (2004), who analyzed the Hα emission, and by West et al. (2011), who included the Ca ii K line, but for stars later than M4. We cannot claim any stronger conclusion because we are investigating a narrow spectral range. Moreover, we note that our subsample of 41 stars is biased towards low activity levels, thus our findings are strongly affected by selection effects.

The dependence on the stellar metallicity is slightly more convincing: more metallic stars tending towards larger ⟨ F_HK ⟩ excess. One possible explanation is that the Ca ii H&K emission increases with the abundance of Ca ii ions in the chromosphere, thus with metallicity at first approximation, as also discussed by Houdebine (2011). Another possibility is that this could be an age effect, as higher metallicities generally correspond to younger ages. This is consistent with the fact that younger stars tend to be more active, as we discuss in Sect. 5.1.

4.2. Ca ii H&K vs. Hα

In the left panel of Fig. 7 we show the emission flux F_Hα radiated in the Hα line plotted against the flux F_HK radiated by the Ca ii H&K doublet. The comparison between Hα and Ca ii H&K fluxes is more scattered than that of Ca ii H and K lines. To avoid cluttering, in the middle panel we show the same diagram for the median of the flux measurements ⟨ F_HK ⟩ and ⟨ F_Hα ⟩, where the plotted uncertainties represent the MAD of the individual flux measurements.

One possible explanation of the scatter in Fig. 7 is that, because they have a larger number of Ca ii ions in the atmosphere, stars with higher metallicities have more efficient Ca ii H&K radiative cooling (i.e., higher fluxes in emission), as found by Houdebine (2011). To test this hypothesis, we split the sample of stars into two halves at [Fe/H] = −0.14, and plot them with different colors in the left panel of Fig. 7; we find that the two subsamples follow the same locus. In particular, we do not find that the F_HK fluxes tend to increase with metallicity at fixed F_Hα. We thus conclude that we do not have enough statistical evidence to state that metallicity affects the Hα vs. Ca ii H&K excess relationship, maybe because the metallicity range spanned by our program stars is not large enough (see Table A.1).

We perform the same analysis splitting the sample of stars in two halves on the basis of their T_eff (Fig. 8). Also in this case we do not find any evidence of correlation, maybe because we are investigating a narrow range in T_eff.

Another interesting feature, visible in the middle panel of Fig. 7, is that the flux-flux relationship is not monotonic. A visual analysis of the plot shows that the flux radiated in the Hα line seems to initially decrease with increasing Ca ii H&K excess flux, then goes below zero at ⟨ F_HK ⟩ ≃ 10^-5 erg cm^-2 s^-1 (i.e., the Hα line deepens), and finally increases with the Ca ii H&K flux. Similar results were obtained by Robinson et al. (1990). The absorption of the Hα line towards negative ⟨ F_Hα ⟩ values is the reason why we plot the Ca ii H&K vs. Hα relationship in linear rather than logarithmic scale.

The non-monotonic trend is supported by statistical tests. Kendall’s correlation test⁴ significantly supports the idea that the two datasets are correlated. Moreover, we naively fit a linear and a quadratic function to the data, and we find that the corrected Akaike Information Criterion⁵ significantly rejects the former compared to the latter.

From the theoretical point of view, this result has been predicted by chromospheric models (Cram & Mullan 1979; Cram & Giampapa 1987; Rutten et al. 1989; Houdebine et al. 1995; Houdebine & Stempels 1997). According to these models, the most feasible scenario is that the chromospheric Ca ii H&K emission lines are collisionally dominated, thus the radiated flux steadily increases with pressure. Conversely, the Hα line is radiation-dominated and the increase of the optical depth initially leads to a deeper absorption profile, until the electron density is high enough to take the Hα line into the collisionally dominated formation regime, leading to the fill-in of the line.

Previous observational studies performed on similar samples of stars (K- and M-type main sequence stars, e.g., Stauffer & Hartmann 1986; Giampapa et al. 1989; Rauscher & Marcy 2006; Walkowicz & Hawley 2009) have found comparable results. In particular, while the deepening of the Hα line is still somehow controversial, it is clear that the Hα is not filled in until a certain Ca ii H&K emission flux is reached. Thus, the Hα line is not a good activity indicator for low to intermediate activity levels. Moreover, Houdebine & Stempels (1997) noted that observational data do not strictly match the models likely owing to surface inhomogeneities, which may also explain the scatter in Fig. 7.

A careful analysis of the left panel of Fig. 7 also shows that the stars follow individual loci with different slopes in the diagram. To better show this effect, in the right panel of Fig. 7 we plot the locus occupied by each single star defined as the RMA best fit (see Sect. 4.1) of its data points.

First of all, we check whether the scatter of each star in the diagram is due to an ongoing long-term activity cycle or to variability on shorter timescales. In this regard, for each star we visually inspect the Ca ii H&K and Hα time series, and in some cases we find hints of increasing or decreasing activity over year-long timescales or even turnarounds (see Fig. 9 for an example). Robertson et al. (2013) found similar evidence analyzing the ~11-yr series of Hα observations of a sample of 93 M-type dwarfs, and concluded that ~10% of M dwarfs show either long-term trends or activity cycles that are longer than one year.

To remove the above-mentioned trends, for each season we reject outliers (generally positive outliers likely due to flares) and we compute the median of the remaining measurements. Then, for each season we subtract the corresponding median flux, thus filtering out any variability longer than one year. We find that long-term variability is generally able to explain up to 40% of the scatter in the emission line fluxes. We thus conclude that ≳60% of the variance in the data is due to daily to weekly to monthly variability.

In support of this conclusion, we find that slopes, correlation coefficients and significances are generally preserved by the correction for long-term trends, which confirms that the displacement of the stars in Fig. 7 is mainly due to correlated shorter term variability. In the following, we always refer to the long-term corrected flux excess.

In Fig. 10 we plot the distribution function of the slopes of the residual F_HK vs. F_Hα relationship for those stars that have slopes with high statistical significance (p-value < 1%). This distribution shows that the slope is always positive, i.e., F_Hα is always an increasing function of F_HK.

Moreover, the distribution function is not symmetric: it has a broad right tail. The Anderson-Darling normality test (Anderson & Darling 1952) supports the non-normality of the sample distribution at the level of ~4 × 10^-3. We thus infer that the distribution of the data is not dominated by measurement errors (which would lead to a Gaussian distribution function), but reflects an intrinsic characteristic of the analyzed sample.

We also test the dependence of the slopes of the F_HK vs. F_Hα relationship on ⟨ F_HK ⟩, T_eff, and [Fe/H], respectively, and we find evidence of decreasing steepness of the flux-flux relationship with increasing T_eff with a confidence level of ~4% (Fig. 11). Stelzer et al. (2013) found a similar trend by analyzing a sample of pre-MS low-mass stars in the 2500–4500 K temperature range. Conversely, Kendall’s correlation test cannot reject the hypothesis that the slope of the flux-flux relationship is independent of the average level of activity or the metallicity of the star. On the other hand, we remark that we cannot claim the independence from the aforementioned parameters as a lack of statistical dependence may be a consequence of the narrow ranges of the investigated stellar parameters.

We compare these results with the model provided by Meunier & Delfosse (2009), who explain the F_HK vs. F_Hα relationship for the Sun. In their work, the authors use spatially resolved observations of the solar surface covering 1.5 activity cycles. They find that the presence of dark filaments affects the correlation between F_HK and F_Hα on short timescales (daily or weekly), sometimes leading to negative correlations depending on the surface coverage of plages and filaments. Analogous phenomena may affect the distribution in Fig. 10 leading to a non-Gaussian shape, and may also explain why we obtain negative slopes in some cases, even though with low statistical significance (see Fig. 7, gray lines in the right panel). In particular, in their model, larger values for the slopes correspond to a scenario in which the chromospheric Hα emission is less affected by the absorption of cold filaments than the rest of the sample is. For a fixed configuration of plages seen in the Ca ii H&K lines on the Sun, the slope of the flux-flux relationship depends on the contrast of the Hα emission in the plages and the absorption by filaments, and their respective filling factors. For a fixed contrast of Hα plages on the solar disk, the slope of the F_HK vs. F_Hα relationship tends to decrease as the contrast of the filaments increases. This is explained by the fact that solar filaments are optically thicker in the Hα than in the Ca ii H&K lines. Thus, the increase in slope with decreasing T_eff (Fig. 11) may indicate that the Hα absorption in equivalent active regions decreases towards later M types. Analyzing the Balmer decrements of the same sample of stars we independently draw the same conclusion (Maldonado et al. 2017), consistently with similar previous studies (Bochanski et al. 2007; Stelzer et al. 2013).

5. Chromospheric variability

5.1. Variability vs. activity vs. age

In Fig. 12 we plot the standard deviation of the residual Ca ii H&K and Hα flux measurements (see Sect. 4.2) of each program star (${\mathit{\sigma }}_{\mathrm{res}}^{\mathrm{HK}}$ and ${\mathit{\sigma }}_{\mathrm{res}}^{\mathrm{H}\mathit{\alpha }}$ respectively) as a function of the median level of Ca ii H&K flux, assumed to monotonically increase with the level of magnetic activity (see Sect. 4). In these diagrams we find a tight correlation (supported by the Spearman’s correlation test) between variability and median excess flux. The linear best fits suggest that the fractional variability with respect to the average Ca ii H&K flux excess (which is represented by the slope of the linear fit) is roughly 10% and that even at the minimum level of activity there is a residual variability of ~0.5 × 10^-5 erg cm^-2 s^-1, which is different from zero at the ~2σ level for the Ca ii H&K lines (see Table 2). If this is a real effect rather than a bias introduced by our analysis, it is evidence that even at minimum activity there is some variability.

We also separate the stars younger and older than ~650 Myr as classified in Maldonado et al. (2017) based on their kinematics (Table A.1). Although kinematics provide only lower limits to stellar ages (Maldonado et al. 2010, and references therein), with this classification we find evidence that the two populations of stars show significantly different distributions in terms of Ca ii H&K emission flux. The Kolmogorov-Smirnov test on the ⟨ F_HK ⟩ we measure confirms this result (top panel in Fig. 12), indicating that younger stars have higher Ca ii H&K fluxes. The same test on ⟨ F_Hα ⟩ is less conclusive probably because the Hα is seen in absorption for a range of ⟨ F_HK ⟩, i.e., the Hα is not a good diagnostic at low and intermediate activity levels (see Sect. 4.2).

The Kolmogorov-Smirnov test on the ${\mathit{\sigma }}_{\mathrm{res}}^{\mathrm{HK}}$ and ${\mathit{\sigma }}_{\mathrm{res}}^{\mathrm{H}\mathit{\alpha }}$ does not provide any firm indication that the amount of variability is age-dependent, i.e., variability is not a good age estimator.

It can be argued that the increase of variability in the Ca ii H&K and Hα residual emission fluxes may lower the degree of correlation in the flux-flux relationship. On the contrary, we find that larger variability generally corresponds to higher levels of Ca ii H&K activity (Fig. 12, middle panel), which in turn are associated with stronger correlation (Pearson’s coefficient r_res closer to 1, smaller p-value_res) of the flux-flux relationship (Fig. 13). Conversely, the quietest stars tend to show high p-value_res, and in general have lower or negative correlation coefficients r_res. The same results were also found by Gomes da Silva et al. (2011) by analyzing the long-term activity of a sample of M dwarfs. This supports the assumption that line emission in the Ca ii H&K doublet and the Hα line are triggered by spatially correlated ARs.

5.2. Variability timescales

To analyze the typical timescales of the variability, we perform a time series analysis using the pooled variance (PV) approach described by Donahue et al. (1997a,b).

The PV is a measure of the average variance in the data over the timescale τ. By construction, the PV computed at τ is expected to be the combination of contributions from a variety of sources, instrumental and/or astrophysical, with timescales shorter than τ. When τ increases, the PV remains constant until the effects from processes with longer timescales become more noticeable, in which case the PV increases with τ. The PV approach is thus suited for time series containing multiple periodic signals with different amplitudes and phases, and is indicated to study the timescales of variability induced by stellar rotation and the life cycle of ARs.

The PV computation needs a large number of data points to be robust (Lanza et al. 2004). In 3.5 yr of observations we have collected up to 100 spectra for a single star. Since this is not a particularly abundant data set, we apply the PV analysis only to the eight stars with more than 90 observations, namely GJ 2, GJ 16, GJ 3942, GJ 625, GJ 4306, GJ 694.2, GJ 3998, and GJ 49, and show their PVs in Figs. 14 to 21.

In the computation of the PV diagrams, we consider the Ca ii H&K and Hα flux excesses corrected for the long-term trend (Sect. 4.2), hence we are confident that any contribution from long-term activity cycles is reduced to a minimum. We use the algorithm described in Donahue et al. (1997a); the only difference is that we use the median and the MAD instead of the mean and the standard deviation, respectively, the former being more robust to the presence of outliers especially for small sample sizes.

We complement the spectroscopic data with the V-band photometric monitoring of the same sample of stars, carried out at the INAF-Catania Astrophysical Observatory in the same epochs of the spectroscopic monitoring using an 80 cm telescope. The complete and detailed analysis of the photometric monitoring will be discussed in a future paper.

In general, for each star we find that the Ca ii H&K, Hα, and V-band PV diagrams look alike. The PV diagrams generally increase at small timescales τ; they reach a plateau at roughly τ ≃10–40 days, and then increase to level off again at τ ≳ 50 days.

Our interpretation of these diagrams is that the first plateau corresponds to the stellar rotation period. This is confirmed by a complementary Lomb-Scargle based analysis⁶ of the Ca ii H&K and Hα fluxes we are carrying out to measure the rotation period of the stars (see Suárez Mascareño et al., in prep.). The increase at longer τ indicates that there is an additional source of variance in the data with timescales longer than the rotation period. The PV then levels off at τ ~ 50–60 days, suggesting that the new source of variance typically has timescales of a few stellar rotations. We exclude that the second flattening is due to cyclic activity for two reasons: dynamo cycles of M dwarfs act on year-long timescales, as indicated by Robertson et al. (2013), and their effects should be greatly reduced by the long-term correction we applied in Sect. 4.2. We thus infer that this source of variance is likely related to the growth and decay of chromospheric ARs, as also suggested by Donahue et al. (1997a,b). Other studies have found that M dwarfs may show variability on timescales on the order of 100 days (Davenport et al. 2015; Robertson et al. 2015; Newton et al. 2016), proposing that it is related to the lifecycle of ARs.

6. Summary

In this study we analyze the spectroscopic database collected in the framework of the HADES project, focusing on the characterization of the variability of Ca ii H&K and Hα lines as chromospheric diagnostics. The database consists of the intensive spectroscopic monitoring performed with the HARPS-N spectrograph at TNG for a sample of 71 low-activity early-type M dwarfs, to search for planets.

We focus on the simultaneous analysis of the flux excess emitted in the Ca ii H&K and Hα lines. For this purpose, we developed a technique used to correct the spectra for instrumental and atmospheric effects by means of synthetic spectra. This leads to the calibration of the spectra to a common absolute flux scale, which enables the measurement of line flux excesses in units of flux at the stellar surface.

Our measurements show that the Ca ii H and K flux excesses are strongly linearly correlated, consistent with previous results (see Martínez-Arnáiz et al. 2011, and references therein). When comparing the Ca ii H&K with the Hα chromospheric line flux we find significantly more scatter, larger than the measurement uncertainties. In the F_HK vs. F_Hα diagram we also find statistical evidence of a non-linear flux-flux relationship. In particular, supported by the results of Houdebine et al. (1995) and Houdebine & Stempels (1997), we argue that the Ca ii H&K emission flux increases monotonically with the stellar activity level, while the Hα line initially goes into absorption and then is filled in by radiative emission processes. This suggests that at very low activity levels the Hα absorption by filaments is more evident, while at higher activity levels the emission by plages dominates.

Searching for the physical origin of the scatter and trends in the F_HK vs. F_Hα diagram, we find that long-term activity cycles on year-long timescales play a minor role in the overall variability. Conversely, short-term activity explains ≳60% of the variance of the flux excess measurements.

After correcting the collected time series for year-long variability (whose analysis is not the scope of the present study), we find marginal statistical evidence that the slope of the F_HK vs. F_Hα relationship increases with decreasing stellar T_eff. This is consistent with a scenario in which cooler stars tend to be less affected by chromospheric filaments as suggested by Meunier & Delfosse (2009), and as we found in Maldonado et al. (2017) with an independent analysis of the Balmer decrements of the same sample of stars.

We also find that the variance of the flux excess is an increasing function of stellar activity, and that even the quietest stars show some degree of variability. We attempt a rough time series analysis using the pooled variance approach suggested by Donahue et al. (1997a,b). We find evidence for rotation periods on the order of 10–40 days and active region lifetime cycles longer than ≃50 days. These findings are in agreement with previous results (Donahue et al. 1997a,b; Reiners et al. 2012; Robertson et al. 2013) and are supported by a more detailed analysis, which we are currently preparing (Suárez Mascareño et al., in prep.).

---

Acknowledgments

G.S. and I.P. acknowledge financial support from “Accordo ASI–INAF” No. 2013-016-R.0 July 9, 2013. J.M. acknowledges support from the Italian Ministry of Education, University, and Research through the PREMIALE WOW 2013 research project under grant “Ricerca di pianeti intorno a stelle di piccola massa”. G.A.P.S. acknowledges support from INAF through the Progetti Premiali funding scheme of the Italian Ministry of Education, University, and Research. J.I.G.H. acknowledges financial support from the Spanish Ministry of Economy and Competitiveness (MINECO) under the 2013 Ramón y Cajal program MINECO RYC-2013-14875, and A.S.M., J.I.G.H., and R.R. also acknowledge financial support from the Spanish ministry project MINECO AYA2014-56359-P. I.R. acknowledges support from the Spanish Ministry of Economy and Competitiveness (MINECO) through grant ESP2014-57495-C2-2-R. This work is based on observations made with the Italian Telescopio Nazionale Galileo (TNG), operated on the island of La Palma by the Fundación Galielo Galilei of the Istituto Nazionale di Astrofisica (INAF) at the Spanish Observatorio del Roque de los Muchachos (ORM) of the Instituto de Astrofísica de Canarias (IAC).

References

Appendix A: Additional table

All Tables

All Figures