## Abstract

Based on the monitoring data of groundwater level, mining data and meteorological data from 1986 to 2015 in Cangzhou City, the change characteristics of depth of groundwater level were studied. On this basis, the nonparametric statistical test method (Mann–Kendall method) was used to study the trend of depth of groundwater level. Besides In addition, the principal component analysis method and grey correlation method were used to study groundwater influence factors. Finally, the multivariable time series model was used to predict the depth of groundwater level. The results showed that in the past 30 years, influenced by human activities and meteorological changes, the groundwater flow field in 1986–2015 was generally move from southwest to northeast, the depth of groundwater level was increasing from the initial 3.26 to 4.06 m, with an annual increase rate of 0.027 m/a; the contribution rate of exploitation factor, precipitation factor, evaporation factor were 40, 20 and 40%; the fitting figure of the observed values and the predicted values were very good, with an average relative error of 7.73%. According to the prediction schemes, when the evaporation increases by 5%, and the agricultural exploitation decreases by 5%, the depth of groundwater will reach 3.39 m; when the evaporation increases by 10%, and the agricultural exploitation decreases by 10%, the depth of groundwater will reach 3.32 m. This study had important reference significance for regional groundwater treatment and rational utilization.

## HIGHLIGHTS

The variation law of groundwater level in Cangzhou is discussed.

The influencing factors of groundwater level are discussed.

This study had important reference significance for regional groundwater treatment and rational utilization.

## INTRODUCTION

As a major water resource, groundwater plays a significant role in water supply. The combined impact of climate change and intensive human activities has caused a substantial decline in groundwater levels, as well as degradation of groundwater quality, groundwater depletion and associated changes in the ecosystems (Wang *et al.* 2018). All over the world, groundwater resources are the main water supply, especially in semi-arid regions, the groundwater is often the only water source, which is vulnerable to contamination, and is prone to depletion (Vries & Simmers 2002). The availability and sustainability of groundwater resources in many basins or plains in the world have been threatened because of continuous groundwater depletion through human activities and climatic stresses, resulting in many environmental geological problems (Brekke *et al.* 2009; Gao *et al.* 2019). According to the new water resources report released by UN, countries around the world waste a lot of water and some 700 million people lack drinking water, and if the current phenomenon continues, about 40 percent of the world's population will lack water by 2030 (Liu 2015a, 2015b), and the exhaustion of groundwater resources will form a groundwater cone of depression, resulting in land subsidence and collapse. Therefore, research on groundwater resources characteristics and evolution is essential.

In recent decades, a wide array of scientific research had been conducted to explore how water resources might respond to climate change. Only recently are researcher recognizing the important role played by groundwater resources in sustaining ecosystems and meeting the demands for drinking water, agriculture and industrial activities, as well as in the mitigation of the impacts of climate change and coupled human activities (Green 2016). Since the 20th century, foreign scholars have studied the dynamic characteristics of groundwater level. Panwar & Chakrapani (2013) discussed the effects of climate change on the groundwater resources in India. Their results showed that the effects of climate change and human activities on the groundwater resources (quality and quantity) and the vulnerability of groundwater system were extremely high. Cerlini *et al.* (2017) indicated climate changes effect on groundwater resources directly, and the data from the European Centre for Medium-Range Weather Forecasts (ECMWF) were compared to local water table measurements, and the correlation between the trend of soil moisture and local water table data was pointed out. Huet *et al.* (2016) compared various approaches for assessing groundwater recharge at a regional scale in the Canadian Shield, and finally found the joint use of several methods to provide a reasonable estimate of recharge was justifiable. Gogolev (2002) used unsaturated zone models for assessing groundwater recharge, and simulated groundwater recharge for a hypothetical homogeneous profile and real profiles of the Waterloo Moraine with deep groundwater level. Healy & Cook (2002) presented a water table fluctuation method which was based on groundwater level data for accurate estimation of groundwater recharge, and other methods that use water levels (mostly based on the Darcy equation) were also described. Smerdon (2017) focused on aspects related to predicting changes to groundwater recharge conditions to study the impact of climate change on groundwater recharge. Patil *et al.* (2020) used Visual-Modflow software to simulate evaluation the groundwater resources in Pingyang Kashgar Basin, and realized model calibration through parameter estimation, so as to study the impact of climate change on the groundwater level. Moustapha *et al.* (2013) determined the trend of time series through the groundwater level and precipitation data of the alluvial aquifer in the northwest plain of Italy, and studied the hydrodynamic characteristics and spatial distribution. Hu *et al.* (2011) used a multiple regression model to analyze the characteristics of karst groundwater in Jinan based on the groundwater monitoring data, and found the results showing obvious seasonality. A wide range of studies on the impact of climate change and human activities on groundwater resources and the variation of groundwater level can be found in this paper (Dong *et al.* 2019).

More recent studies have begun to address the groundwater dynamic in China. China have carried out the national water resources evaluation since the 1980s, and have established a set of relatively mature water resources evaluation methods, and carried out the second national water resources evaluation in the early 21st century. Xu *et al.* (2015) established a Controlled AutoRegressive (CAR) model based on multivariate time series to predict the depth of groundwater in the Hutubi River basin. Guo (2014) analyzed the relationship between precipitation and shallow groundwater depth factors and used a relative analysis method to analyze the forecast of the shallow groundwater level forecasting. Cangzhou City is located on the coast of the Bohai Sea. It is an important open city in the area around the Bohai Sea. It is also one of the areas with serious shortages of water resources in China. With the development of the Beijing–Tianjin–Hebei region, human engineering activities have become increasingly prominent, the industrial structure has been constantly adjusted, and the water resources carrying capacity is obviously insufficient, resulting in a series of environmental geological problems in the area, such as groundwater pollution, land subsidence, and groundwater cone of depression.

Although some studies have improved the general understanding on how a groundwater system may respond to climate change and human activities, studies on the potential impacts have so far made little progress and are poorly understood. The main issues that are involved are: (1) which indicators of climate change and human activities should be selected to evaluate the impacts of climate change and human activities on groundwater systems; (2) how to distinguish the influence caused by climate change and human activities; (3) the unavailability of models for predicting the long-term effects of climate change and human activities on groundwater systems and long historical data.

This paper collected the climate data, depth of groundwater level and exploitation data of 211 monitoring points in Cangzhou from 1986 to 2015. On the one hand, Mann–Kendall (M-K) method and principal component analysis were used to explore the dominant factors affecting the depth of groundwater and its variation trend; on the other hand, the corresponding schemes were considered, and the CAR model was used to predict the depth of groundwater, in order to provide a reference for the evaluation, treatment, development and utilization of groundwater resources in Cangzhou City, and for the sustainable development of social economy provide reference, and provide scientific basis for the similar relevant areas.

