Quantifying the contribution of atmospheric circulation to precipitation variability and changes in the US Great Plains and southwest using self-organizing map–analogue

The Great Plains and southwest regions of the US are highly vulnerable to precipitation-related climate disasters such as droughts and floods. In this study, we propose a self-organizing map–analogue (SOMA) approach to empirically quantify the contribution of atmospheric moist circulation (mid-tropospheric geopotential and column moisture transport) to the regional precipitation anomalies, variability, and multi-decadal changes. Our results indicate that moist circulation contributes significantly to short-term precipitation variability, accounting for 54 %–61 % of the total variance in these regions, though these contributions vary significantly across seasons. As indicated in previous research, Pacific Decadal Oscillation (PDO) is one of the major climate modes influencing the long-term multi-decadal variation in precipitation. By contrasting three multi-decadal periods (1950–1976, 1977–1998, 1999–2021) with shifting PDO phases and linking the phase shift to self-organizing map (SOM) nodes, we found that circulation changes contribute considerably to the multi-decadal changes in precipitation anomaly in terms of the mean and days of dry and wet extremes, especially for the southern Great Plains (GP) and southwest. However, these circulation-induced changes are not totally related to the PDO phase shift (mostly less than half) since internal variability or anthropogenically induced changes in circulation can also be potential contributors. Our approach improves upon flow analogue and SOM-based methods and provides insights into the contribution of atmospheric circulation to regional precipitation anomalies and variability.

1 Introduction

The United States Great Plains (GP) and southwest (SW) are central to national and global agricultural production (Elias et al., 2016; Parton et al., 2015). As the drier half of the continental United States (CONUS) with high variability in precipitation, the GP and SW rely on irrigation much more heavily than the eastern states, and they are known to be susceptible to climate disasters such as extreme droughts, wildfires, storms, and flooding, which are all listed among the top billion-dollar weather and climate disasters (NOAA NCEI, 2023). Understanding the driving factors and mechanisms for the occurrence of precipitation anomalies is crucial for the evaluation and improvement in current prediction systems as well as water resource management.

Past studies have reported various mechanisms that link to precipitation variability in North America, including anthropogenic forcings (Kirchmeier-Young and Zhang, 2020), internal atmospheric variabilities (McKinnon and Deser, 2021), El Niño–Southern Oscillation (ENSO) (Hu and Huang, 2009), Pacific Decadal Oscillation (PDO) (Lehner et al., 2018), Pacific North America (PNA) (Zhuang et al., 2021a), and Atlantic sea surface temperature (Hu and Feng, 2008), each of which can contribute differently for different locations, timescales, seasons, or even individual events. These internal variabilities and large-scale remote forcings can influence the local precipitation either by affecting the atmospheric circulation conditions above or through local and remote land surface feedbacks, both of which control the vertical air motion and moisture supply that finally lead to precipitation.

Many recent studies have argued that the long-term precipitation changes in these regions, at least from the mid-20th century until now, are primarily driven by PDO-induced atmospheric teleconnection to North America instead of anthropogenic warming, as evidenced by the relatively dry periods with negative PDO phase (1940s–1970s and 1999 onwards) when many of the major drought events occurred (Hoerling et al., 2014; Seager et al., 2014; Mankin et al., 2021) and a relatively wet period with positive phase (1980s–1990s). Yet, it is not clear quantitatively how much these long-term changes in precipitation can be explained by the PDO-related circulation variability.

The self-organizing map (SOM) or Kohonen Map (Kohonen, 1990), an unsupervised neural-network-based clustering method, has long been demonstrated as an effective and powerful tool to study dynamics or atmospheric circulation patterns related to precipitation variability, drought, or other atmospheric/oceanic phenomena (Liu and Weisberg, 2011; Zhuang et al., 2020). Different from other clustering methods, such as K means and hierarchical clustering, SOM has the advantage of sorting samples into types (SOM nodes) connected to each other in a “map” that preserves the topological structure of the data, so similar types are close to each other. Based on SOM, Cassano et al. (2007) proposed a quantitative partitioning method (referred to as “C2007” hereafter) to determine the circulation contributions to anomalies or trends of surface parameters such as precipitation. The idea is that the total precipitation anomaly for a certain node K during a period ($P_K^{\prime }$) can be decomposed as the sum of a circulation or dynamic component ($P_{K,\mathrm{dyn}}^{\prime }=\mathrm{\Delta }{f}_K{P}_K^{\prime }$) controlled by the mean frequency change in SOM nodes, a thermodynamic component ($P_{K,\mathrm{the}}^{\prime }=f_K\mathrm{\Delta }{P}_K^{\prime }$) related to the mean state change in precipitation for the same SOM node, and an interaction term ($P_{K,\mathrm{int}}^{\prime }=\mathrm{\Delta }{f}_K\mathrm{\Delta }{P}_K^{\prime }$) related to both. Detailed descriptions of this method can be found in Cassano et al. (2007) and Horton et al. (2015). However, this partitioning method assumes that $P_{\mathrm{dyn}}^{\prime }$ at the timescale of individual sample (e.g., daily) stays the same for circulation patterns classified to the same SOM node. This underlying assumption leads to the result that the determination of $P_{\mathrm{dyn}}^{\prime }$ could be sensitive to the choice of SOM node number. An additional inconvenience related to this method is that the SOM analysis is often conducted in certain seasons as the circulation–precipitation relationship could vary in different seasons.

Flow analogue or dynamic adjustment (e.g., Deser et al., 2016; Horowitz et al., 2022; Jezequel et al., 2018; Lehner et al., 2018; Terray, 2021; Yiou et al., 2007) is another widely used method (referred to as analogue hereafter) to determine the circulation contribution. For a certain day (or other timescales), different from the C2007 method, the constructed flow analogue (Lehner et al., 2018) method searches for a limited number of days (analogue days) with circulation patterns most similar to that of the target day instead of days assigned to the same SOM node in C2007; it then combines precipitation anomaly of these analogue days to form a final analogue anomaly or dynamic component ($P_{\mathrm{dyn}}^{\prime }$), and the residual term which cannot be explained by circulation is defined as the thermodynamic component ($P_{\mathrm{the}}^{\prime }$). Detailed steps and variants of the analogue method can be found in the literature (e.g., Zhuang et al., 2021b). The analogue theoretically better estimates the dynamic contribution by utilizing circulation patterns from the analogue days and accounting for their subtle difference.

In this study, we proposed a modified approach to quantify the multivariate circulation contribution to P^′ in terms of its daily variability, which combines the advantages of both SOM and constructed analogue methods. We further use this approach to quantify the contribution of overall circulation changes and potentially PDO-related circulation changes to the long-term multi-decadal changes in P^′.

