# Measuring the Sun's radial velocity variability due to supergranulation over a magnetic cycle.

Niamh O'sullivan, Suzanne Aigrain, Michael Cretignier, Ben Lakeland, Baptiste Klein, Xavier Dumusque, Nadège Meunier, Sophia Sulis, Megan Bedell, Annelies Mortier, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

Niamh O'sullivan, Suzanne Aigrain, Michael Cretignier, Ben Lakeland, Baptiste Klein, et al. Measuring the Sun's radial velocity variability due to supergranulation over a magnetic cycle. Monthly Notices of the Royal Astronomical Society, 2025, 541 (4), pp.3942-3962. 10.1093/mnras/staf1168. hal-05486803

HAL Id: hal-05486803

https://hal.science/hal-05486803v1

Submitted on 22 Feb 2026

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /></figure>

# Measuring the Sun’s radial velocity variability due to supergranulation over a magnetic cycle.

Niamh K. O’Sullivan,¹* Suzanne Aigrain,¹ Michael Cretignier,¹, Ben Lakeland,² Baptiste Klein,¹, Xavier Dumusque,³ Nadège Meunier,⁴ Sophia Sulis,⁵ Megan Bedell,⁶ Annelies Mortier,²,⁸,⁹

¹Denys Wilkinson Building, Department of Physics, University of Oxford, OX1 3RH, UK

²School of Physics and Astronomy, University of Birmingham, Edgbaston, Birmingham B15 2TT, UK

³Observatoire astronomique de l'Université de Genève, 51 ch. des Maillettes, 1290 Versoix, Switzerland

⁴Univ. Grenoble Alpes, CNRS, IPAG, F-38000 Grenoble, France

⁵Aix Marseille Univ., CNRS, CNES, LAM, 38 rue Frédéric Joliot-Curie, 13388 Marseille, France

⁶Center for Computational Astrophysics, Flatiron Institute, 162, 5th Avenue, New York, NY 10010, USA

⁷Centre for Exoplanet Science / SUPA, School of Physics & Astronomy, University of St Andrews, North Haugh ST ANDREWS, Fife KY16 9SS, UK

⁸Centre for Exoplanets and Habitability, University of Warwick, Coventry CV4 7AL, UK

⁹Department of Physics, University of Warwick, Coventry CV4 7AL, UK

Accepted 2025 July 4. Received 2025 June 30; in original form 2025 April 7

## ABSTRACT

In recent years, supergranulation has emerged as one of the biggest challenges for the detection of Earth-twins in radial velocity planet searches. We used eight years of Sun-as-a-star radial velocity observations from HARPS-N to measure the quiet-Sun’s granulation and supergranulation properties of most of its 11-yr activity cycle, after correcting for the effects of magnetically active regions using two independent methods. In both cases, we observe a clear, order of magnitude variation in the time-scale of the supergranulation component, which is largest at activity minimum and is strongly anticorrelated with the relative Sunspot number. We also explored a range of observational strategies which could be employed to characterize supergranulation in stars other than the Sun, showing that a comparatively long observing campaign of at least 23 nights is required, but that up to 10 stars can be monitored simultaneously in the process. We conclude by discussing plausible explanations for the ‘supergranulation’ cycle.

**Key words:** methods: data analysis – techniques: radial velocities – Sun: granulation.

Over the next decade, extremely precise radial velocity (EPRV) instruments will in principle reach the precision needed to detect Earth-like planets. Since the detection of the first exoplanet around a Sun-like star in 1995 (Mayor & Queloz 1995), over 5800 exoplanets have been confirmed. Most were detected by space-based transit-search missions such as *Kepler* (Borucki et al. 2010) and the *Transiting Exoplanet Survey Satellite* (TESS; Ricker et al. 2014). The PLATO mission (Rauer et al. 2024), due to launch in 2026, will combine a *Kepler*-like collecting area with an extremely wide field-of-view in order to find Earth-like planets around nearby, Sun-like stars, suitable for further characterization. Radial velocity (RV) follow-up is key to confirm transiting planet candidates and to determine the planet’s mass. At the same time, blind EPRV searches employing a very intensive monitoring strategy will aim to detect non-transiting Earth analogues around the same stars (Thompson et al. 2016; Hall et al. 2018; Gupta et al. 2021). The key challenge in both cases is the intrinsic variability of the host star, which typically

overwhelms the planetary signals of interest, and can also mimic them in some cases (e.g. Rajpaul, Aigrain & Roberts 2016; Crass et al. 2021; Meunier 2024).

For the past 15 yr, precise RV spectrographs such as HARPS (Mayor et al. 2003) and HARPS-N (Cosentino et al. 2012) have made the detection of small planets, with RV semi-amplitudes down to around 1 m s⁻¹, a matter of routine (Cretignier et al. 2023; Dalal et al. 2024; Nari et al. 2025), but have struggled to reach below this ‘floor’. Their sensitivity is limited in part by night-to-night calibration errors of order 0.5 m s⁻¹ (Dumusque et al. 2021). This motivated the development of a new generation of ultrastable RV spectrographs equipped with Laser Frequency Combs (LFCs) for wavelength calibration, including EXPRES (Jurgenson et al. 2016), KPF (Gibson et al. 2016), NEID (Schwab et al. 2018), and ESPRESSO (Pepe et al. 2021), which should reach instrument stability of ~20 cm s⁻¹. However, even these state-of-the-art instruments are limited by intrinsic stellar variability, which introduces RV variations up to several m s⁻¹ on a wide range of time-scales. Next-generation RV surveys are thus planning to adopt very intensive monitoring strategies, to characterize these variations

<footnote>*E-mail: niamh.osullivan@physics.ox.ac.uk</footnote>
well enough to disentangle them from planetary signals from Earth analogues (with semi-amplitudes of order $10 \text{ m s}^{-1}$ and periods of several hundred days).

An example of such a survey is the Terra Hunting Experiment (THE; Thompson et al. 2016; Hall et al. 2018), which will use the HARPS3 spectrograph, an upgraded copy of HARPS and HARPS-N that will be installed on the Isaac Newton Telescope (INT) on La Palma. HARPS3 will be equipped with an LFC and will observe bright stars at high Signal-to-Noise Ratio (SNR) to reach a precision of ~30 cm s⁻¹ per epoch, but what sets THE apart is that it will focus on a comparatively small number of target stars (a few dozen), which will be monitored very intensively for the full, 10-yr duration of the survey. To enable this mode of operations, the INT has been refurbished and is fully robotic. While the baseline observing strategy is to observe each of the target stars once per night (subject to visibility constraints), this can be refined further; the optimal strategy depends on the detailed variability properties of the host stars and the strategies used to mitigate that variability Cegla (2019).

## 0.1 SUPERGRANULATION

