Simulation of Seawater Intrusion Area Using Feedforward Neural Network in Longkou, China

: Reliable simulation of seawater intrusion (SI) is necessary for sustainable groundwater utilization. As a powerful tool, feedforward neural network (FNN) was applied to study seawater intrusion area (SIA) ﬂuctuations in Longkou, China. In the present study, changes of groundwater level (GWL) were modeled by FNN Model 1. Then, FNN Model 2 was developed for ﬁtting the relationship between GWL and SIA. Finally, two models were integrated to simulate SIA changes in response to climatic and artiﬁcial factors. The sensitivity analysis of each impact factor was conducted by the “stepwise” method to quantify the relative importance for SIA and GWL. The results from the integrated model indicated that this method could accurately reproduce SIA ﬂuctuations when the Nash–Sutcli ﬀ e e ﬃ ciency coe ﬃ cient was 0.964, the root mean square error was 1.052 km 2 , the correlation coe ﬃ cient was 0.983, and the mean absolute error was 0.782 km 2 . The results of sensitivity analysis prove that precipitation and groundwater pumping for agriculture mainly a ﬀ ect ﬂuctuations of SIA in the study area. It can be concluded that FNN is e ﬀ ectively used for modeling SI ﬂuctuations together with GWL, which can provide enough support for the sustainable management of groundwater resources with consideration of crucial impact factors of seawater intrusion (SI).