The main innovations of this paper are as follows: (1) based on the data of multiple influencing factors and long time series, M-K method, principal component analysis and grey correlation method were used to explore the dominant factors affecting the groundwater depth in Cangzhou and its variation trend; (2) the CAR model was used to predict the groundwater depth, and two schemes were proposed under the condition of policy permission, which can effectively solve problems such as land subsidence.

## STUDY AREA

Cangzhou City is located in the southeast of Hebei Province, China, between 115°6′–117°8′E and 37°4′–38°9′N, with a total area of 14,056 km^{2} (Figure 1) and a population of about 7.74 million. The study area is a semi-arid area, mainly warm temperate continental monsoon climate, with four distinct seasons and significant changes in temperature difference between seasons. The flat terrain and numerous rivers and streams. The Beijing–Hangzhou Grand Canal runs through the study area, with many depressions and the terrain slopes from southwest to northeast. The annual average precipitation is 519.71 mm, mostly from July to September, accounting for 65.67% of the annual precipitation, and the annual average evaporation is 1,573.98 mm. The annual average temperature is 13.49 °C (Liu *et al.* 2014). The terrain in the area is flat and low-lying, slightly inclined from southwest to northeast, and the ground elevation is generally 2–15 m. According to the physiognomy morphological characters and genetic type, the study area is divided into three physiognomy type areas: alluvial plain area, alluvial and coast plain area and marine plain area.

The sedimentary strata in the study area are mainly divided into Quaternary (Q), Tertiary (R), Mesozoic (M_{Z}) and Carboniferous-Permian (C-P) from new to old. The Quaternary aquifer is divided into Holocene series (Q_{4}), Upper Pleistocene series (Q_{3}), Middle Pleistocene series (Q_{2}) and Lower Pleistocene series (Q_{1}). This paper mainly studies the shallow water (the 1st aquifer group), which is equivalent to the bottom of Q_{4} (Holocene series) and the upper part of Q_{3} (Upper Pleistocene series). The lithology of aquifer is fine sand in the west and silt in the east, and the thickness of sand layer is 0–20 m. The relationship between recharge and discharge of shallow groundwater is mainly vertical alternation. The recharge of groundwater is mainly atmospheric precipitation infiltration, river infiltration and irrigation regression. The discharge is mainly artificial extraction and evaporation.

## DATA SOURCES AND RESEARCH METHODS

### Data sources

The monitored data of groundwater logging in Cangzhou City from 1986 to 2015 year were used. There were 211 monitoring points in total, including 201 monitoring points of series logging and 10 monitoring points of long-term logging. The data come from Hebei geological environment monitoring Institute, which collected the data of depth of groundwater level, temperature, precipitation, evaporation, river irrigation and exploitation from 1986 to 2015 (Table 1).

Year . | Temperature (°C) . | Precipitation (mm) . | Evaporation (mm) . | Agricultural exploitation (10^{4} m^{3})
. | Industrial exploitation (10^{4} m^{3})
. | Domestic exploitation (10^{4} m^{3})
. | Depth of groundwater level (m) . |
---|---|---|---|---|---|---|---|

1986 | 12.64 | 397.5 | 1240 | 39,825.65 | 4425.07 | 3982.57 | 3.26 |

1987 | 13.87 | 649.8 | 1488.6 | 44,273 | 4919.22 | 4427.3 | 3.8 |

1988 | 13.13 | 567.4 | 1221.6 | 34,628 | 3847.56 | 3462.8 | 3.86 |

1989 | 13.69 | 372.2 | 1358.2 | 50,608 | 5623.11 | 5060.8 | 3.86 |

1990 | 13.04 | 736.7 | 1065.4 | 33,385.72 | 3709.52 | 3338.57 | 3.87 |

1991 | 13.4 | 629.3 | 1746.9 | 43,994.85 | 4888.32 | 4801.97 | 3.9 |

1992 | 13.59 | 440.6 | 2015.8 | 38,318.56 | 4257.62 | 4264.8 | 3.59 |

1993 | 13.2 | 491.3 | 1687.1 | 44,237.73 | 4915.3 | 4812.18 | 5.22 |

1994 | 14.03 | 623.6 | 1871 | 43,819.72 | 4868.86 | 4760.36 | 4.08 |

1995 | 13.57 | 717.1 | 1866.4 | 38,736.1 | 4304.01 | 4227.9 | 3.31 |

1996 | 13.18 | 538 | 1508.7 | 37,469 | 4163.22 | 3746.9 | 3.43 |

1997 | 14.05 | 346.5 | 1529.9 | 55,002 | 6111.33 | 5500.2 | 5 |

1998 | 14.43 | 551.8 | 1610.8 | 51,442 | 5715.78 | 5144.2 | 5.43 |

1999 | 14.61 | 192.8 | 1671 | 58,651 | 6516.78 | 5865.1 | 6.04 |

2000 | 13.38 | 503.3 | 1505.6 | 50,107 | 5567.44 | 5010.7 | 5.31 |

2001 | 14 | 425.8 | 1834.7 | 52,047.84 | 5783.09 | 5753.06 | 5.85 |

2002 | 13.42 | 329.6 | 1962 | 53,699.92 | 5966.66 | 5967.16 | 6.56 |

2003 | 13.26 | 640.8 | 1775 | 52,347.62 | 5816.4 | 5816.38 | 6.03 |

2004 | 13.06 | 539.8 | 1585.7 | 51,006.7 | 5667.41 | 5667.39 | 5.91 |

2005 | 13.1 | 293.1 | 1796.2 | 48,990.92 | 5443.44 | 5443.42 | 6.46 |

2006 | 13.42 | 539.5 | 1572.4 | 52,955 | 5883.89 | 5295.5 | 5.66 |

2007 | 13.66 | 468.3 | 1687.3 | 52,065 | 5785 | 5206.5 | 5.97 |

2008 | 13.35 | 597 | 1500.6 | 49,840 | 5537.78 | 4984 | 4.99 |

2009 | 13.31 | 676.1 | 1480.1 | 31,239 | 3471 | 3123.9 | 4.83 |

2010 | 12.66 | 534.3 | 1458.3 | 32,040 | 3560 | 3204 | 4.68 |

2011 | 12.36 | 554.9 | 1272.7 | 28,239.87 | 2935.4 | 2715.2 | 4.42 |

2012 | 13.25 | 678.2 | 1550.5 | 27,675.43 | 2384.8 | 2329.2 | 3.93 |

2013 | 13.12 | 569.9 | 1306.1 | 26,258.46 | 2295.07 | 1958.1 | 3.65 |

2014 | 14.72 | 358.3 | 1504.4 | 27,197.83 | 2548.61 | 1498.68 | 4.16 |

