Multivariate statistical air mass classification for the high-alpine observatory at the Zugspitze Mountain, Germany

To assist atmospheric monitoring at high-alpine sites, a statistical approach for distinguishing between the dominant air masses was developed. This approach was based on a principal component analysis using five gas-phase and two meteorological variables. The analysis focused on the Schneefernerhaus site at Zugspitze Mountain, Germany. The investigated year was divided into 2-month periods, for which the analysis was repeated. Using the 33.3 % and 66.6 % percentiles of the first two principal components, nine air mass regimes were defined. These regimes were interpreted with respect to vertical transport and assigned to the BL (recent contact with the boundary layer), UFT/SIN (undisturbed free troposphere or stratospheric intrusion), and HYBRID (influences of both the boundary layer and the free troposphere or ambiguous) air mass classes. The input data were available for 78 % of the investigated year. BL accounted for 31 % of the cases with similar frequencies in all seasons. UFT/SIN comprised 14 % of the cases but was not found from April to July. HYBRID (55 %) mostly exhibited intermediate characteristics, whereby 17 % of the HYBRID class suggested an influence from the marine boundary layer or the lower free troposphere. The statistical approach was compared to a mechanistic approach using the ceilometer-based mixing layer height from a nearby valley site and a detection scheme for thermally induced mountain winds. Due to data gaps, only 25 % of the cases could be classified using the mechanistic approach. Both approaches agreed well, except in the rare cases of thermally induced uplift. The statistical approach is a promising step towards a real-time classification of air masses. Future work is necessary to assess the uncertainty arising from the standardization of real-time data.

1 Introduction

High-alpine observatories such as the Environmental Research Station Schneefernerhaus (UFS) at Zugspitze Mountain, Germany, play an important role in studying changing concentrations of atmospheric constituents such as greenhouse gases, aerosols, and persistent organic pollutants (POPs), which have critical impacts on the climate, environmental integrity, or human health (McClure et al., 2016; Kirchner et al., 2016). In particular, high-alpine observatories frequently offer the opportunity to sample well-mixed air masses of the free troposphere. These air masses are representative for large spatial areas and, thus, are suitable for the determination of large-scale and global trends (Yuan et al., 2019). Therefore, the observational network of the Global Atmosphere Watch Program of the World Meteorological Organization includes many high-alpine observatories. At times, however, high-alpine sites can be affected by local anthropogenic emissions on the mountain or regional emissions if air masses of the boundary layer (BL) are lifted or mixed with the free troposphere via processes such as synoptic lifting (e.g., at fronts), thermally induced anabatic winds, and foehn flows (Zellweger et al., 2003). Foehn winds descend on the lee side of a mountain range but can be associated with air mass lifting on the windward side and with enhanced turbulent mixing in the lee. In relatively rare cases, the air masses at high-alpine sites originate from the stratosphere, e.g., on approximately 6 % of the days at Zugspitze (Stohl et al., 2000).

The BL can be defined as “the part of the troposphere that is directly influenced by the presence of the earth's surface, and responds to surface forcings with a timescale of about an hour or less” (Stull, 1988). Reuten et al. (2007) proposed using a timescale of one diurnal cycle for this definition to explicitly include the residual layer above a stably stratified nocturnal boundary layer and elevated aerosol layers that can result from thermally driven upslope flows (Gohm et al., 2009). For instance, on typical fair weather days in summer, Henne et al. (2004) observed a two-layer structure of the BL in deep Alpine valleys where upslope flows lifted air from a polluted lower layer to a moderately polluted injection layer that reached well above the crest height.

In air pollution studies, the terms mixing layer or mixed layer are more common than BL and have been defined as “the layer adjacent to the ground over which pollutants or any constituents emitted within this layer or entrained into it become vertically dispersed by convection or mechanical turbulence within a time scale of about an hour” (Seibert et al., 2000). Over flat terrain, the mixing layer typically coincides with the BL, at least during daytime. Following Reuten et al. (2007), we use a timescale of one diurnal cycle and consider a residual layer and elevated aerosol layers as parts of the BL and the mixing layer.

Basically, there are two different approaches for air mass classification at high-alpine sites. Mechanistic approaches directly investigate atmospheric transport processes using meteorological measurements or trajectory models, whereas statistical approaches infer the influences of the air layers from the air mass composition at the site and, possibly, meteorological auxiliary data.

Zellweger et al. (2003) used a mixed mechanistic and statistical approach at the Jungfraujoch high-alpine site, Switzerland. They identified three transport processes, namely foehn events, synoptical lifting, and thermally induced uplift, and attributed the other cases to the undisturbed free troposphere. Foehn events and synoptic lifting were detected with mechanistic approaches using standard meteorological measurements and the height of back trajectories, respectively. Thermally induced uplift was determined using a statistical approach based on the diurnal variation of the sum of oxidized nitrogen species (NO_y) or, in the case of data gaps, the aerosol surface area concentration or specific humidity.

In a recent study, Hulin et al. (2019) explored different methods for selecting days with thermally driven circulations influencing the Pic du Midi high-alpine observatory in the French Pyrenees. A statistical method based on the diurnal cycle of specific humidity seemed to be most reliable. Wind profiler measurements above the nearby plain indicated that a plain–mountain circulation with an upper-level return flow played an important role, although this circulation was sometimes not visible in the near-surface in situ wind measurements on the plain.

Other studies used a ground-based lidar or ceilometer that was installed near a high-alpine site at a lower altitude to determine the mixing layer height (MLH) from the vertical aerosol backscatter profile (Gallagher et al., 2012; Ketterer et al., 2014). This method is considered to be a mechanistic approach because the determination of the MLH from vertical profiles of measured quantities requires a priori meteorological knowledge. Recently, Poltera et al. (2017) used a tilted configuration of a ceilometer and demonstrated that Jungfraujoch was rarely embedded in the local convective boundary layer but much more frequently in an injection layer located above, with slightly higher aerosol concentrations compared with the free troposphere.

Some statistical approaches defined a threshold for the concentration of a surface-emitted atmospheric constituent or a ratio of constituents. At the Mount Bachelor Observatory, USA, a seasonal or monthly threshold for the water vapor mixing ratio was used to distinguish between free-tropospheric- and BL-influenced air masses (Ambrose et al., 2011; Zhang and Jaffe, 2017). The threshold was chosen such that the water vapor mixing ratios below this threshold had the same seasonal or monthly mean as the data from National Weather Service soundings that were launched at a lower-elevation site.

Figure 1 Location of the measurement sites in the surrounding of Zugspitze. North is at the top.

For Jungfraujoch, Herrmann et al. (2015) compared a mechanistic approach, which was based on back trajectories and an inventory of carbon monoxide (CO) emissions, with two simple statistical approaches that used a constant threshold of the radon-222 concentration and the ratio of CO to NO_y. The CO∕NO_y threshold appeared to achieve the best distinction between free-tropospheric- and BL-influenced air masses but the authors noted that a single threshold cannot account for varying degrees of BL influence (Herrmann et al., 2015).

Other statistical approaches aimed at selecting baseline (also called background) concentrations of a trace gas at different kinds of remote sites and were based on outlier removal techniques (e.g., Ruckstuhl et al., 2012) or the fact that well-mixed air masses result in a small temporal variability of the trace gas mixing ratio (e.g., Yuan et al., 2018). Baseline concentrations refer to a given species in a well-mixed air mass with minimal influence from anthropogenic impurities with a relatively short lifetime (Calvert, 1990; Yuan et al., 2018), which are conditions associated with the free troposphere at high-alpine sites.