The interplay of convection and magnetic fields in the surface layers of Sun-like stars gives rise to complex structures in their photospheres, which induce RV variability on a wide range of time-scales. Active regions (regions of enhanced surface magnetic flux, containing spots, faculae and plages), give rise to RV variations with amplitudes of several m s⁻¹ (Saar & Donahue 1997; Desort et al. 2007; Meunier, Desort & Lagrange 2010; Lovis et al. 2011; Gomes da Silva et al. 2012) on time-scales ranging from weeks (associated with the star’s rotation rate and the intrinsic evolution of the active regions) to years (associated with cyclic behaviour of the dynamo that powers the star’s large-scale magnetic field). Historically, this stellar activity signal has been the main factor limiting the detection of small-amplitude, long-period planets in RVs. However, the past 15 yr have seen major progress in mitigating activity signals in RV. At the simplest level, decorrelating the RVs against spectroscopic activity indicators (such as $\log R'_{HK} \& H\alpha$) can significantly reduce the activity signals (see Boisse et al. 2009; Dumusque et al. 2011b; Figueira 2013; Holzer et al. 2021, for examples), though such a simplistic approach is generally insufficient to capture their full complexity. More flexible models such as Gaussian Process (GP) regression (Aigrain & Foreman-Mackey 2023) have become increasingly popular to model activity signals in RVs (Haywood et al. 2014), and have proved particularly effective when applied to RVs and activity indicator time-series simultaneously (Aigrain, Pont & Zucker 2012; Rajpaul et al. 2015; Barragán et al. 2022; Delisle et al. 2022). In recent years, new methods have been proposed to model line-shape changes induced by active regions directly in time-series of cross-correlation functions (CCFs; Collier Cameron et al. 2021; de Beurs et al. 2022; Klein et al. 2022), and even spectra and/or line-by-line (LBL) RVs (Jones et al. 2017; Dumusque 2018; Rajpaul, Aigrain & Buchhave 2020; Cretignier, Dumusque & Pepe 2022; Lienhard et al. 2022; Cretignier et al. 2023; Zhao et al. 2024). In the best cases, these methods enable the correction of activity signals down to the sub-m/s level, allowing for the detection of low-mass exoplanets (e.g. Faria et al. 2022).

At the opposite end of the frequency spectrum, stellar oscillations, caused by trapped sound waves in the stellar atmosphere (Goldreich, Murray & Kumar 1994; Kjeldsen & Bedding 1995), induce RV variations on time-scales of a few minutes and with amplitudes of around 1 m s⁻¹ for Sun-like stars. These are typically mitigated by

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Dopplergram of the sun showing the supergranulation pattern. Due to the horizontal nature of supergranulation, the granules are more visible towards the limb of the Sun. (image credit: SOHO/MDI/ESA).</figcaption></figure>

using exposure times that are similar to or longer than the oscillation periods, so that the RV signal is averaged out (Chaplin et al. 2019).

On slightly longer time-scales, granulation is caused by convective upflows resulting in an overall convective blueshift (Roudier et al. 1991; Rimmele et al. 1995; Lefebvre et al. 2008). On the Sun, granulation flows can reach vertical velocities of several hundred m s⁻¹, and vary on time-scales of a few 10 s of minutes within individual granules, which themselves have a typical size of 1 Mm (Meunier 2021). However, the disc-integrated RV signal, which is averaged over 1000 000’s (for a Sun-like star) of granules, is observed to be around 0.4 m s⁻¹ for the Sun (Elsworth et al. 1994; Pallé et al. 1999; Sulis, Mary & Bigot 2020a).

While granulation mainly gives rise to vertical flows in the photosphere, supergranulation (SG), first discovered 60 yr ago (Hart 1954, 1956) is a predominantly horizontal flow. It was first discovered as a horizontal velocity field superimposed on the Sun’s mean solar rotation field near its equator, and confirmed as a disc-wide phenomenon by Leighton, Noyes & Simon (1962) using Doppler images. Supergranulation is a velocity perturbation below of the granulation field, with each supergranule in the Sun measuring 30–35 Mm across (November 1994; Rieutord et al. 2010; Rincon et al. 2017), which share a predominantly horizontal motion, with horizontal speeds of 200–400 m s⁻¹ compared to vertical speeds of 20–30 m s⁻¹ (Rincon & Rieutord 2018). The predominantly horizontal nature of the supergranulation flow is evident in SDO Dopplergram images after subtraction of the rotational motion of the Sun and the mean convective blueshift, where the supergranules disappear near the centre of the solar disk (see Fig. 1 for an example). Despite the slower speeds involved, the large spatial scale of supergranulation compared to granulation leads to larger disc-averaged RV variations. Various studies have estimated the amplitude of supergranulation signals in solar RVs to be in the range 0.5–1 m s⁻¹ using simulations (Meunier et al. 2015; Meunier & Lagrange 2019, 2020) and Sun-as-a-star observations (Al Moulla et al. 2023; Lakeland et al. 2024). This makes supergranulation the next largest stellar variability signal after magnetic activity in the Sun, and a very significant limiting
<table><caption>**Table 1.** List of variations of supergranule sizes and time-scale with magnetic field from the literature. Updated from Meunier et al. (2007). The first 6 references concern studies at a given time (spatial variations), while the last 6 reference studies cover the solar cycle or part of it. FT means Fourier Transform. A plus sign (+) means an increase in cell size or time-scale with increasing activity levels, while a minus sign (−) means the opposite. An equal sign (=) means no variation was found.</caption><thead><tr><th>Reference</th><th>Data</th><th>Method</th><th>Size</th><th>Time-scale</th><th>Comments</th></tr></thead><tbody><tr><td>Sýkora (1970)</td><td>Ca II images</td><td>autocorrelation</td><td>+</td><td></td><td></td></tr><tr><td>Wang (1988)</td><td>magnetograms</td><td>autocorrelation</td><td>+</td><td>—</td><td>Large errorbars</td></tr><tr><td>Wang et al. (1996)</td><td>magnetograms</td><td>autocorrelation</td><td>+</td><td></td><td></td></tr><tr><td>Hagenaar, Schrijver & Title (1997)</td><td>Ca II K images</td><td>segmentation</td><td>=</td><td></td><td></td></tr><tr><td>Raju & Singh (2002)</td><td>Ca II K images</td><td>autocorrelation</td><td>—</td><td></td><td></td></tr><tr><td>Meunier et al. (2007)</td><td>magnetograms</td><td>segmentation</td><td>—</td><td></td><td></td></tr><tr><td>Singh & Bappu (1981)</td><td>Ca II K images</td><td>autocorrelation</td><td>—</td><td></td><td>Via latitude variations</td></tr><tr><td>Muenzer et al. (1989)</td><td>Ca II K images</td><td>segmentation</td><td>+</td><td></td><td></td></tr><tr><td>Kariyappa & Sivaraman (1994)</td><td>Ca II K images</td><td>segmentation</td><td>—</td><td></td><td>Via brightness intensity</td></tr><tr><td>Komm, Howard & Harvey (1995)</td><td>magnetograms</td><td>autocorrelation</td><td>+</td><td></td><td>FWHM of the autocorrelation curves</td></tr><tr><td>Berrilli et al. (1999)</td><td>Ca II K images</td><td>segmentation</td><td>—</td><td></td><td>Over 1 yr only</td></tr><tr><td>De Rosa & Toomre (2004)</td><td>Doppler</td><td>segmentation</td><td>—</td><td>—</td><td>2 time series</td></tr><tr><td>Gizon & Duvall (2003, 2004)</td><td>Doppler</td><td>helioseismology</td><td>=</td><td>—</td><td>Week dependence</td></tr><tr><td>Meunier, Roudier & Rieutord (2008)</td><td>magnetograms</td><td>velocity field divergences</td><td>—</td><td></td><td>1σ detection</td></tr><tr><td>McIntosh et al. (2011)</td><td>magnetograms & Ca II K images</td><td>segmentation</td><td>+</td><td></td><td></td></tr><tr><td>Chatterjee, Mandal & Banerjee (2017)</td><td>Ca II K images</td><td>watershed</td><td>+</td><td></td><td></td></tr><tr><td>Mandal, Chatterjee & Banerjee (2017)</td><td>Ca II K images</td><td>watershed</td><td>+</td><td></td><td></td></tr><tr><td>Rajani et al. (2022)</td><td>Ca II K images</td><td>segmentation</td><td>—</td><td></td><td></td></tr><tr><td>Sowmya et al. (2023)</td><td>Ca II K images</td><td>segmentation</td><td></td><td>+</td><td></td></tr></tbody></table>

factor in the RV detection of Earth-twins (Meunier & Lagrange 2019, 2020). The characteristic lifetime of supergranulation has been widely debated. Time-scales ranging from 0.5 to 2 d have been reported in the literature (see Rincon & Rieutord 2018, for a review).

Convection at the solar surface is strongly coupled to the Sun's magnetic dynamics. For example, Leighton et al. (1962) and Simon & Leighton (1964) found strong correlations between the magnetic field distribution of the quiet-Sun and supergranulation flows. Simon et al. (1988) showed that supergranules are herded by small-scale magnetic structures at their boundaries. Given this link between magnetic field and supergranules, it is natural to ask how the physical properties of the latter vary over the Sun's magnetic cycle.

Several studies have addressed this question using solar observations; we have attempted to summarize their results in Table 1, which is an updated version of a similar table presented in Meunier, Roudier & Tkaczuk (2007). The vast majority of studies focused on the sizes of supergranules, while a small number looked at their lifetimes. All the studies were based on resolved images of the Sun. Despite the relatively similar types of data used, there are obvious disagreements between the results reported in the literature: some studies find supergranule sizes to be positively correlated with the Sun's activity level, others anticorrelated, and the same goes for their lifetimes. Meunier et al. (2007) proposed a possible explanation for these conflicting results: they found that while the size of supergranules is anticorrelated with the magnetic field strength within the supergranulation cells, larger supergranules have a stronger network at their boundary, which can lead to a negative or positive correlation being reported depending on how the magnetic activity level is defined. This highlights the care needed to interpret the results summarized in Table 1.

While spatially resolved studies provide valuable insights into the physical processes governing supergranulation, its impact on RV planet searches is best studied using disc-integrated, Sun-as-a-star RV observations. Lakeland et al. (2024) recently used images from the Solar and Heliospheric Observatory (SDO) to estimate the RV

variations due to active regions, and subtracted them from HARPS-N RVs to study the residual'quiet-Sun' RV variations. They used structure functions, a non-parametric method to quantify the root-mean-square (r.m.s.) variability of an irregularly sampled time-series as a function of time-scale. They found that the r.m.s. of the quiet-Sun RVs is approximately constant over the 8 yr of data they studied, at around 1 m s⁻¹, and that its characteristic time-scale is consistent with supergranulation simulations, indicating that supergranulation is the dominant phenomenon contributing to the residual variations.

The present study aims to revisit this question and expand on this work, with a few important differences. As we have seen when discussing the published results based on spatially resolved images, the way magnetic activity is defined can have a critical effect on our understanding of the connection between supergranulation and magnetism. We therefore consider two different methods to estimate and subtract the RV contribution of the active regions, one based on SDO images (following Lakeland et al. 2024) and one using pixel-level correlations between the spectrum time-series and activity indicators (Cretignier et al. 2021). When using the latter, we consider an additional year of data, up to and including the present activity cycle maximum. Finally, we model the supergranulation-like signal in the quiet-Sun RVs as a Gaussian Process (GP), following the method presented in O'Sullivan & Aigrain (2024), which allows us to measure its variance and characteristic time-scale directly, with robust uncertainties.

The remainder of this paper is structured as follows. Section 1 presents the observations and data reduction steps used in this study as well as the two methods used to extract quiet-Sun RVs, while Section 2 summarizes the GP regression method we use to measure the (super)granulation parameters. We present our results and compare them with the relative sunspot number in Section 3. We then investigate in Section 4 different observing strategies that could be used by EPRV surveys to characterize supergranulation signals in other stars than the Sun. Finally, we summarize our conclusions and discuss the implications of our results in Section 5.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** HARPS-N solar RVs (light blue, top), SDO quiet-Sun RVs (green, middle), and YARARA quiet-Sun RVs (pink, bottom), with arbitrary offset for graphical consideration. The vertical lines delimitate the end of cycle 24, the solar minimum, and the start of cycle 25.</figcaption></figure>

# 1 OBSERVATIONS

This study used observations from the HARPS-N solar telescope and the Solar Dynamics Observatory (SDO). We combine these data sets using different reduction methods to create two'quiet-Sun' RV data sets. The raw observations and the methods to extract the quiet-Sun data are described below.

## 1.1 Raw Data

### 1.1.1 HARPS-N spectra and RVs

The HARPS-N Solar telescope, situated on the Telescopio Nazionale Galileo (TNG) at the Roque de las Muchachos Observatory in La Palma (Spain), has been observing the Sun since 2015 at a 5 min cadence for (on average) 6 hr a day, allowing for a near continuous data set (Cosentino et al. 2012; Dumusque et al. 2015, 2021; Phillips et al. 2016; Collier Cameron et al. 2019). The telescope collects disc-integrated light from the Sun through a 3 arcmin lens that feeds into an integrating sphere, which scrambles the angular information, converting images to a point source. This in turn is fed through an optical fibre into the HARPS-N spectrograph, leading to a sun-as-a-star spectroscopic data set. These observations have very high resolution (R = 115 000), and cover a wide wavelength range, from 383 to 690 nm. The Sun is observed using 5-min exposures, resulting in a typical SNR of 400 per resolution element (at 550 nm). We use HARPS-N solar data taken during the period from 2015 July to 2023 November. As the 5-min cadence matches the characteristic frequency of the Sun's *p*-mode oscillations, these are mostly averaged out, while the granulation with a time-scale greater than five minutes, supergranulation and rotational modulation signals remain, along with instrumental systematics.

The solar spectra collected were processed, and RVs extracted, using version 3.0.1 of the ESPRESSO Data-Reduction Software

(DRS; Pepe et al. 2021) optimized for HARPS-N<fnref n="1" />, with a number of specific adjustments for the solar data as described in Dumusque et al. (2021). The DRS data products used in the present work were the order-merged spectra (S1D) and the RVs extracted by cross-correlation with a digitized mask.

Observations affected by clouds or strong differential extinction across the solar disc were excluded following the prescription described in Klein et al. (2024), which is based on the quality flags and differential extinction corrections derived by Collier Cameron et al. (2019). This selection resulted in a total of 107 602 solar spectra used in our analysis. On average, there are around 50 5-min observations a day, over a 5.3 hr long observing day. We show that this sampling does not affect our results in Appendix A. The HARPS-N solar RVs are shown in the top of Fig. 2. The corresponding General Lomb-Scargle Periodogram (Zechmeister & Kürster 2009) is shown at the top of Fig. 3.

### 1.1.2 SDO/HMI Dopplergrams and magnetograms

The Helioseismic and Magnetic Imager (HMI) onboard SDO has been collecting resolved images of the Sun since 2010, enabling detailed monitoring of the Sun's surface (Pesnell, Thompson & Chamberlin 2012). HMI takes observations in two polarization states, in each of 6 narrow bands around the magnetically sensitive 6176 Å Fe I line. A Gaussian function is fit to the intensities to estimate the local RV at each pixel and construct Dopplergrams, while the polarization information is used to determine the magnetic flux at each pixel and construct magnetograms.

In this work, we used 12 343 SDO/HMI observations, each with a 12-min exposure, taken every 4 hr between 2015 July 29th to 2022 November 12th.

<footnote n="1">¹The ESPRESSO DRS is publicly available on https://www.eso.org/sci/software/pipelines/espresso/?utm_source = chatgpt.com.</footnote>
<figure><img src="image_6.png" /><figcaption>Figure 3. Lomb–Scargle Periodogram of the HARPS-N Solar RVs (light blue), SDO quiet-Sun RVs (green), and YARARA quiet-Sun RVs (pink). The vertical lines indicate the solar rotation period and the first two harmonics (dashed).</figcaption></figure>

## 1.2 Extracting quiet-Sun RVs

This study focuses on the RV variations of the quiet photosphere. To isolate these, we first need to estimate the contribution of magnetically active regions, which otherwise dominate the disc-integrated RV variations, and subtract this contribution from the HARPS-N solar RVs.

We used two different methods to do this, one using SDO data and the other using HARPS-N data alone, allowing us to compare the results of the two approaches in the remainder of the paper.

### 1.2.1 SDO activity correction

We used SDO/HMI Dopplergrams and magnetograms to estimate and subtract the contributions of magnetically active regions to the HARPS-N disc-integrated solar RVs following the methodology developed by Haywood et al. (2016, 2022) and Milbourne et al. (2019). This is the same method used by Lakeland et al. (2024) to study the evolution of the quiet-Sun RV variations over multiyear time-scales, allowing us to perform a direct comparison to the results of that study.

This method involves first identifying which pixels on each HMI epoch belong to active regions by applying a magnetic flux threshold in the magnetogram, and a minimum area requirement (i.e. only active regions above a certain size are considered). We started off using the thresholds defined in Haywood et al. (2016) and Milbourne et al. (2019), which were also used by Lakeland et al. (2024). The magnetic threshold is defined by Haywood et al. (2016) as

$$|B_{r,\text{thresh},ij}| = 24G/\mu_{ij} \qquad (1)$$

where $\mu_{ij}$ is given by

$$\mu_{ij} = \cos \theta_{ij} \qquad (2)$$

where $\theta_{ij}$ is the angle between the line of sight and the outward normal to the feature on the solar surface and the area threshold used was 20 micro-hemispheres.

Once the active regions have been identified, the HMI Dopplergrams are used to estimate their contribution to the overall disc-integrated RVs. Owing to the presence of strong daily systematic

effects in the HMI Dopplergrams, and the very different wavelength range covered by HMI and HARPS-N, this cannot be done by simply summing over the relevant pixels in the Dopplergrams. Instead, the estimation is done using a physically motivated, 2-component model that is empirically calibrated on the HARPS-N RVs. The first component of the model represents the photometric signature of active regions, which results from dark sunspots and bright plage rotating in and out of view as the Sun spins, breaking the symmetry in the disc's rotational profile. The second component corresponds to the localized suppression of the convective blueshift in regions of enhanced magnetic flux (this is the dominant component in the Sun). This procedure was used in this work without modification, and we refer the interested reader to the relevant papers (Haywood et al. 2016, 2022; Milbourne et al. 2019) for a more detailed description.

Following Lakeland et al. (2024), we subtract the SDO-estimated active-Sun RVs from the HARPS-N solar RVs to produce the quiet-Sun SDO RVs used in this study. Linear interpolation was used to bridge the cadence gap between the SDO and HARPS-N observations (4 hr versus 5 min). Lakeland et al. (2024) showed that this does not affect the quiet-Sun results. The SDO quiet-Sun RVs and the corresponding periodogram are shown in Figs 2 and 3.

### 1.2.2 YARARA activity correction

Our second method uses YARARA (Cretignier et al. 2021), a post-processing methodology designed primarily to deliver improved RV precision compared to the DRS data products by correcting some instrumental systematics, as well as the effects of variable telluric absorption and stellar activity at the spectral level. Unlike the SDO activity correction, which makes use of resolved solar images, YARARA is fully data-driven, and relies mainly on the wavelength dependence of the effects being corrected (folding in prior knowledge such as the relevant reference frame.)

The starting point for YARARA post-processing is the time-series of 1D, order-merged spectra produced by the DRS (Dumusque et al. 2021), re-interpolated onto a common wavelength grid. Spectra are daily stacked to increase the SNR of the observations; this is required for a reliable detection and correction of the various effects at the spectrum level. Furthermore, we would not be able to
process more than ~2000 observations due to memory limitations (bearing in mind that each spectrum consists of ~10⁵ pixels). First, the continuum is modelled and the spectra are continuum-normalized using RASSINE (Cretignier et al. 2020b). We then applied corrections for cosmic rays, tellurics, an interference pattern specific to HARPS-N, and ghosts (see Cretignier et al. 2021 for details of each step). The correction of the instrumental point spread function (PSF) variations, first introduced in Stalport et al. (2023) and already applied to solar observations in Klein et al. (2024), was also included to mitigate subtle effects such as the change in PSF following the replacement of the HARPS-N cryostat change, and variations in the projected rotational velocity (v sin i) of the Sun as seen from the Earth (first mentioned in Collier Cameron et al. 2019). Finally, the YARARA activity correction is carried out by linear fitting and subtracting the dependence of the flux at each wavelength on the Ca II H&K S-index, as described in Cretignier et al. (2021). Since the S-index variations are dominated by the contributions of plage/faculae (Cretignier, Pietrow & Aigrain 2024), these are the main types of active regions we expect the YARARA activity correction to account for and contribution from spots and network are therefore still likely at the spectrum level. On the other hand, unlike the SDO correction, no explicit magnetic flux or area threshold is used in this correction. After all these corrections have been applied to the spectra, the RVs are extracted by using the same cross-correlation procedure, with the same G2 mask, as was used in the ESPRESSO DRS. Note that this mask was produced by using a line list tailored to the Sun, where the line centers of the CCF mask were obtained from a parabola fit on the core of the lines as in (Cretignier et al. 2020a, priv. communication). CCFs are computed on the colour corrected spectra following the guidelines of Barragán et al. (2024) in order to remove any airmass dependency or change in atmospheric conditions (Lovis 2007; Cretignier 2022).

Because YARARA is working with daily binned data, the methodology is blind to effects or systematics that act on a time-scale shorter than a day. In order to recreate the native 5-min sampling, we daily centred all the DRS RVs by the daily weighted mean before adding back the YARARA corrected daily binned RVs. In practice, such transformation is similar to a step-wise interpolation of the signals from the daily binned to the 5-min cadence. We slightly improved this extrapolation by assuming that the stellar activity signal should be smooth over a few days baseline. To do so, we interpolated linearly the YARARA activity model for the RVs and compute the difference with the step-wise interpolation. This extra correction vector was applied on the RVs. We refer to the resulting RVs as the YV1 RVs.

The YARARA quiet-Sun RVs and their corresponding periodogram are shown in Figs 2 and 3. The periodograms of the SDO and YARARA quiet-Sun RVs do not show significant peaks at stellar rotation period, indicating that the stellar activity has been effectively removed. We note that there is evidence of activity residuals at around 14 days period (Prot/2) and, to a lesser extent, at around 8 d (Prot/3) in the activity-corrected RVs, weather the activity-correction was done with YARARA or SDO. In Appendix B, we show that these residuals do not have an effect on the supergranulation signal that we model.

## 2 MODELLING THE QUIET-SUN RVS WITH GAUSSIAN PROCESSES

We follow the method described in O'Sullivan & Aigrain (2024) to model the quiet-Sun RVs to measure the standard deviation and characteristic time-scale of the remaining variability signals, and how these change over the Sun's activity cycle. We give a brief description

of the method below; interested readers are referred to O'Sullivan & Aigrain (2024) for a more detailed explanation.

Stellar RV variations are typically modelled as the sum of several stochastic processes, some of which are quasi-periodic (such as p-mode oscillations and rotational modulation of active regions), while others are aperiodic (such as granulation and supergranulation). In this work, we include only the latter, as the 5-min exposure times of the HARPS-N solar data mostly average out the p-modes, while the activity signals have been suppressed in the manner described in Section 1.2.

We use a two-component Gaussian Process (GP) model, where each component is an overdamped harmonic oscillator implemented in the CELERITE2 package (Foreman-Mackey et al. 2017; Foreman-Mackey 2018) as a Simple Harmonic Oscillator (SHO) kernel with a quality factor $Q = 1/\sqrt{2}$. This kernel has a power spectral density (PSD) given by

$$P_{1/\sqrt{2}}(\nu) = \frac{2 S_0}{1 + (\nu/v_0)^4}, \qquad (3)$$

where $v_0$ is the undamped frequency of the oscillator and $S_0$ controls its amplitude. The variance of the resulting process is

$$\sigma^2 = S_0\omega_0 Q = \frac{1}{\sqrt{2}} S_0\omega_0, \qquad (4)$$

where $\omega_0$ is the undamped angular frequency given by $2\pi v_0$. It is also convenient to define the undamped period

$$\tau = 2\pi/\omega_0. \qquad (5)$$

Here, $\tau$ corresponds to the time-scale at which the PSD drops. This is the same definition as used by Al Moulla et al. (2023), however, in the CELERITE2 package this is defined as $\rho$.

The PSD in equation (3) is similar to the 'Harvey law' (Harvey 1985) that is frequently used to model solar and stellar granulation and supergranulation signals (see e.g. Al Moulla et al. (2023)), but asymptotes to a steeper power law (index 4 rather than 2) at high frequencies. However, we note that there is ongoing debate in the literature as to the most appropriate power-law slope to use for convection-related signals (see e.g. Nordlund et al. 1997).

Our model also includes a white noise term added to the diagonal of the covariance matrix to account for the uncertainty due to photon noise and instrumental systematics. Overall, the GP kernel we use to model the quiet-Sun RVs is then

$$k_{\text{quiet}}(\delta_t) = k_{\text{gran}}(\delta_t) + k_{\text{SG}}(\delta_t) + k_n(t_i, t_j), \qquad (6)$$

where $k_{\text{gran}}$ and $k_{\text{SG}}$ are aperiodic SHO kernels:

$$k_{\text{SHO}, Q=1/\sqrt{2}}(\delta_t) = S_0 \omega_0 e^{-\frac{1}{\sqrt{2}}\omega_0 \delta_t} \cos\left(\frac{\pi}{\sqrt{2}}\right), \qquad (7)$$

with $\delta_t \equiv |t_i - t_j|$, representing the granulation and supergranulation components respectively, and $k_n$ is a white noise term:

$$k_n(t_i, t_j) = \delta_{ij}(\sigma_n^2 + \sigma_i^2), \qquad (8)$$

where $\sigma_n$ is the white noise standard deviation and $\sigma_i$ is the formal uncertainty of measurement $i$.

A key goal of this study is to investigate how the granulation and supergranulation characteristics change over the Sun's magnetic cycle. To this end, we split the quiet-Sun RVs into 4 week chunks and model each chunk independently. This also minimizes the impact of any residual activity signals in the quiet-Sun RVs on our results, as the duration of the chunks is similar to the Sun's rotation period. Any 4 week chunks containing fewer than 100 observations (indicative of
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Parameters of the GP fits to the SDO and YARARA quiet-Sun RVs obtained in 4-week chunks over the solar cycle. The parameters of the granulation component, log $S_G$ and log $\omega_G$, are shown in the left hand column, while the supergranulation parameters log $S_{SG}$ and log $\omega_{SG}$ are shown in the right hand column. The SDO and YARARA results are shown in green and pink, respectively. The blue horizontal line and shaded area represent the values and uncertainties reported for these parameters by Al Moulla et al. (2023). The blue vertical lines show the end of solar Cycle 24, the solar Minimum and the start of solar Cycle 25.</figcaption></figure>

significant data gaps) were excluded from the analysis. However, on average the chunks had around 1000 data points.

Following O'Sullivan & Aigrain (2024), we use a Markov Chain Monte Carlo (MCMC) to sample the joint posterior probability distribution over the 5 free parameters of our model, which are the $S_0$ and $\omega_0$ values for the granulation and supergranulation kernels, and the white noise standard deviation $\sigma_w$. We used version 3 of the EMCEE package (Foreman-Mackey et al. 2013, 2019) to perform the MCMC sampling. The walkers are initialized in a tight Gaussian ball (with standard deviation 0.01 dex) around a local optimum found using the MINIMIZE function in SCIPY's OPTIMIZE module. Ln-uniform priors are used for all the parameters, within the interval [-10; 10], in $d^{-1}$ for $\omega$ and $(ms^{-1})^2$ days for $S$. The standard deviation of the white noise $\ln \sigma_n$ is restricted to the interval [-4; 4] $ms^{-1}$. The number of walkers is set to 20 (4 times the number of parameters), and the MCMC chains are run for 50 000 steps. To differentiate between the granulation and supergranulation hyper-parameters, we set the first $\omega_0$ to be larger than the second, corresponding to the granulation signal.

To assess the convergence of the chains and to choose the appropriate burn-in and thinning factors we estimate the auto-correlation length of the chains using EMCEE's built-in functionality. We consider the autocorrelation time estimates to be reliable as long as the longest autocorrelation length estimate across all parameters is $\tau_{max} \le 200$ steps, i.e. less than a 50th of the chain. In this case we say the chains are well-converged. We discard the first 3$\tau_{max}$ samples as part of the

burn-in phase, and thin the remainder of the chains by a factor $\tau_{max}/4$. In the cases where the chains do not converge in 50 000 steps (due to large gaps in the data) we exclude the fit from future analysis.

## 3 RESULTS

The results of our fits for both the SDO and the YARARA quiet-Sun RVs are shown in Fig. 4, where we also compare them with the published values for these parameters reported by Al Moulla et al. (2023), who analysed HARPS and HARPS-N solar data between 2015–2018. For ease of interpretation, we converted the $S_0$ and $\omega_0$ values for each component to the standard deviation $\sigma$ and undamped period $\tau$ of the process, these results are shown in Fig. 5. Corner plots corresponding to the various stages of the solar activity cycle for both the SDO and YARARA results are given in Appendix C.

Three outliers are apparent in both the YARARA and SDO results, with parameter values that differ significantly from the general trend and large uncertainties. These correspond to 4-week chunks that contained at least one large (several day) data gap, despite containing > 100 points. These were discarded from the remainder of the analysis.

In Table 2, we summarize the average values for the whole solar cycle, the end of Cycle 24, the solar minimum, and the beginning of Cycle 25. For the average SDO $\tau_{SG}$ for the beginning of Cycle 25, we remove the anomalous result of a time scale of 8.91 d.
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Same as Fig. 4, but converted to standard deviation $\sigma$ and time-scale $\tau$. The blue horizontal line and shaded area represent the values and uncertainties reported for these parameters by Al Moulla et al. (2023).</figcaption></figure>

<table><caption>**Table 2.** Mean values for $\sigma_G$, $\tau_G$, $\sigma_{SG}$, and $\tau_{SG}$ for the whole solar cycle, the end of Cycle 24, the solar minimum, and the beginning of Cycle 25. The values found by Al Moulla et al. (2023) are also given. We note that the YARARA $\tau_{SG}$ values are systematically shorter than the corresponding SDO counterparts. This may be due to the YARARA post-processing pipelines that correct for a number of different instrumental and telluric effects.</caption><thead><tr><th></th><th></th><th>σ<sub>G</sub> [m/s]</th><th>τ<sub>G</sub> [days]</th><th>σ<sub>SG</sub> [m/s]</th><th>τ<sub>SG</sub> [days]</th></tr></thead><tbody><tr><td rowspan="2">Whole Cycle</td><td>YARARA</td><td>0.329 ± 0.003</td><td>0.031 ± 0.001</td><td>0.752 ± 0.020</td><td>1.053 ± 0.080</td></tr><tr><td>SDO</td><td>0.336 ± 0.004</td><td>0.037 ± 0.002</td><td>0.834 ± 0.020</td><td>1.458 ± 0.136</td></tr><tr><td rowspan="2">End of Cycle 24</td><td>YARARA</td><td>0.324 ± 0.006</td><td>0.034 ± 0.003</td><td>0.794 ± 0.036</td><td>1.020 ± 0.139</td></tr><tr><td>SDO</td><td>0.333 ± 0.006</td><td>0.044 ± 0.006</td><td>0.893 ± 0.033</td><td>1.245 ± 0.149</td></tr><tr><td rowspan="2">Solar Minimum</td><td>YARARA</td><td>0.333 ± 0.005</td><td>0.031 ± 0.001</td><td>0.727 ± 0.030</td><td>1.274 ± 0.125</td></tr><tr><td>SDO</td><td>0.332 ± 0.005</td><td>0.035 ± 0.002</td><td>0.783 ± 0.029</td><td>1.597 ± 0.163</td></tr><tr><td rowspan="2">Beginning of Cycle 25</td><td>YARARA</td><td>0.328 ± 0.006</td><td>0.028 ± 0.002</td><td>0.740 ± 0.037</td><td>0.614 ± 0.071</td></tr><tr><td>SDO</td><td>0.352 ± 0.012</td><td>0.033 ± 0.003</td><td>0.892 ± 0.034</td><td>0.796 ± 0.097</td></tr><tr><td>Al Moulla et al. (2023)</td><td></td><td>0.31</td><td>0.038 ± 0.008</td><td>0.68</td><td>0.554 ± 0.046</td></tr></tbody></table>

## 3.1 4-week chunk granulation and supergranulation results

The granulation parameters appear consistent over the Sun's activity cycle, as expected, see, for example, Sulis et al. (2020b). The values measured from the SDO and YARARA-corrected RVs are approximately consistent with each other. While the time-scale is consistent with the literature values, the amplitude $S_0$ is somewhat smaller. We note that the time-scale found for granulation, on average 44 minutes for the YARARA results, is longer than the lifetime of the granules on the stellar surface. This is because we are looking at the characteristic time-scale, corresponding with the knee of the granulation PSD, which in turn is attributed to the broad dispersion of granulation turnover time-scales (Seleznyov, Solanki & Krivova

2011). 44 min is consistent with the literature values found using Harvey models (e.g. Al Moulla et al. 2023).