2015 | 14.28 | 627.9 | 1546.4 | 24,137.22 | 2088.92 | 1731.1 | 4.06 |

Year . | Temperature (°C) . | Precipitation (mm) . | Evaporation (mm) . | Agricultural exploitation (10^{4} m^{3})
. | Industrial exploitation (10^{4} m^{3})
. | Domestic exploitation (10^{4} m^{3})
. | Depth of groundwater level (m) . |
---|---|---|---|---|---|---|---|

1986 | 12.64 | 397.5 | 1240 | 39,825.65 | 4425.07 | 3982.57 | 3.26 |

1987 | 13.87 | 649.8 | 1488.6 | 44,273 | 4919.22 | 4427.3 | 3.8 |

1988 | 13.13 | 567.4 | 1221.6 | 34,628 | 3847.56 | 3462.8 | 3.86 |

1989 | 13.69 | 372.2 | 1358.2 | 50,608 | 5623.11 | 5060.8 | 3.86 |

1990 | 13.04 | 736.7 | 1065.4 | 33,385.72 | 3709.52 | 3338.57 | 3.87 |

1991 | 13.4 | 629.3 | 1746.9 | 43,994.85 | 4888.32 | 4801.97 | 3.9 |

1992 | 13.59 | 440.6 | 2015.8 | 38,318.56 | 4257.62 | 4264.8 | 3.59 |

1993 | 13.2 | 491.3 | 1687.1 | 44,237.73 | 4915.3 | 4812.18 | 5.22 |

1994 | 14.03 | 623.6 | 1871 | 43,819.72 | 4868.86 | 4760.36 | 4.08 |

1995 | 13.57 | 717.1 | 1866.4 | 38,736.1 | 4304.01 | 4227.9 | 3.31 |

1996 | 13.18 | 538 | 1508.7 | 37,469 | 4163.22 | 3746.9 | 3.43 |

1997 | 14.05 | 346.5 | 1529.9 | 55,002 | 6111.33 | 5500.2 | 5 |

1998 | 14.43 | 551.8 | 1610.8 | 51,442 | 5715.78 | 5144.2 | 5.43 |

1999 | 14.61 | 192.8 | 1671 | 58,651 | 6516.78 | 5865.1 | 6.04 |

2000 | 13.38 | 503.3 | 1505.6 | 50,107 | 5567.44 | 5010.7 | 5.31 |

2001 | 14 | 425.8 | 1834.7 | 52,047.84 | 5783.09 | 5753.06 | 5.85 |

2002 | 13.42 | 329.6 | 1962 | 53,699.92 | 5966.66 | 5967.16 | 6.56 |

2003 | 13.26 | 640.8 | 1775 | 52,347.62 | 5816.4 | 5816.38 | 6.03 |

2004 | 13.06 | 539.8 | 1585.7 | 51,006.7 | 5667.41 | 5667.39 | 5.91 |

2005 | 13.1 | 293.1 | 1796.2 | 48,990.92 | 5443.44 | 5443.42 | 6.46 |

2006 | 13.42 | 539.5 | 1572.4 | 52,955 | 5883.89 | 5295.5 | 5.66 |

2007 | 13.66 | 468.3 | 1687.3 | 52,065 | 5785 | 5206.5 | 5.97 |

2008 | 13.35 | 597 | 1500.6 | 49,840 | 5537.78 | 4984 | 4.99 |

2009 | 13.31 | 676.1 | 1480.1 | 31,239 | 3471 | 3123.9 | 4.83 |

2010 | 12.66 | 534.3 | 1458.3 | 32,040 | 3560 | 3204 | 4.68 |

2011 | 12.36 | 554.9 | 1272.7 | 28,239.87 | 2935.4 | 2715.2 | 4.42 |

2012 | 13.25 | 678.2 | 1550.5 | 27,675.43 | 2384.8 | 2329.2 | 3.93 |

2013 | 13.12 | 569.9 | 1306.1 | 26,258.46 | 2295.07 | 1958.1 | 3.65 |

2014 | 14.72 | 358.3 | 1504.4 | 27,197.83 | 2548.61 | 1498.68 | 4.16 |

2015 | 14.28 | 627.9 | 1546.4 | 24,137.22 | 2088.92 | 1731.1 | 4.06 |

### Methods

#### Mann-Kendall method

The Mann-Kendall method is a kind of nonparametric statistical test, which has the advantage of not following certain distribution and not being interfered by a few outliers. It can well reflect the trend of elements, so it is widely used in the long-term change trend analysis of environmental time series data such as temperature and precipitation.

Among them, the positive and negative of Z indicates the upward or downward trend. At a given level of confidence, if , the original hypothesis is rejected and there is a significant upward or downward trend, such as *α* = 20% (, low confidence level), *α* = 10% (, moderate confidence level), or *α* = 5% (, high confidence level) (Güçlü 2020).

#### Grey correlation analysis method

Grey relation analysis method is the most mature part of grey theory application. It is a method to determine the degree of uncertainty between factors by finding the development law of things or things themselves through data. It is also a method to measure the degree of correlation between factors. Grey correlation analysis is not affected by the sample size and its distribution law, and has the characteristics of small calculation and convenient application. In this paper, the correlation between the influence factors and depth of groundwater level is analyzed by calculating *γ* (Yang *et al.* 2017; Lv & Wu 2020). The main steps are as follows:

Step 1: Determine the comparison sequence and reference sequence between systems.

- Step 2: standardize the data. Then, use the formula to calculate the correlation coefficient :where the is the correlation coefficient between the k-th index of the i-th factor and the reference index;
*ρ*is the resolution, generally 0.5;and are expressed as the two-stage minimum difference and the two-stage maximum difference in i and k, respectively. Step 4: According to the relational order of observation objects, it is sorted to obtain the analysis result, which reflect the advantages and disadvantages of each comparison sequence relative to the reference sequence.

#### Principal component analysis method

Principal component regression analysis method is a mathematical method to reduce the dimension of data. By synthesizing the highly correlated variable information into the principal components with low correlation, the principal components replace the original variable to participate in the regression (Johnson & Wicken 2007; Zhang *et al.* 2009; Liu 2015a, 2015b; Li *et al.* 2018). The main steps are as follows:

- (1)
Standardize the data.

- (2)
Calculate the correlation coefficient matrix. Calculate the covariance matrix of k independent variables, namely the correlation coefficient matrix V.

- (3)
Find the eigenvalues and eigenvectors of correlation matrix R. According to the characteristic equation |V-

