Functional Trait Responses of Macrobenthos to Anthropogenic Pressure in Three Temperate Intertidal Communities

Abstract

With the increasing impact of human activities on marine ecosystems, there is a growing need to assess how the components of marine ecosystems (e.g., macrobenthos) respond to these anthropogenic pressures. In this study, the trait-based approach was used to assess the effects of anthropogenic pressures represented by the area of land-based aquaculture pond (PondArea) and heavy metals on the macrobenthic communities in three intertidal zones (AO, RW and XD) of Laoshan Bay, Shandong Peninsula, China. Compared with RW and XD, AO was under more pressure in terms of the average concentrations of heavy metals and TOC in sediments as well as PondArea. Fuzzy correspondence analysis (FCA) showed that there were significant differences in the composition of functional traits among the three regions (PERMANOVA; p < 0.05). In the highly polluted area, macrobenthic communities exhibited a combination of traits such as relatively short life span, weak mobility, feeding on deposits, and more tolerant to organic matter, while in a less polluted area, they exhibited a combination of traits such as relatively long life span, relatively high mobility, and more sensitive to organic matter. The RDA results showed that the distribution of the trait modalities was significantly affected by heavy metals (Hg and Cd), TOC, PondArea and sampled location. Variation partitioning analysis (VPA) indicated that the shared influence of sediment-related pollution factors and PondArea contributed most to the variance of the functional traits, which implied that human activities directly and/or indirectly lead to changes in functional traits of macrobenthic communities in the intertidal zones.

Introduction

The marine macrobenthos community is a critical component and reliable indicator of the biotic integrity of marine ecosystems, especially the intertidal ecosystems (Teixeira et al., 2009; Piló et al., 2016; Llanos et al., 2020). On the one hand, macrobenthos plays a vital role in maintaining ecosystem functions, such as material cycling in sediments and energy flow in food webs. On the other hand, macrobenthos is relatively sedentary and therefore reflects the ambient conditions of sediments, in which many pollutants (e.g., heavy metals and organic enrichment) are ultimately partitioned (Ryu et al., 2011; Desrosiers et al., 2019).

Heavy metal pollution is one of the most common anthropogenic pressures that impact marine ecosystems (e.g., intertidal zones, coastal waters, and estuaries), which has been documented by many studies throughout the world (Boudaya et al., 2019; Izegaegbe et al., 2020; Wang et al., 2020). Heavy metal contaminants can result in adverse toxic effects on benthic organisms (Dewitt et al., 1996; Dabney et al., 2018), leading to the changes in composition, structure, and ecosystem function of macrobenthic communities (Mucha et al., 2003; Piló et al., 2016; Hu et al., 2019; Roe et al., 2020; Dreujou et al., 2021). For example, in Ria de Aveiro coastal lagoon (Portugal), with the increase of mercury contamination, the total abundance and species richness decreased, and tolerant taxa increased (Nunes et al., 2008); in Incheon North Harbor (Korea) and coastal zone south of Sfax (Tunisia), macrobenthic community gradually changed with the pollution levels, and species diversity decreased with decreased distance from the pollution source (Ryu et al., 2011; Mosbahi et al., 2019). However, most of the studies were conducted in the subtidal zones other than intertidal zones, which are more vulnerable to human activities.

Macrobenthos consists of numerous taxa, and different species have a different tolerance to environmental pressures. For example, polychaetes Capitella capitata and Heteromastus filiformis are naturally tolerant to environmental disturbance, which could live well in a highly organic enrichment and/or heavy metal polluted area (Selck et al., 1999; Ryu et al., 2011; Bae et al., 2017), while some taxa (e.g., polychaete Magelona dakini and amphipods Perioculodeslongimanus) are inherently sensitive to environmental disturbance, and could not survive in such highly polluted zones (de-la-Ossa-Carretero et al., 2012; Ellis et al., 2017). It indicates that each species has evolved a unique survival strategy to adapt to different environmental conditions, even though it may be similar in some ways with other species. When facing loads of contaminants, such as metal(loid)s and organic enrichment or other contaminants gradients, macrobenthos have to make some reactions to resist such adverse environmental conditions. Therefore, macrobenthic responses may reflect different types and levels of pollutant impacts (Ryu et al., 2011; Piló et al., 2016).

There are two major approaches to evaluate the responses of macrobenthic communities to the effect of environmental variables and/or anthropogenic disturbances, namely, the taxonomic- and functional-based approach (Bremner et al., 2003; van der Linden et al., 2012; Henseler et al., 2019). The taxonomic-based approach estimates the impacts of environmental variables and/or anthropogenic pressures by comparing macrobenthos species assemblages in different levels of disturbance, whereas, the functional-based approach compares the macrobenthos functional traits assemblages under different disturbance levels. Although the taxonomic-based method is widely accepted and applied around the world, the functional-based approach has attracted more and more attention in recent years (Beauchard et al., 2017; Lam-Gordillo et al., 2020). This is mainly because that the functional-based method can better describe and understand the responses of macrobenthos community structure and function to multiple stressors than the taxonomic-based approach, which does not consider the changes in the structure of functional biological traits of macrobenthic communities (Bremner et al., 2003; Darr et al., 2014; Llanos et al., 2020).

Biological trait analysis (BTA) is the most applied functional-based approach, which is based on the assumption that phylogenetically independent organisms could have evolved into a similar biological adaptation, leading to functional similarity but taxonomic dissimilarity (Usseglio-Polatera et al., 2000; Dolédec and Statzner, 2008). The biological traits exhibited by organisms in a community provide information about how they behave and respond to stressors, thereby indicating the state of the environment (Bremner et al., 2003). Up to the present, macrobenthos functional traits has been successfully applied to assess the impacts of multiple anthropogenic stressors, such as heavy metals, polycyclic aromatic hydrocarbons (PAH), and sewage pollution (Piló et al., 2016; Bae et al., 2017; Dong et al., 2021; Nunes de Souza et al., 2021).

Land-based pond aquaculture in coastal areas is important for seafood production for many countries in the world. However, long-time land-based pond aquaculture could lead to the increase of heavy metal contents and organic enrichment in sediment (Kolarova and Napiórkowski, 2021). This has become a big challenge in those areas with a long aquaculture history. Therefore, pond aquaculture has also been criticized due to its impact on the environment and ecosystem. The effluent from ponds firstly affects the intertidal zone, with possible negative impacts on the intertidal organisms, and lead to the changes of the structures and functions of macrobenthos community. However, the effects of different pond aquaculture pressure levels on the structure and function of macrobenthic communities in the intertidal zone are less clear.

Around Laoshan Bay, Shandong Peninsula, China, numerous land-based aquaculture ponds cover an area of more than 20 km², and most ponds have been used for aquaculture for more than 10 years. The increase of heavy metal content in Laoshan Bay is mainly due to human activities (Wang et al., 2019). However, studies about the effects of long-time pond aquaculture and heavy metals or organic enrichment on the intertidal macrobenthic community are still limited, especially in China. Previous studies highlighted that the macrobenthos in intertidal zones appeared to be relatively sensitive to metal(loid)s stress than the subtidal zone (e.g., Roe et al., 2020). In addition, in the subtidal zone of Laoshan Bay, the functional approach has revealed that the heavy metals in the sediment resulted in changes in macrobenthos trait composition (Dong et al., 2021). To the best of our knowledge, the composition, structure, function, and drivers of macrobenthic communities in the intertidal zones of Laoshan Bay have not been studied yet. Taking the above situation into consideration, we hypothesized that: (1) the trait composition of macrobenthos communities is different in areas and reflects the strategies that the communities use to cope with the respective environmental stressors; (2) with an increase in anthropogenic pressure, the taxonomic and functional diversity of macrobenthic communities decline; and (3) the pressure of land-based pond aquaculture and sediment contamination (e.g., metal pollution and TOC) are the main drivers of the functional structure of macrobenthic communities.

