Seismic Signatures of Stellar Magnetic Activity -- What Can We Expect from TESS?

TL;DR

This paper explores how TESS data can be used to detect stellar magnetic activity through seismic signatures, deriving a new scaling relation for activity-related frequency shifts and estimating the number of stars where these shifts can be observed.

Contribution

It introduces a new scaling relation for acoustic frequency shifts due to stellar activity and predicts the potential yield of detectable signals with TESS.

Findings

Significant p-mode frequency shifts are expected in hundreds of main-sequence and early subgiant stars.

TESS can potentially detect activity-related seismic signatures in thousands of late subgiant and red giant stars.

A new scaling relation links stellar activity amplitude to seismic frequency shifts.

Abstract

Asteroseismic methods offer a means to investigate stellar activity and activity cycles as well as to identify those properties of stars which are crucial for the operation of stellar dynamos. With data from CoRoT and \textit{Kepler}, signatures of magnetic activity have been found in the seismic properties of a few dozen stars. Now, NASA's Transiting Exoplanet Survey Satellite (TESS) mission offers the possibility to expand this, so far, rather exclusive group of stars. This promises to deliver new insight into the parameters that govern stellar magnetic activity as a function of stellar mass, age, and rotation rate. We derive a new scaling relation for the amplitude of the activity-related acoustic (p-mode) frequency shifts that can be expected over a full stellar cycle. Building on a catalogue of synthetic TESS time series, we use the shifts obtained from this relation and simulate the yield of detectable frequency shifts in an extended TESS mission. We find that, according to our scaling relation, we can expect to find significant p-mode frequency shifts for a couple hundred main-sequence and early subgiant stars and for a few thousand late subgiant and low-luminosity red giant stars.

Tables

Table 1. Table 1 : Predicted number of stars with significant frequency shifts for different spectral types and evolutionary stages. A four year TESS mission extension is assumed.

Spectral type	Region	Number	CC Detections	CC Fraction [%]	PB Detections	PB Fraction [%]
F	I	6676	141	2.1	673	10.1
F	II	1642	24	1.5	196	11.9
G	I	54	30	55.6	44	81.5
G	II	1628	765	47.0	1092	67.1
K	II	2731	1486	54.4	1799	65.9
FG	I	6730	171	2.5	717	10.7
FGK	II	6001	2275	37.9	3087	51.4
total		12731	2446	19.2	3804	29.9

Table 2. Table 2 : Number of stars with N 𝑁 N consecutive sectors of data, number of significant shift detections, and fraction of detections of the total number of stars, for all stars in the sample, the main-sequence stars (region I), and post-main-sequence stars (region II). Shifts measured with the CC and PB methods for a two year TESS mission extension.

Sectors $N$	1	2	3	4	5	6	7	8	9	10	11	12	13
Total number	7930	2601	781	306	9	0	0	40	144	105	114	86	615
CC Detections	1101	595	232	98	3	0	0	12	53	38	44	32	141
CC Fraction [%]	13.9	22.9	29.7	32.0	33.3	0.0	0.0	30.0	36.8	36.2	38.6	37.2	22.9
PB Detections	1749	782	284	114	4	0	0	14	56	45	48	37	172
PB Fraction [%]	22.1	30.1	36.4	37.3	44.4	0.0	0.0	35.0	38.9	42.9	42.1	43.0	28.0
Region I	4251	1391	393	155	3	0	0	20	69	53	50	42	303
CC Detections region I	57	17	8	2	0	0	0	0	2	3	0	1	4
CC Fraction region I [%]	1.3	1.2	2.0	1.3	0.0	0.0	0.0	0.0	2.9	5.7	0.0	2.4	1.3
PB Detections region I	287	104	30	8	0	0	0	1	2	7	2	3	14
PB Fraction region I [%]	6.8	7.5	7.6	5.2	0.0	0.0	0.0	5.0	2.9	13.2	4.0	7.1	4.6
Region II	3679	1210	388	151	6	0	0	20	75	52	64	44	312
CC Detections region II	1044	578	224	96	3	0	0	12	51	35	44	31	137
CC Fraction region II [%]	28.4	47.8	57.7	63.6	50.0	0.0	0.0	60.0	68.0	67.3	68.8	70.5	43.9
PB Detections region II	1462	678	254	106	4	0	0	13	54	38	46	34	158
PB Fraction region II [%]	39.7	56.0	65.5	70.2	66.7	0.0	0.0	65.0	72.0	73.1	71.9	77.3	50.6

Table 3. Table 3 : Number of stars with N 𝑁 N consecutive sectors of data, number of significant shift detections, and fraction of detections of the total number of stars, for all stars in the sample, the main-sequence stars (region I), and post-main-sequence stars (region II). Shifts measured with the CC and PB methods for a four year TESS mission extension.

Sectors $N$	1	2	3	4	5	6	7	8	9	10	11	12	13
Total number	7930	2601	781	306	9	0	0	40	144	105	114	86	615
CC Detections	1139	628	248	102	3	0	0	12	55	39	39	33	148
CC Fraction [%]	14.4	24.1	31.8	33.3	33.3	0.0	0.0	30.0	38.2	37.1	34.2	38.4	24.1
PB Detections	2094	874	311	126	4	0	0	14	65	45	47	40	184
PB Fraction [%]	26.4	33.6	39.8	41.2	44.4	0.0	0.0	35.0	45.1	42.9	41.2	46.5	29.9
Region I	4251	1391	393	155	3	0	0	20	69	53	50	42	303
CC Detections region I	108	30	14	6	0	0	0	0	0	3	1	3	6
CC Fraction region I [%]	2.5	2.2	3.6	3.9	0.0	0.0	0.0	0.0	0.0	5.7	2.0	7.1	2.0
PB Detections region I	473	140	42	20	0	0	0	1	5	7	2	4	23
PB Fraction region I [%]	11.1	10.1	10.7	12.9	0.0	0.0	0.0	5.0	7.2	13.2	4.0	9.5	7.6
Region II	3679	1210	388	151	6	0	0	20	75	52	64	44	312
CC Detections region II	1031	598	234	96	3	0	0	12	55	36	38	30	142
CC Fraction region II [%]	28.0	49.4	60.3	63.6	50.0	0.0	0.0	60.0	73.3	69.2	59.4	68.2	45.5
PB Detections region II	1621	734	269	106	4	0	0	13	60	38	45	36	161
PB Fraction region II [%]	44.1	60.7	69.3	70.2	66.7	0.0	0.0	65.0	80.0	73.1	70.3	81.8	51.6

