Increasing diurnal and seasonal amplitude of atmospheric methane mole fraction in Central Siberia between 2010–2021

Siberia's vast wetlands, permafrost, and boreal forests are significant, but their sources of methane (CH₄) are poorly quantified. Using vertical CH₄ profiles and meteorological data from the ZOtino Tall Tower Observatory (ZOTTO; 60°48^′ N, 89°21^′ E) in Central Siberia, we analyse long-term trends in CH₄ growth rates, seasonal patterns, and diurnal cycles from 2010 to 2021. Our results show a persistent long-term trend in CH₄ mole fractions and an insignificant increasing seasonal cycle amplitude, (2.12 ppb yr⁻¹, p=0.12) along with a pronounced late-summer CH₄ peak. Diurnal analysis reveals a growing late summer (July–October) CH₄ amplitude over the analysed decade (5.29 ppb yr⁻¹, p=0.001), driven by rising nighttime fluxes correlated with late summer soil temperature (R²=0.65, p<0.01), soil moisture (R²=0.36, p=0.032) and with preceding spring snow depth (R²=0.54, p=0.03). Notably high nighttime CH₄ fluxes occurred in 2012 and 2019 mainly due to wildfires. These findings suggest that increasing late-summer CH₄ emissions, primarily from wetlands to the west and southwest of ZOTTO, dominantly contribute to the overall CH₄ rise. Our study underscores the importance of continuous, high-frequency greenhouse gas observations for accurately quantifying regional CH₄ trends.

1 Introduction

Methane (CH₄) is an important greenhouse gas, accounting for approximately 16 % of global greenhouse gas radiative forcing and representing the second largest contributor to current anthropogenic warming (IPCC, 2023). Since 1970, the global mean atmospheric CH₄ mole fraction rose from 1630 to 1774 ppb by 1999. This was followed by a period of stalled growth between 1999 and 2006, after which CH₄ levels increased to 1834 ppb by 2015 and further to 1879 ppb by 2020 (Lan et al., 2022). The renewed increase has been primarily attributed to biogenic sources, particularly wetlands, rather than fossil fuel emissions or changes in atmospheric CH₄ sinks (Basu et al., 2022; Lan et al., 2021; Nisbet et al., 2016, 2019). This rise in wetland emissions is suggested to result from enhanced microbial methane production driven by higher soil moisture content, warmer soil temperature, and extended periods of inundation in tropical and high-latitude regions, all of which promote anaerobic conditions favourable for methanogenesis (Basu et al., 2022; Bridgham et al., 2013). The cause of the 1999–2006 plateau, however, remains unclear, largely due to limited observational data during that period (Nisbet et al., 2019, 2023).

Arctic-Boreal regions are characterised by extensive wetlands, permafrost, and organic-rich soils that act as both sources and sinks of CH₄ (Saunois et al., 2025). Warming in boreal zones has recently been observed to occur three to four times faster than the global average (Rantanen et al., 2022; IPCC, 2023), fuelling significant concerns given the positive feedback between CH₄ emissions and climate warming (Yvon-Durocher et al., 2014; Chang et al., 2021; Zhang et al., 2017; Canadell et al., 2021). Undisturbed and disturbed boreal and temperate peatlands emit 30 and 0.1 Tg CH₄ yr⁻¹, respectively (Olsson et al., 2019), and the northern permafrost region emits 15–38 Tg CH₄ for the 2000–2020 period (Hugelius et al., 2024). While these estimates highlight the significance of boreal CH₄ emissions at a global scale, considerable uncertainty remains regarding their regional responses to long-term environmental changes, resulting in high uncertainty in their estimated contribution to the overall CH₄ budget (Saunois et al., 2025, 2020, 2016; Kuhn et al., 2021). This knowledge gap is primarily due to the boreal forest being sparsely monitored by long-term atmospheric and ecosystem observatories, highlighting the need to leverage existing datasets in this region to refine CH₄ budgets with better large-scale spatial coverage and fine-scale temporal precision, such as the ZOtino Tall Tower Observatory (ZOTTO).

The ZOTTO facility, situated in central Siberia, was established in 2006, and from April 2009 to February 2022, it had been continuously measuring CH₄ at multiple heights up to 301 m along with other atmospheric gases and their isotopic compositions as well as meteorological data (Winderlich et al., 2010; Tran et al., 2024). This long-term monitoring effort makes ZOTTO a valuable atmospheric research station, providing high-time-resolution (half-minute frequency) CH₄ measurements in a key high-latitude (above 50° N) region. An early analysis of ZOTTO data by Winderlich (2012), covering the period 2009–2011, identified pronounced CH₄ mole fraction spikes during mid- to late-summer between July and September. These were suggested to be due to microbial activity in nearby wetlands and episodic emissions from Siberian forest fires during the summer. In light of the globally observed post-2006 increase in atmospheric CH₄ now widely linked to enhanced wetland emissions, there is a compelling motivation to revisit and extend the analysis of Winderlich (2012) with extended CH₄ observations from ZOTTO. To identify the drivers of the mid- to late-summer CH₄ spikes at ZOTTO, the study progresses from seasonal-scale patterns to a detailed examination of diurnal variability and local processes. This approach maximises the scientific value of the ZOTTO observations.

The diurnal variations in atmospheric CH₄ mole fraction are determined by interactions between surface CH₄ emissions originating from wetlands, agriculture, and fossil fuel combustion (Metya et al., 2021) and atmospheric boundary layer (ABL) dynamics including daytime mixing through entrainment and nighttime stratification that traps near-surface emissions. Studies have shown that small-scale ABL processes, such as entrainment from the free troposphere and the daily evolution of boundary layer depth, can strongly influence observed tracer concentrations (Denning et al., 1996; Larson and Volkmer, 2008; Pino et al., 2012; Williams et al., 2011; Schuh and Jacobson, 2023; Faassen et al., 2025). More accurate representations of these processes are essential when interpreting CH₄ observations at diurnal scales or when applying them in high-spatiotemporal resolution models (Yi et al., 2004; Kretschmer et al., 2014; Bonan et al., 2024).

Despite its importance, long-term trends and drivers of CH₄ diurnal cycle have received little attention compared to seasonal or annual CH₄ changes. Previous studies of diurnal CH₄ mole fractions from observational towers have typically focused on characterizing the patterns of the diurnal cycle (e.g., Metya et al., 2021; Mahata et al., 2017) rather than examining long-term changes in these patterns. Moreover, they have not systematically quantified the relative contributions of potential drivers, such as the extent to which observed variations in CH₄ diurnal cycle are controlled by surface processes versus atmospheric dynamics. Other studies have mainly investigated short-term diurnal CH₄ fluxes at the leaf or ecosystem scale under laboratory or field conditions (e.g. Takahashi et al., 2022; Kohl et al., 2023) or have analysed long-term CH₄ fluxes derived from process-based models (Duan et al., 2025).

The goal of this study is, therefore, to investigate long-term CH₄ variability at ZOTTO from 2010 to 2021 across interannual, seasonal, and diurnal timescales. We first analyse interannual and seasonal CH₄ patterns, focusing on the late-summer peak period. We then zoom into the local scale and examine trends and interannual variability in the CH₄ mole fraction diurnal cycle during July–October, when these peaks are most pronounced. Specifically, we quantify changes in CH₄ diurnal amplitude, defined as the difference between daily maximum and minimum mole fractions, to assess how the diurnal cycle has evolved over time. Finally, we combine six-level CH₄ mole fraction and meteorological observed profiles at ZOTTO to calculate local CH₄ fluxes follow the approach of Winderlich (2012), to determine whether the observed trends and variability in the late summer CH₄ diurnal cycle are primarily driven by changes in surface fluxes or by atmospheric boundary-layer processes.

The remainder of this study is structured as follows: Section 2 provides an overview of the ZOTTO site, including its CH₄ and meteorological measurements, and details of the dataset and fundamental concepts employed to examine long-term trends, as well as seasonal and diurnal variations in CH₄. Sections 3 and 4 present the results and discuss the findings, while Sect. 5 summarises the conclusions and acknowledges the limitations of this work.

2 Methods

2.1 ZOTTO Site Description

The 304 m tall ZOTTO tower is located at 60°48^′ N, 89°21^′ E (114 m a.s.l.), approximately 20 km west of the village of Zotino at the Yenisei River. The surrounding area is characterized by gentle hills of 60–130 m a.s.l. covered with light taiga forests (Pinus sylvestris dominated) on lichen-covered sandy soils (Schulze et al., 2002), interspersed by numerous waterlogged old river meanders and bogs. The approximate tree height around the tall tower is 20 m. The nearest airport is 90 km north in Bor (2600 inhabitants), the closest cities are Yeniseysk and Lesosibirsk (20 000 and 61 000 inhabitants) to the south-southeast, more than 300 km away, and Krasnoyarsk (1 million inhabitants) about 600 km south of ZOTTO (Heimann et al., 2014; Kozlova et al., 2008).

The climatic conditions at ZOTTO are characterised by a mean annual temperature of 3.8 °C measuring at 52 m and total annual rainfall of 536 mm measured at the station. There had been a consistent upward trend in temperatures during summer over the 2010–2021 period (Appendix A – Fig. A1). Rainfall was lowest in winter, and peaks in July and August (Appendix A – Fig. A2b). During the 2010–2021 period, climatic anomalies include summer 2012 and 2016, characterised by warmer-than-average temperatures and with reduced precipitation (Appendix A – Fig. A2). Wind patterns at ZOTTO are predominantly south-westerly and westerly reflecting the prevailing regional dynamics (Appendix A – Fig. A3). Winter at ZOTTO in Siberia is also characterised by the presence of the Siberia High (Winderlich, 2012), a persistent high-pressure system that leads to strong temperature inversions, low wind speeds, and limited vertical mixing during the winter in the artic regions (Serreze et al., 1992).

2.2 Data Description

2.2.1 Methane Mole Fraction Observations

Continuous monitoring of atmospheric CH₄ has been conducted at the tall ZOTTO tower since April 2009 (Winderlich et al., 2010). Air is sampled from six inlets located at heights of 4, 52, 92, 156, 227, and 301 m above ground level. CH₄ mole fractions at these heights are measured using an EnviroSense 3000i gas analyser (Picarro Inc., USA) employing the Cavity Ring-Down Spectroscopy (CRDS) technique.

Data were recorded every 30 s from each active sampling line. Each of the six tower levels was sampled for 3 min, discarding the first 1.5 min for stability. Measurements were taken sequentially in an 18 min cycle from the top level. Since only a single gas analyser was available, each sampling line was connected to an 8 L buffer sphere for continuous, synchronised sampling of close-to-the-same air mass at all heights (Winderlich et al., 2010). While one line was analysed, the others were continuously flushed at 150 sccm under 700 mbar. The buffer system integrated data over 37 min to smooth short-term fluctuations, ensuring stable, well-mixed air samples for reliable measurements.

Calibration is achieved using four horizontally stored aluminium tanks. The CH₄ mole fraction in the gas tanks that are currently used at ZOTTO (see also Tran et al., 2024 (supplement, table A1)) were determined in the GasLab of the MPI-BGC Jena and are traceable to scales of the World Meteorological Organization (WMO) maintained in NOAA/ESRL (WMO X2004A for CH₄; Dlugokencky et al., 2005). To monitor measurement accuracy, a target tank is sampled every 200 h for 8 min, interspersed randomly between two calibration cycles. These target measurements are processed in the same manner as ambient air data. Following the calibration procedure described in Winderlich et al. (2010), the CH₄ mole fraction in the target tank was stable at 1946.5 ± 0.2 ppb for the entire period (Appendix B – Fig. B1). A comparison with measurements of the same tank from the Jena GasLab (1946.4 ± 1.4 ppb) revealed a statistically insignificant bias and no discernible long-term stability issue in the measurements. This consistency confirms the reliability of the mole fraction measurements, ensuring their suitability for further analysis.