Materials and Methods

Study Area

Laoshan Bay is the second-largest bay in Qingdao, located on the western coast of the Yellow Sea, with a surface area of approximately 188 km². Although six small rivers flow in this area, these are all seasonal rivers and are dry almost all year round. Even in the rainy season, only a small amount of water flows into Laoshan Bay, which is not enough to affect its salinity (31.80 ± 0.37 PSU) (Wang et al., 2019). Laoshan Bay is a traditional marine culture area, surrounded by many land-based aquaculture ponds. These aquaculture ponds are built not only on the soft bottom coast (mainly on the western and northern side of Laoshan Bay), but also on the hard bottom coast (mainly on the eastern side of Laoshan Bay). Many species, such as crabs, shrimps, and sea cumbers are cultured in these ponds. Long-term mariculture and rapid economic development since 2015 have led to an increase in the environmental pressures and human impacts in this bay (Wang et al., 2019; Dong et al., 2021). In recent years, although many ponds have been removed one after another due to environmental pollution and illegal land use, there are still considerable ponds, both in quantity and area. The impact of long-term pond aquaculture on the surrounding environment is expected to persist for a long time even after removal of the pond. To assess the effects of land-based pond aquaculture on the macrobenthic communities, three intertidal zones of Laoshan Bay representing different levels of pressure were selected (Figure 1). The intertidal zone of Aoshan Bay (AO) in the north of Laoshan Bay is the area under the greatest pressure of pond aquaculture. Within 5 km, there are 17.73 km² of ponds around AO. The distribution and area of these ponds have changed little since 2006. The intertidal zone of Xiaodao Bay (XD) is located in the west of Laoshan Bay, with 2.86 km² ponds around it within 5 km. The area and distribution pattern of these ponds have almost remained unchanged since 2012. The third area is located in the intertidal zone of the estuary where Wenquan River and Daren River meet (RW, in the west of Laoshan Bay). When we sampled, the rainy season had not yet arrived, and the two rivers were in a dry state. Within 5 km, there are only 0.20 km² of ponds around RW since 2012. However, before 2012, there used to be about 7 km² ponds here.

FIGURE 1

Figure 1. Study area with sampled sites in three intertidal zones in Laoshan Bay.

Biological Data

Macrobenthic samples were collected from three intertidal zones of Laoshan Bay at 21 stations of seven transects (two transects with six stations for RW and XD, three transects with nine stations for AO) (Figure 1) during low tide at spring tide, in May 2018, using PVC tube (diameter 10 cm, height 15 cm, 5 replicates per station). The five macrobenthic samples for each station were mixed, washed, and sieved using a 0.5 mm mesh, and then the macrobenthic organisms were picked and preserved in a 4% formaldehyde solution. In the laboratory, all benthic organisms were counted and identified to the lowest taxonomic level, as far as possible, according to species identification guides, using a stereomicroscope (Nikon SM2, Nikon Ltd., Japan). Taxon names were crosschecked against the website WoRMS—World Register of Marine Species http://www.marinespecies.org. Finally, the taxa abundance per station was compiled into the “taxa-by-station” matrix (L table, 50 species × 21 stations).

Environmental Data

Sediment samples were collected simultaneously with macrobenthic samples (one sample per station) using PVC tube (diameter 10 cm, height 10 cm), and then stored in the icebox and transported to the laboratory frozen at –20°C until analysis. The total organic carbon (TOC, %) in the sediment was determined using the potassium dichromate–sulfuric acid (K₂Cr₂O₇–H₂SO₄) oxidization method. Heavy metal (Cr, Pb, Zn, Cd, Cu, and Hg) and metalloid (As) in the sediment were measured following the “Specifications for marine monitoring Part 5: Sediment analysis” (State Bureau of Quality and Technical Supervision of China, 2008). Specifically, 0.5 g freeze-dried sediments were digested with a mixed acid solution (HF–HNO₃–HClO₄) before measurement. Then, the presence of heavy metals (Cd, Cr, Pb, Cu, Zn) was determined using a graphite furnace atomic absorption spectrometer (Shimadzu AA-6300C, Shimadzu Ltd., Japan), whereas the concentration of Hg and As were determined using an atomic fluorescence spectrometer (AFS-930, Titan Ltd., China). To ensure the accuracy of data quality, the blank and the standard reference material GBW07333 (Chinese national standard reference materials) were measured with the same method and procedure as described above for the sediment samples. The area of land-based aquaculture pond (Pond Area) within 5 km of three intertidal zones was measured in Google Earth. The relative location of sampled sites was divided into three groups: low, medium, and high. The low group was far away from the shore, and the high group was near the shore. The relative location of sampled sites was considered as a proxy of the spatial/natural variable (hereafter, referred to as a spatial variable).

Functional Trait

A total of eight traits (body form, adult size, mobility, position, feeding habit, AMBI groups, living habit, and life span) containing 33 categories (Table 1) for their potential ability to indicate environmental disturbance were selected to describe the important functional attributes of the macrobenthic community in the three intertidal zones studied. The traits of these species were collected from a variety of published sources [e.g., species identification guides, scientific papers, and established online databases, such as BIOTIC—Biology Traits Information Catalogue (MarLIN, 2006) and Polytraits (Polytraits Team, 2020)]. According to the affinity of species to the trait categories, they were scored with 0, 1, 2, and 3 using a fuzzy coding approach (Chevenet et al., 1994), with 0 indicating that the species has no affinity to a trait category and 3 indicating a high affinity to a trait category (van der Linden et al., 2016; Hu et al., 2019; D’Alessandro et al., 2020). If a particular trait data for individuals were not currently available, their scores were inferred from those species in the same genus or the same family, and in extreme cases, all the trait categories were scored to 0 (Marchini et al., 2008; Boyé et al., 2019; Dong et al., 2021). The species and their trait scores were then compiled into the “taxa-by-traits” matrix (Q table, 50 species × 33 categories of eight traits) for subsequent analysis.

TABLE 1

Table 1. Definition and labels of the 8 traits (33 categories) used to describe macrobenthos.

Data Analysis and Statistical Tests

Biological traits analysis requires three different matrices: (1) “taxa-by-station” matrix (L table, taxa abundance in each station), (2) “taxa-by-traits” matrix (Q table, fuzzy-coding scores for each trait category of each taxon), and (3) “station-by-traits” matrix (LQ table, 21 station × 33 trait categories), a combination of the L table and Q table. In this study, the value of each trait modality at each station was calculated by multiplying the affinity score for each modality of each species by the relative abundance of the species in each station and then summed over all species (Bremner et al., 2006; Marchini et al., 2008; Paganelli et al., 2012; van der Linden et al., 2012). The formula expresses as follows:

Scores=∑i=1SPi×Ti Scores=∑i=1SPi×Ti

where P_i is the relative abundance of the ith species in a given community, T_i is the score of ith trait modality of species, S is the total number of species in a station. The “station-by-trait” matrix used in this study is an abundance-weighted matrix, which is also called the community-level weighted means of trait values (CWM).