On the other hand, the supergranulation parameters display significant variations over the solar cycle in both the YARARA and SDO results. This is promising, as both methods have limitations that may affect the supergranulation parameters. The SDO results may still include small active regions which fall below the area threshold of 20 micro-hemispheres defined in Milbourne et al. (2019). On the other hand, the YARARA activity correction might remove not only activity, but also some of the supergranulation signals. We recall that the YARARA activity correction consists in a simple linear decorrelation between the measured flux at each wavelength (in the
stellar rest-frame) and the Mount Wilson chromospheric activity index $S_{HK}$. If the supergranulation signal is correlated with the chromospheric index, some of it may be removed in the process. As the two methods are in agreement, and their limitations work in opposite directions, we therefore conclude the effects of the limitations seem minimal, and the cycle we see is the result of a physics effect.

To test for differences in the distributions of the parameters of the granulation and supergranulation models, we employ the use of Kolmogorov–Smirnov tests. We find that KS tests suggest that the supergranulation time-scales found during the End of cycle 24, the minimum, and the beginning of cycle 25 are drawn from different distributions compared to each other, with p-values of $10^{-11}$ or less in both the YARARA and SDO case. We also find that when comparing the parameters found via the two different methods directly that $\sigma_G$, $\tau_G$, and $\tau_{SG}$ are likely drawn from the same distribution. On the other hand, the KS test suggests that $\sigma_{SG}$ comes from a different distribution in the YARARA and SDO case (p-value = 0.0002). We propose that this difference comes from the different systematic attenuation and the differences in activity attenuation between the two methods.

