Source Apportionment and Health Risk Assessment of Zn (II) Ions Based on Gridded Spatial Interpolation and Geochemical Equilibrium Simulation

11 Soil heavy metal pollution had become a global issue involving environmental safety and human health risks. A 12 methodology was explored to quantify the sources of heavy metals in the soils and investigate the spatial distributions 13 of heavy metals by the gridded spatial scale. The case study was implemented in the industrial waste sites in Suzhou 14 city, Jiangsu province. Zinc (Zn) was screened out as the targeted metal (TM) through the potential ecological risk 15 assessment, the species of which was simulated by the geochemical software PHREEQC. The aim of this research 16 was to determine the dominant metal species of TM with potential hazardous health risk to local people to achieve 17 key prevention and pollution control. Herein, according to the morphological evolution of metal species, the activity 18 and concentration of the Zn species was calculated for both carcinogenic and non-carcinogenic health risk assessment. 19 The evaluation of the optimized human health risk demonstrated that the associated health risk of Zn (II) depended 20 predominantly on its metal speciation and was also affected by acidity and soil organic matter. Overall, the optimized 21 carcinogenic and non-carcinogenic risk value of Zn2S3 for adults was 2.01E-04 and for children was 1.31 22 respectively, resulting in corresponding hazardous risk to human, which accounted for high risk level of 61.5% and 23


29
With the development of economy and urbanization, many cities had carried out the adjustment of industrial  contribution rate can be explained and had clear physical significance (Lv et al., 2019). On the other hand, the PMF 125 method did not require the measurement of the source component spectrum, and used error estimates for each 126 individual data point in order to cope with missing and inaccurate data more reasonably (Mamut et al., 2017). The 127 data entered into the program include concentrations and equation-based uncertainties . 128 Compared with the traditional source apportionment methods, the PMF method could weigh all data that can 129 analyze the contribution rate of target variables (Tian et al., 2018). The sample metal concentration data matrix was 130 defined as: 131 (1) Where, the Cx t represented concentration of the t compound at the sampling point x; gxy represented contribution rate 132 of the yth source at the sampling point x; fy t was the mass fraction of the yth source at the metal t; ex t was the deviation 133 of the metal t at sample point x; Objective function Q was defined by Eq. (2) (Duan et al., 2020).

136
To analyze all sources of soil heavy metals in the study area, it was not feasible to rely on the analysis results of 137 sample data alone (Szopka etal., 2013). Through the analysis of the whole area, it was helpful to detect all the sources 138 of heavy metals. With the Ordinary Kriging interpolation method, the values of soil heavy metals can be obtained 139 from all the grids with a small number of samples. For grid A, under the condition that the spatial trend of metal 140 content was unknown constant, the Ordinary Kriging method can predict the metal content K (C(A0, K)) in any non-141 sampled grid as , 142 where, (A0,k) represents the content value of K in adjacent grid elements, A i i =1, 2…n. λ i i =1, 2…n is the 143 8 interpolation weight coefficient. 144 The corresponding equation based on unbiased estimation and minimum variance is as follows: 145 where was the variation function of a certain lag distance and was the Lagrange normal number. 146 According to the above steps, the spatial distribution of content of soil heavy metals in the scale grids can be 147 effectively generated. 148