To identify the differences in functional composition among the three intertidal zone in Laoshan Bay, the “station-by-trait” matrix was ordinated based on Euclidean distances between samples by the fuzzy correspondence analysis (FCA), a correspondence analysis method appropriate for fuzzy coded data (Chevenet et al., 1994). Euclidean distances between samples were calculated from the relative frequencies of abundance-weighted traits in each station (Marchini et al., 2008; van der Linden et al., 2012, 2017). In FCA, the contribution of each trait category to total variation was reflected by correlation ratios, which are expressed as the amount of variance explained by each trait category along the first two FCA axes. The values of correlation ratios greater than 10% (RS > 0.1) can be considered as the traits with the most variable distribution (Conti et al., 2014). Therefore, in the two-dimensional FAC plot, samples close to one another have similar patterns of trait structure (Marchini et al., 2008; Paganelli et al., 2012; van der Linden et al., 2017). Subsequently, the cluster analysis based on the Euclidean distance of the taxon coordinates along the first four FCA axes and Ward’s linkage algorithm was applied to define the functional groups of macrobenthos with a similar combination of the trait modalities. To select the optimal number of groups, the “Phenon line” strategy was used following the previous studies (Usseglio-Polatera et al., 2000; Desrosiers et al., 2019). Taxa belonging to different functional groups are listed in Supplementary Table 1.

To test the significance of differences in the trait composition of different intertidal communities, a permutational multivariate analysis of variance (PERMANOVA) test was performed in the “vegan” package (Oksanen et al., 2019). The indicator traits in each community were determined by the indicator value (IndVal) using the function “multipatt” in the “indicspecies” package (De Cáceres and Legendre, 2009), and only when the statistical significance level of IndVal was <0.05, a trait modality was regarded as the indicator trait modality.

To examine relationships between macrobenthos functional traits and environmental factors, a redundancy analysis (RDA) based on “station-by-trait” matrix was performed as recommended by Kleyer et al. (2012), which showed that it is a more efficient method to assess the patterns of distribution of functional traits across environmental gradients. First, before RDA, the “station-by-trait” matrix was Hellinger-transformed, and the quantitative environment factors were log₁₀(x + 1) transformed and the highly correlated environmental variables (Spearman’s r > 0.7) were removed. Since TOC is an important indicator factor of sediment contamination levels and it has a significant relationship with the Pond Area (Table 2), the interaction between them was also considered in this work. Second, the RDA was constructed and its significance was tested using the function “anova” in the “vegan” package (Oksanen et al., 2019). Environmental variables with a variance inflation factor of more than 10 were excluded from the RDA model. In the final RDA model, only the variables presenting significant relationships (p < 0.05) after a stepwise selection were kept. In this work, we divided all the variables into three groups: sediment-related pollution factors [i.e., metal(loid)s and TOC], spatial factors (i.e., sampled locations), and Pond Area (i.e., direct human-induced factor). To explore the relative contribution of different variable groups, variation partitioning analysis (VPA) was performed using the “varpart” function in the “vegan” package (Oksanen et al., 2019). The adjusted R² of pure and shared contributions of these two data sets was reported because of their impartiality (Peres-Neto et al., 2006). The significance of the pure fraction of each data set was also tested in the “vegan” package using the “anvoa” function (Oksanen et al., 2019).

TABLE 2

Table 2. Spearman correlation between environmental factors.

The taxonomic indices: Shannon-Wiener (H′, log₂), Simpson (D), and Pielou (J′) indices were calculated using the “vegan” package (Oksanen et al., 2019) in R software 4.0.2 (R Core Team, 2020). The functional diversity indices: functional richness (FRic), functional evenness (FEve), and functional divergence (FDiv) were calculated using the “FD” package (Laliberté and Legendre, 2010). Functional redundancy index was estimated by the ratio of Rao’s quadratic entropy to Shannon-Wiener diversity (FRed), with an increase in the ratio indicating a decrease in functional redundancy and vice versa (van der Linden et al., 2012).

All statistical analysis was performed using R software 4.0.2 (R Core Team, 2020). The statistical difference in environmental parameters and diversity indices among different intertidal zones was analyzed using one-way ANOVA, followed by a post-hoc comparison with the Tukey HSD method. Spearman correlation tests were performed in the “corrplot” R package (Wei and Simko, 2017). In addition, the “ggplot2” R package (Wickham, 2016) was used. The map of the study area with the sampling sites was drawn using ArcGIS 10.4.2 (Esri, United States).

Results

Community Structure and Diversity

A total of 50 different taxa were identified in the intertidal zone of Laoshan Bay, of which polychaetes were the most abundant, with a total of 23 species, followed by mollusks (14 species) and crustaceans (10 species), and three species belonged to other taxonomic groups (Supplementary Table 1). Among the three subregions, AO had the highest species richness (35 species), followed by XD and RW with 29 and 20 species, respectively. IndVal method showed that the three intertidal zones have different indicator species. Nephtys californiensis, Anthopleura sp., and Sigalion sp. were the indicator species of AO, Macridiscus multifarius and Diogenes sp. were the indicator species of RW, and Glycera subaenea and Astarte borealis were the indicator species of XD. The Shannon-Wiener (H′), Simpson (D), and Pielou (J′) indices in AO were significantly higher than RW, whereas no significant difference between XD and RW or AO was found. There was no significant difference among the three functional diversity indices except for FDiv in RW and AO (RW higher than AO). In terms of FRed, RW has a higher value than AO and XD; however, there was no significant difference among them (Table 3).

TABLE 3

Table 3. Species and functional diversity indices in three intertidal communities.

Functional Trait Composition

Most of the sampled macrobenthos have a bivalved body form in three intertidal communities followed by vermiform body form (except for RW) (Figure 2). However, in the remaining three modalities of body form, RW has a relatively higher proportion than XD and AO, which showed that the species in RW has a more diversified body form. Most macrobenthos has a small adult size (Size1–5), whereas, the large adult size (size > 10) accounted for some proportion in AO. Burrower (M-Burr) and crawler (M-Craw) were the top two modalities of mobility. The infauna accounted for more proportion than epifauna across three intertidal communities. The macrobenthos showed a diversified feeding habit; suspension feeders (F-Susp) account for the most, followed by deposit feeders (F-Dept) or carnivores (F-Carn). Across the three intertidal zones, most macrobenthos were sensitive species (EG-I) followed by indifferent species (EG-II) and only a small proportion of species were tolerant (EG-III) or opportunistic species (EG-IV and EG-V). Free-living species (H-Free) and burrower dwelling species (H-Bdwel) were more than the tube- (H-Tube) or cave-dwelling (H-Cave) species in XD and RW, whereas, the cave-dwelling (H-Cave) species in AO were more than that in XD and RW. Macrobenthos was characterized by a long life span (Life > 5 year) in this work area, although they displayed relatively diversified modalities of life span.

FIGURE 2

Figure 2. Boxplot of the CWM (Community Weight Mean) values of each compiled functional trait calculated for each intertidal zone. For full labels of traits modalities, see Table 1.

Environmental Conditions

The mean concentration (±SD) of metal(loid)s and TOC in sediment are shown in Figure 3. In general, the mean values of the environmental factors (except Zn) in AO were higher than that in RW and XD, and the RW was slightly higher than XD. Specifically, Cr and TOC in AO were significantly higher than that in RW and XD, Cu and Hg in AO were significantly higher than that in XD, and Zn in RW was significantly higher than that in XD. No significant differences were found in the remaining environmental factors among the three subregions (AO, RW, and XD). Spearman correlation test showed that some metal(loid)s have a significant correlation with each other or with TOC or Pond Area, such as Cr with Hg, TOC, and Pond Area; Hg with Zn and TOC; as with TOC and Pond Area (Table 2). It indicated that these factors may have a similar source and implied that they may be related to the surrounding aquaculture pond.