## 3.2 Short-term scatter in the 4-week chunk results

Another important feature of our results, particularly for the supergranulation component, is the large amount of scatter in the measured parameter values shown in Figs 4 and 5, especially for the SDO results. Leaving aside any long-term trends, the time-scales in particular vary by up to an order of magnitude between consecutive 4-week data chunks, and this scatter is significantly in excess of the formal uncertainties. To check that this was not a fortuitous effect for the specific set of 4 week long chunks used, we shifted all the chunk boundaries by 14 d and reran the analysis, but the results were largely the same. We then repeated our analysis using longer 3-month (12-week) and 1-year chunks, and found that most of the scatter disappeared. The results are shown for the supergranulation component in Fig. 6. Based on our experience using injection-recovery tests on simulated data (O'Sullivan & Aigrain 2024), we believe that the scatter in the 4-week chunk results is real, and is not an artefact of the limited duration of the chunks. We therefore note that caution should be applied when interpreting supergranulation signals, even in observing campaigns lasting a month or more. The residuals in the SDO-corrected case have slightly larger scatter and contain more outliers than in the YARARA-corrected case. A possible explanation for this is that the YARARA post-processing corrects some HARPS-N instrumental systematics, which are not accounted for in the SDO-corrected data, in addition to the different ways in which the variations due to active regions are corrected.

## 3.3 Variation of the supergranulation time-scale over the solar cycle

The supergranulation time-scale $\tau_{SG}$ shows a clear variation over the decade spanned by our data. The variations are approximately in antiphase with the average relative Sunspot number in the corresponding time-chunks, which we downloaded from the Royal Observatory of Belgium's Sunspot Index and Long-term solar Observations (SILSO) portal<fnref n="2" /> (Clette & Lefèvre 2015), and which is also shown in Fig. 6.

<figure><img src="image_9.png" /><figcaption>**Figure 6.** Supergranulation parameters obtained from SDO and YARARA quiet-Sun RVs when fitting longer chunks of 12 weeks (stars) and 1 yr (circles). The colour scheme is the same as for Figs 4 and 5, and the horizontal line and shaded area once again shows the literature values from Al Moulla et al. (2023). Also shown for comparison are the average Sunspot numbers over the same time frame.</figcaption></figure>

To quantify this correlation, we plotted in Fig. 7 the supergranulation parameters measured from 3-month chunks as a function of average relative Sunspot number. We also computed the Spearman's rank correlation between each parameter and the Sunspot number, and report the corresponding *p*-values on Fig. 7. In Appendix D, we show the same for the 4 week and year long chunks. We observe strong anti-correlations (*p*-value below $10^{-2}$) between the time-scale from the SDO-corrected RVs. The correlation is even stronger (*p*-value below $10^{-4}$) for the time-scale derived from YARARA-corrected RVs. Both methods show very weak correlation for the standard deviation. We conclude that the time-scale of the

<footnote n="2">²See https://www sidc be/SILSO /home.</footnote>
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Correlation between the supergranulation parameters obtained from the 12-week (3-month) chunks and the Sun-spot numbers, for the SDO-corrected RVs (left, green) and the YARARA-corrected RVs (right, pink). In each panel, we also report the *p*-values derived from a Spearman's rank correlation analysis, which correspond to the probability of a similarly correlated dataset arising from white Gaussian noise only.</figcaption></figure>

supergranulation component displays a clear variation over the Sun's activity cycle, irrespective of the method used to evaluate it.

## 4 OBSERVING SUPERGRANULATION IN OTHER STARS

Supergranulation has become one of the main issues in EPRV planet searches, but we do not yet know how it scales as a function of stellar parameters. Given that supergranulation is a convective process, we would expect it to decrease in magnitude towards later spectral types, and the limited sample of stars with asteroseismic measurements spanning 5–8 nights analysed by Dumusque et al. (2011a) appears to confirm this. However, magnetic activity level also plays a role, as shown by the present study in the case of the Sun.

A number of observational efforts have been proposed to attempt to characterize supergranulation in stars other than the Sun, but it is very challenging as the signal of interest has characteristic time-scales similar to an Earth day. In this Section, we set out to evaluate the amount of observing time needed to adequately characterize this time-scale and how it depends on the choice of observing strategy.

Our goal is to estimate, for different observing strategies, the minimum duration of an observational campaign that would allow us to measure the supergranulation parameters of each target star to a given precision.

One possible way to do this would be to simulate large numbers of time-series with different parameters and time-sampling, then model them with a GP, sampling the posterior to evaluate the

uncertainties, as we did for the solar data in Section 2. However, this would be computational expensive. While we did this for a few simulations, we noted that it was hard to quantify the success of the simulation, and took significant time to test one strategy. We therefore adopt a Fisher information approach to the problem, following the work of Gupta & Bedell (2024). For data drawn from a multivariate Gaussian distribution with mean vector $\mu$ and covariance matrix $C$, the elements of the Fisher information matrix $F$ are given by:

$$F_{i,j} = \left(\frac{\partial \mu}{\partial \theta_i}\right)^T C^{-1} \left(\frac{\partial \mu}{\partial \theta_j}\right) + \frac{1}{2} tr \left(C^{-1} \frac{\partial C}{\partial \theta_i} C^{-1} \frac{\partial C}{\partial \theta_j}\right). \quad (9)$$

where $\theta$ is a vector of parameters controlling the mean and covariance functions which produce $\mu$ and $C$. The parameter uncertainties are then given by

$$\sigma_{\theta_i}^2 = F_{i,i}. \quad (10)$$

Here, we refer to the Fisher information matrix as $F$, as opposed to $B$, used by Gupta & Bedell (2024), in order to avoid confusion with the magnetic threshold.

In the context of RV surveys, the mean function encodes the planet signals, while the covariance function encodes the noise (photon, instrumental and stellar noise). The goal of Gupta & Bedell (2024) was to assess the impact of observing strategy on the characterization of planets. They therefore worked with a fixed covariance matrix $C$, and focused on the first term in equation (9). In this work, we are interested in the impact of observing strategy on the characterization
<figure><img src="image_11.png" /><figcaption>**Figure 8.** Expected $\sigma$ and $\tau$ fractional uncertainty as a function of number of nights of observation for continuous observations of a single star. The fractional uncertainty is calculated for various $\sigma$ and $\tau$ values corresponding to different phases of the solar cycle. The hyperparameters used are the average values shown in Table 2. The brown horizontal line corresponds to a fractional uncertainty of 0.3.</figcaption></figure>

of a stochastic signal. We therefore assume a trivial mean function, and the Fisher information matrix reduces to equation (9):

$$F_{i,j} = \frac{1}{2} \mathrm{tr} \left( \mathbf{C}^{-1} \frac{\partial \mathbf{C}}{\partial \theta_i} \mathbf{C}^{-1} \frac{\partial \mathbf{C}}{\partial \theta_j} \right), \qquad (11)$$

where the elements of the covariance matrix are set to the sum of two terms given by equation (7): one for granulation, and one for supergranulation (we implicitly assume that activity signals occur on a long time-scale compared to the simulated observations, or are modelled separately). The derivatives of the covariance matrix are given by:

$$\frac{\partial C}{\partial S} = \omega e^{-\frac{\omega \tau}{\sqrt{2}}} \cos \left( \frac{\omega \tau}{\sqrt{2}} - \frac{\pi}{4} \right), \qquad (12)$$

$$\frac{\partial C}{\partial \omega} = \frac{-S e^{-\frac{\omega \tau}{\sqrt{2}}}}{\sqrt{2}} \left( \omega \tau \sin \left( \frac{\omega \tau}{\sqrt{2}} - \frac{\pi}{4} \right) + (\omega \tau - \sqrt{2}) \cos \left( \frac{\omega \tau}{\sqrt{2}} - \frac{\pi}{4} \right) \right), \qquad (13)$$

We also include a 0.3 m s⁻¹ white noise term on the diagonal of the covariance matrix, representative of the expected precision for (e.g.) the HARPS3 instrument. One caveat of the Fisher Information approach is that it assumes that the real granulation and super-granulation signals are perfectly represented by the GP kernels. Of course, this may not be strictly true in practice, so the parameter uncertainties computed according to equation (11) should be treated as approximate.

We started by considering the simplest possible observing strategy, where one star was monitored continuously every night for a number of nights. Each observation consisted of a 5-min block. For simplicity we assume that every night lasts 8 hours, and ignore the impact of weather. Fig. 8 shows how the relative precision on the standard deviation and time-scale of supergranulation improves as the duration of the observing campaign increases, for different assumed values of the parameters (corresponding to the values measured for different parts of the solar cycle in Section 2). Fig. 8 demonstrates that at least a week or two of observations are needed to obtain any meaningful constraints on the supergranulation parameters. The time needed to reach a given precision depends on the time-scale of the signal, but only fairly weakly (a factor ~2 for the different curves shown,

which span an order of magnitude in time-scale). Finally, the relative precision curves fall exponentially, and the benefit of extending the observing campaign beyond one or two e-folding times becomes minimal. We note that while it takes the shortest amount of time to constrain $\sigma_{SG}$ at solar minimum, as we would expect, it takes the longest time to constrain $\tau_{SG}$ during this period. This is due to the increase in $\tau_{SG}$ that we observe at solar minimum.

To test that the Fisher information matrix uncertainties are accurate, we simulated and modelled 5 d of RV observations. The fractional uncertainty we obtained is similar to that calculated using the Fisher information matrix, indicating that the method is reliable.

Next, we considered a more-realistic scenario where the observer cycles randomly between a fixed number of N target stars, still using 5-min exposures, but varying the number M of consecutive exposures of each star taken before switching to the next star. This is intended to explore the trade-off between staying on a given star for longer, which provides better sampling of the granulation signal, versus spacing out the observations of each star in a given night, which improves the sampling of the supergranulation signal. Each time we switch to a new star, we include a pessimistic 5-min overhead for slewing and target acquisition. In actual fact, the overhead time is likely to be less. In this simulation, we also include an additional 0.1 m s⁻¹ white noise to account for the oscillations. Chaplin et al. (2019) showed that with a 5 min exposure time, this would be the noise that oscillations would induce in a Sun-like star. For cooler stars, this would be lower. Having simulated a set of time-stamps for each combination of N and M, we then used the Fisher information matrix to estimate how many nights would be needed to measure the supergranulation time-scale $\tau$ to a fractional uncertainty of 0.3. Since the individual time-stamps are random, we repeated the process 100 times for each observing strategy. As we are trying to see how many stars we can characterize supergranulation in a month, we have an upper limit of 50 d to the test.

A 0.3 fractional uncertainty limit was chosen arbitrarily as from Fig. 8 we see that reducing the fractional uncertainty beyond this limit would require significantly more survey time. Also, this is the level of detection we would need to be able to detect the supergranulation time-scale variance over the solar cycle in the sun. While this level of uncertainty will not tell us the supergranulation time-scales of other stars down to great precision, it will allow us to see if the time-scales
<figure><img src="image_12.png" /><figcaption>Figure 9. Number of nights needed to reach a fractional uncertainty of 0.3 in the supergranulation time-scale $\tau$ as a function of the number of stars observed each night. Different colours correspond to different observational strategies, as described in the text.</figcaption></figure>

varies as we would expect it to, which would be a great improvement on our current understanding of supergranulation.

The results are shown in Fig. 9, where different colours represent different values of $M$. The first take-away point from that figure is that for a given number of stars per night, the time taken to reach the desired precision is systematically lower for $M = 1$ (blue symbols on Fig. 9), i.e. when the target star changes every 10 min, despite the fact that 50 per cent of each night is lost to overheads in this scenario. For example, this rapid switching approach would allow the characterization of the supergranulation signal (albeit only at the modest precision of 30 per cent) to be completed in 23 nights for 10 stars. By contrast, staring at each star for an hour ($M = 11$, shown in purple on Fig. 9), this precision is only reached after 40 nights on a smaller target sample of $N = 6$ stars. The second take-away is that we can get roughly similar characterization results by observing 2 groups of 4 stars for ~10 nights each, compared to 8 stars observed across ~20 nights.

## 5 DISCUSSION AND CONCLUSIONS

We used GP regression to measure the standard deviation and characteristic time-scale of the granulation and supergranulation component of the HARPS-N solar RVs, and their evolution over the Sun's magnetic activity cycle, after removing the contribution of magnetically active regions using two different methods: one using resolved magnetograms and Dopplergrams from SDO, following Lakeland et al. (2024) and references therein, and one using the HARPS-N spectra themselves, following the YARARA methodology (Cretignier et al. 2021).

Our results for granulation were consistent across the two data sets used, and also consistent with literature values after accounting for the 5-min cadence of the observations. Furthermore, the granulation parameters appear stable over the full time-span of our data set, which matches theoretical and observational expectations (Muller et al. 2018; Sulis et al. 2020b).

On the other hand, we have shown that the supergranulation properties, particularly its time-scale, vary by a factor of 3 over the Sun's activity cycle, with peak to peak variations up to a magnitude. Specifically, we observe a longer time-scale at activity minimum compared to the rising and decaying phases of the cycle, and show that the supergranulation time-scale is strongly negatively correlated

with the average relative Sunspot number measured at a given time.

Interestingly, Lakeland et al. (2024), who studied the same SDO-corrected data set using structure functions rather than GPs, did not observe significant changes in the structure functions over the Sun's activity cycle. This may arise from the daily sampling of the data, which makes the characterization of signals on ~ daily time-scales, such as supergranulation, particularly challenging using a non-parametric method such as structure functions. By comparison, the GP regression approach used in the present work makes stronger assumptions about the nature of the signal being modelled. However, as discussed in O'Sullivan & Aigrain (2024), these assumptions are consistent with the common practice in the literature of using Harvey functions to model stochastic signals in stellar RVs.

To enable a direct comparison of our results to those of Lakeland et al. (2024), we analytically computed the structure functions corresponding to the 1 yr results for both the SDO and YARARA quiet-Sun RVs. For a stationary time series exactly described by a known covariance function, $K(\tau)$, as is the case for Gaussian Processes, the corresponding structure function is given by

$$SF(\tau) = 2(\sigma^2 - K(\tau)). \quad (14)$$

These structure functions are shown in Fig. 10. Following Lakeland et al. (2024), we quantify variability as $\sqrt{\frac{1}{2}SF}$.

At $\tau = 10^{-1}$ d, we see that there is a difference between the structure functions calculated at solar minimum and at more activity times; however, by $\tau = 10^0$ d, this correlation is no longer present. The same applies to the structure functions found in Lakeland et al. (2024), indicating that the results are consistent with each other.

Finally, we used a Fisher information approach to compute the number of nights needed to provide useful empirical constraints on supergranulation properties for stars other than the Sun. We showed that a dedicated survey spanning around 20 consecutive nights using an instrument with RV precision of order 30 cm s$^{-1}$, such as HARPS3, can characterize the supergranulation time-scale to a precision of around 30 per cent for up to 10 stars. This would lead to a significant improvement in our understanding of how supergranulation depends on stellar parameters such as temperature, surface gravity and metallicity, as well as activity level. Such a survey is currently being considered as part of the science verification phase of the Terra Hunting Experiment.
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Structure functions corresponding to the YARARA (middle) and SDO (right) hyperparameters calculated using year long chunks. The lines are coloured by the average sun spot number that year. The original results from Lakeland et al. (2024) are shown on the left.</figcaption></figure>

On the other hand, we also saw that the supergranulation time-scales measured for the Sun in 4-week chunks display considerable scatter over and above the formal uncertainties, which disappears when using 12-week chunks. The origin of this scatter needs to be investigated in more detail, as it may affect the interpretability of the results of such a 'supergranulation' survey.

As a sanity check we computed the residuals of the best-fitting GP model in each chunk, divided them by the combined (model + measurement) uncertainties, and checked that the results are consistent with a zero-mean, unit variance Gaussian distribution, which was always the case (irrespective of activity correction method or chunk duration). We also checked if the jitter term (white noise term added to the diagonal of the covariance matrix), which we fit alongside the supergranulation and granulation parameters, varies over the Sun's activity cycle, but it was constant within the uncertainties, with typical values in the range 0.20–0.25 m s⁻¹. The jitter term absorbs any variations not accounted for by the GP model, and this check shows that the GP model's ability to explain the data is quite good, and is not noticeably variable over the solar cycle. This increases our confidence in the fact that the change in the supergranulation time-scale we observe over the solar cycle is not an artefact of the used model (SHO kernel).

## 5.1 What causes the 'supergranulation cycle'?

We have brought to light a clear dependence of the supergranulation time-scale on the magnetic activity state of the Sun. What could be its origin, and does it match theoretical expectations?

As discussed in Section 0.1, there is considerable divergence in the literature on how supergranulation changes as a function of magnetic activity in the Sun. Around half of the papers listed in Table 1 find that the size of supergranulation cells is positively correlated with activity level, while the rest find a negative correlation. The literature on supergranulation cells' lifespan is much sparser, but De Rosa & Toomre (2004) and Sowmya et al. (2023) suggest that larger cells have a longer lifespan. Following this logic, the negative correlation we have evidenced between the supergranulation time-scale in the solar RVs and the Sunspot number would indicate that supergranulation cells tend to be larger at solar minima. In other words, our results would be consistent with the studies which have found a negative, rather than positive, correlation between cell size

and magnetic activity level. It is also interesting to note that Singh & Bappu (1981) and Meunier et al. (2008) suggest that the network magnetic elements have a limiting effect on supergranule sizes. If this is the case, we would expect larger cells at solar minimum, and thus longer time-scales, which is consistent with our results. However, the overall changes in supergranule size over the Sun's activity cycle found by e.g. Meunier et al. (2007, 2008) are too small (around 10 per cent) to explain the factor of two change in time-scale we report in this work. The cause of that time-scale change therefore remains unexplained for the time being.

## 5.2 Impact for RV surveys

The RV signature of supergranulation has a significant impact on the sensitivity of RV surveys to small, temperate planets (Meunier & Lagrange 2020). It introduces correlations between observations taken in a given night which, if not modelled explicitly, constitute a fundamental night-to-night noise floor. If the supergranulation time-scale is well understood, explicitly modelling the covariance of this signal could enhance the sensitivity to planet signals. We intend to quantify this using injection-recovery tests in a forthcoming paper.

However, an important consequence of the results presented in this work is that the RV supergranulation signal cannot be treated as stationary if the observational baselines are on the order of an activity cycle. It is therefore important to (a) continue investigating the correlation between supergranulation properties and activity level, (b) develop strategies to characterize supergranulation in other stars than the Sun, and (c) search for more direct supergranulation indicators in stellar spectra.

We also note that targeting stars which are in the low phase of their activity cycle, in order to minimize the impact of activity signals (as proposed for example by Sairam & Triaud 2022), may inadvertently result in a stronger supergranulation signal. Since this signal cannot be mitigated with current methods, that strategy may not be optimal.

## ACKNOWLEDGEMENTS

We thank the anonymos referee for their valuable comments that improved this manuscript. The authors thank Arvind F. Gupta for fruitful discussions. This publication is part of a project that
has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement no. 865624). NKOS thanks the LSST-DA Data Science Fellowship Programme, which is funded by LSST-DA, the Brinson Foundation, the WoodNext Foundation, and the Research Corporation for Science Advancement Foundation; her participation in the program has benefited this work. AM and BL acknowledge funding from a UKRI Future Leader Fellowship, grant number MR/X033244/1. AM acknowledges funding from an STFC small grant, reference: ST/Y002334/1. XD acknowledges the support from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement SCORE No 851555) and from the Swiss National Science Foundation under the grant SPECTRE (No 200021.215200). This work has been carried out within the framework of the NCCR PlanetS supported by the Swiss National Science Foundation under grants 51NF40_182901 and 51NF40_205606. SS acknowledges support from the 'Programme National de Physique Stellaire' (PNPS) and 'Programme National de Planétologie' (PNP) of CNRS/INSU co-funded by CEA and CNES. ACC acknowledges support from STFC consolidated grant number ST/V000861/1, and UKRI/ERC Synergy Grant EP/Z000181/1 (REVEAL). This publication makes use of The Data & Analysis Center for Exoplanets (DACE), which is a facility based at the University of Geneva (CH) dedicated to extrasolar planets data visualisation, exchange and analysis. DACE is a platform of the Swiss National Centre of Competence in Research (NCCR) PlanetS, federating the Swiss expertise in Exoplanet research. The DACE platform is available at https://dace.unige.ch. Based on observations made with the Italian Telescopio Nazionale Galileo (TNG) operated on the island of La Palma by the Fundación Galileo Galilei of the INAF (Istituto Nazionale di Astrofisica) at the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofísica de Canarias. This work made use of NUMPY (Harris et al. 2020), MATPLOTLIB (Hunter 2007), and PANDAS (pandas development team 2020) libraries. This work made use of ASTROPY<fnref n="3" />: a community-developed core PYTHON package and an ecosystem of tools and resources for astronomy (Astropy Collaboration 2013, 2018, 2022).

