The “Sistema Arrecifal Veracruzano” National Park (SAVNP) in the southwestern Gulf of Mexico, includes 45 fringing and platform coral reef structures (Liaño-Carrera et al., 2019). The SAVNP is located in a terrigenous environment influenced by freshwater discharge from three rivers: the Antigua in the North, the Papaloapan to the South and the Jamapa in the middle (Krutak, 1997). Three weather seasons have been described for the region: dry, rainy and northerly (Gutiérrez de Velasco and Winant, 1996). Water turbidity varies seasonally according to weather conditions. Five coral reefs were selected for this study: three from the northern group and two from the southern group (Figure 1). The northern reefs, known as the “Holandesas” are three submerged banks located in front of the port of Veracruz which we refer to herein as Holandesa 1 (H1), Holandesa 2 (H2), and Holandesa 3 (H3). The southern reefs, commonly known as Santiaguillo and Anegadilla are emergent platform reefs located offshore of Anton Lizardo village (Liddell and Tunnell, 2011).

Hydrographic surveys were carried out during 2015–2017 using a 30 foot Mexican style skiff equipped with an outboard engine and covered a total area of 6.42 km². The boat was equipped with a semi-portable MBES consisting of a 400 kHz SONIC 2020 multibeam echosounder (R2Sonic, Texas, United States) pole-mounted port-side, amidship 0.60 m below the surface, a sound velocity Micro X probe (AML Oceanographic, British Columbia, Canada) mounted at the rear of the transducer, a differential GPS Receiver Vector VS131 (Hemisphere GNSS, Scottsdale, Arizona, United States) with two antennas mounted 1 m apart on top of the pole, an Inertial Navigation System Ekinox-E (SBG systems, Carrières-sur-Seine, France) aligned to the center line of the vessel and used as our central reference point, and the 1PPS box (Hypack/Xylem, Middle-town, Connecticut, United States) installed between the GPS and computer. The box was used to improve time-tagging by reducing the clock sync error to less than 1 ms.

The multibeam echosounder was operated at 300 kHz with an emitted fanned arc of 256 equidistant beams per ping through a 130° swath sector, sampled. Samples were recorded every 15 μs and received using a 60 kHz bandwidth. Surveys were conducted following 180 parallel-line transects at ∼2 km/h average speed. Transects were planned to obtain 100% sector coverage. Sound velocity casts were obtained using an SV-Xchange probe (AML Oceanographic, British Columbia, Canada). Calibration transects were conducted for each survey to determine pitch, roll, and yaw offsets.

Bathymetry data were processed using Hypack hydrographic software (version 2015). Sound velocity, tide, pitch, roll and yaw offsets corrections were applied according to the manual (Hypack, 2015). The noise was manually removed from every transect. The highest estimated resolution was 0.71 m; however, all data were scaled down to generate a standard 2.5 m gridded matrix (x, y, z) and exported to ArcMap software (ESRI, United States). A 3D bathymetry digital terrain model (DTM) was created in raster format with 2.5 m spatial resolution and projected to the UTM zone 14N Transverse Mercator Datum.

Topographic variability analysis was performed using the Benthic Terrain Modeler (BTM) tools for ArcMap, calculated using a 3 × 3 moving window (Walbridge et al., 2018). To describe the geomorphology and structural complexity of the reefs, several key terrain parameters and their variation were measured, as described in Table 1. A Pearson correlation coefficient was computed to assess the relationship between the depth and the slope, depth and ruggedness, as well as depth and rugosity for the H1, H2, and H3 coral reefs.

For the Holandesas reef group, the depth ranged from 15 to 44 m (see Supplementary Material). Three-dimensional models (Figure 2) of those reefs showed three well-defined geomorphologic zones; (1) the leeward side, (2) the reef flat, and (3) the windward side. H1 and H3 present an oval-shape, while H2 is kidney-shaped with a concave bay at the leeward side. H1 is the largest of the northern reef group and has the deepest basement (Table 2).