FIGURE 3

Figure 3. Concentrations of heavy metal(loid) and TOC in sediments from three intertidal zones of Laoshan Bay. *p < 0.05; **p < 0.01; ***p < 0.001.

Fuzzy Correspondence Analysis

The results of FCA showed that the first two axes accounted for 63.32% of the total variability with 37.09% explained by axis 1 and 26.23% by axis 2 (Figure 4). PERMANOVA test revealed significant differences in trait composition among the three subregions communities (p < 0.05) (Table 4). The “IndVal” method revealed that three intertidal communities have distinctive indicator trait modalities (Figure 4). The XD was characterized by indicator trait modalities including infauna, borrower (M-Burr), and borrower dwelling (H-Bwel), whereas epifauna, omnivore, and shrimp/crab-liked body form were associated with RW. For AO, medium adult size (size 5–10), non/semimotile and cave-dwelling (H-Cave) were the indicator trait modalities.

FIGURE 4

Figure 4. Two-dimensional FCA plots weighted by abundance showing the distribution of sites in three intertidal zones of Laoshan Bay, and corresponding indicator trait modalities (p < 0.05) of the three subregions. For full labels of traits modalities, see Table 1.

TABLE 4

Table 4. PERMANOVA test for trait composition of macrobenthic communities in three intertidal zones.

The ordination plots for all traits and their modalities as well as the contribution of each trait to overall variability (i.e., correlation ratios) were shown in Figure 5. The trait AMBI groups, living habits, mobility, and position separated more on the first axis of FCA, whereas the trait body form, feeding habit, and life span separated both for the first and second axis of FCA, which indicated that those three traits correlated with both axes.

FIGURE 5

Figure 5. Results of the fuzzy correspondence analysis (FCA): distribution of the 33 modalities of the eight traits on the factor map. Numbers indicated the correlation ratio with axis 1 (horizontal) and axis 2 (vertical). For full labels of trait modalities, see Table 1.

Functional Groups of Macrobenthos

Macrobenthos sampled from three intertidal zones could be divided into four functional groups, which shared similar trait modalities according to Ward’s method and “Phenon line” strategy (Figure 6A). Functional group 1 (FG1) includes 25 taxa, most of which are polychaetes. FG2 consisted of nine taxa, which share the trait, such as bivalved body form, suspension-feeding habit, long life span, and sensitive to disturbance. FG3 (six taxa) consisted of epifaunal gastropod scavengers, with crawling mobility, free-living habit, long life span, and indifference to disturbance. FG4 included 10 taxa, which share shrimp or crab-liked body form, small adult size, and short life span.

FIGURE 6

Figure 6. Functional groups of macrobenthos. (A) Distribution of the taxa on the first two axes of FCA, taxa were clustered into four groups; (B) the relative abundance of four functional groups species in three intertidal communities. FGs, functional groups; For full names of species, see Supplementary Table 1.

The relative abundance of species from all the functional groups (FG1 to FG4) varied among three intertidal communities (Figure 6B). In AO, the macrobenthic communities were dominated by species from FG1 (51.39%) and FG2 (35.65%), whereas, in RW, the species from FG2 (56.64%) and FG4 (26.55%) were the most abundant. In XD, species from FG2 (66.12%) were the most abundant, followed by species from FG1 (27.27%), whereas the other functional groups (FG3 and FG4) were poorly represented.

Relationships Between Traits Modalities and Environmental Factors

The RDA model revealed the environmental effects on the trait modalities (Figure 7) and explained 43.25% of the total variation in macrobenthos trait modalities, with the first and second axis of RDA accounting for 21.96 and 13.79% of the total variation, respectively. These trait modalities could be divided into three groups (G1–G3) (Figure 7B), based on the cluster analysis of Euclidean distance of the coordinates of trait modalities along the first two axes of RDA (Figure 7C). Heavy metals (Hg and Cd), TOC, Pond Area, and sampled location had significant effects on the distribution of the trait modalities. G2 was located on the positive part of the second RDA axis including trait modalities, such as deposit-feeding habit (F-Dept), vermiform body form (Vermiform), tolerant and second-order opportunistic AMBI groups (EG-III and EG-IV), and cave- and tube-dwelling living habit (H-Cave and H-Tube) (Figure 7B), which are the characteristics of macrobenthos in relatively high polluted sites. It corresponded well with the highest concentrations of Hg, Cd, and TOC in sediment and also the highest area of aquaculture pond (Pond Area) in sites groups of AO. G1 and G3 located on the negative and positive parts of the first RDA axis matched with the site’s group of XD and RW, respectively, with relatively low anthropogenic pressure (e.g., heavy metals and Pond Area). It indicated that the traits modalities in G1 and G3 are more sensitive to anthropogenic pressures than those trait modalities in G2.

FIGURE 7

Figure 7. The results of redundancy analysis (RDA) of macrobenthic communities in three intertidal zones in Laoshan Bay. (A) Environmental variables vs. sites, (B) environmental variables vs. trait modalities, groups referred to panel (C), (C) groups of trait modalities based on the cluster analysis of Euclidean distance of the coordinates of trait modalities along the first two axes of RDA. For full labels of traits modalities, see Table 1.

Variation Partitioning Analysis

Variation partitioning analysis indicated that sediment-related pollution factors (7.93%) have higher pure total explained variance than spatial factors (5.7%) and Pond Area (i.e., direct human-induced factor) (1.14%), whereas, the shared influence of sediment-related pollution factors and Pond Area explained 16.65% of the variance (Figure 8). It can therefore be concluded that sediment-related pollution factors and direct human-induced factors had a greater impact on the distribution of the trait modalities than spatial factors.

FIGURE 8

Figure 8. The results of variation partitioning analysis (VPA) of the distribution of trait modalities.

Discussion

The BTA combined with FCA and variation partitioning technique was used in this study to reveal the structure and function of macrobenthic communities and identify the relationships between environmental (i.e., metals and TOC), spatial (i.e., sampled location), and human-induced (i.e., land-based ponds) predictors structuring benthic meta communities in intertidal systems. The results of this study suggested that the functional-based approach based on the trait-based matrix can efficiently evaluate the effects of anthropogenic pressure, represented by heavy metals and land-based ponds on the intertidal macrobenthic communities in Laoshan Bay, as previously shown in the intertidal zone (Llanos et al., 2020), subtidal zones (Dong et al., 2021), estuary (Piló et al., 2016), and lagoon (Hu et al., 2019). The present results supported our first hypothesis that the trait composition of macrobenthos is not randomly distributed across communities and can reflect their strategies to respond to environmental stressors.

Abiotic Variables

