Rotational Modulation of Spectroscopic Zeeman Signatures in Low-mass Stars

We examined the time-series Habitable-zone Planet Finder (HPF) spectra of GJ 699 and Teegarden's Star. HPF (Mahadevan et al. 2012, 2014) is a high-resolution (R ∼ 55,000), stabilized (Stefansson et al. 2016) spectrograph covering 820–1280 nm. Observing from the 10 m Hobby–Eberly Telescope (HET) at the McDonald Observatory in Texas, HPF is optimized for RV surveys of northern mid-to-late M-dwarf stars. HPF is wavelength calibrated using a custom-built laser frequency comb calibration system (Metcalf et al. 2019) and has a stable instrumental response function (Kanodia et al. 2021), enabling an RV precision of ∼1 m s⁻¹ on-sky (e.g., Lubin et al. 2021).

GJ 699 (M4) and Teegarden's Star (M7) are routinely observed standard stars for HPF, and thus have long and well-sampled baselines. Our sample of the GJ 699 spectra consists of 1095 spectra (corresponding to 126 visits) exceeding a signal-to-noise ratio (S/N) of 100 (at 10700 Å, per reduced spectral pixel) and spanning 2019 April–2021 July. Our sample of Teegarden's Star's spectra contains 155 spectra (82 visits) exceeding the same S/N threshold and spanning 2018 October to 2021 March. These spectra were produced by the standard HPF pipeline HxRGproc, which is described further in Ninan et al. (2018), Kaplan et al. (2019), and Metcalf et al. (2019). We use a synthetic telluric spectrum (Gullikson et al. 2014) to mask out regions that are affected by telluric contamination.

We evaluated the stability of individual absorption line properties (Section 2.2), as well as average line widths (Section 2.3) in these time-series spectra.

We examined the well-studied (Fuhrmeister et al. 2022) 12435.67 Å neutral potassium (K i) absorption line, shown in Figure 1. We compiled line properties from the Vienna Atomic Line Database (VALD; Ryabchikova et al. 2015) (air wavelength: 12432.27 Å; transition: 4p ²P5s ²S_1/2; effective Landé factor g_eff = 1.3; γ_radiative = 7.74836 × 10⁸ rad s⁻¹, γ_Stark = 1.65657 × 10⁸ rad s⁻¹ per 10¹² perturbers per cm³, and γ_vdW = 43.57224 × 10⁸ rad s⁻¹ per 10¹⁶ perturbers per cm³). This line is exceptionally strong and free of both significant telluric contamination and blending with strong neighboring lines. It exhibits the expected Zeeman splitting in the presence of a magnetic field, as shown in Figure 1 by comparison to solar photospheric (Livingston & Wallace 1991; Wallace et al. 1993) and sunspot (Wallace & Livingston 1992; Wallace et al. 1999) spectra and in HPF spectra of AD Leonis (Figure 1), an M dwarf with well-measured ∼kG magnetic fields (Lavail et al. 2018). This behavior is further confirmed in CARMENES spectra of multiple stars by Fuhrmeister et al. (2022), who also establish that this feature is insensitive to chromospheric variability and therefore a powerful indicator of photospheric properties, including magnetism.

Zoom In Zoom Out Reset image size

This K i feature appears in the second-reddest HPF spectral order (order index 26, spanning approximately 12367–12525 Å), which is largely clear of telluric contamination. To track the width of the K i line, we measured FWHM_K directly, using the scipy.signal.peak_widths function (Virtanen et al. 2020). To track the strength of this line, we measured the (pseudo-)equivalent width (EW_K ), defined in the usual way:

where λ₁ and λ₂ are the wavelength limits of the feature (see Figure 1), F_λ is the stellar flux at wavelength λ, and F_c is the continuum flux at λ. We account for fractional pixels at the feature boundaries using a linear interpolation across pixels, and our continuum flux level is defined by the HPF pipeline. In the spectra of GJ 699 and Teegarden's Star, the continuum level is reached at approximately ±2 Å and ±5 Å, respectively.

The K i line is notably broader in the spectra of Teegarden's Star than in the warmer GJ 699; this results from the increasing depth of line formation with decreasing temperature, which reveals the increasing effect of collisional broadening with the H₂ molecule (Fuhrmeister et al. 2022, and references therein). Our simulations indicated that the impact of Zeeman intensification would at either temperature be localized to a region near the line core, so we elect to use uniform narrow boundaries for the EW_K calculation. This helps to isolate potential Zeeman-driven variability in the K i line from other nearby spectral variations; tests with broader windows for Teegarden's star were nonetheless consistent with the results presented below.