*λ*_{i}| = 0, the eigenvalues*λ*_{1},*λ*_{2}…*λ*_{k}, and arranged in order from large to small. At the same time, the eigenvectors*α*_{1},*α*_{2}…*α*_{k}. - (4)
Calculate the contribution rate and determine the principal components. The principal components are constructed by the eigenvector corresponding to the eigenvalue, and the principal components Z

_{i}=*α*_{1}×_{1}+*α*_{2}×_{2}+ …+*α*_{k}x_{k}. The contribution rate of the original index is e_{i}=*λ*_{i}/∑*λ*_{i}, and the cumulative contribution rate is E_{m}= ∑e_{i}reflecting the total contribution rate of the first m principal components. If the cumulative contribution rate of the first p principal components are more than 70%, it shows that the first p principal components basically include all the information of meta indicators. - (5)
Explain the principal components. After determining the main components, explain them with professional knowledge.

- (6)
Principal components were used to replace the original variables for multiple regression. Each principal component was used as a new independent variable for regression analysis of the dependent variable.

#### Construction of multivariate time series model

*et al.*2017). It is assumed that the n-order CAR model is constructed by the time series of m variables in the form ofwhere the are the coefficient, among the m, n are positive integers, , are a time series variable, t is a time series, t > 1.

Modeling steps

- (1)
Make a model factor selection. The choice of model factor is the basis to ensure the prediction accuracy of the model, and the factor of high correlation is selected to study, so that the prediction result is more meaningful.

- (2)
Establish a model. In this paper, the multivariate autoregressive model in data processing system (DPS) is used to obtain the coefficient CAR (n) model for each factor in the study area after parameter estimation and determination of the highest order n. Every time to establish the model, F-test value is used to judge whether the model is ideal.

- (3)
The determination of the real order of the model and its time delay. According to the above method, although a suitable CAR (n) model can be obtained, because the parameters of some variables in the model may be 0, it is necessary to carry out F-test on some coefficients of the model, and then obtain the parameter estimation of the model. If the F-test result is significant, the original CAR (n) model is true; otherwise, the new model is true.

After the above three steps, the experiment will obtain the factors that have a great influence on the system, and establish the CAR (n) model, which will be used to predict Cangzhou City.

## RESULTS AND DISCUSSION

### Dynamic characteristics of groundwater level

#### Spatiotemporal variation of groundwater level

The flow direction of the groundwater flow field in the natural state of the study area was consistent with the terrain slope, that was, it moves from southwest to northeast. In 1990, the groundwater flow field showed the characteristics of flowing from southwest to northeast (Figure 2), the average depth of groundwater was 3.73 m, and the exploitation intensity of the whole thin area around Suning was large. Therefore, the depth of groundwater level around Suning was the largest, generally 8–14 m, and the local areas were greater than 14 m.

In 1995, the characteristics of the groundwater flow field changed greatly from 1990, showing a movement from southeast to west, with an average depth of groundwater level was 4.02 m. The depths of groundwater levels in the most areas were 0–4 and 4–8 m. The depth of groundwater level in Renqiu-Qingxian area was less than 2 m, and the shallowest was 0.7 m. The depth of groundwater level in Cangxian-Suning-Nanpi southeast area was 8–12 m, and the deepest was 14.95 m. Compared with 1990, the range of water level varied. Influenced by the heavy rainfall in 1995, from 1990 to 1995, the shallow groundwater level increased by 0.29 m. After 1995, the groundwater level continued to decline.

The characteristics of groundwater flow field in 2005 were similar to that in 1990, moving from southwest to northeast. The average depth of groundwater level was 6.16 m. From 1995 to 2005, the shallow groundwater level showed a continuous downward trend. In 2005, compared with 1995, the depth of groundwater level in 0–4 m area decreased significantly, and the area of 8–12 m area increased significantly. Compared with 1995, the depth of water level decreased by 3.15 m.

In 2015, the groundwater flow field moved from the middle to the west and from the south to the northeast, with an average groundwater level was 5.55 m and a depth of groundwater was 4.04 m. In the western region, fresh water was relatively developed, the mining volume was large, and the water level was deep. In the east of Renqiu–Hejian west–Suning line, the elevation of shallow groundwater level was lower than the sea level, and the depth of groundwater level was 10–20 m. In 2015, compared with 2005, the groundwater flow direction was basically the same, the depth of groundwater 2–4, 4–6 m division areas was significantly expanded, and other zoning areas were reduced in varying degrees. According to the above analysis, it could be concluded that from 1990 to 2015, the depth of groundwater increased, then decreased, from the initial 3.26 to 4.06 m, with an annual increase rate of 0.027 m/a.

#### Dynamic characteristics of groundwater

According to the data of depth of groundwater level from 1986 to 2015 of 211 groundwater monitoring wells in the study area, the trend of depth of groundwater level in the last 30 years was analyzed by M-K method. The test result showed that the average depth of groundwater was 4.7 m, and the test statistics were S = 94, Z = 1.659. It can be seen that the overall depth of groundwater in the past 30 years showed an upward trend. When *α* = 0.1, . Because | Z | = 1.659 > 1.645, it showed that the depth of groundwater level in Cangzhou City had a significant upward trend in the past 30 years, and passed the significance test with a reliability of 90%.

### Influencing factors and analysis

#### Grey relation analysis

According to statistics, from 1986 to 2015, the total exploitation of shallow aquifer in Cangzhou City was 154.23 × 10^{4} m^{3}, of which agricultural exploitation was 127.4 × 10^{4} m^{3}, accounting for 82.6% of the total exploitation; industrial exploitation was 13.9 × 10^{4} m^{3}, accounting for 9.0% of the total exploitation; domestic exploitation was 12.9 × 10^{4} m^{3}, accounting for 8.4% of the total exploitation. It can be seen that shallow groundwater mining was mainly used for agricultural water.

The depth of groundwater level data of Cang19-1 monitoring well from 1991 to 2015 and the depth of groundwater level sequence from 1986 to 2015 in Cangzhou area were selected for analysis. Firstly, the relationship between the depth of groundwater level and the influencing factors was plotted (Figures 3–5). From Figures 3 and 4, it can be concluded that the depth of groundwater level of Cang19-1 monitoring well and the area were on the rise as a whole, but from 2010 to 2015, the depth of groundwater level fluctuated. After 2008, the amount of artificial mining decreased and the depth of groundwater decreased rapidly. In addition, according to Figure 5, the precipitation had no significant change, and the evaporation was increasing as a whole. The groundwater exploitation was also affected by the rainfall. The higher the precipitation, the less the groundwater was extracted; and vice versa. It can be seen that the reduction of artificial exploitation had a positive effect on groundwater recharge. The depth of groundwater level was positively correlated with evaporation, negatively correlated with precipitation, and results agree well with the theory.