In contrast to mechanistic approaches, statistical air mass classifications are not able to distinguish between different uplift processes but require only local data. So far, however, statistical approaches have only been based on a single constituent or a ratio of constituents, although several atmospheric constituents are typically monitored at high-alpine observatories.

This study proposes a novel statistical approach based on a principal component analysis (PCA) using seven chemical and meteorological variables. This approach is intended for later use in real-time operational mode to enable an automated sampling of ambient air with respect to different air masses using a multi-channel sampling system for the monitoring of POPs (Kirchner et al., 2016). The objectives were to (i) develop a statistical classification scheme for the UFS site at Zugspitze and to (ii) validate this approach, as far as possible, using a mechanistic approach based on ceilometer and standard meteorological measurements.

2 Methods

2.1 Measurement sites

The UFS site (47^∘25^′00^′′ N, 10^∘58^′47^′′ E, 2650 m a.s.l.) is located on a steep south-facing slope, approximately 300 m below the summit of Zugspitze (2962 m a.s.l.), which is the highest mountain in the German Alps and represents the first real barrier for northwesterly advection from the Alpine foreland. At UFS, westerly and easterly wind directions dominate due to the local topography (Risius et al., 2015). Owing to the trace gas and aerosol measurements that take place at the site, the UFS has reduced its emissions of these substances to a minimum. Nevertheless, the measurements can be influenced by local emissions from the surrounding area, which is a highly frequented tourist area all year round. Local emission sources include nearby cable car stations at the Zugspitze Summit (ZSG) and Zugspitzplatt (ZPLT), which is a gently sloping plateau below the UFS site (Fig. 1). Additionally, a skiing area is situated at ZPLT. The large metropolitan area of Munich is also approximately 90 km north of the study site.

Table 1 Chemical measurements used in this study. The given interval is the temporal interval at which the data were provided by the respective institutions. UV denotes ultraviolet light. The symbols for the measured quantities are explained in the main text.

^a Refers to 30 min in the year 2014; ^b see Birmili et al. (2016) and Sun et al. (2018) for details on the scanning mobility particle sizer (SMPS) and the multi-angle absorption photometer (MAAP); ^c see Leuchner et al. (2016) for details; ^d see Steinkopff et al. (2012) for details; ^e German Environment Agency; ^f Ecoclimatology, Technical University of Munich; ^g German Meteorological Service; ^h Experimental Aerosol and Cloud Microphysics, Leibniz Institute for Tropospheric Research.

In addition to UFS, seven weather stations at different altitudes at a maximum horizontal distance of 11 km from UFS were available and were included in this study (Fig. 1). One of these stations is located at ZSG. Another weather station is located on the ZPLT and is surrounded by mountain ridges except to the east where the plateau leads to the narrow, deep Reintal Valley. The Schachen (1830 m a.s.l.), Kreuzalm (1600 m a.s.l.), Felsenkanzel (1250 m a.s.l.), and Brandwiese (900 m a.s.l.) weather stations are part of the KLIMAGRAD project (Schuster et al., 2014). The Schachen site is situated on a plateau above the middle Reintal in close vicinity of some trees. Kreuzalm is situated in a meadow on a mountain saddle. Brandwiese is located in a meadow surrounded by forest, where a tributary valley reaches the Reintal from the west. At Felsenkanzel, the measurements are made on a steep south-facing slope northwest of the town of Garmisch-Partenkirchen (GAP). The GAP site (720 m a.s.l.) is located on the western periphery of the town where a broad west–east oriented valley turns northeast towards the alpine foreland (Fig. 1).

2.2 Instrumentation and data set

The data set that was analyzed in this study spans a period of 1 year from 22 August 2013 to 21 August 2014, which was selected due to high data availability. Table 1 gives an overview of the chemical measurements, associated instruments, measurement principles, and research institutions that provided the data. Most of the atmospheric constituents are measured at UFS, including carbon mono- (CO) and dioxide (CO₂), methane (CH₄), ozone (O₃), the sum of oxidized nitrogen species (NO_y), nitrogen oxides (${\mathrm{NO}}_x=\mathrm{NO}+{\mathrm{NO}}_{\mathrm{2}}$ with NO and NO₂ referring to nitrogen mono- and dioxide, respectively), formaldehyde (HCHO), the ambient particle number size distribution (dN(dlog d_p)⁻¹) for particle diameters (d_p) between 10 and 600 nm, the mass concentration of particulate matter with d_p<10 µm (PM₁₀), and equivalent black carbon (eBC). The radioisotopes beryllium-7 (⁷Be) and radon-222 (²²²Rn) are sampled at ZSG and on the mountain ridge (2825 m a.s.l.) directly above UFS, respectively. At the remaining sites, only meteorological parameters are recorded.

While most of the chemical data were available as 1 min averages, HCHO and the particle size distribution were provided at 10 min intervals and ²²²Rn and ⁷Be were available as 2 h and 12 h averages, respectively (Table 1). The HCHO data contained a large data gap between 15 December 2013 and 15 July 2014 and the ²²²Rn data were only available from 1 January 2014. The ²²²Rn data were downloaded from the World Data Centre for Greenhouse Gases (WDCGG, 2019). Because ⁷Be attaches to aerosol particles, it is measured by collecting the carrier aerosol using a glass-fiber filter. While NO can be measured directly, NO_y and NO₂ have to be converted to NO before detection. This conversion of NO_y and NO₂ is performed using photolysis and gold/CO converters, respectively.

The standard meteorological measurements included air temperature (T), relative humidity (RH), global radiation (R_g), and horizontal wind velocity and direction at all sites. Air pressure (p) was available at the UFS, ZSG, and GAP sites. Additionally, year-round precipitation measurements from an electronic weighing system and a windbreak ring (Sommer Messtechnik, Koblach, Austria) were used at the ZPLT site. At the ZPLT and the four KLIMAGRAD sites, the wind data are measured using propeller anemometers (wind monitor model 05103, Young, Traverse City, USA), whereas at GAP and ZSG, cup anemometers and wind vanes (SK-565 and SK-566, respectively, Thies Clima, Göttingen, Germany) are used. At UFS, an ultrasonic anemometer (model 2-D, Thies Clima, Göttingen, Germany) is used. The standard meteorological data were available as 10 min averages except at UFS (1 min averages).

At GAP, a ceilometer (CHM 15k, Lufft, Fellbach, Germany) provides 15 s averages of up to three aerosol layer heights (ALHs) with a vertical resolution of 15 m using a wavelet algorithm from the manufacturer, which detects strong gradients in the range-corrected attenuated backscatter profile. In the absence of clouds, the signal-to-noise ratio of the ceilometer is typically > 1 up to a height of 4 to 5 km a.g.l. during daytime and up to greater heights at night (Heese et al., 2010).

2.3 Data post-processing and quality control

All time stamps were converted to local standard time. The data were aggregated on a 30 min basis while requiring a data availability of ≥ 66 % in each interval unless stated otherwise. The ²²²Rn and ⁷Be concentrations were interpolated using nearest neighbor interpolation to match them with the 30 min intervals used for the analysis.