For Santiaguillo and Anegadilla reef group, the depth ranged from 2 to 51 m (see Supplementary Material). The reef complex is partially divided by a 0.9 km long channel (Channel-1) oriented northwest to southeast with a width ranging from 0.02 to 0.22 km and depths ranging from 32 to 38 m. Other channels include: Channel-2 which is 0.4 × 0.19 km having depths between 32 to 36 m, Channel-3 is 0.3 × 0.08 km with depths ranging from 33 to 36 m and Channel-4 is 0.7 × 0.07 km with depths ranging from 34 to 43 m. Multiple and well-defined mounds about 0.002–0.008 km² area were identified, away from the main reef within sand beds surrounding the leeward side of Anegadilla (Figure 3A). Spur and groove were the dominant geomorphological features at depths > 30 m, where the reefs show high structural complexity and are continuous between the two reefs (Figure 3B). Our data suggest that Santiaguillo and Anegadilla may constitute a single reef complex, instead of two separate reefs. As such, for the rest of the analyses, these two coral reefs are hereafter referred as the “Santiaguillo and Anegadilla reef complex” (SARC).

The overall results of the estimated terrain parameters for the Holandesas reef group are summarized in Table 2 and the DTM of the topographic variability showing the results of the analysis are presented in Figure 4 in the following order: mean depth, depth standard deviation, slope, aspect, curvature, ruggedness and rugosity. The estimated mean depth ranged from 15 to 42 m and its standard deviation for depth ranged from 0 to 7.2 m. Steepest slopes (Table 2) are found in the 25–43 m depth range in specific areas: northeast windward side for H1, southwest leeward side for H2 and southeast windward side for H3. At the reef’s flat (within 20–25 m depth), H3 (4.8° ± 2.9°) and H1 (4.5° ± 3.1°) showed moderate mean slope. For the H2 reef flat (15–20 m depth), the mean slope was closer to level (3.2° ± 2.8°). The estimated mean slope increased significantly with depth [H1 r_(56,149) = 0.17, p < 0.001, H2 r₍₅₃,₀₃₉₎ = 0.49, p < 0.001, and H3 r_(41,894) = 0.40, p < 0.001; Figure 6A].

The bank walls were mainly southwest and northeast oriented for the well-defined leeward and windward sides, respectively. For H1, the northeast surface area was 0.17 km² (46.8% of the reef) and the southwest surface area was 0.16 km² (46.3% of the reef). For H2, the northeast surface area was 0.18 km² (56.2% of the reef); southwest surface area was 0.09 km² (29%). For H3, the northeast surface area was 0.12 km² (46.2% of the reef); southwest surface area was 0.09 km² (34.6%).

Convex surfaces (curvature range from 100 to 356) delineate the edge of each reef from the Holandesas group. Concave surfaces (curvature range from −100 to −350) revealed notches on the edge of the reefs and depressed areas.

Ruggedness and rugosity analysis revealed highly structurally complex patchy areas that were found at the windward slopes. The H1 showed the maximum ruggedness and H2 showed the highest estimated mean ruggedness among the Holandesas, while H3 presented the lowest ruggedness. Maximum rugosity and estimated mean rugosity were higher for H2, followed for H1 and H3 (Table 2). Estimated mean ruggedness [H1 r₍₅₆,₁₄₀₎ = 0.080, p < 0.001, H2 r₍₅₃,₀₃₉₎ = 0.22, p < 0.001 and H3 r₍₄₁,₈₉₄₎ = 0.14, p < 0.001] and rugosity increased significantly with increasing depth [H1 r₍₅₉,₃₂₂₎ = 0.11, p < 0.001; H2 r₍₅₃,₀₃₉₎ = 0.32, p < 0.001 and H3 r₍₄₁,₈₉₄₎ = 0.33, p < 0.001], showing that deeper areas (up to 43 m) are structurally more complex (Figures 6B,C).