Table 4. Table 4 : Detectability of frequency shifts for different samples of stars. A two year TESS mission extension is assumed.

Spectral type	Region	Number	CC Detections	CC Fraction [%]	PB Detections	PB Fraction [%]
F	I	6676	62	0.9	416	6.2
F	II	1642	6	0.4	49	3.0
G	I	54	32	59.3	42	77.8
G	II	1628	775	47.6	1120	68.8
K	II	2731	1474	54.0	1678	61.4
FG	I	6730	94	1.4	458	6.8
FGK	II	6001	2255	37.6	2847	47.4
total		12731	2349	18.5	3305	26.0

Equations

$$\Delta R_{\text{HK}}^{\prime }\propto {R_{\text{HK}}^{\prime }}^{0.77}.$$

$$\delta \nu \propto \frac{D}{I}\Delta R_{\text{HK}}^{\prime },$$

$$\delta \nu \propto \frac{R^{2.5}{L}^{0.25}}{M^2}\Delta R_{\text{HK}}^{\prime }$$

$$\delta \nu =\frac{\int dVKS}{2I\nu },$$

$$K$$

$$S$$

$$\vec{\xi }={\xi }_r(r){\vec{e}}_r,$$

$$K={\left|\frac{2}{r}{\xi }_r(r)+\frac{\partial {\xi }_r(r)}{\partial r}\right|}^2.$$

$$K\approx \frac{9|{\xi }_r(R){|}^2}{R^2}.$$

$$I=\int dV\rho {\left|\vec{\xi }\right|}^2\approx R^3\overline{\rho }|{\xi }_r(R){|}^2,$$

$$\delta \nu$$

$$=\text{const}\cdot \frac{R}{M\nu }.$$

$$\log \left({\tau }_c/P_{\text{rot}}\right)=-\left(0.324-0.400y+0.283y^2-1.325y^3\right),$$

$${\tau }_c=\frac{\alpha H_p}{v},$$

$$R_{\text{HK, mod}}^{\prime }=\frac{R_{\text{HK}}^{\prime }}{R^2\max \left(\sqrt{\text{Ro}-1},1\right)t_{\text{Age}}^{0.5}}.$$

$$\delta \nu$$

$$\delta {\nu }_{\text{PB}}={\nu }_{\ell ,n}(t)-\overline{{\nu }_{\ell ,n}},$$

$$P_{\text{cyc}}\left[\text{years}\right]=0.5P_{\text{rot}}\left[\text{days}\right]\cdot \max \left(\sqrt{\text{Ro}-1},1\right).$$

$$\delta {\nu }_{\text{input shift}}=\delta \nu \cdot \min \left(\frac{2n}{P_{\text{cyc}}},1\right),$$

$$M=M_{\odot }{\left(\frac{{\nu }_{\text{max}}}{{\nu }_{\text{max, }\odot }}\right)}^3{\left(\frac{\Delta {\nu }_{\odot }}{\Delta \nu }\right)}^4{\left(\frac{T_{\text{eff}}}{T_{\text{eff,}\odot }}\right)}^{\frac{3}{2}},$$

Full text

\helveticabold

1 Keywords:

stellar activity, asteroseismology, stellar activity cycles, TESS, p-mode frequency shifts

2 Introduction

The primary seismic signature of stellar activity and activity cycles is a systematic change in mode frequencies over timescales associated with activity cycles, which are typically of the order of months to decades [1, 2, 3]. Measuring this phenomenon for stars is primarily limited by the baseline length of the observations. For the Sun, however, we have Sun-as-a-star helioseismic observations lasting for decades, covering several eleven-year solar activity cycles. From these data it is found that the frequency of low-degree modes close to the frequency of maximum oscillation power $\nu_{\textrm{max}}$ change by approximately $0.4\,\mathrm{\mu Hz}$ between activity minimum and maximum (e.g., [4, 5]). Furthermore, this shift is strongly correlated with the level of activity along the cycle [4, 6].

The first asteroseismic evidence of variations in stellar activity were found in an F5 V star observed by CoRoT [7, 8], HD 49933. For this star García et al. [9] found that not only did the frequencies vary systematically with time but, as is also observed for the Sun [10, 11, 6], the amplitudes of the modes varied in anti-phase with the frequencies. In a more detailed analysis, Salabert et al. [12] found that the observed frequency shifts are larger than seen in the Sun but the same dependence on mode frequency is observed, i.e., the magnitude of the frequency shift increases with mode frequency (e.g., [13, 14]).

Since then examples of systematic variations in the frequencies of stars observed by Kepler [15, 16] have been found (Salabert et al. [17], Régulo et al. [18], Kiefer et al. [19], Salabert et al. [20], Santos et al. [21]). Indeed, Kiefer et al. [19] found that 23 out of the 24 stars in their sample showed significant and systematic frequency shifts in time, while Santos et al. [21] found quasi-periodic variations on 60% of their sample of 87 solar-type stars. The sample of potential candidates for observing asteroseismic activity cycles is limited by the fact that activity suppresses the amplitude of acoustic oscillations [22]. Since activity also increases with decreasing rotation period [e.g., 23, 24, 25, 26] these stars tend to be more active and consequently are less likely to have high signal-to-noise oscillations. The dependence of the length of stellar activity cycles on surface rotation rate or on Rossby number is complex and is actively researched [e.g., 27, 28, 29]. While it has proven to be difficult to observe complete activity cycles with Kepler, it is still possible to search for systematic variations with time in the properties of p modes. A concise review on the capabilities and prospects of asteroseismology regarding the inference on stellar cycles and activity is also given by Chaplin and Basu [30].