The datasets and methods used in this study are described in Sect. 2. All the main results, including the SOM inferred relationship between large-scale atmospheric circulation and precipitation and SOM–analogue-based quantification of circulation contributions to daily variability and long-term changes in precipitation anomalies, are presented in Sect. 3. Section 4 summarizes the main conclusions of this study and also gives a brief discussion about the implications of this study and limitations of our methodology and results.

2 Data and method

2.1 Data

In previous studies related to SOM and analogue, large-scale circulation is generally represented by mean sea level pressure (SLP) or geopotential height at 500 hPa (Z₅₀₀). Here, we choose Z₅₀₀ over SLP as our experiments have demonstrated analogues derived from Z₅₀₀ show greater similarity in synoptic variability with observed surface anomalies and yield smaller residuals compared to analogues derived from SLP (Zhuang et al., 2021b). Additionally, impacts of anomalous atmospheric moisture transport represented by vertically integrated vapor transport (IVT) on precipitation variability and drought are also suggested in the literature (Zhuang et al., 2020; Pu et al., 2016; Swales et al., 2016). The 3-hourly Z₅₀₀ and IVT data at a $\mathrm{1}{}^{\circ }\times \mathrm{1}{}^{\circ }$ resolution are obtained from the fifth-generation ECMWF atmospheric reanalysis of the global climate (ERA5) (Hersbach et al., 2020). We use daily averages for all months from 1950 to 2021 over the CONUS area (25–55^∘ N, 130–70^∘ W) as input for later analyses. To account for the thermal expansion of the warming atmosphere, we subtract the daily global area-weighted mean Z₅₀₀ from the daily Z₅₀₀ data at each grid point so that there is no linear trend in the Z₅₀₀ data due to warming (Christidis and Stott, 2015; Siler et al., 2019; Zhuang et al., 2021b); however, this procedure does not remove any anthropogenically forced changes in atmospheric circulation patterns. Other studies (e.g., Sippel et al., 2019) also use linear detrending or high-pass filtering to preprocess Z₅₀₀ data at each grid point, assuming that the forced circulation response was smooth and additive; we do not employ these approaches due to our short observational record, and they could remove some of the decadal circulation trend due to internal variability (Zhuang et al., 2021b). To mitigate high-frequency synoptic noise, we employ a simple 5 d moving average filter to both Z₅₀₀ and IVT. Subsequently, we calculate the daily standardized anomaly ($Z_{\mathrm{500}}^{\prime }$ and IVT^′) using the 5 d moving average-filtered data relative to the 1950–1999 climatology. We select the 1950–1999 period as the climatology reference to maximize the utilization of available data for this study; furthermore, this period, which predates the significant warming trend, typically serves as a robust baseline with less climate change impact. In addition, we adopt the running mean approach used in Zhuang et al. (2020) to calculate a smoothed daily climatological mean and standard deviation.

For precipitation, the Climate Prediction Center (CPC) global unified gauge-based analysis of daily precipitation (P) is used. This dataset is at $\mathrm{0.25}{}^{\circ }\times \mathrm{0.25}{}^{\circ }$ resolution and covers the period from 1948 to the present. In this study, we mainly focus on three regions, including the southern GP (SGP; 26–38^∘ N, 105–92^∘ W), the northern GP (NGP; 38–49^∘ N, 105–92^∘ W), and the southwestern US (SW; 31–42^∘ N, 114–105^∘ W). Therefore, area-weighted mean precipitation is calculated over these three regions. A 5 d moving average filter is also applied to precipitation similar to the circulation data. The precipitation anomaly (P^′) is calculated relative to 1950–1999 climatology using the same approach stated above for circulation but without the standardization process as precipitation generally follows a gamma distribution rather than a Gaussian distribution.

Besides circulation and precipitation, other data used in this study include ERA5 2 m dew point temperature (D2), 700 hPa temperature (T700) data, and convective available potential energy (CAPE) as well as the monthly Pacific Decadal Oscillation (PDO) index provided by the National Center for Environmental Information (NCEI). The D2 and T700 data are used to calculate the convective inhibition (CIN) index (CINi = D2 − T700) (Myoung and Nielsen-Gammon, 2010) to represent the lower-atmospheric instability; ERA5 also provides CIN data, but it is not always available around the year due to its original definition from the parcel buoyancy model. Both CAPE and CINi data are processed the same way as circulation data to acquire their standardized anomalies (CAPE^′ and CINi^′) for further analyses. The PDO index is one of the most important Pacific climate variabilities and is often described as an El Niño-like pattern but with long-term persistence of up to several decades (Zhang et al., 1997). The NCEI monthly PDO index is calculated based on the NOAA extended reconstruction of sea surface temperature version 5 (ERSSTv5). A 3-month running mean PDO index is used in this study to smooth out some monthly sea surface temperature (SST) fluctuations and reduce short-term noise. Our whole analysis period 1950–2021 can be divided into three periods, each with a different PDO phase from the previous one: P1, 1950–1976 (negative); P2, 1977–1998 (positive); and P3, 1999–2021 (negative). Long-term change in precipitation in terms of its mean or days of extremes for our following analysis (Sect. 3.3) is defined as the difference between two adjacent periods (P2 − P1 and P3 − P2).

2.2 Method

2.2.1 Multivariate SOM

In this study, the large-scale atmospheric circulation condition is jointly represented by $Z_{\mathrm{500}}^{\prime }$ and IVT^′, as both of them have been demonstrated to be important for the precipitation variability over North America. The data matrices of these two variables are concatenated along the spatial dimension, so the input data matrix has a size of ($\mathrm{2}\times N_{\mathrm{lon}}\times N_{\mathrm{lat}}$, N_t). Specifically for this study, N_lon=61, N_lat=31, and N_t=365 d yr${}^{-\mathrm{1}}\times \mathrm{72}$ years (leap days removed for simplicity), so the input size is (3782, 26 280).

To speed up the SOM training process, we employ a principal component analysis (PCA) preprocessing technique (Zhuang et al., 2020). PCA, or empirical orthogonal functions (EOFs), is a widely used approach for data dimension reduction or extracting spatial modes of variability (Reusch et al., 2005). Here, PCA is used to decompose the high-dimensional circulation data (two variables and thousands of grid points) into the matrix multiplication product of some low-dimensional principal components (PCs) and their corresponding loading vectors (modes). By only retaining the top PCs and modes that explain the majority of the variance in the original data, we can use these PCs as reduced-size input data for SOM. Specifically in this study, we retain the top 26 PCs, which explain 87.3 % of the combined variance in $Z_{\mathrm{500}}^{\prime }$ and IVT^′ (or 90.6 % and 77.4 % of the variance in $Z_{\mathrm{500}}^{\prime }$ and IVT^′, respectively).