The overall results of the estimated terrain parameters for the SARC are summarized in Table 2 and the DTM of the topographic variability showing the results of the analysis are presented in Figure 5 in the following order: mean depth, depth standard deviation, slope, aspect, curvature, ruggedness and rugosity. The estimated mean depth for SARC ranged from 2 to 49 m and its standard deviation ranged from 0 to 7.2 m. The overall estimated mean slope was 10° ± 10.3°. The mean slope was 6.9° ± 5.6° for shallow areas (5–10 m deep), increased up to 18° ± 12.9° between 20–25 m depth and decreased down to 4.4° ± 3.3° between 45 and 49 m depth (Figure 6A). The steepest areas were at the channel walls (68°) and the spur and groove areas ranging from 59° to 35°, within 35–49 m depth, respectively.

Aspect analysis showed that leeward and the windward surfaces are both oriented in a southwest-northeast direction, with geomorphologically well-defined slopes. The northeast surface area was ≈1.88 km², accounting for 47.2% of the total reef complex. The southwest surface area was approximately 0.97 km², accounting for 24.4% of the total reef complex.

The convex surfaces (curvature range 100–356) delineated the channel walls, several mounds dispersed on the sand bed and spur and groove features. Concave surfaces (curvature range from -100 to -350) delineate the edge of the geomorphological features and depressed areas. In shallow areas (2–20 m deep) and the sand bed, the curvature surface ranged from 50 to −50.

The SARC showed the maximum values for slope, ruggedness and rugosity of all reefs (Table 2). Channel slopes and spur and groove zones within the 25–40 m depth range of the SARC had the highest structural complexity, considering both, ruggedness and rugosity. The shallow areas (5–10 m) showed the lowest structural complexity with 0.002 ± 0.004 mean ruggedness and 1.016 ± 0.03 mean rugosity. The highest mean ruggedness (0.012 ± 0.016) was observed for the 30–35 m depth range, while the highest mean rugosity (1.107 ± 0.15) was observed for the 20–25 m depth range. Structural complexity was lower for the 45–49 m depth range, with mean ruggedness of 0.002 ± 0.003 and rugosity of 1.017 ± 0.020 (Figures 6B,C). For the SARC, the middle depths showed the highest structural complexity, while the shallower and deeper areas showed lower structural complexity.

The platform reefs from the SAVPN are characterized by three zones: (1) leeward slope, (2) reef lagoon or reef flat, which shows a crest delimiting the fore reef and (3) windward slope, a surface orientated to the dominant swell (Chávez et al., 2007). For Santiaguillo and Anegadilla, previous descriptions indicated that they are geomorphologically well-formed platform reefs, characterized by relatively steep slopes (Lara et al., 1992), while Santiaguillo is topographically rough with cliffs (Chávez et al., 2007). Gross-scale geomorphological descriptions of Santiaguillo and Anegadilla from Lara et al. (1992) and Chávez et al. (2007) are consistent with our results. Nonetheless, our fine-scale geomorphological description of SARC provided more details than previous studies, demonstrating previously undescribed geomorphological features, including four channels at the windward side with steep slopes (68°) and high structural complexity based on ruggedness (0.26) and rugosity (3.176) and a spur and groove zones that merge at depths > 30 m and sandy flats with patchy mounds at the leeward side (Figure 3). Spur and groove features, shore-normal coral ridges separated by shore-normal patches of sediment (grooves) are typically the most prominent features of fore reefs (Storlazzi et al., 2003), which is consistent with results for the SARC. Fine-scale spur and groove geomorphological characterization is fundamental for the study and understanding of coral reef hydrodynamics processes. Rogers et al. (2013) showed that spur and groove formations drive circulation patterns and are linked to debris transportation processes outside the coral reef.

Definition of the windward and leeward side reef slopes was based on the aspect analysis, which reflects the predominant orientation to the northeast of the reef surface. According to Wilson et al. (2007), such an orientation may be the result of the local hydrodynamics. For the Holandesas reefs, which are bank submerged reefs, the aspect analysis revealed a well-defined slope at both leeward and windward sides deeper than 15 m. Aspect analysis includes the estimation of the slope and curvature in one measure (Walbridge et al., 2018), therefore, as a result of this study, we highlight the relevance of the aspect analysis to be used for defining zones for different types of reefs in an objective way, including emergent or submerged platform, as well as submerged bank coral reefs.