The study was divided into three periods: Flood period (7–9 months), dry season (1–3 months) and the annual year. The depth of groundwater level of flood period (7–9 months), dry season (1–3 months), Cang19-1 monitoring well and the annual year were taken as the comparison sequence, and the temperature, precipitation, evaporation, river irrigation, agricultural exploitation, industrial exploitation and living exploitation were taken as the reference sequence for grey correlation analysis. The results are shown in Table 2.

No. . | Temperature γ_{1}
. | Precipitation γ_{2}
. | Evaporation γ_{3}
. | River irrigation γ_{4}
. | Agricultural exploitation γ_{5}
. | Industrial exploitation γ_{6}
. | Domestic exploitation γ_{7}
. |
---|---|---|---|---|---|---|---|

Flood period (7–9 months) | 0.7299 | 0.5834 | 0.6954 | 0.5881 | 0.8089 | 0.7947 | 0.7889 |

Dry season (1–3 months) | 0.8667 | 0.7461 | 0.8978 | 0.6997 | 0.8561 | 0.8382 | 0.8330 |

Cang19-1 | 0.6606 | 0.6347 | 0.6777 | 0.6416 | 0.7561 | 0.7479 | 0.7261 |

Annual year | 0.7235 | 0.6347 | 0.7533 | 0.6198 | 0.7664 | 0.7499 | 0.7492 |

No. . | Temperature γ_{1}
. | Precipitation γ_{2}
. | Evaporation γ_{3}
. | River irrigation γ_{4}
. | Agricultural exploitation γ_{5}
. | Industrial exploitation γ_{6}
. | Domestic exploitation γ_{7}
. |
---|---|---|---|---|---|---|---|

Flood period (7–9 months) | 0.7299 | 0.5834 | 0.6954 | 0.5881 | 0.8089 | 0.7947 | 0.7889 |

Dry season (1–3 months) | 0.8667 | 0.7461 | 0.8978 | 0.6997 | 0.8561 | 0.8382 | 0.8330 |

Cang19-1 | 0.6606 | 0.6347 | 0.6777 | 0.6416 | 0.7561 | 0.7479 | 0.7261 |

Annual year | 0.7235 | 0.6347 | 0.7533 | 0.6198 | 0.7664 | 0.7499 | 0.7492 |

From Table 2, it can be seen that in the flood period, the correlation coefficient between artificial exploitation and the depth of groundwater level was the largest, which was consistent with the research results of Cang 19-1. And in the dry season, the meteorological conditions were the important factor. In the annual year, the correlation between climate and mining factors was similar. Therefore, both the exploitation factors and the meteorological conditions were the important factors that affecting the depth of groundwater level.

#### Principal component analysis

After the grey correlation analysis, the average temperature (x_{1}), precipitation (x_{2}) and evaporation (x_{3}), river irrigation (x_{4}), agricultural exploitation of groundwater (x_{5}), industrial exploitation (x_{6}) and domestic exploitation (x_{7}) in the study area were taken as the influencing factors, and the depth of groundwater level (y) was taken as the index of groundwater level response to human activities and climate change to study the main factors that affect the dynamic of groundwater level.

##### Flood period (7–9 months)

Process the data of wet season to obtain the correlation coefficient matrix (Table 3).

. | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

x_{1} | 1 | −0.401 | 0.103 | −0.257 | 0.432 | 0.402 | 0.312 |

x_{2} | −0.401 | 1 | −0.035 | 0.424 | −0.476 | −0.464 | −0.440 |

x_{3} | 0.103 | −0.035 | 1 | 0.144 | 0.328 | 0.326 | 0.364 |

x_{4} | −0.257 | 0.424 | 0.144 | 1 | −0.379 | −0.371 | −0.332 |

x_{5} | 0.432 | −0.476 | 0.328 | −0.379 | 1 | 0.993 | 0.969 |

x_{6} | 0.402 | −0.464 | 0.326 | −0.371 | −0.993 | 1 | 0.979 |

x_{7} | 0.312 | −0.440 | 0.364 | −0.332 | −0.969 | 0.979 | 1 |

. | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

x_{1} | 1 | −0.401 | 0.103 | −0.257 | 0.432 | 0.402 | 0.312 |

x_{2} | −0.401 | 1 | −0.035 | 0.424 | −0.476 | −0.464 | −0.440 |

x_{3} | 0.103 | −0.035 | 1 | 0.144 | 0.328 | 0.326 | 0.364 |

x_{4} | −0.257 | 0.424 | 0.144 | 1 | −0.379 | −0.371 | −0.332 |

x_{5} | 0.432 | −0.476 | 0.328 | −0.379 | 1 | 0.993 | 0.969 |

x_{6} | 0.402 | −0.464 | 0.326 | −0.371 | −0.993 | 1 | 0.979 |

x_{7} | 0.312 | −0.440 | 0.364 | −0.332 | −0.969 | 0.979 | 1 |

No. . | Eigenvalue . | Contribution rate . | Accumulating contribution rate % . | Principal eigenvalue . | Extraction rate of variance % . | Cumulative variance extraction rate % . |
---|---|---|---|---|---|---|

1 | 3.793 | 54.190 | 54.190 | 3.793 | 54.190 | 54.190 |

2 | 1.289 | 18.421 | 72.610 | 1.289 | 18.421 | 72.610 |

3 | 0.797 | 11.381 | 83.992 | |||

4 | 0.547 | 8.193 | 92.185 | |||

5 | 0.518 | 7.404 | 99.588 | |||

6 | 0.023 | 0.329 | 99.917 | |||

7 | 0.006 | 0.083 | 100.000 |

No. . | Eigenvalue . | Contribution rate . | Accumulating contribution rate % . | Principal eigenvalue . | Extraction rate of variance % . | Cumulative variance extraction rate % . |
---|---|---|---|---|---|---|

1 | 3.793 | 54.190 | 54.190 | 3.793 | 54.190 | 54.190 |

2 | 1.289 | 18.421 | 72.610 | 1.289 | 18.421 | 72.610 |

3 | 0.797 | 11.381 | 83.992 | |||

4 | 0.547 | 8.193 | 92.185 | |||

5 | 0.518 | 7.404 | 99.588 | |||

6 | 0.023 | 0.329 | 99.917 | |||

7 | 0.006 | 0.083 | 100.000 |

Principal component . | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

1 | 0.143 | −0.166 | 0.092 | −0.131 | 0.255 | 0.253 | 0.247 |

2 | −0.203 | 0.331 | 0.579 | 0.490 | 0.101 | 0.113 | 0.164 |

Principal component . | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

1 | 0.143 | −0.166 | 0.092 | −0.131 | 0.255 | 0.253 | 0.247 |

2 | −0.203 | 0.331 | 0.579 | 0.490 | 0.101 | 0.113 | 0.164 |