Determination of SOM node number and evaluation of the trained SOM quality are often complicated issues. Too few SOM nodes could lead to an underrepresentation of samples, while too many nodes likely result in redundancy and create difficulties for analysis and visual interpretation. In this study, we have tested multiple SOM schemes with various node numbers up to 200. Three metrics are calculated for the trained SOM of each scheme to assess the SOM quality, including the quantization error (QE), topographic error (TE), and combined error (CE). QE represents the average Euclidean distance between each data vector and the SOM node it is assigned to, or the so-called best matching unit (BMU). TE is the proportion of all input data vectors with first and second BMUs not adjacent to each other in the SOM map. CE is a metric that combines the concepts of both QE and TE and was defined by Kaski and Lagus (1996).

SOM training and error metric calculation are performed using the MATLAB SOM Toolbox (Vatanen et al., 2015). Training length is fixed at 100 epochs (50 for rough training and 50 or fine-tuning) for all schemes. Additional experiments show that increasing training length can reduce QE, but the improvement is very marginal after around 100 epochs (not shown). In the SOM Toolbox, after inputting a desired node number, the SOM map shape is automatically determined. During the initialization before the training, we employ the default setting, i.e., linear initialization, where weight vectors are initialized in a linear manner along the subspace defined by the two principal eigenvectors of the input dataset. This choice is made to facilitate the reproducibility of our results. While SOM is generally robust to initialization, slight variations in outcomes may occur when using random initialization, where the weight vectors are initialized with random small values. To assess the impact of initialization on SOM error metrics, we also conduct 100 instances of additional SOM training with random initialization for each node number setting and analyze the range of their error metrics. After training SOMs for all schemes with various node numbers and calculating their quality metrics, the result in Fig. 1 shows that, overall, QE decreases with node number as more nodes mean more detailed classification of circulation patterns, and thus they are more representative of individual samples; TE increases with node number as a more complex network often leads to a larger topographic structure; CE, which combines the concepts of both quantization and topographic errors, decreases with node number first when the number is relatively small and then exhibits a larger fluctuation than QE and TE. A 7×4 SOM scheme is chosen for our following analysis as this scheme has a relatively smaller node number that facilitates the visualization and analysis of the results and also has a lower CE compared to schemes with similar node numbers.

Figure 1 (a) Quantization error (QE), (b) topographic error (TE), and (c) combined error (CE) of SOM schemes with different node numbers (N). The red cross denotes the one (N=28) we selected in this study. Solid lines represent the errors in the SOMs trained with linear initialization, while the shaded areas represent the range of errors within the middle 95 % of the distribution from the 100 SOMs trained with random initialization for each node number setting, and dashed lines are their averages.

2.2.2 Empirical quantification of moist circulation contribution: the SOM–analogue (SOMA) method

As stated in the Introduction, the C2007 dynamic/thermodynamic partition method assumes that the dynamic component $P_{\mathrm{dyn}}^{\prime }$ for days grouped under the same BMU node is determined by their mean P^′. However, this approach can underestimate variability in true $P_{\mathrm{dyn}}^{\prime }$, especially when the node number is small and the circulation pattern of a BMU is less representative of the assigned samples. On the other hand, the flow analogue method (e.g., Yiou et al., 2007; Zhuang et al., 2021b) estimates $P_{\mathrm{dyn}}^{\prime }$ for a given day by identifying historical days with similar atmospheric circulation conditions (analogue days) and combining the P^′ information from these selected analogue days. Specifically, the flow analogue method employed by Zhuang et al. (2021b) follows three main steps. First, for a particular day, analogue days with similar circulation conditions are chosen by minimizing a distance function between the circulation field (e.g., $Z_{\mathrm{500}}^{\prime }$) on the target day and that on historical days within a specific calendar range centered around the target day. Second, the top N analogue patterns with the smallest-distance function are linearly combined to create a “constructed analogue” pattern that resembles the circulation pattern of the target day. Finally, the N coefficients from the linear combination are applied to P^′ values of the N analogue days to obtain $P_{\mathrm{dyn}}^{\prime }$. The constructed flow analogue method generally explains a larger portion of the observed P^′ variance, and its estimate is not sensitive to SOM node number, so it has a higher accuracy in capturing the complex relationship between atmospheric circulation and precipitation anomalies and theoretically provides a better estimate of the circulation contribution compared to the C2007 method. This improvement is primarily due to the flow analogue method considering the differences among circulation patterns in individual analogue days and combining P^′ values from these analogue days. However, flow analogue is less effective than the C2007 method when evaluating the relationship between a specific circulation type and P^′.

In this study, we propose a new method called the SOM–analogue (SOMA), which leverages the strengths of both SOM (C2007) and flow analogue. SOMA aims to be capable of quantifying overall circulation contribution and the contribution of a specific circulation type to P^′ and also to provide a more robust circulation contribution estimate that is less sensitive to parameter choices and explains a large portion of observed P^′ variability. It is worth noting that attempting to separate the actual dynamic and thermodynamic components solely through a circulation clustering approach like SOM can be challenging. Each type of circulation, as represented by SOM nodes, inherently encompasses thermodynamic responses. Therefore, the distinction between “dynamic” and “thermodynamic” components can be ambiguous when using these terms. Consequently, we prefer to refer to the precipitation anomaly influenced by moist circulation patterns involving Z₅₀₀ and IVT as $P_{\mathrm{SOMA}}^{\prime }$ instead of $P_{\mathrm{dyn}}^{\prime }$, and the residual part that cannot be explained by circulation as $P_{\mathrm{RES}}^{\prime }$ instead of $P_{\mathrm{the}}^{\prime }$. This emphasizes that our results regarding circulation contributions are contingent on our chosen set of circulation variables.

The basic idea behind SOMA is to incorporate variability in $P_{\mathrm{SOMA}}^{\prime }$ for days with the same BMU node, similar to the flow analogue method, while adhering to the C2007 framework. In contrast to the flow analogue method, where analogue days are determined by minimizing the Euclidean distance of anomalous circulation patterns, in SOMA, the analogue days are directly selected from the daily samples sharing the same BMU. Additionally, the constructed flow analogue approach relies on a linear combination of precipitation anomaly from analogue days to quantify circulation contribution at a target day j; i.e., $P_{j,\mathrm{SOMA}}^{\prime }=\sum _i{a}_i{P}_i^{\prime }=\sum _i{a}_i{P}_{i,\mathrm{SOMA}}^{\prime }+\sum _i{a}_i{P}_{i,\mathrm{RES}}^{\prime }$, where the coefficients a_i are determined by the linear dependency of the circulation pattern on target day j and those on analogue days i; i.e., $Z_j^{\prime }=\sum _i{a}_i{Z}_i^{\prime }$. However, it is important to note that the analogue day precipitation anomaly $P_i^{\prime }=P_{i,\mathrm{SOMA}}^{\prime }+P_{i,\mathrm{RES}}^{\prime }$ includes both a circulation-contributed component and a residual component. This approach overlooks the potential thermodynamic influences in the precipitation anomaly on analogue days, and the resulting circulation-contributed component $P_{j,\mathrm{SOMA}}^{\prime }$ still retains a thermodynamic residual $\sum _i{a}_i{P}_{i,\mathrm{RES}}^{\prime }$, which may not be negligible with limited samples. To mitigate this issue, we make the assumption that for each group of a similar circulation pattern (same SOM node), $P_{j,\mathrm{SOMA}}^{\prime }=\sum _i{a}_i{P}_{i,\mathrm{SOMA}}^{\prime }$. This differs from the equation form in constructed flow analogue and excludes the term $\sum _i{a}_i{P}_{i,res}^{\prime }$. Following constructed flow analogue, where each circulation pattern is linearly dependent on other circulation patterns within the same SOM node, i.e., $Z_j^{\prime }=\sum _i{a}_i{Z}_i^{\prime }$, we can see that $P_{\mathrm{SOMA}}^{\prime }$ and Z^′ share the same linear form and coefficients. Therefore, $P_{\mathrm{SOMA}}^{\prime }$ and Z^′ can be considered linearly dependent as well, allowing us to treat the calculation of $P_{\mathrm{SOMA}}^{\prime }$ as a regression problem, which is different from the circulation analogue in previous studies (e.g., Deser et al., 2016; Lehner et al., 2018). Detailed steps are outlined below.

