## 2 Materials and Methods
### Methods
#### 2.3.2 Ecosystem Quality Assessment  
In this study, ecosystem quality (\(EQ\)) was comprehensively assessed through three components: ecosystem productivity index (\(EPI\)), ecosystem stability index (\(ESI\)), and ecosystem bearing capacity index (\(EBCI\)). The \(EPI\), \(ESI\), and \(EBCI\) were constructed based on multisource data, such as remote-sensing data and socioeconomic data from 2001 to 2019, and were later assigned weights of 0.40, 0.28, and 0.32, respectively, using the entropy value weighting method [31]. A comprehensive \(EQ\) evaluation model for the URB was established by the weighted summation method [32]. To better reflect the spatial and temporal variation characteristics of \(EQ\) caused by the hydrological process factors and external factors, \(EQ\) was normalized, and the range of values was delineated at (10,100). Its computation is in Equations (4)-(7), which are as follows:  
\[EPI_{t,k}=\left\{\begin{array}{cc}10&G_{t,k}\leq G_{min}\\ 10+(G_{t,k}-G_{min})\times a&G_{min}<G_{t,k}<G_{max}\\ 100&G_{t,k}\geq G_{max}\end{array}\right. \tag{4}\]  
where \(EPI_{t,k}\) is the dimensionless index of productive capacity in year \(t\) of image \(k\), and its larger value indicates a higher level of ecosystem productive capacity; \(G_{t,k}\) is the total GPPin year \(t\) of image \(k\); \(G_{max}\) and \(G_{min}\) are the upper and lower limits of the GPP multi-year average, respectively; and \(a\) is the stretching constant.  
\[ESI_{t,k}=\left\{\begin{array}{cc}10&cv_{t,k}\geq cv_{max}\\ 10+(cv_{max}-cv_{t,k})\times a&cv_{min}<\,cv_{t,k}<cv_{max}\\ 100&cv_{t,k}\leq cv_{min}\\ a=&(100-10)/(cv_{max}-cv_{min})\end{array}\right. \tag{5}\]  
where \(ESI_{t,k}\) is the stability index of the image element \(k\) in year \(t\)--the larger the value, the higher is the stability of the region; \(cv_{t,k}\) denotes the coefficient of variation of the annual mean value of GPP in year \(t\) of the image element \(k\); \(cv_{max}\) and \(cv_{min}\) are the upper and lower limits of the multi-year mean coefficient of variation of GPP, respectively; and \(a\) is the stretching constant.  
\[EBCI_{t,k}=\left\{\begin{array}{cc}10&EHI_{t,k}\leq EHI_{min}\\ 10+(EHI_{t,k}-EHI_{min})\times a&EHI_{min}<EHI_{t,k}<EHI_{max}\\ 100&EHI_{t,k}\geq EHI_{max}\end{array}\right. \tag{6}\]  
\[a= (100-10)/(EHI_{max}-EHI_{min})\]  
where \(EBCI_{t,k}\) is the stability index of the image element \(k\) in year \(t\)--the larger the value, the higher is the stability of the region; \(EHI_{t,k}\) denotes the coefficient of variation of the annual mean value of GPP in year \(t\) of the image element \(k\); \(EHI_{max}\) and \(EHI_{min}\) are the upper and lower limits of the multi-year mean coefficient of variation of GPP, respectively; and \(a\) is the stretching constant.  
\[EQ=\sum_{i=1}^{n}w_{i}\times y_{i} \tag{7}\]  
where \(EQ\) is the \(EQ\) index; \(y_{i}\) represents \(EPI\), \(ESI\), and \(EBCI\) indices; and \(w_{i}\) is the weight of each index.  
The Natural Break clustering method was then used to classify \(EQ\) into the following five levels: (I) 0-45 as the lowest level, (II) 45-50 as the lower level, (III) 50-55 as the medium level, (IV) 55-65 as the higher level, and (V) 65-100 as the highest level, taking into account the URB's actual situation and the degree of influence of each index on the ecosystem.