Starspot rotation rates vs. activity cycle phase: Butterfly diagrams of Kepler stars are unlike the Sun's

TL;DR

This study analyzes starspot rotation rates and activity cycles in Kepler stars, revealing diverse butterfly diagram patterns that differ from the Sun's, influenced by rotation period and temperature, and suggests these may be secondary cycles.

Contribution

It provides the first large-scale analysis of starspot rotation rate variations over activity cycles for Kepler stars, highlighting different butterfly diagram behaviors from the Sun.

Findings

Slow rotators show anti-correlation between rotation rate and activity.

Stars with shorter periods and higher temperatures show positive correlation.

Star activity cycles in Kepler stars are likely shorter and possibly secondary.

Abstract

During the solar magnetic activity cycle the emergence latitudes of sunspots change, leading to the well-known butterfly diagram. This phenomenon is poorly understood for other stars as starspot latitudes are generally unknown. The related changes in starspot rotation rates caused by latitudinal differential rotation can however be measured. Using the set of 3093 Kepler stars with activity cycles identified by Reinhold et al. (2017), we aim to study the temporal change in starspot rotation rates over magnetic activity cycles, and how this relates to the activity level, mean rotation rate, and effective temperature of the star. We measure the photometric variability as a proxy for the magnetic activity and the spot rotation rate in each quarter over the duration of the Kepler mission. We phase-fold these measurements with the cycle period. We perform averages over stars with comparable mean rotation rates and effective temperature at fixed activity-cycle phases. We detect a clear correlation between the variation of activity level and the variation of the starspot rotation rate. The sign and amplitude of this correlation depends on the mean stellar rotation and, to a lesser extent, on the effective temperature. For slowly rotating stars (with periods between 15-28 days) the starspot rotation rates are clearly anti-correlated with the level of activity during the activity cycles. A transition is seen at periods of 10-15 days, where stars with effective temperature above 4200K instead show positive correlation. Our results can be interpreted in terms of a stellar "butterfly diagram", but these appear different from the Sun's as the starspot rotation rates are either in phase or anti-phase with the activity level. Alternatively, the activity cycles seen by Kepler are short (around 2.5 years) and may therefore be secondary cycles, perhaps analogous to the solar quasi-biennial oscillations.

Figures

Tables

Table 1. Table 1: Characteristics of ⟨ Δ ​ A ​ ( t ) / A ¯ ⟩ delimited-⟨⟩ Δ 𝐴 𝑡 ¯ 𝐴 \langle\Delta A(t)/\overline{A}\rangle and ⟨ Δ ​ Ω ​ ( t ) / Ω ¯ ⟩ delimited-⟨⟩ Δ Ω 𝑡 ¯ Ω \langle\Delta\Omega(t)/\overline{\Omega}\rangle in period ranges I, II, III, and IV, including the Sun for reference.

Group	min($⟨\Delta A\left(t\right)/\overline{A}⟩$) [%]	max($⟨\Delta A\left(t\right)/\overline{A}⟩$) [%]	$⟨\Delta A\left(t\right)/\overline{A}⟩$ phase at max.	$⟨\Delta \Omega \left(t\right)/\overline{\Omega }⟩$ amplitude [%]	$⟨\Delta \Omega \left(t\right)/\overline{\Omega }⟩$ phase at min.	$r_{av}$
I, K	$-40.90\pm 4.80$	$52.05\pm 7.39$	$0.45$	$0.59\pm 0.13$	$0.53\pm 0.18$	$-0.50$
I, $T_{\mathrm{eff}}\ge 4200$K	$-32.91\pm 0.93$	$45.02\pm 1.77$	$0.50$	$0.33\pm 0.04$	$0.98\pm 0.13$	$0.76$
II	$-37.09\pm 1.01$	$55.75\pm 2.14$	$0.50$	$0.16\pm 0.05$	$0.57\pm 0.31$	$-0.45$
III	$-41.19\pm 0.89$	$67.92\pm 2.08$	$0.50$	$1.70\pm 0.07$	$0.49\pm 0.04$	$-0.95$
IV	$-49.71\pm 2.27$	$91.10\pm 6.03$	$0.50$	$2.86\pm 0.25$	$0.50\pm 0.09$	$-0.94$
The Sun	$-70.58\pm 5.03$	$77.95\pm 22.88$	$0.50$	$1.24\pm 0.09$	$0.05$	$0.15$

Equations

