Open Access
Issue
A&A
Volume 685, May 2024
Article Number A173
Number of page(s) 14
Section Numerical methods and codes
DOI https://doi.org/10.1051/0004-6361/202347727
Published online 23 May 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1 Introduction

Transmission spectroscopy is a powerful tool for studying the atmospheres of transiting exoplanets. The method is based on the fact that the periodic dimming of a star during the planetary transit has a colour dependence (Seager & Sasselov 2000). If the atmosphere of a planet is more opaque at certain wavelengths due to absorption by atomic or molecular species, then it will result in a larger flux drop (i.e. a deeper transit). By studying these variations in the transit depth at different wavelengths, we can derive the atmospheric properties of the planet. The method has been used to detect atmospheric features such as Rayleigh scattering (Lecavelier Des Etangs et al. 2008; Gibson et al. 2017), clouds (Lendl et al. 2016; Alam et al. 2020), hazes (Pont et al. 2008; Nikolov et al. 2015; Sing et al. 2016), and different atomic and molecular species such as Na (Snellen et al. 2008; Wyttenbach et al. 2015), K (Borsa et al. 2021), Fe (Hoeijmakers et al. 2019; Ehrenreich et al. 2020), He (Allart et al. 2018), H2O, CO2 (JWST Transiting Exoplanet Community Early Release Science Team et al. 2023), and TiO/ VO (Evans et al. 2016) in the atmospheres of exoplanets. Transmission spectroscopy involves measuring small variations in the transit depth. Thus, it is often desirable to study planets with deeper transits as they produce larger amplitude signals. In the case of small Earth-like planets, we are limited to planets revolving around small stars as the transit depth scales with the square of the ratio between the planetary and stellar radii; but these small stars are often magnetically active (Reiners 2012; Shields et al. 2016). This activity manifests itself as surface inhomogeneities in the form of cooler (spots) and hotter (faculae) regions on the photosphere, and these impact the measurement of accurate transit depths and thus the transmission spectra.

The impact of activity features such as spots and faculae has been studied by multiple authors on different planetary systems. Their effect on the transit light curves can be broadly classified into two types: direct and indirect. The direct effect corresponds to the event in which an active region is located on the transit chord during the planetary transit. When the planet occults the active region, it results in anomalies in the transit light curves (Silva 2003; Sing et al. 2011). In the case of spots, it leads to a brief rise in flux that can lead to an underestimation of the transit depth, which can be as high as 10% at wavelengths of ~500 nm (Oshagh et al. 2014). In the case of faculae, the transit depth is overestimated. This effect can be mitigated by using routines that model the occulted active regions along with the transit light curve, such as PyTranSpot (Juvan et al. 2018), SOAP-T (Oshagh et al. 2013), or even Gaussian processes (Bruno et al. 2018).

The indirect effect of active regions on transmission spectra arises from the fact that the measurement of transit depth is differential. Transit depth is measured as the difference between the in-transit and out-of-transit flux of the star. As a transiting planet blocks light from the narrow transit chord of the star, the inho-mogeneities on the stellar surface cause the total disc-integrated flux (i.e. the in-transit or the out-of-transit flux) to differ from the flux blocked from below the transit chord, thus altering the transit depth measurements. In addition, due to spectral differences in the emerging flux, the effect is wavelength-dependent. This phenomenon is commonly referred to as stellar contamination of transmission spectra. In the case of cooler spots, possessing a redder spectrum than the disc-averaged stellar light, we observe a rise in transit depth at bluer wavelengths, which mimics the effect of Rayleigh scattering in the planetary atmospheres. For the highly active transiting system HD 189733, the predicted difference in transit depth between ultraviolet and infrared wavelengths is ~800 ppm (Pont et al. 2008, 2013; McCullough et al. 2014).

To accurately determine the extent of stellar contamination on transmission spectra, accurate estimates of the physical properties of active regions are needed, such as their temperature, size, and location on the stellar surface during the transit. One approach to constraining these properties is to photometrically monitor the star over at least a full rotation period around the transit observation, and use the amplitude and shape of the variability to retrieve the properties of the active regions present on the stellar surface (Pont et al. 2007, 2008). Rackham et al. (2017, 2018) present a forward-modelling approach to link the spot/facula coverage fractions to the amplitude of the variability in the Kepler light curves, and use them to predict the stellar contamination spectra. While allowing an order-of-magnitude evaluation of the effect, this technique is limited in its application to specific systems and observations by the degeneracy in the active region’s size and temperature, and its ignorance of the exact placement of the spots on the stellar disc. An alternative approach by Zellem et al. (2017) utilises the out-of-transit spectrum compared to the spectrum of the star during an inactive phase to calculate the stellar contamination without any prior measurements of the variability or the configuration of active regions. This technique hinges on the precision of the spectra and requires spectra of the quiet star, which are often difficult to secure.

The magnitude of stellar contamination on the transmission spectrum has an epoch dependence due to the evolution of active regions. This evolution, and the variations associated with stellar rotation, pose a significant problem when combining or comparing multi-epoch observations. As transmission spectroscopy becomes increasingly precise with the onset of the James Webb Space Telescope era, there is a need to carefully disentangle the effects of stellar contamination from any effects originating in the atmospheres of exoplanets.

In this work we introduce the Stellar Activity Grid for Exoplanets (SAGE), a tool for constraining and correcting for the epoch-dependent impact of active regions on transmission spectra. The code calculates the stellar contamination from a simulated disc-resolved stellar spectrum with known spot properties. This approach allows us to study additional effects such as limb-darkening (LD), which critically shapes the stellar contamination. It can also simulate the flux variability associated with the stellar rotation for a given spot configuration and can be used in a Markov chain Monte Carlo (MCMC) framework to retrieve spot properties from rotational variability measurements.

This paper is structured as follows: in Sect. 2 we present the astrophysical model that we used to model spots on the stellar photosphere. We also present how the variability and contamination models are calculated in SAGE. In Sect. 3 we discuss the effects of including LD and the projection effects of spots on the strength of stellar contamination. In Sect. 4 we use the variability model to draw inferences about the surface inhomogeneities of WASP-69. Sections 5 and 6 present the discussion, summary, and conclusions of our study.

2 Modelling inhomogeneities of stellar photospheres

Multiple routines, notably SOAP-T (Oshagh et al. 2013), Ksint (Montalto et al. 2014), PRISM (Tregloan-Reed et al. 2013, 2015), StarSim (Herrero et al. 2016), ellc (Maxted 2016), and PyTranSpot (Juvan et al. 2018), are available to model and correct for the ‘direct’ effects of stellar activity on the transit light curves (star spot crossings). With SAGE, we introduce a Python tool that corrects for the ‘indirect’ effects of stellar activity on transmission spectra (stellar contamination). To do so, the code is equipped with two capabilities: (i) it can fit the stellar rotational variability, constraining the distribution and size of active regions on the star, and (ii) based on the derived distribution, it can compute the stellar contamination at any given rotational phase of the star.

With SAGE, the user creates a model star projected on a two-dimensional pixel grid. In the model, the active regions on the stellar surface are parameterised using their number, position, size and temperature. Each pixel in this grid emits a spectrum whose properties are dependent on the local properties of the stellar surface, for example the pixels associated with cooler regions emit a cooler spectrum compared to the pixels associated with hotter regions. The emitted spectra are parameterised using the effective temperature, surface gravity and metallicity of the stellar photosphere. Additionally, stellar effects such as LD, rotational broadening and inclination can be modelled. By rotating this model star and integrating the emitted flux over a desired wavelength band, we obtain the flux variability for a given configuration of active regions.

2.1 Model setup