Surface curvature (concave or convex) is an important geomorphological feature that is liked to ecological processes. Concave surfaces typically accumulate current-carried sediment and foster colonization of sessile organisms (Tokeshi and Arakani, 2012), while convex surfaces favor fish grazing (Bellwood and Choat, 1990). Geomorphological features in coral reefs are related to complex interactions between local oceanographic conditions and ecological functions. Highly complex reef geomorphology is positively correlated to higher biodiversity (Risk, 1972). Analysis of spatial array and complexity of those features can help scientists and decision-makers to identify areas of interest for research, management and conservation, thus, reducing monitoring costs. This is particularly important for developing countries with low budgets, where, having an accurate and detailed bathymetric base map describing the most conspicuous geomorphologic features can be very cost-effective.

Our study is one of the first that describes and characterizes the MCE in the Southern Gulf of Mexico in a systematic way and provides the basis for assessing the ecological and economic implications of including MCE explicitly in management strategies for the region. Results presented herein differ from previous studies regarding reef extension and topographic variability. Our detailed results showed that the platform sizes of the Santiaguillo and Anegadilla reefs have been underestimated in previous studies that were based on aerial photographs and Mexican Naval Charts (Lara et al., 1992), single beam (Ortiz-Lozano et al., 2018), and multibeam (Liaño-Carrera et al., 2019). The differences in reef extension and topography previously reported can be attributed to methodological aspects and the spatial resolution used for data acquisition and interpolation (Kenny et al., 2003; Walbridge et al., 2018).

MBES data offers higher resolution and more accurate maps than the techniques used in previous studies and can thus be valuable for both ecological studies and management including MPA design. High spatial resolution bathymetric maps, such as those produced in our study, are useful for studying ecological, geological and oceanographic processes. According to Kenny et al. (2003), MBES provide larger seafloor coverage with greater precision (≥0.1 m) and it is a useful technique for studying the MCE spatial distribution (Locker et al., 2010). The 3D bathymetric maps and associated physical characterizations can serve as the basis for carefully designed benthic habitat studies in which replicate sampling sites can be selected from within each of the distinctive geomorphological zones. This can be used to validate the predictions of MCE associations with structural complexity in appropriate depth zones. The technological improvement of acoustic methods has allowed the development and use of semi-portable and relatively affordable equipment for geomorphological sea bottom exploration (Kenny et al., 2003). These improvements are critical for both developing countries and remote and poorly studied areas such as those in the southwestern Gulf of Mexico.

Our analyses revealed that the highest structural complexity was shown at mesophotic depths (>30 m), for the Holandesas and SARC, which is a remarkably interesting result that needs further study. Structural complexity has been shown to correlate positively with coral cover (Gardner et al., 2003). We are currently assessing the benthic community structure of MCE in the SAVNP and our preliminary results on the relative benthic cover showed that the general pattern of structure community of deeper areas host MCE communities mainly composed of fleshy macroalgae, turf algae, crustose coralline algae and in less proportion sponges, hermatypic coral and other sessile invertebrates (in preparation).

There is a well-documented decline in shallow coral cover (CC) and structural complexity in the Caribbean region (Gardner et al., 2003; Alvarez-Filip et al., 2009; Medina-Valmaseda et al., 2020). For the SAVNP, there is also a decline in CC that has been documented since 1970s for shallow areas (<20 m), where the coral cover dropped from 50 to 18% (Chávez et al., 2007; Horta-Puga, 2007; Horta-Puga et al., 2015; Pérez-España et al., 2015). Current CC is similar among the northern and southern coral reefs in the SAVNP (Horta-Puga et al., 2015). In this study, mesophotic areas (>30 m depth) showed high structural complexity for both reef groups, considering slope, ruggedness and rugosity (Figure 6). Preliminary results on relative benthic cover for the SARC, show that there is a considerable average of the relative cover of hermatypic coral for mesophotic depths, particularly in spur and groove zones (in preparation). Having the fine-scale geomorphological characterization is important for benthic habitat mapping, according to Sherman et al. (2010); the distribution of well-developed MCE can be in the form of patches (Sherman et al., 2010). According to Bridge et al. (2011), community structure distribution depends on small-scale geomorphology. They demonstrated that distinctive functional ecological groups were consistently associated with specific fine-scale habitat types. For example, autotrophs dominated the submerged reef flat tops while heterotrophs dominated steep slopes.

