Response of regional agricultural water use to the change of climate and plantation structure in the typical agricultural region of China

The study of climate warming is becoming more and more popular, and many observed changes have been unprecedented in previous decades to millennia since the 1950s. The Intergovernmental Panel on Climate Change (IPCC) and the Food and Agriculture Organization of the United Nations (FAO) indicated that the agriculture industry is one of the most vulnerable to climate change. Agriculture is the primary consumer of worldwide water resources (Sun et al. 2013). With the severe shortage of freshwater resources, the shortage of irrigation water has gradually become the norm (Fereres & Auxiliadora 2007). Therefore, reasonable evaluation of the potential impact of climate change on agriculture, especially irrigation water, is of great significance for ensuring food and water resources security.

As a major component of regional and global hydrological cycles, crop evapotranspiration (ET_c) has important implications in the use of agricultural irrigation water. For the impact of climate change on the ET_c of specific crops, researchers have used historical meteorological data to conduct quantitative analysis on the ET_c of certain crops (Liu & Li 2005; Silva et al. 2007). These studies have made a great contribution to exploring the impact of climate change on agricultural water use, while lacked quantitative analysis on specific measures to deal with future climate change. With the progress of modern technology, many researchers have combined crop models such as AquaCrop, SWAP and CERES with traditional experimental research, to optimize irrigation systems (Salemi et al. 2011; Benabdelouahab et al. 2016) and simulate crop growth (Mkhabela & Bullock 2012; Iqbal et al. 2014; Watson et al. 2017) under the conditions of climate change. At present, most researches related to crop models simulated crop water demand or yield based on field experiments, and did not extend simulation predictions at the point scale to large regional scales. It is necessary to apply the crop models to regional large-scale and future research predictions.

As shown in Figure 1, the agricultural water use process is a significant part of the hydrological process. The evapotranspiration and irrigation in it directly participate water resources cycle process. Solar radiation, temperature and other meteorological factors all not only affect crop production but also the reference evapotranspiration (ET₀), which in turn affect the crop water requirement, irrigation water requirement and regional water resource allocation. The increase in irrigation water utilization coefficient will reduce irrigation water consumption. Meanwhile, different crops have different water requirements and growth periods, their needs for irrigation water are also significantly different. So changes in planting structure will also affect regional agricultural water consumption. Therefore, the objective of this research is to use climate data and crop model to assess the effect of climate change on agricultural water use in the future long time series, with the introduction of quantitative changes of data such as plantation structure and irrigation water use efficiency, to put forward more specific mitigation measures on the impact of climate change on agricultural water, and determine a more reasonable proportion of agricultural water use in the water resources system.

Our research consists of three main components: (1) analyzing the spatial and temporal distribution of future changes of major climatic factors (GCM data) in the Huang-Huai-Hai region; (2) calculating the crop water requirement and irrigation water requirement of winter wheat, maize and cotton under the future scenarios and (3) assessing the changes of irrigation water requirement under the condition of changing utilization coefficient of irrigation water and plantation structure in the condition of future climate change.

The Huang-Huai-Hai region is in the north inshore of China (Figure 2). It is mainly the North China Plain formed by the alluvial deposits of the Yellow and the Haihe rivers, as well as the hills in south-central Shandong and the Shandong peninsula. It is a major grain-producing area in China. The study area covers an area of 650,000 km², with the geographic coordinates ranging from 30° to 43° north latitude and 110° to 123° east longitude, covering Hebei, Henan, Shandong and Anhui provinces. The climate type of this region is temperate monsoon.

The AquaCrop model, developed by FAO, is a model designed to deal with the disadvantages in existing crop models such as complexity, lack of transparency and detailed input data. AquaCrop is a water-driven model designed for planning, agricultural management and scenario simulations (Steduto et al. 2008). In this study, we adapted and tested the ability of the FAO-developed AquaCrop model (v3.1) to simulate the yield and ET_c for winter wheat, summer maize and cotton in the Huang-Huai-Hai region.

The data in AquaCrop model includes climate, crop, soil, field and irrigation management parameters.