Stellar magnetic activity also manifests in variability of stellar luminosity [see, e.g., 31, and references therein]. Over the solar cycle the total irradiance is found to change by about 0.1% in temporal correlation with solar activity [32]. This slight increase in luminosity is associated with faculae [32]. In contrast to faculae, spots decrease solar and stellar luminosities [32, 33]. For the Sun, the increase in luminosity caused by faculae dominates over the decrease in luminosity brought about by sunspots, hence the overall increase in total irradiance between cycle minimum and maximum. However, for other stars it is the decrease in luminosity caused by the presence of spots that is the dominant source of variability [34]. Magnetic features have typical life times which are of the order of the rotation period of a star. This has to be taken into account when a meaningful proxy of stellar magnetic activity is to be constructed from the variability of the stellar luminosity. Indeed, Mathur et al. [35] defined such a photometric proxy for stellar activity with the $S_{\text{ph}}$ index. It is given by taking the mean of the standard deviations of segments of the time series each having a length 5 times the rotation period of the star (see also [36] and [37]). As long cadence photometric time series (CoRoT: 512 s, Kepler: 29.4 min, TESS: 30 min) are sufficient for such studies, the number of available light curves for the investigation of stellar magnetic activity through photometry is rather large (CoRoT: $\sim$163 000 targets, Kepler: $\sim$197 000 targets, TESS: $>10^{7}$ targets at end of nominal mission). Studies in which the variability of light curves is explored regarding signatures of activity or activity cycles include, but are not limited to [3, 38, 36, 39, 40]. In this article, we do not consider the variability of light curves, as activity-related changes are not included in the synthetic time series we use (see [41] and Section 3).

With the launch of NASA’s Transiting Exoplanet Survey Satellite [TESS, 42], the opportunity of studying seismic signatures of many more stars is at hand. So far, p-mode frequency shifts have only been detected for a few dozen Kepler targets and one CoRoT target, and as of yet, the dependence of the amplitude of the activity-related frequency shifts (which corresponds to the amplitude of stellar activity cycles) on fundamental stellar parameters is not known. With the large number of stars that will be observed by TESS, asteroseismology can become an important additional tool to learn about stellar cycles. As we will discuss in Section 7, a mission extension of TESS to four or six years is however necessary for this. This paper is aimed at a prediction of how many detections of activity-related frequency shifts we can expect in an extended TESS mission.

3 Synthetic TESS sample

To simulate the yield of seismic signatures of activity that can be expected from TESS, we use the synthetic catalogue of Ball et al. [41], which consists of light curves that realistically mimic the properties of the TESS targets in short-cadence data (120 s). The sample consists of 12 731 stars, covering main-sequence, subgiant, and low-luminosity red giant stars, which were selected from a synthetic Milky Way population simulated with TRILEGAL [43]. The selection criteria are the same as for the real Asteroseismic Target List (ATL; see Schofield et al. [44]) of the TESS Asteroseismic Science Consortium (TASC), which are chosen such as to maximise the yield of detected oscillations and optimise the coverage of different stages of stellar evolution. Ball et al. [41] give a detailed description of their methods to obtain oscillation mode parameters (frequencies, mode damping widths, amplitudes), noise levels, granulation background, and other aspects of the simulated time series. We note that the time series of Ball et al. [41] are provided without noise, which has to be added as $\mathcal{N}\left(0,\sigma_{\text{noise}}\right)\cdot\sqrt{30}$ to the time series. Here, $\sigma_{\text{noise}}$ is the noise amplitude as given in the catalogue of Ball et al. [41] and the factor $\sqrt{30}$ is accounting for the 2-minute sampling of the time series.

Figure 1 shows a Kiel diagram of the synthetic sample with stellar ages given by the colour of the dots. We separated the sample into main-sequence stars (in the following called region I) and post-main-sequence stars (region II) along the Terminal Age Main-Sequence (TAMS), which we defined at a core hydrogen abundance content of $10^{-5}$. Figure 7 shows the Kiel diagram as in Figure 1 but colour coded for the two regions. As can be seen in these figures, the number of stars on the lower main-sequence is rather small. This is due to the low intrinsic mode amplitudes of solar-like oscillations in cool dwarf stars (see, e.g., Kjeldsen and Bedding [45]). Thus, these stars have a low detection probability for their oscillation modes with TESS data [46]. A large number of targets are young, early subgiants and stars towards the end of their main-sequence life with $6000\,\mathrm{K}\lesssim T_{\text{eff}}\lesssim 6900\,\mathrm{K}$. Subgiants ascending to the red giant branch and low-luminosity red giants have large mode amplitudes and high detection probabilities, hence, many of the best asteroseismic targets are in this region of the Kiel diagram.

For a study of the temporal evolution of seismic parameters, it is necessary that the targets are observed for either a long, uninterrupted period of time or that they are observed again after some time. During its nominal two year mission, TESS will observe 13 sectors per hemisphere with an area of four times $24^{\circ}\times 24^{\circ}$ for a duration of about 27.4 days per sector. The sectors partly overlap, chiefly in the region of the ecliptic poles. Thus, some stars will have up to 13 sectors’ worth of data. A detailed description of the TESS observing strategy can be found in Ricker et al. [42]. The time series of Ball et al. [41] reproduce this observing strategy, i.e., their length depends on the position of the respective star in the sky. Here, we shall assume that the TESS mission is extended after its two year nominal duration. We further assume that the observation strategy of the nominal mission will be repeated in possible two year and four year mission extensions, and thus that the same targets are observed again for the same number of sectors after two and four years. However, the scaling relations for p-mode frequency shifts described below are not TESS-specific.

4 Reconsidering the scaling relation for p-mode frequency shifts

There are two models for activity related frequency shifts in the literature: Chaplin et al. [47] simply use the amplitude of the activity cycle $\Delta R^{\prime}_{\text{HK}}$ (see Equation 1) as a proxy for the frequency shift and then scale the values such that the solar frequency shift equals the observed $0.4\,\mathrm{\mu Hz}$. Saar and Brandenburg [48] found that $\Delta R^{\prime}_{\text{HK}}$ scales with the average logarithmic fraction of the stellar luminosity that is emitted in the Calcium II H and K spectral lines $R^{\prime}_{\text{HK}}$:

[TABLE]

For the second model, Metcalfe et al. [49] find that the frequency shifts $\delta\nu$ scale as

[TABLE]

where $I$ is the inertia of the modes and $D$ is the depth of the perturbation that causes the mode frequencies to shift. They further assume that the depth $D$ scales as the pressure scale height $H_{p}$, which is proportional to $L^{0.25}R^{1.5}M^{-1}$ at the stellar photosphere, and that mode inertia scales as $MR^{-1}$, where $R$, $M$, and $L$ are stellar radius, mass, and luminosity respectively. Thus, Metcalfe et al. [49] find

[TABLE]

for the frequency shift amplitude over a full activity cycle.

Here we take a new look at the sensitivity of acoustic mode frequencies to activity. As a starting point we use Equation (1) of Metcalfe et al. [49] to estimate the frequency shift $\delta\nu$:

[TABLE]

where $K$ is a kernel for the sensitivity of p modes to perturbations, $S$ is a source function, and $\nu$ is the frequency of the mode. Following Metcalfe et al. [49], we use

[TABLE]

where $\vec{\xi}$ is the eigenfunction of the considered mode, $A$ is the amplitude of the source function, and the $\delta$ function determines the location of the source $D_{c}$. For now, we set $A=1$ and reintroduce it later via the $R^{\prime}_{\text{HK}}$ index.

In the following, we concentrate on radial modes, which simplifies the calculations considerably. The eigenfunctions for radial modes can be described by

[TABLE]

where $\xi_{r}(r)$ is the radial displacement eigenfunction and $\vec{e}_{r}$ is the unit vector in the radial direction. We do not introduce the spherical harmonic $Y_{0}^{0}$ into the eigenfunction here, as it would only appear as a constant in the following equations.

With this, the kernel function can be written as

[TABLE]

Approximating the derivative $\frac{\partial\xi_{r}(r)}{\partial r}\approx\frac{\xi_{r}(r)}{r}$ and assuming that the radial position $r$ which we are considering is close enough to the surface so that the amplitude of the eigenfunction is similar to that at the photospheric stellar radius $R$ we get

[TABLE]

Thus we have fixed the source to be at $D_{c}=R$. The mode inertia $I$ can be approximated as follows:

[TABLE]

where we used again the eigenfunction of a radial mode, replaced the radial integral $\int dr$ by simply multiplying with the photospheric radius $R$, and replaced the density $\rho$ with the mean density $\overline{\rho}$.

With Equations (8) and (10) we can rewrite Equation (4):

[TABLE]

Throughout this, we assumed that the source is close enough to the photospheric radius so that $D_{c}=R$. All constant factors are absorbed into const. Equation (11) estimates the sensitivity of radial modes to a near-surface perturbation.

The strength of the perturbation and hence the magnitude of the frequency shift is imparted by the strength of the magnetic activity on the star for which a proxy is used. From the synthetic TESS sample, the $R^{\prime}_{\text{HK}}$ index can be estimated with a relation found by Noyes et al. [50] for main-sequence stars:

[TABLE]

where $y=\log\left(R^{\prime}_{\text{HK}}\times 10^{5}\right)$, $\tau_{c}$ is the convective turnover time, and $P_{\text{rot}}$ is the star’s rotation period. Both these values are taken from the models of Ball et al. [41]. The (local) convective turnover time $\tau_{c}$ is calculated half a mixing length $\alpha H_{p}/2$ above the base of the outer convection zone of each star following Gilliland [51], Landin et al. [52]:

[TABLE]

where the solar-calibrated mixing length parameter $\alpha=1.957$ is used, $H_{p}$ is the pressure scale height, and $v$ is the convective velocity.

To take into account that a large fraction of the stars in the synthetic TESS sample are not on the main-sequence (in contrast to the stars considered in Saar and Brandenburg [48] and Noyes et al. [50]), we modify the values of $R^{\prime}_{\text{HK}}$ with three factors: First, we include a factor $R^{-2}$ to account for the increasing surface area over which magnetic activity can be distributed, thus affecting the oscillations less. Second, starting at a Rossby number of $\text{Ro}=\frac{P_{\text{rot}}}{\tau_{c}}=2$, we divide by $\sqrt{\text{Ro}-1}$. This is to accommodate the results of Metcalfe and van Saders [26] who found that, as stars reach a critical Rossby number of $\approx 2$, a transition in the magnetic activity appears to happen: Magnetic braking weakens, rotation periods increase more slowly, and cycle periods grow longer and eventually disappear. Finally, we explicitly include the weakening of magnetic activity with the inverse square root of stellar age $t_{\text{Age}}$ $\left[\text{Gyrs}\right]$ found by Skumanich [53]:

[TABLE]

With this and Equation (11), the expected frequency shift over a complete stellar activity cycle at the frequency of maximum power $\nu_{\text{max}}$ scales as:

[TABLE]

Figure 2 shows the expected p-mode frequency shifts for the synthetic TESS sample between the activity minima and maxima of a full stellar activity cycle. The values are scaled to a solar frequency shift of 0.4 $\mu$Hz at $\nu_{\text{max}}$. Stellar age is given by the colour of the dots. For the Sun, we assume $R^{\prime}_{\text{HK}}=-4.901$ (value from [2]) and do not estimate it with Equation (12). The spike in shifts at 6500 K is caused by a transition of how rotation rates are modelled by [41], cf. their Section 2.5 and their Figure 2.

5 Measuring $\delta\nu$

Two main methods have been developed to determine the change in p-mode frequency as a function of time. The first technique uses the cross-correlation (CC) between power spectra of different epochs of the time series (see, e.g., [54, 55, 18, 19]). The regions around the frequency of maximum oscillation amplitude $\nu_{\text{max}}$ of the power spectra are retained and their cross-correlation function is computed. Typically, only the region in which p-mode peaks can be identified is retained. As the expected frequency shifts between two power spectra of any epoch is only of the order of a few $\mathrm{\mu Hz}$, the complete cross-correlation function does not need to be calculated: restricting the cross-correlation function to a certain lag regime considerably speeds up the computation. The central peak of the cross-correlation function is then fitted to obtain the frequency shift $\delta\nu_{\text{CC}}$ between the two power spectra. To estimate the uncertainty $\sigma$ of the measured frequency shift, we employ the resampling approach described in [18].

The second technique uses the individual mode frequencies (peak-bagging: PB) obtained at different epochs. Here, the frequency shift $\delta\nu_{\text{PB}}$, of a particular mode with harmonic degree $\ell$ and radial order $n$ is given by

[TABLE]

where $\nu_{\ell,n}(t)$ is the mode frequency obtained from a time series observed between $t$ and $t+\delta t$, and $\overline{\nu_{\ell,n}}$ is a reference frequency for that mode. Both methods (CC and PB) require the data to be split into time series of length $\delta t$, where $\delta t$ is some fraction of the full baseline of observations available. Ideally, $\delta t$ will be sufficiently short to resolve changes with time but long enough that the power spectrum has sufficient resolution and signal-to-noise for the modes to be observed clearly. See [21] for a detailed description of the peak-bagging procedure in determining frequency shifts.