We also compared our analysis at ZOTTO site from 2010 to 2022 to the Marine Boundary Layer (MBL) product at 60° N (NOAA GML, 2025). The MBL reference dataset represents well-mixed atmospheric conditions at the same latitude as ZOTTO and is representative for locations situated far from significant anthropogenic and natural sources and sinks.

2.2.2 Meteorological Data

The meteorological measurement system at ZOTTO has been operational since 2007, with its design and functionality described in detail by Winderlich et al. (2010). This section highlights the features relevant to this study.

Wind measurements (in m s⁻¹) are conducted using six 3D sonic anemometers, with 10 Hz frequency, mounted at the same 6 heights as the mole fraction measurements of CH₄ along the tower, and recorded every 30 min, supplemented by air temperature (in °C) and humidity (in %) sensors at all (4, 52, 92, 156, 227, and 301 m a.g.l.) levels and pressure (in hPa) transducers at three heights (4, 92, and 301 m). Pressure values for levels without direct measurements are linearly interpolated. Sensible heat flux (in W m⁻²) is measured with eddy covariance system at all heights. Additional calculated variables include potential temperature (in K), specific humidity (in g kg⁻¹), and vapor pressure deficit (VPD, in kPa), providing a detailed vertical profile of atmospheric conditions.

Vertical profiles of soil temperature (°C) (measured at depths of 2, 4, 6, 8, 16, 32, 64, and 128 cm), soil moisture (%) (at 8, 16, 32, 64, and 128 cm), and precipitation are also recorded every 30 min. These measurements are taken approximately 100 m southeast of the tower at a site within a densely wooded area characterized by sandy soil covered with lichens, representative of typical forest conditions in the region. For this study, we focus on measurements at 32 cm below the ground, which is the depth capturing soil conditions that are relevant to nutrient cycling and microbial processes (Schimel, 1995), both of which influence CH₄ exchange. Snow depth is not directly measured at the ZOTTO site. For this, we use ERA5 reanalysis data (Hersbach et al., 2023b) selecting the grid point closest to the ZOTTO coordinates, located at 60°75^′ N, 89°25^′ E.

To characterise the large-scale atmospheric circulation patterns as well as the atmospheric conditions above the tower height at ZOTTO, we also use ERA5 reanalysis data for variables that are not directly measured at this site, namely daytime boundary layer height (m) (Hersbach et al., 2023b), and horizontal divergence of the wind velocity (s⁻¹) (Hersbach et al., 2023a), extracted from the same ERA5 grid point as above.

2.3 Methane Long-term Trends and Seasonal Signal Processing at ZOTTO

For annual growth rate and seasonal analyses, we used daytime-averaged CH₄ mole fraction measurements (13:00–17:00 LT) collected at a height of 301 m from the ZOTTO tall tower. This time window was selected to ensure sampling of well-mixed boundary layer air, making the data suitable for investigating long-term trends and seasonal variability.

CH₄ mole fractions at continental sites, particularly at locations like ZOTTO, where multiple local sources and sinks influence CH₄ levels, often exhibit an asymmetric annual distribution, characterised by large positive outliers associated with episodic local and regional emission events. To derive a representative background long-term trend, we first aggregated the daytime-averaged CH₄ data into monthly bins and selected the lowest 10 %–30 % of values within each bin. This filtering approach minimises the impact of extreme events, such as the intense wildfire season of 2012, that could otherwise bias estimates of the annual growth rate. The filtered ZOTTO background daytime CH₄ dataset (ZOT_bg) was then processed using the CCGCRV curve-fitting method (Thoning et al., 1989). We employed the Python implementation of CCGCRV, available as a standalone tool from the NOAA CMDL FTP server (https://gml.noaa.gov/aftp/user/thoning/ccgcrv/; last access: 10 January 2025). The curve-fitting configuration included three polynomial terms to capture the long-term trend and four harmonics to represent the seasonal cycle. Cut-off frequencies were set to the default values of 667 d for the long-term component and 80 d for short-term variability. Any data points lying outside 3 times the normalised root mean square deviation from the CCGCRV-derived smoothed curve were iteratively removed (Kozlova et al., 2008). This process eliminated 0.6 % of the ZOT_bg data, ensuring a more accurate representation of the ZOTTO CH₄ long-term trend.

For seasonal cycle analysis, we first detrended the ZOT_bg dataset by subtracting it from the derived long-term trend. Seasonal amplitude was then calculated as the difference between the maximum and minimum of the monthly medians of the detrended ZOT_bg. This seasonal amplitude is used in this study as a diagnostic metric to investigate changes in CH₄ seasonality at ZOTTO.

2.4 Methane Diurnal Signal Conceptual Framework

A common metric used to analyse shifts in the diurnal cycle of a tracer is its diurnal amplitude (Yi et al., 2004; Kretschmer et al., 2014; Bonan et al., 2024). The diurnal amplitude of CH₄ is typically defined as the difference between its maximum nighttime mole fraction (CH_4,max) and minimum daytime well-mixed mole fraction (CH_4,min). The diurnal cycle of atmospheric CH₄ is driven by surface CH₄ sources and sinks, and is modulated by atmospheric boundary layer dynamics, which include enhanced daytime mixing and nighttime stability (Fig. 1). Therefore, the diurnal variability of CH₄ mole fraction results from a combined effect of surface fluxes and atmospheric processes.

Figure 1 Schematic overview of diurnal cycle of the mole fractions of atmospheric CH₄ from the top of a forest canopy (z=ct) to the top of the Atmospheric Boundary Layer (ABL) (z=h) illustrating Eq. (1), adapted from Faassen et al. (2024). The figure illustrates the ABL is characterised by a Convective Boundary Layer (CBL) during daytime and a Nocturnal Boundary Layer (NBL) formation during nighttime. ϕ denotes the CH₄ mole fraction, and the overbar ($\stackrel{\mathrm{‾}}{\mathit{\varphi }}$) represents a 30 min time average. The prime symbol (^′) indicates the deviation of the instantaneous CH₄ mole fraction from its time average. Similarly, w^′ represents the deviation of the instantaneous vertical wind speed from its time-averaged mean ($\stackrel{\mathrm{‾}}{w}$). The “Net Surface flux of CH₄” term (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{NetSurf}}$) refers to the turbulent fluxes from the vegetation layer, up to the top of the canopy (z=ct). The fluxes up to this level depend on terrestrial processes, which contribute to the CH₄ mole fractions observed above the top of the canopy. Entrainment flux at the top of the ABL (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}$) represents the exchange of CH₄ air originated from the residual or free troposphere, into the ABL, and vice versa. The horizontal advection of CH₄ (adv(ϕ)), and the chemical reaction term (S_ϕ) in Eq. (1) are not included in this figure since they are not currently accounted for in this study.

To integrate all the surface and boundary layer dynamics governing the diurnal variations of CH₄ mole fractions, we make use of a time-dependent equation (Eq. 1) inspired on the mixed-layer equations as in Vilà-Guerau de Arellano et al. (2015). Figures 1 and 2 summarises the surface and atmospheric drivers of CH₄ diurnal mole fractions and hence, its diurnal amplitude.

$$\begin{array}{}\text{(1)}& \begin{array}{rl}& {\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}{\int }_{z=\text{ct}}^{z=h}{\displaystyle \frac{\partial \mathit{\varphi }}{\partial t}}\left(z\right)\mathrm{d}z=\\ & \underset{I}{\underbrace{{\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}}}\times \left(\underset{\text{II}}{\underbrace{\stackrel{\mathrm{‾}}{(w^{\prime }{\mathit{\varphi }}^{\prime }{)}_{\text{NetSurf}}}}}-\underset{\text{III}}{\underbrace{\stackrel{\mathrm{‾}}{(w^{\prime }{\mathit{\varphi }}^{\prime }{)}_{\mathrm{e}}}}}\right)-\underset{\text{IV}}{\underbrace{\text{adv}\left(\mathit{\varphi }\right)}}+\underset{\text{V}}{\underbrace{S_{\mathit{\varphi }}}}\end{array}\end{array}$$

In this context, the ϕ refer to CH₄ mole fraction (nmol mol⁻¹), w to vertical wind speed (m s⁻¹), h to the top of the atmospheric boundary layer (ABL) (m a.g.l.), ct to the top of the canopy (m a.g.l.), t to timestep of 30 min (in seconds), adv(ϕ) to horizontal advection (ppb s⁻¹), and S_ϕ to the combination of production and loss of CH₄ from chemical reactions with OH (ppb s⁻¹). The overlines ($\stackrel{\mathrm{‾}}{{}^{}}$) refer to 30 min time averaged values, and the prime (^′) representing the deviations from the mean. All symbols and their corresponding units in this study are provided in Appendix C – Table C1. Equation (1) shows that rate of change in the column-integrated CH₄ mole fraction (left hand side) (i.e. from the top of the canopy (ct) to the top of the atmospheric boundary layer (ABL) (h)) (Fig. 1) depends on:

• I.

  The thickness of the ABL (h−ct) (Fig. 1 and Eq. 1) (in m): The height of the ABL (z=h) exhibits a pronounced diurnal cycle. During the day, the ABL is referred to as the Convective Boundary Layer (CBL), while at night it is termed as the Nocturnal Boundary Layer (NBL) in this study. The thickness of this layer affects the dilution of CH₄ during the day and the accumulation of CH₄ at night (see Sect. D1.1 in Appendix D).

• II.

  The Net Surface flux of CH₄ (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{NetSurf}}$) (ppb m s⁻¹) represents the balance between sources and sinks at the top of the canopy (z=ct) (Fig. 1 and Eq. 1): This flux captures small-scale processes within the canopy, such as CH₄ exchange from vegetation and soil, contributing to the CH₄ mole fractions in the ABL.

• III.

  Entrainment flux at the top of the ABL (z=h) (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}$) (ppb m s⁻¹) (Fig. 1 and Eq. 1): This flux represents the mixing of CH₄ air from above the ABL, either the residual layer during the night or free troposphere during the day, to inside the ABL (see Sect. D1.2 in Appendix D).

• IV.

  The horizontal advection of CH₄ (adv(ϕ)) (ppb s⁻¹) in Eq. (1), representing the meso- and long-range horizontal transport of CH₄ which is currently not accounted for this study similar to Winderlich et al. (2014).

• V.

  The combination of production and loss of CH₄ from chemical reactions with OH (S_ϕ) (ppb s⁻¹), which assumed to be negligible within the diurnal scale due to the slow reaction rate of CH₄ with OH compared to the atmospheric residence time of OH (Patra et al., 2009) and therefore will be omitted from Eq. (1).

Figure 2 Chart illustrating the cause-effect relationships governing the diurnal amplitude of CH₄. Solid arrows indicate positive influences, while dashed arrows represent negative influences. Red and blue highlight key nighttime and daytime processes, respectively, while black denotes processes affecting both daytime and nighttime dynamics. Yellow and pink boxes distinguish atmospheric and surface canopy processes, respectively. Grey text within some boxes references sections and figure where further method details and their corresponding results can be found. Boxes with dashed outline are the variables investigated in this study. The figure illustrates the ABL is characterised by a Convective Boundary Layer (CBL) during daytime and a Nocturnal Boundary Layer (NBL) formation during nighttime. The “Net Surface flux of CH₄” term (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{NetSurf}}$) refers to the fluxes from the vegetation layer, up to the top of the canopy; integrates both turbulent flux (F_Eddy) and storage flux (F_Stor). Entrainment flux at the top of the ABL (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}$) represents the mixing of CH₄ air from above the ABL, to inside the ABL. Δ_{ϕ(ft/rs−cbl/nbl)}: Change in CH₄mole fraction across layers (free troposphere/residual layer vs. CBL/NBL); ${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}>\mathrm{0}/<\mathrm{0}$: Positive/negative surface buoyancy flux, indicating convective/stable boundary layers; $\text{H \& LE}>\mathrm{0}/<\mathrm{0}$: Sensible heat (H) and latent heat (LE) fluxes, with signs indicating net heating or cooling; w_sub(h) represent the mean vertical subsidence velocity at ABL height and $\text{Div}\left(\stackrel{\rightharpoonup }{U_h}\right)$ is the horizontal wind divergence at the ABL height.

