Visualising Daily PM10 Pollution in an Open-Cut Mining Valley of New South Wales, Australia—Part I: Identification of Spatial and Temporal Variation Patterns

Abstract

The Upper Hunter Valley is a major coal mining area containing approximately 40% of the currently identified total coal reserves in New South Wales (NSW), Australia. Due to the ongoing increase in mining activities, PM10 (airborne particles with an aerodynamic diameter of less than 10 micrometres) pollution has become a major air quality concern in local communities. This paper summarises the spatial and temporal variability modes of PM10 pollution in the region, based on long-term multi-site monitoring data and the application of the rotated principal component analysis (RPCA) and wavelet analysis techniques. RPCA identified two distinct air quality clusters/subregions in the valley: one in the west/northwest and the other in the southeast. Wavelet analysis revealed the annual cycle to be the most persistent temporal mode of PM10 variability in both subregions, with intermittent signals also observed at time scales of around 120, 30~90, and under 30 days. How these variation modes are related to the effects of local PM10 emissions and the influence of meteorology at different time scales deserves further attention in future work. The findings will be used in air quality reporting and forecasting in NSW. The methodology and results can also be useful for air quality research in similar regions elsewhere.

1. Introduction

The Upper Hunter Valley is located within the northern end of the Hunter region, around 200 km north of Sydney and 50 km northwest of Newcastle in the State of New South Wales (NSW), Australia (Figure 1) [1]. On a broad scale, the valley is oriented northwest–southeast (NW–SE), approximately 30 km wide and with a terrain elevation of around 300~380 m from the lower (Singleton South) end to the upper (Merriwa) end [2,3]. It is a major coal mining area in NSW, with significant coal mines located between Aberdeen in the north to Bulga in the south, containing approximately 40% of the currently identified total coal reserves in the state [4]. It also has a significant agriculture industry (dairy, beef, horse breeding, and viticulture) and two large coal-fired electricity-generation plants [1,5], with the Liddell Power Station fully decommissioned from April 2022 to April 2023. There are multiple medium-to-small towns in the region. The largest population centres are Singleton (population~24,577) and Muswellbrook (population~16,357), with many small settlements and isolated rural residences scattered throughout the valley.

The prevailing surface winds tend to follow the NW–SE orientation of the valley. The most frequent winds are north-westerlies in winter and south-easterlies in summer, with wind directions less defined in autumn and spring [3,6,7,8,9]. Nocturnal and early-morning down-valley drainage flows, daytime up-valley winds, or south-easterly sea breezes are also observed from time to time [3,8,9,10,11]. Wind strengths vary across the valley, with most locations experiencing annual average wind speeds in the range of around 2~4 m/s. The precipitation in the region is generally low (compared to coastal regions to the east), and varies significantly across years—for example, the annual rainfall ranged from around 549 mm to 853 mm at Singleton and Muswellbrook during 2011–2015 [12]. Higher rainfall tends to occur in summer and early autumn, and lower rainfall occurs in winter and early spring.

Particle pollution is known to have adverse effects on human health and the environment [13,14,15,16,17,18,19,20]. PM10 (airborne particles with an aerodynamic diameter of less than 10 micrometres) pollution is a major air quality concern for local communities in the Upper Hunter Valley [2,21,22,23,24,25]. The main sources of PM10 emissions in the region include local open-cut coal mining activities (e.g., wheel-generated dust, windblown dust from overburden) and surface soil erosions (windblown dust). For example, coal mining contributes around 88% of PM10 emissions in the combined Muswellbrook, Singleton and Upper Hunter local government areas [23]. Coal-fired electricity generation, agriculture, bushfires, prescribed hazard-reduction burnings, and state-wide dust storms also contribute to PM10 pollution in the region [6].

There have been a few studies examining PM10 pollution in the Upper Hunter Valley. Of these, most early investigations were based on data collected (with different types of instruments) at a limited (small) number of locations through short-term campaign monitoring projects. These include those reported in SPCC (1982, 1983) [26,27], Holmes and Associates (1996) [28], Morrison and Nelson (2011) [25], Hibberd et al. (2013) [10], and relevant references therein. For example, SPCC (1982) [26] reported a few observational and modelling studies, which concluded that localised dust pollution from open-cut coal mines and related developments continued to be an issue of concern and that there were unlikely to be serious cumulative, region-wide problems resulting from dust emissions from mines. Holmes and associates (1996) [28] found that the increase in both deposition and concentration levels of dust over 1984–1994 was due to the increase in coal production and the severe drought affecting much of eastern Australia and that the land affected by cumulative effects appeared to be primarily that owned by the coal mines. Holmes (2008) [8] and Hyde et al. (1981) [9] also suggested that local PM10 emissions in the valley can impact air quality in areas away from sources—that is, north-westerlies can transport dust generated in the upper end of the valley to areas near the bottom of the valley or further down over the metropolitan areas of Newcastle.

