Spatial-temporal Distribution of Spawning Grounds in the Western Guangdong Waters Based on Habitat Suitability Index Model

1South China Sea Fisheries Research Institute, Chinese Academy of Fishery Sciences; Guangdong Provincial Key Laboratory of Fishery Ecology and Environment; Scientific Observing and Experimental Station of South China Sea Fishery Resources and Environment, Ministry of Agriculture and Rural Affairs; Guangzhou 510300, PR China. 2College of Marine Science, Shanghai Ocean University, Shanghai 201306, PR China. 3College of Chemistry and Environment, Guangdong Ocean University, Zhanjiang 524088, PR China. Article Information Received 03 March 2020 Revised 01 June 2020 Accepted 29 July 2020 Available online 21 October 2020


INTRODUCTION
S pawning ground is the waters for mating, spawning, hatching and growing of fish, shrimp and shellfish. Therefore, it is an essential space for aquatic creatures living and breeding (Wan et al., 2010).
Fish egg is a crucial stage in the early life of fishery population, which is also the stage most vulnerable to the changes of habitat environment. Fishes concentrate and spawn in an environment suitable for the survival and growth of fish eggs to improve the survival rate of their off springs (Bellier et al., 2010). Previous studies proposed that spatial distribution of the fish eggs were closely related to sea surface temperature, sea surface salinity, depth (Lelièvre et al., 2014) and seafloor relief (Flores et al., 2019) of the spawning ground.
Western Guangdong Waters (WGW) located in the northern South China Sea (SCS), stretches from the estuary of the Pearl River to Hainan Island, and meets the Beibu Gulf in the Qiongzhou Strait (Su et al., 2019). It is an important area for many commercial fishes spawning, feeding, breeding and migration (Yu et al., 2019c). Previous studies centered on biological characteristics of fish eggs and species composition (Emel'yanova and Pavlov, 2012), which provided fundamental information of spawning grounds. However, the response mechanism of spawning grounds to habitat environmental factors was less known.
Habitat suitability index (HSI) was widely applied in studies on fishery forecasts for being able to simulate the organism's response to environmental factors (Eastwood et al., 2001). The HSI model had developed rapidly in fishery research by combining with satellite remote sensing technology (Yu et al., 2019b). For example, HSI was used to forecast the distribution of the fishing ground of Ommastrephes bartramii (Yu et al., 2020), and analyze the optimum fishery habitats in the spawning season (Gillenwater et al., 2006). HSI was also used to analyze the characteristics of its fishing ground of Dosidicus gigas, and explore the relationship between its habitat quality and environmental factors (Yu and Chen, 2018b). This study built a model for fish eggs and environmental factors in the O n l i n e

F i r s t A r t i c l e
WGW based on HSI, explored the impact of environmental changes on the location and size of spawning grounds, analyzed the spatial-temporal changes of fish-egg density, and forecast the location and abundance of spawning grounds in the study waters. Satellite remote sensing offers all-weather, largescale and high-resolution sea surface data, and has already been applied in marine fisheries (Yu et al., 2019a). This study used field investigation data and remote sensing data, built and verified the HSI model based on marine environmental factors. The relationship between major fishes' spawning grounds and environmental factors in the WGW was also explored. The results of this study provide references for protecting the habitats of major fishes in the northern SCS.

Investigation and identification of fish eggs
Fish eggs data came from the spawning ground investigation in the WGW during April to June in 2014 and 2015 respectively, covering the areas of 110°-113°E and 19.5°-22°N (Fig. 1). Fish eggs were sampled by macro-plankton nets, with hauling speed of 1.5 n mile/h. Fish eggs were preserved in 5% formaldehyde solution, for morphological characteristics and counting (Wang et al., 2010). Fish-egg data were calculated as one fishing area according to the spatial resolution (0.25°×0.25°), and the unit of fish-egg density was ind/1000m 3 . Fish eggs data included voyage of fishing boats, investigation time, longitude, latitude, and spawn density.

Satellite remote sensing data
Sea surface temperature (SST), sea surface salinity (SSS) and digital elevation model (DEM) of the seafloor relief were taken from satellite remote sensing (Yu et al., 2018c). SST was from NASA MODIS Aqua (https:// oceandata.sci.gsfc.nasa.gov/), with a temporal resolution of 8 days, a spatial resolution of 4 km. SSS was from Global Ocean Physical Reanalysis Product data of Copernicus Marine Environment Management Service (CMEMS, http://marine.copernicus.eu), with a temporal resolution of the month, and a spatial resolution 1/12°×1/12°. DEM was from Google Earth-elevation data, elevation level was 18, spatial resolution 8.85 m. Depth data was from elevation representation in DEM data.

Spatial autocorrelation analysis
Spatial autocorrelation refers to the following issue: values of an attribute at closer geographical sites are more similar (positive autocorrelation) or more dissimilar (negative autocorrelation) than values at two distant sites (Wang et al., 2016). At the beginning of the study, the global spatial autocorrelation (Cliff and Ord, 1982) analysis was used to determine whether there is a significant relationship between the fish-egg density and geographical factors (Latitude and longitude). Global spatial autocorrelation parameters include Global Moran's I and Getis-Ord General G (Swetnam et al., 2015). They were used to detect global spatial cluster and its variations in the study. Global Moran's I is calculated as follows (Ren et al., 2020): x i is fish-egg density of pixel i, x j is fish-egg density of pixel j, x is average fish-egg density, n is the number, and W ij is the aggregate of all spatial weights. If the pixel i and pixel j are adjacent, the value of corresponding element in the matrix W ij is 1, otherwise it is 0.
To investigate the statistical significance of the Moran's I statistic, Z(I) is calculated as follows: where E(I) is the expected value of I: E(I) = − 1/(n − 1), and VAR(I) is the expected variance of I: VAR(I)=E(I 2 )-E(I 2 ). When a significance level is established, a Moran's I approaching +1 indicates that fish-egg density data is spatially correlative. When the value of Moran's I is close to −1, a discrete data pattern is observed. If the Moran's I value is close to 0 When |Z(I)|> 1.96, P value < 0.05 (Xiao et al., 2018), the null hypothesis is accepted and the fishegg density data are distributed randomly (Getis and Ord, 1992).
Getis-Ord General G is determined as follows (Getis and Ord, 1992): where x i is fish-egg density of pixel i, x j is fish-egg density of pixel j, and W ij is the spatial weight inversely O n l i n e

F i r s t A r t i c l e
correlated with the distance between the two locations. The expectation of G and Z(I) is calculated as follows: In general, if the value of G is greater than E(G), the highvalue data tend to cluster. Otherwise, low-value data tend to cluster. The fish-egg density data in the regions distributes randomly when G is equal to E(G) (Ren et al., 2020).

Suitability index model
Suitability index (SI) model was built accordingly to fish-egg density, SST, SSS and Depth data. This research supposed the waters with the highest fish-egg density was the optimum spawning area, HSI was 1. And the waters with the least fish-egg density (0) was the unsuited spawning area, HSI was 0. Of the spawn data, 80% was applied to build the model, and the left 20% to verify the model. SI of the single-factor habitat was calculated using the equation (Hua et al., 2020).
Where SI i refers to the suitability index in the i th month; Y i is the fish-egg density in the i th month; Y i,max is the maximum fish-egg density being caught in the i th month.

Habitat suitability index model
Using the simple non-linear regression (Fan et al., 2015;Li et al., 2016), the relation model between SI and SST, SSS and Depth could be built, as the equation below.
Where, A and B indicate coefficients of the model, and X the value of environmental factors. HSI was calculated based on the single-factor SI using Maximum Model (MAXM) (Yu et al., 2019b), Minimum Model (MINM) (Van der Lee et al., 2006), Arithmetic Mean Model (AMM) (Yu et al., 2018a) and Geometric Mean Model (GMM) (Tian et al., 2009), in below equations.
Where, SST-SI, SSS-SI and Depth-SI were the SIvalue of SST, SSS and Depth.
The left 20% investigation data were taken into the 4 equations to verify the accuracy of the model. The forecast was considered accurate if the error between theoretical value and actual value was less than 0.4, and inaccurate if not (Tian et al., 2009). HSI distribution was plotted through ArcGIS 10.5 (Wang et al., 2019).

Species identification
According to the morphological characteristics, major fish eggs in the area belong to: Trichiurus haumela, Carangidae, Nemipteras virgatus, Sardinella aurita, Anchoviella commersonii, which are consistence with related researches Zhao and Jia, 2020). The fish-egg density of major fish accounts for more than 90% of the total density. These fish eggs are popular commercial fish eggs and generally distributed in the tropical and warm waters (Yu et al., 2019b).

Global spatial autocorrelation of spawning grounds
In Table I, the Global Moran's I, Getis-Ord General G, their expectations and z-scores were calculated with fish-egg density data in the study area from April to June. The values of Global Moran's I varied from 0.22 to 0.36 and were significant at the 95% confidence level (Z (I) > 1.96, P value<0.05). This indicated that fish-egg density in the study area showed significant spatial correlation from April to June. General G observations are all larger than expectations and were significant at the 95% confidence level (Z (I) > 1.96, P value<0.05). This suggested that fishegg density in the study area tended to be a cluster of high values. Between April and June, the Moran's I decreased from 0.36 to 0.22, demonstrating a significant decrease in the fish-egg density cluster. The difference between estimated and expected General G decreased, indicating that the correlation degree of high-value fish-egg density data decreased. Moran's I peaked at 0.36 in April as did the difference between measured and expected value of General G. In June, the estimated value of Moran's I (0.22) and the difference between estimated and expected General G were the lowest for the study period. The global spatial autocorrelation analysis has shown that the fish-egg density distribution in the study area has significant spatial Spatial-temporal Distribution of Spawning Grounds O n l i n e

F i r s t A r t i c l e
autocorrelation and the high value area has been cluster.

SI of environmental factors
According to the simple linear regression, SI of each environmental factor and fish-egg density in different months were shown in Table II. The significance test (P < 0.05) indicated that the models were quite accurate. Through the fitting of SST, SSS and Depth, the comparison chart between the actual fish-egg density and the fitting result was shown in Figure 2. SI model of fitting factors in each month reflected variations of fish-egg density with the changes of environmental factors. All models showed unimodal distribution, and the simulation results were basically correct.

Optimum value of environmental factors
By obtaining the first derivative of each SI model, making it equal to zero, the value at the peak SI of each environmental factor could be calculated, that is, the optimum value of the spawning ground environmental factor in this month (Li et al., 2016). The results were shown in Table III. As for SST, the suitable range for fish eggs in April was 24.62-25.02°C, the optimum SST was 24.82°C. The suitable SST range for fish eggs in May was 27.96-28.31°C, the optimum was 28.13°C. In June, the suitable SST range for fish eggs was 29.25-29.50°C, the optimum was 29.37°C (Fig. 2, Table III).
As for SSS, the suitable range for fish eggs in April was 32.99-33.19 Practical Salinity Unit (PSU), the optimum SSS was 33.09 PSU. The suitable SSS range for fish eggs in May was 34.67-36.14 PSU, the optimum was 35.40 PSU. In June, the suitable SSS range for fish eggs was 33.39-34.60 PSU, and the optimum was 34.00 PSU (Fig. 2, Table III).
As for Depth, the suitable range for fish eggs in April was 8.55-9.65 m, the optimum was 9.1 m. The suitable Depth range in May was 11.74-13.36 m, the optimum was 12.55 m. The suitable depth range in June was 15.67-24.28 m, and the optimum was 19.97 m (Fig. 2, Table III).

Verification of HSI model
AMM, GMM, MAXM and MINM were applied to establish the HSI model, and data for model verification were imported to the models to get the theoretical HSI value and compare with the actual SI value. The model was considered accurate if the error between the HSI value and actual value was less than 0.4 (Tian et al., 2009), the precision verification results of the models were shown in Table IV. GMM gained the highest precision in April (95.60%) and May (90.91%), and MAXM had the lowest precision in April (15.38%) and May (18.18%) (Table III). In June, MINM had the highest precision (63.64%), and MAXM had the lowest precision (14.14%). Therefore, GMM showed the highest precision in April and May, and MINM had the highest precision in June. According to the verification results, GMM was applied to establish HSI model with SI data in April and May, and MINM was applied to build HSI model with the data in June. The precision of the established HSI model in each month was higher than 60% (Table IV).
HSI forecast results basically complied with the distribution of fish eggs (Fig. 3). Of which, fish-egg density in April was higher (>5000ind/1000m 3 ) in the northwest of the study waters (110.5°-111°E, 21°-21.5°N), HSI in this region was also higher (>0.6). In the northeast of the study waters (111.5°-113°E, 21°-21.5°N), fish-egg density was lower (<2500ind/1000m 3 ), and HSI was also lower (<0.4). The forecast results complied with the actual distribution in these areas. It was noticed that in the waters

O n l i n e F i r s t A r t i c l e
Spatial-temporal Distribution of Spawning Grounds 5 Fig. 2. Comparison between actual fish-egg density and fitting results. around the Qiongzhou Strait (110.5°-111°E, 20°-20.5°N), fish-egg density was higher (>5000ind/1000m 3 ), but HSI in this region was lower (<0.4).
Fish-egg density in the northwest of the study waters (110.5°-111°E, 21°-21.5°N) in May was higher (>5000ind/1000m 3 ), HSI in this region was also higher O n l i n e

Relationship between HSI and investigation data
HSI model can be built according to the non-linear relationship between fish-egg distribution and marine environment, its advantages of forecasting the distribution of spawning grounds accurately (Gillenwater et al., 2006). The spatial autocorrelation analysis showed that the distribution of fish eggs was clustered and the areas with high fish-egg density were clustered (Table I). HSI model can better forecast high-density region. The high-density region (northwest from April to June) and low-density region (northeast in April and June, and the Qiongzhou Strait in June) forecast by the model complied with the actual distribution of fish-egg density (Fig. 3). Previous studies showed that the spatial distribution of high value area forecasted by HSI model was correlation with the area of high stock density (Yu et al., 2018a;Hua et al., 2020), which is consistent with this research. The tendency of high-density fish-egg shifting southeastwards in April and May, and high-density fish-egg shifting southwards in May and June was the same as the results in the HSI model (Fig. 3). Previous studies have shown that the HSI model can accurately forecast spatial position changes (Hua et al., 2020).
However, the regions with actual high fish-egg density (the north region in May, and the Qiongzhou Strait in April and May) were marked as the low-density regions by the model, thus the forecast results varied greatly from the actual distribution of fish-egg density (Fig. 3). As fish eggs could not swim or move spontaneously (Wan et al., 2010), fish eggs floating in water were vulnerable to current. Therefore, complicate seafloor, which slows down the current, enables fish eggs to gather and provide a favorable survival environment for the spawning groups (Yi et al., 2013). Seafloor around the Qiongzhou Strait had greater fluctuations (dense isobath), its large seafloor slopes made complicate underwater environment (Fig. 1b), providing an ideal environment for spawning, resulting in high fish-egg density in this area during April and May. In this research, the seafloor factor, such as underwater topographic slope, was not considered in the HSI model, partly resulted to the deviation between HSI forecast and actual value in the Qiongzhou Strait.
High-value regions in the WGW forecast by HSI model showed the distribution tendency of "higher west and lower east" (Fig. 3), and continuously shifted southwards with the time. This was connected with the ocean current in this area. The previous study showed that there was a westward coastal current all the year-round in the WGW. The current flew westwards and then shifted towards the south nearby the Leizhou Peninsula, and finally entered the Beibu Gulf after flowing westwards through the Qiongzhou Strait (Wang et al., 2010). Fish eggs floating in water in each moth would be influenced by the westward coastal current and drifted to the west, then gathered nearby the Leizhou Peninsula, and finally formed the distribution tendency of higher west and lower east. In addition, the westward coastal current shifted southwards near the Leizhou Peninsula, resulting in the southwards shifting of the high HSI region in each month.

Relationship between the spatial distribution of spawning grounds and marine environment
Water temperature is a key factor influencing quantity and distribution of fish eggs, because it could influence quantity, distribution and population structure by affecting the gonad development and spawning migration of adult fishes (Mendiola et al., 2006). This research obtained the optimum SST for fish eggs in the WGW was 24.82 °C (April), 28.13 °C (May) and 29.37 °C (June), respectively (Table III). Previous researches showed that the optimum water temperature for fish eggs of Trichiurus haumela in the WGW was 25-28 °C (Lin, 1981), which complied with the results of this research. From April to June, suitable temperature for fish eggs increased gradually, but the high HSI area in June achieved the minimum, and fish-egg density in June reduced greatly from that in April and May (Fig. 3). The water temperature had a substantial impact on the process of metamorphosis (Régnier et al., 2018) and hatching speed (Poloczanska et al., 2013). In addition, higher water temperature in June (Fig. 2g) would reduce the survival rate of fish eggs that intolerable to high temperatures. For the fish eggs suitable for higher temperatures, the higher water temperature would accelerate the development of fish eggs, shorten the O n l i n e

F i r s t A r t i c l e
Spatial-temporal Distribution of Spawning Grounds growth period, and promote the speed of hatching, which both led to the decrease of fish-egg density (Mendiola et al., 2006). As a result, water temperature provided a significant reference for defining the quantity and location of spawning grounds in the WGW. Salinity is also an important factor influencing the metabolism of fishes, and fish eggs as the most fragile stage in the fish life are most vulnerable to the environment (Giffard-Mena et al., 2020). According to the SI model, the optimum SSS for fish eggs in the WGW in the northern SCS was 33.09±0.1 PSU (April), 35.40±0.73 PSU (May), and 34.00±0.61 PSU (June) ( Table III). The previous study showed that the Fish eggs in the northern SCS lived in the regions with SSS of 33.24-34.66 PSU (Li et al., 2014), the optimum SSS for fish eggs of Trichiurus haumela was 33.0-34.5 PSU (Lin, 1981), the optimum SSS for fish eggs of Decapterus maruadsi was 32.0-35.0 PSU (Zhang, 1985), and 33-34.4 PSU for fish eggs of Sardinella aurita (Jiang and Lin, 1983), which complied with results of this research. By comparing the actual fishegg density and fitting results of the SI model, it was found that fish-egg density achieved high value (Fig. 2b, 2e, 2h) as SSS was higher (SSS>33). The possible reason was that salinity impacted the embryo development by changing the osmotic pressure of fish eggs, too high or low salinity would reduce the hatching rate by blocking the normal exchange between oosperms and surrounding media, and thus caused the developmental malformation of embryos (Giffard-Mena et al., 2020). Moreover, salinity also influenced the vertical distribution of fish eggs in water. In low-salinity waters, fish eggs were likely to concentrate in piles and lack of oxygen held back the spawn development. On the contrary, in high-salinity waters, fish eggs could suspend or float in water and easily breathe, which would contribute to a higher hatching rate of fish eggs (Nissling et al., 2017).
Fishes preferred to spawn in deeper waters with suitable conditions (Lelièvre et al., 2014). The optimum Depth obtained by the SI model was 9.1 m (April), 12.55 m (May) and 19.97 m (June) ( Table III). Fish eggs of Anchoviella commersonii in the study waters were mainly distributed in water depth < 20 m (Wu, 1989), and fish eggs of Sardinella aurita in water depth around 10 m (Zhang and Huang, 2003), which consisted with this study. Spawning grounds in the WGW were mainly distributed in water depth < 30 m, fish-egg density in shallow waters was obviously higher than that in deep waters (Fig. 2c, 2f,  2i). The breeding interval of fishes in deep water might be prolonged for the harsh environment and less food (annual to perennial), resulting in lower fish-egg density in deep waters than those in shallow waters (Fernandez-Arcaya et al., 2016).

Method selection of HSI model building
In the construction of the HSI model, the most applied methods include Continued Product Model (CPM), MAXM, MINM, AMM, GMM and so on (Vayghan et al., 2013). This study compared HSI models built on the basis of MAXM, MINM, AMM and GMM, the precision of the models varied greatly (Table IV). GMM-based HSI had the highest accuracy (Table IV), a geometrical mean of each SI factor was used to build the model. In the HSI model for Trachurus murphyi in the southeast Pacific, GMM could better define the habitat and predict the habitat area (Li et al., 2016). Using fewer factors in the model building at the same time could improve the accuracy of GMM (Terrell, 1985), in this research GMM has the highest accuracy in April and May. Simultaneous application of 3 SI factors helped improve the accuracy of the GMMbased model, thus GMM was used to build HSI model in April and May. MINM had the second highest accuracy only after GMM (Table IV), and its evaluation of fish-egg density was mostly conservative for using the minimum value of SI factors, thus it was suitable for low-density regions and often applied to evaluate and build the HSI model for fishery conservation areas (Van der Lee et al., 2006). According to the HSI model for Trachurus murphyi in the southeast Pacific Ocean, MINM had higher accuracy in months with lower resource density (Fang et al., 2010). According to the HSI model for Decapterus maruadsi in the north of the SCS, MINM was suitable for model building in spring with fewer resources (Fan et al., 2015). In this research, fish-egg density in June reduced greatly from that in April and May, and low-density spawning regions were widely distributed (Fig. 3c), thus MINM had the highest accuracy in June.
AMM was a frequently-used method of the building model, mostly applied in assessing the quantity of resources by building the model on the basis of SI factors' mean. The forecast results of AMM were less vulnerable to extreme SI values, but different influence degrees of SI factors were not considered, because all SI values were equally treated (Yu et al., 2018a). In the western Pacific, the HSI model based on the AMM could better reflect the distribution of Katsuwonus pelamis under different El Niño events (Yen et al., 2017). AMM showed lower accuracy in this research (Table IV). Using maximum values of SI factors to build HSI model, MAXM always led to too optimistic HSI results, and its accuracy was influenced by being limited by the maximum SI factors (Yu et al., 2019b). MAXM was suitable for forecasting the center of spawning grounds and the months with higher resource density (Fang et al., 2010). To build the HSI model for Trachurus japonicus in the north of the SCS, MAXM was only applicable for the seasons with high resource density (Yan et al., 2018).

O n l i n e F i r s t A r t i c l e
J. Yu et al. In this research, there were maximum values in each month for the extraordinarily high fish-egg density, such extraordinary maximum values in the model building led to over-estimation of fish-egg density. Therefore, MAXM had the lowest accuracy in this research.

CONCLUSIONS
In this study, the HSI model of the distribution of spawn grounds in the WGW was built, using satellite remote sensing and field investigation. The study found that the HSI model based on the GMM and MINM could better reflect the distribution and shifting of spawning grounds, and it could better reflect the relationship between the spawning grounds and the marine environment in the WGW. The distribution of spawning grounds in the WGW was forecasted also. Environmental factors which were acquired easily from remote sensing were used to build the HSI model. In order to improve the accuracy of HSI forecasting, other parameters such as spatial autocorrelation, seafloor relief, current, wind, etc. will be considered in the follow-up research.

O n l i n e F i r s t A r t i c l e
Spatial-temporal Distribution of Spawning Grounds O n l i n e