The results of this study showed that metal(loid)s (except Zn) and TOC concentrations in AO were higher than those in RW and XD, suggesting that the sediments in AO are more polluted than those in RW and XD. According to the guild of marine sediment quality standards of China (State Bureau of Quality and Technical Supervision of China, 2008), in the three intertidal zones, the concentration of Cr, Cd, and Hg exceeded the first grade at 28.75, 66.67, and 100% of sites, respectively, and Cd and Hg concentration exceeded the second grade at one (4.76%) and eight (38.10%) sites. In the subtidal zone of Laoshan Bay, these three metals in the sediment were also reported to exceed the first grade at 18.75–87.5% of the sites (Dong et al., 2021). In addition, several metals (e.g., Hg and Pb) in the seawater of Laoshan Bay were reported by previous studies higher than the primary grade of seawater quality standard of China (Wang et al., 2019). These indicated that Laoshan Bay is polluted by heavy metals, and the heavy metal elements mentioned above should attract our attention. In addition, all metal(loid)s elements (except Zn) and TOC showed positive correlations with the area of land-based cultured pond, which implied that metal(loid)s and TOC are probably related to land-based pond aquaculture (Kolarova and Napiórkowski, 2021). The runoff of land-based cultured pond could result in the accumulation of metal(loid)s and enrichment of organic matter in the surrounding and downstream sediments, which has been documented by many studies (e.g., Mendiguchía et al., 2006; Hong et al., 2020; Kolarova and Napiórkowski, 2021). However, there is no evidence of a direct causal relationship between sediment-related factors (i.e., heavy metals and TOC) and land-based pond area in this study. As reported by a previous study, these metals were also probably related to shipping, ship repair, and coal combustion (Wang et al., 2019).

Taxonomic and Functional Indices

The taxonomic and functional indices showed higher values in higher polluted area AO, except for FDiv, which was contrary to our second hypothesis: “higher disturbance and/or pollution associated with a lower taxonomic and functional diversities” Mouillot et al., 2013). This phenomenon not only appears in our work but also has been documented by previous works, such as in estuaries of the northeastern Brazilian coast and Mondego estuary located on the west central Atlantic coast of Portugal (van der Linden et al., 2016; Nunes de Souza et al., 2021). It can be explained by the intermediate disturbance hypothesis, which stated that species richness is assumed to peak for intermediate disturbance levels, indicating an increase in taxonomic and functional indices from no disturbed to moderately disturbed areas (Mouillot et al., 2013). In our work, we observed that the species richness exhibited a positive relationship with taxonomic and functional indices (except FDiv) (Supplementary Figure 1), which is consistent with the intermediate disturbance hypothesis. The relatively high FDiv value in RW indicated a high different combination of traits by the most abundant species here (Mason et al., 2005; van der Linden et al., 2016; Nunes de Souza et al., 2021), implying that the environment in this area may fluctuate greatly, which is supported by previous studies (e.g., Elliott and Quintino, 2007; van der Linden et al., 2012).

Redundancy plays a key role in maintaining the capacity of an ecosystem to respond to disturbances and changes (De Leo and Levin, 1997). The higher value of FRed represents a lower functional redundancy (van der Linden et al., 2012). Therefore, the FRed can also be considered a potential “early warning” indicator of pollution and/or disturbance (van der Linden et al., 2016; Nunes de Souza et al., 2021). As expected, the lowest value of FRed in XD indicated the relative higher functional redundancy, although no significant indifferences were found among the three intertidal zones (p > 0.05). It indicated that there are more species with similar traits to coexist here, and the disappearance of one or more species will not cause shifts in ecosystem functions, implying that the macrobenthic community in XD has a higher ability to resist the impact of disturbance. What we did not expect is that RW has the highest FRed value, that is, the lowest functional redundancy. One possible reason is that the impact of the long-time pond aquaculture on this area has not been completely eliminated, although these aquaculture ponds have been demolished for 6 years. The second probable reason is that RW is located in the estuarine area where environment factors fluctuate greatly, which may filter out some rare organisms, resulting in the decrease of species richness and diversities (Elliott and Quintino, 2007; van der Linden et al., 2012; Mouillot et al., 2013). FRed measures the amount of trait-dissimilarity among “species” and not among “individuals,” that is, all species are equally important. Therefore, when rare species with rare combinations disappear due to increasing disturbance, FRed will increase (van der Linden et al., 2016).

Traits and Trait-Based Functional Groups

The choice of traits and their modalities are affected by many factors, such as taxa, ecosystem types, and data limitation, and up to now, there is no standard reference method for choosing traits and their modalities in the current scientific community (Siddig et al., 2016; Beauchard et al., 2017). Eight traits were selected to define the functional structure of the macrobenthos in the intertidal zones of Laoshan Bay. These traits included life history, feeding strategy, habitat use, locomotion, and ecological preference, which are commonly used to assess the response of macrobenthic communities to anthropogenic pressures, such as heavy metals, oil spills, and aquaculture by previous studies (Piló et al., 2016; Egres et al., 2019; Lacson et al., 2019; Dong et al., 2021), although the modalities for some traits were different among studies. In general, the closer the phylogenetic relationship is for species, the more similar the traits are. In this work, the trait-based taxa functional groups were observed to be quite homogeneous in their taxonomic structure, suggesting that taxa in the same functional groups have closed phylogenetic relationships (Supplementary Table 1), which is consistent with the results obtained from previous works (Usseglio-Polatera et al., 2000; Desrosiers et al., 2019). Some studies implied that the convergence of traits with taxa phylogeny in some functional groups may have some impact on the trait-based methods (Menezes et al., 2010). On the other hand, some studies indicated that these confounding effects can be alleviated, and have only minor impacts on the functional response of macrobenthos disturbances (Statzner and Bêche, 2010). More studies are needed to explore this issue, especially in marine ecosystems.

Relationship Between Traits and Anthropogenic Pressures

Macrobenthos are composed of numerous organisms, displaying many biological traits, which exhibited different responses to different disturbances. The FCA revealed that the patterns in the spatial distribution of the trait modalities and trait-based taxa functional groups (FGs) were associated with each intertidal zone at different disturbance levels. It indicated that traits are not randomly distributed across macrobenthic communities and the macrobenthic responses can reflect anthropogenic pressures [e.g., metal(loid)s contamination and aquaculture] (Piló et al., 2016; D’Alessandro et al., 2020), confirming our first hypothesis.

In our work, the intertidal zone AO was identified as the most affected area. Pollutants that mainly affect this area are Cr, Cd, and Hg, as well as higher TOC and area of the land-based pond. The macrobenthos in this intertidal zone exhibited a combination of traits leaving them less vulnerable to such stressors as reported by previous studies (Piló et al., 2016; Desrosiers et al., 2019; Hu et al., 2019; Dong et al., 2021). For example, they have a very short life span; they feed on deposits and are more tolerant to organic matter; they have more opportunistic species; they have relatively weak mobility and prefer to live in caves or tubes, which provide protection for organisms from direct exposure to sediment commination and improves their survival possibilities (Hu et al., 2019). In a less polluted area, macrobenthos exhibited a combination of traits, such as a relatively long life span and are more sensitive to organic matter. In addition, they live freely and have relatively high mobility (e.g., swimming and crawling), which is beneficial for macrobenthos to expand their range of activities and increase their access to food, and even help them to escape from contaminated sediments quickly.

Nonetheless, not all trait categories responded to the anthropogenic pressures (e.g., metals contamination and aquaculture) as expected. The small and/or very small adult size was commonly considered as the prevalent traits in the higher polluted area (Piló et al., 2016), whereas, in our work, the medium and large adult size were identified as the dominant traits in the relative highly polluted area and, small and/or very small adult size were the dominant traits in the less polluted area. This phenomenon not only appears alone in this work, but also has been reported by many previous studies, such as in the coast of southern Portugal (Lacson et al., 2019), Cotinga sub-estuary of Paranaguá Bay, Brazil (Gusmao et al., 2016), and Swan Lagoon, China (Hu et al., 2019), as well as in the subtidal zone of Laoshan Bay (Dong et al., 2021). One reason is that the high metal polluted areas are usually associated with high organic matter or TOC in the sediments (Ryu et al., 2011; Hu et al., 2019), which is the important food source for deposit-feeders that are the dominant taxa in high contaminated area. The large food resources may increase the potential for growth of these taxa (Ryu et al., 2011). Secondly, the resistance of a species to disturbances (e.g., metal pollution and organic enrichment) is not determined by a single trait, but by the combined effect of multiple traits. The relative mismatch between the modalities of trait adult size and disturbance levels may be related to unexpected trade-offs among traits (Desrosiers et al., 2019).