SAGE uses a pixelation approach similar to that used in PyTranSpot (Juvan et al. 2018) and PRISM (Tregloan-Reed et al. 2013, 2015). We projected the three-dimensional stellar sphere onto a two-dimensional Cartesian grid of pixels each with the same projected area. This method allows multiple components to be modelled simultaneously, such as spots, faculae, and LD. The size of the pixel grid (n × n) is defined through the star-planet radius ratio together with a user-specified planet size in pixels. A higher planet pixel size provides a better spatial resolution of the inhomogeneities, but at the cost of computation time. A planet radius of 15–50 pixels is recommended when modelling occulted spots in transit light curves and we recommend using these limits for SAGE as well (Juvan et al. 2018; Tregloan-Reed et al. 2015). The choice is also dependent on the resolution of the input spectra; a finer grid will provide better estimates of the effects of LD and rotational broadening relevant at high resolution. However, an increased resolution of spectra and grid comes with a higher required computation time.

We modelled surface inhomogeneities as homogeneous circular regions on the stellar disc. These circular regions deform to ellipses as regions move closer to the limb of the star. This deformation inherently corrects for projection effects that are especially relevant for active regions located near the limb of the star (for example, polar spots). From solar observations, we know that spots and faculae have complicated structures and often appear in groups (Solanki 2003; Pietrow et al. 2020; Mandal et al. 2021). However, we are currently limited in resolving such structures in other stars given the limited information available from spot occultations during transit.

We chose a spherical coordinate system with its origin located at the centre of the stellar sphere. The inhomogeneities on the stellar photosphere are parameterised using their position, size, and temperature or contrast with respect to the immaculate photosphere. These are modelled in SAGE using the following parameters:

  • The longitude (θ) varies between −180° and 180°. The visible stellar disc ranges from −90° to 90° in the counter-clockwise direction.

  • The latitude (ϕ) varies between 90° and −90°. The north pole of the star is at 90° and the south pole is at −90°.

  • The size (α) can be defined either as the angular size on the three-dimensional stellar sphere or in terms of the absolute filling factor. The absolute filling factor is defined as the ratio between the surface area of the active region and the total stellar surface.

  • The contrast (ρspot) is the ratio between the intensities of an active region to the immaculate photosphere, disregarding LD. It is a wavelength-dependent quantity that ranges from 0 to 1 for cooler regions and is >1 for hotter regions.

thumbnail Fig. 1

A model stellar grid and spectral setup in SAGE. Left: Stellar sphere with LD and spots projected onto a two-dimensional Cartesian pixel grid. Right: normalised model spectra used, obtained from the PHOENIX library. For the clear photosphere and for spots, we used temperatures of 4500 and 3500 K, respectively.

2.2 Stellar spectrum

A model spectrum is allocated to each pixel on the stellar grid. Various spectral libraries can be chosen, and we tested both low-to mid-resolution spectral libraries, such as ATLAS9 (Castelli & Kurucz 2003) and PHOENIX (Husser et al. 2013), and high-resolution spectral libraries, such as Coelho (Coelho 2014). To account for the difference in the emergent spectrum between the active regions and clear photosphere, different model spectra are chosen for both clear (FClear,λ) and active regions (FActive,λ). The model spectra are normalised to the peak of the clear spectrum before being embedded in the pixel grid (Fig. 1). This does not alter the shape of the input spectrum but makes the computation more tractable. The calculation of the disc-integrated spectrum at high resolution requires the addition of essential stellar effects such as convective blueshift, which will be part of future work.

The normalised disc-integrated spectrum of the star is calculated by summing the individual flux contributions from each pixel, Fluxλ=iNactiveFactive,λ,i+jNNactiveFclear,λ,jN.${\rm{Flu}}{{\rm{x}}_\lambda } = {{\mathop \sum \nolimits_i^{{N_{{\rm{active}}}}} {F_{{\rm{active}},\lambda }}_{,i} + \mathop \sum \nolimits_j^{N - {N_{{\rm{active}}}}} {F_{{\rm{clear}},\lambda }}_{,j}} \over N}.$(1)

Here, Factive,i is the flux emission from the ith pixel in the active regions and Fclear,j is the flux from the jth pixel of the clear photospheres. N and Nactive are the total numbers of pixels on the stellar disc and the active regions, respectively. The wavelength coverage and resolution of the disc-integrated spectrum are the same as the input spectra. If the spectral resolution of the clear and active photosphere models do not match, then the spectrum with the higher resolution is binned down to the lower-resolution spectrum.

2.3 Rotational broadening and limb-darkening effects

One of the main advantages of using a pixelation approach is the possibility to incorporate effects such as rotational broadening, LD and gravity-darkening in the stellar grid. The magnitude of these effects varies based on the properties of the star. The current version of SAGE is set up to include both, LD and rotational-broadening effects.

In our model, we assume the star is a rigid rotator to calculate the effect of rotational broadening. The model spectrum for each pixel is shifted in wavelength according to its radial velocity, which depends on the pixel location and the rotational velocity of the star. The injected model spectra should be free from the effects of stellar rotation. The rotation of the star in our model is parameterised using the equatorial velocity of the star in km s−1. The resolution of the input spectra puts a natural limit on the velocity resolution, and thus the effects of rotational broadening can only be studied with sufficiently high-resolution input spectra.

To model LD, we chose the quadratic LD law: I(μ)I(1)=1u1(1μ)u2(1μ)2,${{I(\mu )} \over {I(1)}} = 1 - {u_1}(1 - \mu ) - {u_2}{(1 - \mu )^2}{\rm{,}}$(2)

where u1 and u2 are the wavelength-dependent LD coefficients, I(1) is the intensity at the centre of the stellar disc and μ = cos(θ), where θ is the angle between the line of sight and normal to the stellar surface.

In SAGE, the required LD coefficients, u1 and u2, are calculated using the Limb-Darkening Coefficients and Uncertainties (LDCU) package1. LDCU is a modified version of the Python routine by Espinoza & Jordán (2015) that computes the LD coefficients and their corresponding uncertainties using a set of stellar intensity profiles accounting for the uncertainties on the stellar parameters. The stellar intensity profiles are generated based on two libraries of synthetic stellar spectra: ATLAS (Kurucz 1979) and PHOENIX (Husser et al. 2013). The LD coefficients in LDCU are parameterised by the effective temperature (Teff), surface gravity (log 𝑔), metallicity ([M/H]), and turbulent velocity (vtur) of the star.

Limb-darkening can be different for quiet and active regions. From solar observations, Crane et al. (2004) measured a much stronger LD for a spotted region compared to quiet regions at mid-ultraviolet wavelengths (175–210 nm). For faculae on the Sun, surface brightness increases towards the limb (Fontenla et al. 1999; Chapman & McGuire 1977; Meunier et al. 2010). Unfortunately, the LD behaviour of active regions has not been studied for other stars. We therefore chose to use the same LD coefficients for both quiet and active regions, but we note this may not be fully accurate. We therefore included the option to use different sets of LD coefficients for active and quiet regions, should the user prefer to do so.

2.4 Stellar variability model

To compute the stellar flux variability for any specific configuration of active regions, we rotated the stellar grid model at the rotational period of the star. In addition, the user can set the stellar inclination (i), that is, the angle between the stellar rotation axis and the line of sight. We assume the star is a rigid rotator with non-evolving, circular active regions. The disc-integrated spectra are calculated at each rotational phase of the star and light curves are calculated by summing the emergent flux over the desired bandpass. The amplitude of the flux variations depends on the properties of the active regions and the bandpass of the instrument. For example, many Sun-like spots (of angular size ~2°) uniformly distributed across the photosphere will produce little variability. On the other hand, large M-dwarf-type spots (of angular size of ~7–8°) covering the same total area as in the above case, will produce a larger amplitude, as they are less evenly distributed.

Figure 2 show the rotational variability model for a randomly chosen spot configuration on a K5V-type, equator-on star (i = 90°). The light curve was calculated assuming a uniform throughput over a wavelength range of 3000–10 000 Å.