_{1}was the largest, and the absolute value of evaporation and supply coefficient were the largest in Z

_{2}, so Z

_{1}was human factor, Z

_{2}was interpreted as climate factor. These two factors were independent of each other, do not exist collinearity, and can represent most of the information of the original index, so we can do further linear regression analysis. Taking Z

_{1}Z

_{2}as the independent variable and depth of groundwater level y as the dependent variable, the multiple linear regression analysis was made and the equation was obtained:

After calculation, R = 0.846,, its F-test value was 34.070, and the probability of significance was *p* = 0.000 < 0.05, which indicates that the regression was good. Therefore, the results of principal components analysis were representative. According to the influence degree of each factor on the overall average depth of groundwater level of the region, the contribution rate of mining factor was 66.67%, and climate factor was 33.33% respectively. Therefore, in the wet season, mining factor have played a leading role in groundwater dynamics.

##### Dry season (1–3 months)

Process the data of dry season to obtain the correlation coefficient matrix (Table 6).

. | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

x_{1} | 1 | −0.363 | 0.362 | −0.093 | 0.120 | 0.106 | 0.011 |

x_{2} | −0.363 | 1 | −0.403 | 0.661 | −0.198 | −0.176 | −0.103 |

x_{3} | 0.362 | −0.403 | 1 | −0.445 | 0.339 | 0.290 | 0.243 |

x_{4} | −0.093 | 0.661 | −0.445 | 1 | −0.496 | −0.487 | −0.405 |

x_{5} | 0.120 | −0.198 | 0.339 | −0.496 | 1 | 0.993 | 0.969 |

x_{6} | 0.106 | −0.176 | 0.290 | −0.487 | 0.993 | 1 | 0.979 |

x_{7} | 0.011 | −0.103 | 0.243 | −0.405 | 0.969 | 0.979 | 1 |

. | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

x_{1} | 1 | −0.363 | 0.362 | −0.093 | 0.120 | 0.106 | 0.011 |

x_{2} | −0.363 | 1 | −0.403 | 0.661 | −0.198 | −0.176 | −0.103 |

x_{3} | 0.362 | −0.403 | 1 | −0.445 | 0.339 | 0.290 | 0.243 |

x_{4} | −0.093 | 0.661 | −0.445 | 1 | −0.496 | −0.487 | −0.405 |

x_{5} | 0.120 | −0.198 | 0.339 | −0.496 | 1 | 0.993 | 0.969 |

x_{6} | 0.106 | −0.176 | 0.290 | −0.487 | 0.993 | 1 | 0.979 |

x_{7} | 0.011 | −0.103 | 0.243 | −0.405 | 0.969 | 0.979 | 1 |

No. . | Eigenvalue . | Contribution rate . | Accumulating contribution rate % . | Principal eigenvalue . | Extraction rate of variance % . | Cumulative variance extraction rate % . |
---|---|---|---|---|---|---|

1 | 3.595 | 51.361 | 51.361 | 3.595 | 51.361 | 51.361 |

2 | 1.678 | 23.974 | 75.336 | 1.678 | 23.974 | 75.336 |

3 | 0.891 | 12.735 | 88.070 | |||

4 | 0.585 | 8.359 | 96.429 | |||

5 | 0.225 | 3.213 | 99.642 | |||

6 | 0.020 | 0.292 | 99.935 | |||

7 | 0.005 | 0.065 | 100.000 |

No. . | Eigenvalue . | Contribution rate . | Accumulating contribution rate % . | Principal eigenvalue . | Extraction rate of variance % . | Cumulative variance extraction rate % . |
---|---|---|---|---|---|---|

1 | 3.595 | 51.361 | 51.361 | 3.595 | 51.361 | 51.361 |

2 | 1.678 | 23.974 | 75.336 | 1.678 | 23.974 | 75.336 |

3 | 0.891 | 12.735 | 88.070 | |||

4 | 0.585 | 8.359 | 96.429 | |||

5 | 0.225 | 3.213 | 99.642 | |||

6 | 0.020 | 0.292 | 99.935 | |||

7 | 0.005 | 0.065 | 100.000 |

Principal component . | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

1 | 0.070 | −0.131 | 0.150 | −0.198 | 0.261 | 0.258 | 0.245 |

2 | −0.35 | 0.429 | −0.304 | 0.216 | 0.190 | 0.210 | 0.266 |

Principal component . | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

1 | 0.070 | −0.131 | 0.150 | −0.198 | 0.261 | 0.258 | 0.245 |

2 | −0.35 | 0.429 | −0.304 | 0.216 | 0.190 | 0.210 | 0.266 |

_{1}was the largest, and the absolute values of evaporation and climate coefficient were the largest in Z

_{2}, so Z

_{1}was human factor, Z

_{2}was interpreted as climate factor. These two factors were independent of each other, collinearity did not exist, and can represent most of the information of the original index, so we can do further linear regression analysis. Taking Z

_{1}Z

_{2}as the independent variable and depth of groundwater level y as the dependent variable, the multiple linear regression analysis was made and the equation is obtained:

After calculation, R = 0.714, its F-test value was 14.070, and the probability of significance was *p* = 0.000 < 0.05, which indicates that the regression was good. Therefore, the results of principal components analysis were representative. According to the influence degree of each factor on the overall average depth of groundwater level of the region, the contribution rate of mining factor was 45.5%,and climate factor was 54.5% respectively. Therefore, in the dry season, climate factors have played a leading role in groundwater dynamics.

##### The annual year

Process the original data of entire study area to obtain the correlation coefficient matrix (Table 9).

. | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

x_{1} | 1 | −0.282 | 0.360 | −0.133 | 0.261 | 0.211 | 0.119 |

x_{2} | −0.282 | 1 | −0.216 | 0.791 | −0.466 | −0.439 | −0.378 |

x_{3} | 0.360 | −0.216 | 1 | 0.012 | 0.434 | 0.418 | 0.506 |

x_{4} | −0.133 | 0.791 | 0.012 | 1 | −0.379 | −0.371 | −0.331 |

x_{5} | 0.261 | −0.466 | 0.434 | −0.379 | 1 | 0.993 | 0.969 |

x_{6} | 0.211 | −0.439 | 0.418 | −0.371 | 0.993 | 1 | 0.979 |

x_{7} | 0.119 | −0.378 | 0.506 | −0.331 | 0.969 | 0.979 | 1 |

. | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

x_{1} | 1 | −0.282 | 0.360 | −0.133 | 0.261 | 0.211 | 0.119 |

x_{2} | −0.282 | 1 | −0.216 | 0.791 | −0.466 | −0.439 | −0.378 |