Introduction
Groundwater is easily accessible freshwater and an essential resource for domestic, agricultural, and industrial purposes in coastal regions. With the rapid growth of human settlements and social economy, these water demands increase the pressure on aquifer, and overexploitation of groundwater resources has become an unavoidable problem in many coastal regions [1,2]. The excessive and unregular exploitation from the coastal aquifer, hydraulically connected with sea, may cause seawater intrusion (SI) [1,3], which seriously threatens the sustainability of groundwater utilization and ecosystems. SI in different areas of the world has been reported in many researches and has obviously become a global issue [2,[4][5][6][7].
In China, SI was founded as early as the 1960s, and the area was more than 2000 km 2 in 2003 [8]. Longkou, like most coastal areas, has been continuously threatened by SI due to the overexploitation of groundwater resources in the past five decades [8,9]. However, only a handful of previous studies paid attention to it. Wu et al. built a three-dimensional seawater intrusion observation network to provide insight on the characteristics of SI [10]. Xue et al. presented a three-dimensional miscible transport    In 1979, SI was first found in the Longkou. From the 1980s to the early 21st century, the amount of groundwater exploitation has been steadily increased in order to cope with the substantial growth of water resource demand. Gradually, SI developed from patches to strips, which posed a severe threat to sustainable development of the social economy and protection of the ecological environment. In order to gain insight on the SI process, the Longkou government made an effort to monitor GWL and groundwater quality, which provides much support to the present study. Figure 1c displays the distributions of GWL-monitoring wells and groundwater quality-monitoring wells.
In this study, six-monthly groundwater chloride concentration (years 2005-2017) and the data (years 2000-2017), including monthly average GWL, monthly precipitation (P), and monthly average temperature (T), were collected from the Longkou government. Besides, annual groundwater supply for agriculture (GWSA), annual groundwater supply for nonagriculture (GWSNA), and annual surface water supply for agriculture (SWSA) were obtained from the Longkou Water Resources Bulletin. It is noted that GWSNA is composed of domestic and industrial water supply. Then, the datasets of annual water supply for agriculture were processed into a monthly value based on crop structure, water demand of crop growth period, and distribution of monthly precipitation. The monthly GWSNA was calculated by the average distribution of annual value.
In addition, according to the hydrogeological data and management demand, the study area was divided into three zones of water resources management (Figure 1c), named Zone1, Zone2, and Zone3. Figure 2 illustrates the changes of monthly average GWL in the three zones. It is apparent that GWL changes of the three zones are similar in trend, and there are high fluctuations. Based on the evaluation results of the second national water resources investigation in Longkou, the annual average quantity of total recharge is about 21.49 million cubic meters, 17.71 million cubic meters, and 58.15 million cubic meters in Zone1, Zone2, and Zone3, respectively. Recharge sources of the unconfined aquifer are similar in three zones, which mainly include the infiltration of precipitation, infiltration of agricultural water, infiltration of rivers, and lateral inflow from the southern mountainous area. The proportions of four recharge sources in annual average quantity of total recharge are approximately 43%, 42%, 5%, and 10% in Zone1; 43%, 37%, 8%, and 10% in Zone2; and 37%, 40%, 17%, and 6% in Zone3. The discharging sources in Zone1 and Zone2 are primarily composed of groundwater pumping, evaporation, and natural outflow into the sea, in which the annual average quantity of groundwater pumping is about 2.318 million cubic meters, 1.914 million cubic meters, and 5.232 million cubic meters, respectively. However, it is different from Zone1 and Zone2 in that there is no natural outflow into the sea in Zone3, because the underground dam was built 1.2 km from the estuary in 1994 (Figure 1), which completely blocks the hydraulic connection between sea and groundwater system.    The SIA, where the chloride concentration is more than 250 mg/L, was acquired by inverse distance weight interpolation (IDW) of chloride concentration of monitoring wells. It is noted that the area where the chloride concentration is more than 250 mg/L in Zone3 does not need to be considered in the simulative model because of the underground dam. In other words, the GWL drawdown can not make seawater move into land in Zone3. Hence, the SIA referred to the total area of Zone1 and Zone2, which was modeled in this study. As shown in Figure 3, SIA was relieved from 134.4 km 2 in December 2007 to 115.3 km 2 in June 2011, but it showed an increasing trend from 2013 to 2017. Therefore, accurate SIA simulation is very significant for the Longkou government to prevent SI from exacerbating.

Feedforward Neural Network
ANN is distributed parallel processors imitating the structure and operation of the human brain [13]. It seems to be a "black" model or empirical model and has self-organization, self-learning, and strong nonlinear ability. The FNN, which is a typical ANN architecture, has been widely applied to kinds of hydrological research. A common FNN is composed of three separate layers, namely, input, hidden, and output layers, each consisting of processing neurons, with interconnections occurring between neurons in adjacent layers. Figure 4 shows a typical three-layer FNN architecture. The number of neurons in input layers, as well as in output layers, is determined by input variables and targets of the building model. The number of hidden layers is chosen for the best simulation result. Each processing node is interconnected across layers by special nonlinear transfer functions, which means the state of neurons of each layer only has an impact on neurons of the next layers. The distinguishing feature of FNN is the forward propagation of input signals and the backward propagation of the error. During the training stage, the connection weights and bias are adjusted repeatedly until a minimum acceptable error between simulation and observation is achieved. The mathematical expression of output is given by Equation (1).
where i, j, and k refer to the ith input node, the jth hidden node, and the kth output node, respectively. n and m are the number of input nodes and hidden nodes. w and b represent the weight and the bias.
f h and f o are the activation functions of hidden layers and output layers. X and y represent the input and output value.
strong nonlinear ability. The FNN, which is a typical ANN architecture, has been widely applied to kinds of hydrological research. A common FNN is composed of three separate layers, namely, input, hidden, and output layers, each consisting of processing neurons, with interconnections occurring between neurons in adjacent layers. Figure 4 shows a typical three-layer FNN architecture. The number of neurons in input layers, as well as in output layers, is determined by input variables and targets of the building model. The number of hidden layers is chosen for the best simulation result. Each processing node is interconnected across layers by special nonlinear transfer functions, which means the state of neurons of each layer only has an impact on neurons of the next layers. The distinguishing feature of FNN is the forward propagation of input signals and the backward propagation of the error. During the training stage, the connection weights and bias are adjusted repeatedly until a minimum acceptable error between simulation and observation is achieved. The mathematical expression of output is given by Equation (1).
where i , j , and k refer to the i th input node, the j th hidden node, and the k th output node, respectively. n and m are the number of input nodes and hidden nodes. w and b represent the weight and the bias.  A variety of components need to be considered in order to build an appropriate FNN model, including the selection of input variables, the number of hidden layers and nodes, the transfer function of hidden layers and output layers, training algorithm, and criteria of model performance.
The selection of input variables is one of the most important steps in the FNN development. It is noted that not all of the potential input variables will be equally informative because some may be uncorrelated, noisy, or have no significant relationship with the output of models [47]. Commonly, initial input variables are decided based on a priori knowledge about system fundamental properties. Some researches employed a cross-correlation method to reveal linear dependence between two variable time series. Unfortunately, most of the relationship between two variables in the natural groundwater system is nonlinear. In this study, a combination of the priori knowledge and correlation coefficient method was used for determining appropriate input variables.
The activation function is vital to the network, which makes neurons have the ability of perception. The tan-sigmoid function was applied to hidden layers and output layers and is described as Equation (2). A variety of components need to be considered in order to build an appropriate FNN model, including the selection of input variables, the number of hidden layers and nodes, the transfer function of hidden layers and output layers, training algorithm, and criteria of model performance.
The selection of input variables is one of the most important steps in the FNN development. It is noted that not all of the potential input variables will be equally informative because some may be uncorrelated, noisy, or have no significant relationship with the output of models [47]. Commonly, initial input variables are decided based on a priori knowledge about system fundamental properties. Some researches employed a cross-correlation method to reveal linear dependence between two variable time series. Unfortunately, most of the relationship between two variables in the natural groundwater system is nonlinear. In this study, a combination of the priori knowledge and correlation coefficient method was used for determining appropriate input variables.
The activation function is vital to the network, which makes neurons have the ability of perception. The tan-sigmoid function was applied to hidden layers and output layers and is described as Equation (2).
Levenberg-Marquardt algorithm (LM) was applied to train FNN models in the present study, which is a widely popular training algorithm [18]. LM algorithm is a modification of the Newton algorithm for finding an optimal solution about the minimization problem. It has great computational and memory requirements [34], and is faster and less easily trapped in local minima than other optimization algorithms.
In addition, one hidden layer FNN was used, and the number of hidden nodes was also determined by a trial-and-error procedure.
The model performance was evaluated by statistical parameters involving Nash-Sutcliffe efficiency coefficient (E), root mean squared error (RMSE), correlation coefficient (R), and mean absolute error (MAE). The four statistical indicators are defined as follows: Water 2020, 12, 2107 where h ct and h ot are simulated and observed values at time t, h ct and h ot are the average of simulated and observed output series, respectively. The closer that RMSE and MAE are to 0, and E and R are to 1, the better simulations match to observations.

Application for Seawater Intrusion Area Simulation
SI is caused by the comprehensive function of seawater system and coastal groundwater system. Unregulated and excessive exploitation of groundwater is a primary factor in many coastal regions. The GWL changes provide a direct measure of groundwater development, and SIA can display the degree of SI. Besides, critical information about aquifer dynamics is often embedded in the continuously recorded time series. In this paper, changes of SIA together with GWL were simulated based on FNN, and Figure 5 depicts the model structure, in which n represents the number of input variables. The steps of the model are as follows: 1. The FNN Model 1 is built to simulate GWL changes effected by climatic factors and human activities.
2. The FNN Model 2 is developed to simulate the relationship between GWL of three zones and SIA.
3. FNN Model 1 and FNN Model 2 are integrated, which means that the results of the FNN Model 1 are viewed as the input data of the FNN Model 2. The integrated model can help managers understand SIA changes in response to natural and artificial factors.
4. The sensitivity of each factor is studied by the "stepwise" method, in which input variables are deleted one at a time and the corresponding model is developed and validated. Then, the relative importance ranks are obtained by comparison of model performance criteria, and crucial impact factors are identified for SIA and GWL.

Determination of FNN Model 1 Architecture
As aforementioned, the determination of model architecture should be learned from the fundamental conceptual knowledge of the study area. According to the water budget description of the aquifer in Section 2.1, the major influence factors include infiltration of precipitation, infiltration of agricultural water, groundwater exploitation, and evaporation. It is noted that the influence of lateral inflow and outflow into sea for the groundwater system was not considered because the quantity is so small that they have a weak impact on GWL. The infiltration of rivers was likewise ignored in the model due to the long drying-out of rivers. Besides, the evaporation of groundwater system is extremely difficult to quantify over space and time. As an alternative, the monthly average

Determination of FNN Model 1 Architecture
As aforementioned, the determination of model architecture should be learned from the fundamental conceptual knowledge of the study area. According to the water budget description of the aquifer in Section 2.1, the major influence factors include infiltration of precipitation, infiltration of agricultural water, groundwater exploitation, and evaporation. It is noted that the influence of lateral inflow and outflow into sea for the groundwater system was not considered because the quantity is so small that they have a weak impact on GWL. The infiltration of rivers was likewise ignored in the model due to the long drying-out of rivers. Besides, the evaporation of groundwater system is extremely difficult to quantify over space and time. As an alternative, the monthly average temperature was considered as an input variable in this model. In the study area, there is a combination of surface water utilization from water storage projects and groundwater pumping for irrigation. When the total amount of groundwater and surface water for irrigation is constant, the response of groundwater system is obviously different. This is because GWSA is both a discharge and a recharge resource and SWSA is only a recharge resource. According to whether it can feed the aquifer, groundwater exploitation was divided into GWSA and GWSNA. Also, the initial monthly average GWL was considered as a fundamental input.
It should be noted that different time lags of factors for model performance were considered. Different lengths of precipitation, temperature, GWSA, GWSNA, and SWSA varying from 1 to 6 months were tested by the trial-and-error method. It was found that the delay contribution of factors, except precipitation, had no positive effect on model performance, and one time-lag of precipitation had an optimal performance. The autocorrelation coefficients for identifying appropriated lengths of GWL time-lag of three zones from 1 to 6 months are given in Table 1. One time-lag was highly significant with autocorrelation coefficients of 0.975, 0.990, and 0.981 in three zones, respectively. In summary, input variables of FNN Model 1 included monthly precipitation at time t (P(t)), monthly precipitation at time t − 1 (P(t − 1)), monthly average temperature at time t (T(t)), monthly groundwater supply for agriculture at time t (GWSA(t)), monthly groundwater supply for nonagriculture at time t (GWSNA(t)), monthly surface water supply for agriculture at time t (SWSA(t)), monthly average GWL of Zone1 at time t − 1 (GWL 1 (t − 1)), monthly average GWL of Zone2 at time t − 1 (GWL 2 (t − 1)), and monthly average GWL of Zone3 at time t − 1 (GWL 3 (t − 1)). It is noted that t represents the monthly time step of FNN Model 1.
Output variables of FNN Model 1 included monthly average GWL of Zone1 at time t (GWL 1 (t)), monthly average GWL of Zone2 at time t (GWL 2 (t)), and monthly average GWL of Zone3 at time t (GWL 3 (t)). The datasets were separated into two different parts. For the training set, 80% of the collected and preprocessed data were used, and 20% were used for the validation set. The activation function of hidden layers and output layers, the number of hidden layers, and the training algorithm in FNN Model 1 are described in Section 2.2. The number of nodes in the hidden layer was determined by the trial-and-error method. The trial-and-error procedure started with five hidden neurons initially, and the number of hidden neurons was increased to twenty. The optimum number of hidden neurons was found to be ten.  Table 2. It is apparent from this table that during the validation period, E value varied from 0.908 to 0.997, RMSE value varied from 0.133 m to 0.251 m, R value varied from 0.965 to 0.998, and MAE value varied from 0.100 m to 0.184 m. Hence, FNN Model 1 has the ability to model GWL fluctuations with reasonable accuracy, which can help managers obtain more information on the groundwater system to enhance groundwater resource management.    The names of abbreviations are described in Table A1.

Determination of FNN Model 2 Architecture
Seawater intrusion means that the balanced interface between groundwater system and seawater system moves to land. In the study area, the correlation coefficient between six-monthly average sea level and six-monthly SIA was about −0.179, which indicated that variations of sea level have a little influence on SIA fluctuations during the simulated period (years 2005-2017). Therefore, this section pays attention to the relationship between GWL and SIA. In FNN Model 2, the six-monthly average GWL of three zones and initial SIA were considered in the input layer. Correlation coefficients between GWL of each zone and SIA are presented in Table 3. Further, different lengths of GWL varying from 0 to 6 time steps were tested by the trial-and-error method. Altogether, it was found that time lag > 1 of GWL in each zone does not provide a significant contribution to the model performance. In short, input variables of FNN Model 2 included GWL 1 (k), GWL 1 (k − 1), GWL 2 (k), GWL 2 (k − 1), GWL 3 (k), and GWL 3 (k − 1). It should be noted that k represents the six-monthly time step of FNN Model 2. The output variable only included SIA at time k. In addition, datasets were divided into two different parts-80% of the collected data were viewed as the training set; the other 20% for the validation set. In FNN Model 2, the activation function of hidden layers and output layers, the number of hidden layers, and the training algorithm are described in Section 2.2, as well as FNN Model 1. The number of hidden neurons in the network was identified by various trials. The trial-and-error procedure started with two hidden neurons initially, and the number of hidden neurons was increased to twenty. It was found that thirteen neurons were optimal for model performance.  Table 4. It can be seen from the table that during the training, the statistical indicators E, RMSE, R, and MAE were 0.978, 0.877 km 2 , 0.990, and 0.575 km 2 , respectively, and the corresponding parameters were 0.821, 0.493 km 2 , 0.941, and 0.295 km 2 during the validation. The values of statistical indicators implied that the performance of FNN Model 2 was satisfactory during both the training and the validation period. Overall, the FNN Model 2 is capable of an understanding relationship between SIA and GWL. Meanwhile, the results prove that SIA and GWL are closely related, which is in agreement with the previous research [9].

Integrated Model
In this section, FNN Model 1 and FNN Model 2 are integrated, which can help managers study the SIA fluctuations in response to water exploitation activities and climate factors. The results of monthly average GWL (years 2005-2017) from FNN Model 1 were processed into six-monthly average GWL, which was considered as input data for FNN Model 2. Figure 8 displays a comparison between simulation and observation, and statistical indicators of model performance for the computation period are shown, in which E scored for 0.964, RMSE scored for 1.052 km 2 , R scored for 0.983, and MAE scored for 0.782 km 2 . Therefore, the integrated model had an excellent performance based on statistical indicators and graphs, which can be used for accurately predicting SI development.

Sensitivity Analysis for GWL and SIA Simulation
After the FNN models for GWL and SIA simulation were developed, one of the significant objectives in this study was to gain a better understanding of natural and anthropogenic factors for the coastal groundwater system. In other words, the relative importance of each factor for fluctuations of GWL and SIA was analyzed. In this section, the "stepwise" method, in which one input variable was deleted and the corresponding model was developed and validated, is outlined, and statistical indicators of model performance are compared. Based on the input variables set of  The names of abbreviations are described in Table A1.

Integrated Model
In this section, FNN Model 1 and FNN Model 2 are integrated, which can help managers study the SIA fluctuations in response to water exploitation activities and climate factors. The results of monthly average GWL (years 2005-2017) from FNN Model 1 were processed into six-monthly average GWL, which was considered as input data for FNN Model 2. Figure 8 displays a comparison between simulation and observation, and statistical indicators of model performance for the computation period are shown, in which E scored for 0.964, RMSE scored for 1.052 km 2 , R scored for 0.983, and MAE scored for 0.782 km 2 . Therefore, the integrated model had an excellent performance based on statistical indicators and graphs, which can be used for accurately predicting SI development.

Integrated Model
In this section, FNN Model 1 and FNN Model 2 are integrated, which can help managers study the SIA fluctuations in response to water exploitation activities and climate factors. The results of monthly average GWL (years 2005-2017) from FNN Model 1 were processed into six-monthly average GWL, which was considered as input data for FNN Model 2. Figure 8 displays a comparison between simulation and observation, and statistical indicators of model performance for the computation period are shown, in which E scored for 0.964, RMSE scored for 1.052 km 2 , R scored for 0.983, and MAE scored for 0.782 km 2 . Therefore, the integrated model had an excellent performance based on statistical indicators and graphs, which can be used for accurately predicting SI development.

Sensitivity Analysis for GWL and SIA Simulation
After the FNN models for GWL and SIA simulation were developed, one of the significant objectives in this study was to gain a better understanding of natural and anthropogenic factors for the coastal groundwater system. In other words, the relative importance of each factor for fluctuations of GWL and SIA was analyzed. In this section, the "stepwise" method, in which one

Sensitivity Analysis for GWL and SIA Simulation
After the FNN models for GWL and SIA simulation were developed, one of the significant objectives in this study was to gain a better understanding of natural and anthropogenic factors for the coastal groundwater system. In other words, the relative importance of each factor for fluctuations of GWL and SIA was analyzed. In this section, the "stepwise" method, in which one input variable was deleted and the corresponding model was developed and validated, is outlined, and statistical indicators of model performance are compared. Based on the input variables set of FNN Model 1 in Section 3.1, six input scenarios for sensitivity analysis were established, which are presented in Table 5. The absolute relative bias of simulation was displayed by a box plot to make a qualitative analysis. The relative importance rank of each factor was obtained based on the ratio, which is described by Equation (7). A ratio value >1 represents that the elimination of the factor reduces the simulative accuracy. Hence, the larger the ratio, the stronger the influence of the missing factor on the target.

ratio =
MAE without one f actor in input layer MAE with all variables (7) Table 5. Six input scenarios sets for sensitivity analysis.

Sensitivity Analysis for GWL Simulation
The validation period of FNN Model 1 (July 2014-December 2017) was used for sensitivity analysis. Figure 9 illustrates the absolute relative bias distributions of simulation with different scenarios. The relative importance ranks of factors are presented in Table 6. The results of GWL sensitivity analysis showed that natural and anthropogenic variables, except the temperature in Zone1, can improve the capability of FNN Model 1.  No-P(t) P(t) 0.218 1.185 No-T T 0.179 0.973 6 4 No-GWSA GWSA 0.351 1.908 1 5 No-GWSNA GWSNA 0.228 1.239 3  No-SWSA SWSA 0.140 1.373 5 The names of abbreviations are described in Table A1. Figure 9 shows that the median and mean value of absolute relative bias with No-GWSA in all zones were larger than those of the other input scenarios, which demonstrated that GWSA had the greatest impact on the accuracy of GWL simulation among all factors. The same conclusion can be held with the relative importance rank of GWSA. According to the absolute relative bias distributions in Figure 9, the precipitation has more influence on GWL changes than T, GWSNA, and SWSA. However, the relative importance of different time lag of precipitation was not equal in different zones. In Zone1 and Zone2, values of MAE were 0.313 m and 0.221 m, respectively, under the No-P(t − 1) scenario, and were greater than those of the No-P(t) scenario, where values of MAE were 0.218 m and 0.184 m. It was indicated that P(t − 1) is more important than P(t) for GWL simulation in Zone1 and Zone2, but the result was converse in Zone3. It seemed that groundwater depth might influence the recharge rate of precipitation, and the groundwater depth of Zone1 was generally larger than that of the other zones based on the analysis of monitoring data. In addition, it was not surprising that the temperature had minimal influence on model performance, ranked 6th in all zones, and the accuracy of simulation in Zone1 with No-T was even better than that of the scenario of all variables. This is because the groundwater depth is so high that groundwater evaporation is very small, especially in Zone1.
Further structures of water resource utilization can be demonstrated based on results of the sensitivity analysis of artificial factors. As stated by many authors, in semi-arid regions, as surface water resources are scarce, groundwater is heavily pumped to irrigation, which contributes to the rapid deterioration of the groundwater quality [48][49][50]. In the study, this phenomenon has been confirmed by data analysis in Longkou. It is clearly shown from the absolute relative bias distributions of boxplot ( Figure 9) that the contribution of GWSNA for GWL changes was weaker than that of GWSA. It revealed that GWSA constitutes a high proportion in total groundwater consumption, which can be proved by the general survey of groundwater pumping wells in 2011. As expected, the relative importance ranks of SWSA in Zone1, Zone2, and Zone3 were all fifth, which indicated that the contribution of SWSA for GWL fluctuations was the weakest among all artificial factors. The same result is obviously obtained from Figure 9. In Zone1, there are hardly projects of surface water irrigation, and the factor of SWSA has almost no effect on the groundwater system of Zone1. Although the irrigation region of the Wangwu reservoir is located at Zone3, the large reservoir mainly guaranteed domestic and industrial water supply during the simulative period according to the operational data. The SWSA in the study area is mainly composed of the Chijiagou reservoir and the Beixingjia reservoir in Zone2, which provide some surface water for wheat irrigation from March to May and pumping from rivers in flood season. Hence, the quantity of SWSA is very small.

Sensitivity Analysis for SIA Simulation
Based on FNN models with six scenarios, as outlined in Section 3.3.1, sensitivity of natural and anthropogenic factors for SIA was analyzed by the integrated model. The absolute relative bias distributions of SIA simulation are presented in Figure 10, and the relative importance rank of each factor is obtained in Table 7. As expected, factors of P(t), P(t − 1), and GWSA had more contributions to SIA fluctuations, which ranked third, second, and first, respectively. The same conclusion can be held from the qualitative analysis of Figure 10, and GWSA was the closest to SIA among these factors. In addition, the relative importance of T, GWSNA, and SWSA ranked fifth, fourth, and sixth. This is because groundwater depth is so high that the influence of groundwater evaporation can be ignored. Meanwhile, contributions of GWSNA and SWSA for SIA are affected by groundwater pumping distributions. The irrigation region of surface water at piedmont is far from the coastal line, and the amount of surface water for irrigation is also small. Wu et al. presented that the specific distribution of SI has a relation with the strong pumping center [9]. According to the general survey of groundwater pumping wells in 2011, most of domestic and industrial pumping wells are mainly located in central and south parts of the study area. Therefore, it was not surprising that T, GWSNA, and SWSA were weakly correlated to SIA. The results are consistent with the previous research [51]. It is apparent that strengthening groundwater management for agriculture is significant for solving the SI problem in Longkou.
Water 2020, 12, x FOR PEER REVIEW 15 of 18 Meanwhile, contributions of GWSNA and SWSA for SIA are affected by groundwater pumping distributions. The irrigation region of surface water at piedmont is far from the coastal line, and the amount of surface water for irrigation is also small. Wu et al. presented that the specific distribution of SI has a relation with the strong pumping center [9]. According to the general survey of groundwater pumping wells in 2011, most of domestic and industrial pumping wells are mainly located in central and south parts of the study area. Therefore, it was not surprising that T, GWSNA, and SWSA were weakly correlated to SIA. The results are consistent with the previous research [51]. It is apparent that strengthening groundwater management for agriculture is significant for solving the SI problem in Longkou.

Conclusions
In this study, FNN was applied to study SIA and GWL of a shallow coastal aquifer in Longkou, Shandong Province, China. The results from the model indicated that this method could accurately reproduce fluctuations of SIA and GWL. It is concluded that FNN is an excellent choice for estimating  The names of abbreviations are described in Table A1.

Conclusions
In this study, FNN was applied to study SIA and GWL of a shallow coastal aquifer in Longkou, Shandong Province, China. The results from the model indicated that this method could accurately reproduce fluctuations of SIA and GWL. It is concluded that FNN is an excellent choice for estimating SI changes together with GWL in the coastal aquifer, even influenced by highly irregular anthropogenic factors. Besides, this paper demonstrates that SIA can be considered as an alternative to estimate the development of SI, which is easier and more direct to use in water management.
According to sensitivity analysis results, GWSA and precipitation have a strong impact on changes of SI and GWL. However, the relative importance of artificial factors for SIA is different from that of GWL. It depends on the distributions of pumping wells. For example, GWSNA had more influence on GWL, especially in Zone1 and Zone2, but the contribution for SIA changes is weak. Hence, when researchers apply artificial neural network models to the groundwater system, the subzone division of study areas should be considered based on the features of water users and physical properties. Meanwhile, this study can alert the government to capture crucial factors and make effective strategies in SI management.
The modeling results and analysis will help decisionmakers gain insight on the influence of human activities on the groundwater system and promote the sustainable use of groundwater resources. In the future, the present work can be extended to develop artificial neural network models to estimate SIA at a monthly time step and identify control factors in different subregions. This work will provide more reasonable support for SI prevention.
Author Contributions: D.L. conducted the research, performed analyses, and wrote the manuscript. Y.W., E.G., G.W., and H.Z. make helpful discussions and provide technical assistance. Y.X. and W.W. revised the manuscript and given some suggestions about the expression of English scientific papers. All authors have read and agreed to the published version of the manuscript.

Acknowledgments:
We much appreciate the Longkou Water-affair Authority for providing the monitoring data and background information of the study area. The authors would like to thank the editors and reviewers for constructive comments that greatly contributed to improving the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.
Appendix A Table A1. The description of abbreviations in the paper.