$$\delta\nu=\left\{{\begin{array}[]{*{20}c}\alpha&\rm{for}\quad 0.5\,\nu_{rot}<\alpha\\ 0.5\nu_{rot}&\rm{for}\quad 0.5\,\nu_{rot}\geq\alpha\\ \end{array}}\right.$$

Full text

*Aims. *Using the set of 3093 Kepler stars with activity cycles identified by Reinhold et al. (2017), we aim to study the temporal change in starspot rotation rates over magnetic activity cycles, and how this relates to the activity level, the mean rotation rate of the star, and its effective temperature.

*Methods. *We measure the photometric variability as a proxy for the magnetic activity and the spot rotation rate in each quarter over the duration of the Kepler mission. We phase-fold these measurements with the cycle period. To reduce random errors we perform averages over stars with comparable mean rotation rates and effective temperature at fixed activity-cycle phases.

*Results. *We detect a clear correlation between the variation of activity level and the variation of the starspot rotation rate. The sign and amplitude of this correlation depends on the mean stellar rotation and – to a lesser extent – on the effective temperature. For slowly rotating stars (rotation periods between $15-28$ days) the starspot rotation rates are clearly anti-correlated with the level of activity during the activity cycles. A transition is observed around rotation periods of $10-15$ days, where stars with effective temperature above $4200$K instead show positive correlation.

*Conclusions. *Our measurements can be interpreted in terms of a stellar “butterfly diagram”, but these appear different from the Sun’s since the starspot rotation rates are either in phase or anti-phase with the activity level. Alternatively, the activity cycle periods observed by Kepler are short (around $2.5$ years) and may therefore be secondary cycles, perhaps analogous to the solar quasi-biennial oscillations.

1 Introduction

On the Sun, the latitudinal migration of sunspots produces the well-known butterfly diagram (Maunder, 1904). As the solar magnetic activity cycle progresses, the spot coverage of the Sun increases, while the typical emergence latitude moves closer to the equator. After activity maximum the spots continue to migrate further toward the equator, eventually stopping at approximately $8^{\circ}$ latitude, at which point the next 11 year cycle begins again with spots appearing at around $30^{\circ}$ latitudes. The emergence latitudes covered by sunspots span a range of rotation rates, due to the latitudinal differential rotation of the Sun. Observations of spots and active regions are therefore among the methods that have been used to establish the solar latitudinal differential rotation profile (see, e.g., the review by Beck, 2000). The differential rotation profile, combined with the butterfly diagram are clear observable features of the solar magnetic dynamo, and their characteristics therefore act as strong constraints on dynamo models (see, e.g., reviews by Rempel, 2008; Charbonneau, 2010).

Other stars also exhibit magnetic activity, ranging from young fast rotators with irregular magnetic activity through older, slower, rotators that show smooth periodic cycles like the Sun, to in some cases not showing any long term variability at all (Wilson, 1978; Baliunas et al., 1995; Hall et al., 2009). These activity cycles are coupled with a brightening or dimming of the star. Relatively old and slow rotators like the Sun show a brightening of $\sim 0.1\%$ along with increased photometric variability during activity maximum (Froehlich, 1987; Hall et al., 2007). This is in spite of an overall increase in the number of dark spots on the stellar surface, and is attributable to a similarly enhanced number of faculae and plage regions (Foukal & Vernazza, 1979). However, previous studies (Radick et al., 1998; Hall et al., 2009) found that younger and more active stars than the Sun predominantly tend to show the inverse behavior, namely a dimming during activity maximum rather than a brightening. Radick et al. (1998) suggested that this change is due to younger, active stars exhibiting spot-dominated variability while older, less active stars eventually become faculae dominated. Montet et al. (2017) studied photometric variability in 463 solar-like stars (in terms of $\,\log{g}$ and $\,T_{\mathrm{eff}}$) and detected a transition from spot dominated (dimming during variability maximum) to faculae dominated (brightening during variability maximum) photometric variability at a rotation period of $\sim 10-20$ days. This falls neatly in line with efforts to model the variability of the total solar irradiance (TSI), which show that the Sun is predominantly faculae dominated (Chapman, 1987; Shapiro et al., 2016).

As with the Sun, the spot rotation rate of stars can also be expected to vary over the course of their respective activity cycles. While still challenging, measurements of latitudinal differential rotation have now been performed on a variety of other stars. The methods employed are wide ranging and include: the study of power spectra of photometric light curves (e.g. Lanza et al., 1993; Reinhold et al., 2013; Reinhold & Gizon, 2015), time variations of spectral line profiles with Doppler imaging (Donati & Collier Cameron, 1997; Collier Cameron et al., 2002; Barnes et al., 2001) and Zeeman Doppler imaging (Petit et al., 2004; Barnes et al., 2005), changes in starspot rotation rates (Henry et al., 1995; Messina & Guinan, 2002), variability in the magnetically sensitive Ca HK lines (e.g. Donahue et al., 1996), and more recently using asteroseismology (Benomar et al., 2018). However, since activity cycles tend to have periods of several years (Baliunas et al., 1995), only a few studies have been able to trace rotation periods over the course of one or more cycles (Donahue & Baliunas, 1992; Messina & Guinan, 2003).

In this work we extend the sample of stars with both measured activity cycles and rotation periods to include a recently published activity cycle catalog by Reinhold et al. (2017). We then study the change in the spot rotation rate over time and how it relates to the activity cycles of the stars.

2 Stellar activity cycles from Kepler data

Reinhold et al. (2017) use the photometric variability of Kepler light curves as a proxy for magnetic activity. This method is motivated by solar observations showing increased variance of the measured flux during solar activity maximum, and conversely a decrease during activity minimum. In the Sun this enhanced variability stems from an increase in the number of dark spots, bright faculae, and plage regions and their evolution over time. The assumption is that other stars show a similar cyclic change in photometric variability through the activity cycle.

The photometric variability used by Reinhold et al. (2017) was originally defined by Basri et al. (2011) as the $5$th to $95$th percentile range of the relative flux, and is computed for each quarter ($\sim 90$ days) of Kepler observations over the course of the mission lifetime ($\sim 3.4$ years). Reinhold et al. (2017) fit a sine function to the variability measurements, which was tested against a $5\%$ false alarm probability. The catalog was based on the rotation catalog by McQuillan et al. (2014) which contains rotation periods of 34030 stars. Out of these Reinhold et al. (2017) find a total of $3203$ stars that show significant periodicity of their photometric variability. The measured period of the sine function is then assumed to be the period of the activity cycle. The measured cycle periods range from $0.5-6$ years, with a median of $\sim 3$ years.

Figure 1 shows the effective temperature $\,T_{\mathrm{eff}}$ and surface gravity $\,\log{g}$ for these $3203$ stars (adapted from Huber et al. (2014), along with their rotation periods from McQuillan et al. (2014). The sample predominantly consists of cool main-sequence stars, with a range of rotation periods spanning $1-41$ days, and a median of $P_{\mathrm{rot}}\approx 16$ days; the temperature range is from approximately $3200$K to $6900$K with a median of $\,T_{\mathrm{eff}}\approx 5500$K. Based on the trend of slower rotation with increasing age (see, e.g., Barnes & Kim, 2010), this sample of stars is likely on average slightly younger than the Sun.

An inspection of the cycle periods revealed spuriously high numbers of stars with periods of either one year or $185$ days. We checked the light curves of all stars with these cycle periods. For stars with $P_{\mathrm{cyc}}\approx 1$ year the light curves often show discontinuous jumps in the relative flux amplitude every $4$ quarters. This is most likely caused by the quarterly roll of the Kepler spacecraft. Every $4$ quarters any given star will land on the same CCD, thereby being subject to the same instrumental noise sources, thus inducing a $1$ year periodicity in the variability measurements. We removed stars that showed this clear discontinuity in the flux variability, but retained stars that showed a smooth transition between quarters. For stars with $P_{\mathrm{cyc}}\approx 185$ days (77 stars) the cause is less certain, but it may also be related to instrumental effects since it constitutes approximately 2 quarters of observations. An inspection of the light curves did not reveal any remarkable differences compared to the stars in the remainder of the sample. However, we find it unlikely that such a spuriously large number of stars have the same cycle period, and therefore opt to remove them from the sample entirely. In total 110 stars were removed from the sample. This leaves 3093 stars for further analysis.

It should be noted that this sample is likely biased towards stars with significant and coherent photometric variability, as this forms the basis for the rotation measurement method employed by McQuillan et al. (2014). Montet et al. (2017) suggest that stars with rotation periods shorter than $P_{\mathrm{rot}}\approx 20$ are predominantly spot dominated. While this period range constitutes the majority of stars in our sample, we only find an overlap of $83$ stars with Montet et al. (2017). Of these, $80$ stars were defined as spot dominated. Since this only represents a fraction of our total sample we cannot say conclusively if the majority of stars in our sample are spot- or faculae-dominated.

In addition to the above bias, these are all stars with cyclic activity levels. Young, fast rotating stars with strong but incoherently varying activity or old, inactive stars with no variation at all, are therefore underrepresented. The range of cycle periods is between $0.5$ and $6$ years, so exact analogs of the solar activity cycle are also not included.

Searching the catalog of $\log{R\mathrm{{}^{\prime}_{HK}}}$ values collated by Karoff et al. (2016) we find 95 stars overlapping with our sample. These stars have an average $\log{R\mathrm{{}^{\prime}_{HK}}}=-4.569$, while the Sun shows a significantly lower $\log{R\mathrm{{}^{\prime}_{HK}}}=-4.895$. Again, this is not a sufficient overlap to draw conclusions about the average activity level of our sample. However, considering the high amplitude and coherent variability in the majority of our sample, we suggest that these stars are likely on average more active than the Sun.

3 Measuring stellar photometric variability and starspot rotation rate in each quarter

In this work we focus on measuring the amplitude and characteristic frequency of variability induced by magnetic activity, and how it varies over time with the activity cycle. These quantities are measured from the power spectrum111Computed using the Lomb-Scargle method (Lomb, 1976; Scargle, 1982) of the time series of each star from each observation quarter. These measurements are then phase-folded using the cycle period and phase values obtained from Reinhold et al. (2017), to study average changes in spot rotation rate and activity level (see Section 4).

We use the long-cadence Kepler light curves from all quarters of the mission. We use the PDC-msMAP pre-processed version of the light curves to avoid as much contamination of systematic variability as possible. Note that many of the stars in the sample have rotation periods longer than 20 days, beyond which the PDC-msMAP pipeline begins to significantly reduce the amplitude of intrinsic stellar variability (Gilliland et al., 2015).

For each star we start by defining a frequency range around the known average rotation rate $\nu_{rot}$, bounded by $\nu_{1}=\nu_{rot}-\delta\nu$ and $\nu_{2}=\nu_{rot}+\delta\nu$, where $\nu_{rot}$ is taken from McQuillan et al. (2014) (shown in Fig. 2). The half-width $\delta\nu$ of the frequency range is defined as

[TABLE]

For slow rotators, this definition of $\delta\nu$ minimizes the potential contribution of power from the 1st harmonic of the rotation rate. The value of $\alpha$ is set to $0.4\,\mu\mathrm{Hz}$, which is chosen in order to allow for variation of the rotation rate, while also reducing the contribution from instrumental variability on time-scales longer than the stellar rotation period.

We define the spot rotation rate as $\Omega=2\pi\nu_{C}$. Here $\nu_{C}$ is the centroid, which is estimated by the power-weighted sum of the frequency bins in the range from $\nu_{1}$ to $\nu_{2}$. For a single monolithic peak in the power spectrum (corresponding to a single, large spot or group of spots) the centroid will lie close to the frequency of maximum power of the peak. For cases where multiple frequencies show increased power (corresponding to spots at multiple latitudes) the centroid will represent a weighted average frequency of the variability. For stars with strong differential rotation, spots appearing at various latitudes will cause $\nu_{C}$ to vary around $\nu_{rot}$ in accordance with the spot emergence pattern.

The integral of the power $\mathcal{P}$ in this frequency range is taken as a proxy for the amplitude $A$ of the spot variability, such that $A=\left[{\int_{\nu_{1}}^{\nu_{2}}{\mathcal{P}\,d\nu}}\right]^{1/2}.$ This is similar to the photometric range used by Reinhold et al. (2017). However, unless properly filtered in frequency the photometric range is sensitive to additional sources of variability, such as granulation or shot noise. Here, $A$ is akin to the band-pass filtered photometric range. By confining the frequency range to a narrow region around the rotation rate, we ensure that the majority of the variability that contributes to $A$ stems from the rotation modulation.

We perform these measurement for each individual quarter, such that for each star we obtain $A(t_{q})$ and $\Omega(t_{q})$. Where $t_{q}$ is the median time of observation of each quarter $q$. For simplicity we drop the subscript $q$ in the following. The activity level $A(t)$ and $\Omega(t)$ are measured relative to their average values across all quarters, $\overline{A}$ and $\overline{\Omega}$ respectively. This gives us the relative change over time in photometric variability, $\Delta A(t)/\overline{A}$, and spot rotation rate, $\Delta\Omega(t)/\overline{\Omega}$.

Figure 3 shows $\Delta A(t)/\overline{A}$ and $\Delta\Omega(t)/\overline{\Omega}$ for the Sun and one of the stars in the sample: KIC002141150, which has a cycle period of $P_{\mathrm{cyc}}\approx 4$ years. The solar values are measurements of a composite set of TSI222Downloaded from ftp://ftp.pmodwrc.ch/pub/Claus/ISSI_WS2005/ISSI2005a_CF.pdf. PMOD/WRC composite v. 42_65_1709 measurements from multiple spacebased radiometers (Fröhlich, 2006). These observations span almost $42$ years, covering cycles $21$ through $24$.

The variance of the TSI is enhanced by variability from magnetic activity during activity maximum, and thus the integral $A$ that we are measuring is similarly enhanced. Using the Sun as a reference, in the following we therefore define the maximum of $\Delta A(t)/\overline{A}$ as the activity maximum and likewise for the activity minimum. However, it should be noted that we have no information concerning the variation of the chromospheric emission and other magnetic activity indices for the stars in our sample, and so this definition is only based on photometric variability.

For both the Sun and KIC002141150 $\Delta A(t)/\overline{A}$ shows a variation with the cycle period. On the other hand the spot rotation rate $\Delta\Omega(t)/\overline{\Omega}$ measurements for the Sun do not show any apparent variation. Sunspots tend to have lifetimes comparable the rotation period of the Sun (e.g., Petrovay & van Driel-Gesztelyi, 1997), which makes measuring the solar rotation rate from the TSI power spectrum very difficult. In addition, the scatter of spot emergence latitudes during the solar cycle spans up to $\approx 10$ degrees (Hathaway, 2010). This, combined with the short lifetimes makes measuring an average rotation rate challenging. Similarly, measuring the solar differential rotation from the TSI power spectrum is also challenging. In the following we therefore use the relative rotation rate of directly observed sunspots during the same period (shown in blue in the lower left frame of 3). These are computed from sunspot latitudes333RGO and NOAA measurements downloaded from https://solarscience.msfc.nasa.gov/greenwch.shtml and using the solar rotation profile by Snodgrass & Ulrich (1990). In contrast to the centroid values the variation of the directly observed sunspot rotation rates clearly vary with the activity cycle of the Sun.

For KIC002141150, however, the $\Delta\Omega(t)/\overline{\Omega}$ measurements do show variability that is on similar timescales as the cycle period, and anti-correlated with $\Delta A(t)/\overline{A}$, indicating latitudinal differential rotation. Here it is worth emphasizing that other stars in this sample do not necessarily show the same clear trends in $\Delta\Omega(t)/\overline{\Omega}$ as in $\Delta A(t)/\overline{A}$. Indeed, measuring stellar differential rotation using only photometric measurements is often difficult, mainly due the evolution of the spot signal as the spot grows and decays (see, e.g., Aigrain et al., 2015). In the power spectrum the effect of spot evolution is broadening of the distribution of power around the rotation peak. The lifetime of sunspots has indeed been shown to be correlated with the activity cycle (Henwood et al., 2010), however, when averaging over multiple stars we expect this effect to be symmetric around the rotation peak and will therefore not have a systematic effect on $\Delta\Omega(t)/\overline{\Omega}$.

4 Starspot rotation rate versus activity cycle phase

While attempting to measure the variation of the spot rotation rate with the activity cycle for individual stars is typically not feasible, especially if the spot lifetime is short, it is possible to do so when averaging the $\Delta\Omega(t)/\overline{\Omega}$ measurements of multiple stars. We first phase-fold the measurements444These will be made available as online material at the CDS. of all the stars using the cycle periods $P_{\mathrm{cyc}}$ and phases $\phi_{0}$ from (from Reinhold et al., 2017), and then compute the averages $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ in each phase bin. Here $\langle\rangle$ represents an average of quarterly measurements of multiple stars at the same phase of their respective activity cycles. Figure 4 shows $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ as a function of activity cycle phase, where the errors represent the standard errors on the mean. The bin widths in rotation period and effective temperature are $\approx 4$ days and $\approx 400$K respectively.

The mean rotation rate $\overline{\Omega}$ and the effective surface temperature, $\,T_{\mathrm{eff}}$, are thought to be two of the important parameters that influence the dynamo mechanism (Charbonneau, 2010, see, e.g.,). In Fig. 4 we have therefore divided the measurements of $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ into separate rotation and temperature bins, and show the variation of these parameters with cycle phase for each bin.

4.1 Shape of activity cycles

As expected from Fig. 3 the values of $\langle\Delta A(t)/\overline{A}\rangle$ clearly vary with the activity cycle, with generally a sharp peak followed by a shallow minimum. Note that by the construction of the catalog by Reinhold et al. (2017) the general shape of the variation cannot change dramatically, as the detection of the cycles were made using sine fits to the variability. A strong deviation from a symmetric profile might produce a false alarm probability greater than the $5\%$ imposed by Reinhold et al. (2017), and thus would have been removed from the sample. Note that for clarity in the following sections we have phase shifted our measurements of $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ compared to Reinhold et al. (2017), such that the minimum of $\langle\Delta A(t)/\overline{A}\rangle$ (for all the stars in the sample) is at phase $\phi=0$.

It is clear that the amplitude of $\langle\Delta A(t)/\overline{A}\rangle$ changes with the mean rotation period, and only marginally so for changes in temperature. For the fast rotators the change in $\langle\Delta A(t)/\overline{A}\rangle$ over the cycle is much smoother, and lower amplitude. In contrast, the slow rotators have a more peaked structure during activity maximum, with an amplitude relative to the mean of about a factor of 2-3 greater than the fast rotators. Compared to the fast rotators, the slow rotators also appear to have more pronounced excursions to high levels of variability at cycle maximum, whereas the levels of variability at cycle minimum in both cases remain similar.

Temperature appears to have a smaller impact on the variation of the photometric variability. Only the very hottest stars in the sample show any notable difference in $\langle\Delta A(t)/\overline{A}\rangle$. These stars are on average also faster rotators with a mean rotation period of $\approx 10$ day, and so the effect of averaging over many stars with less coherent variability could be the cause of the reduced amplitude that is seen here.

4.2 Detection of changes in starspot rotation rates with cycle phase

The change in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ appears to vary more dramatically with the mean rotation rate and temperature, than for $\langle\Delta A(t)/\overline{A}\rangle$. The slow and cool stars show the largest variation in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$, while the fast and hot stars show very little. The amplitude of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ for slow rotators is comparable to that of the Sun, but has a more sinusoidal form as opposed to the sharp decrease and slow increase seen for sunspot rotation rates. In general, the correlation between the $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ and the activity cycle is negative, that is, the observed spot rotation rate decreases around activity maximum and vice versa during activity minimum. This is different from the variation seen on the Sun (see Fig. 3).

It is important to note that the mean rotation period and the effective surface temperature are strongly correlated (see e.g., Nielsen et al., 2013; McQuillan et al., 2014), with cool stars being on average slower rotators. The observed change in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ with increasing temperature is therefore likely tied to the change in rotation period (see below).

4.3 Average correlation between $\Delta A(t)/\overline{A}$ and $\Delta\Omega(t)/\overline{\Omega}$

While the values of $\langle\Delta A(t)/\overline{A}\rangle$ indicate that a cycle is present at all temperatures and mean rotation periods, there is a stark change in the variation of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ as a function of cycle phase. To investigate this in more detail for our sample, we show in Fig. 5 the averaged Pearson correlation coefficient $\langle r\rangle$ of $\Delta A(t)/\overline{A}$ and $\Delta\Omega(t)/\overline{\Omega}$, for stars of different temperature and average rotation periods. We divide the period range into four ranges of particular interest, which we will focus on in following.

1. I

$P_{\mathrm{rot}}=5$-$11$ days: These stars predominantly show very little correlation between the activity cycle and changes in the spot rotation rate. However, a weak temperature dependence is visible. The hotter ($\,T_{\mathrm{eff}}>5000\mathrm{K}$) stars show a marginal positive correlation on average, while the cool stars are closer to $\langle r\rangle=0$ or marginally negative. 2. II

$P_{\mathrm{rot}}=11$-$15$ days: The stars in this range appear to transition from $\langle r\rangle\approx 0$ toward negative correlation on average. 3. III

$P_{\mathrm{rot}}=15$-$23$ days: These stars begin to show the same anti-correlation as the slower rotators. However, there is very little scatter due to differences in effective temperatures, compared to the faster rotators. 4. IV

$P_{\mathrm{rot}}=23$-$28$ days: Stars in this period range have similar rotation rates as the Sun, but also show the largest degree of anti-correlation of all the stars in our sample. Similar to the fast rotators, $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ and $\langle\Delta A(t)/\overline{A}\rangle$ for the hotter stars in this period range appear to be somewhat less anti-correlated than their cooler counterparts. As was seen in Fig. 4 the hot stars have slightly lower variability in $\langle\Delta A(t)/\overline{A}\rangle$. This may be caused by less coherent activity levels, which may explain the observed decorrelation. However, we note that there are only very few hot stars with these periods (there are no hot stars with periods longer than 25 days in our sample).

Stars with rotation periods longer than $\sim 30$ days appear to show a reduced correlation with the activity cycle, however the sample density in this period range is very low (only $86$ stars). Furthermore, the frequency window used to measure $\Omega$ reaches $\nu=0\,\,\mu\mathrm{Hz}$ at $\nu_{rot}=0.4\,\,\mu\mathrm{Hz}$, corresponding to $P_{\mathrm{rot}}\approx 29$ days. For rotation periods longer than this the measurements of the centroid may be biased if the latitudinal differential rotation is strong. It is therefore unclear if the upward trend starting at $\overline{P}=30$ days is real.

For reference we compute the corresponding correlation coefficient using the spot rotation rates from directly observed sunspots at various latitudes (blue curve in Fig. 3) over the 11 year solar cycle. The correlation coefficient is $\langle r\rangle\approx 0.22$. The reason for the positive correlation in the Sun is due to the migration of sunspot during the solar cycle. Near the start of the cycle sunspots initially appear at high latitudes. Following this the activity rises quickly, which leads to a positive correlation in the first phase of the cycle. Post-maximum, the activity level decreases while the spots continue to migrate toward faster rotating latitudes at toward the equator, thus leading to a negative correlation in the latter part of the cycle.

4.4 Relationship between the averages $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$

Figure 6 shows the change in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ as a function of cycle phase, while Fig. 7 shows it as a function of $\langle\Delta A(t)/\overline{A}\rangle$. Here we also show the correlation of the averages, $r_{av}$, between $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$. Note that in Fig. 6 the activity maximum occurs at phase $\phi\approx 0.5$ for all the period ranges (see Fig. 4). These figures will be discussed in the following.

1. I

$P_{\mathrm{rot}}=5$-$11$ days: This period range is represented in blue in Fig. 7, and the top left frame of Fig. 6. The correlation coefficients in this period range are temperature dependent, and so for clarity in the two figures we omit stars cooler than $4200$K, as the variation of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ for stars cooler than this is the same as in the remaining period ranges. The positive correlation seen for the hotter stars in this period range is clearly visible in Fig. 6, and appear roughly sinusoidal (indicated by the dashed gray sine-curve). 2. II

$P_{\mathrm{rot}}=11$-$15$ days: On average, $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ in this period range does not vary substantially and shows little if any correlation with the activity cycle. Dividing this period range into narrower period bins shows no difference between the stars with periods around $11$ days compared to $15$ days. The decrease in the amplitude of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ is therefore not likely due to averaging two populations in terms of period, with positive or negative correlation respectively. These stars still exhibit a variation in photometric variability, indicating an activity cycle is present. 3. III

$P_{\mathrm{rot}}=15$-$20$ days: These stars begin to show the same anti-correlated profile for $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ as is seen in the slower rotators. However, $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ appears have a more sinusoidal variation as indicated by the sine fit. That is, the transition of spots from rapid to slowly rotating latitudes and vice versa appears to be happening at a similar pace. Moreover, the minimum of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ occurs almost exactly at activity maximum ($\phi\approx 0.5$). Dividing this range into narrower bins in period or separating by temperature shows the same picture of a symmetric variation of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ with cycle phase. 4. IV

$P_{\mathrm{rot}}=22$-$28$ days: These stars show the strongest anti-correlation with the activity cycle, and also show the largest amplitude in both $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$. Most clearly seen in Fig. 7, there is a slight asymmetry in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ between the first and second halves of the activity cycle. During the rise of the cycle, the slope of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ is marginally steeper than during the declining phase of the cycle. Note however, that as with the faster rotators, the minimum in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ occurs at activity maximum.

For comparison we show the solar $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ in the inset in Fig. 7, again based on the 11 year activity cycle. This shows a dramatically different picture of the change in spot rotation rotation with the activity cycle of the Sun, compared to the stars in our sample. The solar $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ variation is also shown in blue in Fig. 6, where we have scaled the differential rotation coefficients from Snodgrass & Ulrich (1990) by a factor of $2.5$ for clarity. This illustrates that not only is the shape of the $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ curve different, but it is also significantly out of phase with the stars in our sample stars. This produces the relation between $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ for the Sun seen in Fig. 7.

The characteristics of the behavior between the activity cycle and the changes in the rotation rate are summarized in Table 5.

5 Discussion

Given the wide range of effects seen among stars of different rotation rates in our sample, we will separately discuss features seen among the fast rotators ($P_{\mathrm{rot}}<15$ days) and slow rotators ($P_{\mathrm{rot}}>15$ days). The discussion of the latter is further broken down based on two separate assumptions: that we are observing the primary magnetic activity cycle, i.e., the equivalent to the solar 11 year cycle; or, given that the cycle periods of these stars is typically only a few years, that we are observing a faster secondary cycle. This division is motivated by the fact that multiple periodicities are seen in both solar and stellar time series, and because the Kepler time series are not long enough to detect stellar cycles longer than about 6 years.

In the following we make the assumption that the differential rotation profiles for the stars are solar-like, in the sense that the equatorial rotation rate is greater than at the poles. It is well established from both theory and numerical models of stellar convection (Ruediger, 1989; Küker & Rüdiger, 2008; Miesch & Toomre, 2009) that stars with rotation periods $\lesssim 20$-$25$ days can establish a rapidly rotating equator compared to the polar regions (i.e. solar-like differential rotation). While strong Lorentz-force feedback from dynamo-generated magnetic fields may suppress the amplitude of the rotational shear in the fastest rotators (particularly for cool stars with deep convection zones), the sense of the rotation gradient should still be equatorward (Brown et al., 2008; Augustson et al., 2012; Guerrero et al., 2013; Gastine et al., 2014). For slowly-rotating stars with rotation periods longer than $\sim 20$ days, the assumption of solar-like differential rotation is less reliable. Global convection simulations suggest that the rotation gradient may reverse to being poleward at the surface when the Rossby number exceeds unity (Gastine et al., 2014; Featherstone & Miesch, 2015; Karak et al., 2015). However, the transition from solar to anti-solar differential rotation in convection simulations is rather abrupt and if it were occurring in our stellar sample, we would expect to see a dramatic change in behavior near a rotation period of about 25-30 days, depending on spectral type. Since we see no indication of such a dramatic transition, we will assume for the remainder of this section that all stars in our sample have a solar-like differential rotation.

5.1 Potential sources of sample contamination

A potential source of noise is if the cycle period found by Reinhold et al. (2017) is not caused by an activity cycle. Instead, the observed variation of the flux such as or similar to that seen for KIC002141150 in Fig. 3, could be caused by a beating effect between two or more spots at slightly different latitudes, and thus with slightly different rotation rates. This is caused by one spot ’lapping’ the other around the star. If in a given light curve, several beat periods are observed, this would be interpreted as a modulation in the photometric variability. With multiple spots on the stellar surface at once, several different spot configurations are possible but they may be roughly divided into two categories: (i) The difference in rotation rate between two active latitudes is small, leading to a long lap-time. If this lap-time is long compared to the typical spot lifetime, many spots will emerge and decay over the course of just a few lap-times. Since spot emergence is generally expected to be random in phase, the beat patterns would be similarly randomly distributed in time. We therefore anticipate the photometric variability to change stochastically, and such cases would likely fail the false-alarm probability test used by Reinhold et al. (2017). (ii) If the spots are long-lived compared to the lap-time, two spots may produce several beat periods. The observed change in the photometric variability would appear very periodic, and may therefore be falsely picked up as a cycle period. However, during one cycle (beat) the observed spot rotation rate as measured by the centroid of the power distribution in the frequency spectrum will not change. While the observed rotation rate might change between cycles (beats), when averaging over multiple cycles (beats) as is done here, the average change in the spot rotation rate will still be close to zero. We therefore expect such cases to contribute only to the number of stars with a very low degree of correlation, and not produce the observed changes in the spot rotation rates as seen in particular for the slow rotators. For the fast rotators however, this may be a source of contamination.

Furthermore, based on inferences from chromospheric and coronal emission, it has been suggested (Berdyugina, 2005; Strassmeier, 2009) that especially fast rotators may, during high activity phases, become completely saturated with spots on the stellar surface. If spots cover a substantial fraction of the stellar surface the photometric variability may decrease, even as magnetic activity increases. For cases such as this the photometric variability is no longer a reliable way of measuring the activity cycle. This may account for the reduction in the amplitude of $\langle r\rangle$ for rotation periods less than 11 days, as shown in Fig. 5.

5.2 Relationship between starspot rotation and activity among fast rotators

The tendency of individual stars among the fast rotators, with a very low degree of correlation between $\Delta A(t)/\overline{A}$ and $\Delta\Omega(t)/\overline{\Omega}$, may in part be caused by contamination effects as discussed above. However, when averaging $\Delta\Omega(t)/\overline{\Omega}$ over multiple stars the correlation with the activity cycle becomes much more readily visible (see Fig. 6 frame I). One interpretation of this change in rotation rate is in terms of latitudinal migration of emerging starspots. If this is the case, then the low correlation for an individual star could be caused by a large scatter in emergence latitudes at a particular phase of the cycle. The mean emergence latitude may still migrate in latitude, but only when observing multiple stars will this become apparent. In contrast to the slower rotators in our sample, the fast rotating stars show a positive correlation between $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$, with primarily the hotter stars showing this tendency. This effect is not unprecedented. Donahue & Baliunas (1992) found that the star $\beta$ Com shows a clear increase of the rotation rate with activity, similar to what is observed here. On the other hand, as shown by Fig. 5 the cool stars ($\,T_{\mathrm{eff}}\lesssim 4200$K) in this period range show the opposite tendency, namely a negative correlation of $r_{av}=-0.50\pm 0.13$ between $\langle\Delta A(t)/\overline{A}\rangle$ and $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$. This suggests that a temperature dependence exists in this period range, between activity level and spot rotation, which is not present in the slower rotators in our sample.

One possible explanation for the low degree of correlation for the individual stars in this period range is that strong Lorentz-force feedback reduces the differential rotation in times of high activity. Since the equatorial region has much more mass and rotational inertia, this could lead to an increase in the rotation rate at mid-latitudes when the activity level is high. This could produce a change in the latitudes at which the spots emerge and a change in their rotation rate.

5.3 Relationship between starspot rotation and activity among slow rotators

The slow rotators in our sample would be expected to be those that show the greatest degree of similarity to the Sun. However, the changes in the spot rotation rate over the cycle show a dramatically different picture. The Sun shows changes in the spot rotation rate with a sharp transition between successive cycles and a large phase-lag between rotation minimum and activity maximum, leading on average to a positive correlation between the spot rotation and activity cycle. The stars in our sample on the other hand show a very symmetric variation in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$, with a minimum at the phase of activity maximum. Assuming once more that the changes in starspot rotation rates are due to migration of active bands over differentially rotating latitudes, we suggest two different hypotheses based on two independent assumptions in order to explain the observations.

5.3.1 Hypothesis #1: Observed cycles are primary activity cycles

The first hypothesis is based on the assumption that we are observing cycles that correspond to the solar $\sim 11$ year cycle. In this case, assuming for the reasons described above that the differential rotation is solar-like, our findings indicate a butterfly diagram where the active regions emerge at low latitudes during the start of a cycle, at higher latitudes towards the peak of a cycle, and then at lower latitudes towards the end of a cycle. Figure 8 shows an illustration of what the spot migration patterns on these stars might look like. Another possibility is that the butterfly diagrams of individual cycles look like those of the Sun (with equatorial propagation over the course of a cycle), but with a greater degree of overlap between successive cycles than the Sun exhibits. In this case at the start of a cycle the spots from the previous cycle could still be determining the inferred rotation rate.

5.3.2 Hypothesis #2: Observed cycles are short secondary cycles

Given the differences in $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ between the Sun and stars, and also considering the relatively short cycle lengths of the stars, an alternative interpretation is that the observed cycles are shorter, secondary cycles and not equivalent to the solar 11 year cycle. Multiple stars have been shown to exhibit two distinct periods in magnetic activity (See et al., 2016; Oláh et al., 2016; Distefano et al., 2017), one which is typically on decadal timescales and another much shorter on timescales of a few years. Since the Kepler mission lifetime limits us to cycle periods on the order of a few years, a bias should be expected toward these shorter cycles.

The Sun also exhibits a secondary cycle, the so-called quasi-biennial oscillation (QBO). The solar QBO is a modulation of the radial magnetic field strength on $\sim 1.5-3$ year timescales. This modulation is strongest at high latitudes, i.e., it affects primarily the top of the solar butterfly-diagram wings (Ulrich & Tran, 2013). In principle, when the QBO suppresses the primary cycle it could therefore cause fewer spots to emerge at high latitudes. On the other hand, when the QBO enhances the primary cycle, sunspot emergence would be enhanced at high latitudes. In terms of an average spot rotation period, this would lead to a modulation of the observed spot rotation rate which is out of phase with the change in activity level, leading to the observed anti-correlation. Figure 8 shows an illustration of the modulation of the spot emergence latitude over a single primary cycle period.

The solar QBO is observable in multiple different proxies including: $10.7$cm radio flux (Valdés-Galicia & Velasco, 2008), X-ray emission (Antalova, 1994), acoustic oscillation frequency modulations (Fletcher et al., 2010), as well as in the TSI and its variance used in this work. However, in our measurements of the solar $\Delta\Omega(t)/\overline{\Omega}$ this modulation is not readily visible. The solar 11 year cycle modulates the measured $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ on the order of $1\%$, and so the QBO should be expected to produce even smaller variations in the spot rotation rate. In contrast, the amplitude of $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ observed for the slow rotators in our sample is on the order of $3\%$. While the stars in our sample are likely much more active than the Sun, this large amplitude of the stellar $\langle\Delta\Omega(t)/\overline{\Omega}\rangle$ is unexpected from a solar-like QBO.

In addition to the QBO’s, there is also the possibility of a non-axisymmetric quadrupole dynamo mode (e.g. Moss et al., 1995). There is some observational support for such a mode existing on the Sun (Berdyugina & Usoskin, 2003) and on other stars (Tuominen et al., 2002; Berdyugina & Järvinen, 2005) in terms of active longitudes in each hemisphere separated by 180∘ in phase. The mode is cyclic, with a flip-flop transition between the two active longitudes every 1.5 to 3 years in the case of the Sun (Berdyugina & Usoskin, 2003), similar to the period of the QBO’s. The interaction between this mode and the primary (axisymmetric) dipole mode of the star could also produce a butterfly diagram similar to that shown in the right panel of Fig. 8.

6 Conclusions

We measure the photometric variability as a proxy for magnetic activity, and star spot rotation rates as a function of time for $3093$ Kepler stars. These stars span a temperature range from $3200-6900$ and with rotation periods from $2-45$ days, and all have previously measured activity cycle periods from Reinhold et al. (2017) ranging from $0.5$ to $6$ years.

We find that the solar-like rotators tend to show a negative correlation between the activity level and rotation rate over the cycle period. For these stars the spots are at their slowest during activity maximum, and vice versa during activity minimum. The fast rotators ($P_{\mathrm{rot}}<20$ days) on the other hand tend toward a decorrelation between activity and rotation. In terms of Rossby number666$Ro=1/(2\Omega_{\mathrm{rot}}\tau_{\mathrm{cz}})$, where the convective turnover time $\tau_{\mathrm{cz}}$ is computed based on Noyes et al. (1984)., $Ro$, this transition occurs between $Ro\approx 0.02$-$0.08$. This decorrelation can be caused by contamination primarily from misidentified cycle periods, or possibly by suppression of differential rotation by Lorentz-force feedback effects. However, in a narrow period range between $5-11$ days, stars with $\,T_{\mathrm{eff}}\gtrsim 4200$K show positive correlation, where the spot rotation maximum occurs during activity maximum. Stars cooler than this show the same negative correlation as the slower rotators. This indicates a temperature dependence on the relation between activity and spot rotation in this period range, something that is not seen for the slower rotators.

For the slow rotators the negative correlation leads to two different hypotheses, under the assumption that the change in rotation rate is caused by spot migration. The first hypothesis is that the spot rotation rates change because of a butterfly-like migration pattern caused by the primary dynamo cycle. The observed change in spot rotation rate would however lead to a dramatically different migration pattern, looking more like a boomerang than a butterfly wing, or would imply that the cycles overlap considerably more than on the Sun. The second hypothesis is that we are observing a shorter secondary cycle. In the Sun this is called the quasi-biennial oscillation, which could in principle modulate the emergence latitudes of spots, and thus produce changes in the spot rotation rate that are anti-correlated with the activity level.

Directly measuring the latitudes of starspots is still not possible for stars of broadly varying spectral types, as the ones in our sample. Measuring the variation of the spot rotation rates on the other hand provides an intriguing look at what phenomena might be present on active stars with short activity cycles.

Acknowledgements.

The authors would like to thank Timo Reinhold, Alexander Shapiro, and Gibor Basri for useful discussions. MBN and LG acknowledge support from NYUAD Institute grant G1502. This paper includes data collected by the Kepler mission. Funding for the Kepler mission is provided by the NASA Science Mission directorate. Data presented in this paper were obtained from the Mikulski Archive for Space Telescopes (MAST). STScI is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. Support for MAST for non-HST data is provided by the NASA Office of Space Science via grant NNX09AF08G and by other grants and contracts.