One of the main limitations when inferring the properties of the stellar surface using the flux variability model is the degeneracy between the different properties of the active regions. The location, number, temperature, and size of the active regions are strongly correlated with each other. This can be remedied to some extent with constraints on the properties of the active regions obtained from other independent techniques. The occultation of an active region by a transiting planet results in the characteristic ‘bumps’ in the transit light curve, which has been proven to be useful in measuring the temperature and size of the active region (e.g. Sing et al. 2011). While adding important information, spot crossings have only been reported for a handful of active stars (Silva-Valio et al. 2010; Désert et al. 2011; Močnik et al. 2016; Morris et al. 2017). Another technique involves simultaneous photometric observations of the star in different filters. The amplitude of the flux variability is wavelength-dependent and by comparing the amplitude in different filters, the temperature of the active regions can be constrained (Henry et al. 1995; Froebrich et al. 2020).

2.5 Stellar contamination of transmission spectra

At any given time, the effect of unocculted active regions on the observed transit depth at different wavelengths can be described as follows: (RpR)λ, observed 2=(RpR)λ, true 21fproj (1Factive ,λFclear ,λ).$\left( {{{{R_p}} \over {{R_ \star }}}} \right)_{\lambda ,{\rm{ observed }}}^2 = {{\left( {{{{R_p}} \over {{R_ \star }}}} \right)_{\lambda ,{\rm{ true }}}^2} \over {1 - {f_{{\rm{proj }}}}\left( {1 - {{{F_{{\rm{active }},\lambda }}} \over {{F_{{\rm{clear }},\lambda }}}}} \right)}}.$(3)

Here, Rp is the radius of the planet, R is the radius of the star, fproj is the total projected surface area of inhomogeneities and Factive,λ and Fclear,λ are the spectrum of the active and clear photosphere, respectively (see also Eq. (7) of Berta et al. 2011, Eq. (1) of Zellem et al. 2017 and Eq. (1) of Rackham et al. 2018). Thus, the measured and true transit depth are related by a wavelength-dependent factor (ϵλ), ϵλ,t=11fproj (t)(1Factive ,λFclearλ).$\left( {{{{R_p}} \over {{R_ \star }}}} \right)_{\lambda ,{\rm{ observed }}}^2 = {{\left( {{{{R_p}} \over {{R_ \star }}}} \right)_{\lambda ,{\rm{ true }}}^2} \over {1 - {f_{{\rm{proj }}}}\left( {1 - {{{F_{{\rm{active }},\lambda }}} \over {{F_{{\rm{clear }},\lambda }}}}} \right)}}.$(4)

The ϵλ,t is often referred to as the stellar contamination factor (Rackham et al. 2018), the stellar activity correction factor (ζ) (Zellem et al. 2017), or the contamination spectrum. The relative change in observed transit depth at any specific wavelength is given by δdF=ϵλ,t1.${\delta _{{\rm{dF}}}} = {_{\lambda ,t}} - 1.$(5)

In SAGE, we calculated ϵλ,t as the ratio between the disc-integrated spectrum of the clear and the spotted star at the time of interest: ϵλ,t= clear stellar spectrum spotted stellar spectrum        =jNFclear ,λ,jiNactiveFactive ,λ,i+jNNactiveFclear ,λ,j.$\matrix{ {{_{\lambda ,t}} = {{{\rm{ clear stellar spectrum}}} \over {{\rm{ spotted stellar spectrum}}}}} \hfill \cr {\,\,\,\,\,\,\,\, = {{\sum\nolimits_j^N {{F_{{\rm{clear }},\lambda ,j}}} } \over {\sum\nolimits_i^{{N_{{\rm{active}}}}} {{F_{{\rm{active }},\lambda ,i}}} + \sum\nolimits_j^{N - {N_{{\rm{active}}}}} {{F_{{\rm{clear }},\lambda }}} ,j}}.} \hfill \cr } $(6)

This method allows us to accurately include projection effects of active regions and any additional atmospheric effects, such as LD. When the above secondary effects are neglected, Eq. (6) transforms into Eq. (4) (see Appendix A).

3 The effects of starspots on transmission spectra

In the previous section, we presented SAGE and the different stellar effects that can be modelled via the pixelation approach. Here, we investigate the impact of projection and LD effects on the extent of stellar contamination in the transmission spectrum. We here only consider dark starspots as the LD behaviour of faculae is little understood at this time. We studied stellar contamination for F5V-, G5V-, K5V-, and M5V-type stars using model input spectra from the PHOENIX spectral database (Husser et al. 2013). Where necessary, spectra at specific temperatures are calculated by linear interpolation of neighbouring grid points. Table 1 contains the adopted stellar parameters for both clear and spotted regions for different spectral types. The temperatures of the clear photospheres are obtained from Pecaut & Mamajek (2013) and the spot temperatures are calculated from Rackham et al. (2018, 2019), which uses the relation Tspot = 0.418 × Tclear + 1620 K for F-, G-, and K-type stars and Tspot = 0.86 × Tclear for M dwarfs. The surface gravity (log 𝑔) is set to 4.5 for F, G, and K stars and to 5.5 for the M dwarf. For all stars, we assumed solar metallicity with no alpha enrichment.

thumbnail Fig. 2

Flux variability and stellar contamination models for different rotation phases of a synthetic star. Top: Stellar surface grid at different rotational phases for an (equator-on) K5V star with ten 2–3° spots randomly distributed across its surface. The effective temperature of the spots and the clear photosphere is 3500 and 4500 K, respectively. Middle: disc-integrated flux calculated between 3000 and 10 000 Å as a function of the stellar rotational phase. Bottom: wavelength-dependent stellar contamination factor (ϵλ) at rotation phases corresponding to those shown in the top panel.

3.1 Position and size of active regions

The extent of stellar contamination on the transmission spectra depends mainly on two properties of the active regions: their projected surface area (a position and size-dependent property) and the surface brightness contrast with respect to the surrounding photosphere. The projected surface area of an active region changes not only with stellar rotation but also due to their evolution over time. We note that stellar rotation alone has no impact on the absolute filling factor of the active regions2.

We studied the dependence of stellar contamination on the projected surface area of spots through simulations with SAGE, assuming spot sizes between 0 and 7° and locations ranging from the centre of the stellar disc to the limb of the star. For this, we assumed a uniformly illuminated stellar disc with no LD and rotational broadening. The simulations were performed for all four spectral types tabulated in Table 1.

As is visible in Fig. 2 (bottom panel), the stellar contamination rises towards bluer wavelengths. To quantify the strength of stellar contamination for different cases, we chose to use the change in observed transit depth at different wavelength bands, hereafter, the contamination offset. We defined the contamination offset (δ) as the difference between the average observed transit depth in between wavelength bands 3000–4000 Å (ϵ0.3–0.4) and 9000–10 000 Å (ϵ0.9–1.0): δ=(RpR)true 2×(ϵ0.30.4ϵ0.91.0).$\delta = \left( {{{{R_p}} \over {{R_ \star }}}} \right)_{{\rm{true }}}^2 \times \left( {{_{0.3 - 0.4}} - {_{0.9 - 1.0}}} \right).$(7)

Here, (Rp/R)true2$\left( {{R_p}/{R_ \star }} \right)_{{\rm{true}}}^2$ is the true transit depth.

In Fig. 3, we present the computed contamination offsets for different stars transited by a hypothetical, bare-rock, planet, that is to say, a planet with no atmosphere and a uniform transit depth of 1% in all wavelengths. In the top panels, we illustrate the offsets calculated for spots located at different locations on the stellar disc. The dependence of stellar contamination on the position of spots is clearly visible. In all simulations, the effect is maximal when the spot is located at the centre of the stellar disc, but decreases significantly when the spot is moved closer to the limb. For the M5V-type star, the maximum offset exceeds 30 ppm. However, we do not measure such high offsets for the earlier spectral types, even though they have cooler spots. Generally, we find a clear anti-correlation between the contamination offset and the effective temperature of the star with the cooler stars having a higher contamination offset than the hotter stars.