We estimated the precision of our FWHM_K and EW_K measurements using a Monte Carlo approach and the HPF pipeline-reported variances for each one-dimensional pixel. By measuring each parameter on 100 different realizations of several different HPF spectra of GJ 699 and Teegarden's Star, we estimated the precision of our FWHM measurement at 10 mÅ (∼1%) and 2 mÅ for the equivalent width (EW) measurement (∼0.5%).

To assess the presence of Zeeman broadening, we also examined the average widths of the absorption lines in the stellar spectra using the differential line-width indicator (dLW) as defined and implemented in the SpEctrum Radial Velocity AnaLyzer (SERVAL) pipeline (Zechmeister et al. 2018). To measure the dLWs for HPF spectra, we used the SERVAL pipeline optimized for HPF spectra (Metcalf et al. 2019; Stefansson et al. 2020). Briefly, the SERVAL algorithm involves the construction of a high-S/N template and matching of individual observations against this template. The dLW results from considering the second (spectral) derivative of the template during the matching process and encodes the combined line-width changes across all lines used in the template. The dLW is analogous to the FWHM of the cross correlation function (CCF) but presents a simpler profile in the context of M-dwarf spectra which have many blended lines. In the calculation of the dLW, we mask out spectral regions impacted by tellurics and sky emission lines, as discussed in Stefansson et al. (2020). We estimated the dLW uncertainty as described in Zechmeister et al. (2018) and found a typical per-spectrum precision of 5 m² s⁻² (5%–10%).

