Assessing Metallic Pollution Using Taxonomic Diversity of Offshore Meiobenthic Copepods

: The Gulf of Gab è s, located on the south-east Tunisian coast, is an important maritime area, with great inﬂuence


Introduction
Chemical pollutants enter the marine ecosystem either directly by their discharge in the ocean or indirectly through contaminated rivers, winds or rain [1]. Furthermore, a

•
Posidonia habitats clearly decline in favor of Caulerpa, as shown by [22], • A very noticeable reduction of the fauna diversity in these meadows (loss of 2/3 of the settled bottom-living macrofauna, in particular the disappearance of two species (Pinna nobilis and Pinctada Radiata), which are well represented in spared sites, such as the Kerkennah islands). • −10 m bottom desiccation from the isobath by fish, • Greater customary appearance of red water [23,24]. Sediments serve as final reservoirs for the pollutants, especially trace metals and hydrocarbons. They differ primarily in their granulometry, composition and by their particles size [27]. Other compounds are known to be carcinogenic and/or mutagenic and belong to the categories of substances that damage the environment. To assess the ecotoxicity of such sedimentary xenobiotics, the characteristics and components of the small food we, meiobenthos, protists and microbs need to evaluated [28]. In addition, the Phosphogypsum contains all of the phosphate ore impurities that are dissolved and released into the seawater. These include iron, chromium and arsenic, the high toxicity of the two latter chemicals being widely known [25]. The phosphogypsum factories in the Gabès discharge acidic liquids and solid discharges into the sea that are heavily loaded with phosphogypsum (600-650 tons/hour on average) and other pollutants (cadmium, fluorine). It was reported that these factories had, since the beginning of their activity in 1972, discharged into the sea water around 90 million tons of phosphogypsum [25]. Due to this high level of contaminants with a negative impact on marine and coastal ecosystems, the Posidonia oceanica meadow's habitat has been significantly reduced [26]. Therefore, the assessment of the Gulf's ecological conditions is crucial for the effective implementation of upcoming environmental monitoring programs.
Sediments serve as final reservoirs for the pollutants, especially trace metals and hydrocarbons. They differ primarily in their granulometry, composition and by their particles size [27]. Other compounds are known to be carcinogenic and/or mutagenic and belong to the categories of substances that damage the environment. To assess the ecotoxicity of such sedimentary xenobiotics, the characteristics and components of the small food we, meiobenthos, protists and microbs need to evaluated [28]. In addition, the hydrocarbon distribution serves as an indicator of identification of the source of the pollutants; indeed [29] reported that organic matter compounds would serve an indicators and tracers for the marine pollution in the area of the Gabès Gulf.
Along with free-living nematodes, the copepods comprise the second major sedimentdwelling taxonomic group of marine meiobenthos [30]. Copepods are widespread in benthic habitats and represent an important part of juveniles fish diets [31]. These crustaceans, dominated by harpacticoids, are very sensitive to various types of environmental stressors and their potential role as reliable ecological indicators was highlighted numerous times in the past (Influence of the Marseille sewer's outflow into the sea, California's continental shelf (Santa Maria Basin), offshore platforms in the Gulf of Mexico, China's Bohai Sea, German Bight North Sea, and a tropical ridge lagoon on Zanzibar Island, Tanzania [32][33][34][35][36][37][38][39][40][41]). The harpacticoids are characterized by small dimensions and a short life span with no planktonic stage, but they are very abundant and diverse [42][43][44][45]. Copepods are widespread and lead a variety of lifestyles. In his initial biogeographic synthesis in 1948, Lang published an inventory of endemic harpacticoids [46], and subsequently [47] examined the eco-ethological information of this widely dispersed group. Harpacticoids are an outstanding predictor of the environmental contamination because they are more susceptible to contaminants compared to nematodes [48,49]. Therefore, harpacticoids have been extensively studied in both the South China Sea and the Baltic Sea [50,51].
Due to the difficulty of taxonomic identification, marine harpacticoids have been included in fewer studies by ecologists compared to nematodes [52][53][54]. Many harpacticoids are cosmopolitan despite their small size and restricted mobility. Benthic harpacticoids enter into the sediments in one of three ways: interstitially, by digging or epibenthically. These modes range from coastal basins to deep abyssal bottoms [55][56][57].
As in many parts of the world, the main focus on Tunisian marine meiofauna was on free-living nematodes, which were commonly used for environmental monitoring programs in the north of the country [58] or in laboratory ecotoxicological experiments [28,59]. The knowledge of the biology and ecology of meiobenthic copepods in the Mediterranean Sea bordering northern Tunisia lags behind that of other taxonomic groups [60][61][62]. On the other hand, compared to planktonic copepods [63][64][65][66][67][68][69][70][71][72][73][74][75][76], the study of meiobenthic copepods is relatively less advanced and less documented. Indeed, most publications do not reach the species level and are limited at the determination of the density [77][78][79][80][81]. This is mainly because of their modest size, the scarcity of qualified experts on taxonomy of meiobenthic copepods and the absence of a global database on the morphology and biometry traits of up-to-date described taxa. Given the taxonomic difficulties mentioned, scientists interested in meiobenthic copepods turned their attention to ecotoxicology, which only required the culture of easily identifiable species [36,37,82,83].
Limited studies were carried out on the exploration of the copepod fauna inhabiting the Mediterranean Tunisian coasts. This group of meiobenthic crustaceans in Tunisia was studied a long time ago by Monard (1935) [60], who was interested in marine harpacticoids in the Salammbô region and in the fishing in Sidi Ahmed (Bizerte). Another study was carried out by Amorri (2004) [62], who investigated the inhabitants of benthic copepods in the Bizerte lagoon. However, the previous studies aimed at analyzing their density, and showed that these inhabitants do not change much over time. The Tunisian faunistic lists of nematodes in Tunisia from 1981 to 2003 and actual ones are strikingly similar [84][85][86][87][88].
This study uses an integrative approach that combines taxonomic and morphometric diversity on the one hand and spatial and temporal diversity on the other. Here, the organo-mineral characterization of sediments from the Golf of Gabès was performed to explain the observed pattern of meiobenthic copepods. Furthermore, the Gulf of Gabès has one of the highest semidiurnal tide levels in the Mediterranean Sea, with a maximum amplitude of about 2 m [89]. This characteristic undoubtedly affects the biodiversity and productivity of planktonic and benthic copepods, resulting in diurnal cycles of dispersal movements into the water and sediment columns, respectively.
The current study addresses this knowledge gap by exploiting the copepod characteristics and developing a checklist of Gulf of Gabès' copepod species that will be used for the assessment of the ecological status of these marine complex habitats. Initially, the measurements of the level of contamination by metals and organometallic compounds of this gulf was carried out, followed by the investigation of the response of the sediment-dwelling copepods to these stressors.

Study Area
The study area was the gulf of Gabès, Tunisia (33 • 45 N-33 • 57 N and 9 • 57 E-10 • 18 E) ( Figure 1). In this area, the currents are regular and the stream comes from the NE or the east, and follows the NW and/or SW directions [90] (Figure 1). Sampling campaigns were carried out seasonally from the winter of 2004 to the autumn of 2005 at nine sites (Table 1). Site S1 was situated in north of the city of Gabès, whereas sites S2, S3 and S4 on the transect R1 are located 0.5 km from the discharge point of the industrial zone of Ghannouch City. Stations S5 and S6 were situated in the south of transect R1 by 2.45 and 13 km, respectively. The sites S7, S8 and S9 were situated on the transect R2, south of Gabès City, within the tourist area. The sampling stations were chosen as such as to cover a pollution gradient that starts from the industrial zone of Gabès City and ends in the gulf. The Gulf of Gabès is known for its low, muddy, sandy beaches and frequent sabkhas on the land side. The El Bibans lagoon (approximately 30,000 acres) is located further south [91,92]. A very high structural similarity was found along the Tunisian Coasts from Cap Bon to the Libyan border. Most of the Tunisian coast is covered by recent Quaternary sediments, some of which extend more than 100 km offshore [93]. The Gulf of Gabès has a wide area of marine abrasion dating from the Pliocene and Miocene periods, where the continental shelf is more extensive than elsewhere, although, this region is frequently covered by marine or continental formations from the recent Quaternary period [94]. Thus, the location of the −200 m isobaths are 250 km from the coastline [95].
Ben Othmen (1973) [94] observed the sedimentological succession at the level of the radial going through the parallel 34 • 30 N while examining the bottoms of the Gulf of Gabès. From −13 to −50 m sand, sand mud and sandy mud sequentially spread across the bottom were observed; from −50 to −80 m sand and mud are gradually spread across the bottom, whereas from −80 to −200 m, sandy-muddy sediment was observed.
Several studies [96][97][98] indicate that the Gulf of Gabès has several unique characteristics, such as the progressive siltation of the area since 1954 and the very fine grain size of the sediment near the seagrass beds. These structures also serve as sediment traps [99].

Sediment Sampling
Sediment samples were taken with a Van Veen grab. In the laboratory, subsamples were taken, dried at 45 • C and used for grain-size measurements [100], whereas others were frozen at −20 • C for chemical analyses of Organic Matter (OM) and lipid contents.

Extraction of Carbohydrates
Following the removal of large shell debris, the sediments were sieved through a 2 mm mesh size sieve before freeze-drying. The free lipids were extracted (flushing, 60%; purge, 120 s; cycles, 2 × 2; total solvent collected, 500 mL) at 80 • C with an accelerated solvent extractor (Dionex ASE 100; Salt Lake City, UT, USA) and a mixture of CH 2 Cl 2 /MeOH (2/1, v/v). The resulting organic matter was fractioned into neutral, polar and acidic compounds using SPE LC-NH 2 cartridges (3 mL/500 mg; Supelco) [101]. According to [101], cartridges were conditioned with 5 mL hexane and samples were loaded at a 1 mL·min −1 rate. The carbohydrates were extracted with 5 mL hexane and alcohols with 5 mL CHCl 3 . Fractions were then eluted on a Cu column to remove elemental sulphur [102].

Trace Metal Analysis
The determination of trace metals was conducted according to the Loring method [103] and the atomic absorption spectrophotometry for sediments [104]. The metallic elements were determined by spectrophotometry Electrothermal Atomic Absorption (SAAE) using a graphite furnace and a color corrector continuous bottom with the Zeeman effect (Varian 220Z).

Fluorine Determination
In order to determine the fluorine content, we applied the method that was designed by [105], utilizing a fluoride selective electrode (DZ/T0167-2006). The accuracy of analysis was checked by comparing the samples with the standard, controlled samples and the errors for all analyzed samples proven to be less than 5%.

Copepod Study
Sediment samples were taken by SCUBA diving (Figure 1). At each sampling site, four hand cores (section of 10 cm 2 , inner diameter of 3.6 cm) were taken. Samples were immediately preserved in 4% buffered formaldehyde until further analysis in the laboratory. Meiobenthic organisms were first extracted with the aid of the levigation-decanting-sieving method of [106], and then stained with Rose Bengal (0.2 g·L −1 ; [107]). The copepods were sorted under a Leica MZ 12.5 stereomicroscope and subsequently transferred to glycerin. The taxonomic identification to species level was done with a Leica DMR microscope, based on the identification keys of [46,57,108,109]. Copepod individuals were measured (body length, maximum width) using a drawing tube attached to the microscope. Their body volumes were calculated using conversion factors given by [110]. The individual copepod biomass [µg, wet weight] was determined using a specific gravity of 1.075 given by [32] and the total copepod biomass was calculated and expressed per 10 cm 2 .

Data Processing
Principal Component Analyses (PCA) was applied to illustrate the pattern of chemical variability within a single factor, providing a comprehensive measure of pollution level (i.e., pollution index). Separate PCA ordinations were done for each season. The pollution index scores, as well as the copepod community-based indices (i.e., species richness, density, and biomass) were checked for Moran's spatial autocorrelation by using correlograms, according to [111][112][113], followed by Bonferroni multiple-comparison corrections [114]. Because data showed significant spatial autocorrelation (see below), conditional autoregressive models [115] were used to investigate the relationship between copepod richness, density and biomass (response variables) and pollution index (explanatory variable). These models assume that the response in each site is influenced not only by the explanatory variables per se, but also by the neighboring sites [116].
Three distance matrices were computed and used: (1) a matrix of geographic distance among sites, (2) a matrix of Euclidian distance in pollution level among sites and (3) a matrix of distance in species composition among sites computed using the Jaccard index. The correlation was made between pollution and species composition matrix, relative to the geographic distance matrix. The significance of the partial Mantel's test was assessed by 10,000 permutations.
Data from the four seasons were treated separately. Spatial autocorrelation analyses and spatial regressions were done in SAM software [117], whereas partial Mantel tests were done in XLSTAT (version 2012; Addinsoft, Paris, France). Throughout the manuscript, means are presented as ± SE.

Abiotic Parameters
In general, there were significant differences across sampling sites in terms of reconstituted oxygen, sediment type and pollutant content ( Table 2). The sediments of radial R1 are mostly silty facies, ranging from fine silts (S2 and S3) to clay silts (S4), whereas sediments at other stations containing sandy granulometry, clearly differ in grain size.  By extracting data from each season separately, a PCA ordination summarized the five measured pollution parameters into a single factor with high eigenvalues (winter: 4.374, spring: 4.657, summer: 4.231, autumn: 4.184), explaining the overall variance from the original datasets (winter: 87%, spring: 93%, summer: 85%, autumn: 84%). This factor was positively correlated with the concentration of pollutants (Table 3), providing a composite measure of pollution level. This pollution index exhibited significant spatial autocorrelation, as shown by the Moran's correlograms ( Figure 2). The sites situated closer were more similar to each another than expected by random distribution of pollution in the studied area.

Copepod Parameters
Thirty-eight copepod species were identified, including five species new to science, comprising three orders, 14 families and 23 genera (Table 4). Three new species belonged to the Ectinosomatidae family, genus Pseudobradya sp. and one to the Canthocamptidae family, genus Stenocaris sp. One species, Apodopsyllus gabesensis (Paramesochridae family), according to Amorri et al., 2010, was discovered and described in the course of this study [61].

Copepod Parameters
Thirty-eight copepod species were identified, including five species new to science, comprising three orders, 14 families and 23 genera (Table 4). Three new species belonged to the Ectinosomatidae family, genus Pseudobradya sp. and one to the Canthocamptidae family, genus Stenocaris sp. One species, Apodopsyllus gabesensis (Paramesochridae family), according to Amorri et al., 2010, was discovered and described in the course of this study [61].  The species richness, density and biomass varied among sites and seasons (Table 5). There was a significant spatial autocorrelation among these parameters, except for the density during winter (Figures 3-5). Considering data from each season separately and using conditional autoregressive models that counted for the spatial autocorrelation, the copepod richness, density and biomass were negatively correlated with the pollution level ( Table 6, Figures 6-8). During each season, the copepod richness, density and biomass decreased in the most polluted sites. The results of the partial Mantel tests showed positive associations between species composition and pollution level, after controlling for the effect of geographic proximity between pairs of sites (winter: r = 0.927, p < 0.0001; spring: r = 0.935, p < 0.0001; summer: r = 0.918, p < 0.0001; autumn: r = 0.937, p < 0.0001). Overall, the results suggest that nearby sites were characterized by similar pollution levels and inhabited by similar copepod communities.         Table 6. Results of regression analyses on species richness, density and biomass (ln-transformed) of copepods (response variables) as functions of pollution index (ln(x + 2)-transformed) as an explanatory variable, using conditional autoregressive models that account for spatial autocorrelation in the data.     During all seasons combined, the Ectinosomatidae family dominated the copepod populations in sites S1, S5, S6, S7, S8 and S9, andthe Harpacticidae family dominated sites S3 and S4, whereas the Canuellidae family dominated in site S2.

Species Richness Density Biomass
Only five families (Ameiridae, Canuellidae, Ectinosomatidae, Harpacticidae and Laophontidae) were common in all sampling sites. The Thalestridae family was encountered in sites S1, S2, S3 and S4, whereas the Paramesochridae family was found only in site S1.
Considering all seasons combined, the most abundant and diverse was the Ectinosomatidae family with seven species (40.78%). The families found during winter were also found in spring, summer and autumn. During summer, just one family of benthic copepods, comprising seven individuals of the new species Apodopsyllus gabesensis [61], were found (Paramesochridae family); this species was also found during autumn.
The species richness increased on both sides of the industrial zone of Ghannouch towards the south. The sediments from sites S2, S3 and S4 comprised the lower taxonomic richness, and the latter site was dominated by Harpacticus gracilis. The species richness, which was the highest in winter, showed variations among sites S2, S3 and S4 and the others going north (S1) or south (S5, S6, S7, S8 and S9). The species richness fluctuated in space and time, following an increasing gradient from sites S2, S3 and S4 towards the south and north. During all seasons combined, the Ectinosomatidae family dominated the copepod populations in sites S1, S5, S6, S7, S8 and S9, andthe Harpacticidae family dominated sites S3 and S4, whereas the Canuellidae family dominated in site S2.
Only five families (Ameiridae, Canuellidae, Ectinosomatidae, Harpacticidae and Laophontidae) were common in all sampling sites. The Thalestridae family was encountered in sites S1, S2, S3 and S4, whereas the Paramesochridae family was found only in site S1.
Considering all seasons combined, the most abundant and diverse was the Ectinosomatidae family with seven species (40.78%). The families found during winter were also found in spring, summer and autumn. During summer, just one family of benthic copepods, comprising seven individuals of the new species Apodopsyllus gabesensis [61], were found (Paramesochridae family); this species was also found during autumn.
The species richness increased on both sides of the industrial zone of Ghannouch towards the south. The sediments from sites S2, S3 and S4 comprised the lower taxonomic richness, and the latter site was dominated by Harpacticus gracilis. The species richness, which was the highest in winter, showed variations among sites S2, S3 and S4 and the others going north (S1) or south (S5, S6, S7, S8 and S9). The species richness fluctuated in space and time, following an increasing gradient from sites S2, S3 and S4 towards the south and north.

Discussion
The current study provided the first large-scale investigation of benthic copepods' spatial distribution in south-east of Tunisia, as well as information on their biomass, species diversity and response to environmental influences, with a main focus on harpacticoids.

Variation of Environmental Variables
In general, the stations that were investigated in the Gulf of Gabès showed extremely heterogeneous surface sediments. Additionally, the obtained results show a typical dissolved oxygen mass. According to [118], the oversaturation limit is 8.8 mg·L −1 . In our study, the average value of the oxygen concentration measured at the bottom of the water resulted in a level close to the seawater's natural level (6.8 mg·L −1 ) [118].
The hydrocarbon content in the sediments sampled from the Gulf of Gabès was moderate-high [119,120]. These values were close to those measured on the coast of the city of Sfax [121,122] but higher compared to those reported from other coastal areas (Table 7). The high concentration of hydrocarbons from sites S2, S3 and S4 have a mixed origin: biogenic sources with low molecular weight n-alkanes (C14-C19) and degraded organic matter and oil contamination highlighted by the presence of MCUs [130][131][132][133]. The hydrocarbons from the coast of Gabès are the result of a mixture of natural and anthropogenic sources. The former input could have its origins from plant debris, carried by inland waters and discharged in coastal areas.
High levels of iron, zinc, cadmium and fluorine were found at sites S2, S3 and S4 (the radial R1) in all seasons. These results reflect the industrial activity in the Gabès region. The fluorine is an important component of phosphogypsum, an end-product of the Ghannouch industrial complex [134]. The sediments analysis from the Gulf of Gabès showed that the contamination levels above the tolerated threshold (i.e., lowest threshold effect level) was likely to cause mutations and the death of certain marine organisms.
The quantities of fluorine, iron, zinc and cadmium in the sediments of R1 were quite similar to those reported at the Ghannouch industrial complex's discharge region by Sarbaji et al. in 1993 [135]. The contents we measured were significantly correlated among them and with the sediments particle size (r > 0.70 and r > 0.90, respectively). The ability of fine silt to trap organic matter compared to the coarse sediment was associated with the higher levels of contaminants, such as trace metals [136].
The monitoring of the different sedimentary parameters (MOT, total hydrocarbons, n-alkanes, fatty acids, phosphorus (P), fluorine (F), Fe, Zn, Cd and Cu) along the infralittoral zone of the city of Gabès highlighted the impact of industrial discharges, without underestimating the importance of wind direction and peri-littoral currents in the distribution of pollutants in the sea [135]. In addition to the chemical contamination of the seawater, phosphogypsum deposits alter the ecology of the local marine environment and, according to [22], destroy the Posidonia meadow and reduce the coast's halieutic wealth [135]. The contaminants' concentrations were considered high compared to the Tunisian discharge standards for the public ocean domain [135]. In fact, the main compo-nents of the discharges coming from the Ghannouch industrial complex-especially those of phosphogypsum-included phosphorus and fluorine [137].
Finally, the distribution of the different sediment chemicals found in the Gulf of Gabes indicated a highly polluted zone near to the Ghannouch industrial complex's discharge (Radial R1). They are one of the most important sources of pollution in the Gulf's littoral and infralittoral zones.
The area near the northern coast of the discharge spot resulted in less contamination, and when heading south, the impact of the discharges also decreases (stations S6, S7, S8 and S9). The direction of the prevailing winds and currents plays a crucial role in the spread of contaminants into the sea [135]. In addition, several previous studies in the field [135,138], and those carried out in the other regions of Mediterranean [139], indicated that the self-cleaning capacity of seawater helps to mitigate the absorption of pollution as it moves from the coast to the open sea.
Metals originate from rock, soil and sediment weathering, but also from human activities. The measured trace metals from the current study were iron, cadmium and zinc, the dominant but equally among the most dangerous types of pollutants in marine habitats [140,141]. The current study showed the presence of high concentrations of Fe, Zn, Cu, Cd and fluorine that exceeded international standards (Table 8) [142,143]. The results agree with those obtained in similar habitats by [144] and in accordance with those of [135]. Table 8. Bibliographic comparison of sedimentary trace metals in the Gulf of Gabès with internationals norms [142] (LTEL: lowest threshold effect level) and fluorine [143]. Fluorides are also released into the environment (air, soil, water) by various human activities, leading to an increase beyond the natural background of fluorides in waters. Their concentrations can reach more than 100 times the natural content of inorganic fluorides in surface and ground waters, creating an ecological risk for aquatic organisms [145].

LTEL Gulf of Gabès (Current Work)
Aluminum smelting, steel, phosphate, fertilizer production, glass and enamel production, brick and ceramic production, glues and adhesives, fluoride-based pesticides and fluoridation of some drinking water supplies are examples of anthropogenic sources of inorganic fluorides [146][147][148][149]. The manufacturing of phosphate fertilizers, phosphoric acid and aluminum metallurgy are the main sources of fluoride emissions [150,151]. Humans who consume fluorides may become acutely or chronically poisoned [152][153][154].
The digestive and skeletal systems are the main systems affected long term by fluoride exposure [155][156][157].

Variation of Meiobenthic Copepod Quantitative Traits from the Gulf of Gabès
Overall, despite a large fluctuation in numbers across sites and seasons, the mean density of copepods in the Gulf of Gabès was generally low. The density values measured in the current study (Table 9) are similar to those of [57], where the highest values were recorded in intertidal sediments (1,000,000 ind./m 2) . In the same vein, the total biomass of copepods closely followed that of mean density. The values recorded (4.59 ± 0.18-5317.4 ± 573.14 (µg/10 cm 2 ) for the the Gulf of Gabès were higher compared to those from the lagoon of Bizerte (6.61 µg/10 cm 2 -155.90 µg/10 cm 2 ) or the Bizerte Bay (Rimel station: 165.99 µg/10 cm 2 ) [62]. The biomass values are similar to those measured by [32] in Anse de Cortiou before the installation of the Marseille wastewater treatment plant (63.4-8261 µg/10 cm 2 ). These biomass values are also comparable to those measured in the Bay of Bizerte, where the mean individual biomasses of copepods ranged between 0.94 ± 1.12 µg-6.53 ± 2.11 µg [62]. Furthermore, biomasses between 1 and 7 µg were reported by Wieser (1960), Guille and Soyer (1968) and Stripp (1969) (in [32]). Table 9. Comparison of the average numbers of copepods recorded in the Gulf of Gabès with data from polluted to very polluted environments in the literature.
Stenocaris sp. (Canthocamptidae) and Apodopsyllus gabesensis (Paramesochridae) had the lowest biomass because of their small body size, which is an adaptation to the minute dimensions in their interstitial microhabitats [54,109,170,172].

Taxonomic Diversity of Meiobenthic Copepods from Gulf of Gabès
The species richness was the highest in sites S1, S5, S6, S7, S8 and S9, which is comparable to those recorded in more clean areas from the Bizerte lagoon. An increasing gradient in species richness was detected from stations S2, S3 and S4 towards the south and north. In turn, the copepod communities from the polluted sites were species-poor, indicating that pollution of sediments with metals and hydrocarbons eliminated the most sensitive taxa. The sites S1, S4, S5, S6, S7, S8 and S9 were characterized by sandy sediments and hosted the most diverse copepod communities, dominated by four species: Pseudobradya sp.1, Halectinosoma aff. itoi, Halectinosoma herdmani and Pseudoameiropsis sp. It was proven that the abundance increases with the increment of the grain size [173,174]. H. herdmani, on the other hand, is a psammophilous species that prefers medium-sized sand. [47]. This species was dominant at sites S5, S6, S7, S8 and S9 on the Gulf of Gabès. However, [60] reported this species as very rare at Salammbô. The species Phyllopodopsyllus berrieri, which was dominant in site S1 during summer, is actually a cosmopolite species (i.e., Scotland, Norway, Rovigno, Algiers, Bermuda, Ceylon (Malaysia), Antarctica [171]).
In site S4, only two species were found throughout the year. In sites S2 and S3, characterized by medium silt sediments and S4 with silty clay sediment, and which were situated very close to the effluent discharge of the Ghannouch industrial complex, the copepod fauna was dominated by the species Harpacticus gracilis, Harpacticus flexus, Brianola stebleri and Halectinosoma curticorne. The former species was previously reported from the polluted and euryhaline waters of the lagoon of Tunis and La Goulette [60]. H. flexus was also found in euryhaline, muddy and sandy environments, with little disturbance in the La Goulette lagoon of Tunis, the Bay of Algiers, and the Adriatic and Banyuls seas [47,60]. This species prefers sandy-muddy environments with no algae [57]. B. stebleri is usually found in small numbers in coastal zones (Lagune de Tunis, [61], with a preference for muddy, euryhaline and more or less polluted biotopes. This species was previously reported in the Lagoon of Tunis before its development at La Goulette, Castiglione, Sète and Roscoff [60,171,175]. This copepod was replaced in spring by H. curticorne in site S2 (polluted clay sediment). The latter species prefers sandy-muddy sediments and was reported in high numbers in soiled muds (Roscoff, Salammbô, Algiers [60,171]. Co-occurring species with low abundances in sites S2, S3 and S4 were Enhydrosoma propinquum, Enhydrosoma sordidum, Heterolaophonte stroemii brevicaudata, Paralaophonte congenera, Delavalia tethysensis, Dactylopusia tisboides, Ameira parvula, Paralaophonte brevirostris and Phyllopodopsyllus sp. These species are known to prefer polluted muddy environments (i.e., La Goulette, Roscoff, Banyuls, Sète, Genoa, Aegean Sea, Sardinia, Rovigno, Dalmatia, Suez Canal Red Sea Port vente and le Vivier) and are resilient to wide variations in salinity [60,171,175]. These species were reported by [60] from the Mediterranean Sea, and are very common in Salammbô, La Goulette, Tunis Lake, Castiglione and in the port of Algiers.
Ameira scotti was the only common species in the surveyed sites from the Gulf of Gabès. This copepod species was previously reported from other confined habitats (i.e., Sète, Banyuls, Suez Canal, Castiglione, etc.) and dominant in muddy sites polluted with oil, such as in La Goulette and in the lagoon of Tunis and the Port of Algiers [60,171,175].
In the pristine sites, three species (i.e., Pseudobradya sp.1, Halectinosoma aff. itoi and Halectinosoma herdmani) were dominant, mainly in the north site S1, but also in the southern sites with sandy bottoms. An explanation for this pattern could be that the Ectinosomatidae family comprises cosmopolite, euryhaline and eurythermal species [57,109].
In addition to anthropogenic factors, the taxonomic composition of marine copepod assemblages can equally be explained through stochastic factors, such as the spatial variability of environmental conditions, the specific ecological zonation of harpacticoids, which are linked to the nature and availability of trophic resources, local mechanical agitation, interspecific competition and reproductive activity [57], whereas the sediment type (fine or coarse) should not be overlooked. In fact, it has been demonstrated that sediment type plays an important role for meiobenthic copepods in terms of population dynamics, density and species [41,176]. The presence of a specific substrate directly affects the vertical and horizontal distribution of species within sediments [56]. Harpacticoids, which are known to be intolerant to anaerobic environments, accumulate in the oxygen-rich layers of the sediments and are restrained in the upper layers of muddy sediments, but they are occasionally found in depths of 50-100 cm [57]. Additionally, copepod populations can migrate vertically, either entirely or partially, depending on the tides, temperature, dissolved oxygen concentration and other environmental factors.
According to [60], high levels of pollution lead to scarce harpacticoid fauna, despite their well-known resilience and increased taxonomic richness. In time, the taxonomic composition of such communities can strongly fluctuate, potentially related to the development cycles of the dominant species from the Gulf of Gabès, but also due to environmental stochasticity. Such patterns are nonetheless not easily discernible from anthropogenic influences, such as the fine sediment fraction that favors the trapping of several metals such as zinc, iron and cadmium [136,177]. Statistical analyses showed an antagonism between density, biomass and species richness of the copepod assemblages and sedimentary parameters, such as the concentrations of fluorine, carbohydrates and trace metals (i.e., Fe, Cd, Zn).
The sensitivity of copepods to certain contaminants was proven, such as hydrocarbons [178], high MOT load [179], trace metals [180,181] or mixed pollution [182], along with the structure of the facies (i.e., fine or coarse). The latter aspect was proven to play an important part in meiobenthic copepods' dynamics, densities and taxonomic composition [176]. The example of the copepod Cletodes tuberculatus is instrumental, and is known to be sensitive to changes in sediment structure and associated organic carbon content [33]. Moreover, Ref. [34] observed that harpacticoids dwelling in hydrocarbon-contaminated sediments in the Gulf of Mexico decreased their diversity. In the same context, Ref. [37] also demonstrated that the harpacticoids are reliable bioindicators of pollution with hydrocarbons. According to [35], several species of ectinosomatid copepods were sensitive to sediment contamination with hydrocarbons. Similarly, Ref. [36] demonstrated that the harpacticoid copepod Tigriopus brevicornis is a good pollution indicator [38] and reported that early generations of littoral Microarthridion and Amphiascus tenuiremis were affected by trace mixed sediment pollution.
In conclusion, our findings support the necessity of including benthic copepod populations in metal pollution biomonitoring programs in both Mediterranean and Tunisian marine environments.

Conclusions
Marine coastal sediments usually comprise the final collector of trace metals and hydrocarbons. In the current survey, the impact of anthropogenic pollution on the distribution and taxonomic composition of the meiobentic copepods from the littoral area of the City of Gabès (Tunisia), in the southern Mediterranean Sea was investigated for the first time. The results showed that the pollution index exhibited significant spatial autocorrelations, as shown by Moran's correlograms. The sites situated closer to each other were more similar than what was expected by pure random distribution of pollution. The species richness, density and biomass of marine copepods were negatively related to pollution levels, a pattern consistent throughout all seasons. High sediment loads with hydrocarbons and fluorine were the main abiotic factors controlling the spatio-temporal distribution of copepod communities, closely followed by trace metals such as iron, cadmium and zinc.
The sites with the least anthropogenic impact were those situated the furthest from the industrial complex of Ghannouch, characterized by low sediments concentrations of fluorine, iron, cadmium, zinc as well as hydrocarbons. Consequently, these zones hosted the most diverse copepod communities, whereas the sampling sites S2, S3 and S4, which were the most disturbed, were characterized by species-poor and low-abundance associations. The current survey also highlights the fact that benthic copepods may serve as excellent indicators in biomonitoring programs of marine coastal ecosystems.  Data Availability Statement: The data are not to be shared due to restrictions, e.g., privacy and regulation.