This study presents evidence of previously undescribed structural complexity and geomorphological continuity in deep SARC reefs (Figure 3). Our new data may have implications for the management of the SAVNP. In 1992, SAVNP was declared a National Park with a protected area of 532.38 km² (Diario Oficial de la Federación [DOF], 1992). The general boundary of the SAVNP was modified in 2012 and is currently covering 655.16 km². Core zone boundaries for both, Santiaguillo and Blanca reefs, cover 11.14 km² (Diario Oficial de la Federación [DOF], 2012; Figure 1). Activities permitted in the core zone include only monitoring, scientific research and maritime equipment installation for navigation. According to the management plan (Programa de Manejo Parque Nacional Sistema Arrecifal veracruzano, 2017) the Santiaguillo core zone includes 7.12 km² of reef area. However, when overlaid onto our bathymetric model using the published coordinates, the polygon covered only 1.84 km², equal to only 25.8% of the SARC and includes flat sandy areas with low structural complexity. Based on our results we propose expanding the core zone to 12.2 km² to include all of the areas with high structural complexity and associated MCE (Figure 1). Marine reserves can enhance the recovery of corals both within and outside the boundaries (Mumby and Harborne, 2010), promote resilience and connectivity by expanding the range of an MPA to include deep reefs (Abelson et al., 2016) since submerged reefs or deep habitats may contribute significantly to larval production (Thomas et al., 2015).

High structural complexity has been related to the provision of ecosystem services, supporting, local culture, tourism and fisheries. MCE provides essential habitat for fish and other species to spawn, shelter, feed and grow to maturity, particularly to important commercial fish (Bejarano et al., 2014; Holstein D. et al., 2016). The more structurally complex are MCE, the greater abundances of fish they have, as they can shelter under steep overhangs and in caves and crevices (Bejarano et al., 2014). In the SAVNP important ecosystem services have been reported: (1) regulation: protection to the coast, protection from erosion processes, (2) support: habitat and biodiversity maintenance, shelter; (3) cultural: recreation, tourism-recreational activities, and academic and research activities; and (4) Provision: fisheries (Arceo et al., 2010; Reyna-González et al., 2014). The structural complexity is an important element that divers take into consideration when selecting diving sites for recreational purposes (Williams and Polunin, 2000). Hence, the highest structural complexity areas are potential sites for recreational diving. However, appropriate regulation is needed for the management of these areas.

We support improved understanding of the MCE and their inclusion in studies, monitoring programs and management strategies (Rooney et al., 2010). The Flower Garden Banks National Marine Sanctuary (FGBNMS) in northwestern Gulf of Mexico, and the Papahānaumokuākea Marine National Monument (PMNM) in northwestern Hawaiian Islands (Papahānaumokuākea Marine National Monument, 2011), are some examples where MCE were included in specific management plans for conservation. FGBNMS includes three separate undersea features, the East and West Flower Garden Banks were designated a National Marine Sanctuary in 1992 and Stetson Bank was added in 1996. Fourteen additional reefs and banks were recently added (Office of National Marine Sanctuaries, 2020), in large part to provide additional protection for sensitive underwater features and marine habitats associated with continental shelf-edge reefs and banks in the northwestern Gulf of Mexico including MCE habitat, as Habitat Areas of Particular Concern (HAPC). HAPCs are high priority areas for conservation, management, or research because they are rare, sensitive, stressed by development, or important to ecosystem function, and they help prioritize and focus conservation effort (U.S. Department of Commerce National Oceanic and Atmospheric Administration, 2012).