2.3.1 Atmospheric constituents

The chemical measurements were quality controlled by the research institutions that provided the data. Invalid data, which resulted from calibration, instrument repair, or power failure, were discarded.

Due to measurement uncertainties, NO_y occasionally exhibited a lower mixing ratio than NO_x. If NO_y was lower than 75 % of NO_x and NO_x was at least 0.03 ppb, both quantities were treated as artifacts and were discarded. If NO_y still remained lower than NO_x after averaging on a 30 min basis, both quantities were assumed to be equal and were replaced by the mean of NO_x and NO_y. Negative eBC concentrations were treated in a similar way. Values below −0.05 µg m⁻³ were discarded. If eBC still remained negative after averaging on a 30 min basis, the concentration was set to zero.

The processing of the particle number size distributions is described in Birmili et al. (2016) and includes, among others, a multiple charge inversion and corrections for particle losses. Following Herrmann et al. (2015), the number concentration of accumulation-mode particles (N₉₀) was used as an indicator of BL air masses and approximated by including particle diameters between 90 and 600 nm.

Data of most trace gases and PM₁₀ were flagged manually with respect to locally polluted air masses. These air masses are a special kind of BL air mass because they do not result from uplift processes but from human activities on the mountain. Local emissions were evident from short but pronounced peaks in trace gas and aerosol concentrations and were most frequently observed for NO_x and NO_y. The flag for local NO_x emissions was approximately reproduced by selecting data with a 30 min NO_x standard deviation (${\mathit{\sigma }}_{{\mathrm{NO}}_x}$) of > 0.4 ppb (not shown). According to this criterion, local pollution events affected 9 % of the NO_x data but half of these events lasted no longer than one time interval, i.e., 30 min. Consequently, local pollution events are generally so short that they cannot be predicted on a 30 min basis in real-time applications. Owing to this, in addition to the fact that the classification scheme was developed for later use in real-time operational mode, the air mass classification does not account for local pollution events.

2.3.2 Standard meteorological data

The wind data at ZSG and GAP were discarded for wind velocities below the starting threshold of the wind vane of 0.3 m s⁻¹ for a displacement of 90^∘ (Löffler, 2012). For the wind vane of the propeller anemometers at the ZPLT and the KLIMAGRAD sites, the starting threshold was 1.1 m s⁻¹ (Young, 2018). To achieve a trade-off between high data quality and high data availability, the wind data of these anemometers were discarded for wind velocities of < 0.7 m s⁻¹.

The standard meteorological data were aggregated on a 30 min basis using the vector mean for wind direction, the sum for precipitation, and the arithmetic mean for the other variables. Specific humidity q (g kg⁻¹) and virtual potential temperature θ_v (K) were calculated for each site as follows: (Foken, 2008),

$$\begin{array}{}\text{(1)}& {\displaystyle }& {\displaystyle }q=\mathrm{0.622}{\displaystyle \frac{e}{p-\mathrm{0.378}e}}\text{(2)}& {\displaystyle }& {\displaystyle }{\mathit{\theta }}_{\text{v}}=(\mathrm{1}+\mathrm{0.608}q)T{\left({\displaystyle \frac{p_{\mathrm{0}}}{p}}\right)}^{R_{\text{L}}/c_{\text{p}}},\end{array}$$

where e (hPa) is vapor pressure, p (hPa) is air pressure, p₀ is 1000 hPa, T (K) is air temperature, R_L (J kg⁻¹ K⁻¹) is the gas constant of dry air, and c_p (J kg⁻¹ K⁻¹) is specific heat capacity of dry air. The vapor pressure was calculated from relative humidity (RH) and T using an approximation of the Clausius–Clapeyron equation (e.g., Wallace and Hobbs, 2006):

$$\begin{array}{}\text{(3)}& e_{\text{s}}\left(T\right)=e_{\text{s}}\left(T_{\mathrm{0}}\right)\exp \left\{{\displaystyle \frac{L_{\text{v}}}{R_{\text{v}}}}\left({\displaystyle \frac{\mathrm{1}}{T_{\mathrm{0}}}}-{\displaystyle \frac{\mathrm{1}}{T}}\right)\right\},\end{array}$$

where e_s (hPa) is saturation vapor pressure, T₀ (K) is a reference temperature with known e_s, L_v (J kg⁻¹) is latent heat of evaporation, and R_v (J kg⁻¹ K⁻¹) is the gas constant for water vapor. Equation (3) assumes that the specific volume of liquid water is negligible compared with that of water vapor and that the latent heat of vaporization is independent of temperature, which is expected to cause negligible errors with regards to the purpose of this study.

Because p was only measured at the ZSG, UFS, and GAP sites, it was calculated for the other sites using the hypsometric equation:

$$\begin{array}{}\text{(4)}& p=p_{\text{ref}}\exp \left\{-{\displaystyle \frac{g_{\mathrm{0}}(Z-Z_{\text{ref}})}{R_{\text{L}}{\stackrel{\mathrm{‾}}{T}}_{\text{v}}}}\right\},\end{array}$$

where p (hPa) and p_ref (hPa) are the air pressures at the site of interest and a reference site, respectively; Z (m) and Z_ref (m) are the geopotential heights that were approximated using the altitudes of the respective sites; g₀ is 9.8 m s⁻²; R_L (J kg⁻¹ K⁻¹) is the gas constant of dry air; and ${\stackrel{\mathrm{‾}}{T}}_{\text{v}}$ (K) is the mean virtual temperature for the air layer between the two sites under consideration. p was computed in two steps. First, p was approximated using the dry-bulb temperature instead of ${\stackrel{\mathrm{‾}}{T}}_{\text{v}}$. Second, ${\stackrel{\mathrm{‾}}{T}}_{\text{v}}$ was calculated with the approximated p and used to recalculate p. The UFS site was used as the reference site to calculate p at the next lower site (ZPLT), which was then used as the reference site to calculate p at the next lower site (Schachen). This procedure was continued until Kreuzalm. Similarly, GAP was used as the reference site to determine p at the next higher site and this procedure was continued until Schachen. At Kreuzalm and Schachen, the two p estimates based on the next lower and next higher sites were averaged.

In order characterize the static stability of the valley atmosphere with a single quantity, the range of the pseudo-vertical profile of virtual potential temperature, Δθ_v (K), was defined as follows (considering all stations except Schachen):

