## 2 Study Area and Datasets
### Datasets  
#### 2.2.1 GNSS Vertical Displacements  
We collect daily GNSS station time series and the corresponding standard deviations (STDs) released by the Nevada Geodetic Laboratory (NGL) in the International GNSS Service 2014 reference frame from 1 January 2018 to 31 August 2022 [37], in order to be consistent with the GFO observation span. The GNSS time series were processed with GipsyX software-1.0 [38]; details of the GNSS data analysis strategies are available on the NGL website ([[http://geodesy.unr.edu/gps/ngl.acn.txt](http://geodesy.unr.edu/gps/ngl.acn.txt)]([http://geodesy.unr.edu/gps/ngl.acn.txt](http://geodesy.unr.edu/gps/ngl.acn.txt)), accessed on 18 December 2022). We use the GNSS stations selected by [PERSON] et al. [39]. Following their recommendations, we remove the GNSS sites affected by volcanic activity in California, as well as those influenced by solid Earth's porous responses due to groundwater withdrawal in California and oil extraction in the Gulf Coast, which are the two most prominent aquifer systems across the CONUS [40]. We remove the GNSS stations whose displacements are out of phase with the nearby stations [27; 41; 42] and the stations with observation durations of less than 2 years. Finally, 1437 GNSS stations are employed to study the TWS changes across the CONUS. The distribution and continuity of GNSS stations and total number of GNSS stations recorded each day are shown in Figure 2b, c. The area and number of GNSS stations of the 18 watersheds and the CONUS are listed in Table 1. Figure 3 shows the flow chart of the GNSS data processing process. \begin{table}
\begin{tabular}{c c c c} \hline \hline
**HUC-2** & **Watershed** & **Area (km2)** & \(n\) \\ \hline
11 & Arkansas–White–Red & 646,504 & 42 \\
12 & Texas–Gulf & 467,444 & 34 \\
13 & Rio Grande & 354,447 & 2 \\
14 & Upper Colorado & 307,390 & 5 \\
15 & Lower Colorado & 390,390 & 3 \\
16 & Great Basin & 402,552 & 58 \\
17 & Pacific Northwest & 779,423 & 227 \\
18 & California & 477,853 & 265 \\  & CONUS & 8,270,645 & 1315 \\ \hline \hline \end{tabular}
\end{table} 
Table 1: _Cont._  
Figure 2: Study region and basic information of GNSS stations. (**a**) Overview of the geographical area across the CONUS. The inset map in the bottom left shows the global view of the CONUS. (**b**) Distribution of GNSS stations; the color inside the circle indicates the continuity of each GNSS station. (**c**) Total number of GNSS stations recorded each day. Surface displacements induced by non-linear variations, such as non-tidal atmospheric loading (NTAL) and non-tidal oceanic loading (NTOL), are unrelated to hydrological loading and need to be eliminated from the GNSS time series [43; 44]. To isolate hydrological loading, we remove the vertical displacements due to NTAL and NTOL at each GNSS station. NTAL and NTOL can be downloaded from the Earth System Modeling group at Deutsches GeoForschungsZentrum [45]. We fit the vertical displacements caused by NTAL and NTOL with Equation (1) using the least squares method, as shown in Figure 4; the amplitudes of the annual terms in vertical displacements induced by NTAL and NTOL are 0\(\sim\)1.5 mm and 0\(\sim\)0.7 mm, respectively. \[\mathrm{y}(t)=A+Bt+Ccos(wt+\theta_{1})+Dcos(2 wt+\theta_{2})+\varepsilon \tag{1}\]  
where \(t\) represents the epoch with a daily or monthly resolution; \(y\) represents vertical displacement or TWS; \(A\) represents the offset; \(B\) represents the linear trend, \(C,D,\ \theta_{1}\), and \(\theta_{2}\) represent the amplitudes and phases of annual and semi-annual components; and \(\varepsilon\) represents a residual.