Terms I and III in Eq. (1) describe key atmospheric boundary layer dynamics that influence CH₄ variations from the canopy top to the ABL, primarily driven by ABL height and entrainment flux at the ABL top. A more detailed derivation of Eq. (1), along with explanations of each atmospheric driver terms in Fig. 2 (yellow-coloured boxes), is provided in Appendix D1. The key concepts driving the diurnal amplitude of CH₄ in this study (Eq. 1, excluding Terms IV and V) are assumed to apply under high-pressure atmospheric circulation systems. High-pressure systems are generally associated with subsidence and stable, calm weather, which limits horizontal advection. This assumption is important, as horizontal advection (Term IV in Eq. 1) is not explicitly included in our analysis.

2.5 Methane Diurnal Signal Processing at ZOTTO

To calculate the diurnal amplitude of CH₄ at ZOTTO, we used hourly CH₄ mole fraction measurements at 52 m, which is the closed available measurement height above the forest canopy. Measurements at this level are well-suited to capturing short-term diurnal variations, where surface-atmosphere exchange processes are most active and pronounced.

To investigate the drivers behind shifts in CH₄ diurnal amplitude over the study period at ZOTTO, our analysis focused exclusively on days influenced by high-pressure systems. High-pressure conditions were identified using ERA5 geopotential height data at 550 hPa (Hersbach et al., 2023a), selecting periods where geopotential height exceeded the 90th percentile of the distribution of each season following the approach of Marín et al. (2022).

2.5.1 Atmospheric Boundary Layer Height

The ABL is distinguished between the daytime CBL and the nighttime NBL (Figs. 1 and 2).

To identify the CBL height, we used ERA5 reanalysis data (Hersbach et al., 2023b). Daily values were averaged between 12:00 and 15:00 LT to capture the peak convective period while avoiding the sunrise and sunset transitional periods. By examining interannual variations in the summer CBL height, we assess whether daytime dilution effects on CH₄ mole fractions have strengthened or weakened over time.

To determine the NBL height, we utilised the in-situ vertical potential temperature (θ) observations from our 300 m measurement tower, as ERA5 reanalysis data have been observed to overestimated nighttime conditions (Sinclair et al., 2022). We applied a least-squares regression fit to the nighttime vertical θ gradients using the following equation (Oncley et al., 1996; Frenzen and Vogel, 2001; Johansson et al., 2001):

$$\begin{array}{}\text{(2)}& X\left(S\right)=A+B\times \ln \left(S\right)+C\times \ln (S{)}^{\mathrm{2}}\end{array}$$

Where X(S) is the fitted function of the variable S, which in this case represents θ (in K). For nighttime data at ZOTTO, this fitting was applied to the vertical profile from 4 heights above the boreal canopy (52, 92, 156 and 227 m) (red line in the Appendix E – Fig. E1). Data from 4 m (within the canopy) and 301 m (residing in the residual layer during nighttime (see explanation in Appendix A – Fig. A4)) were excluded. Fits with an R² value greater than 0.7 were retained. This process eliminated 14.5 % of the vertical potential temperature data at ZOTTO during the July to September months for the 2010–2021 period. The first derivative of the fitted curve X(S), representing the temperature lapse rate ($\frac{\partial \mathit{\theta }}{\partial z}$), was computed and normalised (X_norm) as in Eq. (3):

$$\begin{array}{}\text{(3)}& X_{\text{norm}}={\displaystyle \frac{X^{\prime }\left(S\right)}{X_{\text{max}}}}\end{array}$$

where X^′(S) (K m⁻¹) is the first derivative of X(S) and X_max is the maximum value of X^′(S). The NBL height was identified as the altitude where the X_norm curve (blue line Appendix E – Fig. E1) decayed to its e-folding value ($\sim \mathrm{1}/e$ of the maximum) (Stull, 1988). We analysed temporal variations in these derived NBL heights, restricting the dataset to timestamps between 00:00–04:00 LT, period when the vertical θ profile at ZOTTO is the most pronounced.

Since both the vertical profiles of θ and CH₄ provide insight into the structure and evolution of the nighttime atmospheric column, with the CH₄ profile typically mirroring the temperature profile in the opposite direction (Appendix A – Fig. A4 and Fig. 3), we applied the same regression methodology to the CH₄ vertical profiles (Eqs. 2 and 3) to compare the trends in NBL height derived from both variables over time. The NBL heights derived from CH₄ and θ vertical profiles exhibit similarities (see Appendix E), further validating the reliability of the least-squares regression fit approach.

Figure 3 The climatological late summer (JASO) vertical profile of CH₄ mole fractions every 2 h at ZOTTO (2010–2021 average). The shaded area represents the calculated storage flux (F_Stor) between consecutive time steps and tower heights. The dashed green line denotes the canopy height at ZOTTO, while the dashed brown and grey lines indicate the 52 and 72 m levels, where the F_Stor and F_Eddy terms in Eq. (4) are calculated, respectively.

By examining interannual variations in the NBL height, where a shallower NBL normally related to stronger thermodynamic stable stratification and a higher NBL indicates weaker stability, we can assess whether the nighttime stability leading to accumulation of near-surface CH₄ mole fractions have strengthened or weakened over time.

Entrainment Effect

As entrainment is negligible during the nighttime (see Sect. D1.2 in Appendix D), we focus on quantifying the daytime entrainment rate over the 2010–2021 period. We examined trends in CBL growth ($\frac{\partial h}{\partial t}$) (m s⁻¹), subsidence velocity (w_sub(h)) (m s⁻¹), and CH₄ mole fraction differences between the CBL and its overlaying layer – the free troposphere (FT) (Δ_ϕ(ft−cbl)) (nmol mol⁻¹). As shown in Fig. 2, these factors collectively contribute to entrainment strength. Given that Δ_ϕ(ft−cbl) is influenced by the stability of the previous night (Fig. 2 and See Sect. D1.2.3 in Appendix D), which is already assessed as NBL height as in the Atmospheric Boundary Layer Height section above, we primarily focused on CBL growth and subsidence velocity to assess the daytime entrainment effects.

We used positive sensible heat flux, a proxy for buoyancy turbulence, to assess whether mixing strength and CBL growth rate have changed over time. Higher values indicate increased turbulence, enhancing CBL growth and entrainment. Sensible heat flux data from 52 m at ZOTTO (where flux measurements peak just above the canopy) were compared with ERA5 surface sensible heat flux values to determine temporal trends in daytime CBL growth.

In this study, we analysed long-term variations in subsidence by examining positive divergence data from ERA5 at the 750 hPa pressure level over the course of 2010–2021 period. The 750 hPa level resides in the mid-troposphere, where large-scale vertical motions are the most prominent (Stepanyuk et al., 2017). This level, therefore, effectively captures the dynamics of vertical air movement and its influence on atmospheric stability.

2.5.2 Surface Processes

This section gives an overview of how to assess net CH₄ surface fluxes (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{NetSurf}}$) at the top of the canopy (z=ct) (Fig. 1 and pink-coloured boxes in Fig. 2) using the vertical CH₄ profile at ZOTTO following the approach of Winderlich (2012). This can be calculated using the following terms in Eq. (4) adapting from Finnigan (1999), Yi et al. (2000), and Feigenwinter et al. (2004):

$$\begin{array}{}\text{(4)}& {\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{NetSurf}}={\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{ct}}+{\int }_{z=\mathrm{0}}^{z=\text{ct}}{\displaystyle \frac{\partial \stackrel{\mathrm{‾}}{\mathit{\varphi }}}{\partial t}}\mathrm{d}z\end{array}$$

The sign convention used here gives positive ${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{NetSurf}}$ for ecosystem emissions, where a positive flux term (i.e., source) corresponds to transport out of the control volume (Feigenwinter et al., 2004). In Eq. (4) on the right hand side, the first term represents the turbulent vertical flux (F_Eddy) measured at the top of the canopy (ct) (ppb m s⁻¹), while the second term represents the storage of CH₄ (F_Stor) (ppb m s⁻¹), which is the temporal dynamics of CH₄ in the air column below the F_Eddy measurement height, not influenced by turbulence, as calculated by the integral. Each of these two terms is discussed in detail in the subsequent subsections.

Inferred Turbulent Vertical Flux

Direct eddy flux measurements of CH₄ are not feasible at ZOTTO tall tower, due to the low measurement frequency of the CH₄ analyser (0.2 Hz), long tubing lengths (up to 320 m), and the use of buffer volumes with extended mixing times (∼40 min, corresponding to ∼0.0004 Hz). Instead, we apply the modified Bowen ratio approach (Businger, 1986; Winderlich et al., 2014), in which the turbulent vertical flux at the top of the canopy (F_Eddy) in Eq. (4) can be written as:

$$\begin{array}{}\text{(5)}& F_{\text{Eddy}}={\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\text{ct}}\cong {\displaystyle \frac{H_{\left(\text{ct}\right)}}{{\mathit{\rho }}_{\left(\text{ct}\right)}\cdot C_p}}{\displaystyle \frac{\partial \mathit{\varphi }/\partial z}{\partial \mathit{\theta }/\partial z}}\end{array}$$