$$\begin{array}{}\text{(5)}& \mathrm{\Delta }{\mathit{\theta }}_{\text{v}}=\left\{\begin{array}{ll}{\mathit{\theta }}_{\text{v}}^{\max }-{\mathit{\theta }}_{\text{v}}^{\min },& \text{if}{z}_{{\mathit{\theta }}_{\text{v}}^{\min }}\le z_{{\mathit{\theta }}_{\text{v}}^{\max }}\\ {\mathit{\theta }}_{\text{v}}^{\min }-{\mathit{\theta }}_{\text{v}}^{\max },& \text{if}{z}_{{\mathit{\theta }}_{\text{v}}^{\min }}>z_{{\mathit{\theta }}_{\text{v}}^{\max }},\end{array}\right.\end{array}$$

where ${\mathit{\theta }}_{\text{v}}^{\min }$ (K) and ${\mathit{\theta }}_{\text{v}}^{\max }$ (K) are the minimum and maximum virtual potential temperature of the pseudo-vertical profile, respectively, and $z_{{\mathit{\theta }}_{\text{v}}^{\min }}$ (m) and $z_{{\mathit{\theta }}_{\text{v}}^{\max }}$ (m) are the associated altitudes. Schachen was excluded from this calculation because θ_v was particularly high at Schachen if the global radiation was high and the wind velocity was low, which suggested that the naturally aspirated thermometer was affected by radiative errors (not shown). ${\mathit{\theta }}_{\text{v}}^{\min }$ was mostly found at the low-elevation sites GAP and Brandwiese, which applied to 58 % and 37 % of the cases, respectively. ${\mathit{\theta }}_{\text{v}}^{\max }$ mostly occurred at the high-elevation sites ZSG, UFS, or ZPLT, which applied to 68 %, 18 %, and 13 % of the cases, respectively. The values of Δθ_v ranged between −5 and +31 K, where positive values indicate stable conditions (see Fig. S1a in the Supplement).

Figure 2 Post-processing of the ceilometer-based aerosol layer heights (ALHs), which were used to determine the mixing layer height (MLH) at GAP, illustrated for the period from 17 to 21 July 2014. Outliers have been removed from the MLH. Points of the 30 min median MLH are connected by a line if separated by a maximum of 3 h.

2.3.3 Ceilometer data

To determine the MLH at GAP from the ALH retrievals of the ceilometer, several post-processing steps were necessary (Fig. 2). ALHs below 600 m a.g.l., corresponding to 1320 m a.s.l., were regarded as artifacts and were discarded because of high measurement uncertainties in this part of the backscatter profile resulting from an incomplete overlap of the laser beam and the field of view of the ceilometer (Flentje et al., 2010). The backscatter signal only reflects aerosol concentrations below clouds. Thus, ALHs that were greater than or equal to the cloud base height, which was measured by the same instrument, were discarded. For this reason, it was not possible to determine the correct MLH if clouds were present at the top of or within the mixing layer. These first two post-processing steps strongly reduced the data availability from 89 % to 45 % of the time.

If two or three ALHs were found at the same time, the following procedure was applied. For each ALH, the number of ALHs in its vicinity, i.e., within centered intervals of 300 m and 30 min in the vertical and temporal dimensions, respectively, was determined. The ALH with the maximum number of ALHs in its vicinity was attributed to the MLH, while the other ALHs were discarded. If each of the ALHs had less than three ALHs in its vicinity, a missing value was registered. The MLH attribution, as described above, was based on the idea that the MLH only varies gradually; thus, similar MLHs are expected shortly before and after the time under consideration. In particular situations, in which elevated aerosol layers enter or leave the field of view of the ceilometer, this assumption may not be valid.

From the resulting time series at 15 s intervals, outliers were removed iteratively as follows. From each pair of subsequent data points with an absolute rate of change exceeding 240 m min⁻¹, one data point was removed so that the standard deviation of the MLH in the preceding 15 min and the following 15 min was most strongly reduced. This procedure was repeated until the absolute rate of change was a maximum of 240 m min⁻¹, corresponding to four vertical averaging intervals (4×15 m) per time interval (15 s).

The MLHs were aggregated using the 30 min median while requiring a data availability of ≥ 50 % and a standard deviation of ≤ 170 m, corresponding to the 90 % percentile of all 30 min standard deviations. From the aggregated time series, outliers were removed in the same way as described above but with a maximum allowed absolute rate of change of 400 m h⁻¹. Finally, the MLH was available 34 % of the time.

2.4 Statistical classification approach

2.4.1 Definition of the air mass classes

The statistical analysis aimed at distinguishing the following three classes of air masses: (a)  BL – air masses with recent contact with the BL, characterized by increased concentrations of surface-emitted atmospheric constituents suggesting a recent uplift with a dominating influence of regional or local sources; (b) UFT/SIN – air masses of the undisturbed free troposphere or stratospheric intrusions, characterized by very low concentrations of surface-emitted constituents suggesting recent subsidence or horizontal advection; and (c) HYBRID – air masses that are influenced by both the BL and the free troposphere or that exhibit ambiguous characteristics. The HYBRID class was anticipated, for example, if the air mass had been exported from the BL to the free troposphere at a significant horizontal distance from the UFS and had been mixed with the free troposphere during transport, resulting in intermediate concentrations of surface-emitted constituents.

Stratospheric intrusions are characterized by high ⁷Be and O₃ concentrations and a low humidity and were observed at Zugspitze on approximately 5 % of the days (Stohl et al., 2000). Because of this low frequency and the coarse temporal resolution of ⁷Be of 12 h, stratospheric intrusions were not classified individually but were attributed to the same class as air masses of the undisturbed free troposphere, which exhibit similar characteristics.

2.4.2 Principal component analysis

A PCA allows for the reduction of the number of dimensions of a data set while maintaining as much variance as possible. Principal components (PCs) are uncorrelated standardized linear combinations of the original variables (Mardia et al., 1979). If the variance of the original variables is mainly caused by shifts between the air mass classes, the first few PCs will be suitable indicators for air mass classification. The PCA was performed separately for 2-month periods (i) to account for seasonally changing relationships between the variables and (ii) to largely eliminate the influence of the seasonal variability of atmospheric constituents because the seasonality reflects not only the frequency of BL-influenced air masses but also other factors such as the source and sink strength (e.g., residential heating in winter, photosynthetic CO₂ removal in summer), chemical reactivity, and deposition. The original variables were standardized by subtracting the arithmetic mean and dividing by the standard deviation of the respective 2-month period in order to avoid biased results due to different units and ranges of the data. The subtraction of the 2-month mean roughly removes the seasonal cycle.

Figure 3 (a) Scheme for statistical air mass classification at UFS using the first (PC1) and second principal components (PC2) of a 6- or 7-dimensional data set. The mapping of air mass regimes to air mass classes is only shown for February and March and was different for other 2-month periods. BL denotes recent contact with the boundary layer, UFT/SIN includes the undisturbed free troposphere and stratospheric intrusions, and HYBRID reflects influences from both BL and UFT/SIN and includes the LRMD (long-range transport of mineral dust) and MBL/UFT (influence of the marine boundary layer or UFT) subclasses. (b) Numbering of the air mass regimes.

Only variables with an expected unambiguous link to vertical transport processes and a high data availability were used as input variables for the PCA. These variables include the CO, CH₄, CO₂, and O₃ gases, and water vapor. Among the meteorological variables, the air pressure at GAP and Δθ_v were considered to be most suitable for the PCA because they are physically linked to uplift (low pressure, low static stability) and subsidence processes (high pressure, high static stability). The remaining variables were only used for validation purposes. Aerosol measurements were excluded from the PCA because low aerosol concentrations do not necessarily indicate UFT/SIN air masses and can also result from wet deposition during the uplift of BL air masses. ²²²Rn and HCHO were only used for validation purposes due to low data availability. NO_y and NO_x were excluded from the PCA because their variability was particularly strongly influenced by local emissions, which are not indicative of uplift processes. O₃ was only used as a PCA input variable in the winter half year because high O₃ mixing ratios are not only caused by a subsidence of UFT/SIN air masses but can also result from photochemical O₃ production in BL air masses, which typically contain high precursor concentrations. In the winter half year, photochemical O₃ production was assumed to play a minor role at high-alpine sites because of the generally low solar irradiance and the weak thermally induced uplift. This argument is supported by the trajectory residence time statistics of Kaiser et al. (2007), which demonstrated that the O₃ mixing ratio at European high-alpine sites was generally lower in winter if the air mass originated from lower altitudes, whereas in summer, this was often not the case.

The 2-month periods were defined as December to January, February to March, and so on, such that the winter and summer half years included the months with the lowest and highest solar forcing, respectively. The PCs were computed from the following five gas-phase variables and two meteorological variables, where O₃ was only included in the 2-month periods of the winter half year (October–March):

$$\begin{array}{}\text{(6)}& \begin{array}{rl}{\mathrm{PC}}_i=& a_{\mathrm{1}i}\widehat{\left[\mathrm{CO}\right]}+a_{\mathrm{2}i}\widehat{\left[{\mathrm{CH}}_{\mathrm{4}}\right]}+a_{\mathrm{3}i}\widehat{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+a_{\mathrm{4}i}\widehat{\left[{\mathrm{O}}_{\mathrm{3}}\right]}\\ & +a_{\mathrm{5}i}\widehat{q}+a_{\mathrm{6}i}\widehat{\mathrm{\Delta }{\mathit{\theta }}_{\text{v}}}+a_{\mathrm{7}i}\widehat{p_{\mathrm{GAP}}}.\end{array}\end{array}$$

Here, PC_i represents the scores of the ith PC; $\widehat{\left[\mathrm{CO}\right]},\widehat{\left[{\mathrm{CH}}_{\mathrm{4}}\right]},\widehat{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$, and $\widehat{\left[{\mathrm{O}}_{\mathrm{3}}\right]}$ are the standardized mixing ratios of the respective trace gases at UFS; $\widehat{q}$ is the standardized specific humidity at UFS; $\widehat{\mathrm{\Delta }{\mathit{\theta }}_{\text{v}}}$ is the standardized range of the pseudo-vertical profile of virtual potential temperature; $\widehat{p_{\mathrm{GAP}}}$ is the standardized air pressure at GAP; and $a_{\mathrm{1}i},\mathrm{\dots },a_{\mathrm{7}i}$ are the loadings of the ith PC. All quantities in Eq. (6) are dimensionless due to standardization. The loadings represent standardized eigenvectors of the correlation matrix of the original variables, while the associated eigenvalues correspond to the variances of the PC scores (Mardia et al., 1979).

2.4.3 Isolating regimes and classes of air masses

In each 2-month period, the air masses were divided into nine regimes using the 33.3 % and 66.6 % percentiles of the first two PCs as thresholds (Fig. 3). The division of each of the two PCs into three equally large subsets follows the same logic as the definition of the three air mass regimes. ML and UFT/SIN air masses are expected on opposite edges of the distributions of the PC1 and PC2 scores while a mixture of these air masses may be found in the center of the distributions. The nine regimes were interpreted with respect to vertical transport and attributed to the three air mass classes, BL, UFT/SIN, and HYBRID, by visually comparing summary statistics of the PCA input variables between the regimes. In the winter half year, the regimes were assigned to the BL class if CO, CH₄, CO₂, and q were relatively high and Δθ_v, p_GAP and O₃ were relatively low compared with the other regimes. The UFT/SIN class was assigned in the opposite case, and the remaining regimes were assigned to the HYBRID class. Apart from CO₂ and O₃, the same criteria were used in the summer half year. Due to photosynthesis, CO₂ was required to be relatively low for BL air masses compared with other air masses during the summer half year. O₃ was not used as an input variable in the summer half year.

This subjective mapping of regimes to classes of air masses was validated using summary statistics of the remaining chemical (e.g., NO_y, NO_x, ²²²Rn) and standard meteorological measurements (e.g., precipitation, relative humidity). Long-range transport of mineral dust (LRMD) was regarded as a subclass of HYBRID and was identified using the following three criteria: (i) a high PM₁₀ concentration of ≥13 µg m⁻³, (ii) a relatively low 30 min standard deviation of the PM₁₀ of ≤ 0.2 PM₁₀ to avoid the attribution of local pollution events, and (iii) an air mass regime that does not suggest a current uplift of BL air masses. The PM₁₀ threshold of 13 µg m⁻³ was motivated by a local minimum in the PM₁₀ histogram for non-BL air masses in June–July (not shown) and the plausibility was checked by visual inspection of the whole PM₁₀ time series. Ambiguous air mass regimes, which carried the fingerprint of the marine boundary layer or the lower free troposphere (among others, very low CO and CH₄ mixing ratios but high q and low O₃ mixing ratios) were attributed to MBL/UFT, another subclass of HYBRID, unless the criteria for LRMD were fulfilled (Fig. 3). Balzani Lööv et al. (2008) demonstrated that characteristics of the marine boundary layer could be observed in the Alps if the air mass was lifted above the Atlantic Ocean and transported to the Alps within several days.

Figure 4 Wind rose plots showing the frequency of counts by wind direction and velocity (v) at the GAP (a, c) and ZPLT (b, d) sites during day- (a, b) and nighttime (c, d) in summer (June to August). Daytime was defined by a global radiation of > 5 W m⁻² at a minimum of one of the sites used in this study; all other cases were treated as nighttime. The plots were used to define wind direction sectors for thermally induced anabatic and katabatic winds. The percentage in brackets specifies the data availability.

2.5 Mechanistic classification approach

In search of criteria for thermally induced anabatic and katabatic winds that influence the UFS site, the wind directions of all sites were investigated. The wind direction at UFS was not useful because upslope and downslope flows were not clearly visible but appeared to be superimposed by synoptically driven winds due to the relatively exposed location of the UFS site (not shown). The wind directions only exhibited pronounced diurnal patterns that indicated thermally induced mountain winds at GAP, Felsenkanzel, and ZPLT, especially in summer (Fig. 4). However, the wind data from Felsenkanzel were not used because of problems with shifts in the sensor orientation. Up-valley and down-valley winds were characterized by northeasterly and southwesterly wind direction sectors at GAP, respectively, and by easterly and west-southwesterly wind direction sectors at ZPLT, respectively (Fig. 4). The wind patterns at ZPLT were consistent with the study of Gantner et al. (2003), which reported a thermal circulation above the same plateau with an easterly inflow of a few hundred meters depth at daytime and a westerly outflow at night under fair weather conditions in summer. Modeling suggested that the inflow ascended the narrow, deep Reintal Valley before reaching ZPLT.

Table 2 Criteria for thermally induced anabatic and katabatic winds occurring simultaneously at GAP and ZPLT. φ is wind direction, v is wind velocity, and Δθ_v is the range of the pseudo-vertical profile of virtual potential temperature.

Figure 5 Percentage of explained variance as a function of the principal component (PC) number for each 2-month period: (a) winter half year and (b) summer half year.

The wind direction sectors at ZPLT and GAP were used as a criterion for anabatic and katabatic winds (Table 2). Additionally, anabatic and katabatic winds were restricted to cases with a weak static stability (Δθ_v<8 K) and a strong static stability (Δθ_v≥8 K), respectively, which was justified by the intersect between the frequency distributions of Δθ_v for the wind direction sectors described above (see Fig. S1b in the Supplement). For katabatic winds, a low wind velocity of <3 m s⁻¹ was required at ZPLT and GAP to exclude cases with a strong synoptic forcing. Additionally, both anabatic and katabatic winds were required to persist for at least 1 h, i.e., two time intervals (Table 2).

The detection of these winds and the ceilometer-based MLH at GAP were combined using a mechanistic approach that primarily accounted for thermally induced vertical transport. For this purpose, the MLH_GAP was discarded if clouds were detected at a height of < 4 km a.s.l. because low-level clouds increase the uncertainty of the MLH retrieval. The mechanistic approach distinguished between the following three conditions: (a) anabatic winds occur or the UFS site is below MLH_GAP, which suggests BL air masses; (b) katabatic winds occur and the UFS site is above MLH_GAP, which suggests UFT/SIN or HYBRID air masses; (c) the UFS site is above MLH_GAP and the winds are not thermally induced, which also suggests UFT/SIN or HYBRID air masses under the assumption of identical MLHs at UFS and GAP. This approach was limited by the low availability of the MLH (22 % after discarding periods with low-level clouds) but still allowed for a partial verification of the statistical approach.

3 Results and discussion

3.1 Principal components and their seasonal dependence

The PCA converts several input variables into the same number of PCs while the variance and, thus, the importance of the PCs is highest for the first few PCs. In the winter half year, the percentage of explained variance of the PCs decreased strongly from PC1 to PC2 and more slightly towards higher-order PCs. In the summer half year, the explained variance decreased approximately linearly with increasing PC number (Fig. 5). The first two PCs explained a total of 60 % to 72 % and 53 % to 58 % of the variance in the winter half year and summer half year, respectively. Hence, most of the variance was maintained when reducing the number of considered PCs from six or seven to two, especially in the winter half year.

The loadings of the first two PCs were similar among the 2-month periods of the winter half year but different and more variable in the summer half year (Fig. 6). Although the loadings differed between the 2-month periods, PC1 was always a meaningful indicator for air mass classification, which will be explained in the following. PC2, as well as higher-order PCs, did not always allow for an unambiguous interpretation.

According to Eq. (6), variables with a high absolute loading determine the PC scores, i.e., the linear combinations, to a large extent. If, for example, an original variable is much higher than its 2-month mean value and its loading on PC1 is strongly negative then the variable will strongly contribute to a large negative score of PC1. In the winter half year, large negative scores of PC1 primarily reflected high CO, CO₂, and CH₄ mixing ratios and rather low values of Δθ_v (i.e., rather low static stability) and the air pressure at GAP (p_GAP), which suggested an uplift of BL air masses, while large positive scores of PC1 reflected the opposite characteristics suggesting a subsidence of UFT/SIN air masses (Fig. 6). PC2 primarily separated air masses with low q, high O₃, and rather high CH₄ mixing ratios from air masses with the opposite characteristics in the winter half year. On its own, PC2 was a less reliable indicator for vertical transport because low q and high O₃ mixing ratios suggested a subsidence of UFT/SIN air masses while high CH₄ mixing ratios suggested an uplift of BL air masses. Possibly, PC2 reflects not only vertical but also horizontal concentration gradients.

Figure 6 Loadings of the input variables on the first (PC1) and second principal component (PC2) for each 2-month period. O₃ was only used as an input variable in the winter half year (black arrows).

In the summer half year, the loading of CO₂ on PC1 had the opposite sign compared with the loadings of CH₄ and CO on PC1. This observation suggests that low CO₂ mixing ratios were generally indicative of BL air masses during the vegetation period due to CO₂ removal by photosynthesis. Additionally, p_GAP contributed more strongly to the first two PCs than in the winter half year and was anticorrelated with q. Δθ_v exhibited very small absolute loadings on the first two PCs in the April–May and June–July periods, indicating a poor correlation with the other variables.

In April and May, an uplift of BL air masses was indicated by large negative scores of PC1, which primarily reflected high q, low CO₂ mixing ratios, and low p_GAP, and by large negative scores of PC2, which primarily reflected high CO and CH₄ mixing ratios.

In June and July, PC1 separated cases with high CH₄, low CO₂, and low p_GAP, which indicated an uplift of BL air masses, from cases with opposite characteristics, which indicated a subsidence of UFT/SIN air masses. In the same period, PC2 could not be interpreted unambiguously with respect to vertical transport because high CO mixing ratios, which suggested BL air masses, but also low q and high p_GAP, which suggested UFT/SIN air masses, contributed to large negative scores of PC2.

In August and September, PC1 separated cases with high CO and CH₄ and low CO₂ and Δθ_v (i.e., low static stability), which was typical for BL air masses, from cases with the opposite characteristics, which was typical for UFT/SIN air masses. Again, PC2 was a less reliable air mass indicator because low p_GAP and high q, which suggested an uplift of BL air masses, but also high CO₂ mixing ratios, which were typical for UFT/SIN air masses, increased the scores of PC2 in August and September (Fig. 6).

3.2 Air mass regimes in February and March

The interpretation of the air mass regimes, which were confined by the 33.3 % and 66.6 % percentiles of the first two PCs, differed among the 2-month periods except for the December–January and February–March periods. As an example, the period February–March is discussed in detail using summary statistics. A case study illustrating the classification results in the measured time series is shown in Sect. S3 in the Supplement. The joint data availability of the PCA input variables was 87 % of the time in February–March.

For the regimes I, II, and III, the summary statistics of all seven PCA input variables were consistent with a recent uplift of BL air masses (Fig. 7a–g): compared with the other regimes, CO, CH₄, and CO₂ generally exhibited high mixing ratios and q was also relatively high while O₃ and p_GAP were low or intermediate. Δθ_v indicated a weakly stable stratification with a median of approximately 8 K for regimes I to III.

Figure 7 Summary statistics of PCA input variables (a–g) and validation variables (h–r) for the nine statistical air mass regimes in February and March 2014. The boxes show the quartiles and the box width is proportional to the number of data points; each whisker is limited to a length of 1.5 times the interquartile range. Total precipitation (precip_ZPLT) is depicted in a bar plot. Unless labeled differently, the data were measured at UFS except for Δθ_v, which represents the maximum difference among the sites. Symbols are explained in Sect. 2.2 and 2.3. The data availability is specified as fraction of time and is shown above the plots. For NO_x and NO_y, local pollution events were excluded. For ⁷Be_ZSG, N₉₀, eBC, and PM₁₀, cases with precipitation (precip_ZPLT>0) were excluded. The ordinate axis was limited to m±3 MAD, where m and MAD are the median and median absolute deviation of individual regimes, respectively.

The other chemical and meteorological measurements supported the interpretation of regimes I to III (Fig. 7h–r). Relatively high RH (medians between 80 % and 90 %) and total precipitation values were consistent with a current uplift of air masses. For NO_y and NO_x, local pollution events (${\mathit{\sigma }}_{{\mathrm{NO}}_x}>\mathrm{0.4}$ ppb) were excluded, while for ⁷Be, N₉₀, eBC, and PM₁₀, cases with precipitation at ZPLT were excluded from the summary statistics because local emissions or wet deposition obscure the chemical signature resulting from vertical transport processes. Nevertheless, the data for these variables were available 75 % to 78 % of the time. Regimes I to III exhibited high NO_y, NO_x, and ²²²Rn and low ⁷Be concentrations compared with the other regimes, which emphasizes the strong influence of the BL. The high NO_x mixing ratios suggest a strong influence of combustion processes. High ²²²Rn concentrations are typical for the continental BL because of natural emissions from the ice-free land surface (Griffiths et al., 2014). ⁷Be concentrations are generally low in the BL due to the formation by cosmic rays in the stratosphere and upper troposphere (Stohl et al., 2000). The aerosol concentrations, N₉₀, eBC, and PM₁₀ were somewhat elevated for regimes I to III but not as strongly as the gas-phase measurements. This finding could be due to wet deposition during the uplift of the air masses. At ZSG, southeasterly to southerly winds with varying velocities were frequently observed for the three regimes under consideration (see Fig. S2a, d, g in the Supplement). The observations suggest that the influence of BL air masses was mostly caused by low-pressure systems or south foehn events in February and March.

Regimes IV, V, and VI exhibited intermediate CO, CH₄, and CO₂ mixing ratios and an intermediate p_GAP (Fig. 7a–g). For these air masses, vertical transport processes were not inferable but the intermediate mixing ratios suggest an influence from both the BL and the free troposphere. Thus, the three regimes under consideration were attributed to the HYBRID air mass class, which was in line with intermediate NO_y, NO_x, and ²²²Rn concentrations (Fig. 7l–n).

Regime VII was generally characterized by the lowest CO, CH₄, and CO₂ mixing ratios and a strongly stable stratification (median Δθ_v of 15 K). However, low O₃ mixing ratios, high q, and an intermediate and strongly variable p_GAP indicated that the air masses did not originate from the upper troposphere or stratosphere but from the lower troposphere or the marine boundary layer (Fig. 7a–g). This interpretation was in line with intermediate ⁷Be and ²²²Rn concentrations and low NO_y, NO_x, N₉₀, eBC, and PM₁₀ concentrations (Fig. 7l–r). Intermediate RH with a median of 54 % corroborated some recent influence of the free troposphere. The air temperature was high compared with the other regimes (Fig. 7i), the wind direction at ZSG almost always had a southerly component, and the wind velocity at ZSG showed the highest mean wind speed of 10.22 m s⁻¹ among the regimes (see Fig. S2c in the Supplement). The summary statistics did not allow for a final conclusion on the influence of the marine boundary layer to be drawn; thus, regime VII was attributed to the HYBRID air mass class and the MBL/UFT subclass unless the criteria for LRMD were fulfilled.

Balzani Lööv et al. (2008) clustered back trajectories for Jungfraujoch, another mountain site in the Alps, and demonstrated that some air masses originated from the marine boundary layer above the tropical Atlantic Ocean and were transported to Jungfraujoch inside the free troposphere within 5 to 15 days, resulting in a rather high RH (average of 67 %), high temperatures, and low O₃, CO, NO_y, and NO_x mixing ratios. Air masses with similar trajectories may also reach Zugspitze and represent, at least partly, regime VII.

For regimes VIII and IX, the CO, CH₄, and CO₂ mixing ratios were generally low but not as low as for regime VII. O₃ was intermediate for regime VIII and high for regime IX. In combination with a strongly positive Δθ_v (median of approximately 15 K), i.e., a strongly stable stratification, and a generally high p_GAP, these observations indicated a subsidence of UFT/SIN air masses for the two regimes under consideration unless the criteria for LRMD were fulfilled (Fig. 7a–g). This interpretation was supported by a generally low RH (medians of 45 % and 25 %), a total precipitation of almost zero, low NO_y, NO_x, N₉₀, eBC, and PM₁₀ concentrations, and high ⁷Be concentrations compared with the other regimes (Fig. 7h–r).

For all regimes except III, the box plots of global radiation (R_g) were similar with medians close to zero (Fig. 7k), indicating that the air mass characteristics were generally independent of day- and nighttime in February and March.

Figure 8 Relative frequency of the air mass classes as fraction of the classified cases in each 2-month period. The joint data availability is given as fraction of time above the plot. LRMD (long-range transport of mineral dust) represents certain cases within the regimes that were not attributed to BL (recent contact with the boundary layer). LRMD and MBL/UFT (influence of the marine boundary layer or UFT) are part of HYBRID (influences of both the boundary layer and the free troposphere or ambiguous). UFT/SIN denotes undisturbed free troposphere or stratospheric intrusion.

3.3 Seasonal frequencies of the air mass classes

In the entire year, the data availability limited the percentage of classified cases to 78 % of the time. Data gaps predominantly occurred in the June–July and October–November periods (Fig. 8). On average, BL, UFT/SIN, and HYBRID air masses accounted for 31 %, 14 %, and 55 % of the classified cases, respectively. The percentage of UFT/SIN air masses is in reasonable agreement with the study of Yuan et al. (2019), in which 13.6 % of the long-term CO₂ data from UFS were selected as baseline concentrations using an univariate statistical approach called “adaptive diurnal minimum variation selection”. For earlier CO₂ data from ZSG and a pedestrian tunnel, approximately 70 m above the UFS, the percentage of baseline concentrations was 19.5 % and 9.9 %, respectively, according to the same study.

The BL air mass class was attributable to three of the nine air mass regimes in all 2-month periods except June–July and October–November, when the BL class comprised two regimes (Fig. 8). This finding implies that air masses are lifted up to UFS with similar frequencies in all seasons although the underlying processes may vary. However, the frequency of the air mass classes was constrained to some extent by using the 33.3 % and 66.6 % percentiles of the first two PCs as thresholds between the air mass regimes. These percentiles were somewhat arbitrary but allowed for a distinction between more unambiguous air masses at both edges of the distribution and more ambiguous air masses in the center of the distribution.

Figure 9 Comparison between the statistical and mechanistic classifications of air masses. The mechanistic classes were based on the boundary layer height at GAP (MLH_GAP) and thermally induced anabatic and katabatic winds (Sect. 2.5) and accounted for 3 %, 5 %, and 17 % of the statistically classified cases (from left to right). The thick black lines show the percentages considered as an agreement between the approaches. MLH_GAP was discarded if clouds were detected below 4 $\mathrm{km}\mathrm{a}.\mathrm{s}.\mathrm{l}.$. LRMD and MBL/UFT are subclasses of HYBRID.

UFT/SIN air masses were attributable to two regimes in the winter half year but only to one regime in the August–September period and to none of the regimes in the April–May and June–July periods. UFT/SIN air masses were only evident if the summary statistics of an air mass regime indicated a dominating influence of subsidence. From April to July, subsidence appeared to be very rare at UFS, such that UFT/SIN air masses could not be isolated by splitting the data at the 33.3 % and 66.6 % percentiles of the first two PCs. This result can be explained by the fact that subsidence in high-pressure systems is counteracted by thermally induced uplift during daytime, especially in the summer half year. However, the data availability in June–July was only 56 %, which calls the long-term representativeness for this period into question. Other statistical classifications at the Jungfraujoch and Mount Bachelor Observatory high-alpine sites indicated that air masses of the undisturbed free troposphere occur year-round but less frequently in summer than in winter (Herrmann et al., 2015; Ambrose et al., 2011). In principle, these air masses do not have to descend but can also be advected horizontally to high-alpine sites. Such cases may be included in the HYBRID class and most likely in the MBL/UFT subclass.

HYBRID air masses including MBL/UFT comprised four to five regimes from August to March and even six and seven regimes in the April–May and June–July periods, respectively. In most of the 2-month periods, MBL/UFT air masses accounted for one regime (Fig. 8). Seven LRMD events were identified during the whole year, predominantly from April to July, which is in line with the study of Flentje et al. (2015) that found 5 to 15 Saharan dust events per year at the Hohenpeissenberg mountain site, approximately 40 km north of Zugspitze.

3.4 Comparing the statistical and mechanistic classifications

The ceilometer-based MLH at GAP exhibited a pronounced diurnal cycle with a maximum in the late afternoon or early evening on days with a high total global radiation (not shown). Little diurnal variation in the MLH was observed for days with a low total global radiation. Thermally induced anabatic winds were most frequently observed during summer and daytime, whereas thermally induced katabatic winds occurred most frequently at night and in the first few hours after sunrise (not shown). The patterns described above were expected and confirmed the plausibility of the mechanistic approach.

A total of 25 % of the statistically classified cases were evaluated using the mechanistic approach. In 3 % of the statistically classified cases, the mechanistic approach suggested an influence of the local BL on the UFS, mostly due to thermally induced anabatic winds. It should be noted that the mechanistic approach is not able to detect synoptic uplift processes because low-level clouds and precipitation result in data gaps of the MLH_GAP. In 5 % of the statistically classified cases, katabatic winds and a low MLH_GAP compared with the UFS level suggested non-BL air masses. In another 17 % of the statistically classified cases, the winds were not thermally induced and the MLH_GAP suggested non-BL air masses.

When the mechanistic approach suggested BL air masses, the statistical approach agreed in 56 % of the cases and yielded HYBRID air masses in almost all other cases (Fig. 9). This finding indicates poor performance of the statistical approach in the case of thermally induced uplift. However, thermally induced uplift only accounts for a small fraction of the time according to the mechanistic approach. Other uplift processes including non-local uplift were not identifiable with the applied mechanistic approach. LRMD was found in 3 % of the cases, which were mechanistically classified as being influenced by the BL, likely due to an overestimation of the ceilometer-based MLH_GAP if the advected dust layer merged with the mixing layer.

In the presence of katabatic winds and a low MLH_GAP compared with the UFS level, the UFT/SIN and HYBRID statistical classes including LRMD and MBL/UFT were both considered to agree with the mechanistic approach (illustrated by the thick black line in Fig. 9). Thus, the level of agreement was 91 % for katabatic winds and a relatively low MLH_GAP, whereby UFT/SIN and HYBRID air masses accounted for approximately half of these cases, respectively. When the MLH_GAP was lower than the UFS level and the winds were not thermally induced, the level of agreement between the classifications was a bit lower (86 %) with UFT/SIN and HYBRID air masses accounting for 34 % and 52 %, respectively.

Uncertainties in the mechanistic approach mainly arise from the assumption of identical MLHs at the UFS and GAP sites and the lack of information on non-local uplift followed by horizontal advection of air masses. Uncertainties in the statistical approach primarily result from the lack of objective thresholds between the air mass classes and a varying significance of sources and sinks of the gases used. In view of these uncertainties, a reasonable agreement between the approaches was found when the mechanistic approach indicated non-BL air masses.

3.5 Implications for real-time operational mode

To classify air masses in real-time, the scores of the first two PCs can be approximated using the loadings that were calculated in this study. Thus, the same thresholds described here can be applied to determine the regime and the class of the current air mass. To standardize the original variables, the arithmetic mean and standard deviation of the current 2-month period are required but these statistics are only known at the end of the 2-month period. Due to the long-term trends of CO₂, CH₄, and CO and the interannual variability of all variables, it is advisable not to use the 2-month means and standard deviations from the year investigated in the present study and to estimate these statistics from recent multiyear records. To make use of the known part of the current 2-month period, the mean (μ) and standard deviation (σ) could be updated regularly using weighted averages,

$$\begin{array}{}\text{(7)}& {\displaystyle }& {\displaystyle }\mathit{\mu }=(\mathrm{1}-f){\mathit{\mu }}_{\text{l}}+f{\mathit{\mu }}_{\text{c}}\text{(8)}& {\displaystyle }& {\displaystyle }\mathit{\sigma }=(\mathrm{1}-f){\mathit{\sigma }}_{\text{l}}+f{\mathit{\sigma }}_{\text{c}},\end{array}$$

where μ_l and σ_l are the 2-month mean and standard deviation, respectively, estimated from the long-term record; μ_c and σ_c are the mean and standard deviation of the known part of the current 2-month period; and f is the time fraction that is known in the current 2-month period. In real-time operational mode, additional uncertainties can arise from uncorrected measurement artifacts.

4 Conclusions

In this study, a novel statistical approach was developed to distinguish between the BL (recent contact with the boundary layer), UFT/SIN (undisturbed free troposphere or stratospheric intrusion), and HYBRID (influences of both the boundary layer and the free troposphere or ambiguous characteristics) air mass classes at Schneefernerhaus at Zugspitze. A main purpose of the classification scheme is later use in real-time operational mode. The scheme was based on the first two principal components, which were calculated from five gas-phase (CO, CH₄, CO₂, O₃, q) and two meteorological (Δθ_v, p_GAP) variables but the principal components and their interpretation differed between the six 2-month periods. Additionally, the PM₁₀ concentration and its 30 min standard deviation were needed to identify long-range transport of mineral dust. In retrospect, local pollution events on the mountain were evident from a high 30 min standard deviation of NO_x approximately 9 % of the time. But these events are too short to be predictable on a 30 min basis in real-time applications; thus, they are neglected by the classification scheme.

While the first principal component was a suitable air mass indicator throughout the year, the second principal component left room for ambiguities. BL air masses were detected in all seasons with similar frequencies (average of 31 % of the classified cases). UFT/SIN air masses were predominantly found in the winter half year (23 % of the classified cases from October to March) but subsidence was so rare from April to July that the UFT/SIN class was not determinable in these months. HYBRID air masses (average of 55 % of the classified cases) mostly exhibited intermediate characteristics leaving room for ambiguities. For 17 % of HYBRID air masses, it remained unclear whether the air mass originated from the lower undisturbed free troposphere or the marine boundary layer. To achieve a distinction between these two areas of origin, future work with trajectory models is necessary. A total of 5 % of the HYBRID air masses was explained by long-range transport of mineral dust.

Independent chemical and standard meteorological measurements such as NO_y, ²²²Rn, and precipitation were in line with the statistical approach. A mechanistic classification based on ceilometer and standard meteorological measurements was feasible in 25 % of the statistically classified cases and predominantly suggested non-BL air masses in these cases, which was in good agreement with the statistical approach. In the rare cases of thermally induced uplift, however, the statistical approach often misclassified BL air masses as HYBRID.

In principle, the statistical classification scheme can be used in real-time operational mode if the unknown arithmetic means and standard deviations of the observational variables in the current 2-month period can be estimated with sufficient accuracy. Future work should test the real-time applicability and quantify the uncertainty arising from the use of approximated means and standard deviations for standardization.

The framework of the presented statistical classification might also be useful at other high-alpine sites because it is based on common measurements.