Most of the targets in the synthetic TESS sample have only few sectors’ worth of data (one sector: 7930 stars, two sectors: 2601 stars, three sectors: 781 stars, four sectors or more: 1419; see also Tables 2 and 3). We thus calculate the power spectra for the complete synthetic time series of each star without separating them into shorter segments: here, $\delta t$ is equal to the length of the time series.

In general, the two methods CC and PB are found to produce results in good agreement with each other but the uncertainties associated with the cross-correlation technique tend to be systematically larger than the frequency shifts determined directly from the mode frequencies (see [55, 21] for comparisons of the frequency shift uncertainty obtained with the two methods). However, the peak-bagging method requires good resolution and signal-to-noise for the mode frequencies to be determined accurately and precisely and so there are instances where the cross-correlation method is preferable. The PB technique is also much more reliant on input parameters like initial guesses for mode frequencies and widths than the CC method. We do not apply the PB method here, but assume that the frequency shifts measured with both methods are in agreement and that the uncertainties of the PB method are smaller by a factor of two than those obtained with the CC method.

We implemented the CC technique to measure the frequency shifts in the following way: the power spectrum is calculated from the synthetic time series with an oversampling factor of 8. This power spectrum is then shifted in frequency by the amount that can be expected over the fraction of the stellar cycle that is covered in the two (four) years that pass between the TESS nominal mission and the second observation of any given star (see Section 3). We estimate the length of the activity cycles by

[TABLE]

We chose this rather crude estimation as a compromise between the short and long activity cycle branches that can be found in a $P_{\text{cyc}}$-$P_{\text{rot}}$ diagram like that of Böhm-Vitense [23] or Metcalfe and van Saders [26]. The lengthening of activity cycles above $\text{Ro}\approx 2$ found by Metcalfe and van Saders [26] is again modelled into this with the factor $\max\left(\sqrt{\text{Ro}-1},1\right)$. Thus, the input frequency shift is given by

[TABLE]

where the factor 2 accounts for the fact that the full cycle frequency shift is observed over a half-cycle, between activity minimum and maximum, $n=2$ or 4 is the number of years between observations, depending on the length of the TESS mission extension, and $\delta\nu$ is calculated with Equation (15). After shifting the power spectrum, we return to only the natural frequencies without oversampling. The oversampling was only introduced to be able to shift by frequencies smaller than the natural frequency resolution, which is $\approx 0.422\,\mathrm{\mu Hz}$ for stars with only one sector of data.

The regions around the frequency of maximum oscillation amplitude $\nu_{\text{max}}$ of the unshifted and shifted power spectra are then cross-correlated. We retain a width of $\pm 5\Delta\nu$ around $\nu_{\text{max}}$, where $\Delta\nu$ is the large frequency separation. Both $\Delta\nu$ and $\nu_{\text{max}}$ are given in the catalogue of Ball et al. [41]. To obtain the frequency shift between the two power spectra, we fit the cross-correlation function with a Lorentzian profile.

6 Influence of assumptions on the predicted shift amplitudes

There are several assumptions made and approximations applied in the derivation of the scaling relation given by Equation (15). The first assumption is that the magnitude of activity-related frequency shifts should be proportional to the strength of stellar magnetic activity. This is certainly justified by the fact that for the Sun there is a clear positive correlation between the p-mode frequency shifts and proxies for magnetic activity as shown by, e.g., [6]. Also, in a follow-up study to [21], in which they investigate the dependence of the observed frequency shifts on stellar fundamental parameters, Santos et al. (submit.) find that for a sample of 30 solar-like stars, the amplitude of the frequency shifts decreases with decreasing $R^{\prime}_{\text{HK}}$. Hence, the $R^{\prime}_{\text{HK}}$ activity index, which is available for a relatively large number of stars (see, e.g., [1] and the database of [56]) and is employed in the scaling relation of [55], is a good starting point for stars similar to the Sun. In principle, spectroscopic activity indices other than $R^{\prime}_{\text{HK}}$ can be used in any of the frequency shift scaling relations we discuss here; [49], for example, use the Mg II emission index $i_{\text{Mg II}}$. However, the scaling of the amplitude of the index over a complete activity cycle, as is done for $R^{\prime}_{\text{HK}}$ in Equation (1), has to be adjusted accordingly.

To account not only for the strength of the activity, as is done in the Chaplin et al. scaling, but also for the physics of the oscillation modes, [49] derived the sensitivity of mode frequencies to near surface perturbations. We note that the sensitivity of modes we find with Equation (11) is in fact very similar to that found by [49] and quoted by [57] in their Equation (14). Inserting the scaling relation of [45] $\nu_{\text{max}}\propto MR^{-1}T_{\text{eff}}^{-0.5}$ into Equation (11) and using the Stefan-Boltzmann law $\sqrt{T_{\text{eff}}}\propto\left(R^{-2}L\right)^{1/8}$, gives $\delta\nu\approx\frac{R^{2.75}L^{0.125}}{M^{2}}$. All quantities are assumed to be normalised to solar values in this discussion. The difference of our Equation (11) to the mode sensitivity factor of [49] is thus only a factor $R^{-0.25}L^{0.125}\propto\sqrt{T_{\text{eff}}}$, which is in the range 0.91–1.14 for the sample of stars investigated here, i.e., is of order unity. Thus, the approximations we made in the derivation of Equation (11) — approximation of the radial integral by simply multiplying with the stellar radius $R$, approximation of the derivatives, localising the perturbation to the surface — yield essentially the same result as the calculations of [49].

The three scaling relations are most easily distinguished for stars reaching the lower red giant branch, which is where the disagreement is largest. In Figures 8 and 9 we show the full-cycle frequency shifts of the synthetic TESS sample for the Chaplin et al. and Metcalfe et al. scalings, respectively. The Chaplin et al. [47] relation, simply using the $R^{\prime}_{\text{HK}}$ estimation of [50] for all stars from main-sequence to the lower red giant branch, gives maximal full-cycle shift values of $\approx 3\,\mathrm{\mu Hz}$. With estimated cycle lengths of at least 10 years, no shifts larger than $\approx 1.2\,\mathrm{\mu Hz}$ should be found in TESS data with a two year mission extension. This value increases to $\approx 5\,\mathrm{\mu Hz}$ for the new scaling relation presented in this article and to even larger values for the scaling relation of [49].