## DATA AVAILABILITY

This work makes use of the HARPS-N solar RVs, which are available at https://doi.org/10.82180/dace-h4s8lp7c and will be described in Dumusque et al. (in preparation). The SDO/HMI images are publicly available at https://sdo.gsfc.nasa.gov/data/

## REFERENCES

* Aigrain S., Foreman-Mackey D., 2023, A&A, 61, 329

* Aigrain S., Pont F., Zucker S., 2012, MNRAS, 419, 3147

* Al Moulla K., Dumusque X., Figueira P., Lo Curto G., Santos N. C., Wildi F., 2023, A&A, 669, A39

* Astropy Collaboration, 2013, A&A, 558, A33

* Astropy Collaboration, 2018, AJ, 156, 123

* Astropy Collaboration, 2022, ApJ, 935, 167

* Barragán O., Aigrain S., Rajpaul V. M., Zicher N., 2022, MNRAS, 509, 866

* Barragán O. et al., 2024, MNRAS, 531, 4275

* Berrilli F., Ermolli I., Florio A., Pietropaolo E., 1999, A&A, 344, 965

* de Beurs Z. L. et al., 2022, AJ, 164, 49

* Boisse I. et al., 2009, A&A, 495, 959

* Borucki W. J. et al., 2010, Science, 327, 977