149
Ecological risk index was used to assess environmental sensitivity and the pollution degree of soil heavy metals 150 (Kusin et al., 2018).The standardized heavy metal developed by Hakanson was used for evaluation basis 151 (Hakanson, 1980). It was calculated by the following equation for the potential ecological risk index of an individual 152 Where c represented the measured concentration of heavy metal and c0 was the background value of heavy metal 155 in soil. Er i represented the environmental risk index of heavy metal i and Tr i represented the toxicity response 156 coefficient of heavy metal i, which mainly reflected the toxicity level of heavy metals and the sensitivity of the 157 environment to heavy metal pollution. 158 The potential ecological risk index (RI), based on heavy metal concentration and toxicity, was conducted to 159 calculate the integrated potential environmental risk for the total hazard heavy metal (Zhu et  The toxicity response coefficient for As, Cd, Cr, Ni, Zn and Pb were 10, 30, 2, 6, 1, and 5, respectively (Guo et al., 163 2010a). 164

165
PHREEQC was a computer software which was used to calculate a variety of low temperature hydrogeochemical 166 reactions. It could offer a more comprehensive range of databases based on geochemical model. Compared with 167 experimental results, it indicated that the results of PHREEQC were scientific (Masindi et al., 2018). PHREEQC can 168 calculate the species of substance formation and the saturation index and simulate the geochemical inversion process. 169 It was also capable of investigating geochemical reactions in a variety of natural and human-affected environments, 170 including acid mine water discharge, radioactive waste isolation, pollutant transport, nutrient enrichment, natural and 171 artificial aquifer restoration, aquifer reserve recovery, potable water treatment, laboratory experiments, and regional 172 aquifer systems. The software of PHREEQC had a powerful thermodynamics database for input and operation. It can 173 be very convenient to obtain the molar concentration and the activity of the component, pH, pE, saturation index, etc 174 (Khoshgoftarmanesh et al., 2006). 175

176
The human health risk model, which includes an analysis of the average daily dose (ADD), as well as estimation 177 of carcinogenic and non-carcinogenic risks as determined by the US Environmental Protection Agency, had been 178 proved successful and has been successfully adopted globally (Tapia-Gatica et al., 2019). Among the exposures of 179 toxic heavy metal elements to human body through ingestion, inhalation, and dermal contact, the pathway of ingestion 180 was considered in this study to estimate the ADD of targeted metal, owing to that it was the most significant pathway 181 for toxic heavy metals (Giri and Singh, 2015;. 182 In general, the average daily dose for human health risk assessment was always estimated in terms of total 183 concentration of toxic element, but this tended to ignore the bioavailability and absorption of toxic elements. Heavy 184 metals in nature exsited in a variety of ionic forms, and not every kind of ion form with toxicity. Optimized human 185 10 health risk assessment through weight analysis could avoid imperfect and uncertainty of traditional risk assessment 186 model. Therefore, the ADD of TM species used for risk assessment were modified according to the following formula : 187 Where, Cp,q was the concentration of the q ionic speciation in p heavy metal element (mg .L −1 ), Mp,q was the molar 188 concentration of the q ionic speciation in p heavy metal element (mol. L −1 ), mp represented the relative atomic mass 189 of p heavy metal, np,q was the number of TM from the q ionic speciation in p heavy metal element (i.e. ZnSO4, np,q = 190 1); wp,q was the weight value of the q ionic speciation in p heavy metal element; Ap, q was the activity of the q ionic 191 speciation in p heavy metal element; Cp,q' was the modified concentration of the q ionic speciation in p heavy metal 192 element (mg. L −1 ); rp, q represented the weight assignment of the q ionic speciation in p heavy metal element.
According to the actual daily situation of the local resident population, the average annual exposure frequency (EF) 11 of the area was used for 365 days. RfD is the reference dose (mg/kg . d); IR was the water infestation rate (L.day −1 ) 199 and BW was the body weight (kg); AT reprensented the averaging time (day); ED is the exposure duration, wher 6 200 and 30 years was used for children and adults, respectively ( The study area of Suzhou city was located in the middle of the Yangtze river delta, the southeast of Jiangsu province, 211 located at 119°55'E to 121°20'E, 30°47''N to 32°02'N, with the total area of 8657.32 km 2 . The city was low-lying and 212 flat, with many rivers and lakes. Most of the water surface of taihu lake were in suzhou. Rivers, lakes and beaches 213 account for 36.6% of the city's land area. Suzhou was a subtropical monsoon maritime climate, with an average 214 temperature of 17.8℃ and precipitation of 1369.2mm in 2018. The prevailing wind direction was southeast wind. In 215 the shallow layer, the clay soil with slight deformation and high strength was mainly grey, with compact texture. 216 The study area focused on an abandoned dye factory, which was surrounded by four land uses, such as agricultural, 217 12 urban and industrial land, and a mixed zone. The research dyestuff factory was mainly engaged in the production and 218 sales of neutral dyestuff, cationic dyestuff and reactive dyestuff. The factory may have discharged a large amount of 219 heavy metals before being relocated and abandoned. According to the principles of Technical Specifications for Soil 220 Environmental Monitoring, 30 topsoil samples were randomly collected from the target center of the 1km grid in the 221 region on sunny days in September 2019, and the distribution of sampling points were shown in Fig. 2. 222

Collection and analysis 223
The soil samples were taken back to the laboratory and naturally dried and ground crushing, first through a 20 224 mesh sieve for pH analysis, then used for the determination of physical and chemical properties (Papa et al., 2010). 225 The processed soil samples were dissolved in the mixture acid solution (HNO3-HF-HClO4) at a high temperature of 226 150℃ for 4h (Bryanin et al., 2019). Soil pH was determined with a ratio of 2:5 (w/v) soil/water mixture using a pH  respectively. This phenomenon preliminarily indicated that the two heavy metal elements were not naturally enriched 248 in the study area, which may be caused by anthropogenic factors (Jiang et al., 2016). The relatively high coefficient 249 of variation (CV) of heavy metal element were Cd, Pb and Zn, with the CV values of 38.70%, 28.56% and 53.09%, 250 14 respectively. Among them, Zn and Cd were higher, indicating that they were considered to be affected by more factors. 251 The relatively low total CVs for Ni and Cr were mainly due to their absence of external disturbances, an assumption 252 based on previous studies (Onianwa and Fakayode, 2000). The higer values of As distributed in in the farm soil ( Fig.  253 3), indicating that the overuse of fertilizer and pesticide promoted the accumulation of As. Moreover, the research 254 considered that the larger values of Cd were affected by multivariate factors, owing to its broad distributions. Cr was 255 mainly concentrated in the northwest and central areas for bare land and dye factory, suggesting that the existence of     characterized for the total concentrations of As (53.8%). As was used largely in the manufacture of pesticides, but 287 also in pigments, dyes and paints (Zheng et al., 2003). 288 Comparatively speaking, industrial activities occupied the largest contribution to soil heavy metals in the study 289 area; followed by the human activities such as waste emission. The third largest contribution was made by agricultural 290 activities, which were related to the farmland in the study area. Owing to the frequent utilization of vehicle and fuel 291 near the roads, traffic activities accounted for the fourth largest contribution and the followed one was atmospheric 292 deposition. The results were consistent with the current situation of the whole study area, with abandoned factories 293 developed transportation around and frequent agricultural activities, indicating the close relationship between the dye 294 factory and its surrounding soil.   index of various heavy metals showed that the other areas in the study area were at low ecological risk level. In 334 summary, Er i values for the studied heavy metals were in the decreased order of Zn > Cd > As > Pb > Cr > Ni. But 335 on the whole, the industrial areas with serious pollution had strong ecological risk. The enrichment of Cd and Zn was 336 obvious in bare land and industrious area. However, the risk level of these two elements was significantly different, 337 which was mainly due to the large difference in toxicity coefficient between the two elements. Nevertheless, the 338 accumulation concentration of Zn heavy metal element was higher, which ultimately led to the higher ecological risk 339 index. The source analysis showed that the heavy metal of Zn in the whole area had an impact on the main pollution 340 sources around the abandoned dye plant. Therefore, it was necessary to strengthen the remediation of Zn element in 341 the study area to further prevent the heavy metal pollution in soil. On the basis of the methodology explored in this 342 paper, Zn was screened out as the TM for further study by means of potential ecological risk assessment. Further 343 research on TM to evaluate the health risks of its metal species and acquaintance with the metal speciation of TM 344 with higher risk levels can achieve targeted pollution control effectively and accurately, which could implement 345 environmental management policies. 346 Table 2 347 The relationship between ecological risk index and classification of heavy metal environmental pollution.

Parameters
Sampling points (n=30)     Table S1. Subsequently, the non-carcinogenic and carcinogenic risk of 390 Zn (II) species were evaluated by optimized ADD of Zn (II) species. Moreover, the non-carcinogenic and 391 carcinogenic risk based on non-modified and modified ADD for adults and children were shown in Table 5. Modified 392 ADD values for metal species were calculated according to Eq. (11) and ADD values for both adults and children 393 decreased in terms of activity of metal speciation, which was consistent with previous reports that metal absorbed by 394 organisms was primarily determined by the dissolved portion and that free ion activity was associated with these 395 processes. The modified ADD values of Zn (II) species were all higher than those of non-modified for both adults 396 and children. Due to the fluidity and resulting bioavailability of neutral ion, the ADD values of Zn(OH)2, ZnS and 397 ZnCO3, were much higher than those of negatively corresponding charged ions, such as Zn(OH)3, Zn2OH 3+ , Zn4S6 4-, 398 Zn(CO3)2 2etc. (Kotas and Stasicka, 2000;Accornero et al.,2010). 399 In terms of carcinogenic and non-carcinogenic health risk (Table 5), the HI values of Zn (II) species decreased in 400 the order as followed: Zn2S3 2-> ZnS > Zn 2+ > Zn4S6 4-> ZnCO3 > ZnHCO3 + > ZnSO4 > ZnF + > ZnOH + > Zn(CO3)2 2-> 401 Zn(OH)2 > ZnCl + > Zn(SO4)2 2-> ZnNO3 + > ZnCl2 > Zn(OH)3 -> Zn(NO3)2 > Zn2(OH) 3+ > ZnCl3 -> Zn(OH)4 2-> ZnCl4 2-. 402 Optimized HQ value of Zn2S3 2for children was 1.31 surpassing 1 compared with the non-optimized HQ value of 403 0.40, which exhibited non-carcinogenic risk. In addition, optimized HQ value of Zn2S3 2for adults was 0.94 404 approaching to 1. Thus special attention should be paid to this speciation because the potential toxicity of heavy metal 405 species may vary with the accumulation of concentration and change of environment factors. Accordingly, the 406 optimized CR value of Zn2S3 2for adults was 2.01E-04 and the level of risks surpassing 1E-04 are viewed as 407 unacceptable, therefore it was considered to pose significant hazardous effects on humans. However, the optimized 408 CR value of Zn2S3 2for children was 1.36E-05, in the range of 10 −6 −10 −4 , which was in a tolerable level associated 409 with the exposure and real environment (Granero and Domingo, 2003). Meanwhile, the HQ values of Zn4S6 4after 410 optimization were changed from 1.65E-06 and 1.11E-06 to 2.83E-05 and 1.74E-05 for adults and children, 411 respectively, as shown in Fig .7a. The CR values of Zn 2+ after optimization were changed from 5.17E-07 and 8.53E-412 24 07 and to 4.70E-05 and 1.87E-05 for adults and children, respectively. It was worth noting that the weighted average 413 optimized data were able to assess potential non-carcinogenic risks that were ignored by classical assessment methods. 414 Compared with the traditional evaluation model, using the activity of simulated metal species as a factor to optimize 415 human health risk model was more specific and reliable. The carcinogenic risk for ZnS underwent an acceptable 416 degree for adults after optimization with the CR values of 4.10E-05. However, before modification, the non-417 carcinogenic risk values of ZnS for adults was 3.62E-06, less than 10 -4 , which further underscores the need to modify 418 average daily doses to obtain accurate and effective human health risk assessment, as shown in Fig .7b. Moreover,  419 the rest of the metal species of Zn (II) ion were all beyond the non-carcinogenic and carcinogenic risk level. The 420 results indicated that Zn2S3 2-, ZnS, Zn4S6 4and Zn 2+ should be converted to other forms to reduce adverse effects on 421

humans. 422
Traditional risk assessment methods tended to ignore, overestimate or underestimate regional risk level when 423 conducting potential risk assessment. Therefore, when evaluating the relationship between the bioavailability of 424 various metals and health risks, attention should be paid to the absorbability of organisms to different metal ion forms. 425 Hence, it was necessary to assess the level of carcinogenic and non-carcinogenic risk of each metal species for 426 regional management and rehabilitation. A prerequisite for reducing the level of health risks in the region was the 427 screening of high-risk metal forms, which should be converted to non-toxic forms or microbial decomposition for 428 targeted environmental management.   In this paper, Zn2S3 2-, ZnS, Zn 2+ and Zn4S6 4was evaluated as the primary metal species by optimized human health 435 risk assessment. However, this did not mean that other low-risk ion forms can be ignored because the total health risk 436 value of zinc element was shared by each metal species. Changes in the environment, such as pH and alkalinity, will 437 affect ion activity, thus affecting ion concentration and weight distribution, and eventually led to changes in the level 438 of health risks. Therefore, the governance of different regions should be combined with the local actual situation to 439 carry out specific analysis and research. 440

Uncertainty and limitation of the risk assessment 441
Uncertainty estimation was an important part of human health risk assessment. Monte Carlo simulation was used 442 to verify that the simulated values were close to those calculated based on human health risk assessment. The results 443 showed that the health risk and output of each factor were stable. However, SOM and pH values can strongly affect 444 soil quality and metal distribution and fraction, further triggering changes in soil bioavailability, exposure dose and 445 toxicity of heavy metals, and finally causing health risks of different degrees. At the same time, the metal speciation 446 may affected by the interaction of metals with humus in soil or soil solution. Other influencing factors, such as sand 447 and clay, should also be considered. In addition, due to the limitation of conditions, the number of sampling points 448 27 was limited, the analysis types of pollutants were not comprehensive enough, and there was a lack of historical 449 monitoring data of site pollution. Although there was some uncertainty, this study can come up with an unbiased 450 effective risk index, which could screen out targeted metal by classical ecological risk assessment and optimize the 451 health risk assessment based on the simulation of morphological evolution of certain heavy metal species. 452

453
In this study, an optimized methodology for the contribution of metal speciation to health risk was developed by 454 combining chemical balance modeling with activity correction for average daily doses through exposure pathways. 455 Primarily, the ecological health risk levels and main pollution sources of six heavy metals (As, Cd, Cr, Ni, Pb and 456 Zn) in the region were studied. Source apportionment based on PMF spatial scales was conducted and the spatial 457 distribution of Zn and Cd with high enrichment were mainly in the center of abandoned industrial area, which were 458 polluted by discharge of Industrial waste water. The optimization simulation of the target heavy metals with high 459 ecological risk levels was carried out to establish the optimal health risk model. The results showed that the ecological 460 risk level of zinc was high, which need to be paid more attention, and heavy metals were mainly concentrated near 461 abandoned factories, which required the government to carry out policy management for soil remediation. The 462 carcinogenic and non-carcinogenic risk assessment suggested that the health risks of the four metal species of Zn2S3 2-, 463 ZnS, Zn 2+ and Zn4S6 4need to be paid attention, which can be converted into other harmless substances for effective 464 degradation and emission. The addition of all forms also indicated carcinogenic and non-carcinogenic risks of 465 metallic zinc in this region. Moreover, the harm of heavy metals to human health was mainly caused by the forms of 466 absorbable ions, rather than the total concentration. Traditional risk models may overestimate or underestimate the 467 actual risk level of heavy metals. In addition, soil pH was a key factor affecting soil heavy metal content, and organic 468 matter and soil pH value had significant effects on inorganic zinc binding state. Therefore, changing soil conditions 469 can effectively transform metal speciation with high toxicity. 470 The method was accurate and reliable and the evaluation results can provide decision-making basis for heavy metal 471 28 screening and discharge as well as targeted remediation for site management. Soil pH was found to be a crucial factor 472 affecting soil Zn fraction, and carbonate-bound Zn can be significantly affected by both organic matter and pH of 473 soils. 474