1. After PCA preprocessing, the original circulation anomaly matrix ${\mathbf{Y}}_{N_{\mathrm{p}}\times N_t}$ can be written as ${\mathbf{Y}}_{N_{\mathrm{p}}\times N_t}={\mathbf{V}}_{N_{\mathrm{p}}\times N_e}\cdot {\mathbf{X}}_{N_t\times N_e}^T$, where N_p is the number of grid points (61×31) multiplied by the number of circulation variables (2), N_t the number of days, and N_e the number of all PCs. By only retaining the top K PCs, ${\mathbf{Y}}_{N_p\times N_t}\approx {\mathbf{V}}_{N_p\times K}\cdot {\mathbf{X}}_{N_t\times K^T}$.

2. The PC matrix ${\mathbf{X}}_{N_t\times K}$ is used as input to train a SOM. For each daily sample x_t(1×K) at time t, we can find its corresponding BMU; i.e., BMU(x_t)=B_n, where B_n is one of the SOM nodes, and n (n=1, 2, 3, …, N) can be determined by minimizing the Euclidean distance between x_t and B_n. The circulation pattern for node n can then be recovered as ${\mathbf{V}}_{N_p\times K}\cdot {\mathbf{B}}_n^T$.

3. For any day t₀, assuming its BMU is node n, we find all days t_i (i=1, 2, …, M) that have this node n as their BMU and are within the 91 d calendar window centered on the day t₀ but not in the same year as t₀, to account for the seasonal dependence of the circulation–precipitation relationship; i.e., $\mathrm{BMU}\left({\mathbf{x}}_{t_i}\right)=\mathrm{BMU}\left({\mathbf{x}}_{t_{\mathrm{0}}}\right)$, $\left|\mathrm{DOY}\right(t_i)-\mathrm{DOY}(t_{\mathrm{0}}\left)\right|\le \mathrm{45}$ d, Year(t_i)≠Year(t₀).

4. Solve the regression problem: ${\mathbf{P}}_{\mathrm{1}\times M}={\mathbf{U}}_{\mathrm{1}\times K}\cdot {\mathbf{X}}_{M\times K}^T$, where P_1×M is the precipitation vectors for days t_i (i=1, 2, …, M) from the previous step and X_M×K is input PC samples for these M days. The regression coefficient vector U_1×K is then used to calculate the dynamic precipitation anomaly for the day t₀, $P_{\mathrm{SOMA},t_{\mathrm{0}}}^{\prime }={\mathbf{U}}_{\mathrm{1}\times K}\cdot {\mathbf{x}}_{t_{\mathrm{0}}}^T$.

5. Repeat steps (3) and (4) for all days until all $P_{\mathrm{SOMA}}^{\prime }$ are calculated.

To account for the sampling uncertainty, we apply a bootstrap technique in step (3) where we create multiple datasets of X_M×K by randomly sampling ${\mathbf{x}}_{t_i}$ with a replacement. After bootstrapping and repeating steps (3) to (5) 1000 times, the 1000 calculated $P_{\mathrm{SOMA},t_{\mathrm{0}}}^{\prime }$ are then averaged as the final result.

Figure 2 Composite standardized $Z_{\mathrm{500}}^{\prime }$ (contours; solid and dashed lines are for positive and negative values, respectively) and standardized IVT^′ (shading) for each node in the 7×4 SOM. Mean frequency ($\stackrel{\mathrm{‾}}{f}$) and seasonality (f_range, i.e., the difference between maximum and minimum seasonal frequencies as shown in Fig. 3; a larger value indicates stronger seasonal variation) of each node are labeled in each panel.

3 Result

3.1 Atmospheric circulation pattern clustering and link to precipitation

After the 7×4 SOM for $Z_{\mathrm{500}}^{\prime }$ and IVT^′ is trained and each daily sample is assigned to its BMU node, circulation patterns with the same BMU node are averaged as the general representation of circulation conditions for the particular node. The composite circulation anomaly patterns for all SOM nodes are shown in Fig. 2. Similar to our earlier study (Zhuang et al., 2020), the SOM map shows a continuum of anomalous circulation and moisture transport patterns over the CONUS. Regions with large $Z_{\mathrm{500}}^{\prime }$ gradients often correspond to large IVT^′ due to the geostrophic balance relationship. For example, when an anomalous geopotential low is centered over the western or central US and a geopotential high over the eastern or southeastern US (e.g., nodes A1, B1, C1, C2, D1, D2, and D3), enhanced moisture transport from the Gulf of Mexico represented by the large positive IVT^′ is often observed in between the geopotential low and high with an elevated geopotential gradient.

As the SOM is trained with daily samples of all seasons and considering the potential seasonal variation in occurrence frequency and circulation–precipitation relationship, we calculate the frequency and average P^′ for each SOM node during every 91 d calendar window during all years. The resulting seasonal variations in node frequency and mean P^′ are shown in Figs. 3 and 4, respectively. In Fig. 3, some nodes exhibit very clear seasonal variation with the range of frequency being as large as 4 % (mean frequency and seasonal frequency ranges for each node are also available in Fig. 2); they either peak in cold seasons, such as C1 and D1, with the anomalous geopotential high centered more on the southeastern US and strong moisture transport to the SGP, or peak in the warm seasons, such as C4 and D4, where the high is centered more on the northeastern US, so the moisture transport is weaker and more towards the north. Some other nodes show a relatively stable frequency throughout the year, with a range of less than 1 %, such as nodes C6, C7, and D7, where the CONUS is mostly dominated by an extensive geopotential high.