The bottom panels of Fig. 3 illustrate the dependence of the contamination offset on the spot size. As expected, contamination offsets are largest for large spots located close to the centre of the stellar disc. Again, M5 V stars show the highest contamination, with offsets exceeding 60 ppm for spots with an angular size of 7°. Such large spots have been detected in M dwarfs through occultations by a transiting planet (Almenara et al. 2022). These simulations isolate the dependence of the variations of the stellar contamination with the projected surface of the active regions, as we assumed a non-limb darkened star emitting equally at all distances from the centre.

thumbnail Fig. 3

Variation in contamination offset with location and size of active regions for stars of different spectral types. Top: Stellar contamination offset for a single spot of 5° angular size, placed at different locations on the stellar disc for stars of different spectral types. Bottom: Stellar contamination offset for the same stars with spot sizes ranging from 1 to 7° for different longitudes. LD was not included in these simulations.

Table 1

Model spectrum parameters used in the presented simulations.

3.2 Limb-darkening

A more physically accurate simulation includes the effect of LD (i.e. the drop in intensity between the centre of the stellar disc and its limb). We thus recreated our simulations incorporating LD in the pixelated stellar grid. To do so, we adjusted the intensity of every pixel on the stellar grid following a quadratic LD law. The coefficients are calculated at each wavelength using LDCU with the stellar parameters Teff and log(𝑔) as tabulated in Table 1.

In Fig. 4, we illustrate the effect of LD on the stellar contamination spectra using as an example a 5° spot at longitudes of 0°, 30°, 50°, and 70° for the four spectral types considered. The contamination spectra show a diversity of slopes and features for different spectral types and spot positions. For all spectral types, the LD model predicts higher contamination than the non-LD case when the spot is located closer to the centre of the stellar disc (μ = 1 and μ = 0.866). Indeed, the more realistic LD simulations show a ~1.8 × larger contamination offset (calculated using Eq. (7)), compared to the simplified non-LD case. On the other hand, in the case of spots being located closer to the limb of the star, the non-LD model overestimates the contamination effect. Here, the offset is up to ~0.5 × reduced compared to the non-LD model. This behaviour is easy to understand: due to LD, the bright disc centre contributes more strongly to the total flux, enhancing the impact of spots located there. Conversely, the impact of spots near the darker limb is reduced.

Limb-darkening also has an impact on the shape of the contamination spectrum due to its wavelength dependence. For spots close to the disc centre, the characteristic rise in the contamination spectrum at lower wavelengths is visible in both the LD model and the non-LD model. However, for spots located very close to the limb of the star, LD considerably flattens the contamination spectra, most notably in the wavelength regime of 4500–7500 Å. In addition, at wavelengths <4500 Å, spots at the stellar limb produce a positive slope in the spectra for the F5V-, G5V-, and K5V-type stars.

We could again compute the contamination offset δ (similar to Sect. 3.1) and study its variation with the position and size of the active regions. The results are summarised in Fig. 5. For spots located closer to the centre of the stellar disc, the amplitude of the stellar contamination offset increases by roughly 30 ppm with the addition of LD. For the M-type star, spots with an angular size of 7° can introduce transit depth offsets exceeding 120 ppm. Contrary, the contamination offset reduces strongly for spots closer to the limb with the addition of LD.

thumbnail Fig. 4

Effect of limb-darkening on stellar contamination factor and observed transit depth for stars of different spectral types. Top: illustration of the stellar disc with a single spot at four different positions, i.e. at different stellar rotation phases. Bottom: corresponding stellar contamination spectra for F5V, G5V, K5V, and M5V stars with LD (solid) and without LD (dashed). The secondary y-axis represents the change in transit depth due to stellar contamination assuming a planet with a uniform transit depth of 1% at all wavelengths.

thumbnail Fig. 5

Same as Fig. 3, but including LD.

4 Stellar surface retrievals

4.1 Model setup

SAGE can be used to retrieve the properties and positions of active regions from observations of the rotationally modulated flux of the star. To do so, we used the capacity of SAGE to predict the photometric variability produced by a given spot configuration by rotating a stellar surface model as described in Sect. 2.4. Here, SAGE includes the full stellar atmospheric effects described above, most importantly the LD. Furthermore, the inclination of the star and its impact on the flux variability is also accounted for.

To derive the properties of the stellar surface inhomo-geneities, we created model predictions from SAGE and used an MCMC framework to evaluate these against the observed data. The fitted parameters are: the stellar rotation period; the stellar inclination (i.e. the angle between the sky plane and the stellar rotation axis); the latitude, longitude, size, and contrast for each active region; and a jitter term to account for any additional photometric uncertainties.

We recommend imposing uniform distributions as priors for all parameters, unless prior information is available about the distribution of active regions. The width of the priors (for latitude and longitude) can be set to cover the entire stellar surface or restricted to remove the possibility of polar spots.

One of the main limitations in accurately constraining stellar surfaces using the variability models is that the retrievals are performed under the assumption that the observed flux variability of a star is due to a limited number of large active regions. Thus, this technique is unable to capture the variability associated with many small inhomogeneities, randomly distributed across the stellar surface. However, it does provide accurate estimates of the relative coverage of the active regions over the different rotation phases of the star.

4.2 Application to WASP-69

We used the setup described above to study the flux variability of the exoplanet host WASP-69 in data collected by the Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2014). The K5-type star has an effective temperature of 4715 ± 50 K and hosts a puffy Saturn-mass planet (WASP-69b, 0.26 MJ, 1.06 RJ), orbiting with a period of 3.87d (Anderson et al. 2014).

Anderson et al. (2014) detected emission features in the CaH + K line cores, indicative of activity, and estimated an activity index of log RHK −4.54, which puts the age of the star at ~0.8 Gyr (Mamajek & Hillenbrand 2008). The stellar rotational period calculated from light curve modulation is 23.07 ± 0.16 d, broadly agreeing with that derived from the spectroscopic rotational broadening method, <18.7 ± 3.5 d. This suggests that the star is seen equator-on.

The atmosphere of the planet WASP-69b has been studied using the transmission spectroscopy technique in both low and high resolution, which has led to robust detection of multiple chemical species such as CH4, NH3, CO, C2H2, and H2O (Guilluy et al. 2022). The presence of Rayleigh scattering by molecular hydrogen in the upper atmosphere has also been detected (Murgas et al. 2020). However, unocculted spots on the stellar surface are speculated to be a possible alternative explanation for the bluewards rise in transit depth.

Here, we used the recently acquired observations of WASP-69 from TESS Sector 55 (UT 2022-August-05 to UT 2022-September-01). We used the 20 s cadence light curves generated by the Science Processing and Operations Center (SPOC; Jenkins et al. 2016) and employed the lightkurve package (Lightkurve Collaboration 2018) to calculate the flux using Simple Aperture Photometry (SAP). The flux outliers are removed using a running median method. The points above 3×MAD3 are removed along with the 6 transits of WASP-69b. In Sector 55, the median background flux increased by a factor of 4 between 2805 and 2807 TESS barycentric Julian date (TBJD) due to stray light from the Moon, and we removed all observations taken during this time window. A clear offset is present in the light curve at 2806 TBJD. We calculated the average flux over a 2-day window before and after the offset and shifted the data taken before 2806 TBJD accordingly to remove this offset. The processed light curve shows a periodicity in flux with a period of ~22.1 d, in good agreement with the previously determined stellar rotation period. As the total science data are collected over a period of 26.26 days, we obtain a little over one complete rotation of the star. In addition, we find no evidence of stellar flares in the data. For the sake of computational speed, we phase-folded the data (with Prot = 22.1 d) and binned them into 28 uniformly spaced time intervals before progressing with the stellar surface retrievals.

