Asteroids and Ophiuroids Associated With Sponge Aggregations as a Key to Marine Habitats. A Comparative Analysis Between Avilés Canyons System and El Cachucho, Marine Protected Area

This study analyzes the fauna composition of the community of brittle and sea stars associated with sponge aggregations located in Avilés Canyons System and El Cachucho, Marine Protected Area (MPA). Diverse sampling methods were used depending on bottom morphology, such as rock dredges and specific samplers for sedimentary bottoms, mainly beam trawl models. These banks are made up of sponge and coral species that build a very appropriate substrate for the proliferation of benthic species, which together create Vulnerable Marine Ecosystems that are highly relevant for management and conservation. Among these benthic species, echinoderms are of great interest due to their value as indicators of good habitat. In total, 1261 specimens were collected (934 brittle stars and 327 starfishes), belonging to 42 species (28 ophiuroids and 14 asteroids) from INDEMARES AVILÉS, ECOMARG, and SponGES project surveys. Specimens were distributed among four sponge aggregations (F: fields) that were considered according to the sponge records obtained in the same stations (36). These fields were defined and named based on the five most common sponge species: Aphrocallistes beatrix and Regadrella phoenix (F1: Avilés Canyon); Pheronema carpenteri (F2: Intraslope basin of Le Danois Bank); Asconema setubalense (F3: Le Danois Bank); and Neoschrammeniella aff. bowerbankii (F4: Corbiro Canyon). Faunistic results show that Ophiactis abyssicola (55.55% occurrence), Brisinga endecacnemos, Ophiolycus purpureus, and Peltaster placenta (33.33%) were the most frequently found species in F1; Psilaster andromeda (80%), Pseudarchaster parelii (60%), and Nymphaster arenatus (46.67%) in F2; Ophiura carnea (71.43%) and Ophiacantha smitti (42.86%) in F3; and Ophiacantha densa, and Henricia caudani (100%) in F4. The asteroid and ophiuroid species collected seem to be composed of four different communities that fit to areas with particular morphological and biological features, related to the presence of species specialized in the use of the resources they can find there. In general, the abiotic factor controlling this structure is depth. This assemblage structure, which favors the dissimilarity between the canyons and the bank, is not so clear, since the deepest stations are located on the intraslope basin of El Cachucho, Marine Protected Area (MPA), therefore, using it a priori could lead to misunderstandings. Once the structure of the echinoderm community was known, we compared the expected and obtained results to analyze evidence which should prove the existence of any association between echinoderms and sponges, which enabled us to refute the incongruous hypothesis.