* Cegla H. M., 2019, Geosciences, 9, 114

* Chaplin W. J., Cegla H. M., Watson C. A., Davies G. R., Ball W. H., 2019, AJ, 157, 163

* Chatterjee S., Mandal S., Banerjee D., 2017, ApJ, 841, 70

* Clette F., Lefèvre L., 2015, SILSO Sunspot Number V2.0, https://doi.org/10.24414/qnza-ac80

* Collier Cameron A. et al., 2019, MNRAS, 487, 1082

* Collier Cameron A. et al., 2021, MNRAS, 505, 1699

* Cosentino R. et al., 2012, in McLean I. S., Ramsay S. K., Takami H., eds, Proc. SPIE Conf. Ser. Vol. 8446, Ground-based and Airborne Instrumentation for Astronomy IV. SPIE, Bellingham, p. 84461V

* Crass J. et al., 2021, preprint (arXiv:2107.14291)

* Cretignier M., 2022, PhD thesis, University of Geneva, Switzerland

* Cretignier M., Dumusque X., Allart R., Pepe F., Lovis C., 2020a, A&A, 633, A76

* Cretignier M., Francfort J., Dumusque X., Allart R., Pepe F., 2020b, A&A, 640, A42

* Cretignier M., Dumusque X., Hara N. C., Pepe F., 2021, A&A, 653, A43