where H_(ct) and ρ_(ct) are the sensible heat flux (W m⁻²) and air density (g m⁻³) at the top of the canopy. Since direct measurements at the exact canopy height at ZOTTO (∼28 m) are unavailable, we use the sensible heat flux measured using eddy-covariance system and the density at 52 m, which is the closest available measurement height above the canopy. C_p is specific heat capacity constant (C_p=1.00467 $\mathrm{J}{\mathrm{g}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$). The mole fraction and temperature gradients between two adjacent heights (52 and 92 m) are used to compute the turbulent fluxes at the intermediate level 72 m. This derived F_Eddy at 72 m represents the turbulent flux at the top of the canopy.

Storage Flux

The storage term F_Stor in Eq. (4) represents the temporal changes in column-integrated CH₄ mole fraction below 72 m, the height at which F_Eddy is calculated in the above section. For illustration, the diurnal development of the CH₄ profile along the tower is given in Fig. 3. The F_Stor can be visualised through the shaded trapezoidal areas between the half-hourly time steps t_i and t_i+1, and the two different tower heights below 72 m (i.e., z₁=4 m and z₂=52 m) (see Fig. 3, grey-shaded area). The F_Stor term in Eq. (4) can be expanded as in Winderlich et al. (2014):

$$\begin{array}{}\text{(6)}& \begin{array}{rl}{F}_{\text{Stor}}& (t_i,z_{\text{ct}})={\int }_{z=\mathrm{0}}^{z=\text{ct}}{\displaystyle \frac{\partial \stackrel{\mathrm{‾}}{\mathit{\varphi }}}{\partial t}}\mathrm{d}z\\ \cong & {\displaystyle \frac{\frac{\mathrm{1}}{\mathrm{2}}\left(\right({\mathit{\varphi }}_{\mathrm{1}}\left(t_{i+\mathrm{1}}\right)-{\mathit{\varphi }}_{\mathrm{1}}\left(t_i\right))+({\mathit{\varphi }}_{\mathrm{2}}\left(t_{i+\mathrm{1}}\right)-{\mathit{\varphi }}_{\mathrm{2}}\left(t_i\right)\left)\right)}{t_{i+\mathrm{1}}-t_i}}\\ & \times (z_{\mathrm{1}}-z_{\mathrm{2}})\end{array}\end{array}$$

The ϕ₁(t_i) and ϕ₂(t_i) represent CH₄ mole fraction at 4 and 52 m at the time step t_i, respectively. Storage fluxes are calculated up to 52 m in this study, which is the highest available measurement below the F_Eddy estmation at 72 m.

During high-pressure systems, the downward subsidence velocity can lower the effective height at which the storage flux is calculated up to. This downward movement of air means that the storage flux calculations need to be adjusted to reflect this displacement. Equation (6), can be expanded to below:

$$\begin{array}{}\text{(7)}& \begin{array}{rl}{F}_{\text{Stor}}(t_i,z_{\text{ct}})\cong & {\displaystyle \frac{\frac{\mathrm{1}}{\mathrm{2}}\left(\right({\mathit{\varphi }}_{\mathrm{1}}\left(t_{i+\mathrm{1}}\right)-{\mathit{\varphi }}_{\mathrm{1}}\left(t_i\right))+({\mathit{\varphi }}_{\mathrm{2}}\left(t_{i+\mathrm{1}}\right)-{\mathit{\varphi }}_{\mathrm{2}}\left(t_i\right)\left)\right)}{t_{i+\mathrm{1}}-t_i}}\\ & \times (z_{\mathrm{1}}-z_{\mathrm{2}})\times {\displaystyle \frac{z_{\mathrm{2}}-w_{\text{sub}\left(z_{\mathrm{2}}\right)}\cdot (t_{i+\mathrm{1}}-t_i)}{z_{\mathrm{2}}}}\end{array}\end{array}$$

Where $w_{\text{sub}\left(z_{\mathrm{2}}\right)}$ is the vertical wind component at 52 m. We do not derive this vertical wind component from 3D anemometer measurements at ZOTTO due to its sensitivity to sensor misalignment, where parts of the horizontal wind components are inadvertently reallocated into small parts of the vertical component, causing significant errors (Winderlich, 2012). To address this challenge, horizontal divergence from 3-hourly short-term forecast fields from the operational archive of the European Centre for Medium-Range Weather Forecasts (ECMWF, http://www.ecmwf.int/, last access: 18 November 2025) has been used to derive vertical mean wind speed $w_{\text{sub}\left(z_{\mathrm{2}}\right)}$ at 52 m at ZOTTO.

In summary, the net surface flux at the top of the canopy (in Eq. 4) will be represented as the sum of F_Stor up to 52 m and F_Eddy at 72 m in this study. We focus the calculation of this net surface flux on the nighttime period (00:00–04:00 LT). This is because the Bowen ratio method used to estimate the turbulent flux (F_Eddy in Eq. 5), which depends on vertical gradients of CH₄ and θ, becomes less reliable during the day. Strong daytime mixing reduces these vertical gradients, introducing significant noise into F_Eddy and compromising the accuracy of the total net surface flux. By restricting our analysis to the nighttime period (00:00–04:00 LT), when the vertical gradients are the most pronounced, we ensure more reliable F_Eddy and hence net surface flux estimates.

During nighttime at ZOTTO, the F_Stor becomes the dominant component of the total net surface flux, contributing approximately 60 %–75 % (Appendix F – Fig. F1), surpassing F_Eddy. A similar pattern was previously reported for CO₂ net surface flux at Missouri Ozark by Yang et al. (2007) and at ZOTTO by Winderlich (2012).

2.6 Statistics

To analyse trends in variables over the 2010–2021 period, we applied the Theil-Sen regression method (Theil, 1992; Sen, 1968), a robust non-parametric approach known for its resistance to outliers. To investigate potential drivers of any observed significant trends, we performed orthogonal regression, which accounts for errors in both dependent and independent variables, providing a more reliable assessment of relationships between variables.

3 Results

3.1 Long-term Trend and Growth Rate

The CH₄ mole fraction recorded at ZOTTO is generally higher than the MBL by around 50 ppb (Fig. 4a). ZOTTO also displays a slightly more pronounced increase in CH₄ levels as its trend diverges from MBL over time (Fig. 4b).

Figure 4 (a) Daytime CH₄ mole fractions at 301 m a.g.l. from the ZOTTO tower: 13:00–17:00 LT averaged daytime data (ZOT, grey) and filtered-background daytime data (ZOT_bg, red), shown alongside bi-weekly CH₄ measurements from the marine boundary layer at 60° N (MBL, blue); (b) Long-term CH₄ trends of ZOT_bg and MBL data derived from the CCGCRV curve-fitting method (Thoning et al., 1989), normalized to their respective 2010 baseline values.

The average annual CH₄ growth rate at ZOTTO from 2010 to 2021 is slightly higher (9.85 ± 7.1 ppb yr⁻¹) than MBL (9.05 ± 5.5 ppb yr⁻¹). The growth rate of ZOTTO also reflects more interannual variability compared to the MBL, indicated by the larger standard deviations (Table 1) due to local and regional sources. The annual growth rates in Table 1 show that the CH₄ growth rate at ZOTTO peaked in 2014 at 21.22 ppb yr⁻¹, and in 2016 at 12.22 ppb yr⁻¹ followed by an acceleration in 2019 to 2021. The MBL growth rates show similar temporal patterns, with notable increases in 2014, 2016, and 2019–2021.

Table 1 Annual growth rate of ZOTTO, and MBL in ppb yr⁻¹ derived from the first derivative of the trendlines in Fig. 4b (without normalising to their respective 2010 baseline values).

^∗ Years indicate strong growth rate.

3.2 Seasonal Cycle

There is a consistent shape of the seasonal cycle across the two datasets, with higher mole fractions in the colder months (winter) and lower mole fractions in the warmer months (spring and summer) (Fig. 5). At ZOTTO, the seasonal fluctuations of CH₄ are more pronounced compared to the MBL data, with several clear peaks during the cold months. While the shape of the seasonal cycle at ZOTTO shows variability from year to year, it remains relatively constant at MBL across the years. A slight time lag is observed between the seasonal CH₄ minima at ZOTTO and in the marine boundary layer (MBL). There is a shift in the seasonal cycle phase between the two datasets, with the MBL phase occurring one or 2 months later than that of ZOTTO. This lag likely results from the atmospheric transport of CH₄-depleted air masses originating over the continents and moving toward the ocean. A similar pattern has been reported in the ZOTTO CO₂ record (Tran et al., 2024) when compared with the MBL reference data, supporting the influence of large-scale air mass transport on the observed seasonal timing.

Figure 5 The yearly seasonal cycles of background-filtered daytime CH₄ at ZOTTO (ZOT_bg) and biweekly marine boundary layer CH₄ at 60° N (MBL), shown after removing long-term trends (i.e., subtracting the trend components presented in Fig. 4b from the data in Fig. 4a). The line plots with circle markers represent the monthly medians. The darker shaded boxes indicate the interquartile range (IQR). The lighter shaded boxes extend from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{1.5}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{1.5}\times \text{IQR}$, where Q₁ and Q₃ are the 25th and 75th percentiles, respectively. Coloured ticks on the x axis highlight: dark red (ZOT) and dark blue (MBL) for winter maxima between December of the previous year and March of the current year (star indicates the maximum belong to the next year), light red (ZOT) and light blue (MBL) for spring minima, and yellow for the late summer peak at ZOT.

The seasonal amplitude of CH₄ at ZOTTO is consistently larger than that observed in the MBL dataset. The seasonal amplitude is calculated as the difference between the winter maximum and the spring minimum median values (in Fig. 5). Although both datasets exhibit increasing trends in seasonal amplitude over the period 2010–2021 (2.12 ppb yr⁻¹, p=0.12 and 0.49 ppb yr⁻¹, p=0.09 for ZOT and MBL respectively), neither trend is statistically significant (Fig. 6). The increase in the ZOTTO seasonal amplitude (Fig. 6), while not statistically significant, appears to be mainly driven by an increasing winter CH₄ mole fraction maximum (1.42 ppb yr⁻¹, p=0.17) accompanied by a slight decrease in the CH₄ spring minimum (−0.69 ppb yr⁻¹, p=0.32) (Appendix G – Fig. G1).

Figure 6 Time series of the CH₄ seasonal cycle amplitude (square markers) for detrended background-filtered daytime CH₄ at ZOTTO (ZOT_bg) and biweekly marine boundary layer CH₄ at 60° N (MBL). Seasonal amplitude is calculated as the difference between the winter maximum and the seasonal minimum median values derived from Fig. 5. For ZOTTO, the amplitude is defined as the difference between the December–February maximum and the May–July trough; for MBL, it is calculated as the December–February maximum minus the July–August trough. The Theil-Sen regression trend is depicted by the solid line, with the 95 % confidence interval of the trend shown as dashed curves. The p value indicates whether the slope of the regression is significantly different (at 0.05 level) from zero.

Notably, ZOTTO displays a secondary peak in late summer during the late summer (July–October) period, whereas for the MBL dataset this feature is less pronounced, occurring later in the season, concurrent with the winter peak, around November–January (Fig. 5). The amplitude of this late summer peak at ZOTTO, calculated as the difference between the late summer maximum (between July–October) and the seasonal minimum (during May–July period) median values (in Fig. 5), shows a significant increasing trend at 0.05 level (1.35 ppb yr⁻¹, p=0.02) along with notable interannual variability (Fig. 7). The enhanced late-summer peak amplitude (Fig. 7) is primarily attributed to the strong significant increase in the late-summer CH₄ maximum (1.35 ppb yr⁻¹, p=0.02), rather than to a slight decrease in the springtime minimum (−0.69 ppb yr⁻¹, p=0.32) (Appendix G – Fig. G1). To further explore the potential factors that might contribute to this observed increase in the late summer peak unique to ZOTTO, we will focus our analysis of variations and trends in diurnal amplitude over the years during the late summer months (July to October) in the next section.

Figure 7 Time series of the CH₄ seasonal late summer peak amplitude (circle markers) for detrended background-filtered daytime CH₄ at ZOTTO (ZOT_bg). The amplitude here is calculated as the difference between the late summer maximum (between July–October) and the seasonal minimum (during May–July) median values in Fig. 5. The Theil-Sen regression trend is depicted by the solid line, with the 95 % confidence interval of the trend shown shaded areas. The p value indicates whether the slope of the regression is significantly different (at 0.05 level) from zero.

3.3 Diurnal Cycle

CH₄ measurements in the summer from the ZOTTO site exhibit a distinct diurnal cycle, characterized by peak mole fractions around 06:00 LT, followed by a sharp decline around 07:00 LT (Appendix G – Fig. G2). Lower values persist throughout the day until approximately 18:00 LT. Seasonally, the diurnal cycle is most pronounced during the warmer months, i.e. spring, summer, and autumn, while it remains minimal in winter (Appendix G – Fig. G3).

Between 2010 and 2021, the late summer (July to October averaged) diurnal amplitude increased significantly at a rate of 5.29 ppb yr⁻¹ (p=0.001) (Fig. 8a). Both daytime and nighttime CH₄ mole fractions significantly increased at the 0.01 level over this period (Fig. 8b and c), driven largely by the long-term trend observed in Fig. 4c. However, the increase was more pronounced at night (17.10 ppb yr⁻¹, p<0.001) compared to the daytime (11.81 ppb yr⁻¹, p<0.001), emphasising the dominant role of nighttime CH₄ mole fractions in driving the observed rise in summer diurnal amplitude, which may in turns contributed to the statistically significant rise in the late-summer CH₄ peak observed at ZOTTO (Fig. 7). When the long-term anthropogenic trend derived from the full timeseries at the 52 m level (including both daytime and nighttime) is removed, the influence of nighttime CH₄ becomes even more evident (Appendix G – Fig. G4). After detrending, the diurnal amplitude continues to show the same significant increase, which is driven solely by the rise in nighttime CH₄ mole fractions (6.19 ppb yr⁻¹, p=0.008), while daytime mole fractions show no significant trend (Appendix G – Fig. G4). This further suggesting that the increase of the summer diurnal amplitude is primarily due to increased nighttime CH₄ mole fractions.

Figure 8 Time series of yearly late summer (JASO): (a) averaged CH₄ diurnal cycle amplitude; (b) its daytime (10:00–16:00 LT averaged) CH₄ mole fraction, and (c) its nighttime (00:00–04:00 LT averaged) CH₄ mole fraction (right) (circle markers) at ZOTTO using 52 m a.g.l. data. The Theil-Sen regression trend is depicted by the solid line, with the 95 % confidence interval of the trend shown as a dashed line. The p value indicates whether the slope of the regression is significantly different (at 0.05 level) from zero.

To better understand the drivers behind the observed increase trend in the summer CH₄ diurnal amplitude, we focused our analysis on summer days occurring under high-pressure conditions, when the fundamental concepts underlying the potential drivers of CH₄ diurnal amplitude (Presented in Fig. 2) are assumed to hold (see Sect. 2.4). Appendix H – Fig. H1 summarises the number of high-pressure days identified for each summer month over the 11 year period (2010–2021), based on the filtering criteria described in Sect. 2.5. The following analysis is restricted to these high-pressure cases to ensure consistency in atmospheric conditions.

Our results indicate that increasing nighttime CH₄ surface fluxes during summer are the primary driver of the observed rise in nighttime CH₄ mole fractions and the associated increase in summer diurnal amplitude from 2010 to 2021, rather than changes in the NBL dynamics. Orthogonal regression analysis confirms a significant positive relationship between the diurnal amplitude and nighttime net CH₄ surface flux (R²=0.67, p<0.001), while other potential atmospheric dynamic drivers (dashed-outline boxes in Fig. 2) show no such correlation (Fig. 9). At ZOTTO, the inferred nighttime net surface CH₄ fluxes during summer range from −0.05 to 0.2 ppb m s⁻¹, with the highest values occurring in August (Fig. 10). July and August show an increasing trend in nighttime CH₄ flux, with a statistically significant rise in August at 0.05 level (p=0.016). No significant trends are observed in the other atmospheric drivers of the diurnal amplitude (dashed-outline yellow boxes in Fig. 2). A detailed analysis of the interannual variations in these atmospheric drivers over the 2010–2021 period is provided in in Appendix I.

Figure 9 Relationship between July–October averaged summer CH₄ diurnal amplitude and its potential drivers (dashed-outline boxes in Fig. 2): July–October averaged CBL height, NBL height, divergence at 750 hPa, cumulative sensible heat flux at 52 m a.g.l. and nighttime net CH₄ surface flux. Shaded areas represent 95 % confidence intervals.

Figure 10 Box and whisker plot of yearly nighttime (00:00–04:00 LT) Net CH₄ flux for each summer month. The box denotes the interquartile range (IQR), showing the median with a thick black line. The whiskers range from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{1.5}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{1.5}\times \text{IQR}$, where Q₁ and Q₃ are the 25th and 75th percentiles, respectively. The blue line is the monthly mean. The Theil-Sen trend slope for the mean and its p value are denoted on the top left corner of each plot. Numbers above each box indicate the sample size or the number of available days (based on the number of high-pressure days and the availability of vertical profile (4, 52 and 92 m a.g.l.) of CH₄, potential temperature and sensible heat flux at 52 m) for analysis in that month.

Given a clear increasing trend observed in the nighttime net CH₄ surface flux, we further investigate the potential environmental drivers of the increase in this surface flux. The relationship between July–October averaged nighttime net CH₄ surface flux and various environmental variables (Fig. 11) indicates strong positive correlations with July–October averaged soil temperature (R²=0.65, p<0.01), soil moisture at 32 cm below ground (R²=0.36, p=0.032), and with preceding February–May averaged snow depth (R²=0.54, p=0.029). This result suggests that warmer late-summer soil conditions, higher soil moisture and thicker spring snow cover are associated with increased late-summer CH₄ emissions.

Figure 11 Relationship between July–October averaged nighttime net CH₄ Flux and July–October averaged precipitation, soil temperature at the depth of 32 cm, soil moisture at the depth of 32 cm, air temperature measured at 52 m a.g.l., Vapour Pressure Deficit (VPD), and preceding February–May averaged snow depth. Shaded areas represent 95 % confidence intervals. The colour gradient in the data points indicates temporal trends, with more recent years (darker blue) tending toward higher temperatures, lower soil moisture, and increased VPD. The x axis of VPD is reversed to maintain consistency in the direction of positive correlations across the plots, since lower VPD values represent drier conditions.

4 Discussion

4.1 Long-term Trend and Growth Rate

The observed persistent increase in the CH₄ growth rates at ZOTTO are consistent with the global trend reported by NOAA and other long-term monitoring networks (Lan et al., 2021, 2022). Peaks in CH₄ growth at ZOTTO, particularly in 2014, 2016, and 2019–2021, correspond to global trends of increasing CH₄ levels, which began around 2014, steepened after 2018, and further intensified in 2020 (Worden et al., 2017; Nisbet et al., 2019). The global CH₄ mole fraction increases in 2014 and 2020 have been largely attributed to reductions in atmospheric OH radical concentrations (Zhang et al., 2021). In 2020, this OH-driven effect was likely amplified by decreases in anthropogenic NO_x emissions and the associated reduction in free-tropospheric ozone during the COVID-19 lockdowns (Peng et al., 2022). High global growth rates from 2016 to 2020 were additionally driven by strong emissions from boreal wetlands in Eurasia (Yuan et al., 2024; Zhang et al., 2021).

A comparison of CH₄ mole fraction time series between the inland tall tower at ZOTTO and the MBL reference at 60° N shows consistently higher CH₄ concentrations at ZOTTO. However, only in some of the peak years highlighted in Table 1 (i.e., 2014, 2016, and 2020) does ZOTTO exhibit annual growth rates exceeding those of the MBL, suggesting that, in addition to the global baseline increase, ZOTTO may have been influenced by additional regional sources in those years.

4.2 Seasonal Cycle

ZOTTO exhibits a pronounced seasonal cycle in atmospheric CH₄ mole fractions, characterised by maxima during December–January and minima in May–June. This seasonal pattern is consistent with earlier short-term observations at ZOTTO (2009–2012) (Winderlich, 2012). Comparable sinusoidal seasonal cycles have been documented at other monitoring sites, including urban locations such as tall tower Białystok in Poland (Popa et al., 2010), and Guro and Nowon in Seoul, South Korea (Ahmed et al., 2015), as well as remote inland stations such as Fraserdale in Ontario, Canada (Worthy et al., 1998), and Ulaan Uul in Mongolia (Kim et al., 2015).

The seasonal minima in CH₄ mole fractions observed at ZOTTO during May–June period is likely driven by a combination of enhanced OH-driven atmospheric oxidation and hydrological constraints on CH₄ production. OH radicals typically peak in abundance during late spring to early summer, thereby intensifying CH₄ oxidation and contributing to lower atmospheric CH₄ levels during this period. Elevated water table elevation (WTE), defined as the depth below which the soil is saturated, may further suppress CH₄ emissions during the spring. A high WTE in early spring often reflects substantial snow accumulation from the preceding winter, which acts as an insulating layer that limits soil freezing (Granberg et al., 1999). Under such conditions, methanotrophic communities can remain active throughout the winter (Einola et al., 2007; Trotsenko and Khmelenina, 2005), oxidizing CH₄ produced during early thaw and thereby reducing net emissions (Feng et al., 2020). Moreover, because water has a higher heat capacity and latent heat of fusion compared to air, soils under high WTE conditions warm more slowly in spring. This delayed warming postpones the onset of microbial methanogenesis, further limiting CH₄ production during the spring-summer seasonal transition (Feng et al., 2020). Collectively, these processes contribute to the pronounced CH₄ trough observed at ZOTTO in the May–June timeframe.

Wintertime CH₄ peaks at ZOTTO could also be influenced by biomass burning and anthropogenic emissions (Saunois et al., 2025). While fossil fuel activities may contribute to winter CH₄ spikes at ZOTTO, their influence is likely limited due to the distance of the station from major oil and gas production sites (>300 km) (Winderlich, 2012). Former studies suggest that this might be too far to notably modify the CH₄ mole fraction at ZOTTO. Tohjima et al. (1996) discovered CH₄ mole fraction of 2900 ppb in only 150 m altitude above an oil production site, while the signal has already been vented away in 250 m altitude. Moreover, natural gas emissions sum up to 1 % to 10 % of the overall wetland emissions only during the 1999–2003 period (Tarasova et al., 2009) and no significant increase in CO level has been detected during the winter to suggest substantial burning processes (Kozlova et al., 2008). Elevated winter CH₄ levels at ZOTTO are most likely driven by meteorological conditions influenced by the Siberian High, a persistent high-pressure system that leads to strong temperature inversions, low wind speeds, and limited vertical mixing during the winter in the artic regions (Serreze et al., 1992). These conditions trap CH₄ near the surface, contributing to episodic of CH₄ enrichment during the winter (Winderlich, 2012). This winter phenomenon has been witnessed on subcontinental scale in Western Siberia, when CH₄ enrichments above 2000 ppb occurred during high pressure situations in combination with temperature inversions, with temperatures below −20 °C (Sasakawa et al., 2010).

Among the years with high CH₄ growth rates (highlighted in Table 1), the seasonal amplitude at ZOTTO remained relatively small in 2014, 2016, and 2020. During these years, elevated CH₄ mole fractions were distributed relatively uniformly across seasons, resulting in consistent increases throughout the year. This pattern enhanced the annual mean CH₄ mole fractions while leaving the winter-spring amplitude unchanged. Notably, these are the same years in which ZOTTO growth rates exceeded those of the MBL. In contrast, 2019 and 2021 exhibited both strong growth rates and enhanced seasonal amplitudes, particularly driven by elevated winter CH₄ mole fractions. The underlying causes of these patterns remain uncertain. To address this, future studies should employ atmospheric inverse modelling, which integrate CH₄ observations with atmospheric transport models to estimate fluxes. Analysing these inverted fluxes would help identify the regional sources responsible for the enhanced CH₄ levels that increased the ZOTTO growth rate in 2014, 2016 and 2019–2021, as well as clarify the drivers behind the elevated winter CH₄ observed in 2019 and 2021.

We observed an increasing, though statistically insignificant, trend in the seasonal amplitude of background CH₄ mole fractions at ZOTTO over the 2010–2021 period. At the Waliguan station (WLG), an insignificant upward trend was also reported, where the seasonal amplitude rose from 15.1 to 23.7 ppb between 1994 and 2019 (Liu et al., 2021). In contrast, Dowd et al. (2023) and Liu et al. (2025) reported a significant decreasing trend in the seasonal CH₄ amplitude across high-latitude Northern Hemisphere sites since the 1980s. This discrepancy may be partly attributed to the longer observational periods used in Dowd et al. and Liu et al. (2025) studies. Additionally, their analyses were based on marine or remote background sites, which are less influenced by local emissions. In contrast, ZOTTO (and WLG) is located inland and is subject to regional influences, including variable contributions from wetlands, fires and to some extend also fossil fuel sources. These regional emissions may have contributed to enhanced seasonal amplitude observed at ZOTTO over the study period.

A notable and distinguishing feature of the CH₄ seasonal cycle at ZOTTO is the presence of a secondary peak in late summer (between July–October). This distinct secondary CH₄ peak observed at ZOTTO is primarily attributed to increased emissions from wetlands located to the west of the station. This pattern is consistent with observations from other boreal wetland systems, including sites in Western Siberia (Sasakawa et al., 2010) and Canadian boreal regions, where late summer CH₄ maxima have been linked to enhanced microbial activity under persistently anaerobic conditions and elevated temperatures in late summer (Pickett-Heaps et al., 2011). Additionally, a late summer decline in OH reactivity over boreal forests may contribute to elevated CH₄ levels. Measurements at the SMEAR II station in Hyytiälä, Finland, indicate that OH reactivity in boreal forest peaks in spring but decreases in late summer (Nölscher et al., 2012), potentially reducing CH₄ oxidation and allowing for more CH₄ in late summer. The late summer peaks at ZOTTO show a significant increasing trend over 2010–2021 period, suggesting there is potential enhancing wetland activity over the years.

4.3 Diurnal Cycle

At ZOTTO, CH₄ mole fractions follow a marked diurnal cycle, with higher nighttime and lower daytime mole fractions. After sunset, rapid canopy cooling creates a stable NBL that traps CH₄ near the surface. By day, surface warming enhances vertical mixing, dispersing CH₄. The CH₄ diurnal cycle is most pronounced in warmer months due to larger diurnal temperature ranges, which amplify both daytime mixing and nighttime inversions, increasing discrepancy in CH₄ mole fraction between day and night.

We examine the 2010–2021 trend and interannual variations in the diurnal amplitude at ZOTTO focusing on late summer months (July–October) to further explore the potential factors contributing to the increased amplitude of the unique late-summer seasonal peak observed at the site. We observed a significant increase in the late summer diurnal amplitude of CH₄ at ZOTTO from 2010 to 2021, primarily driven by the significant rise in the nighttime CH₄ maxima. Our analysis of high-pressure system cases revealed no observed trends in synoptic or local atmospheric processes over the 11 years, either during the day or at night, that would suggest changes in boundary layer and synoptic dynamics as the main drivers of the increasing diurnal amplitude trend. Instead, there is a strong significant positive correlation between nighttime surface flux and the CH₄ diurnal amplitude. This relationship is expected, as our flux estimates are derived from CH₄ mole fraction measurements; increases in nighttime CH₄ mole fraction naturally led to higher inferred nighttime surface fluxes. For a more independent assessment of surface fluxes, eddy covariance measurements would be preferable. Nonetheless, our findings underscore that the increase in summer CH₄ diurnal amplitude at ZOTTO is mainly attributable to changes in surface emission processes, rather than shifts in atmospheric boundary layer structure or synoptic conditions over the study period.

There is also a significant increase in nighttime net CH₄ surface flux, particularly in August over the study period, indicating an intensification of late-summer emissions. This could potentially contribute to the increasing in the late-summer seasonal peak at ZOTTO. Similar to our finding at ZOTTO, Rößger et al. (2022) also observed pronounced seasonal CH₄ flux peaks in both July and August in the North Siberian Lena River Delta tundra site (72.37° N, 126.50° E) (2002–2019). However, they found that long-term increases in CH₄ emissions were limited to the early summer months (June–July), with no significant upward trend in August fluxes, despite August being the period of maximum emissions. The main differences in the August CH₄ flux trends between our study at ZOTTO and Rößger et al. (2022) likely stem from differences in the long-term trends of environmental drivers at each site, particularly soil temperature. Rößger et al. (2022) attributed the stability of August fluxes at the Lena River Delta to relatively insignificant small increases in August soil temperature over their study period. In contrast, the increase in late summer nighttime net surface CH₄ fluxes observed at ZOTTO is significantly positively correlated with rising soil temperature, soil moisture during the late summer, and snow depth in the preceding spring during the 2010–2021 period.

The strong correlations between summer nighttime net surface CH₄ flux and soil temperature and soil moisture reinforce the well-established relationship between microbial CH₄ production and environmental conditions (Basu et al., 2022; Bridgham et al., 2013). In contrast, the weak correlation with precipitation suggests that short-term rainfall events have a limited influence on CH₄ variability. These findings align with studies from other boreal and wetland-dominated regions, where CH₄ emissions peak in late summer due to sustained high soil temperature and moisture leading to high microbial activity. For instance, Bohn et al. (2015) observed increasing late-summer CH₄ emissions in Siberian peatlands, driven by persistent anaerobic conditions and enhanced methanogenesis. Our results also reveal a strong positive relationship between spring snow depth and late-summer CH₄ fluxes, indicating that deeper snow in the preceding winter-spring period could enhance CH₄ emissions during late growing season. This finding is consistent with Kivimäki et al. (2025), using satellite observations, identified snow depth as a key driver of the variations in CH₄ seasonality. We hypothesise that thicker winter snowpacks act as an insulating layer, maintaining warmer subsurface temperatures that promote CH₄ production during winter. During spring, snowmelt of larger snowpacks lead to stronger increases soil moisture, creating and maintaining anaerobic conditions that persist longer throughout the late growing season. The flat topography, impermeable subsurface layers, and poor drainage characteristic of western Siberia further enhance water retention, while the warmer soil temperatures in July–October promote the persistence of methanogenesis that sustains elevated CH₄ emissions in the late summer. These findings underscore the importance of assessing the effects of environmental drivers not just for isolated snapshots in time but also considering their interactions over seasonal timescales.

High net surface CH₄ fluxes recorded in June and July 2012, as well as July and August 2019, coincided with major wildfire events, which have been associated with increased emissions of CO, CO₂, and PM_2.5 aerosols in Siberia (Tran et al., 2024; Mokhov and Sitnov, 2022; Bondur et al., 2020). The 2012 wildfire season, one of the most severe in the decade, saw more than 17 000 wildfires detected in July and August alone. Satellite data indicated approximately 29 000 fire sources with a total fire radiative power of ∼3 TW across a region from 50–75° N, 60–140° E in July 2019 (Bondur et al., 2020). Such extreme events significantly contribute to regional CH₄ variability, both directly through biomass burning and indirectly by altering wetland hydrology and soil organic matter decomposition.

In this study, the analysis of the diurnal CH₄ net surface fluxes is limited to nighttime, as our method, based on vertical CH₄ and potential temperature gradients, becomes unreliable during the day due to strong mixing, which minimises vertical gradients and hinders accurate flux estimation. Future studies could apply other alternative methods, such as the Monin-Obukhov Similarity Theory (MOST) approach (Physick and Garratt, 1995), or to derive more accurate estimates of daytime net CH₄ surface fluxes. Another limitation of this study is that the flux analysis is constrained to the summer months, as meteorological sensors at ZOTTO are highly susceptible to icing during colder months, leading to inaccurate or incomplete meteorological measurements, which are required for the net surface flux estimation. Consequently, this seasonal limitation restricts our ability to assess the complete annual cycle of the net CH₄ flux, particularly to investigate whether there is a trend in wintertime fluxes contributing to the observed increasing in the seasonal amplitude of CH₄ at ZOTTO. An eddy covariance system, which provides continuous and more direct measurements of surface-atmosphere exchange, could help overcome both the daytime and seasonal limitations by enabling more accurate, year-round flux estimates independent of vertical gradient assumptions.

Further studies are needed to accurately attribute the sources responsible for the observed increasing nighttime summer CH₄ fluxes. The prevailing summer wind patterns at ZOTTO primarily originate from the west and southwest (Appendix B – Fig. B1), where extensive inland marshes dominate the landscape (Zhang et al., 2023). While this suggests that enhanced wetland activity is a major driver of increased CH₄ emissions, wind direction analysis alone does not provide precise source attribution.

To overcome this limitation, future studies could incorporate CH₄ isotopic analyses from flask samples collected at the 301 m level of the ZOTTO tower, combined with atmospheric inverted modelling. This integrated approach enables more precise attribution of CH₄ signals to distinct emission types, such as wetlands, pyrogenic sources, and fossil fuel combustion, and improves our ability to quantify their relative contributions. Ultimately, this methodology could further enhance understanding of the seasonal and spatial dynamics of CH₄ sources at ZOTTO and provide a stronger basis for evaluating and refining CH₄ emission inventories across the region.

5 Conclusions

We investigate the temporal variability of CH₄ in Central Siberia across annual, seasonal and diurnal scales by utilising the 2010–2021 ZOTTO continuous dataset. This study provides new evidence that warming over the last decade has enhanced CH₄ surface net fluxes in Central Siberia. We demonstrate that the observed enhancement of the late-summer (July–October) diurnal amplitude (5.29 ppb yr⁻¹, p=0.001) is driven by an increase in nighttime surface fluxes and not by changes in atmospheric dynamics. This nighttime surface flux is positively correlated with soil temperature and soil moisture (R²=0.65, p<0.01; and R²=0.36, p=0.032 respectively) during July–October, and with February–May snow depth (R²=0.54, p=0.029). The increase in net CH₄ surface flux during the summer months suggests a growing contribution from wetlands. Episodes of high CH₄ growth rates at ZOTTO are observed in 2012 and 2019 mainly due to wildfires, in 2016–2021 mainly due to increase wetland activity and in 2014 and 2020 mainly due to reduction in OH concentration.

The seasonal analyses reveals an insignificant upward trend of CH₄ winter peak-spring trough amplitude at ZOTTO over the past decade driven by an increase in winter peak. The underlying causes of this increase remain uncertain and warrant further investigation through atmospheric inversion and isotopic studies. A significant increase in the distinct late summer CH₄ peak in August at ZOTTO underscores the enhancement of regional wetland emissions. The persistent rise in CH₄ mole fractions at ZOTTO reflects global trends, underscoring the sustained impact of biogenic emissions, especially from wetlands, which are increasingly active due to rising global temperatures.

Our study highlights the importance of CH₄ surface fluxes in driving diurnal variations, but the need to include in the analysis the effects of the dynamics of the ABL. In this context, advancing our understanding of soil microbial activity through direct measurements could improve estimates of ecosystem surface fluxes at both daily and seasonal scales. Extending these findings to the regional scale, combined with inverse modelling and isotopic analysis, will enhance our ability to accurately attribute CH₄ sources and sinks at ZOTTO. Continued monitoring and improved modelling efforts are critical for refining our understanding of CH₄ variability and assessing its implications for future climate feedback.

Appendix A: ZOTTO Meteorology Analysis

Figure A1 Timeseries of yearly temperature measured at 52 m a.g.l. at ZOTTO for each season over the 2010–2021 period.

Figure A2 (a) Temperature measured at 52 m a.g.l. (in °C) and (b) precipitation measured at 2 m a.g.l. (mm) for the period 2009–2021 at ZOTTO.

Figure A3 Wind rose showing the average wind speed distribution at 300 m a.g.l. at ZOTTO from 2010 to 2021, categorized by 12 horizontal wind directions across the four seasons.

Figure A4 Climatological (2010–2021) late summer (JASO) vertical profile of potential temperature at ZOTTO. The dashed green line denotes the canopy height at ZOTTO, while the dashed brown and grey lines indicate the 52 and 72 m levels, where the F_Stor and F_Eddy terms in Eq. (4) in the main text are calculated, respectively. At night, the vertical gradient in potential temperature between the 227 and 301 m levels is minimal, indicating that the 301 m level is already within the residual layer during nighttime conditions.

Appendix B: Target Tank Time Series

Figure B1 Target tank time series (coloured line represents the mean ± standard deviation, grey is laboratory standard ± error).

Appendix C: Formula Symbols and Units

Table C1 Formula Symbol and Units.

Appendix D: Details explanations of the Atmospheric Drivers of CH₄ Diurnal Amplitude

D1 Governing Equations for the Atmospheric Drivers of CH₄ Diurnal Amplitude

Terms I and III of Eq. (1) represent key atmospheric dynamics with distinct daytime and nighttime characteristics that influence both CH_4,min and CH_4,max and will be discussed below.

D1.1 Atmospheric Boundary Layer Height – Term I

Regarding Term I in Eq. (1), the diurnal dynamics of the atmospheric boundary layer height (ABL) (h) are driven by the surface buoyancy flux introduced into the ABL. This process is represented by the potential temperature variable (θ), solved using three additional equations as described in Appendix D2 and more detailed in Vilà-Guerau de Arellano et al. (2015). Note that we are not using Appendix D2 to solve for the ABL height in this study but to solely explain the processes. In short, the effect of h on CH₄ is the following. In the morning, solar heating destabilizes the atmospheric column, warming the surface and causing less dense air to rise. This results in an upward transfer of surface heat flux (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}>\mathrm{0}$), driving turbulent convective motions. These turbulent processes increase the height of the CBL, expanding the atmospheric volume available for CH₄ dilution, decreasing CH₄ mole fractions.

At night, radiative cooling of the ground creates a temperature gradient where heat flows downward (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}<\mathrm{0}$) from the warmer air to the cooler surface. This cooling stabilizes the lower atmospheric layers, leading to the formation of a stratified NBL, typically ranging between 100 and 300 m above ground level (Kubiak and Zimnoch, 2022). The NBL traps surface-emitted CH₄, limiting its vertical dispersion and promoting nighttime CH₄ accumulation, increasing CH₄ mole fractions.

