key: cord-0868686-sm5zxha3
authors: Jeong, Dajeong; Yoo, Changhyun; Yeh, Sang-Wook; Yoon, Jin-Ho; Lee, Daegyun; Lee, Jae-Bum; Choi, Jin-Young
title: Statistical Seasonal Forecasting of Winter and Spring PM(2.5) Concentrations Over the Korean Peninsula
date: 2022-03-28
journal: Asia Pac J Atmos Sci
DOI: 10.1007/s13143-022-00275-4
sha: 8186c6301f0257696a5d895cbb5e355fc673ed66
doc_id: 868686
cord_uid: sm5zxha3

Concentrations of fine particulate matter smaller than 2.5 μm in diameter (PM(2.5)) over the Korean Peninsula experience year-to-year variations due to interannual variation in climate conditions. This study develops a multiple linear regression model based on slowly varying boundary conditions to predict winter and spring PM(2.5) concentrations at 1–3-month lead times. Nation-wide observations of Korea, which began in 2015, is extended back to 2005 using the local Seoul government’s observations, constructing a long-term dataset covering the 2005–2019 period. Using the forward selection stepwise regression approach, we identify sea surface temperature (SST), soil moisture, and 2-m air temperature as predictors for the model, while rejecting sea ice concentration and snow depth due to weak correlations with seasonal PM(2.5) concentrations. For the wintertime (December–January–February, DJF), the model based on SSTs over the equatorial Atlantic and soil moisture over the eastern Europe along with the linear PM(2.5) concentration trend generates a 3-month forecasts that shows a 0.69 correlation with observations. For the springtime (March–April–May, MAM), the accuracy of the model using SSTs over North Pacific and 2-m air temperature over East Asia increases to 0.75. Additionally, we find a linear relationship between the seasonal mean PM(2.5) concentration and an extreme metric, i.e., seasonal number of high PM(2.5) concentration days. SUPPLEMENTARY INFORMATION: The online version contains supplementary material available at 10.1007/s13143-022-00275-4.

High concentrations of fine particulate matter smaller than 2.5 μm in diameter (PM 2.5 ) are a major environmental contributor to cardiovascular and respiratory diseases in East Asian countries (West et al. 2016; Song et al. 2017) . Annual mean concentrations of PM 2.5 over the Korean Peninsula exceed 25 µg m − 3 (Lee 2014; Bae et al. 2020) , more than double of the World Health Organization's (WHO's) air quality guideline value of 10 µg m − 3 annual mean of (Chae et al. 2021 ). On certain extreme days (e.g., March 3-5, 2019), the daily mean PM 2.5 over the Seoul metropolitan area reach approximately 110 µg m − 3 (Lee et al. 2019 ). This record also substantially exceeds the WHO's daily mean guideline of 25 µg m − 3 . To minimize citizens' exposure to high PM 2.5 concentrations, the Korea Ministry of Environment (KME) issues advisories (warnings) when the 2-hour mean concentration exceeds or is expected to exceed 75 (150) µg m − 3 .

The National Institute of Environmental Research (NIER) of KME conducts operational PM 2.5 concentration forecasts for lead times up to a week. Along with observation and meteorological field forecasts, their primary forecasting tool is the Community Multiscale Air Quality System (CMAQ) developed by the U.S. Environmental Protection Agency (Byun and Schere 2006) . A recent report indicated their daily forecasts have an accuracy as high as 85% when Korean Meteorological Society measured based on four categories: good (0-15 µg m − 3 ), moderate (16-35 µg m − 3 ), unhealthy (36-75 µg m − 3 ), and very unhealthy (76 µg m − 3 and above) (NIER 2019). However, demand to extend their forecast out to seasonal time scale, i.e., 1-3 months ahead, has been growing. Information regarding the air quality of upcoming seasons can play an important role in planning emission controls and bolstering vulnerable populations' resilience against harmfully high concentration events.

Many operational centers provide seasonal climate predictions, such as temperature and precipitation forecasts. Such predictions are understood as boundary condition problems because a few months in the future, errors in initial atmospheric conditions become too substantial to provide useful information. Instead, these predictions rely on the sensitivity of the atmosphere to its slowly varying boundary conditions, such as sea surface temperature and soil wetness. In general, operational centers utilize two types of forecast guidance for their seasonal climate predictions: dynamical model forecast guidance and statistical model forecast guidance. Dynamical models provide numerical solutions for physical climate processes, which has become increasingly accurate as a result of advances in ensemble forecasting techniques and initializations of land and ocean. Statistical models, on the other hand, are based on empirical relationships between predictors and predictands. El Niño Southern Oscillation (ENSO) is, for example, one of the most important predictors of the globe's seasonal climate; as a result, many institutions carefully monitor and predict it (e.g., https:// iri. colum bia. edu/ our-exper tise/ clima te/ forec asts/ enso/ curre nt/). One clear advantage of seasonal statistical models is that they have comparable accuracy to dynamical models and are much less expensive to develop and maintain.