To limit the degeneracy between temperature and size of active regions, we fixed the contrast of the spots in our simulations. We set the clear photospheric temperature to 4715 K and, using the linear relation between the clear and spotted photosphere discussed in Sect. 3, we set the temperature of spots to 3590 K. We convolved the PHOENIX spectra of both clear and spotted photospheres with the TESS bandpass, finding a spot contrast (ρspot) of 0.26. We used LDCU to compute quadratic LD coefficients assuming the temperature of the clear photosphere, and we kept these coefficients fixed during the retrieval. The system is likely seen equator-on, evidenced by the fact that the stellar rotation period computed from the V sin(i) (<18.7 ± 3.5 d) and the rotation period derived from the TESS light curve (22.1 d) agree. We, therefore, fixed the stellar inclination to 90°.

We assumed uniform priors for latitude, longitude and size of each spot, restricting the latitudes to the northern hemisphere of the star to avoid bimodal posteriors due to the symmetry between the north and south. The spot size was limited to be between 1° and 10°. This range covers small Sun-like spots to large M-dwarf-type spots. To determine the number of spots, we first ran a retrieval with a single spot and then iteratively increased the number of spots from 1 to 5. Each time, we also calculated the absolute filling factor of the active regions. The two-spot model has the lowest Bayesian information criterion, and we thus deem it the most appropriate. We also find that increasing the number of spots to model the light curve does not have a significant impact on the absolute filling factor of the active regions.

In Fig. 6, we present the phase-folded TESS data together with 100 variability models randomly selected from the posteriors of our MCMC analysis. We also show a model of the stellar surface corresponding to the peak of the posterior distribution. We present the posteriors of our retrieval in Fig. B.1 and list the derived parameters together with their uncertainties in Table 2. While a range of spot latitudes is possible given the available data, our analysis strongly rejects the possibility of polar spots on WASP-69. From the posterior distribution, the correlation between the spot size and latitude is clearly visible; the spots at higher latitudes are larger and vice versa. This relation preserves the projected filling factor of the active regions, which modulate both the stellar flux and contamination.

We measure the absolute filling factor of the active regions to be fs, absolute =0.770.19+0.29 % ${f_{s,{\rm{absolute}}}} = 0.77_{ - 0.19}^{ + 0.29}{\rm{\% }}$. Assuming a uniform radius of planet b (Rp = 11.85 R) at different wavelengths, we also constrained the stellar contamination in the optical regime and its variation with the different rotation phases of the star as shown in Fig. 6.

To probe the impact of the i = 90° assumption, we also performed a separate analysis with Gaussian priors on the inclination (µ = 90°, σ = 15°), truncated at <90° to avoid degeneracies due to the symmetry of the problem. Here, the σ of 15° is calculated by taking the lower limit on inclination derived using V sin(i) of 2.2 ± 0.4 km sec−1, a stellar radius of 0.813 ± 0.028 R and Prot of 22.1 days (Anderson et al. 2014). We used uniform priors on the entire stellar surface and spot sizes between 1° and 10°. The posteriors of the analysis are presented in Fig. B.2. The retrievals qualitatively reproduce those found for i set to 90°, namely two spots at mid-latitudes separated by 85° in longitude. The posterior for the inclination peaks at 68°, but extends out to include solutions with 90°. Due to the inclination, the symmetry between the north and south is broken, and the derived spot sizes depend on the hemisphere on which the spot is located. While we adopted the 90° inclination for WASP-69 due to the limited quality and quantity of data available, this assumption might need revision once more detailed observations are available.

Table 2

Fitted spot parameters and their priors for the stellar surface retrievals of WASP-69.

5 Discussion

In this paper, we have used a stellar surface pixelation approach to investigate the two-fold impact of active regions on exoplanet observations. First (see Sect. 3), we studied their impact on transmission spectra, including a thorough treatment of projection effects and LD. Second (see Sect. 4), we illustrated the photometric effect of spots by rotating our model and devised a retrieval approach to constrain the properties of active regions from photometric observations.

5.1 Advantages and novelty

For all but the most quiet stars, the stellar surface, and most importantly the quantity and position of active regions vary with time, leading to varying stellar contamination in transmission spectra. This makes it difficult to combine multi-epoch transit observations to enhance detection of atmospheric features. One of the main advantages of SAGE is its ability to use the stellar surface map, obtained using photometric time series observations, to directly constrain the stellar contamination at different rotation phases of the star. This feature is currently not included in tools like Starry (Luger et al. 2019) and SOAP (Boisse et al. 2012), which are mainly used to create stellar surface maps from photometric observations. In addition, SAGE also calculates the stellar contamination models using a numerical approach accounting for LD, which is more precise than the purely analytical models.

An important result of our study is the impact of LD on stellar contamination. As an active region moves closer to the limb of the star, the overall surface brightness of the stellar photosphere decreases, and so does its impact on the measurements. In the case of spots located close to the centre of the stellar disc, the non-LD model underestimates the stellar contamination. In contrast, for spots closer to the limb, the non-LD model overestimates the stellar contamination. This result has been tested with multiple simulated stars of different spectral types with different spot temperatures.

As LD is a chromatic effect, the changing surface brightness of an active region as it approaches the limb also has a wavelength-dependent effect on the observed transit depths. The LD effect is stronger at shorter wavelengths, impacting the shape of the contamination spectrum. In Fig. 4, we presented the diversity in contamination spectra for F5V, G5V, K5V, and M5V stars with spots at different distances from the centre of the stellar disc. In all cases, as a spot moves close to the limb of the star, the slope of the contamination spectrum flattens significantly for the wavelength range 4500–7500Å. At shorter wavelengths of 30004500 Å, the slope can even invert for F-, G-, and K-type stars, leading to a decrease in contamination at short wavelengths. This is contrary to non-LD models, which predict a monotonous rise in stellar contamination towards shorter wavelengths, an effect that can be misinterpreted as atmospheric Rayleigh scattering in exoplanets (Oshagh et al. 2014; Murgas et al. 2020). Thus, the inclusion of LD in stellar contamination calculations is vital for disentangling planetary and stellar effects on the transmission spectrum.

thumbnail Fig. 6

Flux variability models and the expected variability in contamination spectrum of WASP-69. Top: model stellar surface map for variability models associated with the peak of the posteriors. Middle: TESS data of WASP-69, phase-folded on the stellar rotation period. The grey points are the individual measurements, and the black points are the binned measurements. The red curves show 100 variability models randomly selected from the posterior distribution of the fit. Bottom: observed radius of planet b at different rotational phases of the star and different wavelengths.

5.2 Further development and caveats

Although SAGE takes more stellar atmospheric effects into account compared to the existing routines, it is still limited by the complexity of stellar atmospheres, and we have made some simplifying assumptions in this study.

One of the central assumptions of the model is that the star is rotating as a rigid body. This is not the case, as observed clearly for the Sun, which is known to show differential rotation (Balthasar et al. 1986; Beck 2000), rotating significantly faster at the equator than at the poles. Not including for differential rotation can affect our retrievals of stellar surface properties from photometric flux variability, as the contributions of equatorial and polar spots will have different timescales. Similarly, we also assumed that active regions are not evolving within the investigated timescales. This assumption works well for fast rotators where the stellar rotation period is much smaller than the evolutionary timescales of active regions such as AU Mic (Szabó et al. 2021), or for stars with long-lived and stable spots that last a few stellar rotations, for example HD 189733 and α Cen B (Dumusque et al. 2014).

Stellar surface maps are also affected by an important degeneracy between the position, size and temperature of active regions. The flux variability due to a hot, large spot (or a cluster of small spots) on the stellar surface can be reproduced equally well with a smaller, cooler spot. The potential presence of multiple stellar spots with different temperatures and sizes further complicates the scenario. To break this degeneracy, an independent measurement of the size or temperature of the active regions is needed. This may be attempted with multi-colour observations, or transit observations that show short-term brightening due to the starspot crossings. Another complexity lies in the LD behaviour of active regions. While we assumed the same LD coefficients for active regions as for the clear photosphere, their LD behaviour might differ in reality.