As can be seen in Figure 9, the full-cycle frequency shifts obtained by employing the scaling relation of [49] would lead to very large values of over 50 $\mu$Hz for stars reaching the lower red giant branch. This is why we looked for reasonable adjustments. As the sample of stars on which [50] based their results is exclusively comprised of main-sequence stars, we modified the $R^{\prime}_{\text{HK}}$ index with the factors discussed above. This aims to capture some of the evolution that can naively be expected in stellar magnetic activity. Here, the decrease of activity with $t_{\text{Age}}^{-0.5}$ found by [53] is again extrapolated from the main-sequence up to the lower red giant branch. The inclusion of this factor must be understood as an attempt to circumvent our ignorance of how exactly activity evolves after the main-sequence. It is completely possible that activity does indeed decrease faster than just proportional to $t_{\text{Age}}^{-0.5}$ or even stops completely at a certain evolutionary stage. The same applies for the factor $\max\left(\sqrt{\text{Ro}-1},1\right)$.

Contemporaneous measurements of the $R^{\prime}_{\text{HK}}$ index and photometric time series for asteroseismology would facilitate a better differentiation between these scaling relations. If indeed $\delta\nu\propto A\cdot\left(R^{\prime}_{\text{HK}}\right)^{B}$, with $A,B$ factors dependent on stellar fundamental parameters regulating mode sensitivity to perturbations ($A$) and the amplitude of activity over a full cycle ($B$), then measurements of frequency shifts and activity level can be used to learn about mode physics as well as stellar dynamos. We note that measurements of stellar rotation periods are essential for the prediction of the $R^{\prime}_{\text{HK}}$ index with Equation (12) or with the modified expression (14) should direct spectroscopic measurements of stellar activity not be available.

The estimation of the cycle period as described by Equation (17) is obviously very crude. If cycle periods are underestimated with this, then the resulting shifts we can expect to observe decrease and our predictions for $\delta\nu$ are an upper limit in this case. More, reliably measured cycle periods are needed to replace this estimation with a more rigorous prediction. We remark however, that Equation (17) does include the latest findings concerning the evolution of cycle periods by [26] through the lengthening of cycles with the factor $\max\left(\sqrt{\text{Ro}-1},1\right)$.

7 Potential of TESS

We evaluated the frequency shifts obtained from four settings: the TESS mission is either extended by two years or by four years and the frequency shifts are determined by the cross-correlation method (CC) or the peak-bagging method (PB), where the PB shift results are simply taken as those of the CC method but with half the uncertainty as described in Section 5. The predicted number of detections of significant frequency shifts for different subsets of the cohort of TESS ATL short-cadence stars is given in Table 1 for a four year mission extension and Table 4 for a two year mission extension. We define significant frequency shifts as those which are at least $1\sigma>0$, and within $3\sigma$ of the input frequency shift $\delta\nu_{\text{input shift}}$, where $\sigma$ is the uncertainty of the measured frequency shift. Tightening the second condition to $2\sigma$ or $1\sigma$ would eliminate some of the more precisely but less accurately measured frequency shifts.

We separated the spectral types according to their effective temperature with boundaries at $T_{\text{eff}}=5200\,\mathrm{K}$ between G and K and at $T_{\text{eff}}=6000\,\mathrm{K}$ between F and G. There are no K-type main-sequence stars in the sample. The separation into region I and II is discussed in Section 3. Table 1 gives the number of stars for each subset.

The suppression of p-mode amplitudes by activity has not been taken into account in our simulations. However, only stars whose oscillations have a low detection probability and small intrinsic mode amplitudes are likely to be affected by this. These are mainly late-type dwarf stars whose number in the mock sample is rather small, thus the overall number of expected detections is not strongly affected by this. As can be seen in Tables 2 and 3 the fraction of stars with detectable shifts increases with the length of the time series as frequency resolution increases. However, even for stars which are observed for only one sector, it should be possible to detect shifts at a significant level, especially for stars with long-lived modes.

With the CC method we were able to measure significant frequency shifts in 171 stars in region I and in 2275 stars in region II for a four year TESS mission extension. The large difference between the two regions reflects the increased amplitude and longer mode life times of solar-like oscillations for stars in region II. Also, the much smaller mode widths of cooler stars allow one to measure even smaller shifts, especially with the CC method. Figure 3 shows a Kiel diagram of the stars with significant shifts measured with the CC method for the four year extension. By comparison of Figure 1 with Figure 3 it can be seen that stars with detectable frequency shifts are largely found at effective temperatures below $\lesssim$6500 K. Note the different colour tables used in these two figures. For region I, the fraction of F-type stars for which we predict significant shifts is only 2.1% (141 out of 6676 stars) and only 1.5% (24 out of 1642) for region II F stars. This reflects the predicted smaller full cycle frequency shifts and large mode widths of F stars. The fraction of stars with predicted detectable frequency shifts increases for G stars: 55.6% (30 out of 54) for region I and 47.0% (765 out of 1628) for region II. For the G dwarfs of region I this increase is mainly due to the smaller mode line widths compared to F dwarfs. Furthermore, for region II G-type stars the detectability of shifts is enhanced by the larger predicted full cycle shifts and the larger mode amplitudes. A large fraction of K-type region II stars are predicted to have detectable shifts with 54.4% (1486 out of 2731). With the smaller uncertainties which we assumed for the PB method, these numbers increase for all subsets of stars, adding up to a total of 3804 out of 12731 stars (29.9%) for the complete cohort. Table 4 gives the predicted number of detections for a 2 year TESS extension.

We also tested for the number of false positives caused by the realisation noise. For this we set the input shift to zero for all stars and proceeded as described above. We found no false positives detections with the PB method for either a two or a four year extension even when ommiting the condition that the measured shift be within $3\sigma$ of the input shift value.

In Figures 4 and 5 the significant frequency shifts obtained with the cross-correlation method for a four year TESS extension are plotted for three subsets of stars (rows from top to bottom: all stars, region I, region II) as a function of effective temperature (left column of panels) and rotation period (right column of panels). The colour of the dots gives stellar age in Figure 4 and stellar mass in Figure 5. We capped the ordinate in these plots at $5\,\mathrm{\mu Hz}$. Higher frequency shifts tend to have large uncertainties and including them would only distort the presentation.