Previous studies have found that various forms of interannual climate variability, such as ENSO and East Asian winter monsoons, strongly influence PM 2.5 concentrations over East Asia (Shimadera et al. 2013; Jia et al. 2015; Yang et al. 2016; Jeong and Park 2017; Lang et al. 2017; Yu et al. 2019; Kumar et al. 2021) . This suggests that the use of a statistical model to predict seasonal PM 2.5 over the Korean Peninsula could be successful. For example, circulation changes associated with winter monsoons can result in anomalous PM 2.5 transport that modulates surface-level concentrations (Jeong and Park, 2017) . Applying this information regarding wintertime monsoons and Siberian Highs in a simple linear regression model generates a skillful PM 2.5 concentration prediction for the same season with a correlation coefficient of approximately 0.72 (Jeong et al. 2021 ). However, due to the insufficiency of long-term surface PM 2.5 observations, previous studies have relied heavily on either proxy data of PM 2.5 , such as aerosol optical depth from satellite observations (e.g., Jia et al. 2015) , or model-generated data, such as atmospheric chemical transport model simulations guided by atmospheric reanalysis data (Jeong et al. 2021) . Also, while some studies have sought to generate short-term PM 2.5 predictions (e.g., Lee et al. 2013; Xiao et al. 2020; Chae et al. 2021) , models have not been developed yet for long lead times of 1-3 months and for different seasons, meaning our understanding of PM 2.5 seasonal prediction remain limited. Lastly, extreme metrics, such as number of high concentration days, have only recently been investigated for fine particulate matter smaller than 10 μm in diameter, PM 10 (e.g., Cho et al. 2021; Ku et al. 2021) , and connections between seasonal mean PM 2.5 concentrations and high concentration days are not examined.

Recognizing these gaps in existing research, this study seeks to make statistical seasonal predictions of PM 2.5 concentrations over the Korean Peninsula based on atmospheric boundary conditions, which is the first attempt to the best of our knowledge. We use long-term in situ observation of PM 2.5 to construct our model and verify its predictions. To understand the characteristics of the data, we explore seasonal distributions and connections between seasonal mean PM 2.5 concentrations and high concentration days. This is followed by examination on boundary conditions by land and ocean as potential sources of prediction skill. Our analyses verify predictions for lead times from 1 to 3 months during winter and spring, respectively. The remainder of the study is organized as follows: Section 2 describes the data and model formulation; Section 3 presents the results of correlation analyses and model verifications; and Section 4 explains and discusses the conclusions of the study.