INTRODUCTION
The benthic fauna of the Avilés Canyon System and the El Cachucho, Marine Protected Area (MPA) have been studied from a systematic point of view (mostly in Altuna, 2013;Altuna and Ríos, 2014;Manjón-Cabeza et al., 2014b;García-Guillén et al., 2018;Taboada et al., 2019), and under a general ecosystem approach (Sánchez et al., 2008(Sánchez et al., , 2009(Sánchez et al., , 2014a, providing a complete analysis of habitat scene for management and conservation. These studies were very valuable for these areas to be considered a vulnerable ecosystem and, as such, a protected area by the European Union (92/43/CEE) (Sánchez et al., 2014a(Sánchez et al., ,b, 2017Punzón et al., 2016;Rodríguez-Basalo et al., 2019a). However, specific level studies are scarce, from the biological and ecological approach.
Echinoderms, together with sponges and corals, constitute the most important groups, in both biomass and/or abundance (Sánchez et al., 2008(Sánchez et al., , 2009(Sánchez et al., , 2014a(Sánchez et al., , 2017Manjón-Cabeza et al., 2014b;García-Guillén et al., 2018), as well as specific richness of the deep seabed. The importance of considering the association between these three groups is considerable, and all these previous studies provide us with a unique opportunity to carry out different kinds of analysis. Moreover, many of these bottoms are made up of sponge and coral species that build a very appropriate substrate for the proliferation of other benthic species. Among these benthic species, echinoderms are of high interest, mostly because of their inherent needs in order to survive (Murillo et al., 2012;Manjón-Cabeza et al., 2014a,b;Andrino-Abelaira, 2015;Gómez-Delgado, 2015;Murillo, 2015;Palma-Sevilla, 2015;Hurtado-García, 2016;Moya, 2016;Mah, 2020), or they are found to be associated to specific communities, leading them to be indicators of good habitat.
Recently, Ríos et al. (2020) studied the community composition and characterization of sponge aggregations in the Cantabrian Sea, showing that these aggregations constitute structuring habitats (3D communities) that function as support for different types of benthic communities (Peña-Cantero and Manjón-Cabeza, 2014, among others) and that certain benthic groups seem to have an intrinsic association with these sponge aggregations.
Based on these results, it should be possible to find an echinoderm reliance, preference, or association with these aggregations, but precise studies on the subject were not carried out, which is the main objective of the present study.
For this purpose, four sampling areas were chosen, two in the canyons and two on El Cachucho, Marine Protected Area (MPA), as well as two groups of echinoderms, starfishes, and brittle stars, because these echinoderm taxa tend to be indicators of different morphological and granulometry bottom features and of the assemblage of benthic communities.
According to this objective, to ascertain whether there exists any link between asteroids and ophiuroids and sponge aggregations on sea beds, we proposed the following hypotheses ( Figure 1): 1 Although sponge aggregations present different specific compositions, there is not enough evidence to consider that the asteroid/ophiuroid community in canyons is different from on El Cachucho, Marine Protected Area (MPA) (Figures 1A,C). 2 These echinoderm species distributions will enable us to define four different echinoderm assemblages related to each sponge field ( Figure 1B). 3 There are two different taxocoenosis, one in the canyons, and another on El Cachucho, Marine Protected Area (MPA) ( Figure 1D).

Study Area and Description of Work Fields
The study area is located in the northern continental margin of the Iberian Peninsula in the Cantabrian Sea (Figure 2A), specifically the Avilés Canyon Systems (ACS: Figure 2B) and  El Cachucho, Marine Protected Area (MPA) (LDB: Figure 2C), which present peculiar geomorphological and habitat features described in several previous papers (Ballesteros et al., 2006;Sánchez et al., 2008Sánchez et al., , 2009Van Rooij et al., 2010). These features together with station locations were analyzed and georeferenced by ArcGIS 10.7. This work focuses on two canyon heads, Avilés Canyon (AC) and the Corbiro Canyon (CC) (Figure 2B), and Le Danois Bank and its intraslope basin (Figure 2C), based on sponge aggregation types of settlements described by different authors (Prado et al., 2019;Rodríguez-Basalo et al., 2019b;Ríos et al., 2020). Four sponge fields were established according to sponge species: F1 (AC), Aphrocallistes Beatrix Gray, 1858, and Regadrella phoenix Schmidt, 1880; F2  These four species were studied based on previous knowledge of the study area (García-Alegre et al., 2014;Sánchez et al., 2014a) and the criteria by which vulnerable marine ecosystems such as sponge grounds are considered (Hogg et al., 2010;Maldonado et al., 2016): they support high biodiversity of other species, are fragile and unlikely to recover from trawl damage, and are limited to discrete areas with suitable environmental conditions. Sampling was carried out using different trawl gear and surveys, as described in previous studies (Sánchez et al., 2008;Rapp, 2019). A total of 36 stations ( Table 1) meeting the aforedescribed criteria requirements were selected.

Material
The biological material consisted of 1261 specimens: 934 ophiuroids and 327 asteroids. Specimens were photographed and conserved in ethanol and identified based on their morphological characteristics using specialized literature (Mortensen, 1927(Mortensen, , 1933Lieberkind, 1935;Paterson, 1985;Clark and Downey, 1992;Southward and Campbell, 2006), and the appropriate protocols for their visualization (light microscopy or SEM).

Data Analysis
Echinoderm and sponge occurrence frequencies were calculated (to analyze the general faunal composition). Echinoderm species were classified into four categories according to their frequency across stations, which is a surrogate for evaluating their importance in the community: the most common species (50% of stations), very common species (between 25 and 50%), common species (between 25 and 10%), and rare or accidental species (<5%) (Mora, 1980;Manjón-Cabeza and García Raso, 1994;Manjón-Cabeza and Ramos, 2003, among others). To investigate the structure of echinoderms, assemblage similarities, related to sponge aggregations, were computed by a hierarchical cluster analysis (classification) using the UPGMA agglomerative algorithm (Sneath and Sokal, 1973; RMACOQUI ver. 1.0 software Olivero et al., 2011;RStudio Ver. 0.99.473) made on the similarity matrix of the Baroni-Urbani coefficients calculated from presence/absence data (Baroni-Urbani and Buser, 1976). The robustness of each cluster was supported by a test of biological significance of the boundaries between echinoderm assemblages. Strong and weak boundaries were defined between assemblages following (McCoy et al., 1986, P < 0.001). A strong boundary separates two significantly different clusters (red node number in Figure 7). A weak boundary (green asterisk in Figure 7) measures the homogeneity of species distribution between stations. When boundaries are not significant, it means that species are randomly distributed. Boundary analysis followed Olivero et al. (1998Olivero et al. ( , 2011. Stations were identified using Canonical Correspondence Analysis (CCA) computed from the presence/absence matrix and based on the eigenvalues of χ 2 distances between all data points (Ter Braak and Prentice, 1988;Hennebert and Lees, 1991;Legendre and Legendre, 1998), using PAST (paleontological statistics, ver. 3.25 computer program (Hammer et al., 2001).
Three analyses were performed: (CCA1) only with% sponge occurrence as the biotic variable; (CCA2) with all noncorrelated abiotic variables, in this case only depth and granulometry (latitude and longitude were discarded) and biotic variable (% sponge occurrence); (CCA3) only abiotic variables (granulometry and depth). These were used to define ordination axes on which echinoderm data (with both stations and specimens) were plotted. Environmental variables were plotted as well as correlations with ordination axes.

General Faunal Composition
The faunal composition of the study area presented 42 echinoderm species (28 ophiuroids and 14 starfishes) and 21 sponge species (Table 2 and Figures 3-6). Ophiolycus purpureus, Ophiophrixus spinosus, and Ophiotreta valenciennesi were new records for the area and as such will be included in the "Echinodermata Spanish Check List" (2020 in press, update of Ministerio de Agricultura y Pesca, Alimentación y Medio Ambiente, 2017).

Work Field Faunal Features
Sponge (S) and echinoderm (E) composition of each field are shown in Figure 3. F1 (Aphrocallistes and Regadrella aggregation) presented nine sponge species, the most frequent of which were Aphrocallistes beatrix (77.78%) followed by Regadrella phoenix (22.22%) and Pachastrella ovisternata (22.22%) (Figures 3A-S). Twenty-three echinoderm species were recorded, 12 of which were exclusive ( Figure 3B): the most commonly occurring echinoderm species was Ophiactis abyssicola (55.5%) (Figures 4A-E). F2 (Pheronema aggregation) was represented by Pheronema carpenteri with 92% of sponge species occurrence (Figure 4A-S), while the echinoderm community consisted of 24 species, nine of which were only present in F2 (Figure 3B). In this case the most common species were Psilaster andromeda with 75% and Pseudarchaster parelii with 56.25% occurrence, respectively. F3 (Asconema aggregation): Two sponges, Asconema setubalense, Podospongia loveni, established a new field with 40% occurrence. Regarding echinoderms, 12 species were distinguished, three of which were exclusive. Echinoderm composition was 14 species, 10 of which were starfish, and four of which were very rare brittle stars (in terms of occurrence). The species featured were Ophiura carnea with 55.56% and Ophiacantha smitti with 33.33%. F4 (Neoschrammeniella aggregation): the most frequent sponges (50%) were Neoschrammeniella aff. bowerbankii, Pachastrella ovisternata, and Geodia pachydermata. Six echinoderm species made up this field, three of which were exclusive. Ophiacantha densa and Henricia caudani, were recorded in all stations of this field.

Classification Analysis
Presence-absence species matrix, and occurrence percentage used for data analysis (Table 3). Cluster results display a clear discontinuity between different station groups (Figure 7, see dot colors), revealing the existence of four distinctive assemblages (G1, G2, G3, and G4), divided by strong boundaries (red nodes), while three stations did not match up with any other.
G2 contains four stations from F1 and two stations from F4. In this case there was no evidence of homogeneity, but it was a consolidated group (node 31; P < 0.00001) as well as G1. This group was mainly composed of Henricia caudani and Ophiactis balli ( Figure 6E) (both 50% occurrence), followed by Ceramaster grenadensis, Ophiacantha densa, and Peltaster placenta (33.33%).
G1 also presents a mix of stations from F2 (2) and F3 (6). It was made up of a group with a high dissimilarity with the rest, consisting exclusively of ophiuroid species. The most frequent species was Ophiura carnea.

Reliability of Setting to Preset Fields
Percentage fit of the different cluster groups to the working fields are: G1, 75% of F3 stations; G2, 33.33% of F4 (considering G2 is the only group with F4 stations); G3, 94.33% of F2; and G4 representing 100% of F1stations.

Ordination Analysis
Results from CCA analysis are shown in Table 4 (Figures 8, 9). CCA1 was carried out only with sponge frequency as a biotic variable and echinoderms were ordered according to these, which was not significant (Table 4). However, when granulometry was taken into account (CCA2) (Figure 8), the CCA results became significant despite the very low% explanation. In the case of CCA3, only granulometry and depth significance were taken into account, showing the highest significance (Table 4 and Figure 9).

Echinoderm Assemblages and Control of Their Environmental Variables
The echinoderm community seems to be composed of four different communities that fit to areas with particular morphological and biological features, related to the presence of species specialized in the use of the resources they can find there (Sánchez et al., 2008;Ríos et al., 2020).
In general, the abiotic factor that mainly controls this community structure is depth. In fact, it is very frequent in echinoderm assemblage studies (Manjón-Cabeza and Ramos, 2003;Moya, 2016). This assemblage structure, which favors the dissimilarity between the canyons and El Cachucho, Marine Protected Area (MPA), is not so clear, since the deepest stations are located on the intraslope basin of the bank, therefore, its use a priori could lead to misunderstandings (Figures 8, 9).
Group G4 would represent a community associated to deep, hard bottoms covered by coarse sands. The stations fit perfectly at the head of the Avilés Canyon (Figures 2, 7). The taxa making up this community are suspension feeder species, such as brisingid species like Brisinga endecacnemos and Novodina pandina (Downey, 1986;Clark and Downey, 1992), which take advantage of the pedestals offered by the rock outcrops, or coral patches of Madrepora oculata Linnaeus, 1758; and Desmophyllum pertusum; Linnaeus, 1758 (Sánchez et al., 2014b). Coral aggregations are used as support by some ophiacanthids such as Ophiochondrus armatus, or species of the Genus Ophiacantha. However, Ophiactis abyssicola, like the rest of species in the Ophiactis genus, lives associated with bottoms that have cavities available, such as oscula sponges (Schejter et al., 2012;Sivadas et al., 2014;Çinar et al., 2019), little holes in stones or associated with dead corals or rest of calcareous algae (rhodoliths), which provides them with shelter (Gofas et al., 2014;Manjón-Cabeza et al., 2014c;Palma-Sevilla, 2015).
The community closest to G4 is G3, which is found on the intraslope basin of the Le Danois bank. This affinity is mainly due to species richness (Figures 3, 7). Depth seems to be the abiotic environmental factor controlling these two communities (Figures 7, 8), although other variables should be taken into account, such as the slope. In fact, steep areas could favor the settlement of structuring species (3D) such as Pheronema carpenteri, which would determine the echinoderm community (Cristobo et al., 2010;Sánchez et al., 2010Sánchez et al., , 2014a. In contrast with the rest of the communities studied, the species making up this one are mainly Asteroids, like Psilaster

Species
Ocurrence (%) A11DR07 A11DR11 A410DR07 A410DR08 A710DR01 A710DR06 A710DR09 A710DR10 A710DR12 E3St2 E3V01 E3V03 E3V07 E4St1 E4V03 E4V08 E4V10  Species E8G02 E8G03 E8V03 E8V06 E8V09 E9G09 E9V01 E9V02 E9V03 E9V10 S17BT03 S17BT04 S17BT09 S17BT10 S17BT11 S17BT12 S17DR04 S17DR11 S17DR16 andromeda, Pseudarchaster parelii, Nymphaster arenatus, and Zoroaster fulgens, and others of lesser occurrence such as Plutonaster bifrons and Pontaster tenuispinus. This community is the most homogeneous one, and the one with the greatest specific richness. In this case, the abiotic factor mainly affecting species composition is the presence of fine sands with a high content of organic matter, preferred by sand burrow species such as those belonging to the Genus Amphiura (Sánchez et al., 2008). On the other hand, G1 is represented by stations located mostly on Le Danois Bank which has fine sand bottoms with Asconema setubalense. These features are very well defined in previous publications where the characteristic habitats of the area are described, and a large occurrence of Callogorgia veticillata (Pallás, 1766) is attributed to the upper zone of the bank (Sánchez et al., 2008(Sánchez et al., , 2017. These features, indeed, explain the presence of species as diverse as Ophiura carnea, Ophiomyces grandis, Ophiocten affinis, and Ophiothamnus affinis, which live on sandy bottoms; and species of the Genus Ophiacantha and Ophiothrix (Granja-Fernández et al., 2014) that have a preference for corals, especially gorgonians.

Amphiura filiformis
G2 group does not make much biological sense and its stations seem to be a consequence of the scarcity of stations from F4 in the Corbiro Canyon (only two). Another reason that could explain this artifact would be due to stations from F1, associated would have a similar sponge species contents, and Aphrocallistes beatrix were not as frequent as in the rest of stations from Avilés Canyon or in the other way round, Neoschrammeniella aff. bowerbankii, does not represent any echinoderm association. Therefore, the sea star and brittle star community of Corbiro Canyon should be more profusely studied in the near future.

Does There Exist a Real Association of Echinoderms With Sponge Aggregations?
Once the structure of the echinoderm community was known, we were able to compare the expected and obtained results in order to analyze the evidence which should prove the existence of any association between echinoderms and sponges, which enable us to refute the incongruous hypothesis.
In this way, the results obtained do not conform to any of the proposed hypotheses (Figures 1, 7), and the reasons that would explain this issue are developed below.
(1) Although station fit is quite high in G1, G3, and G4 clusters, the G2 cluster has a very low percentage station affiliation.
FIGURE 7 | Cluster resulting from Echinoderm species classification (Baroni-Urbani index). Noted group color and dot colors. They are related to the cluster (Baroni-Urbani similarity coefficient) significant group, G1: purple; G2: black; G3: green; G4: blue. Node number in red illustrates strong boundaries segregating significant (P < 0.001) clusters (or groups), whereas green asterisks denote where weak boundaries (P < 0.001) were found, measuring the homogeneity of species distribution between stations included in these clusters or group. No node number shows non-significant boundaries (P > 0.001), in these cases species are randomly distributed (following Olivero et al., 1998Olivero et al., , 2011. Sector diagrams show species occurrence percentage related to each significant cluster. (2) Asteroid/ophiuroid community assemblages do not fit the sponge species composition.
Given the high percentage of adjustment that some of the fields present, it is possible that these small imbalances can be explained, since fields were delimited on the basis of main sponge species, although this occurrence may vary between stations (Figures 4, 8), and then some of them did not fit the field as we expected.
On the other hand, CCA using only an environmental variable set is more significant than using it in combination with occurrence of sponges. There are two ways to address this question: (1) sponges are not a very good biotic factor to control the echinoderm community; (2) echinoderms depend more on other bottom types (for instance related to granulometry, Figure 9), so the variable data set should be improved.
Finally, these results enable us to infer that the association of asteroids and ophiuroids with sponge aggregations is conditioned to environmental factors, like granulometry, which control fields such as habitat. Sponge species composition, or the structure they FIGURE 8 | Ordination studies. Canonical Correspondence Analysis (CCA2). Physical environmental variables and biotic variable (most relevant sponges): Granulometry characteristics and depth. F1, 2, 3, 4 described in Figure 2. Statistical parameters and null hypothesis shown in Table 4. FIGURE 9 | Ordination studies. Canonical Correspondence Analysis (CCA3). Physical environmental variables and granulometry characteristics and depth. F1, F2, F3, F4 described in Figure 2. CCA3: variables: sponge occurrence%, depth and granulometry. Statistical parameters and null hypothesis shown in Table 4. provide, would not be the main reason for explaining echinoderm assemblage structure.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
MM-C and PR conceived and designed the study and wrote the manuscript. PR, TI, AR-B, FS, and JC collected the specimens and the pre-identification major taxa. PR and JC identified the sponges. MM-C, PR, AM-R, and LG-G identified the echinoderms, analyzed the data, prepared the figures and tables, reviewed drafts of the manuscript, and helped to writing the manuscript. MM-C, FS, and JC acquired the funding. All authors contributed to the article and approved the submitted version.