By comparing panel (a) of Figure 4 with Figure 2 it can be seen that for older stars, as stellar cycles grow longer with stellar age, only a fraction of the predicted shift over a full activity cycle can be detected with data that are separated by only four years. As expected, younger and faster rotating stars tend to show larger detectable shifts. The older the stars become (green and blue data points), the smaller the shifts that can be expected to be observed after two years become. This can also be seen in panel (b), where the measured shifts are plotted as a function of rotation period: as stars age, their rotation period becomes longer, the level of activity decreases, and the frequency shifts become smaller. Obviously, this by design for our model for the full-cycle frequency shifts (Equation 15) with the factors $t_{\text{age}}^{-0.5}$ and the decreasing rotation period with age as it is modelled in the synthetic TESS sample (cf. Section 2.5 and Figure 2 of [41]). It is noteworthy that for more evolved stars, which have smaller mode damping widths, shifts are detected more precisely and smaller shifts are detectable at a significant level. The lack of stars with a large frequency shift at short rotation periods in panel (b) of Figures 4 and 5 is due to the fact that most of these young, faster-rotating stars in the synthetic TESS sample are hot stars ($T_{\text{eff}}\gtrsim 6500\,\mathrm{K}$) with broader p-mode peaks for which shifts are harder to measure precisely. Also, the full cycle frequency shifts of these stars are rather small as can be seen from Figure 2.

The discrete structure in the measured frequency shifts in Figures 4 and 5 is due to the large number of stars with only one sector of observation (7930 stars). The frequency resolution of one sector’s worth of observation is 0.422 $\mu$Hz. As the cross-crorrelation function is computed from periodograms with this resolution, the fitted shift value is likely to be a multiple of this value.

For the stars in region I (panels (c) and (d)) there is no clear tendency in frequency shift with effective temperature or rotation period. This is partly due to the way that rotation periods are modelled in the synthetic sample. For stars with $T_{\text{eff}}>6500\,\mathrm{K}$, most of which are in region I, the rotation period is drawn from a normal distribution with mean 5 days and a standard deviation of 2.1 days, i.e., there is no dependency on age or effective temperature. Overall, only a small percentage of targets with $T_{\text{eff}}>6500\,\mathrm{K}$ have predicted detections. The same explanation for this as given above applies. Also, stellar mass introduces a selection bias here: due to their shorter evolutionary time scales, the more massive stars in this sample do not have the time to brake their rotation and thus stars with $M\gtrsim 1.8\,\mathrm{M_{\odot}}$ for which detection of frequency shifts are predicted, are all at rotation periods $P_{\text{rot}}\lesssim 20\,\mathrm{d}$ for stars in region I (panel d) and at $P_{\text{rot}}\lesssim 40\,\mathrm{d}$ for stars in region II (panel f). Similarly, only lower mass stars evolve slowly enough for rotation period to increase and the activity to decrease: note the accumulation of old, low-mass stars at $T_{\text{eff}}\approx 5000\,\mathrm{K}$ in panels (e), which corresponds to the base of the red giant branch seen in Figure 1.

7.1 Comparison to real observations

To test the agreement between our predictions and real observations we concentrate on main-sequence stars as there are currently no measured frequency shifts of more evolved stars. We use the results of [21] available online at this link. This is the largest available record on activity related frequency shifts to date. [21] employ a PB approach on segments of Kepler time series with a length of 90 days which overlap by 45 days. We calculate the frequency shift for each star in their sample by taking the difference between the minimum and maximum value of the mean frequency shift, which combines their results for modes with harmonic degree $l=0$ and 1.

The left panel of Figure 6 shows the shifts of [21] in blue and the shifts of the region I stars from this work in red as a function of effective temperature. As can be seen, our predictions cover essentially the same range of shifts as the measurements. The baseline of the Kepler data used by [21] is typically around three years, which is shorter than the baseline we assumed for our simulations. Thus, when assuming that the Kepler data do not cover a full activity cycle, extrapolating onto the same baseline would give somewhat larger shift values for the measured frequency shifts compared to our predictions. However, our sample is dominated by stars closer to the TAMS as can be seen in the right hand panel of Figure 6, which shows a Kiel diagram of the stars of [21] in blue and our region I stars in red. This biases the predicted shifts more towards lower values, as these stars tend to have slowed down their rotation compared to their younger counterparts of similar mass and effective temperature. Restricting our sample further would however leave us with only a very small number of stars. We will have to wait for data of a hopefully extended TESS mission to be able to truly test our scaling relation Equation (15) and those of [47] and [49] on a larger sample of stars.

8 Overestimation of global stellar parameters

The full-cycle frequency shifts as shown in Figure 2 have typical values of $6\,\mathrm{\mu Hz}$ for a low-luminosity red giant. For such a star, let $\nu_{\text{max}}=324\,\mathrm{\mu Hz}$, $\Delta\nu=23.3\,\mathrm{\mu Hz}$, $T_{\text{eff}}=5058\,\mathrm{K}$, $M=1.064M_{\odot}$ (target 00197 in the sample of [41]). Using the scaling relation (see, e.g., Kjeldsen and Bedding [45])

[TABLE]

the mass value given above is recovered. The solar reference values are those used in [46], which are $\nu_{\text{max, }\odot}=3090\,\mathrm{\mu Hz}$, $\Delta\nu_{\odot}=135.1\,\mathrm{\mu Hz}$, $T_{\text{eff,}\odot}=5777\,\mathrm{K}$. If $\nu_{\text{max}}$ is shifted by $6\,\mathrm{\mu Hz}$ to higher frequencies, a mass of $M_{\text{shift}}=1.12M_{\odot}$ is obtained from the scaling relation Thus, if $\nu_{\text{max}}$ were measured during times of higher activity, stellar mass would be overestimated by $0.06\,\mathrm{M_{\odot}}$ with this scaling relation. The statistical uncertainties of asteroseismic stellar masses obtained from the standard scaling relations are typically of the order of 10% [58, 59]. Thus, an activity-related shift to higher frequencies would introduce a systematic error and lead to an overestimation of mass and radius when using the standard scaling relations. Indeed, [60] found in their study of 10 eclipsing binary systems with at least one oscillating red giant that the dynamical masses were consistently smaller than those obtained with the asteroseismic scaling relations. The same arguments apply for stellar radii. We note that for more evolved red giants, a shift to higher frequencies of only $1\,\mathrm{\mu Hz}$ would lead to an even larger overestimation: Considering the values of [60] for the star KIC 5786154 ($\nu_{\text{max}}=29.75\,\mathrm{\mu Hz}$, $\Delta\nu=3.523\,\mathrm{\mu Hz}$, $T_{\text{eff}}=4747\,\mathrm{K}$) and using the solar reference values and scaling relation as given above, we find a mass for this star of $M=1.44M_{\odot}$. With $\nu_{\text{max}}=30.75\,\mathrm{\mu Hz}$, which incorporates the frequency shift predicted here, the scaling relation gives $M=1.58M_{\odot}$. This example shows that it might be necessary to calibrate the global seismic parameters $\nu_{\text{max}}$ and $\Delta\nu$ for the level of activity of the star when the scaling relations are to be used for an estimation of mass and radius of more evolved stars. A paper on activity-related frequency shifts in evolved Kepler stars is currently in preparation. As we will show there, frequency shifts of a few tenths of $\mu$Hz are not uncommon on the red giant branch.