As expected, the sediment-related pollution factors (e.g., metals and TOC) and the area of land-based pond aquaculture accounted for most variance of trait composition of macrobenthic communities, which confirmed that our third hypothesis, i.e., that they are the main drivers of the functional structure of macrobenthic communities. It also indicated that the macrobenthos is relatively sensitive to pond aquaculture and local heavy metal pollution, especially the contamination in sediment, as reported by previous studies (Ryu et al., 2011; Piló et al., 2016; Hu et al., 2019; Roe et al., 2020; Dong et al., 2021). However, it is difficult to further distinguish which pollution-related factors contributed more to the variance of the functional traits of macrobenthic communities (Piló et al., 2016; Hu et al., 2019; Dong et al., 2021). The accumulation of heavy metals and enrichment of organic matter in sediments commonly co-occurred with intensive land-based pond aquaculture (Kolarova and Napiórkowski, 2021). The intertidal zones are in high-energy environments and are affected by many factors and/or their combination. For more accuracy, at least parts of the trait composition of macrobenthic communities were caused by anthropogenic pressures, such as heavy metals and land-based pond aquaculture.

Conclusion

With the increase in the world’s population and the demand for aquatic products, it can be predicted that aquaculture ponds and/or facilities will be intensified in the future (Schubel and Thompson, 2019; Costello et al., 2020). This will further aggravate the impact of human activities on the marine ecosystem. Therefore, knowing the dynamics of the benthic communities subjected to the effect of anthropogenic disturbances is essential to assess the ecological status and protect the ecosystem. In this work, we investigated the effects of the anthropogenic disturbance represented by land-based aquaculture pond area (direct human impact) and sediment-related pollution factors (indirect human impact) on macrobenthic communities by using the trait-based approach in three intertidal zones of Laoshan Bay. The results confirmed that anthropogenic disturbance would lead to the change in functional traits of macrobenthic communities and highlighted that a combined effect of land-based pond area and sediment-related pollution factors on macrobenthic communities. Moreover, a combination of traits (e.g., relatively short life span, feeding on deposits, and more tolerant to organic matter) was identified in a highly disturbed area, which can be used as an early warning indicator of the ecosystem and those under similar disturbances. However, we should note that the intertidal ecosystems are in high-energy environments and are affected by many factors. Although the results of this work indicated that the Pond Area and sediment-related pollution factors, as well as a spatial factor were the drivers of macrobenthic communities, the effects of other factors (e.g., sediment grain size) on macrobenthic communities cannot be excluded. In addition, based on an independent field survey, the results may be biased to more accurately understand the impact of human stress macrobenthic communities, long-term, and multiple sampling surveys are recommended.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

Author Contributions

J-YD and XZ contributed to the conception and design of the study. J-YD and LZ gathered the data in the field and organized the database. J-YD and XS performed the statistical analysis. J-YD wrote the first draft of the manuscript. XZ and XY provided guidance and reviewed drafts of the manuscript. All authors contributed to the revision of the work.

Funding

This work was supported by funds from the National Key R&D Program of China (Grant No. 2019YFD0901303), the Province Key Research and Development Program of Zhejiang (Grant No. 2021C02047), and the National Natural Science Foundation of China (Grant No. 41906115).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank those who provided help during sample collection.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmars.2021.756814/full#supplementary-material

References

Bae, H., Lee, J. H., Song, S. J., Park, J., Kwon, B. O., Hong, S., et al. (2017). Impacts of environmental and anthropogenic stresses on macrozoobenthic communities in Jinhae Bay. Korea. Chemosphere 171, 681–691. doi: 10.1016/j.chemosphere.2016.12.112

Beauchard, O., Veríssimo, H., Queirós, A. M., and Herman, P. M. J. (2017). The use of multiple biological traits in marine community ecology and its potential in ecological indicator development. Ecol. Indic. 76, 81–96. doi: 10.1016/j.ecolind.2017.01.011

Boudaya, L., Mosbahi, N., Dauvin, J. C., and Neifar, L. (2019). Structure of the benthic macrofauna of an anthropogenic influenced area: Skhira bay (gulf of gabès, central mediterranean sea). Environ. Sci. Pollut. Res. 26, 13522–13538. doi: 10.1007/s11356-019-04809-8

Boyé, A., Thiebaut, E., Grall, J., Legendre, P., Broudin, C., Houbin, C., et al. (2019). Trait-based approach to monitoring marine benthic data along 500 km of coastline. Divers. Distrib. 25, 1879–1896. doi: 10.1111/ddi.12987

Bremner, J., Rogers, S. I., and Frid, C. L. J. (2003). Assessing functional diversity in marine benthic ecosystems: A comparison of approaches. Mar. Ecol. Prog. Ser. 254, 11–25. doi: 10.3354/meps254011

Bremner, J., Rogers, S. I., and Frid, C. L. J. (2006). Matching biological traits to environmental conditions in marine benthic ecosystems. J. Mar. Syst. 60, 302–316. doi: 10.1016/j.jmarsys.2006.02.004

Chevenet, F., Dolédec, S., and Chessel, D. (1994). A fuzzy coding approach for the analysis of long-term ecological data. Freshw. Biol. 31, 295–309. doi: 10.1111/j.1365-2427.1994.tb01742.x

Conti, L., Schmidt-Kloiber, A., Grenouillet, G., and Graf, W. (2014). A trait-based approach to assess the vulnerability of European aquatic insects to climate change. Hydrobiologia 721, 297–315. doi: 10.1007/s10750-013-1690-7

Costello, C., Cao, L., Gelcich, S., Cisneros-Mata, M., Free, C. M., Froehlich, H. E., et al. (2020). The future of food from the sea. Nature 588, 95–100. doi: 10.1038/s41586-020-2616-y

Dabney, B. L., Clements, W. H., Williamson, J. L., and Ranville, J. F. (2018). Influence of metal contamination and sediment deposition on benthic invertebrate colonization at the North Fork Clear Creek Superfund Site, Colorado, USA. Environ. Sci. Technol. 52, 7072–7080. doi: 10.1021/acs.est.7b06556

D’Alessandro, M., Porporato, E. M. D., Esposito, V., Giacobbe, S., Deidun, A., Nasi, F., et al. (2020). Common patterns of functional and biotic indices in response to multiple stressors in marine harbours ecosystems. Environ. Pollut. 259:113959. doi: 10.1016/j.envpol.2020.113959

Darr, A., Gogina, M., and Zettler, M. L. (2014). Functional changes in benthic communities along a salinity gradient- a western Baltic case study. J. Sea Res. 85, 315–324. doi: 10.1016/j.seares.2013.06.003

De Cáceres, M., and Legendre, P. (2009). Associations between species and groups of sites: indices and statistical inference. Ecology 90, 3566–3574. doi: 10.1890/08-1823.1

De Leo, G. A., and Levin, S. (1997). The multifaceted aspects of ecosystem integrity. Conserv. Ecol. 1:3.