D1.2 Entrainment flux – Term III

The entrainment flux (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}$) in Term II in Eq. (1) could be written as in Vilà-Guerau de Arellano et al. (2015) and expressed as below:

$$\begin{array}{}\text{(D1)}& {\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}=-\left(\left(\underset{\text{III.1}}{\underbrace{{\displaystyle \frac{\partial h}{\partial t}}}}-\underset{\text{III.2}}{\underbrace{w_{\text{sub}\left(h\right)}}}\right)\times \underset{\text{III.3}}{\underbrace{{\mathrm{\Delta }}_{\mathit{\varphi }(\text{ft/rs}-\text{cbl/nbl})}}}\right)\end{array}$$

According to Eq. (D1), the entrainment flux (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}$) depends on: the growth rate of the ABL height ($\frac{\partial h}{\partial t}$, Term III.1) which is dependent on the surface buoyancy flux; the large-scale vertical subsidence velocity (w_sub(h), Term III.2); and the difference in CH₄ mole fraction between the CBL (or NBL) and the overlying layer (Δ_{ϕ(ft/rs−cbl/nbl)}, Term III.3) – either the free troposphere (FT) during the day or the residual layer (RS) at night (Fig. 1 in the main text). For the latest we assume that this jump occurs in an infinitesimal layer (zero-order approach) (Driedonks and Tennekes, 1981) Higher entrainment flux rates introduce more CH₄-depleted air from the overlying layer into the integrated CH₄ column, reducing the overall CH₄ mole fraction. This exerts a negative impact on both CH_4,max and CH_4,min.