x_{3} | 0.360 | −0.216 | 1 | 0.012 | 0.434 | 0.418 | 0.506 |

x_{4} | −0.133 | 0.791 | 0.012 | 1 | −0.379 | −0.371 | −0.331 |

x_{5} | 0.261 | −0.466 | 0.434 | −0.379 | 1 | 0.993 | 0.969 |

x_{6} | 0.211 | −0.439 | 0.418 | −0.371 | 0.993 | 1 | 0.979 |

x_{7} | 0.119 | −0.378 | 0.506 | −0.331 | 0.969 | 0.979 | 1 |

No. . | Eigenvalue . | Contribution rate . | Accumulating contribution rate % . | Principal eigenvalue . | Extraction rate of variance % . | Cumulative variance extraction rate % . |
---|---|---|---|---|---|---|

1 | 3.817 | 54.533 | 54.533 | 3.817 | 54.533 | 54.533 |

2 | 1.367 | 19.534 | 74.067 | 1.367 | 19.534 | 74.067 |

3 | 1.084 | 15.479 | 89.546 | 1.084 | 15.479 | 89.546 |

4 | 0.542 | 7.741 | 97.288 | |||

5 | 0.178 | 2.541 | 99.829 | |||

6 | 0.008 | 0.108 | 99.937 | |||

7 | 0.004 | 0.063 | 100.000 |

No. . | Eigenvalue . | Contribution rate . | Accumulating contribution rate % . | Principal eigenvalue . | Extraction rate of variance % . | Cumulative variance extraction rate % . |
---|---|---|---|---|---|---|

1 | 3.817 | 54.533 | 54.533 | 3.817 | 54.533 | 54.533 |

2 | 1.367 | 19.534 | 74.067 | 1.367 | 19.534 | 74.067 |

3 | 1.084 | 15.479 | 89.546 | 1.084 | 15.479 | 89.546 |

4 | 0.542 | 7.741 | 97.288 | |||

5 | 0.178 | 2.541 | 99.829 | |||

6 | 0.008 | 0.108 | 99.937 | |||

7 | 0.004 | 0.063 | 100.000 |

Principal component . | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

1 | 0.094 | −0.173 | 0.143 | −0.147 | 0.249 | 0.246 | 0.241 |

2 | 0.004 | 0.471 | 0.337 | 0.564 | 0.132 | 0.142 | 0.199 |

3 | 0.781 | −0.144 | 0.404 | 0.048 | −0.159 | −0.207 | −0.244 |

Principal component . | x_{1}
. | x_{2}
. | x_{3}
. | x_{4}
. | x_{5}
. | x_{6}
. | x_{7}
. |
---|---|---|---|---|---|---|---|

1 | 0.094 | −0.173 | 0.143 | −0.147 | 0.249 | 0.246 | 0.241 |

2 | 0.004 | 0.471 | 0.337 | 0.564 | 0.132 | 0.142 | 0.199 |

3 | 0.781 | −0.144 | 0.404 | 0.048 | −0.159 | −0.207 | −0.244 |

_{1}was the largest, therefore Z

_{1}was interpreted as human factor; the absolute value of Z

_{2}the absolute value of correlation coefficient of precipitation and river irrigation were largest, so Z

_{2}was interpreted as precipitation factor; and the absolute value of the coefficient of temperature and evaporation in the principal component Z

_{3}were the largest, so Z

_{3}was interpreted as evaporation factor. These three factors were independent of each other, do not exist collinearity, and can represent most of the information of the original index, so we can do further linear regression analysis. Taking Z

_{1}Z

_{2}Z

_{3}as the independent variables and depth of groundwater level y as the dependent variable, the multiple linear regression analysis was made and the equation was obtained

After calculation, R = 0.701, R^{2} = 491, its F-test value was 8.352, and the probability of significance was *P* = 0.001 < 0.05, which indicates that the regression was good. Therefore, the results of principal components analysis were representative. According to the influence degree of each factor on the overall average depth of groundwater level of the region, the contribution rate of mining factor, precipitation factor, evaporation factor were 40, 20 and 40% respectively. It can be seen that in the past 30 years, climate factors have played a leading role in groundwater dynamics.

#### Discussion

In flood period, the correlation coefficient of mining factor was the largest, followed by climate factor, which was consistent with the result of principal components analysis: mining factor was dominant. In the dry season, the correlation coefficient of climate factors was the largest, and the principal components showed that the mining factor accounts for 45.5%, and the climate factor accounts for 54.4%, so the climate factor was dominant. In the whole region, the correlation between climate and mining factors was similar. In the principal components, mining components account for 40%, precipitation components account for 20%, and evaporation components account for 40%. Therefore, climate factors played a dominant role, which was consistent with the conclusion of dry season. Besides the exploitation was mainly concentrated in the deep groundwater (III aquifer), so the groundwater was greatly affected by climate factors. Human activities also play an important role. For example, the quantitative analysis of water resources variation for Shiyang River Basin in recent 50 years shows that exploitation factors play a leading role in water resources change, and its effect accounts for 58.48%. Climate factors also play an important role, it accelerates the decrease of water resources (Li *et al.* 2007, 2020).

### The establishment and application of the CAR model

#### Modeling

_{(α=0.05)}= 3.028. F < F

_{(α=0.05)}, explained that the model was suitable. The expression of the CAR model after excluding nonsignificant items was as follows:where the t was the time series number, t > 1.

#### Model validation

In order to verify the prediction results of the CAR model, 28 years of observation data of depth of groundwater level from 1988 to 2015 were used to calculate, and the results were shown in Table 12.

Year . | Observed value/m . | Predicted value/m . | Relative error/% . | Year . | Observed value/m . | Predicted value/m . | Relative error/% . |
---|---|---|---|---|---|---|---|

1988 | 3.86 | 3.80 | 1.62 | 2002 | 6.56 | 6.18 | 5.74 |

1989 | 3.86 | 4.39 | −13.60 | 2003 | 6.03 | 6.30 | −4.51 |

1990 | 3.87 | 3.51 | 9.36 | 2004 | 5.91 | 5.92 | −0.21 |

1991 | 3.9 | 3.77 | 3.22 | 2005 | 6.46 | 5.94 | 7.99 |

1992 | 3.59 | 4.30 | −19.81 | 2006 | 5.66 | 6.26 | −10.69 |

1993 | 5.22 | 4.42 | 15.31 | 2007 | 5.97 | 5.67 | 4.99 |

1994 | 4.08 | 5.04 | −23.43 | 2008 | 4.99 | 5.72 | −14.66 |

1995 | 3.31 | 4.19 | −26.48 | 2009 | 4.83 | 4.31 | 10.82 |

1996 | 3.43 | 3.88 | −13.03 | 2010 | 4.68 | 4.39 | 6.28 |

