*2.2. Methods*

The linear regression method was used to analyze the inclination of climatic elements over a long time scale. To evaluate inclination, we established the unary linear regression equation *y* = *a* + *bx* between the SD (*y*) and time series (*x*, year), where *a* is the regression constant and *b* is the regression coefficient—namely, the inclination rate. Both *a* and *b* can be calculated using the least square method. *b* > 0 indicates that *y* increases with an increase in time *x*; on the contrary, *b* < 0 means *y* decreases as time *x* increases.

A moving average was introduced to investigate the SD trends in this study. For a climate element series, its moving average series can be depicted using Equation (6):

$$\overline{y}\_j = \frac{1}{k} \sum\_{i=1}^k y\_{i+j-1} \qquad (j = 1, \ 2, \ \dots, n-k+1) \tag{6}$$

where *n* stands for the number of years, *yj* is the *j*th moving average of SD, *k* is the step length and *i*th is the sequence of the step length. From the curve of the moving average, we can see whether the climate variable is increasing or decreasing, which helps determine the turning point of the climate variable.

We also analyzed the trend coefficient *ryx*, which is defined as a correlation coefficient between the climatic element series and the time series, and can be expressed as in Equation (7) [10]:

$$\sigma\_{y\chi} = \frac{\sum\_{i=1}^{n} (y\_i - \overline{y})(i - \overline{x})}{\sqrt{\sum\_{i=1}^{n} (y\_i - \overline{y})^2 \sum\_{i=1}^{n} (i - \overline{x})^2}} \tag{7}$$

where *n* denotes the number of years, *i* is the year sequence, *yi* is the element value in the *i*th year, *y* represents the mean of all the element values in *n* years and *x* is equal to (*n* + 1)/2. If *ryx* > 0, the element increases during *n* years, whereas if *ryx* < 0, the trend in the element declines during *n* years; *ryx* = 0 means no change. Furthermore, correlation and fitting analyses were used in our study.

#### **3. Results and Discussion**

#### *3.1. Trends in Sunshine Duration from 1981 to 2020*

Figure 2 shows the spatial distribution of the trend coefficients of SD for the ten stations over China from 1981 to 2020. SD trend coefficients are presented in Table 2. There was a clear decline in SD at seven stations, but not in Kunming, Guangzhou and Shenyang located in Southwest, South and Northeast China, respectively. Owing to its environment and high altitude, Kunming showed a significant increase in SD, with a trend coefficient greater than 0.5. Guangzhou, thanks to its being an inshore region in southeast China and its advanced technology and tertiary industry, had a positive trend coefficient of 0.21. The trend coefficient of Shenyang was 0.03; hence, the increase in SD in this region was minimal. In contrast, stations in the other representative cities in China showed a decreasing trend in SD. The Beijing, Shanghai and Wuhan stations had trend coefficients lower than −0.5. Notably, in Beijing and Shanghai, the trend coefficients were both less than −0.6, which were attributed to their larger population densities and levels of anthropogenic pollution due to urbanization.

**Figure 2.** Trend coefficients of SD from 1981 to 2020 for the ten selected stations.


**Table 2.** Statistical summary of SD.

Through a regression analysis of the daily average SD, a ranking of the cities with decreasing SD trends can be obtained: Lanzhou < Chengdu < Urumqi < Wuhan < Harbin < Beijing < Shanghai, ranging from −0.03 h d−<sup>1</sup> per decade to −0.36 h d−<sup>1</sup> per decade. Increasing trends in SD were found in Kunming, Guangzhou and Shenyang, with the biggest increasing trend of 0.38 h d−<sup>1</sup> per decade in Kunming. The regression equations are listed in Table 2.

Setting the step length *k* as 5 years, the moving average curves of the annual daily mean SD of the ten stations are shown in Figure 3. Except for obvious upward trends in Kunming and Guangzhou and a weak upward trend in Shenyang, downward trends in the SD of the other seven stations can be observed from Figure 3—which is consistent with the results of the linear analysis and trend analysis. Moreover, except for the turning points of Urumqi and Harbin in 2015 and Shenyang in 2006, inflection points in the SD of most stations could be identified at around 2010.

**Figure 3.** Moving average of the daily mean sunshine duration per year.

Figure 4 shows the seasonal trends in SD from 1981 to 2020. Consistently, the changes in SD in the studied cities exhibited a similar variation. The overall SD decreased in 1981– 2010 and started increasing after 2010. This fall–rise trend might be due to environmental protection measures taken by the Chinese government over the past decade. In Shanghai, Wuhan and Chengdu, SD trends seasonally varied in fall–rise and rise–fall patterns. This was primarily due to the climate characteristics of spring and winter in these areas, which are not conducive to the diffusion of pollutants, further reducing the solar radiation reaching the ground, and consequently reducing SD.

**Figure 4.** Seasonal trends in SD from 1981 to 2020 for the ten selected stations. (**a**) Seasonal trends in SD for Beijing; (**b**) Seasonal trends in SD for Shenyang; (**c**) Seasonal trends in SD for Harbin; (**d**) Seasonal trends in SD for Shanghai; (**e**) Seasonal trends in SD for Wuhan; (**f**) Seasonal trends in SD for Guangzhou; (**g**) Seasonal trends in SD for Chengdu; (**h**) Seasonal trends in SD for Kunming; (**i**) Seasonal trends in SD for Urumqi; (**j**) Seasonal trends in SD for Lanzhou.