To account for the distinct atmospheric dynamics between day and night, each term in Eq. (D1) is applied differently for nighttime and daytime, as detailed in the following sections.

Growth Rate of the ABL – Term III.1

The growth rate of the boundary layer height ($\frac{\partial h}{\partial t}$) is closely linked to the surface buoyancy flux (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}$), as discussed earlier similarly. This term is more pronounced during the daytime, as strong turbulence and convection drive rapid changes in CBL height, while weaker nighttime turbulence results in slower changes in NBL height. A stronger $\frac{\partial h}{\partial t}$ amplifies the entrainment flux, making entrainment more pronounced during the daytime than nighttime.

Vertical Subsidence Velocity – Term III.2

The vertical subsidence velocity represents large-scale downward motion in the atmosphere, primarily driven by synoptic-scale conditions. Using the mass conservation equation assuming incompressibility, we represent the vertical subsidence velocity at ABL height (w_sub(h)) as:

$$\begin{array}{}\text{(D2)}& w_{\text{sub}\left(h\right)}=-\text{Div}\left(\stackrel{\rightharpoonup }{U_h}\right)\cdot h\end{array}$$

Where $\text{Div}\left(\stackrel{\rightharpoonup }{U_h}\right)$ is the horizontal wind divergence. While subsidence velocities are generally small (rarely exceeding a few cm s⁻¹) and the same magnitude as the entrainment velocity, they can significantly influence mass conservation and the growth of the CBL/NBL (Stull, 1988), and consequently, the entrainment flux.