1997 | 5 | 4.48 | 10.32 | 2011 | 4.42 | 4.13 | 6.47 |

1998 | 5.43 | 5.12 | 5.76 | 2012 | 3.93 | 3.68 | 6.26 |

1999 | 6.04 | 6.09 | −0.77 | 2013 | 3.65 | 3.65 | −0.08 |

2000 | 5.31 | 5.87 | −10.58 | 2014 | 4.16 | 3.59 | 13.69 |

2001 | 5.85 | 5.46 | 6.68 | 2015 | 4.06 | 3.63 | 10.59 |

Year . | Observed value/m . | Predicted value/m . | Relative error/% . | Year . | Observed value/m . | Predicted value/m . | Relative error/% . |
---|---|---|---|---|---|---|---|

1988 | 3.86 | 3.80 | 1.62 | 2002 | 6.56 | 6.18 | 5.74 |

1989 | 3.86 | 4.39 | −13.60 | 2003 | 6.03 | 6.30 | −4.51 |

1990 | 3.87 | 3.51 | 9.36 | 2004 | 5.91 | 5.92 | −0.21 |

1991 | 3.9 | 3.77 | 3.22 | 2005 | 6.46 | 5.94 | 7.99 |

1992 | 3.59 | 4.30 | −19.81 | 2006 | 5.66 | 6.26 | −10.69 |

1993 | 5.22 | 4.42 | 15.31 | 2007 | 5.97 | 5.67 | 4.99 |

1994 | 4.08 | 5.04 | −23.43 | 2008 | 4.99 | 5.72 | −14.66 |

1995 | 3.31 | 4.19 | −26.48 | 2009 | 4.83 | 4.31 | 10.82 |

1996 | 3.43 | 3.88 | −13.03 | 2010 | 4.68 | 4.39 | 6.28 |

1997 | 5 | 4.48 | 10.32 | 2011 | 4.42 | 4.13 | 6.47 |

1998 | 5.43 | 5.12 | 5.76 | 2012 | 3.93 | 3.68 | 6.26 |

1999 | 6.04 | 6.09 | −0.77 | 2013 | 3.65 | 3.65 | −0.08 |

2000 | 5.31 | 5.87 | −10.58 | 2014 | 4.16 | 3.59 | 13.69 |

2001 | 5.85 | 5.46 | 6.68 | 2015 | 4.06 | 3.63 | 10.59 |

It can be seen from Table 12, the relative errors in 1992, 1994 and 1995 were around 20%, while the relative errors in other years were very small, with an average relative error of 7.73%. As can be seen from Figure 6, the fitting figure of the observed values and the predicted values were very good. The overall trend was consistent with that of the fitting value, and the deviation degree was very small. Therefore the CAR model can be used to predict the groundwater depth.

#### Model application

According to the analysis of the historical meteorological data of Cangzhou City, in the past 30 years, the average temperature had increased by 0.55 °C every 10 years. According to the future development of the region and the change of water shortage, 2015 was taken as the base year, and the following two schemes were proposed to predict the future change trend of groundwater (Table 13).

Plan . | Annual precipitation/mm . | Annual evaporation/mm . | Agricultural production/10^{4} m^{3}
. | Depth of groundwater level/m . |
---|---|---|---|---|

2015 | 627.9 | 1,546.4 | 24,137 | 4.06 |

Plan1 | 627.9 | 1,623.72 | 22,930 | 3.39 |

Plan 2 | 627.9 | 1,701.04 | 21,723 | 3.32 |

Plan . | Annual precipitation/mm . | Annual evaporation/mm . | Agricultural production/10^{4} m^{3}
. | Depth of groundwater level/m . |
---|---|---|---|---|

2015 | 627.9 | 1,546.4 | 24,137 | 4.06 |

Plan1 | 627.9 | 1,623.72 | 22,930 | 3.39 |

Plan 2 | 627.9 | 1,701.04 | 21,723 | 3.32 |

According to the meteorological data analysis of Cangzhou, the precipitation had no significant change. From Table 12, under local government exploitation control conditions, we can see that, when the annual evaporation increases by 5%, and the agricultural exploitation decreases by 5%, the depth of groundwater will reach 3.39 m; when the annual evaporation increases by 10%, and the agricultural exploitation decreases by 10%, the depth of groundwater will reach 3.32 m. According to the actual situation of Cangzhou City, in a period of time in the future, properly reducing the amount of agricultural exploitation can reduce the depth of groundwater level, which had positive significance to the local groundwater resources.

## CONCLUSIONS

The study comprehensively analyzes the impacts of climate change and human activities on the groundwater system in the Cangzhou, China. It was found that there were serious impacts of climate change and human activities on the groundwater system. From the analyses, the following three conclusions can be drawn:

- (1)
Based on the 30-year long time series data of depth of groundwater level in Cangzhou City, the M-K method was adopted to determine the depth of groundwater, which was generally increasing for many years. The results showed that the maximum correlation coefficient between agricultural exploitation and depth of groundwater level was 0.7664, followed by evaporation.

- (2)
Based on the data of temperature, precipitation, evaporation, river irrigation and exploitation from 1986 to 2015 in Cangzhou City, using the method of principal components regression analysis, the change trend of each factor and its influence on depth of groundwater level were quantitatively studied, the contribution rate of exploitation factor, precipitation factor, evaporation factor were 40, 20 and 40%. It showed that in the past 30 years, climate factors had played a leading role in the whole research. Although there was uncertainty, it had certain research significance. The results showed that the trend of temperature and global warming was the same.

- (3)
The CAR model was applied to establish a prediction model of depth of groundwater level in Cangzhou City. The fitting figure of the observed values and the predicted values were very good, which showed that the prediction accuracy of the model was high and the effect was good. When the annual evaporation increases by 5%, and the agricultural exploitation decreases by 5%, the depth of groundwater will reach 3.39 m; when the annual evaporation increases by 10%, and the agricultural exploitation decreases by 10%, the depth of groundwater will reach 3.32 m. This study was of great significance to the evaluation, treatment and rational development and utilization of water resources in Cangzhou City.

## ACKNOWLEDGEMENTS

This work was financially supported by the National Natural Science Youth Fund Project (42002251); China’ Post-doctoral Science Fund (2018M631874); Scientific Research Projects of the Higher University in Hebei (ZD2019082) Youth Foundation of Hebei Province Department (QN2017026); Natural Science funds Project in Hebei Province (D2018403040 and D2020403022); and Hebei Key Laboratory of Geological Resources and Environmental Monitoring and Protection Fund (JCYKT201901).

## DATA AVAILABILITY STATEMENT

All relevant data are included in the paper or its Supplementary Information.