9 Discussion and Conclusion

To estimate the yield of detections of activity-related p-mode frequency shifts in an extended TESS mission, we first derived a scaling relation for full-cycle shifts (Equation 15) for stars from the main-sequence up to low luminosity red giants. We used a catalogue of simulated light curves (Ball et al. [41]) and input frequency shifts which can be expected to occur within the two (four) years which would pass between a first observation of a given star and a repeated observation by TESS. We then measured these frequency shifts with a cross-correlation technique and found that we can expect to find a couple hundred main-sequence and early subgiant stars and a few thousand more evolved subgiant and low-luminosity red giant stars with activity-related p-mode frequency shifts in an extended TESS mission according to our scaling relation.

Our scaling relation does not account for the varying sensitivity of modes of different harmonic degree: for the Sun it is observed that modes of higher harmonic degree experience larger shifts over the solar cycle [e.g., 5]. For several Kepler stars, larger shifts have been measured for $l=1$ modes than for $l=0$ modes by [21]. A dedicated study of how the frequency shifts of different harmonic degrees can influence the estimation of fundamental stellar parameters is presented in this collection by Pérez Hernández et al. [61].

Up until now there is no confirmed detection of frequency shifts in subgiant and low-luminosity red giant stars. This may either be because their activity cycles are too long to be detectable with the available time series of Kepler and CoRoT, they no longer have activity cycles, or nobody has looked for this phenomenon in these stars yet. We also note that it is entirely possible that the scaling relation we have derived here is not valid for stars that have evolved far off the main-sequence even with the adjustments we have included to account for suppression of activity, see Section 6.

A dedicated search for activity-related frequency shifts (or frequency offsets in the case of very long cycles, which could also be attributed to activity see Section 8), should therefore be carried out with the existing Kepler data sets. Such results can then also be used to refine the estimates we presented here. Given our scaling relation, Equation (15) and the large number of confirmed oscillating subgiants and red giants (e.g., [62, 63]) such a study could also help to decide between the scaling relations for the p-mode frequency shifts by Chaplin et al. [47], Metcalfe et al. [49], and the one presented here.

So far, only a few dozen solar-like stars, most of them Kepler targets, have measured p-mode frequency shifts. An extended sample of seismic detections of stellar activity and activity cycles will also be useful to determine the fundamental parameters which govern stellar dynamos and their evolution as stars age. [64] have shown that it is possible to measure stellar differential rotation with asteroseismology. This can be connected with photometric measurements of activity and the temporal variation of p-mode frequencies to learn more about the intricate connection of rotation, its change through stellar evolution, and stellar activity cycles.

Carrying out searches for and investigations of stellar magnetic activity with seismic measures (mode frequencies, and amplitudes) are greatly enhanced in their reliability by contemporaneous spectroscopic measurements of magnetic activity. As, e.g., [65] have shown, it is feasible to employ seismology together with ground-based observation to study stellar cycles and activity in detail. Stars near the ecliptic poles, which are observed almost continuously during TESS’s observation (615 stars in the synthetic sample have 13 sectors of data; a similar number is expected for the real observations) can be used for equivalent studies to fill the $P_{\text{rot}}$-$P_{\text{cyc}}$ diagram. Long-term ground-based spectroscopic measurements of the level of activity of these targets would be expedient for such studies. The rationale to optimally select targets for such campaigns was laid out by Karoff et al. [57, 66]. Furthermore, the investigation of stellar variability can strengthen the seismic detection of magnetic activity as, e.g., [9, 17] have shown.

Our predictions for frequency shifts measured after four years for stars on the main-sequence and early subgiant stars (region I) agree qualitatively with those found by [19] and Santos et al. (submit.) for Kepler stars of similar evolutionary states. The reason the decrease of shifts with decreasing effective temperature seen by Santos et al. (submit.) is less clear in our simulations is most likely due to a difference in the underlying sample of stars. We will have to await TESS data from a hopefully extended mission, calculate models of the observed stars, measure the frequency shifts, and then compare predictions from our scaling relation and those of [47] and [49] to these observed shifts. As also mentioned above, measurements of a spectroscopic activity proxy like $R^{\prime}_{\text{HK}}$ and reliable measurements of the stellar rotation periods and age would facilitate this considerably.

We found that it should be possible to reliably measure even small average frequency shifts between two TESS observations each with a length of only one month. This should also enable the search for photometric and asteroseismic signatures of short activity cycles in those young, fast rotating stars which have many sectors of data in TESS observations (see [3, 67] for such studies with Kepler data). Obviously, long and uninterrupted photometric time series with high stability are optimal for asteroseismology in general and specifically for the search for seismic signatures of activity cycles. Even though most of the targets of TESS only have few consecutive sectors of observation, time series provided by TESS, especially by an extended mission, can be used for such searches. In the future, PLATO [68] will pick up the torch, observe two target fields for 2-3 years each, and provide data for asteroseismology. These data will make another major contribution in the seismic sounding of stellar activity cycles.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Author Contributions

RK led the work on this article and wrote most of the text. AMB contributed to the writing of Sections 2, 5, and the overall refinement of the writing. WHB generated and provided the dataset which is utilised in this article, computed additional necessary quantities, and helped in the writing of the relevant sections.

Funding

This work was funded by STFC consolidated grant ST/L000733/1. WHB acknowledges support from the UK Science and Technology Facilities Council (STFC).

Acknowledgements

RK & AMB acknowledge the support of the STFC consolidated grant ST/L000733/1.

Data Availability Statement

The datasets used for this study can be found in this repository: Link.

Appendix A Additional plots and tables

[FIGURE:]

[FIGURE:]