During the daytime, subsidence slows the growth of the CBL by introducing downward motion, which counters the upward expansion driven by surface heating and turbulence. This downward motion brings warmer, drier air from the free troposphere into the CBL, stabilizing the atmosphere and weakening convective activity. At the same time, the temperature and moisture contrast between the warm, dry overlying air and the CBL air enhances entrainment at the top of the CBL.

At night, divergence associated with subsidence laterally transport cold air masses generated by the longwave radiative cooling at NBL air causing the NBL to not grow as rapidly as would otherwise be expected (Carlson and Stull, 1986). This additional stabilization further suppresses turbulence and reduces vertical mixing within the NBL. While nighttime entrainment is minimal, subsidence could still contribute to the transport of air from the residual layer above, influencing the temperature and composition of the NBL (Carlson and Stull, 1986).

The Difference in CH₄ Mole Fraction between the CBL (or NBL) and the Overlying Layer – Term III.3

The entrainment flux of CH₄ (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\varphi }}^{\prime }\right)}}_{\mathrm{e}}$) is also influenced by the difference in CH₄ mole fractions between the daytime (nighttime) CBL (NBL) and the layer above it, which is the FT (RS) (Δ_{ϕ(ft/rs−cbl/nbl)}). This difference evolves over time and can be expressed as:

$$\begin{array}{}\text{(D3)}& \begin{array}{rl}{\displaystyle \frac{\partial {\mathrm{\Delta }}_{\mathit{\varphi }(\text{ft/rs}-\text{cbl/nbl})}}{\partial t}}=& {\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{ft/rs}\right)}}{\partial t}}\\ & -{\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}{\int }_{z=\text{ct}}^{z=h}{\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{cbl/nbl}\right)}}{\partial t}}\mathrm{d}z\end{array}\end{array}$$

Here, ϕ_(ft/rs) represents the CH₄ mole fraction in the FT during the daytime or in the RS during the nighttime. The last term on the right-hand side of Eq. (D3) refers to the averaged column integrated CH₄ mole fraction from the canopy top to the top of the CBL during the day or the NBL during the night.

During the daytime, strong convective turbulence in the CBL leads to a well-mixed CH₄ distribution, making its mole fraction independent of height. This is evident from minimal vertical gradients in CH₄ mole fraction, as shown in Fig. 3 in the main text and Fig. G2 in Appendix G (i.e. P3 period at ZOTTO). Under such conditions, the Eq. (D3) could be simplified to:

$$\begin{array}{}\text{(D4)}& {\displaystyle \frac{\partial {\mathrm{\Delta }}_{\mathit{\varphi }(\text{ft}-\text{cbl})}}{\partial t}}={\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{ft}\right)}}{\partial t}}-{\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{cbl}\right)}}{\partial t}}=\mathit{\gamma }{\displaystyle \frac{\partial h}{\partial t}}-{\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{cbl}\right)}}{\partial t}}\end{array}$$

Here, the rate of change of the free tropospheric CH₄ mole fraction over time ($\frac{{\mathit{\varphi }}_{\left(\text{ft}\right)}}{\partial t}$) is proportional to the product of the tropospheric CH₄ lapse rate (γ) and the boundary layer height growth ($\frac{\partial h}{\partial t}$). Assuming a zero tropospheric CH₄ lapse rate (γ=0) as in Faassen et al. (2024), integration yields:

$$\begin{array}{}\text{(D5)}& {\mathrm{\Delta }}_{\mathit{\varphi }(\text{ft}-\text{bl})}\left(t_i\right)={\mathrm{\Delta }}_{\mathit{\varphi }}\left(t_{\mathrm{0}}\right)+\left({\mathit{\varphi }}_{\text{cbl}}\right(t_{\mathrm{0}})-{\mathit{\varphi }}_{\text{cbl}}(t_i\left)\right)\end{array}$$

Where ϕ_cbl(t₀) and Δ_ϕ(t₀) represent the initial CH₄ mole fraction and the initial difference in CH₄ mole fraction between the ABL and the layer above it just before sunrise (e.g., P2 period in Fig. G2 in Appendix G). (ϕ_cbl(t₀)−ϕ_cbl(t_i)) represents the change in total CH₄ mole fraction in the CBL over time. The initial difference (Δ_ϕ(t₀)) is influenced by nighttime stability, which will be discussed below.

During nighttime, the is a clear vertical gradient of CH₄ mole fraction (Fig. 3 in the main text and Fig. G2 in Appendix G, i.e., P1 period at ZOTTO), suggesting integrated column of CH₄ depends on height. Equation (D3) during nighttime could be written as:

$$\begin{array}{}\text{(D6)}& {\displaystyle \frac{\partial {\mathrm{\Delta }}_{\mathit{\varphi }(\text{rs}-\text{nbl})}}{\partial t}}={\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{rs}\right)}}{\partial t}}-{\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}{\int }_{z=\text{ct}}^{z=h}{\displaystyle \frac{\partial {\mathit{\varphi }}_{\left(\text{nbl}\right)}}{\partial t}}\mathrm{d}z\end{array}$$

During the night, the CH₄ mole fraction in the RS remains almost constant on time ($\frac{\partial {\mathit{\varphi }}_{\left(\text{rs}\right)}}{\partial t}=\mathrm{0}$) because the RS is largely decoupled from surface processes (Stull, 1988). In absence of sources and sinks of CH₄, the mole fraction of CH₄ remains almost constant. Essentially, the RS “stores” the composition of the previous daytime mixed layer. As a result, the difference in CH₄ mole fractions between the RS and NBL depends primarily on the rate of averaged CH₄ accumulation from the top of the canopy to NBL height (${\int }_{z=\text{ct}}^{z=h}\frac{\partial {\mathit{\varphi }}_{\left(\text{nbl}\right)}}{\partial t}$). Since this term is relatively constant overnight (as seen in P1 period in Fig. D1 and in Winderlich, 2012), $\frac{\partial {\mathrm{\Delta }}_{\mathit{\varphi }(\text{rs}-\text{nbl})}}{\partial t}$ remains minimal, resulting in limited nighttime entrainment.

Figure D1 Climatological (2010–2021) late summer (JASO) diurnal cycle of the column-integrated of CH₄ mole fraction from 52 to 300 m, representing the layer above the canopy to the top of the NBL during P1. The 52–300 m column-integrated flux is calculated using a method similar to Eq. (6) in the main text but applied to the 52–300 m layer. This calculation accounts for temporal changes in the 52–300 m column-integrated CH₄ mole fraction and the influence of vertical subsidence velocity at 300 m.

However, this dynamic changes significantly at sunrise. As surface heating begins, convection resumes, and the turbulent eddies lead to well-mixed conditions and reconnecting it with the RL. This process leads to a sharp reduction in the averaged CH₄ accumulation from the top of the canopy to NBL height (as seen in P2 period in Fig. D1), as CH₄ accumulated near the surface during the night is rapidly mixed into the expanding convective boundary layer (CBL). Consequently, there is a dramatic increase in the entrainment flux at sunrise.

The thermal stability of the nighttime atmospheric column plays a key role in this process. A more stable column leads to a lower NBL height, which increases the near-surface CH₄ mole fraction and enhances the storage flux overnight as well as larger CH₄ difference between the RS and the NBL (Δ_ϕ(t₀)). At sunrise, the greater storage flux results in a larger reduction when convection begins, as well as a higher Δ_ϕ(t₀) right before sunrise, driving stronger entrainment flux into the CBL.

D2 Governing Equations for the evolution of the Atmospheric Boundary Layer (ABL)

The dynamic evolution of the ABL is solely driven by the heat introduced in the ABL, represented by the virtual potential temperature (θ) variable.

To describe the evolution of the ABL potential temperature (θ) and the discontinuity jump of potential temperature (Δθ) at the inversion (i.e. top of the ABL), we solve three fundamental equations as described in Vilà-Guerau de Arellan et al. (2015) and in Driedonks and Tennekes (1981). These equations result from the vertical integration of the conservation equations for θ within the canopy-top-to-ABL-top layer and the entrainment zone.

$$\begin{array}{}\text{(D7)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}& {\int }_{z=\text{ct}}^{z=h}{\displaystyle \frac{\partial \mathit{\theta }}{\partial t}}\left(z\right)\mathrm{d}z=\\ & {\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}\times \left({\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}-{\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{e}}\right)-\text{adv}\left(\mathit{\theta }\right)\end{array}\text{(D8)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{\partial \mathrm{\Delta }\mathit{\theta }}{\partial t}}& ={\mathit{\gamma }}_{\mathit{\theta }}\left({\displaystyle \frac{\partial h}{\partial t}}-w_{\mathrm{s}}\right)-{\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}{\int }_{z=\text{ct}}^{z=h}{\displaystyle \frac{\partial \mathit{\theta }}{\partial t}}\left(z\right)\mathrm{d}z\\ & ={\mathit{\gamma }}_{\mathit{\theta }}\cdot w_{\mathrm{e}}-{\displaystyle \frac{\mathrm{1}}{h-\text{ct}}}{\int }_{z=\text{ct}}^{z=h}{\displaystyle \frac{\partial \mathit{\theta }}{\partial t}}\left(z\right)\mathrm{d}z\end{array}\end{array}$$

