2.1 Selection and evaluation of thresholds
Table 1 lists the clustering results based on our initially-set thresholds. After four-iterations, a total of 42 GPS stations are included in 10 clusters, while the remaining 16 stations are considered to be different with neighboring stations, thus they are not included in any cluster. Specifically, in the first iteration, the thresholds of three criteria matrices were set as 0.5 mm, 150 km and 10 days. In this case, only one cluster of level-one similarity is obtained, including three GPS stations. The maximal differential WRMS of this cluster is 2.02 mm. These three stations acquired from the first iteration were then removed and the thresholds were enlarged to perform the second iteration. At this stage, four clusters of level-two similarity are obtained, containing 18 GPS stations. In the third and fourth iterations, three and two clusters are obtained, respectively. With increasing thresholds, the maximal differential WRMS also increased, and the annual signal consistency of GPS stations within a cluster decreased correspondingly.
Table 1. Clustering results derived from initially-set thresholds
Similarity Level
|
Thresholds (mm/km/day)
|
No. of Stations
|
No. of Clusters
|
Maximal Differential WRMS (mm)
|
1
|
0.5/150/10
|
3
|
1
|
2.02
|
2
|
1.0/300/20
|
18
|
4
|
2.82
|
3
|
1.5/450/30
|
15
|
3
|
2.97
|
4
|
2.0/600/40
|
6
|
2
|
2.93
|
Sum
|
|
42
|
10
|
2.79
|
The selection of thresholds would affect the number of clusters at all levels of similarity and also the number of GPS stations. Therefore, the reason of those 16 stations that were not being classified into clusters probably may be due to too much strict thresholds of one criterion or more. To verify whether the initially-set thresholds are reasonable, the step lengths of these three criteria thresholds were increased, respectively, and the corresponding clustering results are listed in Table 2.
Table 2. Clustering results derived from several groups of thresholds
Thresholds Step Length
|
Similarity Level
|
Thresholds (mm/km/day)
|
No. of Stations
|
No. of Clusters
|
Maximal Differential WRMS (mm)
|
Group One
(WRMS: +0.5 mm)
|
1
|
1.0/150/10
|
3
|
1
|
2.02
|
2
|
2.0/300/20
|
20
|
4
|
3.26
|
3
|
3.0/450/30
|
20
|
4
|
3.41
|
4
|
4.0/600/40
|
3
|
1
|
2.33
|
Sum
|
|
46
|
10
|
3.10
|
Group Two
(Distance: +50 km)
|
1
|
0.5/200/10
|
7
|
2
|
2.67
|
2
|
1.0/400/20
|
19
|
4
|
2.85
|
3
|
1.5/600/30
|
12
|
2
|
3.07
|
4
|
2.0/800/40
|
8
|
2
|
3.01
|
Sum
|
|
46
|
10
|
2.89
|
Group Three
(Phase: +5 days)
|
1
|
0.5/150/15
|
3
|
1
|
2.02
|
2
|
1.0/300/30
|
23
|
5
|
2.96
|
3
|
1.5/450/45
|
16
|
4
|
2.66
|
4
|
2.0/600/60
|
7
|
2
|
3.11
|
Sum
|
|
49
|
12
|
2.81
|
From group one and two in Table 2, when step length of WRMS increased by 0.5 mm and interstation distance increased by 50 km, the total number of clusters remained the same, while the number of clustered GPS stations increased by four compared with results in Table 1. Meanwhile, the maximal WRMS increased from 2.73 mm to 3.10 and 2.89 mm, respectively. For group three, if step length of annual phase difference increased by 5 days, the total number of clusters and the contained GPS stations obviously increased, with only slight change in the maximal differential WRMS, indicating that the initially-set thresholds of annual phase difference can be strict, and the results of group three in Table 2 can be appropriate. It not only guarantees a sufficient number of clusters and clustered GPS stations, but also keeps the clustering accuracy steady compared with clustering results in Table 1. Hence, we take group three as the final clustering results for later analysis.
2.2 Clustering results
Fig. 2 shows the average annual signals of GPS stations in 12 clusters from group three, namely C1-C12. In general, the annual signals in each cluster have high consistency, and most GPS stations in Southwest China have obvious annual sinusoidal signals, which further demonstrates that application of the third criterion matrix, that is, the annual phase difference, helps to identify GPS stations with similar motion features. Moreover, higher-level clusters contain more GPS stations with obvious sinusoidal signals, while lower-level clusters contain less.
Specifically, C1 is the only level-one similarity cluster, with the highest similarity level and the strongest consistency of annual signal among internal GPS stations. C2-C6 are five level-two clusters, with a lower similarity level than C1, but also show strong consistency among annual signals of their own internal stations. C7-C10 represent four level-three clusters, and the annual signals of internal stations tend to be consistent overall, although there are a few stations show inconsistency at some certain period of a year. C11 and C12 indicate two level-four clusters, with the lowest similarity level, and the annual signal consistency of their internal stations is the weakest among all 12 clusters.
Fig. 3 depicts the geographic distribution of 12 clusters. In general, clusters with high-level similarity are concentrated in the middle of the region, while those with low-level similarity are mostly distributed in the edge. Combined with Fig. 2, it can be found that most GPS stations located in southern Sichuan and Yunnan Province show obvious annual sinusoidal signals, while cluster C8 in the Sichuan Basin and C4 in the Western Sichuan Plateau do not.
After four-iteration process, 12 clusters are generated containing 49 GPS stations, while the remaining 9 stations represented as black dots in Fig. 3 are preliminarily considered as stations with local abnormal signals. These stations are mainly distributed at the marginal area. Generally, the number of adjacent stations of a marginal station is less than that of an internal station, so there is less opportunity for the marginal station to form a cluster with at least two similar surrounding stations. This is probably the reason why stations assumed as with abnormal signals are more distributed in the marginal area.
Fig. 3 also shows that the distribution patterns of all 12 clusters are different. Some clusters show an approximate linear shape, such as C1, C3 and C10, while some present a nest-like distribution, like C2, C5 and C7. These different distribution patterns are partly due to the insufficient density of GPS stations, because there are more clusters showing as nest-like shape when involved stations are dense enough. However, based on the current density of GPS stations in this region, surface deformation characteristics can also be displayed in Fig. 3 to some extent.
2.3 Annual signal analysis at GPS stations
Based on the above clustering results, we further analyze the annual motion characteristics of the GPS stations in Southwest China, as shown in Fig. 4. Among the four southwestern provinces, the annual amplitudes of three GPS stations in Guizhou Province are the smallest, while those in the remaining three provinces have much larger annual amplitudes, with average annual amplitudes of GPS stations in Sichuan, Yunnan and Chongqing as 4.6, 6.2 and 5.4 mm, respectively. Stations with the largest annual amplitudes are mainly concentrated in central and western Yunnan, which is consistent with existing research results (Sheng et al., 2014; Jiang et al., 2017).
For the annual phase, the study area has shown obvious systematic changes. From western Yunnan along the Hengduan Mountains to the north of Western Sichuan Plateau, the annual phase of GPS stations gradually decrease from April to February. However, those of GPS stations in eastern Sichuan, namely the Sichuan Basin and Chongqing, are increased from May to July. Interestingly, for several GPS stations along the Yangtze River, starting from eastern station CQWZ in Chongqing, passing upstream through stations CQCS, LUZH, and SCJU, the annual phases gradually decrease from July to April. We suppose the potential reason may be the diverse hydrological conditions along the Yangtze River.
In addition, we observe that for clusters of high similarity levels (such as C1 to C6), their internal GPS stations have similar annual amplitudes and phases, while for clusters with lower similarity levels (such as C7 and C12), the annual phase of internal stations can be within thresholds of the corresponding level, but their annual amplitudes are quite different. The remaining 9 stations are mainly distributed at the edge of the area. Some stations have obvious abnormal annual phases, such as SCGY, SCMX and YNGM marked in Fig. 4. Such stations will be mainly analyzed in section 3.1.
To fully and truly detect the characteristics of land surface deformation, it is necessary to integrate other geodetic techniques, such as GRACE and SLM. The employed GRACE Level-2 GSM data were the latest RL06 version released by the Center for Space Research (CSR) at the University of Texas (Bettadpur, 2007; Gao et al., 2019). The time span was 163 months from April 2002 to June 2017, after excluding missing months. Compared with the last version RL05, RL06 data adopts some new background gravity field models and improved processing methods. The north-south stripes have been significantly reduced compared with RL05 data (Save et al., 2016; Save, 2019), among which quadratic polynomial fitting method (Swenson and Wahr, 2006) and the Gaussian filtering with a radius of 300 km were used to reduce stripe errors. Since the original C20 terms of the SH coefficients have large uncertainty, it needs to be replaced by estimates from satellite laser ranging (SLR, Cheng et al., 2013). In addition, replacement of the first-degree SH coefficients is also necessary, and we used the replacement values provided by Swenson et al. (2008). In order to weaken the influence of high-frequency errors, all SH coefficients were truncated to the 60th order/degree. In this study, we apply the software GRAMAT developed by Feng (2019).
Several organizations around the world, such as NASA, GFZ and EOST, have established different SLM models based on geophysical observation data and Earth models, and have provided a variety of land surface loading products (Li et al., 2020). In this study, the loading grid product provided by GFZ (Dill and Dobslaw, 2013) was used to calculate the displacement series generated by non-tidal atmospheric and oceanic loading, as well as hydrological loading at CMONOC stations (Yuan et al., 2018; Wu et al., 2019). Considering that the temporal resolution of GRACE data is one month, the single-day series of GPS and SLM were also averaged to obtain monthly displacement series. Correlation coefficients and RMS reduction rates among the three types of monthly displacement series, including GPS, GRACE and SLM, are shown in Fig. 5 and Fig. 6.
As shown from Fig. 5, the surface displacement series obtained by three types of techniques have strong correlations, with average correlation coefficients as 0.62, 0.70 and 0.87 for GPS/GRACE, GPS/SLM and GRACE/SLM, respectively, at 58 GPS stations. In comparison, the correlation coefficients for GPS/GRACE and GPS/SLM are lower than that for GRACE/SLM. This is because GPS is more sensitive to surface deformation caused by local loading effects, and could continuously monitor millimeter-level changes at GPS stations. Specifically, GPS belongs to a kind of point measurement mode (Nahmani et al., 2012; Wei et al., 2015). However, GRACE uses a surface measurement mode with spatial scale of ~400 km, while the spatial resolution of SLM provided by GFZ is also close to that of GRACE. Such local changes will not be recognized by GRACE and SLM as effectively as GPS (Wang et al., 2014), thus causing the sensitivity difference of these techniques in detecting land surface displacement.
From Fig. 6, we observe that for all 58 GPS stations, after deducting GRACE and SLM series from GPS series, the RMS of GPS time series decrease at 51 and 50 stations, respectively, with average RMS reduction rates of 20.0% and 21.8%. Stations with better GRACE correction effect are distributed in western part of the region, while the distribution of seven stations with poor correction effect is dispersive (Fig. 6a). From Fig. 6b, stations with better SLM correction effect are distributed in the Hengduan Mountains, while those in the Sichuan Basin and some along the Yangtze River display poor correction effect. Comparing Fig. 6a with 6b, the RMS reduction rates of four common stations GZFG, SCJU, SCMX and YNGM are negative, while these stations all show positive values in Fig. 6c. Therefore, it is reasonable to believe that the GPS monthly series of these four stations probably contain some abnormal local signals, which could be detected by joint employment of GPS, GRACE and SLM.
To illustrate the ability of our improved clustering algorithm to detect regional anomaly signals at GPS stations, Table 3 lists the average values of correlation coefficients and RMS reduction rates at all 58 stations, 49 stations of the clustering results (C1-C12), and 26 stations of high similarity level clusters (C1-C6). Before executing clustering algorithm, the average correlation coefficients of GPS/GRACE and GPS/SLM at all 58 stations are 0.62 and 0.70, respectively. After clustering process, the corresponding values for 49 stations in 12 clusters increase to 0.65 and 0.73. In particular, for 6 clusters of similarity level one and two (26 stations in total), the corresponding values increase to 0.67 and 0.77. Similarly, before clustering, the average RMS reduction rates of GPS-GRACE and GPS-SLM series at all 58 stations are 20.0% and 21.8%, respectively. After clustering, the corresponding values for 49 stations in 12 clusters increase to 23.0% and 24.2%, while for 26 stations in cluster C1-C6, the values increased to 26.6% and 33.0%, respectively.
Table 3. Average correlation coefficients and RMS reduction rates at all 58 stations, 49 stations of all 12 clusters, and 26 stations of 6 clusters (C1-C6).
Stations
|
Correlation Coefficient
|
RMS reduction rate (%)
|
GPS/
GRACE
|
GPS/
SLM
|
GRACE/
SLM
|
GPS-GRACE
|
GPS-SLM
|
GRACE-SLM
|
All 58
|
0.62
|
0.70
|
0.87
|
20.0
|
21.8
|
35.6
|
49 in C1-C12
|
0.65
|
0.73
|
0.86
|
23.0
|
24.2
|
34.5
|
26 in C1-C6
|
0.67
|
0.77
|
0.86
|
26.6
|
33.0
|
36.1
|