To study the dLW behavior as a function of wavelength, we also calculated the dLW for individual spectral orders. Estimating the dLW per order enabled 10 independent measurements (corresponding to HPF's most telluric-free orders) of dLW spanning 840–1250 nm.

Teegarden's Star has also been monitored extensively with CARMENES (Zechmeister et al. 2019), which has visible (R ∼ 94,600, 5200–9600 Å) and near-IR (NIR; R ∼ 80,400, 9600–17100 Å) channels. The public data used here consist of 244 measurements in the NIR channel and 245 in the visible channel of the dLW, spanning an observing baseline from 2016 to 2019 and with typical S/N of 58 per extracted pixel around 746 nm.

To measure the periodicity of modulations present in our spectroscopic indicators, we employed the generalized Lomb–Scargle (GLS) periodogram (Zechmeister & Kürster 2009), as implemented in astropy (Astropy Collaboration et al. 2013, 2018). We searched frequencies corresponding to periods of 3–500 days, using a frequency grid sufficiently dense to provide 10 samples across a given periodogram peak. We estimated the false-alarm probability using the method described in Baluev (2008), following similar recent work (Lafarga et al. 2021). ²⁰

To assess whether the observed variations were consistent with variable Zeeman splitting, we simulated the K i line using the NICOLE (Socas-Navarro et al. 2015) program, with plane-parallel Model Atmospheres in Radiative and Convective Scheme (MARCS) atmosphere (Gustafsson et al. 2008) inputs corresponding approximately to GJ 699 and Teegarden's Star (3300 K, [Fe/H] = 0.0 and 2900 K, [Fe/H] = −0.25, respectively, and , no microturbulence) and line parameters extracted from the Vienna Atomic Line Database (VALD; given above).

For each star, we manually adjusted the K abundance to approximately match the depth of the observed line profile. For Teegarden's Star, we used [K] = 4.92 ²¹ , and for GJ 699, we used [K] = 4.62. These values are generally consistent with the old age and subsolar metallicity estimates for both stars (Muirhead et al. 2012; Zechmeister et al. 2019), though we emphasize that our goal in this work is not precise abundance determination but to look at the temporal variability of the line. We set the macroturbulent velocity to 1 km s⁻¹ for both stars, in order to approximately match the observed line widths. For simplicity, we assumed local thermodynamic equilibrium (LTE) in these simulations. We caution that departures from LTE are expected for other K i lines in cool stars (e.g., Reggiani et al. 2019) and may impact the comparison between simulated and observed spectra, although Reggiani et al. (2019) found minimal departures from LTE for the low temperatures and high log(g) values associated with M dwarfs.

We implemented a radially oriented two-component magnetic field strength distribution, described by a magnetic field strength B and a filling fraction f (Kochukhov 2021), which are uniform for all simulated lines of sight.

We disk-integrated over the visible hemisphere of the star by first simulating in NICOLE a series of 512 one-dimensional spectra along a line from disk center to disk edge. These spectra were evenly spaced in a projected radius from the disk center, and account for the effect of limb darkening as well as the changing angle between the magnetic field and the line of sight. To account for stellar rotational broadening, we applied appropriate Doppler shifts to each line of sight assuming rigid-body rotation perpendicular to the line-of-sight and rotational periods corresponding to the measured period for each star (see Section 3.2). We then formed the disk-integrated line profile by integrating over the rotating and limb-darkened lines of sight.

The time-series measurements of FWHM_K , EW_K , and dLW in GJ 699 and Teegarden's Star are shown in Figure 2 and reveal significant, similar, long-term variations over >2 yr. FWHM_K varies at the 50–100 mÅ level (∼10%), EW_K varies at the 5–10 mÅ level (0.5%–3%), while the dLWs vary at ≳100 m² s⁻². The amplitudes of these variations exceed the estimated per-spectrum uncertainties.

Zoom In Zoom Out Reset image size

Figure 3 shows the dLW measurements in detail. CARMENES data dramatically extend the dLW baseline for Teegarden's Star. CARMENES data from both the visible and NIR channels show a long-term periodicity, which is consistent with the variations detected with HPF. This is discussed in more detail in Section 4.

Zoom In Zoom Out Reset image size

We also examined the dLW measurements at different wavelengths by comparing the HPF and CARMENES dLW values (for Teegarden's Star), as well as the per-order dLW measured in the HPF spectra (Figure 3). A trend of increasing variability amplitude with increasing wavelength is apparent for both stars in the per-order dLW values and between instruments for Teegarden's Star. As discussed below, this increasing amplitude with wavelength is expected of Zeeman broadening (Reiners 2012).

To gain insight into the behavior of the K i line, we also examined the line-shape change between epochs of low and high EW_K , as shown in Figure 4. This figure confirms that higher EW_K correlates with a broader line profile (as suggested in Figure 2) and also with a shallower line core.

Zoom In Zoom Out Reset image size

We also explored whether the variation observed in EW_K and FWHM_K is consistent with what might be produced by rotational modulation of magnetic features on the surfaces of GJ 699 and Teegarden's Star, using the simulations described in Section 2.6. Figure 4 shows the resulting variation of the line-shape parameters when the magnetic field strength and filling factor of our two-component model are varied. We note that the EW_K and FWHM_K are dependent on the precise choice of continuum level, which cannot be defined equivalently in the observed and simulated spectra, because the latter has only a single source of opacity. We therefore include an offset for each of these quantities to align the simulated and measured values, and we emphasize that the relevant comparison here is the range of variation.

For both stars, we observe that the range of variation and the measured line shapes themselves can be approximately reproduced by magnetic field strength and filling factor variations in the simulated spectra (Figure 4). We note that our two-component model is likely a significant oversimplification of the true magnetic field configuration, and that there are numerous degeneracies between the line profile measurements and the line parameters adopted for these simulations which complicate the reliable determination of a single magnetic field and filling factor for a given spectrum. ²² Nonetheless, we conclude that the range of variations observed in our line profile measurements can at least be plausibly explained by variations in our magnetic field model.

The GLS periodogram results for our measurements, including the periods of the most significant peaks, are shown in Figure 5. For GJ 699, the most recent observing season appears qualitatively different from the previous two; this change may be related to a recently uncovered detector bias voltage drift. For completeness, we retain these data points in our data set, but for our period analysis we conservatively use only the data through 2020.

Zoom In Zoom Out Reset image size

For Teegarden's Star, the spectroscopic indicators are unanimous in the detection of a strong periodicity. We estimate the period by taking the mean and standard deviation of the periodogram peaks from the five time series (HPF dLW, Carmenes-NIR and VIS dLW, HPF FWHM_K , and HPF EW_K ), finding 99.6 ± 1.4 days, which we identify with the rotation period.

For GJ 699, we found strong peaks in the periodograms near the suspected 145 ± 15 days rotation period of this star (Toledo-Padrón et al. 2019). Among the HPF time series, the mean and standard deviation of these peak values are 155 days and 2 days, respectively. We note that the time baseline of our data set is smaller than that of Toledo-Padrón et al. (2019), who find a forest of periodogram peaks for GJ 699 at 145 ± 15 days, which may be related to differential rotation. The periodicity detected in our data set is more tightly localized around 155 days; this may be an artifact of the shorter time baseline. Another significant peak in our periodograms around 300 days may be a harmonic of the rotation period or due to a longer-term activity cycle.

We argue that the observed variations are consistent with a changing magnetic field imprinting on the spectra via the Zeeman effect. The indicators (FWHM_K , EW_K , and dLW) change in phase with each other (Figure 2), and at levels consistent with our simple model of a changing surface magnetic field. Further, the line-width changes are stronger at longer wavelengths (Figure 3), as is expected for Zeeman broadening, in which the wavelength separation of the split line components (Δλ) from the unsplit line center (λ₀) behaves as

where g_eff is the effective Landé factor and B is the magnetic field (Kochukhov 2021). Finally, the line-shape change of the K i line is consistent with simple simulations of the impact of a changing magnetic field (Figure 4).

We note, however, the potential contribution of other effects on the stellar surface. For example, stellar rotation that brings regions of varying magnetic field strength into and out of view may also modulate the appearance of associated bright/hot regions (faculae) or dark/cool regions (spots). Spot and faculae sizes and filling fractions are not well constrained for these M dwarfs (e.g., Rackham et al. 2018), and we cannot rule out that the combined effects of these regions may be the origin of some of the rotational modulation we observe. However, we do note that in simple (single-temperature) simulations, the K i line becomes both deeper and broader with decreasing temperature, rather than the stronger/shallower behavior shown in Figure 4. This favors the Zeeman explanation over a simple temperature contrast.

A fully self-consistent description of the effect of magnetic fields on the stellar surface would account for the impact of the magnetic fields on the gas dynamics via magnetohydrodynamic simulation. Indeed, the resulting changes of line positions and profile asymmetries are expected to be an important pathway for stellar activity to print through to RV measurements (Meunier et al. 2017), although the dynamics may be somewhat different for M dwarfs than for FGK stars (Liebing et al. 2021). With only the K i line examined in detail, we cannot measure the low-level line asymmetries which accompany changes in the visible convective pattern and so cannot constrain if or how this may impact our measurements. However, we do note that Toledo-Padrón et al. (2019) examined the bisector behavior in GJ 699 and found no detectable signal of changing line asymmetry at their level of measurement precision.

Finally, it is possible that the measured line strengths are also impacted by variable chromospheric emission, which is thought to drive some of the variability in line-by-line studies of correlation with established activity indicators (e.g., Wise et al. 2018). However, Fuhrmeister et al. (2022) considered correlations between the established chromospheric Hα activity indicator and the K i line for a range of M dwarfs and found that it exhibits no measurable chromospheric component (in contrast to K I lines at shorter wavelengths). Further, pure chromospheric "filling-in" of the line core would not be expected to exhibit the tight correlation observed between EW_K and FWHM_K , nor would we expect it to be so widespread across the spectrum as indicated by the dLW measurements.

The strength of the rotational modulations stands in contrast to the low photometric variability in both GJ 699 and Teegarden's Star. Careful analysis revealed no rotation period signal in the photometry of Teegarden's Star by Zechmeister et al. (2019), and the rotation period was only marginally detected in some subsets of photometric data used by Toledo-Padrón et al. (2019). The clear spectroscopic detection of rotation in both stars echos the findings of Lafarga et al. (2021), who showed that the rotation period is frequently recoverable in the cross-correlation function parameters, chromospheric indicators, and dLWs of M dwarfs. These findings indicate that spectroscopic monitoring may often provide the best handle on the rotation period for these types of stars.

The need to disentangle stellar activity from exoplanetary signals in RV data is well recognized and has motivated careful statistical (e.g., Davis et al. 2017) and line-by-line (e.g., Dumusque 2018; Wise et al. 2018) analyses to tease these effects apart. To the extent that our observed modulations are driven by the Zeeman effect, this provides a compelling justification for spectroscopic monitoring of RV targets (especially fully convective M dwarfs) at infrared wavelengths, where the Zeeman impact is largest. Together with the results of Lafarga et al. (2021), the K i line studied here appears to provide a powerful indicator of the changing photospheric magnetic field on these stars, which may be expected to perform well as a tracer of rotational modulation and as an RV activity indicator.

Recently, Klein et al. (2021) used Zeeman Doppler imaging to recover the magnetic field geometry of Proxima Centauri, an old, slowly rotating, fully convective star similar to GJ 699 and Teegarden's Star (albeit more active). They additionally found a significant rotational modulation in the Stokes I spectra of Proxima Cen and discuss the value of contrasting polarized and unpolarized indicators for diagnosing the properties of the dynamos of these fully convective stars. Based on the significance of the unpolarized Zeeman signatures presented here, we suggest that long-term spectropolarimetric and spectroscopic monitoring of GJ 699 and Teegarden's Star may prove similarly revealing about the (differential) rotation, dynamo mechanisms, and activity cycles of fully convective stars.