* Cretignier M., Dumusque X., Pepe F., 2022, A&A, 659, A68

* Cretignier M., Dumusque X., Aigrain S., Pepe F., 2023, A&A, 678, A2

* Cretignier M., Pietrow A. G. M., Aigrain S., 2024, MNRAS, 527, 2940

* Dalal S. et al., 2024, MNRAS, 534, 2410

* De Rosa M. L., Toomre J., 2004, ApJ, 616, 1242

* Delisle J. B., Unger N., Hara N. C., Ségransan D., 2022, A&A, 659, A182

* Desort M., Lagrange A. M., Galland F., Udry S., Mayor M., 2007, A&A, 473, 983

* Dumusque X., 2018, A&A, 620, A47

* Dumusque X., Udry S., Lovis C., Santos N. C., Monteiro M. J. P. F. G., 2011a, A&A, 525, A140

* Dumusque X. et al., 2011b, A&A, 535, A55

* Dumusque X. et al., 2015, ApJ, 814, L21

* Dumusque X. et al., 2021, A&A, 648, A103

* Elsworth Y., Howe R., Isaak G. R., McLeod C. P., Miller B. A., New R., Speake C. C., Wheeler S. J., 1994, MNRAS, 269, 529

* Faria J. P. et al., 2022, A&A, 658, A115

* Figueira P., 2013, in Chavez M., Bertone E., Vega O., De la Luz V., eds, ASP Conf. Ser. Vol. 472, New Quests in Stellar Astrophysics III: A Panchromatic View of Solar-Like Stars, With and Without Planets. Astron. Soc. Pac. San Francisco, p. 137

* Foreman-Mackey D., 2018, Res. Notes Am. Astron. Soc., 2, 31

* Foreman-Mackey D., Hogg D. W., Lang D., Goodman J., 2013, PASP, 125, 306

* Foreman-Mackey D., Agol E., Ambikasaran S., Angus R., 2017, AJ, 154, 220

* Foreman-Mackey D. et al., 2019, J. Open Source Softw., 4, 1864

* Gibson S. R., Howard A. W., Marcy G. W., Edelstein J., Wishnow E. H., Poppett C. L., 2016, in Evans C. J., Simard L., Takami H. eds, SPIE Conf. Ser. Vol. 9908, Ground-based and Airborne Instrumentation for Astronomy VI. SPIE, Bellingham, p. 990870

* Gizon L., Duvall T. L. J., 2003, in Sawaya-Lacoste H., ed., ESA Special Publication Vol. 517, GONG + 2002. Local and Global Helioseismology: the Present and Future. ESA, p. 43

* Gizon L., Duvall T. L., 2004, in Stepanov A. V., Benevolenskaya E. E., Kosovichев A. G., eds, Proc. IAU Symp. Vol. 223, Multi-Wavelength Investigations of Solar Activity. Cambridge Univ. Press, Cambridge, p. 41

* Goldreich P., Murray N., Kumar P., 1994, ApJ, 424, 466

* Gomes da Silva J., Santos N. C., Bonfils X., Delfosse X., Forveille T., Udry S., Dumusque X., Lovis C., 2012, A&A, 541, A9

* Gupta A. F., Bedell M., 2024, AJ, 168, 29

* Gupta A. F. et al., 2021, AJ, 161, 130

* Hagenaar H. J., Schrijver C. J., Title A. M., 1997, ApJ, 481, 988

* Hall R. D., Thompson S. J., Handley W., Queloz D., 2018, MNRAS, 479, 2968

* Harris C. R. et al., 2020, Nature, 585, 357

<footnote n="3">³http://www.astropy.org</footnote>
Hart A. B., 1954, MNRAS, 114, 17

Hart A. B., 1956, MNRAS, 116, 38

Harvey J., 1985, in Rolfe E., Battrick B. eds, ESA Special Publication Vol. 235, Future Missions in Solar, Heliospheric and Space Plasma Physics. ESA, p. 199

Haywood R. D. et al., 2014, MNRAS, 443, 2517

Haywood R. D. et al., 2016, MNRAS, 457, 3637

Haywood R. D. et al., 2022, ApJ, 935, 6

Holzer P. H., Cisewski-Kehe J., Zhao L., Ford E. B., Gilbertson C., Fischer D. A., 2021, AJ, 161, 272

Hunter J. D., 2007, Comput. Sci. Eng., 9, 90

Jones D. E., Pennington D. C., Ford E. B., Wolpert R. L., Loredo T. J., Gilbertson C., Dumusque X., 2017, preprint (arXiv:1711.01318)

Jurgenson C., Fischer D., McCracken T., Sawyer D., Szymkowiak A., Davis A., Muller G., Santoro F., 2016, in Evans C. J., Simard L., Takami H., eds, Proc. SPIE Conf. Ser. Vol. 9908, Ground-based and Airborne Instrumentation for Astronomy VI. SPIE, Bellingham, p. 99086T

Kariyappa R., Sivaraman K. R., 1994, Sol. Phys., 152, 139

Kjeldsen H., Bedding T. R., 1995, A&A, 293, 87

Klein B. et al., 2022, MNRAS, 512, 5067

Klein B. et al., 2024, MNRAS, 531, 4238

Komm R. W., Howard R. F., Harvey J. W., 1995, Sol. Phys., 158, 213

Lakeland B. S. et al., 2024, MNRAS, 527, 7681

Lefebvre S., García R. A., Jiménez-Reyes S. J., Turck-Chièze S., Mathur S., 2008, A&A, 490, 1143

Leighton R. B., Noyes R. W., Simon G. W., 1962, ApJ, 135, 474

Lienhard F., Mortier A., Buchhave L., Collier Cameron A., López-Morales M., Sozzetti A., Watson C. A., Cosentino R., 2022, MNRAS, 513, 5328

Lovis C., 2007, PhD thesis, University of Geneva, Switzerland

Lovis C. et al., 2011, preprint (arXiv:1107.5325)

Mandal S., Chatterjee S., Banerjee D., 2017, ApJ, 844, 24

Mayor M., Queloz D., 1995, Nature, 378, 355

Mayor M. et al., 2003, The Messenger, 114, 20

McIntosh S. W., Leamon R. J., Hock R. A., Rast M. P., Ulrich R. K., 2011, ApJ, 730, L3

Meunier N., 2021, in Bigot L., Bouvier J., Lebreton Y., Chiavassa A., Lèbre A., eds, Star-Planet Interactions, Evry Schatzman School 2019. EDP Sci. Proc., p. 22

Meunier N., 2024, Comptes Rendus Physique, 24, 140

Meunier N., Lagrange A. M., 2019, A&A, 625, L6

Meunier N., Lagrange A. M., 2020, A&A, 642, A157

Meunier N., Roudier T., Tkaczuk R., 2007, A&A, 466, 1123

Meunier N., Roudier T., Rieutord M., 2008, A&A, 488, 1109

Meunier N., Desort M., Lagrange A. M., 2010, A&A, 512, A39