In partnership with the Upper Hunter coal and power industries, the NSW Government commissioned the Upper Hunter Air Quality Monitoring Network (UHAQMN) during 2010–2012 (Figure 1), to provide the community, industry, and government with reliable and up-to-date information on air quality within the valley [12]. Pollution (including PM10) data from the network are reported as air quality categories (AQCs) in near real-time and in quarterly or annual data reports on the NSW government website (https://www.airquality.nsw.gov.au; accessed on 10 December 2023). Multi-year PM10 data summaries were available in two main reports, i.e., OEH (2017) [12] for data in 2011–2015 and DPE (2022) [6] for data in 2011–2021. Three main findings are worth highlighting: (1) the annual PM10 concentrations in the UHAQMN were observed to be amongst the highest across the NSW Air Quality Monitoring Network (NSW AQMN, with over 90 stations across NSW) [29]; (2) PM10 levels at some sites can exceed the national benchmark (i.e., Australian standard for PM10) [30] from time to time, and vary significantly from year to year due to impacts of draught conditions and occurrence of hazard reduction burnings or bushfires; and (3) emissions from open-cut mines can lead to elevated PM10 pollution typically in the lower (southeast) end of the valley, in particular under the north-westerly winds.

Jiang (2017) [3] performed a detailed analysis of the 2012–2015 (4-year) PM10 data from the UHAQMN in an air quality management campaign project. The results showed that: (1) poor air quality days generally occurred in spring and secondarily summer and autumn, with winter and February/March having relatively good air quality; (2) there were more poor air quality days in 2012–2013 but significantly reduced number of events in 2014–2015; and (3) the correlations between high PM10 pollution and individual meteorological variables were complex and non-linear, varying with time and location. The author also proposed the existence of two air quality clusters, i.e., the south–east (SE) and west–northwest (WNW) air quality subregions in the valley. Drawing up the project findings, the author indicated that the application of more sophisticated (holistic) analytic methods such as pattern recognition techniques may provide an increased understanding of PM10 pollution in the study region.

Globally, there are relatively few studies on air quality issues in rural valley environments, with much of the air quality literature primarily focusing on urban air pollution problems. Of those few studies, most work was undertaken with data on a small number of sites and based on correlation analysis, e.g., on spatial and/or temporal variations of PM10 by Mohd Shafie et al. (2022) and NPS (2023) [31,32], and correlations between PM10 levels and selected atmospheric parameters by Mannis (1988), Giri et al. (2008), Fortelli et al. (2016), Reisen et al. (2017), Czernecki et al. (2017) and Quimbayo-Duarte (2021) [33,34,35,36,37,38]. In summary, elevated PM10 pollution in those valleys is associated with (prolonged) dry conditions (low rainfall and humidity), low winds, thermal inversions (low mixing heights), and/or the influence of high-pressure systems. To date, to the best of our knowledge, there is little or no research in the literature on the topic of examining the spatial and temporal variability modes of PM10 pollution.

This project expands on Jiang (2017) [3], based on the 11-year (2012–2022) PM10 data from 14 stations in the UHAQMN, by applying advanced analytic methods to holistically examine (1) the spatial-temporal variation patterns of PM10 pollution in the Upper Hunter Valley and (2) how elevated pollution is related to local- and synoptic-scale meteorological configurations in the region. This paper is focused on presenting the investigation results of (1), i.e., the spatial and temporal variation modes from the long-term multi-site PM10 data. The investigation is unique in at least three aspects: (a) the air quality subregions initially proposed in Jiang (2017) [3] are verified and analysed based on a longer dataset using rotated principal component analysis (RPCA); (b) the temporal modes of PM10 pollution were identified for the air quality subregions through wavelet analysis and were illustrated via heat map visualisations; and (c) the impact of exceptional events including dust storms and bushfires on PM10 pollution were also examined in some depth for the subregions. The results on the linkage between PM10 pollution and local and synoptic meteorological features will be reported separately in a companion paper of this text. The findings will used for air quality forecasting in NSW, and the methods and results can also be useful for air quality research in similar regions elsewhere.

2. Data

2.1. Air Quality Data

The 14 monitoring stations in the UHAQMN (Figure 1) can be grouped into station types that serve different purposes [39], as described in

Table 1. Air quality station details and PM10 data used in this study. Source: adapted from OEH (2013) [39].

Station Type	Station Purpose	Station Name *	Total Number of Days with Valid Daily PM10 Data in 2012–2022	Total Number of Days with Invalid or Missing Daily PM10 Data in 2012–2022
Larger population centre	Monitoring air quality in the larger population centres	Aberdeen	3984	34
		Muswellbrook	3968	50
		Singleton	3978	40
Smaller population centre	Monitoring air quality in the smaller communities	Bulga	3967	51
		Camberwell	3970	48
		Jerrys Plains	3950	68
		Maison Dieu	3953	65
		Warkworth	3940	78
		Wybong	3962	56
Diagnostic	Providing diagnostic data that helps to diagnose the likely sources and movement of particles across the region as a whole; they do not provide information about air quality at population centres	Mt Thorley	3876	142
		Muswellbrook NW	3986	32
		Singleton NW	3978	40
Background	Providing background data at the top (Merriwa) and bottom (Singleton South) ends of the valley	Merriwa	3862	156
		Singleton S	3951	67

* NW: Northwest; S: South; there is a total of 130 exceptional event days (Section 2.2).

. Air quality data for these stations were obtained from the NSW Air Quality Data System (AQDS). These included (1) daily (24 h) average PM10 concentrations (μg/m³) for each monitoring station and (2) information on whether and when PM10 measurements at any of these stations were significantly impacted by exceptional events (further details in Section 2.2), such as hazard-reduction burnings (HRBs), bushfires, or widespread dust storms. The data were checked and confirmed to be of high quality, with missing values generally less than 2% at individual stations, except that there was up to 3.9% of missing data at the Merriwa and Mt Thorley stations.

2.2. Definition of Exceptional Event Day, Normal Day, and Poor Air Quality Day

Of the up to 3986 days of valid daily measurements in 2012–2022, there were a total of 130 days when air quality in the Upper Hunter Valley was significantly impacted by air emissions from bushfires, planned hazard-reduction burning (HRB), and/or continental-scale dust storms, with daily PM10 levels above the 24 hr average national benchmark level of 50 µg/m³ [30] at one or more stations in the UHAQMN. These 130 days are referred to as “exceptional event days” in this text. All other days, which were not significantly impacted by bushfires, HRBs, or continental-scale dust storms, are referred to as “normal days”, or interchangeably, “non-exceptional event days”. Most of the exceptional event days occurred during the 2019–2020 spring–summer bushfires and the widespread dust storm events in 2018, both associated with impacts of prior prolonged droughts across Australia [6,40,41,42].

These definitions led to two derived daily PM10 datasets for 14 monitoring stations, i.e., the normal-day dataset (excluding data for exceptional event days listwise across all stations) and the exceptional-event-day dataset. The normal-day dataset reflects the air quality conditions mainly associated with PM10 emissions within the valley, primarily from open-cut mining activities and secondarily from soil erosion, power generation, and agriculture activities. The exceptional-event-day dataset is related to the significant air quality impacts due to emissions from (local or remote) exceptional events such as vegetation fires or dust storms. For easy reference, the original PM10 dataset is referred to as an all-day or full dataset (i.e., including data from both above datasets), depicting the full air quality conditions associated with effects from both local (i.e., within-valley) and broad-scale PM10 emissions in the study region.

To facilitate our analysis, a day was marked as a “poor air quality day”, also referred to as “exceedance day” or “poor air pollution day”, if the PM10 concentration on the day was above the national benchmark level. This definition led to counts of poor air quality days for each monitoring station. The terms “PM10 pollution” and “air pollution” are used in an exchangeable manner in this text.

3. Methods

The investigation was conducted in four steps to examine the spatial and temporal variability modes of daily PM10 pollution. The first step was to describe the general properties of PM10 pollution in the Upper Hunter Valley, based on boxplots and general summary statistics for the long-term (2012–2022) PM10 measurements from 14 stations in the UHAQMN (Section 2.1). This sets out the context for the subsequent analyses and interpretation of results from other steps. The second step was to identify the spatial co-variation structure (i.e., spatial regionalisation) of PM10 pollution by applying RPCA (rotated principal component analysis) to the multi-site PM10 data in 2012–2022 (Section 3.1). This analysis verifies the air quality clustering (subregions) previously identified by Jiang (2017) [3] from a dataset with a shorter time period (2012–2015). In the third step, wavelet analysis was applied to the principal component (PC) time series obtained from RPCA to determine the dominant temporal modes of PM10 pollution in each subregion, as well as how those modes changed over time (Section 3.2). Then, in the fourth step, heat maps were used to illustrate the variation patterns of average PM10 levels and the number of exceedance days in each subregion on the annual and interannual time scales.

Throughout steps 1 to 4, the effects of exceptional events (vegetation fires and dust storms) on the findings were examined by repeating the relevant analyses separately for the all-day and normal-day datasets, and where appropriate, the dataset for exceptional event days only (see definitions in Section 2.2). Due to the space limit, this text is focused on presenting findings from the normal-day dataset, highlighting the spatial and temporal variability of PM10 pollution primarily associated with dust emissions from open-cut mining activities and soil erosion (within the valley), which are of great concern in local communities. The findings from the analyses on the all-day dataset are briefly summarised in a comparative manner, with details given as Supplementary Materials to support our interpretations of results. The implementation of RPCA and wavelet analysis are described next.

3.1. Rotated Principal Component Analysis (RPCA)

Principal component analysis (PCA) can be used for data reduction, variation mode identification, or feature classification [43]. In mode identification or feature classification, a rotation technique is often applied to PCs for the easy interpretation of results [44]. Here, PCA was applied to the correlation matrix of the daily PM10 concentration time series (2012–2022) for 14 air quality stations to examine the structure of inter-site covariations, i.e., the spatial regionalisation of PM10 pollution in the valley. Missing data points were not used, i.e., ignored listwise across stations, when calculating the correlation matrix to suppress the potential for introducing new noise due to a data imputation process (note: PCA was applied differently in Section 3.2 for the purpose of wavelet analysis). The retention of the number of principal components (PCs) was decided through a scree-test, following Jiang (2011) [45] by the combined use of the Cattell (1966) [46] and North et al. (1982) [47] methods. A plot of initial eigenvalue and sampling error by PC number was used for this purpose, keeping in mind that degenerated multiplets (i.e., unrotated PCs) should not be separated from one another and as much variance in the dataset as possible should be explained by the PCs retained. Varimax (orthogonal) rotation was then applied to the retained PCs to facilitate physically meaningful interpretations while preserving the linear independence between PCs [48].

In addition, an obliquely rotated PCA (without the orthogonality constraint for PC rotation) [45] was performed on the same data to confirm if the linearly independent PCs from the Varimax rotation could approximate the true (underlying) “simple structure” in the data [49]. The RPCA was repeated for the daily PM10 data with or without exceptional event days included, i.e., all-day dataset vs. normal-day dataset, to test the stability of the RPCA findings.

3.2. Wavelet Analysis

Wavelet analysis can simultaneously generate a representation of a signal (stationary or non-stationary) in the time and frequency domains, thereby allowing access to localised information about the signal [50]. Wavelet analysis was performed to identify the dominant temporal variability modes of PM10 pollution if existing in the study region. Wavelet analysis requires that no missing data exist in the time series. Hence, an RPCA was performed on a processed daily PM10 dataset, where the missing records (if any) at each station were filled with the median of the PM10 measurements in the 2012–2022 period for that station, resulting in PC time series (scores) without missing values but otherwise (almost) identical to those from the RPCA described in Section 3.1. A wavelet analysis was then applied to these PC scores to determine the dominant modes of co-variability in PM10 pollution for each subregion, and how those modes changed over time.

We followed the procedure outlined in Torrence and Compo (1998) [51], using the Morlet wavelet with the following parameter configurations: sampling rate/resolution = 1 (day); frequency resolution = 0.25, i.e., four suboctaves (voices) per octave; lower and upper Fourier periods (scales) for wavelet decomposition set to 2 days and 2048 days (~5.6 years), respectively. The lower and upper Fourier periods were chosen to allow for variability modes on the sub-weekly, annual (seasonality) and multi-year (5–6 years) scales, as well as anything in between being considered in this analysis. We assumed a red-noise process with the lag-1 serial correlation for each PC time series when testing the significance of wavelet spectrum power at the 0.05 level for a Chi-square test.

We repeated the above analysis process on the PCs derived from the RPCA for the full dataset and the normal-day dataset separately (note: for RPCA on the normal-day dataset, both the missing data and the readings for exceptional event days were replaced with the median of PM10 measurements for normal days by individual stations). This comparison was to assess the (potential) impact of the inclusion of high particle measurements in exceptional event days on the wavelet analysis results. The wavelet analysis on the PCs for the full dataset should identify all signals (temporal variability modes) if existing in the time series. In contrast, the wavelet analysis on the PCs for the normal-day dataset may identify the existence of regular signals (if any) on non-exceptional days.

4. Results and discussions

4.1. General Description of Daily PM10 Pollution and Impacts of Exceptional Events

Box plots and summary statistics for the normal-day dataset (i.e., exceptional event days excluded) are given in Figure 2 and Figure 3 to illustrate the distributional properties of the long-term (2012–2022) PM10 concentration data from 14 stations in the UHAMQN. In general, relatively higher daily PM10 levels were recorded at four (direct) source-impacted locations, i.e., the Camberwell, Maison Dieu, Mt Thorley and Singleton NW stations (Figure 2), which are relatively close to open-cut mining sites in the southeast (lower end) of the valley (Figure 1). These stations had larger variabilities (indicated by larger interquartile ranges, i.e., box lengths) and more outlier/extreme values (indicated by dots and asterisks). PM10 levels at the Warkworth station were comparably high. In contrast, Merriwa (background station near the top of the valley), and secondarily, Wybong, Jerrys Plains, Bulga and Aberdeen, recorded relatively lower and less variable PM10 pollution, as indicated by the lower box positioning and relatively fewer outlier/extreme values. PM10 levels at the two largest population centres, Singleton and Muswellbrook, appeared in between the above two cases, with extreme values also observed at these locations. Minimum daily PM10 levels were generally similar across 14 stations, ranging from 1.5 to 3.5 µg/m³ (Figure 3). Overall, these results are generally consistent with the distributional properties reported by Jiang (2017) [3] for a shorter (normal days in 2012–2015) PM10 dataset. Holmes and Associates (1996) [28] noted that the land affected by cumulative effects appeared to be primarily owned by the coal mines. Essentially, as is the case for other regions [52,53], the significant inter-site variability in PM10 pollution reflects the combined effects of changes in local emissions (primarily from open-cut mining and soil erosion) and environmental factors such as meteorological conditions on different time scales.

An analysis was also performed for the all-day (full) dataset (including exceptional event days), as illustrated in Figure 3 and Figure S1 in Supplementary Materials. The distributional properties of PM10 data were found to be very similar to those described above for the normal-day dataset. The exception is that the full dataset showed slightly higher means, medians, and standard deviations, significantly (around two or more times) higher maximum values (Figure 3), and many more days with outlier or extreme values (Figure S1). These differences between datasets were expected, essentially reflecting the significant impact of exceptional events on regional air quality, such as the widespread dust storms in 2018 and the spring–summer bushfires across large areas of NSW in 2019–2020 [41,42,54], as is further illustrated on the interannual scale in Section 4.3 and Section 4.4.

We also compared the summary statistics between the normal-day and exceptional-event-day datasets (Figure 4). The mean and median of PM10 concentrations for exceptional event days were above the national benchmark (50 µg/m³) at all stations, equivalent to an around 2~3 times increase over those for normal days. The minimum PM10 pollution for exceptional event days reached around half of the national benchmark level, or an around 5~10 times increase over that for normal days. The maximum and standard deviation values for exceptional event days were also around 1~4 and 2~7 times higher, respectively, compared to normal days. Overall, the four source-impacted stations were most significantly affected, with Merriwa recording the highest maximum PM10 level due to the impact of a widespread dust storm on 11 January 2020 [55]. This finding indicates the accumulated (combined) impacts on air quality from local (within-valley; primarily associated with open-cut mining and soil erosion), remote (mainly continental dust storms) and incidental (HRBs or bushfires) particle emissions.

4.2. Spatial Pattern—Identification of Two Air Quality Subregions

The RPCA for the normal-day dataset (Section 3.1) led to the retention of two leading principal components (PC1 and PC2). The two Varimax rotated PCs explain around 88% of the total variance in the dataset: 45% by PC1 and 43% by PC2. Figure 5 shows the loadings of two PCs, equivalent to Pearson correlations between each PC time series and daily PM10 concentrations at individual stations. The PC loadings (correlations) identify two distinct air quality clusters: (1) the WNW subregion, with high loadings on PC1 for stations in the northern and western parts of the valley, and (2) the SE subregion, with high loadings on PC2 for stations in the south-eastern part of the valley. Hence, the variability of PM10 pollution in the valley is summarised into two linearly independent (dimensionless) time series, i.e., PC1 and PC2 scores, despite the significant inter-site difference in PM10 concentrations (as shown in the previous section).

RPCA on the all-day dataset (including exceptional event days) could produce similar results, except that the two leading PCs are reversed in order and two clusters show slightly reduced separations (Figure S2 in Supplementary Materials). Obliquely rotated PCA [45] on the relevant datasets generated very similar (almost identical) results, confirming the existence of a strong underlying simple structure (clusters) [49] in the daily PM10 dataset, i.e., the existence of two distinct air quality subregions in the valley.

These results are highly consistent with Jiang (2017) [3], who initially proposed the two-subregion property of PM10 pollution from a shorter (4-year) normal-day dataset for the study region. Hence, based on the longer (11-year) dataset, the present analysis has verified the stability and robustness of two air quality subregions in the Upper Hunter Valley, despite the impact of (a small number of) exceptional events.

The division of the valley into two air quality subregions can somehow be expected, primarily due to the valley’s NW–SE-oriented slope terrain and the prevalence of valley-following airflows in the lower boundary layer. Previous work revealed that the most frequent winds are north-westerlies and south-easterlies in the valley [6,7,8,12]. North-westerly winds can blow dust generated in the upper (NW) part of the valley south-eastward (downslope), contributing to elevated PM10 concentrations in areas near the bottom end (i.e., the SE subregion) of the valley [9]. In contrast, south-easterly winds may transport pollutants from the lower valley north-westward (upslope), resulting in elevated PM10 pollution over the upper end (i.e., the WNW subregion) of the valley. Some recent observational case studies [6,12] also indicated that emissions from open-cut mines can lead to elevated PM10 pollution, typically in the lower (southeast) part of the valley, in particular under north-westerly winds.

4.3. Temporal Pattern—Identification of Key Variation Modes

The previous section shows that the variability of PM10 pollution in the study region can be expressed with two linearly independent PC time series—that is, PC1 and PC2 scores represent the distinct covariational features of PM10 pollution at stations in the WNW and SE subregions. Wavelet analysis was used to identify the dominant variability modes in the PC time series and how those modes changed over time (Section 3.2). To reveal the existence of temporal variability patterns on normal days (i.e., PM10 levels were not significantly impacted by exceptional events), Figure 6 and Figure 7 illustrate the wavelet analysis results for PC1 (Figure 6a; WNW subregion) and PC2 (Figure 7a; SE subregion) scores from the RPCA of the normal-day dataset. It is most distinguishable that the wavelet power spectrum peaks near the 360-day time scale (annual mode, significant at the 0.05 level for a Chi-square test) persistently across all years for both subregions (Figure 6b and Figure 7b). The higher variance occurred in spring to summer for the WNW subregion (Figure 6c) but in winter to spring for the SE subregion (Figure 7c).

The annual variability mode in PM10 pollution can be linked to the seasonal changes in weather and climatic conditions that influence pollutant (PM10) emissions and dispersion conditions in the valley. Higher rainfall (hence, lower dust generation and higher wet deposition) tends to occur in summer and early autumn, and lower rainfall (hence higher dust generation and lower wet deposition) occurs in winter and early spring [12]. The most frequent winds in the valley are north-westerlies in winter and south-easterlies in summer, with wind directions less defined in autumn and spring (Section 1) [3,8]. Consequently, the prevailing north-westerly winds and higher PM10 emissions in winter and spring provide a high potential for transporting or accumulating PM10 pollution over the SE subregion. In contrast, the prevailing south-easterly flows or sea breezes in summer tend to transport pollutants north-westward, potentially resulting in elevated PM10 pollution in the upper part of the valley (i.e., in the WNW subregion). This seasonality in PM10 pollution is further illustrated in Section 4.4 for two air quality subregions.

There are also intermittent wavelet power peaks at around 120 days, 30~90 days, and shorter (less than 30 days) time scales (Figure 6b and Figure 7b), indicating the signals of triannual, intraseasonal, and shorter-term variation modes in the PC time series. These variation modes appeared stronger in the SE subregion (PC2 time series) than in the WNW subregion (PC1 time series). The signal strengths changed dramatically across the years, manifesting phase differences between the two subregions. For example, the signal for the variability mode of around 120 days in the WNW subregion was more intense in summertime across years 2013–2015 and 2017–2018 but a lot weaker during some periods in 2021–2022 (Figure 7b). In comparison, this variability mode in the SE subregion showed relatively higher power during wintertime in 2012, 2013 and 2018 but lower power during some periods in 2016 and 2019–2022. Of note is the high wavelet power at around 5–6 years for both subregions (Figure 6b and Figure 7b), implying a signal at the interannual scale in the PC time series. However, most of the spectra at this time scale fall within the cone of influence (COI) region (under the curved line), which is within the uncertainty range associated with the edge effect in Fourier transforms [51].

The variability modes at times scales of 30–90 days and 120 days in PM10 pollution are not yet readily understood. Many studies have examined the signal of intraseasonal oscillations at time scales of 30–90 or 30–60 days in tropical atmospheric or oceanic variables, which can often be linked to Madden–Julian Oscillation (MJO) in the tropics and extra-tropical teleconnection patterns such as El Niño-Southern Oscillation (ENSO) [56,57,58]. Such oscillations are also investigated extensively in atmospheric variables associated with monsoon activities [59,60] and to a lesser extent in those at higher latitudes [61,62]. For example, Rimbu et al. (2013) [62] revealed two intraseasonal variability patterns in synoptic observations (temperature, wind speed, sea level pressure) at a high-latitude Antarctic research station, Neumayer (70° S, 8° W). The dominant pattern manifests with out-of-phase variations of temperature and wind speed with sea level pressure at time scales of around 40 and 80 days, which can be related to tropical intraseasonal oscillations via large-scale eastward-propagating atmospheric circulation wave-trains. In contrast, the second pattern was characterised by the in-phase variations between temperature, wind, and sea level pressure at time scales of around 35, 60, and 120 days, which can be attributed to the Southern Annular Mode/Antarctic Oscillation (SAM/AAO). Drawing upon these studies, one may speculate whether the PM10 variability modes at time scales of 30–90 and 120 days (and potentially at the multi-year time scale) could be related to the influence of large-scale climate drivers such as MJO, ENSO and SAM, which are known to modulate the weather and climate in Australia. This aspect deserves further attention in future work.

In addition, the shorter-term (under 30 days) variability in PM10 pollution is attributable to day-to-day changes in emission conditions and the effects of local and synoptic weather variability, as is demonstrated in some previous studies for other regions, e.g., Jiang et al. (2014) [52] for Auckland and Jiang et al. (2017) [53] for Sydney. This aspect will be discussed in a companion paper of this text.

A wavelet analysis was also performed on the two leading PCs from the RPCA on the all-day PM10 dataset (including measurements for exceptional event days; Section 3.2), with results shown in Figures S3 and S4 in Supplementary Materials (note that the order of two PCs was reversed). Overall, the wavelet power of the PC time series could identify generally similar modes to those described above. However, clearly, there are significant distortions (increases) in the wavelet power peaks during mid-2017 to mid-2020 for the SE subregion (PC1 times series; Figure S3) and during mid-2019 to mid-2020 for the WNW subregion (PC2 time series; Figure S4). The distortions appear more intense in the SE subregion than the WNW subregion, essentially reflecting the significant impact of high PM10 measurements during the widespread dust storms in 2018 and the unprecedented Black Summer bushfires in 2019–2020 across large areas of NSW. In addition, the wavelet powers (signals) at time scales of around 120 days and 30–90 days appear less distinguishable if compared to those on the PC time series derived from the normal-day dataset, especially for the WNW subregion. Hence, the inclusion or exclusion of exceptional-event-day measurements in the wavelet analysis can lead to appreciable differences in the wavelet transform results, somehow indicating a limitation of the wavelet analysis technique for analysing time series containing extreme values.

It is worth noting that we also performed a set of wavelet analyses directly on the daily PM10 concentration time series at individual stations (in comparison with the PC scores), separately for the all-day and normal-day datasets (missing data gaps at each station were filled with the median value of the records for that station). To illustrate, Figures S5 to S8 (in Supplementary Materials) show the wavelet transform results for PM10 time series at the Singleton (in the SE subregion) and Muswellbrook (in the WNW subregion) stations, for measurements on all days and normal days only, respectively. Overall, the wavelet power distributions are very similar between the PM10 time series and the PC scores for the relevant subregions, except that the wavelet transform of PCs tends to further highlight the persistence of the (common) annual variability modes—for example, Figure S6 on PM10 data from the Singleton site compared to Figure 7 on the PC time series for the SE subregion, and Figure S8 on PM10 data from the Muswellbrook site compared to Figure 6 on the PC time series for the WNW subregion. Hence, the station-specific wavelet analyses confirm the subregional temporal variability patterns identified in the PC times series from RPCAs, both for the normal-day and all-day datasets.

4.4. Illustrating the Annual and Interannual Variability in PM10 Pollution

4.4.1. Mean PC Scores by Year and Month for Each Subregion

The knowledge of the spatial and temporal variability patterns (Section 4.2 and Section 4.3) can be used to facilitate the summarisation of PM10 data with increased clarity. To demonstrate this utility, we use PC scores to visualise the general annual (seasonality) and interannual variation patterns of PM10 pollution for two subregions in the same (one) framework. Figure 8 shows the heat maps of mean PC1 and PC2 scores for normal days (exceptional event days excluded) by month and year for each subregion. The analysis of PC scores for the all-day (full) dataset makes no significant difference to the findings reported here.

It is visually clear that mean PC scores identify very distinct annual variation patterns (seasonality) in PM10 pollution between two subregions. In the WNW subregion (Figure 8a), positive PC1 mean scores (hence, high mean PM10 levels) occurred mostly in warmer months (late spring to summer; highest in November to February), but negative scores (low mean PM10 levels) in cooler months (May to September; lowest in June). In the SE subregion, however, the seasonal variability pattern appears complex (Figure 8b). Positive mean PC2 scores (thus higher mean PM10 levels) tended to occur in some cooler months (in particular, May and July-October) and negative scores (thus lower PM10 levels) in some warmer months (November-April) in most years. As expected, the mean scores for November and December 2019 are also high, reflecting the broad-scale higher mean PM10 pollution associated with the impacts of the unprecedented 2019–2020 spring–summer bushfires in NSW [42].

The different patterns in mean PC scores for the two subregions may be broadly associated with the seasonal variations in meteorological conditions, as is discussed in Section 4.2. The most frequent winds in the valley are north-westerlies in winter and south-easterlies in summer, with wind directions less defined in spring and autumn [6,12]. Higher rainfall (thus lower dust generation and higher wet deposition) tends to occur in summer and early autumn, and lower rainfall (higher dust generation and lower web deposition) in winter and early spring. Hence, during winter and spring, the prevailing north-westerly winds can blow dust generated in the upper (NW) part of the valley south-eastward (downslope), contributing to elevated PM10 concentrations in areas near the bottom end (i.e., the SE subregion) of the valley. In contrast, during summer, the prevailing south-easterly winds may transport pollutants from the lower valley north-westward (upslope), resulting in elevated PM10 pollution at the upper end (i.e., the WNW subregion) of the valley.

The interannual changes in PC scores appear similar between two subregions (Figure 8). Higher values occurred in 2018–2019 and secondarily 2012–2013 (this variability signal is slightly weaker in the WNW subregion), but lower values in 2021–2022 and secondarily 2015–2016, with 2022 being the cleanest year on record [63]. Of note is that June has generally negative mean scores across all years in both subregions, indicating generally better air quality for this time of the year.

4.4.2. Mean PM10 Levels and Total Number of Poor Air Quality Days

In this section, the annual and interannual variability patterns in mean PM10 concentrations and total number of poor air quality days are illustrated at the station level for two air quality subregions (Figure 9 and Figure 10). The results confirm the above identified variability patterns in PC scores, but with greater details on PM10 pollution at individual locations/stations.

In the WNW subregion (Figure 9a), the mean PM10 levels and number of poor air quality days were generally higher in warmer months (i.e., late spring to summer, highest in November to January) but lower in cooler months (lowest in June). In the SE subregion (Figure 9b), however, two statistics were generally higher in August-November (late winter to spring) and May (highest in September and October), but lower in February, March and June (lowest in June), with near average values in December to January. As for the subregional mean PC scores (Figure 8), the significant difference in mean PM10 concentrations and the number of poor air quality days for two subregions (Figure 9) may also be associated with seasonal variability in meteorological conditions. In winter and spring, the lower rainfall condition provides a potential for higher dust generation but lower wet deposition, and the prevailing north-westerly winds can blow dust generated in the upper (NW) part of the valley south-eastward (downslope), contributing to elevated PM10 pollution at stations near the bottom end (the SE subregion) of the valley. In contrast, during summer the prevailing south-easterly winds may transport air pollutants from the lower valley north-westward (upslope), leading to elevated PM10 pollution at stations near the upper/western end (i.e., the WNW subregion) of the valley. This aspect will be further investigated in a companion paper of this text.

Broadly, the interannual variation patterns are similar between the two subregions, with higher values in 2018–2019 and 2012–2013 (this variability signal is slightly weaker in the WNW subregion) but lower values in 2021–2022 and 2015–2016 (Figure 10). Consistently, Jiang (2017) [3] also noted higher PM10 pollution in 2012 and 2013 than in other years during 2012–2015. OEH (2014) [40] and DPE (2019, 2020) [41,42] suggested that the (prior) higher-than-average temperature and prolonged drought conditions across NSW (broadly, most of the Australian continent) had contributed to the poorer air quality observed in years 2012–2013 and 2018–2019. DPE (2022, 2023) [6,63] noted that the cooler and wetter climate conditions contributed to improved air quality in 2021–2022.

Overall, the SE subregion has generally higher mean PM10 pollution and more exceedance days (highest at the four source-impacted sites) than the WNW subregion (as expected), as is consistent with local air quality experience [6,28]. Of the three larger population centre sites, higher PM10 pollution occurred at Muswellbrook and Singleton. In contrast, Aberdeen experienced relatively better air quality, with the mean PM10 levels and total number of exceedance days comparable to Merriwa (background site). Also of note is that June appears to be the cleanest across all months, recording the lowest mean pollution and the least number of poor air quality days across almost all stations.

The distribution patterns identified for the normal-day dataset were further confirmed through a similar analysis on the all-day dataset (Figures S9 and S10 in Supplementary Materials). In comparison, as expected, the mean PM10 levels and poor air quality days for the all-day dataset show a significant increase for November-January and 2018–2020, primarily associated with the broad-scale impact of dust storms in 2018 and the bushfires in the 2019–2020 spring–summer months [42].

5. Summary and Conclusions

The present study has examined the spatial-temporal variability patterns of PM10 pollution in the Upper Hunter Valley, based on long-term (2012–2022) multi-site air quality data from the UHAQMN and the application of RPCA (rotated principal component analysis) and wavelet analysis techniques. The impact of exceptional events has also been examined in some depth, with some results given as Supplementary Materials to support the interpretations in this text. The main findings are summarised below:

• The RPCA identified two air quality subregions in the Upper Hunter Valley: one in the west/northwest (WNW) part and the other in the southeast (SE) part of the valley. The inclusion or exclusion of exceptional-event-day measurements in the RPCA analysis does not make a significant difference to the results. This finding verifies the previous work by Jiang (2017) [3], confirming the spatial regionalisation property of PM10 pollution in the valley despite of the significant impacts of extreme events such as widespread dust storms or vegetation fires. This suggests that the identified spatial variability modes in PM10 pollution are primarily associated with changes in local (within-valley) PM10 emissions and the influence of meteorological conditions in the Upper Hunter Valley. Hence, it is possible to characterise the air quality variability in the valley with two linearly independent and dimensionless PC-score time series, or alternatively, PM10 time series from a subset of (representative) monitoring stations from the subregions.
• Wavelet analysis has identified the annual cycle (and potentially interannual variability) to be the most persistent temporal variation mode(s) of PM10 pollution in two subregions, in particular on normal days when air quality is mainly affected by the within-valley emissions sources such as open-cut mining and soil erosion, rather than significantly impacted by exceptional events. There were also intermittent signals at time scales of around 120 days (triannual), 30~90 days (intraseasonal), and under 30 days, with the mode intensities changing dramatically across time. These intermittent signals in PM10 pollution are not yet readily understood and deserve further attention in future research.
• The inclusion or exclusion of exceptional-event-day measurements in the wavelet analysis can lead to appreciable differences in the wavelet transform results. With the (extreme) exceptional event data removed, the wavelet analysis helped to extract the regular signal patterns (prevailing on normal days), with a loss of information at time points of extreme measurements. In contrast, if the exceptional event data were included, the wavelet power spectra were significantly dominated (distorted) by extreme measurement values—in this case, the wavelet analysis somehow showed a reduced ability to identify the (normal-day) temporal variability patterns in the time series. This finding indicates that cautions need to be taken when applying the wavelet transform method in air quality research, where air pollutant measurements can be (unavoidably) very irregular or extreme, primarily due to the impacts of incidental pollutant emissions on the region of interest.
• The knowledge of the temporal and special variability modes can be used to facilitate the summarisation of PM10 data with increased clarity, as has been demonstrated for illustrating the seasonal and interannual variability patterns for two air quality subregions. The seasonal variation patterns differed between subregions: higher pollution occurred in warmer months (summer in particular) in the WNW subregion but in late winter to spring in the SE subregion. The interannual variation patterns were broadly similar for the two subregions, with higher PM10 pollution in 2018–2019 and 2012–2013 but lower PM10 pollution in other years, whereas 2022 observed the lowest PM10 pollution on record.
• Relatively higher daily PM10 levels were recorded in the southeast of the valley and were highest at four direct source-impacted locations, which are relatively close to the open-cut mining sites scattered near the bottom end of the valley. The stations in the WNW subregion had generally lower PM10 pollution, with one background site (Merriwa) and one larger population site (Aberdeen) recording relatively good air quality. The other two larger population centre sites, Singleton and Muswellbrook, experienced generally moderate PM10 pollution levels when compared to other stations.
• The exceptional events, including widespread dust storms in 2018 and bushfires in 2019–2020, had significant impacts on regional air quality, resulting in significant increases in minimum, mean and maximum PM10 levels, as well as the number of poor air quality days, across all monitoring stations in the valley. The impact was most significant in the SE subregion, in particular at the four source-impacted stations, indicating the combined impacts from local, remote, and/or incidental emissions. This suggests a need for the enhanced management of open-cut mining activities, even during exceptional event days.
• The division of two air quality subregions in the valley, as well as the significant differences in mean and elevated pollution levels between the two subregions, may be associated with the interactions between NW–SE-oriented valley-slope terrains and the seasonal variability in meteorological conditions. In winter and spring, the lower rainfall conditions provide a potential for higher dust generation (but lower wet deposition), and the prevailing (lower-boundary-layer) north-westerly winds blow dust generated in the upper part of the valley south-eastward, contributing to elevated PM10 pollution at locations near the bottom end (the SE subregion) of the valley. In contrast, during summer, the prevailing south-easterly winds may transport air pollutants from the lower valley north-westward (upslope), resulting in elevated PM10 pollution at stations near the upper/western end, i.e., the WNW subregion of the valley.

In conclusion, the present study has illustrated a quantitative, holistic approach to examining the spatial-temporal variability patterns of PM10 pollution in the Upper Hunter Valley. The findings from this study will be used to improve air quality reporting and forecasting in NSW. Future work may point to a few directions. For example, one area to investigate is the underlying mechanisms associated with the PM10 variation modes at time scales of 30–90 days, 120 days, and potentially multiple years, in particular when even a longer dataset becomes available. The present analysis suggests that the division of the valley into two air quality subregions can be linked to the valley’s NW–SE-oriented slope terrain and the prevalence of valley-following airflows in the lower boundary layer. Previous studies showed that both local-, synoptic-, and large-scale climatic conditions could affect local air quality in a region [52,53]. Research in this aspect is ongoing, with details to be presented in a companion (Part II) paper of this text. The methodology and results can also be useful for air quality research in similar regions elsewhere.