Figure 3 Seasonal variation in frequency for SOM nodes. For each node in each panel, the frequency at a certain calendar day represents the node frequency during a 91 d calendar window centered on that calendar day for all years (1950–2021).

Figure 4 Seasonal variation in mean P^′ averaged over each region (blue for SGP, red for NGP, and yellow for SW) for each SOM node. For each node in each panel, the mean P^′ at a certain calendar day represents the P^′ of that node during a 91 d calendar window centered on that calendar day for all years (1950–2021). The shading represents the interquartile range of all daily P^′ used to calculate the mean.

Similarly, precipitation also shows apparent seasonal variations for certain nodes (Fig. 4). For example, node D1, located at the lower-left corner of the SOM, brings the strongest positive IVT^′ from the Gulf of Mexico and the Gulf of California of all nodes, and it is one of the major wet nodes related to positive P^′ in all three regions during the cold season. However, during the warm season, node D1 is more of a neutral or weakly dry node as it is linked to a small negative P^′ in the SGP; for the NGP, this node is still a wet node and the average P^′ is larger than that during the cold season. This can in part be explained by the seasonal difference in precipitation mechanisms. In the warm season, precipitation is mostly convection-driven, reducing the impact of large-scale vertical air motion associated with certain circulation patterns. Additionally, the composite CAPE^′ and CINi^′ map for node D1 (Fig. S1 in the Supplement) shows relatively smaller CAPE^′ and a more stable lower troposphere (negative CINi^′) compared to other nodes in the warm season (May to July), suppressing convective development and thus limiting precipitation regardless of the strong moisture transport, whereas in the cold season (November to January), CAPE^′ and CINi^′ shows larger positive values associated with more precipitation compared to other nodes. However, it is important to note that the above comparison is within the same season, as node D1 in the warm season still corresponds to higher traditional CAPE and CINi values (not anomalies), as well as less total precipitation, compared to the cold season. Similar situations occur to many other nodes, either exhibiting opposite signs of mean P^′ in different seasons or clear seasonal variations with the same sign. Overall, these seasonal fluctuations in the circulation–precipitation relationship can be attributed to factors influencing the types and mechanisms of precipitation, which can vary seasonally, such as atmospheric stability and thermodynamic feedback (Myoung and Nielsen-Gammon, 2010), or to potential slight sampling bias in different seasons, as evidenced by the contrast between composite circulation maps for the same node in different seasons (Fig. S1).

Figure 5 (a) Seasonal variation in the percentage of P^′ variance explained by $P_{\mathrm{SOMA}}^{\prime }$ ($R_{\mathrm{SOMA}}^{\mathrm{2}}$) for the three regions; each data point is calculated within a 91 d calendar window. Thin dashed lines represent the value for all year round. (b) Same as (a) but for seasonal variation in the percentage of moist circulation (SOMA) contribution determined by fitting a simple linear regression model between $P_{\mathrm{SOMA}}^{\prime }$ and P^′ during a 91 d calendar window ($P_{\mathrm{SOMA}}^{\prime }=\mathit{\beta }{P}^{\prime }$).

3.2 Moist circulation contributions to precipitation

The above analyses show that strong seasonal variations could exist for either frequency or mean P^′ in certain nodes due to different types of precipitation and factors other than large-scale circulation and moisture transport that contribute to precipitation. Therefore, assuming node-mean P^′ as dynamic components calculated using the C2007 method may result in an incorrect representation of dynamic contribution and underestimation of its variability.

The second part of our analysis is to employ the SOMA method described in the Method section to estimate the contribution of dynamics and thermodynamics to P^′. We decompose the daily average P of certain regions (SGP, NGP, and SW) into the sum of a moist circulation-contributed component $P_{\mathrm{SOMA}}^{\prime }$, which is explained by the mid-tropospheric circulation and moisture transport, and a residual component $P_{\mathrm{RES}}^{\prime }$, which ideally should be independent of circulation condition and related to thermodynamic processes such as land surface feedbacks. Further analysis shows that composite seasonal cycles of $P_{\mathrm{SOMA}}^{\prime }$ for different nodes are similar to that of P^′ (Fig. 4) but with a narrower interquartile range (not shown). In contrast, $P_{\mathrm{RES}}^{\prime }$ does not show a clear season cycle for any node (not shown). This indicates the majority of precipitation variation at daily to synoptic scales can be mainly explained by circulation, and thermodynamic residual does not have a preference for certain nodes in any season. The $P_{\mathrm{SOMA}}^{\prime }$ calculated by our SOMA method overall is highly correlated to P^′ in all seasons, indicating that circulation and moisture transport conditions can explain a large portion of precipitation variability.

Figure 5a shows the percentage of the variance of P^′ explained by the $P_{\mathrm{SOMA}}^{\prime }$ ($R_{\mathrm{SOMA}}^{\mathrm{2}}$). NGP, which is located more to the north and more susceptible to frontal weather systems, has the largest mean $R_{\mathrm{SOMA}}^{\mathrm{2}}$ of 60.5 % for all year round, while SW and SGP have slightly smaller values of 56.3 % and 53.6 %, respectively. These three regions also show very clear seasonal variations in $R_{\mathrm{SOMA}}^{\mathrm{2}}$ as well: regions located more to the south, like the SGP and SW, show a maximum of ∼60 % and ∼70 % in winter (December to January) when precipitation is more affected by frontal systems and a minimum of ∼45 % and ∼40 % in summer (July to August) when precipitation is mostly convectively driven and surface thermodynamics have more influence. On the other hand, the more northern region NGP shows continuously higher $R_{\mathrm{SOMA}}^{\mathrm{2}}$ in the warm season (March to October); its seasonal variation in $R_{\mathrm{SOMA}}^{\mathrm{2}}$ shows two local minima (one in August, ∼55 %, and the other in January, ∼50 %) and two local maxima (one in April–May, ∼65 %, and the other in October–November, ∼64 %).

In addition to daily precipitation variability, we fit a simple linear regression model $P_{\mathrm{SOMA}}^{\prime }=\mathit{\beta }{P}^{\prime }$ for each 91 d calendar window in each region to estimate the relative percentage dynamic contribution (β) to intensity of the precipitation anomaly. Figure 5b shows the seasonal variation in β for all samples at a daily scale for each region. Overall, the percentage of SOMA contributions in the SGP and SW shows a similar seasonal cycle as the $R_{\mathrm{SOMA}}^{\mathrm{2}}$, with the highest β (∼69 % and ∼73 %) in the cold season (December to January) and the lowest (∼56 % and ∼50 %) in the warm season (May to August), while the NGP has a mean β of ∼68 % with a smaller range of seasonal variation (∼65 % to ∼75 %).