Climate Data: The climate data of AquaCrop model involves daily maximum (T_max) and minimum air temperatures (T_min), precipitation (P), ET₀ and average annual carbon dioxide concentration. We selected 1997–2017 as the historical study period, selected the climate data (mean, maximum and minimum air temperatures, precipitation, relative humidity, sunshine hours and wind speed) of 89 stations in Huang-Huai-Hai from China Meteorological Data Network (http://cdc.nmic.cn/home.do). The period of 2030–2089 was selected as the future study period, and the average of 30 climate models at 89 sites in the dataset was used as future climate data including T_max, T_min, radiation (R_n) and P. ET₀ was calculated according to the FAO Irrigation and Drainage Paper 56 (Allen et al. 1998). The applied climate model information was shown in Supplementary Table S1. The NWAI-WG statistical downscaling method was from Liu & Zuo (2012), data were from (http://stdown.agrivy.com/##).

Crop Parameter: Winter wheat, summer maize and cotton were chosen as the main crop. The actual statistical yields of the three crops from 1997 to 2017 were obtained from the China Rural Statistical Yearbook (NBSC 2019), which was used to revise the crop parameters required by the model. We refer to the reference manual of AquaCrop Version 6.0/6.1 (Raes et al. 2018) for the three crops and the Study of Contour Map of Water Requirement for Major Crops in China (1993) for the K_cTr,x values of the four provinces. It is shown in Supplementary Table S2. Meanwhile, the sowing and harvest dates of the crops to be input in the model are shown in Supplementary Table S2. Growing period information was from Zheng (2015).

Soil Data: The soil parameter was from Harmonized World Soil Database (HWSD; FAO 2012). We used the HWSD-Viewer to examine and find soil data for the corresponding meteorological station in this raster soil map. The soil water content parameter was from Du et al. (2011). The soil texture is shown in Supplementary Figure S1.

When the soil water content is sufficient, ET_a is equal to crop water requirement (ET_c). Therefore, the irrigation method in this research was set as sufficient irrigation. Since this paper sets sufficient water supply conditions, ET_c is used to represent crop water requirement.

In this study, in order to ensure the accuracy of the output of AquaCrop, we used the Monte Carlo method (MCM) to correct the most sensitive crop parameters based on the yield data of the three crops from 1997 to 2017. Data from 1997 to 2007 were used as the calibration of the model, while data from 2008 to 2017 were used as the validation of the model.

In order to analyze the applicability of GCM downscaling data in the study area, we used the measured data from 1997 to 2017 to do Kruskal–Wallis Rank (KW) test on the average values of 30 Global Climate Model (GCM) climate data, which was used to test for whether the distribution has the same mean variance (Stallings 1995; Hu 2019).

The Mann-Kendall trend test is a rank-based non-parametric approach, which has been widely used in hydrometeorological time series (Xu et al. 2006). This study used the MK trend test to examine the annual variation trend of meteorological factors, ET₀ and the ET_c of the major crops in the study area.

The amount of irrigation water from water sources is the sum of the net irrigation water use and the lost water. The irrigation water utilization coefficient (η) is an indicator of the utilization of irrigation water (Li & Liu 2011).

Grey prediction is a method to predict the known but incomplete information systems and also to predict some range of changes in the grey process associated with time. And it is to do a correlation analysis of the components in the system to predict the future development trend (Cao 2007). Based on the known plantation structure of the three main crops in the study area during 1979–2017, the future plantation structure of these crops was predicted by Grey prediction. The details of this method are given in Supplementary material. In this study, the total regional irrigation water requirement was estimated by using the known or the predicted plantation area multiplied by the calculated irrigation water requirement per unit area.

The statistical magnitude (the value of KW) represented the significance level: when the significance level was lower than 0.05, there would be a significant difference between the two distributions, otherwise no significance. We have done the KW test for T_max, T_min and precipitation at all the research sites, and the results show that (Supplementary Figure S2) most values of T_max and T_min for most sites had no significant difference. Therefore, we can judge that using the average value of various GCM data to predict future climate conditions has some credibility. This indicates that the GCM data we applied can better predict the future climate change and then help to predict the impact of climate change on agricultural water use on this basis. Although the significance level of precipitation in most sites was lower than 0.05, the analysis and comparison between the statistical and predicted precipitation data shows that the spatial distribution and numerical range of the two groups of the data are very close (Figure 3). In this study, the average value of 30 climate models was temporarily used to predict future precipitation changes.

Table 1 shows the validation results of the AquaCrop model. After calibration of the model parameters with the statistical data from 1997 to 2007, the yield from 2008 to 2017 was simulated by the model, and the error analysis between the simulated data and the statistical data was performed. The sites in the list were representative sites of the four provinces. By comprehensively comparing the three error parameters, it can be found that the simulation results of wheat and cotton in Anhui and Shandong provinces were excellent, which were good in Henan and Hebei provinces. The simulation results of wheat and cotton in the provinces were good. Thus, we consider that the AquaCrop model could simulate the future crop water requirement.

NRMSE≤10%, the simulation is considered as excellent; 10%≤NRMSE≤20%, the simulation is considered as good; 20%≤NRMSE≤30%, the simulation is considered as fair; NRMSE≥30%, the simulation is considered as poor.

Figure 3 shows that under the RCP 4.5 scenario, the three climatic factors T_max, T_min and P had the same spatial distribution characteristics, and their values increased from north to south, which presented the distribution characteristics of ‘high in the south and low in the north’, while R_n showed an opposite characteristic. With the evolution of time, T_max, T_min and P in the future period all showed different increasing degrees compared with that in the historical period.

As seen from Figure 4, under the RCP 8.5 scenario, the spatial distribution of the four climatic factors was basically the same as that of the RCP 4.5 scenario, and the climate scenario thus had little effect on the spatial distribution of meteorological factors.

From the value of the MK trend test in Figure 5, it can be observed that in both climate scenarios, the four climatic factors all increased to different degrees and the rising trend of T_max, T_min and P under RCP 8.5 was more significant than under RCP 4.5. Under the RCP 4.5 scenario, R_n rose in the entire study area, whereas it showed regional differences under RCP 8.5. For example, the upward trend in the south was clear, while in most parts of Hebei province, there were no such significant trends. In particular, a slight downward trend can be seen in the northernmost part of Hebei.

As observed from Supplementary Figure S3, the largest value of ET₀ was in Shandong, which showed a decreasing trend from the middle to both sides. In terms of ET_c for the three crops, cotton required more water and had a smaller range of ET_c, varying between 600 and 860 mm. Between 1997 and 2017, the variation range of ET_c-wheat was large, in Anhui province, where the temperature was relatively high, the value of ET_c-wheat was about 400 mm; while in Hebei, to the north of the study area, the ET_c-wheat reached as high as 700 mm. The value of ET_c-maize was between 370 and 850 mm, which was the largest in the southwest of Henan and north of Hebei, exceeding 700 mm, whereas it was less than 520 mm in Anhui in the 2050s. Supplementary Figure S4 shows that in both climate scenarios, the spatial distribution of ET_c for different crops was almost the same, except that it increased more significantly under RCP 8.5.

From the values of the MK tests of ET₀ and ET_c in Figure 6, it can easily be seen that the ET_c-wheat in most of the areas showed increasing trends under both scenarios, which was more significant under the RCP 8.5 scenario, the increasing trend was most pronounced in the southern and middle parts of the study area. There were a few sites showing a downward trend in the northernmost region of Hebei. For maize, the rising trend of the whole study area was quite obvious, the value of Z_c was all around 4–5, and the trend difference between the two climate scenarios was not significant. For cotton, the increasing trend of ET_c-cotton in the RCP 8.5 scenario was slower than that in the RCP 4.5 scenario, and the change of ET_c-cotton was not particularly large.

Based on the simulated ET_c, we calculated the irrigation water requirement of these three crops. Irrigation directly affects the allocation of regional water resources, and then affects the regional hydrological processes. Figure 7 shows the calculation results of irrigation water requirement (IRR) of the three main crops. IR_wheat was the largest, followed by IR_cotton and IR_maize. For the future three eras, the spatial distribution of IR_wheat was higher in the north and lower in the south. For example, Hebei province had the largest water requirement for irrigation, between 550 and 700 mm in both the 2050 and 2080s, whereas Anhui had a very small water requirement for wheat irrigation, less than 100 mm. The IR_maize were between 0 and 468 mm, which has the largest value in the northwest of Hebei province and increased year by year. The IR_cotton varied between 0 and 624 mm, which was largest in northern Hebei, then Shandong and most of Henan. In Anhui province, all the IRR values of the three crops were the least, and there were even some places that need no irrigation such as the southernmost part of Anhui.

Under the RCP 8.5 scenario (Figure 8), the spatial distribution of the values of IRR of the three crops was basically similar to that of RCP 4.5. However, compared with the IRR in the historical years, the future change degree under RCP 8.5 was different from that under RCP 4.5. Under RCP 8.5, IR_wheat increased even more significantly in the northern part of the study area. For crops like maize with less irrigation water, under RCP 4.5, the areas with irrigation water requirement between 0 and 100 mm account for the majority of Anhui province; while under RCP 8.5, the area in Anhui with less water requirement was significantly reduced. In Shandong, Hebei and a small part of Henan, different from the above two crops, under RCP 8.5, the future IR_cotton would be smaller than that under the scenario of RCP 4.5. In Anhui, the future value of IR_cotton was almost the same as that under the scenario of RCP 4.5, both of which were higher than that in the historical years. Therefore, the increasing trend of IR_cotton under RCP 8.5 was not as obvious as that under RCP 4.5.

In this study, the irrigation coefficient in 2016 was taken as the irrigation water utilization coefficient (η) in the current year. On this basis, it was assumed that the η in future years would be increased by 2, 5 and 10% to different degrees. As shown in Figure 9, under the condition of the four values of η, the total irrigation water of the three crops all increased with time in the study region. Under the RCP 8.5 scenario, more irrigation water would be needed in most of the future eras than under the RCP 4.5, except for the 2030s. Assuming that the future η will be increased, the total amount of irrigation water in the study region would be reduced, and the higher the η, the less the water used for total irrigation. For example, under the RCP 8.5 scenario, assuming that the η remained unchanged, the total irrigation water in the study area would be up to 139 Gm³ in the 2080s, whereas assuming a 10% increase in η, the total irrigation water would reach 126 Gm³ in the 2080s.

This study used the Grey prediction method to make three reasonable predictions of the future plantation areas of wheat, maize and cotton in each province, the plantation proportions of these three crops were also calculated (the predicted results are shown in Supplementary Figure S4). In the circular graph (Figure 10(c)), the plantation proportions of the three crops in the current year (2001–2010) and the predicted years 2030, 2050 and 2080 were shown from inside to outside. It can be seen that the proportion of wheat and cotton in the three crops slightly decreased, whereas the proportion of maize significantly increased.

As can be seen from Figure 10, under the current value of η, if the future plantation structure did not change, the total amount of irrigation water in the Huang-Huai-Hai region would increase year by year with the future climate change. For example, the average total amount of irrigation water each year in the 2080s was approximately 137 Gm³, which was only 117 Gm³ from 2001 to 2010; moreover, the upward trend was more obvious under the scenario of RCP 8.5. From the perspective of the development trend of plantation structure, if the plantation area of the three crops in the future developed in line with the predicted development trend, the total amount of irrigation water would change differently. For instance, during the 2030s, the total irrigation water in the Huang-Huai-Hai region was about 124 Gm³, which was 2.4 Gm³ less than the total irrigation water when the plantation structure was unchanged. While during the 2080s, the total amount of irrigation water decreased the most. Compared with the predicted total amount of irrigation water with an unchanged plantation structure, the total amount of irrigation water in both climate scenarios decreased by about 10 Gm³.

Under the two climate scenarios, ET₀ in the study area showed a clear upward trend, closely related to the change of climatic factors. ET₀ is determined by the climatic factors, which provides energy for the vaporization process and the migration of water vapor from the evaporating surface (Allen et al. 1998). To research the effects of different climatic factors on ET₀, the sensitivity of ET₀ to the changes of major climatic factors and the contribution rate of climatic factors to the changes of ET₀ were analyzed (Figure 11).

It can be observed that ET₀ was most sensitive to the changes in R_n, followed by T_max, T_min and P. The sensitivity coefficient of P in the entire study area and T_min in the northern area were negative, which indicated that the two climatic factors were negatively correlated with ET₀. Liu et al. (2010) studied the change characteristics of ET₀ and its influencing factors in the North China Plain, they found that R_n accounted for a large proportion of the factors that influenced ET₀, which is consistent with the results of this article. The reason is that the process depends on how much energy is needed to vaporize soil water, and solar radiation is the biggest source of energy, turning large amounts of liquid water into vapor (Allen et al. 1998). This study also found that ET₀ was highly sensitive to changes in temperature. Allen et al. (1998) found that a rise in temperatures would increase crop water requirements, whereas a reduction in crop yields would result when temperatures constrained crop development. Although precipitation is not directly involved in the calculation of ET₀, it is closely related to ET₀, because ET₀ is contrary to the meteorological conditions required by precipitation weather. The thicker the clouds, the shorter the sunshine hours, the higher the air humidity, and the slower the air flow, the more conducive to the formation of precipitation, but not to evapotranspiration (Liu et al. 2010).

The contribution rate of climatic factors is determined by the sensitivity of ET₀ to climatic factor changes and the range of climatic factor changes (Li et al. 2017). Therefore, the spatial distribution of the contribution value is different from that of the sensitivity. In the southern part of the study area, R_n contributes the most to ET₀. In the northern part of the study area, T_max contributes the most to ET₀, while ET₀ was most sensitive to the changes of R_n in this area. The reason could be that the interannual variation trend of the value of R_n was not pronounced, R_n thus contributed less to ET₀. Similarly, even if the region's ET₀ is not the most sensitive to T_max, T_max still contributes the most to ET₀ in the region because of T_max's prominent annual variation. Si et al. (2017) had the same results as us: in their study on the actual evapotranspiration in the Hetao region, they pointed that the net solar radiation was the most sensitive factor on ET₀, while the dominant factor which leads to the increase of ET₀ was air temperature. In our study, we can infer that the increase of R_n is the main reason for the increase of future ET₀ in the southern part of the Huang-Huai-Hai region, while the rise of temperature is the main reason for the increase of ET₀ in the northern part. When studying reference evapotranspiration sensitivity coefficients to climate factors in the Huang-Huai-Hai area, Yang et al. (2013) found that solar radiation was the most sensitive and primarily controlling variable for the variation of ET₀ in summer maize season, and higher sensitive coefficient values of ET₀ to solar radiation and temperature were detected in the east part and southwest part of 3H plain, respectively. Such regional differences may be due to differences in other climate factors such as wind speed and relative humidity in different regions. Gong et al. (2006) revealed that the large spatial variability of the sensitivity coefficients of all the climatic variables in the middle and tower regions of the basin was to a large extent determined by the distinct wind-speed patterns in different regions.

ET₀ in the study area showed a slight downward trend in the historical scenario, but a significant upward trend in the future long time series. In the study of effects of climate change on reference crop evapotranspiration, Yao (2009) pointed out that the major food-producing countries such as China and India would have a clear upward trend in ET₀ in the future period, which is consistent with our findings. In the prediction of crop water requirement ET_c, Farahani et al. (2009) evaluated the ET_c-cotton with the AquaCrop model, the deviation for ET_c was between 2.1 and 10.2%, which was within reasonable bounds. Heng et al. (2009) presented total measured and simulated ET the deviations of which were in the range of −1.23 and −8.4%. These studies confirmed that the model could acceptably simulate the value of ET_c and could be used to study water balance then further estimate water consumption. However, most of these studies were based on the field test, in our research, we tried to use AquaCrop to simulate the future crop water demand on a large regional scale and achieved good results as well. The results of the range of ET_c-wheat, ET_c-maize, ET_c-cotton calculated by the AquaCrop model are similar to those of the three crops by Luo & Liu (2020), and our simulated yield output was verified with statistical data, therefore, the AquaCrop model can be used to predict general trends in crop water requirements in the future. ET_c is largely determined by the change of ET₀ due to the narrow variation range of crop coefficients. Therefore, the temporal variation trend of ET_c of the three crops is basically consistent with that of ET₀, a significant upward trend. The upward trend of ET_c in the RCP 8.5 scenario is more pronounced than that under the RCP 4.5 scenario. Lu et al. (2019) pointed that greenhouse gas emissions are relatively low under the RCP 4.5 scenario, so that the effect of climate change will be slighter than that of RCP 8.5. In this way, meteorological factors such as temperature, and net radiation, which contribute a higher rate to ET₀ change less, and the ET_c variation trend under RCP 4.5 is not as obvious as that under RCP 8.5.

The IRR values of the three crops in the study area showed different rising trends. For example, the rising trend of IR_wheat was more significant than that of ET_c-wheat, because the increasing trend of precipitation during the growth period of wheat was less than that of ET_c-wheat. However, ET_c-maize showed a slight upward trend, and the irrigation water requirement of maize thus also showed an upward trend due to no significant change of its effective precipitation during the growth period. During the growth period of cotton, the average effective precipitation varied between 500 and 550 mm, and the increase of ET_c-cotton is not apparent. Therefore, the change of IR_cotton is the least among the three crops. When studying the cotton irrigation water requirement in North China Plain from 1965 to 2016, Yang et al. (2021) has reached a similar conclusion with us that the average ET_c of cotton is around 573 mm. However, due to the high value of effective precipitation in the growth period of cotton, only 25% of cotton's water requirement comes from irrigation on average.

The irrigation water utilization coefficient is one of the characteristics of irrigation efficiency (Feng 2013). With the improvement of science and technology, the water use efficiency in agriculture will be gradually improved in the future, and the irrigation water utilization coefficient for agriculture will also change. Therefore, it is necessary to explore the effect of improving the irrigation water utilization coefficient on mitigating its adverse effects under the conditions of climate change. The above results suggested that irrigation water requirement would decrease with the increase of η. As Deng et al. (2006) said in their research, improving irrigation efficiency is essential, because it is related not only to the use of water resources but also to agricultural production in the region.

Changes in plantation patterns will also affect the amount of agricultural water. According to the analysis of irrigation water use of the three crops above, the growth periods of maize and cotton are in the rainy season, and the irrigation water consumption is thus relatively small. Maize and cotton both belong to water-saving crops in this region. In contrast, the effective precipitation in the growth period of wheat is less, and more irrigation water is thus needed. Therefore, the expansion of maize plantation area and the decrease of wheat plantation area in the future could be the most important reason for the decrease of the total irrigation water for the three crops in this region. In the study of Ma et al. (2011), higher irrigation water requirements were found in counties with relatively high cultivated areas of wheat and vegetables while counties with relatively high cultivated areas of cotton and soybean had low irrigation requirements in North China Plain. Therefore, if the variation of plantation structure is similar to the prediction in our study, it will reduce the agricultural water use in the study region, thus alleviating the regional water shortage.

Future changes in climate, irrigation efficiency and plantation structure will affect agricultural water use to different degrees. Therefore, corresponding measures should be taken to mitigate the adverse effects according to the regional characteristics. In view of the impact of climate change on agricultural water use, some studies also put forward some targeted strategies. Sun et al. (2018) found that the Loess Plateau in Northern Shaanxi would become warmer and more humid in the future; they suggested that administrators should appropriately reduce the irrigation water and focus on the water requirement based on ecological principles. Joyce et al. (2011) discovered that increasing temperature would lead to the increase of agriculture water requirement, which would also lead to additional regional stresses on the water system. The Huang-Huai-Hai region is short of water and the precipitation varies greatly. As the main grain-producing area in China, to fully use the limited water resources, the region's precipitation should be fully used and the plantation structure should be appropriately adjusted. The combined adjustment of irrigation efficiency and plantation structure will bring greater water-saving benefits in the agricultural production. On the premise of ensuring national food security, reducing water-intensive crops such as wheat and expanding the plantation of water-saving crops such as maize and cotton should be considered (Qu et al. 2003).

The three climatic factors (T_max, T_min and P) increased in different degrees under the RCP 4.5 and RCP 8.5 scenarios and were more obvious under the RCP 8.5 scenario. The R_n with climate shows a slight increase trend with the climate change.

The ET_c of wheat had the most obvious rising trend with climate change. However, the variation tendency of IRR was not as obvious as that of ET_c due to the different trends of the effective precipitation during the growth period.

The increase of the value of η will reduce the irrigation water requirement of the region and greatly improve the water efficiency of agriculture. Under RCP 4.5, by the 2080s, the 10% increase of η in the value would reduce the irrigation water requirement of 12.0 Gm³ in the area. This water-saving effect is more obvious under RCP 8.5.

Due to the future expansion of plantation area of maize, which requires less water, the irrigation water requirement will be relatively reduced in the future. The combination of adjustment of these irrigation efficiencies and plantation structure will bring greater water-saving benefits to agricultural production.

Kaixuan Wang calculated the data and wrote the manuscript; Jiahui Wang and Yingjie Li analyzed the data; Xinyu Qi and Chong Li performed the statistical analysis; Fei Gao and Xinyu Hu revised the manuscript; Shikun Sun designed research; all authors read and reviewed the manuscript.

This work is jointly supported by the National Natural Science Foundation of China (51979230; 52122903), Fok Ying-Tong Education Foundation (171113) and Science Fund for Distinguished Young Scholars of Shaanxi Province (2021JC-20).

There is no conflict of interest.

Data cannot be made publicly available; readers should contact the corresponding author for details.