Meunier N., Lagrange A. M., Borgniet S., Rieutord M., 2015, A&A, 583, A118

Milbourne T. W. et al., 2019, ApJ, 874, 107

Muenzer H., Schroeter E. H., Woehl H., Hanslmeier A., 1989, A&A, 213, 431

Muller R., Hanslmeier A., Utz D., Ichimoto K., 2018, A&A, 616, A87

Nari N. et al., 2025, A&A, 693, A297

Nordlund A., Spruit H. C., Ludwig H. G., Trampedach R., 1997, A&A, 328, 229

November L. J., 1994, Sol. Phys., 154, 1

O'Sullivan N. K., Aigrain S., 2024, MNRAS, 531, 4181

Pallé P. L., Roca Cortés T., Jiménez A., GOLF Team, Virgo Team, 1999, in Gimenez A., Guinan E. F., Montesinos B., eds, ASP Conf. Ser. Vol. 173, Stellar Structure: Theory and Test of Connective Energy Transport. Astron. Soc. Pac., San Francisco, p. 297

Pepe F. et al., 2021, A&A, 645, A96

Pesnell W. D., Thompson B. J., Chamberlin P. C., 2012, Sol. Phys., 275, 3

Phillips D. F. et al., 2016, in Navarro R., Burge J. H., eds, Proc. SPIE Conf. Ser. Vol. 9912, Advances in Optical and Mechanical Technologies for Telescopes and Instrumentation II. SPIE, Bellingham, p. 99126Z

Rajani G., Sowmya G. M., Paniveni U., Srikanth R., 2022, Res. Astron. Astrophys., 22, 045006

Rajpaul V., Aigrain S., Osborne M. A., Reece S., Roberts S., 2015, MNRAS, 452, 2269

Rajpaul V., Aigrain S., Roberts S., 2016, MNRAS, 456, L6

Rajpaul V. M., Aigrain S., Buchhave L. A., 2020, MNRAS, 492, 3960

Raju K. P., Singh J., 2002, Sol. Phys., 207, 11

Rauer H. et al., 2024, Exp. Astron., 59, 26

Ricker G. R. et al., 2014, in Oschmann Jacobus M. J., Clampin M., Fazio G. G., MacEwen H. A., eds, Proc. SPIE Conf. Ser. Vol. 9143, Space Telescopes and Instrumentation 2014: Optical, Infrared, and Millimeter Wave. SPIE, Bellingham, p. 914320

Rieutord M., Roudier T., Rincon F., Malherbe J. M., Meunier N., Berger T., Frank Z., 2010, A&A, 512, A4

Rimmele T. R., Goode P. R., Harold E., Stebbins R. T., 1995, ApJ, 444, L119

Rincon F., Rieutord M., 2018, Living Rev. Solar Phys., 15, 6

Rincon F., Roudier T., Schekochihin A. A., Rieutord M., 2017, A&A, 599, A69

Roudier T., Vigneau J., Espagnet O., Muller R., Mein P., Malherbe J. M., 1991, A&A, 248, 245

Saar S. H., Donahue R. A., 1997, ApJ, 485, 319

Sairam L., Triaud A. H. M. J., 2022, MNRAS, 514, 2259

Schwab C. et al., 2018, in Evans C. J., Simard L., Takami H., eds, Proc. SPIE Conf. Ser. Vol. 10702, Ground-based and Airborne Instrumentation for Astronomy VII. SPIE, Bellingham, p. 1070271

Seleznyy A. D., Solanki S. K., Krivova N. A., 2011, A&A, 532, A108

Simon G. W., Leighton R. B., 1964, ApJ, 140, 1120

Simon G. W., Title A. M., Topka K. P., Tarbell T. D., Shine R. A., Ferguson S. H., Zirin H., SOUP Team, 1988, ApJ, 327, 964

Singh J., Bappu M. K. V., 1981, Sol. Phys., 71, 161

Sowmya G. M., Rajani G., Paniveni U., Srikanth R., 2023, Astrophys. Bull., 78, 606

Stalport M. et al., 2023, A&A, 678, A90

Sulis S., Mary D., Bigot L., 2020a, A&A, 635, A146

Sulis S., Lendl M., Hofmeister S., Veronig A., Fossati L., Cubillos P., Van Grootel V., 2020b, A&A, 636, A70

Sýkora J., 1970, Sol. Phys., 13, 292

Thompson S. J. et al., 2016, in Evans C. J., Simard L., Takami H., eds, Proc. SPIE Conf. Ser. Vol. 9908, Ground-based and Airborne Instrumentation for Astronomy VI. SPIE, Bellingham, p. 99086F

Wang H., 1988, Sol. Phys., 117, 343

Wang H., Tang F., Zirin H., Wang J., 1996, Sol. Phys., 165, 223

Zechmeister M., Kürster M., 2009, A&A, 496, 577

Zhao Y. et al., 2024, A&A, 687, A281

pandas development team T., 2020, pandas-dev/pandas: Pandas, https://doi.org/10.5281/zenodo.3509134

## APPENDIX A: SIMULATIONS WITH REDUCED COVERAGE

To test whether the fact that our solar RV time series consist of sequences of 5.3 hr of continuous observations, compared to the day-long time-scale of supergranulation, has an effect on the parameters we derived we simulate a month-long data set with observations throughout the whole day. We simulated a granulation, supergranulation and white noise signals, and modelled them with the GP method. We subsample this data set down to 18, 12, and 6 h d⁻¹ and model each subsampled data set. The corner plots of these tests are shown in Fig. A1.

These results show that the accuracy of the parameters are not degraded with decreased daily sampling. We do note, however, that the precision degrades with decreasing sampling. There also seems to be some systematic increase in the deviations of the peak values of the parameter distribution compared to the true values, especially for the granulation parameters, however, these are always within the precision of the measurements.
<figure><img src="image_14.png" /><figcaption>**Figure A1.** MCMC posterior distribution plots for simulated data with a baseline of 28 d. Observations were simulated for 24 hr (black), 18 hr (blue), 12 hr (green), and 6 hr (cyan) in a 24 hr period. The 1D posterior distributions for each parameter, marginalized over all the other parameters, are shown by the histograms in the diagonal panels, with the true parameter values indicated by the blue lines. The 2D posteriors are shown in the off-diagonal panels, with 1, 2, and $3\sigma$ contours.</figcaption></figure>

## APPENDIX B: SIMULATIONS WITH ACTIVITY

Here, we show that activity residuals do not affect the convective signal that we observe. We simulated a 1-yr long time-series including activity residuals (using 1 yr's worth of HARPS-N timestamps and the parameters measured from the observed RVs for that yearly chunk). The activity residuals were simulated using a quasi-periodic celerite SHO term with $Q > 1$ at the rotation period. The activity signal included was at a period of 13.5 d, to correspond with the activity residuals observed in the periodogram seen in Fig. 3. We then modelled the data with a GP consisting of two SHO aperiodic

kernels and a white noise term (with a prior of [-10; 10] m s⁻¹. Fig. B1 shows the result. We see that the inclusion of the activity signal in our simulation does not change the accuracy of the granulation and supergranulation signals found, and the uncertainties are unchanged compared to the case where we do not inject activity residuals. We therefore do not expect these residuals to have a strong effect on the super-granulation properties measured using our GP models. Our explanation for this is that the GP model containing only SG and granulation terms has enough flexibility to absorb longer-term variations, but they do not dominate the variance or characteristic
<figure><img src="image_15.png" /><figcaption>**Figure B1.** MCMC posterior distribution plots for the granulation and supergranulation components. The data modelled were simulated with the addition of an activity signal with period = 13.5 hr. The 1D posterior distributions for each parameter, marginalized over all the other parameters, are shown by the histograms in the diagonal panels, with the true parameter values indicated by the blue lines. The 2D posteriors are shown in the off-diagonal panels, with 1, 2, and 3σ contours.</figcaption></figure>

time-scale estimates. We note that the white noise term is badly constrained, this is most likely the effect of the activity signal not being modelled.

corresponding to solar minimum shows diverging results compared to the end of cycle 24 and the start of cycle 25, as is expected.

## APPENDIX C: CORNER PLOTS

Here, we show the corner plots for three chunks representative of the solar cycle for both the SDO (Fig. C1) and YARARA (Fig. C2) quiet Sun RVs. We see that in both cases, the results for the three chunks are the same for granulation, while for supergranulation, the chunk

## APPENDIX D: CORRELATION PLOTS FOR 4 WEEK AND 1 YEAR LONG CHUNKS

Here, we present the correlation plots, like the one shown in Fig. 7, and the corresponding p-values, for the 4 week (Figure D1) and 1 yr (Figure D2) long chunks.
<figure><img src="image_16.png" /><figcaption>**Figure C1.** MCMC posterior distribution plots for SDO quiet sun results for 3 sample chunks. The results for the chunk 2016.47–2016.55 are in blue, 2018.77–2018.84 in green, and 2022.45–2022.52 in cyan. The 1D posterior distributions for each parameter, marginalized over all the other parameters, are shown by the histograms in the diagonal panels. The 2D posteriors are shown in the off-diagonal panels, with 1, 2, and 3σ contours.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure C2.** MCMC posterior distribution plots for YARARA quiet sun results for 3 sample chunks. The results for the chunk 2016.47–2016.55 are in blue, 2018.77–2018.84 in green, and 2022.45–2022.52 in cyan. The 1D posterior distributions for each parameter, marginalized over all the other parameters, are shown by the histograms in the diagonal panels. The 2-D posteriors are shown in the off-diagonal panels, with 1, 2, and 3σ contours.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>**Figure D1.** Correlation between the supergranulation parameters obtained from the 4-week chunks and the Sun-spot numbers, for the SDO-corrected RVs (left, green) and the YARARA-corrected RVs (right, pink). In each panel, we also report the *p*-values derived from a Spearman’s rank correlation analysis, which correspond to the probability of a similarly correlated data set arising from white Gaussian noise only.</figcaption></figure>
<figure><img src="image_19.png" /></figure>

**Figure D2.** Correlation between the supergranulation parameters obtained from the year long chunks and the Sun-spot numbers, for the SDO-corrected RVs (left, green) and the YARARA-corrected RVs (right, pink). In each panel we also report the *p*-values derived from a Spearman's rank correlation analysis, which correspond to the probability of a similarly correlated dataset arising from white Gaussian noise only.

This paper has been typeset from a TEX/LATEX file prepared by the author.