Figure 6 further shows the historical evolution of 12-month moving average P^′ and its moisture circulation contributed (SOMA) and residual (RES) components from 1950 to 2021 to demonstrate how they contribute to the longer-term variability in precipitation. Overall, $P_{\mathrm{SOMA}}^{\prime }$ matches better with P^′ than $P_{\mathrm{RES}}^{\prime }$ does for all three regions, even for some P^′ extremes during major drought and flooding events. Percentage variance explained by the SOMA component is 87.3 %, 82.1 %, and 82.8 % for SGP, NGP, and SW, respectively, which are considerably higher than that explained by the residual (62.6 %, 77.2 %, and 65.5 %). This suggests atmospheric circulation and moisture transport are the main controls for precipitation variations at an interannual scale as well.

Figure 6 Twelve-month moving average P^′ (assigned to the last month of a 12-month period) for the three regions. Black, blue, and red lines represent the total anomaly P^′, the moist-circulation-contributed component $P_{\mathrm{SOMA}}^{\prime }$, and the residual component $P_{\mathrm{RES}}^{\prime }$, i.e., $P^{\prime }-P_{\mathrm{SOMA}}^{\prime }$.

To further demonstrate what type of atmospheric circulation patterns are most related to regionally dry and wet anomalies, Fig. 7 ranks the nodes according to their mean seasonally cumulative $P_{\mathrm{SOMA}}^{\prime }$ and $P_{\mathrm{RES}}^{\prime }$ during four seasons, including January to March (JFM), April to June (AMJ), July to September (JAS), and October to December (OND). As these seasonal contributions are averaged over all years, many nodes show a clear tendency to be dynamically related to either dry or wet anomalies during a particular season in a region, but almost all of them have negligible residual components. Consistently through all three regions, the top four nodes that contribute most to dry anomalies in the cold season (JFM and OND) are A7, B7, C7, and D7, all located on the right side of the SOM map, featuring either a west high-east low or high dominant geopotential pattern over the CONUS with reduced IVT^′ in both SGP and SW. On the other hand, for the wet anomalies, top contributing nodes are mostly located on the left or bottom side of the SOM map, featuring a west low-east high geopotential pattern favoring increased IVT^′.

Figure 7 Average moist circulation contribution (SOMA; a–c) and residual (RES; d–f) to P^′ from each node in each 3-month season (January–March, April–June, July–September, and October–December) over each region (SGP, NGP, SW). The colors of nodes are assigned in such a manner that adjacent nodes exhibit greater similarity in color compared to nodes that are farther apart.

In the warm season (AMJ and JAS), top dry and wet nodes all change a bit for the three regions. For example, for the SGP, besides the high dominant pattern D6 and 7, a few west low-east high patterns (C3, D2, D3) with increased IVT^′ also contribute considerably to $P_{\mathrm{SOMA}}^{\prime }$ due to the increased lower-atmospheric stability related (CIN^′) to them. For the SW, A1 becomes the top contributing node in its monsoon season due to the weak geopotential gradient there and the strong geopotential high on the west coast, all of which reduce the main moisture source, including that from the Gulf of Mexico and the Gulf of California.

3.3 Long-term changes in precipitation

The results from the previous section suggest that dynamic factors, specifically circulation and moisture transport, account for the majority of variance in precipitation anomalies from daily to interannual variability. Our next objective is to investigate whether there are long-term changes in precipitation and, if so, how these changes can be explained by circulation changes, which are mainly related to internal variability in the atmosphere and oceanic forcings. The PDO has been demonstrated as a main oceanic forcing that controls the multi-decadal variability in precipitation over North America, as the PDO can generate atmospheric circulation patterns that favor dry or wet conditions depending on its phases (Lehner et al., 2018). Our period of analysis, spanning from 1950 to 2021, encompasses three periods with shifting PDO phases, including P1 (1950–1976, mostly negative), P2 (1977–1998, mostly positive), and P3 (1999–2021, mostly negative).

Figure 8 Mean P^′, $P_{\mathrm{SOMA}}^{\prime }$, and $P_{\mathrm{SOMA}\text{-}\mathrm{PDO}}^{\prime }$ during three periods with shifting PDO phases (P1: 1950–1976; P2: 1977–1998; P3: 1999–2021) in the three regions (a) SGP, (b) NGP, and (c) SW. Error bars in (a) and (b) represent the 95 % confidence intervals of the mean P^′ and $P_{\mathrm{SOMA}}^{\prime }$.

Figure 8 (first column in each panel) shows the period mean P^′ for different regions and contrasts them in three periods with different PDO phases. We can see that, overall, all three regions show a significant increase in P^′ from P1 to P2, followed by a decrease from P2 to P3. The SOMA component $P_{\mathrm{SOMA}}^{\prime }$ (Fig. 8, second column in each panel), though it generally exhibits smaller-amplitude anomalies, captures the multi-decadal changes in the P^′ in the three periods (from 14 % to 82 % for all seasons, Table 2), especially for the SW (75 % and 82 % of the changes in P1 − P2 and P2 − P3).

Figure 9 Frequency-weighted $P_{\mathrm{SOMA}}^{\prime }$ difference between positive and negative PDO phases (monthly PDO index > 0.5 and $<-\mathrm{0.5}$) for different nodes in different months. For each row or month, nodes marked with a circle (with a plus sign) indicate that the frequency differences are significant positive or negative at the 0.10 (0.05) level, and they can be considered PDO-related nodes for that month.