Finally, we currently do not include high-resolution effects such as the convective blueshift and variations in line shapes, limiting the applicability of the code to low- and medium-resolution observations. While updates are foreseen in future releases, we expect these effects to be minimal for the results presented in this paper.

Using the pre-tabulated models of stellar emission spectra, we also assumed that the emergent flux from each photospheric pixel is the same as that of a disc-integrated spectrum of a star with the same temperature, metallicity and gravity. This assumption is commonly used in this field (Pont et al. 2008; McCullough et al. 2014; Rackham et al. 2017), but we do acknowledge the need for more sophisticated spectral models.

6 Summary and conclusions

In this paper we examine the effect of stellar activity on exo-planet transmission spectra using realistic models of stellar photospheres implemented through a pixelation approach. In particular, we investigate the effect of the spot position and stellar LD on the contamination spectrum at optical wavelengths (3000–10 000Å). The developed tool, SAGE, not only provides the contamination spectra needed to correct exoplanet transmission spectra, but also allows stellar surface maps to be derived from photometric stellar rotation curves. Our key findings are summarised as follows:

  • 1.

    The extent of stellar contamination depends on the effective temperature of the star. From simulations of F5V, G5V, K5V, and M5V stars, we detect a clear negative correlation between the temperature of the star and the amplitude of the stellar contamination effect.

  • 2.

    The shape and amplitude of the stellar contamination depends strongly on the location of the spot. Due to projection effects, spots near the centre produce more pronounced effects.

  • 3.

    Stellar LD further exacerbates this projection effect. For spots located close to the centre of the stellar disc (µ ≳ 0.6), the stellar contamination effect is enhanced due to LD, but for spots located close to the limb (µ ≲ 0.6), the contamination decreases.

  • 4.

    Limb-darkening also affects the shape of the contamination spectrum due to its colour dependence, especially for active regions located close to the limb of the star. For spots near the limb, we find a flattened contamination spectrum between 4500 and 7500 Å and even a negative slope for <4500 Å for the simulated F-, G-, and K-type stars. Thus, the monotonous rise in transit depth associated with unocculted spots depends on the location of spots.

  • 5.

    We also used SAGE to model the flux variability of WASP-69 in TESS observations. We find that the ~22.1 d periodicity in the flux can be modelled as an effect of active regions on the stellar surface. We estimate that the spots cover ~1 % of the total stellar surface.

To conclude, the purpose of the developed tool is to calculate precise stellar contamination spectra, accounting for realistic and time-varying spot configurations, and couple these computations with the photometric flux variability that encodes the distribution of spots on the star. The increasingly detailed understanding of stellar effects on transmission spectra will enable a meaningful comparison of multi-epoch transmission spectra, opening the door towards studies of variability in exoplanet atmospheres.

Acknowledgements

H.C. acknowledges the fruitful discussions with Xavier Dumusque, Yinan Zhao, Alex Pietrow and Matthew Battley. H.C., M.L., B.A. and D.P. further acknowledge the support of the Swiss National Science Foundation under grant number PCEFP2_194576. 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. D.P. further acknowledges funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (project Four Aces; grant agreement No 724427). A.D. and D.P. acknowledge financial support from the Swiss National Science Foundation (SNSF) for project 200021_200726. This paper includes data collected by the TESS mission. Funding for the TESS mission is provided by the NASA’s Science Mission Directorate. The computations in this paper were performed at the University of Geneva on the “Yggdrasil” HPC cluster.

Appendix A Equivalence of stellar contamination calculations in the case of no secondary effects on the stellar disc

Under the assumption of a planet transiting a uniformly illuminated stellar disc, the measured transit depth is given by δtransit ,λ=(RpR)λ2=Fout ,λFin,λFout ,λ.${\delta _{{\rm{transit}},\lambda }} = \left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2 = {{{F_{{\rm{out}},\lambda }} - {F_{{\rm{in}},\lambda }}} \over {{F_{{\rm{out}},\lambda }}}}.$(A.1)

Here, δtransit,λ is the flux drop associated with the planetary transit, Fout,λ and Fin,λ are the spectra of the star measured before and during the transit.

In the case of a planet transiting an active host star with FClear,λ is the spectrum of the clear photosphere and FActive,λ is the spectrum of an active region. The out-of-transit spectrum of the star (Fout,λ) under the assumption of a stellar disc with no LD effect is Fout ,λ,t=[ 1fproj (t) ]Fclear ,λ+fproj (t)Factive ,λ.${F_{{\rm{out}},\lambda ,{\rm{t}}}} = \left[ {1 - {f_{{\rm{proj}}}}(t)} \right]{F_{{\rm{clear}},\lambda }} + {f_{{\rm{proj}}}}(t){F_{{\rm{active}},\lambda }}.$

Here, fproj(t) is the projected covering fraction of the active regions with respect of an observer at any particular time. Assuming a planet with no atmosphere and a true transit depth of (RpR)λ2$\left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2$, the in-transit spectrum (Fin, λ) of the star is given by Fin ,λ,t=[ 1fproj (t)(RpR)λ2 ]Fclear ,λ+fproj (t)Factive ,λ.${F_{{\rm{in}},\lambda ,{\rm{t}}}} = \left[ {1 - {f_{{\rm{proj}}}}(t) - \left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2} \right]{F_{{\rm{clear}},\lambda }} + {f_{{\rm{proj}}}}(t){F_{{\rm{active}},\lambda }}.$

Thus, the observed transit-depth (Rp^R)λ,t2$\left( {{{\widehat {{R_p}}} \over {{R_ \star }}}} \right)_{\lambda ,t}^2$ is given by 1(Rp^R)λ,t2=Fin ,λ,tFout ,λ,t                       =[ 1fproj (t)(RpR)λ2 ]Fclear ,λ+fproj (t)Factive ,λ[ 1fproj (t) ]Fclear ,λ+fproj (t)Factive ,λ                       =1fproj (t)+fproj (t)Factive ,λFclear, λ(RpR)λ21fproj (t)+fproj (t)Factive ,λFclear, ,λ                     =1(RpR)λ21fproj (t)+fproj (t)Factive ,λFclear ,λ(Rp^R)λ,t2=(RpR)λ21fproj (t)(1Factive ,λFclear, ,λ).$\eqalign{ & 1 - \left( {{{\widehat {{R_p}}} \over {{R_ \star }}}} \right)_{\lambda ,t}^2 = {{{F_{{\rm{in}},\lambda ,{\rm{t}}}}} \over {{F_{{\rm{out}},\lambda ,{\rm{t}}}}}} \cr & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, = {{\left[ {1 - {f_{{\rm{proj}}}}(t) - \left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2} \right]{F_{{\rm{clear}},\lambda }} + {f_{{\rm{proj}}}}(t){F_{{\rm{active}},\lambda }}} \over {\left[ {1 - {f_{{\rm{proj}}}}(t)} \right]{F_{{\rm{clear}},\lambda }} + {f_{{\rm{proj}}}}(t){F_{{\rm{active}},\lambda }}}} \cr & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, = {{1 - {f_{{\rm{proj}}}}(t) + {f_{{\rm{proj}}}}(t){{{F_{{\rm{active}},\lambda }}} \over {{F_{{\rm{clear,}}\lambda }}}} - \left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2} \over {1 - {f_{{\rm{proj}}}}(t) + {f_{{\rm{proj}}}}(t){{{F_{{\rm{active}},\lambda }}} \over {{F_{{\rm{clear,}},\lambda }}}}}} \cr & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, = 1 - {{\left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2} \over {1 - {f_{{\rm{proj}}}}(t) + {f_{{\rm{proj}}}}(t){{{F_{{\rm{active}},\lambda }}} \over {{F_{{\rm{clear}},\lambda }}}}}} \cr & \left( {{{\widehat {{R_p}}} \over {{R_ \star }}}} \right)_{\lambda ,t}^2 = {{\left( {{{{R_p}} \over {{R_ \star }}}} \right)_\lambda ^2} \over {1 - {f_{{\rm{proj}}}}(t)\left( {1 - {{{F_{{\rm{active}},\lambda }}} \over {{F_{{\rm{clear,}},\lambda }}}}} \right)}}.\, \cr} $