Equation (D7) shows that within the top of the canopy height to the ABL height, the tendency term of θ on the left hand sides depends on the vertical turbulent flux difference between the surface heat flux (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}$) and entrainment zone heat flux ${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{e}}$, the horizontal advection adv(θ) which is current not accounted for.

The evolution of the discontinuity or jump value Δθ at the entrainment zone, see Eq. (D8), is a function of the tendency value at the residual layer or free troposphere (first term right-hand side) and the evolution of the canopy-top-to-ABL-top value (2nd term r.h.s). Above the jump, the profile of θ in the layer above the ABL is dependent on the vertical gradient (γ_θ) and on the mean subsidence vertical velocity (w_s). This velocity is normally opposite to the boundary layer growth ($\frac{\partial h}{\partial t}$), i.e., subsidence (w_s<0). We assume that the w_e is a function of the entraiment flux and the jump in the virtual potential temperature in the inversion layer. This assumption is known as zero-order closure, and it was first suggested by Lilly (1968). It is expressed mathematically by:

$$\begin{array}{}\text{(D9)}& w_{\mathrm{e}}={\displaystyle \frac{\partial h}{\partial t}}-w_{\mathrm{s}}=-{\displaystyle \frac{{\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{e}}}{\mathrm{\Delta }\mathit{\theta }}}\end{array}$$

Equation (D9) assumes that the inversion is characterised by a sharp discontinuity (Driedonks and Tennekes, 1981). Under conditions of weak inversion, it is convenient to include explicitly the inversion depth requiring a modification of Eq. (D9) (Kim et al., 2006). For θ, Eqs. (D7)–(D9) contain seven variables: h, $\frac{\mathrm{1}}{h-\text{ct}}{\int }_{z=\text{ct}}^{z=h}\mathit{\theta }\left(z\right)\mathrm{d}z$, Δθ, ${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}$, ${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{e}}$, γ_θ and w_s. The first three are solved by the system Eqs. (D7)–(D9) the other four need to be prescribed or calculated using additional equations or closure assumptions. The heat surface fluxes (${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}$) are either prescribed based on field measurements or calculated using a coupled land-surface scheme. The subsidence velocity (w_s) and the potential temperature lapse rate in the ABL overlaying layer (γ_θ) depend on the atmosphere at large scales. These upper boundary conditions are thus obtained either from large-scale models or by a radiosounding taken in the early morning hours of the ABL development. In consequence, to close the set of Eqs. (D7)–(D9), we still need to relate the entrainment of heat flux to the surface flux. We assume the following relation, ${\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{e}}=-\mathit{\beta }{\stackrel{\mathrm{‾}}{\left(w^{\prime }{\mathit{\theta }}^{\prime }\right)}}_{\mathrm{s}}$, where β represents an additional percentage of entrainment of warm air into the ABL. Here, it needs to be mentioned that the β value can increase, depending on the contribution of shear in the ABL development (Angevine et al., 1998; Pino et al., 2003; Conzemius and Fedorovich, 2006).

Appendix E: Nocturnal Boundary Layer Height Estimation

This section visual the method to estimate the Nocturnal Boundary Layer (NBL) height from the regression fit (Eqs. 2 and 3 in the main text) applying to the vertical gradient of CH₄ and potential temperature in Fig. E1. In some years, summer NBL heights based on vertical temperature profiles are unavailable (Fig. E3) due to meteorological instrument malfunctions at specific heights, resulting in incomplete data for constructing full vertical profiles.

The NBL heights derived from vertical CH₄ concentrations and potential temperature exhibit similar ranges, generally falling between 100–150 m (Figs. E2 and E3), and show comparable interannual variability over the 2010–2021 period. However, there is no strong 1:1 correlation between the NBL derived by the two parameters as shown in Fig. E4. This discrepancy may result from a time lag in the vertical profile development between potential temperature and CH₄. The nighttime vertical CH₄ profile decreases with height (Fig. 3), mirroring the pattern of potential temperature (Fig. A4). The nighttime vertical CH₄ profile stability persists from 00:00 to 08:00 LT, peaking at 08:00 LT and beginning to weaken around 10:00 LT, later than the temperature profile. This observed delay between potential temperature and measurement gases has also been observed in the 213 m tower in Cabauw (CBW: 51°97^′ N, 04°93^′ E, 0 m a.s.l.) in the Netherlands (Casso-Torralba et al., 2008). This could be caused by the larger difference of the CH₄ mole fraction between the ABL and the free troposphere compared to potential temperature (Casso-Torralba et al., 2008).

Figure E1 Vertical profile of potential temperature during the nighttime on 21 June 2010, presented as an example day for estimating the nocturnal boundary layer height. The curve fit and normalised lapse rate are calculated using Eqs. (2) and (3) respectively in the main text.

Figure E2 Box-and-whisker plot of the yearly nocturnal boundary layer (NBL) height (averaged from 00:00 to 04:00 LT) at ZOTTO, derived from CH₄ vertical profiles for each late summer month (JASO) by applying Eqs. (2) and (3) in the main text to CH₄. The box denotes the interquartile range (IQR), showing the median with a thick black line. The whiskers range from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{15}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{15}\times \text{IQR}$, with Q₁ and Q₃ being the 25th and 75th percentiles, respectively. The blue line is the monthly mean. Numbers above each box indicate the sample size or the number of available days (based on the number of high-pressure days and fits (Eq. 2) with an R² value greater than 0.7) for analysis in that month.

Figure E3 Box and whisker plot of yearly nocturnal boundary layer (NBL) height (00:00–04:00 LT averaged) height derived from potential temperature vertical profile criteria at ZOTTO for each late summer month (JASO) by applying Eqs. (2) and (3) in the main text to potential temperature. The box denotes the interquartile range (IQR), showing the median with a thick black line. The whiskers range from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{15}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{15}\times \text{IQR}$, with Q₁ and Q₃ being the 25th and 75th percentiles, respectively. The blue line is the monthly mean. Numbers above each box indicate the sample size or the number of available days for analysis (based on the number of high-pressure days and fits (Eq. 2) with an R² value greater than 0.7) in that month.

Figure E4 Comparison of the monthly average summer NBL height derived from the potential temperature vertical profile and the CH₄ vertical profile using Eqs. (2) and (3).

Appendix F: Climatology Nighttime Column Integrated CH₄ flux up to 52 m at ZOTTO

Figure F1 Climatological (2010–2021) late summer (JASO) nighttime (00:00–04:00 LT) CH₄ net surface flux up to 52 m. The individual components of Eq. (4) in the main text are plotted in (a) and their contributions to the total flux net surface flux (in percentage) are plotted in (b).

Appendix G: Additional Annual Growth Rate, Seasonal and Diurnal Analyses for CH₄ mole fractions at ZOTTO

Figure G1 Time series of the CH₄ spring minimum, winter maximum and late summer maximum for detrended background-filtered daytime CH₄ at ZOTTO (ZOTbg). The Theil-Sen regression trend is depicted by the solid line, with the 95 % confidence interval of the trend shown as dashed lines. The p value indicates whether the slope of the regression is significantly different (at 0.05 level) from zero.

Figure G2 Climatological (2010–2021) late summer (JASO) diurnal cycle of CH₄ mole fractions measured at six different heights at ZOTTO. The shadings represent different time periods: P1 corresponds to the nighttime period (17:00 to 07:00 LT (KRAT, UTC+7)) the following day), P2 marks the onset of mixing around sunrise (07:00 to 11:00 LT), and P3 represents the well-mixed period (11:00 to 17:00 LT).

Figure G3 Yearly diurnal cycle of CH₄ from 2010 to 2021 at ZOTTO 52 m a.g.l. (i.e. above the canopy), in different seasons using non-detrended data. The shaded colours show 95 % confidence interval of the mean.

Figure G4 Time series of yearly late summer (JASO): (a) averaged CH₄ diurnal cycle amplitude; (b) its daytime (10:00–16:00 LT averaged) CH₄ mole fraction, and (c) its nighttime (00:00–04:00 LT averaged) CH₄ mole fraction (right) (circle markers) at ZOTTO using detrended 52 m a.g.l. data. The Theil-Sen regression trend is depicted by the solid line, with the 95 % confidence interval of the trend shown as a dashed line. The p value indicates whether the slope of the regression is significantly different (at 0.05 level) from zero.

Appendix H: Number of High-pressure Days for Each Month from 2010–2021

Figure H1 Number of high-pressure days identified using ERA5 geopotential height data at 550 hPa, selecting periods where geopotential height exceeded the 90th percentile of yearly summer distribution.

Appendix I: Trend Analysis for Atmospheric Drivers of CH₄ Diurnal Amplitude in Summer Months

The interannual variations in atmospheric process drivers influencing the summer CH₄ diurnal amplitude (dashed yellow boxes in Fig. 2 in the main text) over the 2010–2021 are analysed in detail in this section.

Interannual variations are observed in the heights of both the Convective Boundary Layer (CBL) and the Nocturnal Boundary Layer (NBL), but no significant long-term trends are detected over the study period (Figs. I1 and E3). At ZOTTO, the 12:00–16:00 LT averaged CBL height typically reaches approximately 1500 m, while the NBL height, averaged from 00:00 to 04:00 LT, generally ranges between 100 and 150 m.

Figure I1 Box and whisker plot of yearly daytime (12:00–16:00 LT averaged) convective boundary layer (CBL) height from ERA5 for each late summer month (JASO). The box denotes the interquartile range (IQR), showing the median with a thick black line. The whiskers range from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{15}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{15}\times \text{IQR}$, with Q₁ and Q₃ being the 25th and 75th percentiles. The red line is the monthly-mean. Numbers above each box indicate the sample size or the number of available days for analysis in that month.

The cumulative daytime sensible heat flux at 52 m from ZOTTO and the ERA5 surface heat flux exhibit similar interannual variability, with no significant long-term trend detected (Fig. I2). However, notable month-to-month differences in magnitude are observed between the two datasets. In July and October, both datasets align closely, showing similar median values, interquartile ranges, and interannual variability. In August and September, larger discrepancies occur, with ERA5 underestimating the observed flux, while observations display greater variability and higher extreme values. A distinct anomaly was observed in 2012 in August and September, when both datasets showed daytime sensible heat flux reached 1400–2000 W m⁻².

The divergence at 750 hPa from ERA5 show no significant trends and minimal variation over the 2010–2021 period, indicating the absence of a long-term change in synoptic-scale subsidence over the ZOTTO region (Fig. I3).

Figure I2 Box and whisker plot of yearly cumulative positive sensible heat flux of ERA5 and ZOTTO observational measurement at 52 m a.g.l. for each late summer month (JASO). The box denotes the interquartile range (IQR), showing the median with a thick black line. The whiskers range from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{15}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{15}\times \text{IQR}$, with Q₁ and Q₃ being the 25th and 75th percentiles, respectively. The line plots are the monthly mean. Numbers above each box indicate the sample size or the number of available days for analysis in that month (based on the number of high-pressure days).

Figure I3 Box and whisker plot of yearly divergence at 750 hPa from ERA5 for each late summer month (JASO). The box denotes the interquartile range (IQR), showing the median with a thick black line. The whiskers range from ${\mathrm{Q}}_{\mathrm{1}}-\mathrm{15}\times \text{IQR}$ to ${\mathrm{Q}}_{\mathrm{3}}+\mathrm{15}\times \text{IQR}$, with Q₁ and Q₃ being the 25th and 75th percentiles, respectively. The blue line is the monthly mean. Numbers above each box indicate the sample size or the number of available days for analysis in that month (based on the number of high-pressure days).