Study area

Four sites in the GC northwest region of Mexico were studied. Three of the sites share similar arid environmental conditions, with high evaporation rates (>1,750 mm·yr-1) that greatly exceed precipitation (^{López-Medellín and Ezcurra 2012}, Table 1). Additionally, they lack permanent freshwater contributions from rivers, and the mangrove trees consequently tend to be dwarfed; the water depth in the coastal mangrove fringe is relatively shallow, with a semi-mixed tidal regime. The most common mangrove species are Rhizophora mangle, Avicennia germinans, and Laguncularia racemosa (^{Payan-Alejo 2012}, ^{Padilla-Serrato et al. 2016}).

Bahía de La Paz (Baja California Sur, BCS) is the largest protected coastal system on the eastern side of the Baja California Peninsula. It is located in the southern region of the Central GC (sensu^{Brusca et al. 2005}) between 24º07′ and 24º21′ N and 110º17′ and 110º40′ W. Bahía de la Paz has several mangrove patches, with different substrate conditions. The bay has a surface area of 45 km2 and a maximum depth of 450 m, although the depth of the selected mangroves is shallow (2 m). Annual rainfall is around 200 mm. Santa María-La Reforma (Ramsar designation date 02 February 2009, No. 2025) is found on the western coast of the state of Sinaloa (SIN) in the region of the Southern GC. It is located between 24º43′ and 25º15′ N and 107º55′ and 108º26′ W, with a surface area of 40 km2 and depths ranging from 0.5 to 32.0 m. This location is protected by sandy barriers formed by the action of the ocean currents, and annual rainfall there is 650 mm. Las Guásimas (Ramsar designation date 02 February 2008, No. 1790) is a coastal lagoon located on the eastern coast of the GC in the southern region of the state of Sonora (SON) between 27º49′ and 27º55′ N and 110º29′ and 110º45′ W. This is a shallow lagoon (0.7 m average depth) with a surface area of 37 km2 and it is protected by 2 sand barriers that separate the water body from the sea, with 250 mm of annual rainfall (Fig. 1).

Figure 1 Map of the locations within the Gulf of California. Baja California Sur (BCS), Sinaloa (SIN), Sonora (SON), and Nayarit (NAY). 

The fourth site is an estuarine location and was selected to test the first hypothesis. The differences between the 3 target arid localities and this site provide the opportunity to compare if the functional structure of fish is driven by more elusive conditions in tropical environments (e.g., greater presence of freshwater, lower evaporation rates), and how this could change the functional and taxonomic diversity of the fish community. Marismas Nacionales is one of the most important mangrove sites on the Pacific coast of Mexico because of its diversity and productivity (Ramsar designation date 22 June 1995, No. 732; ^{CONABIO 2009}). It is located to the south and outside of the GC region in the state of Nayarit (NAY), between 21º30′ and 23º51′ N and 105º14′ and 106º01′ W. This is a subhumid site with estuarine conditions that receives a supply of freshwater from different rivers year-round, and an annual rainfall of 1,500 mm.

Species composition

In BCS, fish composition was determined from monthly field samples (2010-2011) in 4 sites within 2 localities using a beach seine (50 × 3 m; mesh size: 1 cm). At each sampling site, 2 hauls were collected during high tide. For the other localities, species lists were extracted from the scientific lit erature for SIN (^{Payan-Alejo 2012}), SON (^{Padilla-Serrato et al. 2016}), and NAY (^{CONABIO 2009}). All surveys from the literature followed systematic seasonal sampling lasting a yearlong and used seine nets and in some cases fish trawling to increase species detection effort. All extractions were done near the mangroves, an aspect that limits the maximum depth at which sampling can be performed. In the end, the overall sampling effort portrayed the ichthyofauna of each site in a similar way, allowing for analysis based on the presence and absence of species per location.

Analysis of similarity (ANOSIM) was used to test for differences in species composition similarity among sites and between sites (BCS, SIN, SON, and NAY). ANOSIM is a nonparametric test that uses a rank dissimilarity (Jaccard distance), where random iterations of the relationship between 2 objects are calculated and the resulting values range from -1 to 1; values closer to 1 indicate statistical difference between groups, while values closer to 0 indicate that the relationships are no different from a random construction. We separated 2 groups, the arid locations in one (BCS, SIN, and SON) and NAY in the second. Six random re-sam ples from NAY were required in order to generate enough data to create the second group, from which we ran 999 permutations of the analysis (^{Clarke 1993}).

Functional data