The stellar contamination factor between the measured transit depth and the true transit depth of the planet is ϵλ,t=11fproj (t)(1Factive ,λFclear ,λ).${_{\lambda ,t}} = {1 \over {1 - {f_{{\rm{proj}}}}(t)\left( {1 - {{{F_{{\rm{active}},\lambda }}} \over {{F_{{\rm{clear}},\lambda }}}}} \right)}}.$(A.2)

Rearranging equation A.2, we get ϵλ,t=11fproj (t)+fproj (t)Factive ,λFclear ,λ        =Fclear ,λ[ 1fproj (t) ]Fclear ,λ+fproj (t)Factive ,λ.        =Fclear ,λFout ,λ        = Clear stellar spectrum  Stellar spectrum at epoch of interest $\matrix{ {{_{\lambda ,t}} = {1 \over {1 - {f_{{\rm{proj}}}}(t) + {f_{{\rm{proj}}}}(t){{{F_{{\rm{active}},\lambda }}} \over {{F_{{\rm{clear}},\lambda }}}}}}} \hfill \cr {\,\,\,\,\,\,\,\, = {{{F_{{\rm{clear}},\lambda }}} \over {\left[ {1 - {f_{{\rm{proj}}}}(t)} \right]{F_{{\rm{clear}},\lambda }} + {f_{{\rm{proj}}}}(t){F_{{\rm{active}},\lambda }}}}.} \hfill \cr {\,\,\,\,\,\,\,\, = {{{F_{{\rm{clear}},\lambda }}} \over {{F_{{\rm{out}},\lambda }}}}} \hfill \cr {\,\,\,\,\,\,\,\, = {{{\rm{Clear stellar spectrum}}} \over {{\rm{Stellar spectrum at epoch of interest}}}}} \hfill \cr } $(A.3)

This equation can also be useful for constraining the contamination effect without any prior knowledge about the properties of the active regions like its covering fraction or temperature during the planetary transit. However, that requires obtaining a stellar spectrum when no surface inhomogeneities are present. This is quite difficult as some level of inhomogeneities will always be present even during the minima of the magnetic activity of the star.

Appendix B Corner plot from the stellar surface retrievals

thumbnail Fig. B.1

Posterior distribution of the spot parameters from the stellar surface retrievals of WASP-69 using a two-spot model.

thumbnail Fig. B.2

Same as Fig. B.1, but with stellar inclination (i) as a free parameter.