We use hourly PM 2.5 concentration data observed at 25 stations in Seoul, South Korea, from 2005 to 2018, provided by Seoul Research Institute of Public Health and Environment (SRIPHE; https:// data. seoul. go. kr/ dataL ist/ OA-15526/S/ 1/ datas etView. do). This is the longest period of PM 2.5 concentration observation in South Korea. We find that across the stations approximately from 40.1 to 62.8% of the hourly observations are missing. Unrealistic values are masked when they are either less than 0 µg m − 3 or greater than 300 µg m − 3 , removing overall less than 0.1% of the entire data. We obtain a daily time series by taking averages for each day and over all the stations (black line in Fig. 1) .

We also employ another hourly PM 2.5 concentration dataset by KME (NIER 2019). Collection of this data began in 2015, and the dataset therefore covers a much shorter time period than the one collected by SRIPHE. However, the KME dataset has extensive spatial coverage over South Korea, is quality controlled, and is publicly available after a few months delay (https:// www. airko rea. or. kr). Measurements from 9 stations (Seoul, Busan, Daegu, Daejeon, Gwangju, Gangneung, Taean, Chungju, and Andong) are selected to cover temporally extensive length of the 2015-2019 period and spatially the Korean Peninsula. At each station, between 46.2 and 23.4% of data are missing, except Taean station that misses 76.9% of its data. As for the PM 2.5 of SRIPHE, we construct a daily time series without missing values by averaging the hourly concentration from KME over days and stations (red line in Fig. 1 ).

To benefit from both datasets, we construct a long-term observation by appending the KME dataset (2015-2019) to the SRIPHE dataset (2005) (2006) (2007) (2008) (2009) (2010) (2011) (2012) (2013) (2014) , dropping SRIPHE data of the overlapping 2015-2018 period. Note that between 2015 and 2018, the two datasets exhibit an extremely high correlation: correlation coefficient (R) = 0.88, which exceeds the 99% statistical significance level. Seasonal averages for DJF and MAM will further ensure that the differences between the two data are smoothed out. We find that the overall conclusions of this study remain the same when we use only the SRIPHE data for 2005-2018 (not shown). Therefore, in this study, we present results for 14 seasons based on the merged data during DJF of 2005 /2006 -2018 and MAM of 2006 We define the number of high PM 2.5 concentration days by counting the total number of days for a season that daily mean concentrations exceeds 36 µg m − 3 . In terms of air quality category of Korea, this includes both unhealthy (36-75 µg m − 3 ) or very unhealthy (76 µg m − 3 and above) air quality categories (NIER 2019). Also, the definition of high PM 2.5 concentration days resembles the definitions on high PM 10 concentration days in previous studies (Yun and Yoo 2019; Ku et al. 2021 ).

As potential predictors of seasonal PM 2.5 concentration, we examine five variables that provide boundary conditions of atmospheric circulation: sea surface temperature (SST), soil moisture, and 2-m air temperature (T 2m ), along with sea ice concentration and snow depth. For SST, we employ monthly data obtained from Met Office Hadley Centre (Rayner et al. 2003) . The data has a horizontal resolution of 1° and is fill valued when the value of the SST data is lower than − 1.8 °C, indicating formation of sea ice. For soil moisture, T 2m , and snow depth, we use the monthly European Centre for Medium-Range Weather Forecasts (ECMWF) fifth generation of atmospheric reanalysis (ERA5; Hersbach et al. 2020) dataset. We integrate four volumetric soil water layers in the reanalysis, which ranges approximately 0-3 m, to represent soil moisture. The horizontal resolution of the datasets was originally at 0.25° but is interpolated to 1.5° using bilinear interpolation. For sea ice concentration, we take the data derived using measurements from the Scanning Multichannel Microwave Radiometer (SMMR) and from the Special Sensor Microwave/ Imager (SSM/I) sensors through bootstrap algorithm (Comiso 2017) . The sea ice data is at 25 km resolution on the polar stereographic grid.

For all the climatological variables, we first apply the seasonal averages. We then compute anomalies by removing the seasonal cycles, which are defined as the calendar month means. Lastly, we subtract the linear trends defined for each calendar month from the anomalies to contain variability of interannual time scale or longer without trends. 

Multiple linear regression model is built for predictions of seasonal mean PM 2.5 concentrations and numbers of high PM 2.5 concentration days during DJF and MAM at lead times of 1-3 months. To identify the predictors, we apply the forward selection stepwise regression approach, which is described in Section 3 (see also Myoung et al. 2020) .

Note that the regression model is built after the linear trends are removed respectively from predictors and predictands for each season. Also note that here we limit the maximum number of predictors to be equal to or less than three, including the linear trend of the PM 2.5 concentration. This helps us avoid overfitting the model, given the short time period of observation. As a result, all models turn out to have three predictors, which make the equation of the model:

where ŷ is the predictands, such as DJF mean PM 2.5 concentration or number of high PM 2.5 concentration days during DJF; and x 1 and x 2 are two predictors based on boundary conditions (Section 3.2) with their lagged regression coefficients, b 1 and b 2 . The model includes the linear trend of the PM 2.5 concentration b t , which is defined for each season, and time x t , as well as the intercept term, b 0 . To validate the prediction skill of the regression model, we apply a leave-one-year-out cross-validation approach. This means we construct the regression model for a certain year based on the data from all other years and verify the model using the data for the chosen year. We repeat this process each year for each model. Three categories can be defined for categorical forecast using lower and upper climatological terciles of the predictand distribution. Because only 14 seasonal values are available, we obtain the terciles from 30-day running averaged data for each season, instead of seasonal values, an approach we discuss further in the following section.

We first explore the probabilistic distributions of PM 2.5 concentrations by seasons to determine the characteristics of our chosen predictand variable. A histogram of daily averaged PM 2.5 concentrations reveals overall positive skewness of the data (Fig. 2 ). Many days of low and moderate PM 2.5 concentrations are clustered around the left tail of the distribution while some extremely high PM 2.5 concentrations can be seen

in the right tail. Compared to the annual case (Fig. 2) , the skewness is alleviated for the DJF and MAM distributions (Fig. 2) , which is due to reductions in relative frequency of occurrence of low PM 2.5 concentrations.

In contrast, we find seasonally averaged PM 2.5 concentrations are expected to be much more symmetrical than the daily means. Because we only have data for 14 years, daily PM 2.5 concentrations are 90-day running averaged, instead of being seasonally averaged to examine possible Gaussianity in the seasonal mean field, assuming that all daily concentrations are independent from one another. Given that the decorrelation time scale of the PM 2.5 concentration over the Korean Peninsula is only a few days, the assumption seems reasonable. The distributions of the 90-day running average concentrations are bell shaped for all seasons, despite some greater-than-normal kurtosis for the DJF distributions (Fig. 2) .

Similarly, 90-day running days can be used to estimate underlying distributions of high concentration days instead of seasonal number of days (Fig. 2) . Positive skewness is apparent in the all-season case (Fig. 2) . However, for DJF and MAM, the skewness is alleviated, and the distributions are roughly normal distributions with no pronounced kurtosis (Fig. 2) . MAM show a larger number of high concentration days than DJF.

Next, we examine seasonally averaged PM 2.5 concentrations and seasonal high PM 2.5 days. Note that we use these metrics as predictands in our statistical forecasting models. Interestingly, a linear relationship can be seen between the quantities, and more strongly for MAM (Fig. 3) . The Pearson correlation between the metrics for MAM is 0.95, suggesting that the same statistical model may be used to predict both predictands, particularly for MAM (Fig. 3) , because the number of daily extreme events primarily determines the seasonal mean value. For DJF, the correlation between the metrics is slightly less, at 0.71. Despite some linearity, we observe a range of high PM 2.5 days near the seasonal mean PM 2.5 value of approximately 30 µg m − 3 (Fig. 3) , suggesting that each predictand may have a different predictor. In terms of location of the predictands, the means of the seasonal values for DJF are 29.91 µg m − 3 and 8.17 days, respectively. These seasonal values are slightly higher than those from MAM (28.02 µg m − 3 and 7.36 days). For the spread, however, MAM exhibits a greater range of values in both mean concentration and number of days (see also Fig. 2 ).

We begin selecting the models' first predictors using SST data. As described earlier, we build respective models for DJF and MAM, first targeting seasonal mean PM 2.5 concentration predictions out to 3 months. Next, we select Korean Meteorological Society predictors by exploring soil moisture, T 2m , sea ice concentration, and snow depth through stepwise regression approach.

In search of potential predictors in SST, we examine lagged correlation maps between seasonal mean concentrations and SSTs north of 30°S (Fig. 4) . Note that in the maps, the three-month averaged SSTs are lagged for the PM 2.5 concentrations of DJF and MAM by 1-3 months. For DJF PM 2.5 concentrations, significant values can be observed over the North Atlantic, most pronouncedly over the equator (Fig. 4) , with a hint of tripole-like SST mode in the North Atlantic including the area south of Greenland. The North Atlantic is known to influence East Asia through changes Fig. 2 Probability distributions of (a-c) daily mean and (d-f) 30-day running mean PM 2.5 concentrations, along with 90-day running total number of high PM 2.5 days. We obtain the distributions for (left panels) all seasons, (middle panels) DJFs, and (right panels) MAMs. The overlaid lines are normal distributions obtained using the means and standard deviations of the underlying data in atmospheric circulation that exhibit poleward and then equatorward wave trains (e.g., Ok et al. 2017; Myoung et al. 2020 ). In addition, due to the slowly varying nature of the ocean, the correlation patterns remain largely unchanged from September to November (i.e., 3 months ahead) through October-December (2 months ahead) and November to January (1 month ahead), suggesting that the equatorial Atlantic SST exerts a significant influence at seasonal time scale.

Area averaged at the three regions with the largest correlation coefficients, we select the DJF SSTs as individual predictors in our regression models (Table 1) . The model built with the equatorial Atlantic SST (10°S-10°N, 30°W-10°E; green boxed region in Fig. 4) shows a 0.52 correlation with detrended DJF PM 2.5 concentration. The same calculation using the western North Pacific (30°N-65°N, 110°E-150°E) and North Atlantic (50°N-70°N, 20°W-30°E) SSTs show smaller correlations (i.e., 0.36 and 0.24, respectively). Moreover, the strongly negative correlation we find between the equatorial Atlantic and western North Pacific SSTs suggests the predictors may be Table 1 Pearson correlation between the detrended seasonal mean PM 2.5 concentration and the simple linear regression model built with a selected SST predictor

We select three SST predictors for (left) DJF and (right) MAM, respectively, based on Fig. 4 . Note that models are built no lag between the predictor and predictand but leave-one-year-out cross-validated. Models with the highest correlation are shown from the top, and the correlations that exceed the 95% confidence level of the two-tailed Student t-test are asterisked. Also shown are the correlations between a predictor and the most skillful predictor Fig. 4 ) as a first predictor in our models. As with the DJF models, we make this decision based on correlation maps (Fig. 4) and simple linear regression model results (Table 1 ). Note that through ENSO, the North Pacific SST is likely to be connected with the equatorial eastern Pacific SST (e.g., mode 1 of global SST in Fig. 2 of Messié and Chavez 2011) , which also exhibits significant correlations. Our analysis shows significantly negative correlations in the North Pacific throughout the lags for MAM (Fig. 4) . When built with this North Pacific predictor, the model demonstrates an astonishing accuracy, with R = 0.77 (Table 1) . SSTs over the equatorial Atlantic (20°N-30°N, 80°W-40°W) and equatorial eastern Pacific (20°S-10°N, 140°W-50°W) also prove high skills (R = 0.64 and 0.47, respectively) in MAM concentration predictions, but we do not select them due to high cross correlations with the North Pacific SST (-0.57 and − 0.71, respectively).

After establishing the first predictors based on SSTs, we examine lagged relationships between seasonal mean PM 2.5 concentration and soil moisture (Fig. 5) . Based on the statistical significance of lagged correlation values for each season, two potential regions are selected. During DJF, soil moistures over the eastern Europe (45°N-55°N, 40°E-50°E; green boxed region in Fig. 5 ) and the Gobi Desert (40°N-50°N, 90°E-100°E) show large correlation coefficients (Fig. 5) . We examine these two candidate Fig. 4 , except that the lagged correlations with seasonal mean soil moisture are shown Note that we linearly remove the SST predictor from the detrended PM 2.5 concentration before constructing and verifying the models for the chosen soil moisture predictor. Also shown are the correlations between the selected SST predictor (see main text and regions using the stepwise procedure; after linearly removing the equatorial Atlantic SST predictor from the PM 2.5 data, regression models are built and verified using the soil moistures ( Table 2 ). The resulting analysis highlights the benefit of using the eastern Europe soil moisture as a second predictor. The correlation is R = 0.82 for the DJF predictions of the detrended and equatorial Atlantic SST removed PM 2.5 concentration. The eastern Europe soil moisture exhibits a low correlation with the equatorial Atlantic SST, suggesting weak dependency between the predictors. Repeating the process for MAM, we find that soil moisture over the eastern Europe (45°-60°N, 55°E-70°E) and southern China (20°N-30°N, 90°E-110°E ) are good potential sources of predictability (Fig. 5) . The correlations obtained from the models using these regions are reasonably high, R = 0.61 and 0.57, respectively (Table 2) . However, in contrast to the DJF case, both soil moisture predictors are highly correlated with the North Pacific SST (R = -0.53 and − 0.38, respectively), making them less ideal as predictors.

Including the linear trend of each season, we have so far identified three and two predictors for DJF and MAM PM 2.5 concentration prediction models, respectively. To complement the model, we now examine T 2m over the land as potential predictors. For MAM, high correlations can be seen over almost all of the Eurasian continent (Fig. 6) , which we divide into three subdomains: East Asia (25°N-45°N, 110°E-120°E; green boxed region in Fig. 6) , Siberia (45°N-65°N, 60°E-110°E) , and western Europe (35°-50°N, 10°W-20°E). Based on its extremely high accuracy and negligible correlation with the SST predictor, the T 2m over East Fig. 6 As in Fig. 4 , except that the lagged correlations with seasonal mean T 2m are shown Asia is selected as the final predictor for the MAM model (Table 3) . Despite their reasonable accuracy and moderate correlations with the SST predictors, we leave aside the other two regions to limit the total number of predictors in our models. For DJF, a high correlation emerges over the eastern Europe (45°N-55°N, 40°E-50°E), but this region is already used for soil moisture. The predictor selection process is reiterated for sea ice concentration and snow depth. However, the latter two do not make to our final selection because our analyses identify no large areas with significant correlations (Figs. S1-S2 ). This may be surprising given that some previous studies have shown these two factors to be good predictors of springtime rainfall over the Korean Peninsula (Myoung et al. 2020) , which is an important wet removal mechanism of PM 2.5 concentrations. Moreover, there are studies that find evidence of linkage between the first mode wintertime PM 10 concentrations over the Korean Peninsula and Arctic sea ice (Kim et al. 2019; Lee et al. 2020) . Although our results do not rule out the impact of sea ice and snow depth on PM 2.5 concentrations, we speculate that physical and chemical processes other than precipitation may play an important role in PM 2.5 concentrations. Also, the influence of sea ice and snow depth may be pronounced only during certain months of winter or spring, instead of the entire season (Kim et al. 2019) .

To examine whether the circulation patterns modulated by the chosen predictors resemble the circulation anomalies associated with seasonal PM 2.5 concentrations, we first inspect the interannual relationship identified using detrended PM 2.5 concentrations with detrended seasonal mean geopotential height at 500-hPa (Fig. 7) . This is compared with the 2-month lagged correlation maps between October-November-December (OND) predictors and Fig. 7 Correlation coefficients using seasonal-mean detrended 500-hPa geopotential height anomalies and the PM 2.5 concentrations of the Korean Peninsula for (a) DJF and (b) MAM. We obtain (c, e) lagged correlation coefficients between DJF detrended 500-hPa geopotential height anomalies and October-December (OND) SST anomalies averaged over equatorial Atlantic (10°S-10°N, 30°W-10°E) and soil moisture anomalies averaged over (45°N-55°N, 40°E-50°E). For (d, f), lagged correlation coefficients are computed using MAM detrended 500-hPa geopotential height anomalies and January-March (JFM) SST anomalies averaged over North Pacific (35°N-45°N, 160°E-180° ) and T 2m anomalies averaged over East Asia (25°N-45°N, 110°E-120°E) . Regions with correlation values that exceed the 95% confidence level of the two-tailed Student t-test are hatched DJF height field (Fig. 7) and also between January-February-March (JFM) predictors and MAM height field (Fig. 7) . During DJF, we find zonally elongated areas of positive coefficients over Mongolia and in the east of the Korean Peninsula of R near 0.3 (Fig. 7) . The coefficient values are not statistically significant by the 95% confidence level, which weakens their robustness. Nonetheless, this points to a decrease in the climatological wind speed as well as advection of PM 2.5 or relevant substances by anomalous southerly during high PM 2.5 concentration. The high anomaly over the Korean Peninsula also indicates a weakening of vertical ventilation due to lowered atmospheric boundary layer and increased static stability (e.g., Fig. 4 in Cho et al. 2021) . The SST and soil moisture predictors of DJF likewise exhibit anomalously high coefficients over Mongolia and in the east of the Korean Peninsula, respectively (Fig. 7) . Pattern correlation between the circulation anomalies associated with the PM 2.5 concentrations and equatorial Atlantic SST over 30°S and poleward is R = 0.63. However, the value is weak (R = 0.16), for the eastern Europe soil moisture-understanding the mechanisms of this factor will require further study.

During MAM, anomalously high correlations over Mongolia and China and a tripole structure over the North Pacific are pronounced in the correlation maps between the PM 2.5 concentration and the 500-hPa geopotential height (Fig. 7) . The high anomaly over China forms a northwesterly to the Korean Peninsula, while the tripole indicates a weakening of the climatological jet stream. The MAM correlation map does not resemble that of DJF, suggesting that different mechanisms may work during the spring season, and further studies will be required to figure out the mechanisms for MAM PM 2.5 concentrations. Nonetheless, based on the importance of local T 2m in the MAM predictions, we speculate that chemical processes may play an important role in MAM PM 2.5 concentrations. The lagged correlation maps for the SST and T 2m predictors of MAM align closely with the pattern correlations of 0.76 and 0.73, respectively, over poleward of 30°S (Fig. 7) .

To evaluate the performance of the chosen predictors, we first examine observed and forecast seasonal mean PM 2.5 concentration time series (Fig. 8) . As discussed earlier, equatorial Atlantic SST, eastern Europe soil moisture, and the linear trend are selected for the DJF regression models (Fig. 8a) , and North Pacific SST, East Asia T 2m , and the linear trend for the MAM models (Fig. 8) . We find overall skills of the DJF forecasts were shown as R = 0.87, 0.84, and 0.69 in terms of correlations for lead times of 1-3 months, indicating that the measures of the accuracy gradually decrease with lead times. For the MAM models, the correlations are R = 0.93, 0.76, and 0.75, slightly higher for lag 1 and 3 months than the skills for the DJF models. Examining the time series, we find declining trends in the PM 2.5 concentrations of both seasons, indicating that setting the linear trend as a predictor has a positive effect. Interannual variation is more pronounced during MAM than DJF, and this is well captured by the predictors.

A scatterplot of the observed and forecast concentrations can reveal over-or under-forecasting biases and refinements of forecast resolutions (Fig. 9) . We use the scatterplot instead of a calibration function due to the very small sample size. Overall, the scattered values line up well with the diagonals, suggesting the forecasts have good accuracy, as already demonstrated by the correlations. The fact that the values are not one sided from the diagonals indicates that the forecasts have no pronounced bias. However, the DJF forecasts hint at some conditional bias, as the regression lines (colored lines) exhibit slopes less than 1 (Fig. 9 ). This indicates that the DJF models are underconfident and tend to predict smaller (greater) concentration values when higher (lower) seasonal concentrations are observed. The MAM models show better calibration than the DJF models despite featuring greater spread from the diagonal.

Because most operational climate forecasts provide categorical information, we reevaluate deterministic predictands as categorical forecasts. Three categories, i.e., below normal, normal, and above normal seasonal concentrations, are defined using the upper and lower terciles obtained from the 90-day running averaged daily PM 2.5 concentrations of each season (dashed black horizontal lines in Fig. 8 ). We adopt this approach because the terciles from 14 observed seasonal values alone are not likely to be reliable. We use the Heidke skill score (HSS; Eq. 8.22 in Wilks 2011), a common seasonal climate forecast measure, to quantify the accuracy of the three categorical forecasting. HSSs are 77.42%, 29.50%, and 29.50% for DJF and 78.63%, 100%, and 46.97% for MAM forecasts. The abrupt decline in the DJF skill is caused by the small interannual variation during this season (Figs. 2 and 3) , which makes the normal category very narrow (between the horizontal dashed lines in Fig. 8a) . Nevertheless, given that 100% HSS signifies a perfect forecast accuracy and 0% HSS represents the accuracy that can be expected from a random forecast, these are skillful predictions.

We constructed prediction models for seasonal mean PM 2.5 concentrations and seasonal number of high PM 2.5 concentration days using multiple linear regression models based on slowly varying atmospheric boundary conditions and long-term trends in the predictands, targeting the winter and spring seasons over the Korean Peninsula at 1-3-months ahead. To build the statistical models, we extended the KME's PM 2.5 observations from 2015 to 2019 by appending 2005-2014 observations from the SRIPHE. For this period, we found significant empirical relationships were found between the predictands and sea surface temperature, soil moisture, and 2-m air temperature. Skillful predictions were possible for seasonal mean PM 2.5 concentrations; the correlations between the observations and the forecasts at 3-months ahead reached 0.69 and 0.75 for the December-February (DJF) and March-May (MAM) seasons, respectively. We found slightly lower accuracy for the models of seasonal number of high PM 2.5 concentration days; the correlations were 0.64 and 0.58 for DJF and MAM, respectively.

As predictors, we chose equatorial Atlantic SST, eastern Europe soil moisture, and the linear term trend for DJF mean concentrations and North Pacific SST, East Asian 2-m temperature, and the linear trend for MAM mean concentrations. Our analysis show that anomalous circulation patterns accompanied with the predictors, which may influence the PM 2.5 concentrations through changes in atmospheric transport and ventilation. However, these processes were not demonstrated in this study. A recent clustering analysis-based study revealed that five different circulation patterns are responsible for the occurrence of extreme PM 10 concentration events over the Korean Peninsula (Ku et al. 2021) . We speculate that Fig. 9 A scatterplot of the observed and forecast seasonal mean PM 2.5 concentrations for (a) DJF and (b) . The open red circles indicate forecasts at 1-month lead time, while the blue crosses and green dots represent the 2-and 3-month lead times, respectively. The black diagonal line denotes a 1:1 line between the observations and the forecasts. The colored dashed lines represent the regression lines of the two variables (red, blue, and green for 1-3-month lead times) PM 2.5 concentration is also governed by multiple patterns of meteorological conditions via its physical transport and chemical processes. Analyses of daily evolution are required because the processes may be obscured in seasonal mean fields.

Our models benefit from the downward trend in seasonal mean PM 2.5 concentrations (black lines in Fig. 8) . Enhancement in forecast accuracy by the linear trend can be quantified if the model is built without the trend (not shown). For DJF, the accuracy of the forecast without the trend begins at 0.66 for 1-month lead and drops to 0.62 and 0.40 for 2-3 months leads. The values for MAM are 0.82, 0.62, and 0.60 for 1-3-month forecasts without the trend. By comparing the values, we find that an increase in skills of roughly 0.1-0.3 owes to the existence of linear trends, which are approximately − 0.24 µg m − 3 year − 1 and − 0.33 µg m − 3 year − 1 in DJF and MAM PM 2.5 concentrations, respectively.

We may apply the same procedure to configure prediction models for number of high PM 2.5 days. It is inherently more challenging to predict extreme metrics than seasonal mean values because the former targets the tails of distributions. However, the near-linear relationship between seasonal mean concentrations and numbers of high concentration days, particularly during MAM (Fig. 3b) , indicates a potential of the model to skillful predictions of the number of high PM 2.5 days. For DJF, a reselection of predictors can be expected given the weak linear relationship between the two predictands (Fig. 3a) . Nonetheless, the number of high PM 2.5 days in DJF is relatively well spread out unlike the DJF mean PM 2.5 concentration clustered near its climatological mean (Fig. 3a) .

The immense impact of Coronavirus disease 2019 on regional air quality could undermine the accuracy of statistical model-based predictions. This is because the models do not directly take emission information into account, although one may argue that some of the emission information is concealed in the PM 2.5 concentration and high PM 2.5 day trends. Aerosol optical depth, a proxy for PM concentration, declined between February and March in 2020, even with stagnant atmospheric conditions, indicating a strong pause in anthropogenic emissions (Koo et al. 2020 ). Observations of PM 2.5 concentrations indeed showed unprecedented seasonal mean values of 23.64 µg m − 3 and 17.98 µg m − 3 over the Korean Peninsula during DJF 2019/2020 and MAM 2020, respectively. Similarly, the mean value was even lower at 21.80 µg m − 3 for DJF 2020/2021, while official records were not yet available for MAM 2021 during the time of this research. Statistical models may be trained directly using emission information. However, ensuring emission inventories up to date could take several years, making the use of emission information as a predictor not possible soon. This is also an obstacle for atmospheric chemistry models that require emission data as an input to initialize their simulations.

Previous climate prediction studies have outlined ways that regression model can be improved. One recommended improvement is to include uncertainty-related information when switching a deterministic forecast to a probabilistic forecast (Chang et al. 2021) . Specifically, the variance of a forecast field can be estimated using the errors in the forecast, which then can be used as a parameter to build a probabilistic distribution of the forecast field. This allows the model to provide probabilities for each category. Another way to improve statistical models is to bridge with numerical forecast. Kim et al. (2021) demonstrated how a numerical forecast of slowly varying atmospheric mode of variability, such as Arctic Oscillation, can be integrated into a statistical model to extend its surface air temperature forecast skill. Similarly, numerical forecasts of atmospheric blocking can be utilized to predict PM 10 concentrations several days ahead ).

The online version contains supplementary material available at https:// doi. org/ 10. 1007/ s13143-022-00275-4.

Long-range transport influence on key chemical components of PM2.5 in the Seoul Metropolitan Area

Review of the governing equations, computational algorithms, and other Components of the Models-3 Community Multiscale Air Quality (CMAQ) modeling system

PM10 and PM2.5 real-time prediction models using an interpolated convolutional neural network

Evaluation of subseasonal impacts of the MJO/BSISO in the East Asian extended summer

Cold-season atmospheric conditions associated with sudden changes in PM10 concentration over Seoul

Bootstrap sea ice concentrations from nimbus-7 SMMR and DMSP SSM/I-SSMIS, Version 3. NASA National Snow and Ice Data

The ERA5 global reanalysis

Statistical predictability of wintertime PM2.5 concentrations over East Asia using simple linear regression

Winter monsoon variability and its impact on aerosol concentrations in East Asia

A new indicator on the impact of large-scale circulation on wintertime particulate matter pollution over China

Possible Link Between Arctic Sea Ice and January PM10 Concentrations in South Korea

Classification of wintertime atmospheric teleconnection patterns in the Northern Hemisphere

The implication of the air quality pattern in South Korea after the COVID-19 outbreak

Classification of large-scale circulation patterns and their spatio-temporal variability during High-PM10 events over the Korean Peninsula

Contributions of international sources to PM2.5 in South Korea

Trends of PM2.5 and chemical composition in Beijing

Analysis of a severe PM2.5 episode in the Seoul Metropolitan area in South Korea from 27

PM10 data assimilation over south Korea to Asian dust forecasting model with the optimal interpolation method

Dominance of large-scale atmospheric circulations in long-term variations of winter PM10 concentrations over East Asia

An analysis on the concentration characteristics of PM2.5 in Seoul

Global modes of sea surface temperature variability in relation to regional climate indices

Long-lead predictions of warm season droughts in South Korea Using North Atlantic SST

How does the SST variability over the western North Atlantic Ocean control Arctic warming over the Barents-Kara Seas?

Global analyses of sea surface temperature, sea ice, and night marine air temperature since the late nineteenth century

Analysis of summertime atmospheric transport of fine particulate matter in Northeast Asia

Predictability of PM2.5 in Seoul based on atmospheric blocking forecasts using the NCEP global forecast system. Atmos. Environ

Health burden attributable to ambient PM2.5 in China

What we breathe impacts our health: improving understanding of the link between air pollution and health

Statistical methods in the atmospheric sciences

An improved deep learning model for predicting daily PM2.5 concentration

Increase in winter haze over eastern China in recent decades: Roles of variations in meteorological parameters and anthropogenic emissions

Impacts of different types and intensities of El Niño events on winter aerosols over China

The effects of spring and winter blocking on PM10 concentration in Korea

Acknowledgements This research was supported by the National Institute of Environment Research of Korea through NIER-2020-01-02-023 and the National Research Foundation of Korea through grants NRF-2018R1A6A1A08025520 and NRF-2019R1C1C1003161.

Competing Interest The authors declare no competing financial interests or personal relationships that influence the work reported in this paper.