Representative photographs of each species were obtained from specialized fish webpages (e.g., FishBase, https://www.fishbase.de/; Naturalista, https://www.naturalista.mx/). Five species were not considered in the analysis because available photographs did not match the quality criteria to be selected. Using ImageJ software (National Institutes of Health, Bethesda, Maryland, USA; http://imagej.nih.gov/ij/), 19 homologous points were digitized in each photo (Fig. 2). From such points, 13 measurements (heights and lengths) were calculated for each species (Table 2). Finally, a total of 11 functional traits were obtained, which were derived from the proportions between measurements (Table 3). The traits were selected from those proposed by ^{Villéger et al. (2010)} and ^{Soares et al. (2013)} to best represent the strategies of fish species in relation to food acquisition, habitat, and locomotion. Correlation tables were used to identify and exclude those variables that were highly correlated with one another (Pearson coefficients >0.80), to avoid overestimation of the process that is represented in the functional attribute. From the remaining variables, a principal component analysis was used to identify those that were most significant in explaining the trait variation in the regional species pool (^{Houseman and Gross 2011}). As long as trait selec tion is representative of species performance (^{Villéger et al. 2017}), functional diversity indices are more reliable measures describing communities using even a small number of samples (^{Bady et al. 2005}).

Figure 2 Set of homologous points measured in every species. 

Some authors have suggested modifying the measurement procedure to reduce the variability in the measurements of the most extreme morphologies (e.g., flatfishes) (^{Toussaint et al. 2016}). However, we decided not to do so since it is important to consider differences in morphology for a better understanding of the ecological processes studied (^{Toussaint et al. 2016}).

Functional space

To create a multidimensional functional space from the traits selected, an affinity matrix using the Euclidean distance was calculated. Then, with the resulting matrix, a principal coordinate analysis (PCoA) was performed and results from the species coordinates were interpreted for the first 5 axes (^{Buisson et al. 2012}).

Functional indices

The functionality of the fish community was characterized using the functional richness (FRic, ^{Villéger et al. 2010}), functional dispersion (FDis, ^{Laliberté and Legendre 2010}), and functional redundancy (FRed, ^{Bello et al. 2007}). Using the R statistical language (v.3.6.1, R Core Team 2019), the first 2 indices were calculated with the multidimFD package (^{Villéger et al. 2017}) and functional redundancy was cal culated with the SYNCSA package (^{Debastiani and Pillar 2012}).

To test our first hypothesis, functional richness was used to refer to the ecological niche occupied by species within a community (^{Mason et al. 2005}). Functional dispersion is the weighted mean distance in the multidimensional trait space of individual species to the centroids (^{Laliberté and Legendre 2010}), and it works as a complementary analysis to understand the distance of the traits from the center of the functional space. Both indices allowed us to create the base line of the traits present in the community and to understand how the assembly process could be shaped.

Environmental and biotic interactions

With the aim of relating these traits to environmental conditions such as salinity, maximum and minimum temperature, type of substrate, precipitation rate, evaporation rate, precipitation deficit (precipitation rate - evaporation rate; mm·yr-1), and mangrove coverage area, distance-based redundancy analysis with Jaccard distance was used. The analysis uses PCoA eigenvalues into a redundancy analysis, allowing us to identify the relationship between multiple variables x and y in a non-symmetrical way. The analysis was processed with the function capscale of the package vegan in R (v.3.6.1, R Core Team 2019). Additionally, we described the environmental characteristics of each site from scientific publications and specialized fish webpages to corroborate this relationship (e.g., ^{Castro-Aguirre et al. 1999}; FishBase, https://www.fishbase.de/).

Specialized traits

In order to verify our second hypothesis and detect the specialized traits in the communities, functional redundancy analysis was selected to represent the maximum saturation of species with similar traits that are found in a given community; a zero value indicated that all the species are functionally different, in contrast to a value of one, which indicated that all the species shared similar traits in the community (^{Ricotta et al. 2016}).

To complement the detection of the specialized traits, a taxonomic distinctness analysis was used to add another perspective into the differences in the community structure. As less taxa are found, it is most likely that the community would be represented by dominant traits, pointing out to a specialized community. This technique calculates the average distance based on the species phylogenetic tree of the locality, pairing all the species together (^{Clarke and Warwick 1999}) based on their taxonomic relationships (species, genus, family, order, and class). The nature of this technique has proven effective in overcoming the limitation of comparing sites with different sampling efforts, by a res ampling process of the “regional” pool of the species and the mathematical pairing of the analysis. An a posteriori test was used to detect the statistical similarity in the taxo nomic structure of the sites, through the respective indices calculated from random subsets of communities to create a probability funnel, in which the localities are positioned; the ones outside the funnel are considered statistically different (^{Clarke and Warwick 1999}). The analyses were performed using PRIMER software (v.6, ^{Clarke and Gorley 2006}).

Species richness

A total of 123 species were identified in the arid mangrove habitats in the GC region in northwest Mexico. Specifically, 69 species were found in Las Guásimas (SON), 58 in Santa María-La Reforma (SIN), and 50 in La Paz Bay (BCS). We found 12 shared species within 9 families across these localities: Haemulidae (3 spp.), Lutjanidae (2 spp.), Ephippidae, Gerreidae, Gobiidae, Mugilidae, Paralichthyidae, Synodontidae, and Tetraodontidae (one species each). For each pair of localities in the gulf, at least 8 to 9 species were shared (SIN-SON, SON-BCS, and BCS-SIN). Marismas Nacionales (NAY) is represented by 111 species, with 10 species in common with the species pool of the GC sites. It shares species from 8 out of 9 families with the exception of Gobiidae, with 59 exclusive species. ANOSIM detected differences between the arid mangrove sites and the external community (R = 0.951, P = 0.007).

Functionality

The variables that defined the grouping of the functional space in the community were the position of the mouth, aspect of the caudal peduncle, and position of the eye. The explained variance for the 4 localities was 95% in the first 5 ordination axes.

The results pertinent to the first hypothesis were that func tional richness and dispersion presented the highest values in SON because of a higher diversity of species traits related to the presence of flatfishes in the community (Fig. 3). There was a large clustering of species from the order Perciformes in the middle of the functional space for all localities, where the attributes tended to be redundant because of their similar morphologies (perch-shape). Outside the cluster of redundancy, the nearest points are related to benthic species (Gobiidae, Tetraodontidae), while the farthest points are related to the morphology of flatfishes. The NAY locality, on the other hand, presented the lowest functional values for richness and dispersion because of its estuarine condition.

Figure 3 Principal coordinate analysis for functional dispersion (FDis) and functional richness (FRic) in each locality: Baja California Sur (BCS), Sinaloa (SIN), Sonora (SON), and Nayarit (NAY). The functional space of the mangrove fishes in the Gulf of California represents the species distributed along the first and second axes using 11 functional traits. The clear polygon represents the functional space of the total pool of species at all localities, while the smaller polygon (shaded) represents the functional space at each locality. The explained variability in the axis is 60.0% for BCS, 57.9% for NAY, 52.3% for SON, and 49.3% for SIN. In general, the northern localities had higher values in both indices than the southern locality (NAY). 

The distance-based redundancy analysis detected 3 separate groups among the sites analyzed and the environmental factors. In the first axis, there was a salinity gradient with the arid locations closely associated with higher values, BCS and SIN presented similar salinity concentrations, and NAY had lower salinity. The variable related to the surface area describes how NAY (with the highest richness) was related to a larger land extension of mangrove vegetation. Finally, in SON the community was governed by higher evaporation rates (Fig. 4).

Figure 4 Distribution-based redundancy analysis ordination diagram of northwest region with environmental variables and species scores for Baja California Sur (BCS), Sinaloa (SIN), Sonora (SON), and Nayarit (NAY). 

In the description of the milieu of overall species (179 species in the 4 localities), 56% were defined as strictly marine, 32% had a marine-estuarine affinity, 10% were present from marine to freshwater habitats (e.g., Achirus mazatlanus, Eucinostomus currani, Gerres cinereus, Mugil curema, and Mugil cephalus), and 2% were peripheral with continuous migration from the sea into the river (e.g., Colpichthys hubbsi and Colpichthys regis). Only one introduced species of freshwater origin was found in NAY (Oreochromis aureus). Regarding the species with wide salinity tolerance, 6, 8, 11, and 16 species were found in the SON, BCS, SIN, and NAY sites, respectively. In the case of the species that occurred at all sites (11 spp.), those with marine-estuarine affinity were the most common (5 spp.), followed by those that were strictly marine (3 spp.), and finally those with an affinity to all environments (3 spp.).

The results related to the second hypothesis showed that functional redundancy had similar results for all localities with overlap in the community’s functional space, with an increasing gradient starting from the northern to the southern localities (SON: 0.586; BCS: 0.595; SIN: 0.600; and NAY: 0.641). The taxonomic diversity showed significant statis tical variation in the taxonomic composition of the localities within the GC because SIN was outside the probabilistic funnel. However, if we compare the position of NAY, there is a higher level of variance within this location (Fig. 5). Additionally, this analysis showed that the organizational levels with the greatest probability of change were found at the tax onomic level of order.

Figure 5 The probability funnel of the taxonomy shows the northwest localities of Mexico, and the ones that are inside the funnel had statistically similar taxonomic composition. Baja California Sur (BCS), Sinaloa (SIN), Sonora (SON), and Nayarit (NAY). The values with respect to the y-axis (delta +) show the taxon level of variation of the community (order level). 

Composition of fish communities

A comparison of the diversity of fish species in mangroves from the GC with that of other mangroves with arid condi tions worldwide (e.g., ^{Shahraki et al. 2016}) indicated that our localities (>55 spp. on average) presented almost twice the richness of species than that at other sites (±30 spp.). A possible explanation is that the GC presents temporal and spatial changes in the climatic and environmental conditions, promoting the annual and interannual signals that enable the entrance of fish species into the gulf (^{Escalante et al. 2013}). Consequently, it is common to find species with tropical and temperate affinities (^{González-Acosta et al. 2005}).

On the other hand, when comparing the total richness of the GC sites against NAY (111 spp.), we found reduced richness of species in northwest Mexico. For example, in NAY, the presence of marine and brackish species in the families Centropomidae, Eleotridae, and Sciaenidae was common, as was the presence of some freshwater euryhaline species like the those in the family Cichlidae. In the case of marine demersal fish, some species of Haemuidae (Anisotremus interruptus and Genyatremus dovii), Lutjanidae (Lutjanus inermis and Lutjanus jordani), and Carangidae (Carangoides otrynter) were absent in the northwest.

The functional approach

The first 2 axes in the functional space of the mangrove fish community explained a high percentage of the overall variance since almost all the species belonged to Perciformes. The traits with the higher functional explanation power were the position of the mouth and position of the eye on the first axis, both related to the feeding process. Meanwhile, the aspect of the caudal peduncle better described the second axis and is related to locomotion. Similar results were found by ^{Soares et al. (2013)} for fish communities in Brazil. However, the results of the explained variance appear to be magnified by the percentage of species present from the families Lutjanidae, Haemulidae, Sciaenidae, Serranidae, and Gerreidae (fishes with similar morphology, from now on referred as generalists). This has been reported for other fish fauna that have experienced high radiation events but had less significant differences in morphological features (^{Toussaint et al. 2016}).

Other fish adaptations in the results showed sessile ben thic species (flatfishes, the farthest species from the cluster of generalist species) on the first axis. These species are characterized by relatively thick peduncles, indicating the use of short burst movements. Further, the mouth and eyes are higher with respect to the size of the head (although this aspect is more related to an ontogenetic process than a feeding one, ^{Stoner and Ottmar 2003}). Another important group that we could identify is that of pelagic fish (present in the upper positive zone from the generalist species), which have thin peduncles, as they are more compact for constant movements. Additionally, the origin of their mouths is in the lower part of the head and far towards the back of the body, allowing the greater opening of the jaw and also generating more suction to capture prey (^{Webb 1984}). From the upper positive part of the PCoA to the lower negative area (axis 2), we could mainly differentiate a change in locomotion strategies from pelagic species that can be characterized by a continuous, sustained swimming speed to those with less locomotion capacity (^{Sambilay 1990}).

The sites in the GC showed lower functional redundancy than NAY. Consequently, the values of the functional indices were higher in northwest locations because of the diversifi cation of traits (Fig. 3). Similar results had been found by other authors, where localities at higher latitudes presented a higher functional diversity of traits than localities closer to the equator (^{Stuart-Smith et al. 2013}, ^{Toussaint et al. 2016}). Considering the geographical proximity of our localities, we corroborated what has been proposed by other authors, which is that these functional indices allow for more sensitive detection of small disturbances (^{Villéger et al. 2017}). Therefore, the communities in northwest Mexico had a greater diver sity of feeding methods and locomotion traits than those at the other sites, due more to the environmental affinity of the species and oceanic conditions than to the heterogeneity of resources (^{Shahraki et al. 2014}). Consequently, the mangrove fish community’s richness was structured by a latitudinal gradient (regardless of how small it was), and the commu nity’s functional composition was governed by the hydroclimatic conditions. We identified different traits in relation to seasonality. In the warm season, ^{Payan-Alejo (2012)} and this study found predators in arid habitat conditions usually related to migratory swimmers, such as members of the Carangidae and Scombridae families. In both studies, adults enter the environment in order to feed on juveniles of some benthopelagic species (e.g., Mugilidae, Gerreidae); we therefore consider them as opportunistic predators.

In the cold season, all authors have reported fewer organisms. However, flatfishes of different families (e.g., Achiridae and Paralichthydae) are found as the temperature in the BCS sites usually drops down to 10 ºC because of the shallow habitats less than 1.5 m depth. Most of the species migrate to deeper waters where the temperature is more stable. In this context, the mangrove is available with low competition. Therefore, the flatfish, being a species that prefers cold temperatures, capitalizes on this habitat (^{Wilber et al. 2013}). Thus, we do not differentiate 2 community structures but rather temporal functional traits that appear/rise from temperature-related conditions. Additionally, some studies have reported that the substratum is also an important distribution factor for flatfishes, particularly that sediment size is related to their habitat preference (^{Spinner et al. 2016}). In this context, it is possible that the high richness found in the localities of BCS and SON (^{Padilla-Serrano et al. 2016}) is related to the presence of sediments of different sizes in the mangrove substratum (a mixture of lime and small-grained sand).

Fish species and environment

From an environmental perspective, fish richness in NAY is separated based on its estuarine conditions, such as greater contribution of freshwater and generally lower temperature. Also, as reported by ^{Henriques et al. (2017)}, habitat extension is related to a higher diversity of species. On the other hand, the localities in the northwest were divided by the evaporation rate, where SON had significantly higher values that separated it from SIN and BCS. From these results, we are able to suggest that salinity, rather than temperature, is a more influential factor in the redundancy of functional traits for the fish in mangroves, even considering the low contribu tion of freshwater in the arid zones.

Given the hypersaline characteristic of some of the sites analyzed, most of the species had a strictly marine and marine-estuarine affinity. In addition, the presence of periph eral species such as silversides had been related to a wide salinity range (^{Veale et al. 2014}). On the other hand, most of the species that had some freshwater affinity were found only in environments with regular freshwater contributions (SIN and NAY), whereas some fish like mojarras and sanddabs were also found in areas with the lowest freshwater input (SON and BCS, ^{Vega-Cendejas and De Santillana 2004}). Although the species composition from arid localities within the GC differ almost 50% from that of NAY, the environmental affinity of the ichthyofauna was dominated mainly by marine species. The absence of some pelagic and freshwater species in the GC may be related to the collection method (seine nets) and to the location of the sampling sites (brackish and coastal sections).

Taxonomic diversity

In the case of the taxonomic structure of the fish commu nity, we identified significant changes at the order level. In other habitats, similar results have been found at this taxo nomic resolution (^{Mueller et al. 2013}), suggesting a highly stable habitat. However, it seems more convenient to consider stability at the family level to generate a probabilistic list of families present in these environments. The most diverse and frequent families in the area are in agreement with those reported by ^{Sheaves (2017)} as the most likely families to be found in mangrove habitats. However, according to the arid conditions of the mangroves present in the GC, the degree of importance of the families was different. Particularly, the most frequent families are those with a high tolerance for environmental stress, such as Mugilidae, Gerreidae, and Lutjanidae (^{Mora and Ospína 2001}).

The differences in the taxonomic structure were greater when introducing NAY into the analysis. In the probability funnels, SIN was significantly different from the other loca tions in the northwest (Fig. 5). However, when comparing the distances that separate the points of SIN and NAY from those for the rest of the localities, SIN behaves as a transi tive area. The taxonomic structure of the mangrove fish of the northwest gradually changes with respect to the more tropical localities in a latitudinal gradient, an aspect iden tified in other regions for which these techniques were used (^{Henriques et al. 2017}). Both localities, SIN and NAY, show an increase in species of the Carangidae, Centropomidae, and Sciaenidae families. Contrasting this, NAY presents a greater number of species with estuarine affinity that are absent from the northern section of the GC, which increases the difference in the taxonomic composition of the arid mangroves with respect to the more southern localities.

Future studies must increase the number of sites covered (particularly sites that currently lack information) and evaluate environmental differences in the assembly process of the fish community. This will aid in understanding whether it is the latitudinal position that governs these functional and taxonomic changes or if it is the suboptimal characteristics of the environment that force a higher rate of diversification of species traits (^{Villéger et al. 2017}).