To further determine how much of the dynamic contribution to the multi-decadal changes in P^′ can be linked to PDO changes, we first investigate whether the moist circulation and precipitation anomalies associated with each node may be modulated by PDO phases. We assess whether the difference in frequency-weighted $P_{\mathrm{SOMA}}^{\prime }$ between positive and negative PDO phases is statistically significant through a Monte Carlo test. Detailed steps are as follows. (a) For node k in month m, calculate node frequency for the positive PDO phase (monthly PDO index > 0.5) and the negative phase (monthly PDO index < −0.5) for all years, i.e., $f_{\mathrm{PDO}+}(k,m)$ and $f_{\mathrm{PDO}-}(k,m)$, as well as average SOMA-contributed P^′, i.e., $P_{\mathrm{SOMA},\mathrm{PDO}+}^{\prime }(k,m)$ and $P_{\mathrm{SOMA},\mathrm{PDO}-}^{\prime }(k,m)$. (b) The frequency-weighted $P_{\mathrm{SOMA}}^{\prime }$ difference between positive and negative PDO phases for node k in month m can be calculated as ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }(k,m)=f_{\mathrm{PDO}+}(k,m)\times P_{\mathrm{SOMA},\mathrm{PDO}+}^{\prime }(k,m)-f_{\mathrm{PDO}-}(k,m)\times P_{\mathrm{SOMA},\mathrm{PDO}-}^{\prime }(k,m)$, and the overall PDO impact for node k can be represented by summing this difference across all months, i.e., ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }\left(k\right)=\sum _{m=\mathrm{1}}^{\mathrm{12}}\left[{\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }\right(k,m\left)\right]$. (c) Randomly shuffle the sequence of years for the PDO index data and recalculate the abovementioned statistics. (d) Repeat step (c) many times (10 000 here) if ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }(k,m)$ is larger than 90 % of all simulated values, and then for node k in month m, $P_{\mathrm{SOMA}}^{\prime }$ is significantly larger during positive PDO phase than negative phase; in contrast, if ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }(k,m)$ is smaller than 90 % of all simulated values, then $P_{\mathrm{SOMA}}^{\prime }$ is considered significantly smaller during positive PDO phase. (e) Similar conclusions can be drawn for node k overall by assessing ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }\left(k\right)$ instead of ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }(k,m)$. The resulting nodes statistically related to PDO (for each month and overall) are shown in Fig. 9. Additionally, the selection of a 0.5 threshold for PDO phases is mainly based on two considerations: (1) it has been commonly used in some prior studies (e.g., Hu and Guan, 2018; Kiem et al., 2003); (2) days or months categorized as positive or negative PDO phases using the 0.5 threshold constitute about two-thirds of the total samples (66.2 %), striking a balance between inclusivity and specificity. Using a smaller threshold, such as 0.25, would result in a much higher percentage (81.8 %) of samples categorized as positive or negative PDO-related, leading to results that are less representative of the true PDO impact due to over-inclusiveness. Conversely, a larger threshold, like 1.0, would yield fewer (40.6 %) samples, although it is still feasible, which could introduce greater sampling uncertainty due to limited data availability. Results obtained using the 1.0 threshold can be found in the Supplement (Figs. S2–S4) for reference; overall, the 0.5 and 1.0 thresholds produce slightly different results, but these differences do not impact our subsequent discussion.

A significant ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }$ value could be due to differences in node frequency, node-mean $P_{\mathrm{SOMA}}^{\prime }$, or both between the two PDO phases. Our result in Fig. 9 indicates that, for most nodes, PDO phases can have very different influences on precipitation through moist circulation conditions identified in this study, for different seasons and regions. For example, node A1, characterized by geopotential high anomalies along the west and east coasts and low anomalies over the central US, mirrors a Pacific North America (PNA)-like pattern, previously linked to Pacific variability and its impact on precipitation in the western US and GP (e.g., Ciancarelli et al., 2014; Zhuang et al., 2021a). Interestingly, node A1 does not display a significantly higher frequency during the positive PDO phase compared to the negative phase (even lower in January to April; Fig. 10). However, it stands out as one of the nodes exhibiting the most substantial contrast in $P_{\mathrm{SOMA}}^{\prime }$ between the two PDO phases over the SGP. This contrast arises because the composite circulation pattern for node A1 during the positive PDO phase exhibits stronger negative $Z_{\mathrm{500}}^{\prime }$ over the central US, which can promote upward motion and is generally linked to reduced convection inhibition (Fig. S5), thereby leading to larger $P_{\mathrm{SOMA}}^{\prime }$.

Figure 10 Same as Fig. 9 but for frequency difference.

Table 1 Long-term mean precipitation anomaly changes across three periods with shifting PDO phases (P1: 1950–1976, PDO-; P2: 1977–1998, PDO+; P3: 1999–2021, PDO−) and the total and PDO-related moist circulation (SOMA) contribution.

Overall, when analyzing node frequency alone (Fig. 10), there are considerably more PDO-related nodes in the cold season from October to March (∼13.5 nodes per month) than in the warm season from April to September (∼6.8 nodes per month). This aligns with previous studies indicating that oceanic forcing has more influence on North American precipitation variability in the cold season than in the warm season, largely due to increased atmospheric wave train activities excited by Pacific SST anomalies (Ciancarelli et al., 2014; Ropelewski and Halpert, 1986). However, when we account for the potential PDO impact on specific circulation patterns and corresponding within each node, the seasonal contrast in ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }$ (Fig. 9) is reduced; 6.5, 9.3, and 10.3 nodes per month in the cold season and 5.7, 6.5, and 9.0 nodes per month in the warm season are considered PDO-related for the three regions.

Once the PDO-related nodes are identified, we accumulate $P_{\mathrm{SOMA}}^{\prime }$ only for days associated with significant ${\mathrm{\Delta }}_{\mathrm{PDO}}{P}_{\mathrm{SOMA}}^{\prime }$ as PDO-related circulation contributed P^′ ($P_{\mathrm{SOMA}\text{-}\mathrm{PDO}}^{\prime }$), and the result is shown in Fig. 8 (third column in each panel). We can see that although the PDO-related nodes we identified only constitute a small fraction of all nodes (about 22 %, 28 %, and 35 % for SGP, NGP, and SW, respectively), $P_{\mathrm{SOMA}\text{-}\mathrm{PDO}}^{\prime }$ shows very similar multi-decadal changes as the total $P_{\mathrm{SOMA}}^{\prime }$. Table 1 summarizes the total and PDO-related dynamic contribution to the multi-decadal changes in P^′. The dynamic contributions are highest for the SW (75 % to 82 %), with about one-third being PDO-related (32 % to 39 %). For the other two regions (SGP and NGP) with smaller dynamic contributions, $P_{\mathrm{SOMA}\text{-}\mathrm{PDO}}^{\prime }$ also explains a relatively smaller portion of the total P^′ change (mostly below 30 %).

Table 2 Same as Table 1 but for change in extreme dry days (P^′ below 10th percentile) and wet days (P^′ above 90th percentile).

In addition to the change in mean P^′, the change in P^′ extremes is also of interest because these extremes are more impactful. Figure 11 shows the probability density plot of percentile precipitation anomalies for P^′ and $P_{\mathrm{SOMA}}^{\prime }$ in the three periods. During P1 → P2, all three regions show a reduced probability of dry extremes and an increased probability of wet extremes. For the SGP and NGP, this decrease in dry extremes was more prominent in the warm season than in the cold season; the increase in wet extremes, though it also varies seasonally, shows a larger change than dry extremes and can be found in all seasons. In the second period change P2 → P3, as the PDO shifted back to the negative phase, all three regions showed an increase in dry extremes and a decrease in wet extremes, but it did not recover to the level as in P1, i.e., less dry extremes and more wet extremes than in P1. This is potentially because anthropogenic warming has increased the moisture in the atmosphere; thus the recent P3 period favored wetter conditions and wet extremes more than the P1 period with the same negative PDO phase.

Figure 11 Kernel density plot of percentile P^′ and $P_{\mathrm{SOMA}}^{\prime }$ during the three periods with shifting PDO phases (P1: 1950–1976; P2: 1977–1998; P3: 1999–2021) in three regions: (a) SGP, (b) NGP, and (c) SW. Solid and dashed lines represent percentile P^′ and $P_{\mathrm{SOMA}}^{\prime }$, respectively.