References

  1. Alam, M. K., López-Morales, M., Nikolov, N., et al. 2020, AJ, 160, 51 [NASA ADS] [CrossRef] [Google Scholar]
  2. Allart, R., Bourrier, V., Lovis, C., et al. 2018, Science, 362, 1384 [Google Scholar]
  3. Almenara, J. M., Bonfils, X., Forveille, T., et al. 2022, A&A, 667, L11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Anderson, D. R., Collier Cameron, A., Delrez, L., et al. 2014, MNRAS, 445, 1114 [NASA ADS] [CrossRef] [Google Scholar]
  5. Balthasar, H., Vazquez, M., & Woehl, H. 1986, A&A, 155, 87 [Google Scholar]
  6. Beck, J. G. 2000, Sol. Phys., 191, 47 [NASA ADS] [CrossRef] [Google Scholar]
  7. Berta, Z. K., Charbonneau, D., Bean, J., et al. 2011, ApJ, 736, 12 [NASA ADS] [CrossRef] [Google Scholar]
  8. Boisse, I., Bonfils, X., & Santos, N. C. 2012, A&A, 545, A109 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Borsa, F., Allart, R., Casasayas-Barris, N., et al. 2021, A&A, 645, A24 [EDP Sciences] [Google Scholar]
  10. Bruno, G., Lewis, N. K., Stevenson, K. B., et al. 2018, AJ, 156, 124 [NASA ADS] [CrossRef] [Google Scholar]
  11. Castelli, F., & Kurucz, R. L. 2003, IAU Symp., 210, A20 [Google Scholar]
  12. Chapman, G. A., & McGuire, T. E. 1977, ApJ, 217, 657 [CrossRef] [Google Scholar]
  13. Coelho, P. R. T. 2014, MNRAS, 440, 1027 [Google Scholar]
  14. Crane, P. C., Floyd, L. E., Cook, J. W., et al. 2004, A&A, 419, 735 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Désert, J.-M., Charbonneau, D., Demory, B.-O., et al. 2011, ApJS, 197, 14 [Google Scholar]
  16. Dumusque, X., Boisse, I., & Santos, N. C. 2014, ApJ, 796, 132 [NASA ADS] [CrossRef] [Google Scholar]
  17. Ehrenreich, D., Lovis, C., Allart, R., et al. 2020, Nature, 580, 597 [Google Scholar]
  18. Espinoza, N., & Jordán, A. 2015, MNRAS, 450, 1879 [Google Scholar]
  19. Evans, T. M., Sing, D. K., Wakeford, H. R., et al. 2016, ApJ, 822, L4 [NASA ADS] [CrossRef] [Google Scholar]
  20. Fontenla, J., White, O. R., Fox, P. A., Avrett, E. H., & Kurucz, R. L. 1999, ApJ, 518, 480 [NASA ADS] [CrossRef] [Google Scholar]
  21. Froebrich, D., Scholz, A., Eislöffel, J., & Stecklum, B. 2020, MNRAS, 497, 4602 [NASA ADS] [CrossRef] [Google Scholar]
  22. Gibson, N. P., Nikolov, N., Sing, D. K., et al. 2017, MNRAS, 467, 4591 [NASA ADS] [CrossRef] [Google Scholar]
  23. Guilluy, G., Giacobbe, P., Carleo, I., et al. 2022, A&A, 665, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Henry, G. W., Eaton, J. A., Hamer, J., & Hall, D. S. 1995, ApJS, 97, 513 [CrossRef] [Google Scholar]
  25. Herrero, E., Ribas, I., Jordi, C., et al. 2016, A&A, 586, A131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Hoeijmakers, H. J., Ehrenreich, D., Kitzmann, D., et al. 2019, A&A, 627, A165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Husser, T. O., Wende-von Berg, S., Dreizler, S., et al. 2013, A&A, 553, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Jenkins, J. M., Twicken, J. D., McCauliff, S., et al. 2016, SPIE Conf. Ser., 9913, 99133E [Google Scholar]
  29. Juvan, I. G., Lendl, M., Cubillos, P. E., et al. 2018, A&A, 610, A15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. JWST Transiting Exoplanet Community Early Release Science Team, Ahrer, E.-M., Alderson, L., et al. 2023, Nature, 614, 649 [NASA ADS] [CrossRef] [Google Scholar]
  31. Kurucz, R. L. 1979, ApJS, 40, 1 [NASA ADS] [CrossRef] [Google Scholar]
  32. Lecavelier Des Etangs, A., Vidal-Madjar, A., Désert, J. M., & Sing, D. 2008, A&A, 485, 865 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Lendl, M., Delrez, L., Gillon, M., et al. 2016, A&A, 587, A67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Lightkurve Collaboration (Cardoso, J. V. d. M., et al.) 2018, Astrophysics Source Code Library [record ascl:1812.013] [Google Scholar]
  35. Luger, R., Agol, E., Foreman-Mackey, D., et al. 2019, AJ, 157, 64 [Google Scholar]
  36. Mamajek, E. E., & Hillenbrand, L. A. 2008, ApJ, 687, 1264 [Google Scholar]
  37. Mandal, S., Krivova, N. A., Cameron, R., & Solanki, S. K. 2021, A&A, 652, A9 [Google Scholar]
  38. Maxted, P. F. L. 2016, A&A, 591, A111 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. McCullough, P. R., Crouzet, N., Deming, D., & Madhusudhan, N. 2014, ApJ, 791, 55 [NASA ADS] [CrossRef] [Google Scholar]
  40. Meunier, N., Desort, M., & Lagrange, A. M. 2010, A&A, 512, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Močnik, T., Clark, B. J. M., Anderson, D. R., Hellier, C., & Brown, D. J. A. 2016, AJ, 151, 150 [CrossRef] [Google Scholar]
  42. Montalto, M., Boué, G., Oshagh, M., et al. 2014, MNRAS, 444, 1721 [Google Scholar]
  43. Morris, B. M., Hebb, L., Davenport, J. R. A., Rohn, G., & Hawley, S. L. 2017, ApJ, 846, 99 [NASA ADS] [CrossRef] [Google Scholar]
  44. Murgas, F., Chen, G., Nortmann, L., Palle, E., & Nowak, G. 2020, A&A, 641, A158 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  45. Nikolov, N., Sing, D. K., Burrows, A. S., et al. 2015, MNRAS, 447, 463 [NASA ADS] [CrossRef] [Google Scholar]
  46. Oshagh, M., Boisse, I., Boué, G., et al. 2013, A&A, 549, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Oshagh, M., Santos, N. C., Ehrenreich, D., et al. 2014, A&A, 568, A99 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Pecaut, M. J., & Mamajek, E. E. 2013, ApJS, 208, 9 [Google Scholar]
  49. Pietrow, A. G. M., Kiselman, D., de la Cruz Rodríguez, J., et al. 2020, A&A, 644, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Pont, F., Gilliland, R. L., Moutou, C., et al. 2007, A&A, 476, 1347 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Pont, F., Knutson, H., Gilliland, R. L., Moutou, C., & Charbonneau, D. 2008, MNRAS, 385, 109 [Google Scholar]
  52. Pont, F., Sing, D. K., Gibson, N. P., et al. 2013, MNRAS, 432, 2917 [NASA ADS] [CrossRef] [Google Scholar]
  53. Rackham, B., Espinoza, N., Apai, D., et al. 2017, ApJ, 834, 151 [Google Scholar]
  54. Rackham, B. V., Apai, D., & Giampapa, M. S. 2018, ApJ, 853, 122 [Google Scholar]
  55. Rackham, B. V., Apai, D., & Giampapa, M. S. 2019, AJ, 157, 96 [Google Scholar]
  56. Reiners, A. 2012, Liv. Rev. Sol. Phys., 9, 1 [Google Scholar]
  57. Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2014, SPIE Conf. Ser., 9143, 914320 [Google Scholar]
  58. Seager, S., & Sasselov, D. D. 2000, ApJ, 537, 916 [Google Scholar]
  59. Shields, A. L., Ballard, S., & Johnson, J. A. 2016, Phys. Rep., 663, 1 [NASA ADS] [CrossRef] [Google Scholar]
  60. Silva, A. V. R. 2003, ApJ, 585, L147 [Google Scholar]
  61. Silva-Valio, A., Lanza, A. F., Alonso, R., & Barge, P. 2010, A&A, 510, A25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Sing, D. K., Pont, F., Aigrain, S., et al. 2011, MNRAS, 416, 1443 [NASA ADS] [CrossRef] [Google Scholar]
  63. Sing, D. K., Fortney, J. J., Nikolov, N., et al. 2016, Nature, 529, 59 [Google Scholar]
  64. Snellen, I. A. G., Albrecht, S., de Mooij, E. J. W., & Le Poole, R. S. 2008, A&A, 487, 357 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Solanki, S. K. 2003, A&A Rev., 11, 153 [NASA ADS] [CrossRef] [Google Scholar]
  66. Szabó, G. M., Gandolfi, D., Brandeker, A., et al. 2021, A&A, 654, A159 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Tregloan-Reed, J., Southworth, J., & Tappert, C. 2013, MNRAS, 428, 3671 [Google Scholar]
  68. Tregloan-Reed, J., Southworth, J., Burgdorf, M., et al. 2015, MNRAS, 450, 1760 [Google Scholar]
  69. Wyttenbach, A., Ehrenreich, D., Lovis, C., Udry, S., & Pepe, F. 2015, A&A, 577, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Zellem, R. T., Swain, M. R., Roudier, G., et al. 2017, ApJ, 844, 27 [NASA ADS] [CrossRef] [Google Scholar]

2

As defined in Sect. 2.1, the absolute filling factor refers to the ratio of the total surface area of the active regions to the surface area of the star, while the projected filling factor is the ratio of the projected surface area of active regions to the surface area of the stellar disc.

3

Mean absolute deviation from the median for each 15-point window

All Tables

Table 1

Model spectrum parameters used in the presented simulations.

Table 2

Fitted spot parameters and their priors for the stellar surface retrievals of WASP-69.

All Figures

thumbnail Fig. 1

A model stellar grid and spectral setup in SAGE. Left: Stellar sphere with LD and spots projected onto a two-dimensional Cartesian pixel grid. Right: normalised model spectra used, obtained from the PHOENIX library. For the clear photosphere and for spots, we used temperatures of 4500 and 3500 K, respectively.

In the text
thumbnail Fig. 2

Flux variability and stellar contamination models for different rotation phases of a synthetic star. Top: Stellar surface grid at different rotational phases for an (equator-on) K5V star with ten 2–3° spots randomly distributed across its surface. The effective temperature of the spots and the clear photosphere is 3500 and 4500 K, respectively. Middle: disc-integrated flux calculated between 3000 and 10 000 Å as a function of the stellar rotational phase. Bottom: wavelength-dependent stellar contamination factor (ϵλ) at rotation phases corresponding to those shown in the top panel.

In the text
thumbnail Fig. 3

Variation in contamination offset with location and size of active regions for stars of different spectral types. Top: Stellar contamination offset for a single spot of 5° angular size, placed at different locations on the stellar disc for stars of different spectral types. Bottom: Stellar contamination offset for the same stars with spot sizes ranging from 1 to 7° for different longitudes. LD was not included in these simulations.

In the text
thumbnail Fig. 4

Effect of limb-darkening on stellar contamination factor and observed transit depth for stars of different spectral types. Top: illustration of the stellar disc with a single spot at four different positions, i.e. at different stellar rotation phases. Bottom: corresponding stellar contamination spectra for F5V, G5V, K5V, and M5V stars with LD (solid) and without LD (dashed). The secondary y-axis represents the change in transit depth due to stellar contamination assuming a planet with a uniform transit depth of 1% at all wavelengths.

In the text
thumbnail Fig. 5

Same as Fig. 3, but including LD.

In the text
thumbnail Fig. 6

Flux variability models and the expected variability in contamination spectrum of WASP-69. Top: model stellar surface map for variability models associated with the peak of the posteriors. Middle: TESS data of WASP-69, phase-folded on the stellar rotation period. The grey points are the individual measurements, and the black points are the binned measurements. The red curves show 100 variability models randomly selected from the posterior distribution of the fit. Bottom: observed radius of planet b at different rotational phases of the star and different wavelengths.

In the text
thumbnail Fig. B.1

Posterior distribution of the spot parameters from the stellar surface retrievals of WASP-69 using a two-spot model.

In the text
thumbnail Fig. B.2

Same as Fig. B.1, but with stellar inclination (i) as a free parameter.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.