de-la-Ossa-Carretero, J. A., Del-Pilar-Ruso, Y., Giménez-Casalduero, F., Sánchez-Lizaso, J. L., and Dauvin, J.-C. (2012). Sensitivity of amphipods to sewage pollution. Estuar. Coast. Shelf Sci. 96, 129–138. doi: 10.1016/j.ecss.2011.10.020

Desrosiers, M., Usseglio-Polatera, P., Archaimbault, V., Larras, F., Méthot, G., and Pinel-Alloul, B. (2019). Assessing anthropogenic pressure in the St. Lawrence River using traits of benthic macroinvertebrates. Sci. Total Environ. 649, 233–246. doi: 10.1016/j.scitotenv.2018.08.267

Dewitt, T. H., Swartz, R. C., Hansen, D. J., McGovern, D., and Berry, W. J. (1996). Bioavailability and chronic toxicity of cadmium in sediment to the estuarine amphipod Leptocheirus plumulosus. Environ. Toxicol. Chem. 15, 2095–2101. doi: 10.1002/etc.5620151205

Dolédec, S., and Statzner, B. (2008). Invertebrate traits for the biomonitoring of large European rivers: An assessment of specific types of human impact. Freshw. Biol. 53, 617–634. doi: 10.1111/j.1365-2427.2007.01924.x

Dong, J.-Y., Zhao, L., Sun, X., Hu, C., Wang, Y., Li, W. T., et al. (2021). Response of macrobenthic communities to heavy metal pollution in Laoshan Bay, China: A trait-based method. Mar. Pollut. Bull. 167:112292. doi: 10.1016/j.marpolbul.2021.112292

Dreujou, E., Desroy, N., Carrière, J., Tréau, de Coeli, L., McKindsey, C. W., et al. (2021). Determining the ecological status of benthic coastal communities: A case in an Anthropized Sub-Arctic Area. Front. Mar. Sci. 8:121–133. doi: 10.3389/fmars.2021.637546

Egres, A. G., Hatje, V., Gallucci, F., Machado, M. E., and Barros, F. (2019). Effects of an experimental oil spill on the structure and function of benthic assemblages with different history of exposure to oil perturbation. Mar. Environ. Res. 152:104822. doi: 10.1016/j.marenvres.2019.104822

Elliott, M., and Quintino, V. (2007). The estuarine quality paradox, environmental homeostasis and the difficulty of detecting anthropogenic stress in naturally stressed areas. Mar. Pollut. Bull. 54, 640–645. doi: 10.1016/j.marpolbul.2007.02.003

Ellis, J. I., Clark, D., Atalah, J., Jiang, W., Taiapa, C., Patterson, M., et al. (2017). Multiple stressor effects on marine infauna: Responses of estuarine taxa and functional traits to sedimentation, nutrient and metal loading. Sci. Rep. 7, 1–16. doi: 10.1038/s41598-017-12323-5

Gusmao, J. B., Brauko, K. M., Eriksson, B. K., and Lana, P. C. (2016). Functional diversity of macrobenthic assemblages decreases in response to sewage discharges. Ecol. Indic. 66, 65–75. doi: 10.1016/j.ecolind.2016.01.003

Henseler, C., Nordström, M. C., Törnroos, A., Snickars, M., Pecuchet, L., Lindegren, M., et al. (2019). Coastal habitats and their importance for the diversity of benthic communities: A species- and trait-based approach. Estuar. Coast. Shelf Sci. 226:106272. doi: 10.1016/j.ecss.2019.106272

Hong, A. H., Hargan, K. E., Williams, B., Nuangsaeng, B., Siriwong, S., Tassawad, P., et al. (2020). Examining molluscs as bioindicators of shrimp aquaculture effluent contamination in a southeast Asian mangrove. Ecol. Indic. 115:106365. doi: 10.1016/j.ecolind.2020.106365

Hu, C., Dong, J., Gao, L., Yang, X., Wang, Z., and Zhang, X. (2019). Macrobenthos functional trait responses to heavy metal pollution gradients in a temperate lagoon. Environ. Pollut. 253, 1107–1116. doi: 10.1016/j.envpol.2019.06.117

Izegaegbe, J. I., Vivier, L., and Mzimela, H. M. (2020). Trace metal contamination in sediment in the Mhlathuze Estuary, northern KwaZulu-Natal, South Africa: effects on the macrobenthic community. Environ. Monit. Assess. 192:401. doi: 10.1007/s10661-020-08352-9

Kleyer, M., Dray, S., Bello, F., Lepš, J., Pakeman, R. J., Strauss, B., et al. (2012). Assessing species and community functional responses to environmental gradients: which multivariate methods? J. Veg. Sci. 23, 805–821. doi: 10.1111/j.1654-1103.2012.01402.x

Kolarova, N., and Napiórkowski, P. (2021). Trace elements in aquatic environment. Origin, distribution, assessment and toxicity effect for the aquatic biota. Ecohydrol. Hydrobiol. 2021:2. doi: 10.1016/j.ecohyd.2021.02.002

Lacson, A. Z., Piló, D., Pereira, F., Carvalho, A. N., Cúrdia, J., Caetano, M., et al. (2019). A multimetric approach to evaluate offshore mussel aquaculture effects on the taxonomical and functional diversity of macrobenthic communities. Mar. Environ. Res. 151:104774. doi: 10.1016/j.marenvres.2019.104774

Laliberté, E., and Legendre, P. (2010). A distance-based framework for measuring functional diversity from multiple traits. Ecology 91, 299–305. doi: 10.1890/08-2244.1

Lam-Gordillo, O., Baring, R., and Dittmann, S. (2020). Ecosystem functioning and functional approaches on marine macrobenthic fauna: a research synthesis towards a global consensus. Ecol. Indic. 115:106379. doi: 10.1016/j.ecolind.2020.106379

Llanos, E. N., Saracho Bottero, M. A., Jaubet, M. L., Elías, R., and Garaffo, G. V. (2020). Functional diversity in the intertidal macrobenthic community at sewage-affected shores from Southwestern Atlantic. Mar. Pollut. Bull. 157:111365. doi: 10.1016/j.marpolbul.2020.111365

Marchini, A., Munari, C., and Mistri, M. (2008). Functions and ecological status of eight Italian lagoons examined using biological traits analysis (BTA). Mar. Pollut. Bull. 56, 1076–1085. doi: 10.1016/j.marpolbul.2008.03.027

MarLIN (2006). “BIOTIC - Biological Traits Information Catalogue,” in Marine Life Information Network (Plymouth: Marine Biological Association of the United Kingdom). Available online at: www.marlin.ac.uk/biotic (accessed May 13, 2021).

Mason, N. W. H., Mouillot, D., Lee, W. G., and Wilson, J. B. (2005). Functional richness, functional evenness and functional divergence: The primary components of functional diversity. Oikos 111, 112–118. doi: 10.1111/j.0030-1299.2005.13886.x

Mendiguchía, C., Moreno, C., Mánuel-Vez, M. P., and García-Vargas, M. (2006). Preliminary investigation on the enrichment of heavy metals in marine sediments originated from intensive aquaculture effluents. Aquaculture 254, 317–325. doi: 10.1016/j.aquaculture.2005.10.049

Menezes, S., Baird, D. J., and Soares, A. M. V. M. (2010). Beyond taxonomy: a review of macroinvertebrate trait-based community descriptors as tools for freshwater biomonitoring. J. Appl. Ecol. 47, 711–719. doi: 10.1111/j.1365-2664.2010.01819.x