The moist circulation contribution $P_{\mathrm{SOMA}}^{\prime }$ has relatively flat curves at the drier side of the distribution (< 10th percentile; defined as dry extremes), inconsistent with the P^′ distribution, as dry extremes tend to be accompanied by stronger thermodynamic feedback; on the wetter side of the distribution (>90th percentile; defined as wet extremes), a higher consistency between P^′ and $P_{\mathrm{SOMA}}^{\prime }$ can be observed. Table 2 summarizes these observed changes in dry and wet extreme days and the corresponding total and PDO-related dynamic contribution. For the SGP and NGP, there is a ∼5 d yr⁻¹ decrease in dry extremes and a ∼10 d yr⁻¹ increase in wet extremes during the earlier PDO shift (P1 → P2, negative to positive); both the total and PDO-related circulation contribution accounts for over 90 % to the change in the SGP, while the percentages are much smaller for the NGP. However, during the recent PDO shift (P2 → P3, positive to negative), the change in both dry and wet extremes in the SGP and NGP have a considerably smaller amplitude compared to that during P1 → P2, which could be related to the anthropogenic warming-induced increase in atmospheric moisture. SW has a greater change in extremes (−24 and +17 d yr⁻¹ for dry and wet in P1 → P2, +5 and −14 d yr⁻¹ in P2 → P3) compared to the other two regions. Moist circulation overall explains a significant portion of the changes, and most of the circulation contribution can be potentially related to PDO phase shift, especially for wet extremes; a lower circulation contribution for dry extremes (both overall and PDO-related circulation) may be due to the larger influence of thermodynamic feedback in the development of a dry event.

4 Summary and discussion

SOM has been demonstrated as an effective way of identifying and visualizing large-scale circulation patterns related to various surface anomalies, such as precipitation and temperature (Liu and Weisberg, 2011). However, the SOM-based approach to quantify the contribution of anomalous circulation patterns on surface anomalies is not as reliable. On the other hand, flow analogue (including many of its variants) or the so-called dynamical adjustment technique has long been used as an important empirical approach to acquire or remove circulation anomaly variability, though it is not an effective way in terms of demonstrating what types of circulation are responsible for specific surface anomalies or quantifying the contribution from a specific group of circulation. Building upon the previous studies related to these two methods, we develop a circulation contribution quantification approach based on both the two above-mentioned methods and use it to investigate the large-scale atmospheric circulation and moisture transport pattern related to the dry and wet precipitation anomalies over the US Great Plains (SGP and NGP) and the southwest (SW), as well as the moist circulation contribution to these anomalies in terms of their daily to interannual variability and long-term multi-decadal changes.

The first part of our results is a traditional SOM clustering analysis based on the large-scale atmospheric moist circulation anomaly pattern, including mid-tropospheric geopotential height ($Z_{\mathrm{500}}^{\prime }$) and moisture transport (IVT^′). We showed that many of these SOM nodes, i.e., types of circulation, have distinctive seasonal variations in terms of their frequencies of occurrence and regional precipitation anomalies related to them. This change in the circulation–precipitation relationship is in part due to the seasonal difference in the main mechanisms producing the precipitation. For example, warm-season precipitation is mainly convectively driven, so it is more related to convective energy and atmospheric stability; therefore nodes favoring wet conditions in the cold season which feature intense moisture transport (e.g., node D1, Fig. 4) could be linked to dry or neutral precipitation anomalies due to the less desirable convective statistics related to the same circulation condition.

We then went on to develop a SOM–analogue approach and used the trained SOM to estimate the moist circulation contribution to precipitation anomalies over the SGP, NGP, and SW regions. Overall, the moist circulation condition represented by Z₅₀₀ and IVT explains the majority of the variability and amplitude of P^′, with a higher explained ratio in the cold season than in the warm season and a higher one over the northern region (NGP) than the southern regions (SGP and SW). We also showed that the main nodes responsible for dry or wet anomalies could differ largely depending on seasons and regions of interest.

The last part of the study further uses the partition results to estimate how much circulation can explain the multi-decadal shift in dry and wet mean P^′ and extremes between 1950–1976, 1977–1998, and 1999–2021 and how much the PDO-related circulation can explain them given the overall distinctive PDO phase in each of the three periods. Our result suggests circulation explains the shift in mean P^′ best in the SW (75 %–82 %), with the PDO-related circulation occupying a large portion (32 %–39 %); while circulation may account for a greater proportion of variability and amplitude in daily P^′ over the NGP compared to the other two more southern regions, it contributes less to the multi-decadal change in mean P^′. The result is somewhat similar for the change in wet extreme days between the periods, with circulation and PDO-related circulation explaining a higher portion over the SW than over the SGP and NGP; circulation contribution drops significantly for change in dry extreme days than that of wet extreme days, indicating that thermodynamic factors such as land surface feedbacks and anthropogenic warming may play a more important role in the increase in dry extremes during the past decades.

This work demonstrates the effectiveness of the proposed SOM–analogue in estimating circulation contribution to precipitation anomalies. This approach has the advantages of the C2007 SOM-based method, i.e., easy visualization of the circulation–precipitation anomaly relationship and quantification of contribution from different types of circulation, yet is less sensitive to the choice of node numbers. It also has a comparable performance as the constructed flow analogue method (Zhuang et al., 2021b) in terms of $R_{\mathrm{SOMA}}^{\mathrm{2}}$ (not shown), but $P_{\mathrm{SOMA}}^{\prime }$ calculated with SOM–analogue does not contain a thermodynamic residual like that in the analogue method. However, this is still a statistical method with certain assumptions, so limitations still exist. For example, it could become less reliable when fewer samples are presented; both the analogue and SOM–analogue approaches all assume a constant circulation–precipitation anomaly relationship which might not hold true in a changing climate. Furthermore, the identified PDO-related SOM nodes and P^′ are likely also modulated by tropical SST variability. This is due to the well-established understanding that ENSO and PDO can generate similar atmospheric and oceanic anomaly patterns (e.g., Hu and Huang, 2009). Additionally, other internal climate variability modes, such as the Atlantic Multidecadal Oscillation (AMO; e.g., Hu et al., 2011), North Atlantic Oscillation (NAO; e.g., Whan and Zwiers, 2017), and Interdecadal Pacific Oscillation (IPO; e.g., Dai, 2013), can also influence these patterns. The SOM circulation patterns defined by $Z_{\mathrm{500}}^{\prime }$ and IVT^′ simply represent a regional-scale manifestation of larger-scale variability simultaneously influenced by multiple internal climate modes. Attempting to isolate the individual contribution of these modes using statistical methods with limited data, such as SOM or SOMA, can be a complex challenge. Therefore, our results related to PDO contribution serve as a preliminary starting point which demonstrates the combined contribution of the SOM node patterns statistically linked to the PDO. To gain a more detailed quantification of the PDO's influence excluding the effects of other climate variability modes, further research integrating both observational data and climate model output is needed.