Uncertainty of hourly-average concentration values derived from non-continuous measurements

Continental greenhouse gas monitoring networks extensively use tall towers for higher spatial representativeness. In most cases, several intakes are built along the tower to give information also on the vertical concentration profile of the components considered. Typically, a single gas analyzer is used, and the intake points are sequentially connected to the instrument. It involves that the continuous concentration signal is only sampled for discrete short periods at each intake point, which does not allow for a perfect reconstruction of the original concentration variation. It increases the uncertainty of the calculated hourly averages usually used by the atmospheric transport and budget models. The purpose of the study is to give the data users an impression of the potential magnitude of this kind of uncertainty, as well as how it depends on the number of intakes sampled, on the length of the sampling period at each intake, on the season, and on the time of the day. It presents how much improvement can be achieved using linear or spline interpolation between the measurement periods instead of the simple arithmetic averaging of the available measurements. Although the results presented here may be site-specific, the study calls attention to the potentially rather heterogeneous spatial and temporal distribution of the uncertainty of the hourly-average concentration values derived from tall-tower measurements applying sequential sampling.

1 Introduction

Continental greenhouse gas (GHG) monitoring networks extensively use tall towers for measurements for obtaining higher spatial representativeness (Tans, 1991; Bakwin et al., 1995; Wofsy and Harriss, 2002; Vermeulen, 2007; Gerbig et al., 2009; Kadygrov et al., 2015; Oney et al., 2015; White et al., 2019; ICOS RI, 2020). In most cases, taking advantage of the tower, several intakes are mounted along the tower to get information on the vertical distribution of the components considered. Typically, a single gas analyzer is used, and the intake lines are connected sequentially to the instrument (Bakwin et al., 1998; Haszpra et al., 2001; Vermeulen, 2007; Thompson et al., 2009; Popa et al., 2010; Winderlich et al., 2010; Vermeulen et al., 2011; Andrews et al., 2014; Satar et al., 2016; Conil et al., 2019; ICOS RI, 2020). Using a single instrument is not only cheaper but also avoids any potential scale differences among the instruments, which could inevitably occur even at frequent synchronization. In the case of continuous in situ measurements, international databases, such as the World Data Centre for Greenhouse Gases of the World Meteorological Organization (https://gaw.kishou.go.jp/, last access: 7 May 2021), ICOS Carbon Portal (https://www.icos-cp.eu/icos-carbon-portal, last access: 7 May 2021), or ObsPack that is maintained by the National Oceanic and Atmospheric Administration, USA (https://www.esrl.noaa.gov/gmd/ccgg/obspack/, last access: 7 May 2021), store and disseminate hourly-average concentration values, and these values are used in the different mathematical models (atmospheric budget models, inverse transport models, etc. – see e.g., Ciais et al., 2010; Rivier et al., 2010; Bergamaschi et al., 2015; Diallo et al., 2017; Shirai et al., 2017; Bergamaschi et al., 2018; Lin et al., 2018; Chevallier et al., 2019). As long as an instrument receives air from a given intake continuously, the hourly-average concentration is the temporal integral of the instrument's signal for the given hour. However, if the instrument switches between the different intakes, only short samples of the continuous concentration signal are available from each intake for the estimation of the hourly average, which increases the uncertainty of the hourly-average concentration values reported.

Both the systematic and random errors of the measurements decrease the reliability of the results of the atmospheric models. Systematic errors such as the scale bias may distort the source/sink distributions calculated by the models and may result in false emission values, while random errors increase the uncertainty of the calculated values. Instrument noise, scale instability, and other processes may cause random errors in the measurements. In this paper, we focus on the random error caused by the non-continuous sampling of the continuous concentration signal, which increases the uncertainty of the calculated hourly-average concentrations.

According to the Nyquist–Shannon sampling theorem, a fluctuation having a higher frequency than half of the sampling one cannot be reconstructed from the data recorded. At continental carbon dioxide monitoring stations surrounded by vegetation, significant short-term concentration changes are rather common, especially during the growing season due to the activity of the plants and the dynamic processes of the atmosphere. Losing the high-frequency part of the spectrum may introduce significant uncertainty into the hourly-average concentration values especially in certain seasons of the year and in certain periods of the day, which adds to the common instrument noise and scale uncertainty. While the reduction of the instrument noise would require a longer signal integration time, it would also result in the loss of high-frequency concentration fluctuation. Significant concentration changes, like those during the morning transition periods (when the low-level inversion breaks up) or during a frontal passage (Pal et al., 2020), cannot be followed properly by discrete sampling either.

For the reduction of the uncertainty derived from episodic sampling, both physical and mathematical methods can be applied. The physical method involves a buffer volume in the sampling line, which physically integrates the high-frequency concentration fluctuation (see e.g., Winderlich et al., 2010). A properly designed weighted averaging scheme can significantly reduce the uncertainty of the hourly-average values (Cescatti et al., 2016). However, in this case, it is difficult to determine the exact start of the integration period for which the calculated hourly average is characteristic, and spike detection (El Yazidi et al., 2018) cannot be applied for quality assurance.

In the present paper, we analyze the uncertainty of the hourly-average carbon dioxide concentration values originating from the discrete episodic signal sampling for the case when no buffer volume is applied in the measuring system. We use measurement data from a mid-continental monitoring site to show how much uncertainty is introduced by the episodic sampling and how it varies in time.

2 Methodology and data

For the study, the continuous in situ carbon dioxide concentration measurements carried out at Hegyhátsál tall-tower GHG monitoring site (46^∘57^′ N, 16^∘39^′ E; 248 m above the sea level, NOAA/WMO GAW code: HUN) at 82 m above the ground are used. The tower is located in a European mid-continental rural environment, at a low elevation above sea level, in the temperate climate zone, and it is surrounded by natural and agricultural vegetation (Haszpra et al., 2001). At such a site, the vegetation and the dynamic processes of the atmosphere may generate significant short-term changes in the concentration of carbon dioxide in the planetary boundary layer sampled by a tower. At the Hegyhátsál tall-tower GHG monitoring site, a Picarro G2301 cavity ring-down spectrometer (CRDS) analyzer (Picarro, Inc., Santa Clara, California, USA) is operated for carbon dioxide (and methane) measurements. The instrument provides the dry mole fraction, also referred to as concentration in this paper, at a 5 s temporal resolution. For the study, the data measured in 2018 were selected. This high temporal resolution data series is considered continuous from which the true hourly-average concentrations are calculated for each hour of the year. Different subsamples of this series were generated to simulate the multilevel measurements when the analyzer cycles through the monitoring elevations sampling each of them only for a short time.

In this study, we suppose 2, 3, 4, and 5 measurement elevations and study a wide range of sampling periods. The length of the sampling period was determined so that each measurement elevation is sampled equal times within the hour (Fig. 1). So, the maximum sampling period is 1800 s for two intakes, and it is 720 s for five intakes. In such extreme situations, only one sampling would be performed at each measurement elevation in each hour.

Figure 1 Schematic of an example measurement protocol with three intakes and 600 s sampling time, providing two complete cycles within an hour. Vertical dashed lines represent the 1 min flush time after switching the intakes and before the signal integration. The signal integration period is indicated by thick black line segments.

The minimum length of the sampling period is technically determined: after switching the intakes, a certain time is needed for flushing the tubing and instrument before valid measurements can be performed. The length of this time depends on the design of the monitoring system. In our case, the intake lines from the tower are permanently ventilated, and only a relatively small volume (short tubes, selector valve, Nafion drier tube, etc. – see Haszpra et al., 2001, for general design) has to be flushed at each switch before the air to be measured arrives at the measuring cell of the analyzer. Our system is operated at a flow rate of ∼220 mL min⁻¹ of air, and at this flow rate it takes <10 s for the sample air to reach the measurement cell of the analyzer after switching the intake. Theoretically, flushing the measuring cell follows an exponential model where the air from intake 2 gradually replaces the air remaining in the measuring cell from intake 1 at continuous mixing. Also theoretically, the perfect flushing needs infinite time in this way. In practice, the difference between the true concentration of air from intake 2 and the concentration in the measuring cell quickly disappears in the range of the noise of the instrument. Our experience with standard gases shows that the deviation of the concentration in the measuring cell from the true value falls below 0.1 µmol mol⁻¹ within 5–20 s if the concentration difference is below 20 µmol mol⁻¹. More precise determination of the response time is not possible due to the coarse temporal resolution of the readings (5 s), which may introduce an error up to 10 s immediately. In the case of higher concentration difference, the response time increases reaching 35–45 s at 70 µmol mol⁻¹ difference, which indicates some sort of a memory effect in the system.

Another question is the length of the signal integration after the acceptable flushing to reduce the effect of the instrument noise in the measured concentration. Laboratory tests with standard gases may help the determination of the optimum averaging time (Yver Kwok et al., 2015). However, laboratory tests do not help much if the concentration in the sample gas changes during the averaging period as it happens during real-world atmospheric measurements. Long averaging time involves less frequent measurement at each intake and increased uncertainty of the calculated hourly averages. On the other hand, at an averaging time that is too short, the instrument noise may dominate the result. For the determination of a reasonable trade-off between these conflicting requirements, a linear regression model was constructed in which the independent variable was the time, while the dependent one was the mole fraction measured. As long as the instrument noise dominates above the natural variability the slope of the linear regression line fitted to the data series does not differ statistically significantly from zero (null hypothesis). We calculated the slopes of the linear regression lines as the function of the length of the involved data series. Using the F test at a significance level of 5 % for the null hypothesis mentioned above, it turned out that even a 40 s long data series shows a statistically significant linear trend in as much as 38 % of the cases, which increases with the length of the data series involved. Taking into account the above experience, we selected the shortest sampling time for this study as 100 s, including 60 s flushing and 40 s signal averaging for the cases when 2, 3, or 4 intakes were assumed. For five intakes, the shortest sampling time was 120 s, including 60 s averaging to provide an equal number of sampling periods at each intake in an hour taking into account the 5 s temporal resolution of the instrument readings.

The dynamic processes in the planetary boundary layer (PBL) significantly differ from those above. Mixing of the cases when the PBL is sampled and when the air sample is taken from the free troposphere or from the nighttime residual layer above would lead to a hardly interpretable result. Therefore, in the present study, we selected only those periods when the measurement elevation (82 m above the ground) was well within the planetary boundary layer; i.e., the depth of the PBL was at least 120 m at both the beginning and end of the hour considered. The PBL depth data were taken from the ERA5 reanalysis dataset of the European Centre for Medium-Range Weather Forecasts (Copernicus Climate Change Service, 2017) with 1 h temporal resolution. Using the criterion mentioned above, we could also avoid periods when the top of the PBL is around the measurement elevation, and its fluctuation causes extreme variability in the measured carbon dioxide concentration. These periods are also difficult to handle in the atmospheric models.

Figure 2 Distribution of the deviations of the hourly averages estimated as the arithmetic average of the available measurements from the true hourly averages for three intakes with 120 s (a) and 400 s (b) sampling periods. Scale parameters of the fitted Cauchy distribution are 0.045 and 0.095, respectively.

In addition to the arithmetic averaging of the data available in the given hour, more sophisticated methods were also applied for the estimation of the hourly averages. In this study, the linear interpolation and the cubic spline interpolation were tested on how much they can reduce the uncertainty of the hourly-average concentrations relative to the common arithmetic averaging. The timestamps of the data for the linear and spline interpolations were the middle of the sampling periods, i.e., not considering the flush time. For cubic spline interpolation, the SPLINE_P procedure of IDL 6.3 (Interactive Data Language, Harris Geospatial Solutions) was applied. Using the cubic spline interpolation, the concentration was reconstructed at 5 s resolution and integrated through the given hour.

3 Results and discussion

Typically, the hourly-average concentration is calculated as the arithmetic mean of the average concentrations of the short periods sampled. The resulting value randomly deviates from the true hourly average due to the unsampled periods. The probability distribution of the deviations is far from Gaussian as it can be seen in Fig. 2. The rather peaked distribution can be approximated by a Cauchy distribution fairly well. As it is reasonable to suppose that the mean deviation is zero, the density of the Cauchy distribution can be written in a simple form:

$$\begin{array}{}\text{(1)}& f(x,\mathit{\gamma })={\displaystyle \frac{\mathit{\gamma }}{\mathit{\pi }\left({\mathit{\gamma }}^{\mathrm{2}}+x^{\mathrm{2}}\right)}},\end{array}$$

where γ is the so-called scale parameter. It gives half of the interquartile range of the distribution. As the distribution is symmetric to zero, γ also gives the median of the absolute values of the deviations. The probability distribution of the deviations has a rather long tail. Therefore, for practical purposes, it is reasonable to know the probability of high deviations, which might significantly increase the uncertainty of the results of atmospheric models. In this paper, the 90-percentile values are used to characterize the extreme values, which still have a non-negligible probability. The Cauchy distribution was fitted to the empirical distribution using the MPFITEXPR function written in IDL (Markwardt, 2009), applying Levenberg–Marquardt least-squares fit.

With increasing length of the sampling period, the uncertainty of the calculated hourly average is also increasing due to the fewer measurements per intake. The increased uncertainty is indicated by the increased γ value. Figure 2 shows the probability distribution of the deviation from the true value for three intakes and two sampling periods: 120 and 400 s. In the first case (120 s), a certain intake is sampled 10 times during an hour, while in the latter case (400 s) it is sampled only three times. The difference is remarkable. The full set of distributions can be found in the Supplement.

It is fairly obvious that there is little difference among the results of arithmetic averaging, linear interpolation, and spline interpolation if the sampling period is short and a given intake is sampled frequently, such that only short periods are missed by the measurements. The advantage of the more sophisticated methods appears when the sampling period is longer and fewer samples are available within the hour. The linear interpolation can better estimate the concentration at the beginning and the end of the hour as it also uses measurements from the previous and subsequent hours. Cubic spline interpolation follows the temporal course of the concentration as closely as possible given the limited number of data.

Figure 3 visualizes the difference between the methods and their results in an example of a typical summer morning. Figure 4 shows the γ values (equivalent to half of the interquartile range or the median of the absolute deviations) and the 90-percentile values of the absolute deviations for the different number of intakes as a function of the sampling time and averaging method (arithmetic average, linear interpolation, cubic spline interpolation) based on the whole dataset.

Figure 3 Difference between the hourly averages estimated by arithmetic averaging, linear interpolation, and spline interpolation compared to the true value in a typical summer morning hour (2 May 2018, 07:00–08:00 LST (local standard time, UTC+1). Gray dots are the measurement data with 5 s temporal resolution, while the thick red lines give the sampling periods and their average concentrations supposing three intakes and 400 s sampling time (60 s flushing + 340 s signal integration at each intake). The height of the planetary boundary layer changed from 399 to 555 m from 07:00 to 08:00 LST.

Figure 4 The median (p50) and the 90-percentile (p90) values of the absolute deviation from the true hourly average as a function of sampling time and number of intakes calculated from the fitted Cauchy distribution.

The activity of the vegetation, uptake, and release of carbon dioxide and the dynamics of the atmosphere dispersing it in the air have significant seasonal and diurnal variations, which also influence the uncertainty of the calculated hourly averages. To get an insight into the temporal variation of the uncertainty, we grouped the calculated deviations from the true values by month and by the time of the day with hourly resolution. At such a resolution, the available data are insufficient in number for a reliable estimation of the scale parameter of the Cauchy distribution. Instead, we calculated the approximate median and 90-percentile of the absolute deviations based on the empirical distributions. The results are rather qualitative than quantitative, but they indicate when modelers should be aware of potentially high uncertainty of the data and what the magnitude can be. Figure 5 shows a color-coded temporal distribution map of the rate of the uncertainty for three intakes and two sampling periods: 120 and 600 s, giving 10 and 2 measurement cycles per hour, respectively. For other sampling periods, numbers of intakes, and averaging methods, the data are presented in the Supplement.

Figure 5 The temporal variation of the 90-percentile of the absolute deviations of the hourly averages calculated as the arithmetic average of the available measurements from the true values for three intakes and for 120 s (a) and 600 s (b) sampling times.

Our analysis attributes quantitative results to the qualitatively more or less expectable findings. As it can be seen in Fig. 4, the shorter the sampling period (the higher the number of the sampling periods within the hour), the lower the uncertainty of the calculated hourly averages, and the more sophisticated methods result in reduced uncertainty. While both the systematic errors derived from, for example, scale bias and the random noise of the gas analyzers are relatively stable within a short time, the uncertainty in the hourly averages caused by the discrete sampling of the continuous signal shows a remarkable temporal variation (Fig. 5). It is lower during the afternoon hours when the atmosphere is fairly well mixed but may be large during the morning transition periods due to the sudden concentration changes in the planetary boundary layer monitored, which cannot be properly followed by episodic sampling. Due to the seasonal variation in both the activity of the vegetation and the dynamics of the atmosphere, the uncertainties are higher in summer than in winter.

In the case of arithmetic averaging, the mean of the absolute values of deviations from the true hourly averages remain below 0.3 µmol mol⁻¹ during the winter (December–February) afternoon hours (the 90-percentile values are <0.8 µmol mol⁻¹) supposing at least two measurement cycles per hour. It is somewhat higher in summer (June–August) but still less than 0.5 µmol mol⁻¹ (the 90-percentile values are <0.9 µmol mol⁻¹). Although the mixing of the atmosphere is more vigorous during summer afternoons, the intensive carbon dioxide uptake by the vegetation causes higher concentration fluctuation in the atmosphere than it experienced in winter.

The uncertainty of the hourly averages derived from episodic samplings is much higher during the morning transition period when the concentration may change by several tens of µmol mol⁻¹ within an hour. The uncertainty is especially high in summer when the diurnal amplitude of the carbon dioxide concentration is the highest due to the intensive photosynthesis/respiration of the vegetation. If at least two measurement cycles per hour are assumed again, the mean uncertainty remains below 0.4 µmol mol⁻¹ during the morning transition hours in winter (the 90-percentile values are <1.5 µmol mol⁻¹) but it may reach 1.7 µmol mol⁻¹ in summer (the 90-percentile values may be as high as 5 µmol mol⁻¹). During these periods, the sampling frequency is critical. The higher the sampling frequency, the better the arithmetic mean mirrors the concentration course and the lower the uncertainty of the estimated hourly average becomes. At three cycles per hour, the mean uncertainty does not exceed 1.2 µmol mol⁻¹ in the data series studied (the 90-percentile values are <3 µmol mol⁻¹), which is further reduced to 0.9 µmol mol⁻¹ (the 90-percentile values are <2.5 µmol mol⁻¹) when at least four cycles were assumed. The evening collapse of the convective boundary layer causes a relatively fast increase in the concentration in the boundary layer but the rate of the change is lower than during the morning build-up; therefore, the resulting uncertainty in the estimated hourly-average concentrations is also lower.

Figure 6 The seasonal variation of the median (p50) and the 90-percentile (p90) values of the absolute deviations of the hourly averages calculated using the different methods (arithmetic averaging, linear interpolation, spline interpolation of the available data) in the case of three intakes and for 120 s (upper panels) and 600 s (lower panels) sampling times. The left panels characterize the morning transition periods (06:00–08:00 LST), while the right panels show the corresponding values for the early afternoon hours (13:00–15:00 LST).

When the sampling time is short, each intake is sampled relatively frequently; therefore, the more sophisticated averaging methods do not give significant improvements in the uncertainty of the calculated hourly averages. Their advantage appears when the sampling frequency is low and/or the concentration changes are significant within the hour, because the interpolation methods follow the actual temporal course of the concentration better than the short-term averages. Although the temporal distribution of the uncertainty does not change, the high values (summer morning periods) are significantly reduced as can be seen in Fig. 6. Theoretically, the cubic spline interpolation follows the course of the concentration more faithfully than the linear interpolation; however, a significant difference cannot be seen in their performance under the conditions studied here.

The full set of the uncertainty maps (different number of intakes, sampling periods, methods) can be found in the Supplement.

It should be emphasized that the numerical results presented here may be highly site-specific. They depend on the signal variability, which in general depends on the height of the sampling elevation above the ground, on the geographical location and environment of the monitoring site, and partly on instrument setup. Sources and sinks of carbon dioxide are located at the surface, and the fluctuation generated by the surface processes gradually attenuates with height (Stull, 1999). The relative role of the high-frequency part of the spectra is reduced at higher elevations, and so the hourly averages can be estimated with lower uncertainty. Under our conditions, the intra-hour trend is the dominant term determining the uncertainty of the hourly-average concentrations, and it might also be true for other sites. The reason is fairly obvious: in the case of a significant trend within the hour, the start time of the first measurement period significantly determines the calculated hourly average. Depending on the start time, the highest or lowest values at the beginning or end of the hour are missed from the averaging. To give an impression of the temporal variation at Hegyhátsál tall-tower site, Fig. 7 shows the diurnal variation of the concentration for July when the diurnal amplitude is the highest. The uncertainty of the hourly averages presented in Fig. 5 is in synchrony with the diurnal variation of the rate of concentration change.

Figure 7 The average diurnal variation of carbon dioxide concentration at Hegyhátsál tall-tower site (82 m above the ground) in July, when the daily amplitude is the highest. The gray zone indicates the 1σ range.

Going farther from the active vegetation, on the top of a mountain, seashore or island, in a desert, or a poorly vegetated region, the concentration fluctuation is also lower than at a mid-continental site in the temperate zone. It means the atmospheric carbon dioxide models face a complex spatiotemporal structure of measurement uncertainty when incorporating data from tall towers applying sequential sampling. The uncertainty of the hourly data available in the databases may be magnitudes higher in certain regions and times than that of the continuous measurements.

4 Conclusions

The accuracy of atmospheric carbon dioxide budget calculations and source/sink allocations essentially depends on the accuracy of the measurements performed at the monitoring stations. Uncertainties derived from scale transfer, scale inconsistency, scale drift, or instrument noise may be assessed (Vermeulen, 2016) and reduced with careful work, measurement intercomparisons, and applying high-quality instrumentation. Our analysis has shown that the uncertainty derived from the non-continuous sampling at the tall-tower sites may be significantly higher than the other terms of the measurement uncertainty.

In the case of continental monitoring sites, the present-day atmospheric inversion models typically use only the early afternoon measurements, the uncertainty of which is the lowest. However, it also means that ∼80 % of the measurements are not used. The progress in the representation of atmospheric dynamics in the models may make it possible to use data from a wider time window (e.g., from late morning till evening), making more measurement data useful. However, the wider time window also means that data with higher uncertainty also penetrate the model.

There is a good reason to suppose that the results presented in this paper are site-specific, and the uncertainties of the hourly-average values in the databases vary in both time and space. Handling this variable uncertainty in the models and assessing its consequences on the results may be challenging. The data users should be aware of the spatially and temporally variable uncertainty of the measurement data they use. Metadata on sampling frequency and integration time, as well as more uncertainty studies, may help their work. It might also be reasonable to store the data with high temporal resolution (e.g., minute values) in public databases.

Because of the presumable site-specificity of the results, no general recommendation can be given for the measurement strategy. The trade-off between the high-frequency sampling to follow the atmospheric changes precisely and the long integration time to get proper flushing and to reduce the instrument noise must be found locally knowing the local environmental conditions and the characteristics of the given monitoring system.