Mosbahi, N., Serbaji, M. M., Pezy, J. P., Neifar, L., and Dauvin, J. C. (2019). Response of benthic macrofauna to multiple anthropogenic pressures in the shallow coastal zone south of Sfax (Tunisia, central Mediterranean Sea). Environ. Pollut. 253, 474–487. doi: 10.1016/j.envpol.2019.06.080

Mouillot, D., Graham, N. A. J., Villéger, S., Mason, N. W. H., and Bellwood, D. R. (2013). A functional approach reveals community responses to disturbances. Trends Ecol. Evol. 28, 167–177. doi: 10.1016/j.tree.2012.10.004

Mucha, A. P., Vasconcelos, M. T. S., and Bordalo, A. A. (2003). Macrobenthic community in the Douro estuary: relations with trace metals and natural sediment characteristics. Environ. Pollut. 121, 169–180. doi: 10.1016/S0269-7491(02)00229-4

Nunes, M., Coelho, J. P., Cardoso, P. G., Pereira, M. E., Duarte, A. C., and Pardal, M. A. (2008). The macrobenthic community along a mercury contamination in a temperate estuarine system (Ria de Aveiro. Portugal). Sci. Total Environ. 405, 186–194. doi: 10.1016/j.scitotenv.2008.07.009

Nunes de Souza, J. V., Martins Vieira, M. L., De Assis, J. E., Zanardi - Lamardo, E., Gomes, P. B., and Botelho de Souza, J. R. (2021). Responses of functional traits of macrobenthic communities on the presence of Polycyclic Aromatic Hydrocarbons contamination in three tropical estuaries. Estuar. Coast. Shelf Sci. 250:107105. doi: 10.1016/j.ecss.2020.107105

Oksanen, J., Blanchet, F. G., Friendly, M., Kindt, R., Legendre, P., McGlinn, D., et al. (2019). vegan: Community Ecology Package. Available online at: https://CRAN.R-project.org/package=vegan (accessed date 2020-11-28).

Paganelli, D., Marchini, A., and Occhipinti-Ambrogi, A. (2012). Functional structure of marine benthic assemblages using Biological Traits Analysis (BTA): A study along the Emilia-Romagna coastline (Italy, North-West Adriatic Sea). Estuar. Coast. Shelf Sci. 96, 245–256. doi: 10.1016/j.ecss.2011.11.014

Peres-Neto, P. R., Legendre, P., Dray, S., and Borcard, D. (2006). Variation partitioning of species data matrices: estimation and comparison of fractions. Ecology 87, 2614–2625. doi: 10.1890/0012-9658200687[2614:vposdm]2.0.co;2

Piló, D., Ben-Hamadou, R., Pereira, F., Carriço, A., Pereira, P., Corzo, A., et al. (2016). How functional traits of estuarine macrobenthic assemblages respond to metal contamination? Ecol. Indic. 71, 645–659. doi: 10.1016/j.ecolind.2016.07.019

Polytraits Team (2020). “Polytraits: A database on biological traits of polychaetes. Lifewatch greece,” in Hellenic Centre for Marine Research Available online at: http://polytraits.lifewatchgreece.eu (aaccessed May 13, 2021).

R Core Team (2020). R: A Language and Environment for Statistical Computing. Vienna: R Foundation for Statistical Computing.

Roe, R. A. L., Tran, T. K. A., Schreider, M. J., and MacFarlane, G. R. (2020). Assessment of the effects of sediment-associated metals and metalloids on mangrove macroinvertebrate assemblages. Water. Air. Soil Pollut. 231:352. doi: 10.1007/s11270-020-04731-7

Ryu, J., Seong, J., Kang, S., Kang, D., Lee, C., and Koh, C. (2011). The impact of heavy metal pollution gradients in sediments on benthic macrofauna at population and community levels. Environ. Pollut. 159, 2622–2629. doi: 10.1016/j.envpol.2011.05.034

Schubel, J. R., and Thompson, K. (2019). Farming the sea: The only way to meet humanity’s future food needs. GeoHealth 3, 238–244. doi: 10.1029/2019GH000204

Selck, H., Decho, A. W., and Forbes, V. E. (1999). Effects of chronic metal exposure and sediment organic matter on digestive absorption efficiency of cadmium by the deposit-feeding polychaete Capitella species I. Environ. Toxicol. Chem. 18, 1289–1297. doi: 10.1002/etc.5620180631

Siddig, A. A. H., Ellison, A. M., Ochs, A., Villar-Leeman, C., and Lau, M. K. (2016). How do ecologists select and use indicator species to monitor ecological change? Insights from 14 years of publication in Ecological Indicators. Ecol. Indic. 60, 223–230. doi: 10.1016/j.ecolind.2015.06.036

State Bureau of Quality and Technical Supervision of China (2008). “Part 5: Sediment Analysis” National Standard of China for the Specification for Marine Monitoring (GB 17378.5-2007). Beijing: Standards Press of China.

Statzner, B., and Bêche, L. A. (2010). Can biological invertebrate traits resolve effects of multiple stressors on running water ecosystems? Freshw. Biol. 55, 80–119. doi: 10.1111/j.1365-2427.2009.02369.x

Teixeira, H., Neto, J. M., Patrício, J., Veríssimo, H., Pinto, R., Salas, F., et al. (2009). Quality assessment of benthic macroinvertebrates under the scope of WFD using BAT, the benthic assessment tool. Mar. Pollut. Bull. 58, 1477–1486. doi: 10.1016/j.marpolbul.2009.06.006

Usseglio-Polatera, P., Bournaud, M., Richouxy, P., and Tachety, H. (2000). Biological and ecological traits of benthic freshwater macroinvertebrates: relationships and definition of groups with similar traits. Freshw. Biol. 43, 175–205. doi: 10.1046/j.1365-2427.2000.00535.x

van der Linden, P., Marchini, A., Dolbeth, M., Patrício, J., Veríssimo, H., and Marques, J. C. (2016). The performance of trait-based indices in an estuarine environment. Ecol. Indic. 61, 378–389. doi: 10.1016/j.ecolind.2015.09.039

van der Linden, P., Marchini, A., Smith, C. J., Dolbeth, M., Simone, L. R. L., Marques, J. C., et al. (2017). Functional changes in polychaete and mollusc communities in two tropical estuaries. Estuar. Coast. Shelf Sci. 187, 62–73. doi: 10.1016/j.ecss.2016.12.019

van der Linden, P., Patrício, J., Marchini, A., Cid, N., Neto, J. M., and Marques, J. C. (2012). A biological trait approach to assess the functional composition of subtidal benthic communities in an estuarine ecosystem. Ecol. Indic. 20, 121–133. doi: 10.1016/j.ecolind.2012.02.004

Wang, X., Fu, R., Li, H., Zhang, Y., Lu, M., Xiao, K., et al. (2020). Heavy metal contamination in surface sediments: A comprehensive, large-scale evaluation for the Bohai Sea. China. Environ. Pollut. 260:113986. doi: 10.1016/j.envpol.2020.113986

Wang, X., Liu, L., Zhao, L., Xu, H., and Zhang, X. (2019). Assessment of dissolved heavy metals in the Laoshan Bay. China. Mar. Pollut. Bull. 149:110608. doi: 10.1016/j.marpolbul.2019.110608

Wei, T., and Simko, V. (2017). R Package “Corrplot”: Visualization of a Correlation Matrix. Version 0.84. Available online at: https://github.com/taiyun/corrplot

Wickham, H. (2016). ggplot2: Elegant Graphics for Data Analysis. New York, NY: Springer-Verlag.