Anthropogenic impact is negatively related to coral health in Sicily (Mediterranean Sea)

Shallow-water marine organisms are among the first to suffer from combined effects of natural and anthropogenic drivers. The orange coral Astroides calycularis is a shallow-water bioconstructor species endemic to the Mediterranean Sea. Although raising conservation interest, also given its special position within the Dendrophylliidae, information about the threats to its health is scant. We investigated the health status of A. calycularis at five locations in northwestern Sicily along a gradient of cumulative human impact and the most probable origin of the threats to this species, including anthropogenic land-based and sea-based threats. Cumulative human impact appeared inversely related to the performance of A. calycularis at population, colony, and polyp levels. Sea-based human impacts appeared among the most likely causes of the variation observed. The reduction in polyp length can limit the reproductive performance of A. calycularis, while the decrease of percent cover and colony area is expected to impair its peculiar feeding behaviour by limiting the exploitable dimensional range of prey and, ultimately, reef functioning. This endangered habitat-forming species appeared susceptible to anthropogenic pressures, suggesting the need to re-assess its vulnerability status. Creating microprotected areas with specific restrictions to sea-based human impacts could be the best practice preserve these bioconstructions.


Methods
Study area and assessment of impacts. The study was conducted at five locations lying approximately at the same latitude (between 37°N and 38°N) on the north-western coast of Sicily, western Mediterranean. The study area is characterized by a gradient of increasing urbanization moving eastward from the Egadi Archipelago to the Gulf of Palermo (Fig. 1). Five locations tens of kilometres apart, characterized by increasing levels of anthropogenic pressure were randomly selected among those hosting well-developed populations of A. calycularis after a preliminary survey (Fig. 1). Locations, according to a scale from minimum to maximum pressure, included: Marettimo (Loc 1) and Favignana (Loc 2) Islands (both inside the Egadi Islands marine protected area), the Zingaro coast (Loc 3) (whose land part is a nature reserve), Capo Gallo At each location, human impact was assessed adapting the method proposed by Micheli et al. 24 to the specific case study. Among the agents of change reported by these authors, the following were considered as potentially impacting A. calycularis: human population, fertilizers, pesticides, urban runoff, shipping, invasive species, artisanal fishing and four types of commercial fishing (destructive demersal, demersal with low by-catch, demersal with high by-catch, pelagic). In order to obtain the matrix of pressures of the above mentioned agents for each location, the value of each agent of change at each location was extracted using ArcGis 9.3 from the georeferenced layers, 1 km 2 resolution, developed by Micheli et al. 24 . The raw matrix thus obtained was then standardized to make each response variable (agent of change) vary within the same scale, from 0 to 1. This was achieved by dividing the scores of the agent at the five locations by the maximum value obtained for such given agent across the five studied locations. The obtained standardized matrix allowed for relative comparisons of the importance of the different agents across the study area. Subsequently, the potential effect of each agent of change on A. calycularis was elicited from experts on the biology of the species and the study area. Experts, who were not aware of the experimental design and data, were asked to rank the importance of each agent of change along a scale of discrete values between 0 (not important at all) and 10 (very important). For each agent, the average score was obtained after reducing the variability of judgments, by excluding a posteriori the minimum and maximum value (outliers) for each variable (if two or more experts expressed the same outlier judgment, only one of the values was excluded).
The estimated impact of each agent of change on the studied species was calculated by multiplying the standardized value of the pressure exerted by a given agent at a given location (from the GIS layers by Micheli et al. 24 ) by its average importance score across the elicited experts. For example, if the standardized value of urban runoff at location 1 was 0.5, and the average expert judgment for that agent was 7, the impact of urban runoff on the species was 0.5 * 7 (see Supplementary Table S6). The cumulative human impact (Ic) for each location was obtained by summing up the estimated impact values for the different agents of change. Ic was partitioned into sea-based (I S ) and land-based (I L ) human impacts, so that I C = I L + I S at any given location. I S was obtained as described above considering shipping, invasive species, artisanal fishing and four types of commercial fishing. I L was obtained considering human population, fertilizers, pesticides, and urban runoff. The obtained I C , I L and I S values are shown in Supplementary Table S2. I C and I L increased along the west-east direction ( Fig. 1; see Supplementary  Table S2), while I S increased along the west-east direction from Loc 1 to Loc 4, while Loc 5 showed an I S value lower than Loc 3 and Loc 4, but higher than Loc 1 and Loc 2.
The five locations were chosen within the narrowest latitudinal range possible in order to reduce potential effects of natural environmental variation on the performance of this stenotherm species. Nonetheless, we also tested the potential correlation with average sea surface temperature (SST AVG) and average annual variation of sea surface temperature (SST VAR), i.e. average difference between maximum and minimum annual temperatures. Five-year temperature data (of the five years preceding the sampling: 2009 to 2013) were extracted from the NASA online application Giovanni (https://giovanni.gsfc.nasa.gov/giovanni/). SST AVG increased along the west-east direction form 19 Biological response variables. Fieldwork was carried out between May and July 2014. In order to assess the health status of the species along the cumulative impact gradient and to check it at multiple levels of biological organization, three response variables were considered. The percent cover, considered an important indicator of health status in coral reef studies 75,76 , was analyzed as a proxy of health at the population level. Colony area and polyp length were selected as proxies of health at the colony and polyp levels, as these biometric parameters have been used in biometric, reproductive biology, and population dynamic studies of several temperate coral species 55,77 as a measure of size, which is a good proxy of health in long-lived marine invertebrates like A. calycularis 55,77,78 . percent cover. The percent cover of A. calycularis was quantified in three sites hundreds of m apart within each location, at depths between 1 and 3 meters. Sites were chosen in order to represent the optimum condition for the population in the analyzed location (i.e., where the species appeared more abundant). At each site, percent cover was determined by taking six photos, which acted as replicates, on vertical/sub-vertical walls in order to avoid any confounding effect of substrate slope on the response variables. The camera was equipped with two laser pointers set at a fixed 20 cm distance, which produced a reference scale for each picture (Fig. 2). The digitizing software Image J was used to calibrate each photo setting the reference distance given by the two pointers and to measure the area occupied by A. calycularis colonies inside 25 × 25 cm 2 quadrats by tracing their outline with a hand-controlled mouse on the digital image. colony and polyp biometry. Three to five colonies of A. calycularis were collected in the same sites used for the analysis of percent cover. The sampled colonies were randomly chosen among those settled within previously photographed frames. In each location, samples were collected from the three sites, except for Loc 1 and Loc 2 where only two sites were sampled for logistic reasons. The collected colonies were kept immersed in 10% sodium hypochlorite for 10 days to remove the soft tissue, then dried at 50 °C for 4 days. Skeletons were observed under a binocular microscope to remove fragments of substratum and calcareous deposits produced by epibionts. The length (L C : major axis of the colony) and width (W C : minor axis of the colony) were measured with a calliper in order to obtain the colony area that was calculated using the formula for the ellipse = π * * A (L W ) 4 C C . Colony area was chosen as the main biometric parameter for colonies as it is considered a more accurate and representative measure of colony size than colony length [79][80][81][82] .
The number, length (i.e., maximum axis of the oral disc) and width (i.e., minimum axis of the oral disc) of all polyps was determined for each colony 55 . Polyp length was selected as a proxy of polyp size. The distribution of polyp size classes at each location was then determined by considering 1 mm as the threshold between a size class and the subsequent one 55 .
Statistical analyses. The proportions of variation of (1) cover, (2) colony area, and (3) polyp length explained by I C , I L , I S , SST AVG and SST VAR were checked by non-parametric multiple regression analyses using the DISTLM-forward procedure (distance-based multivariate analysis for fitting a linear model using forward selection of factors). The analyses were performed using the software PRIMER v6 including the PERMANOVA + add-on package 83 .
Differences in the frequency distribution of polyp size classes among locations were tested by means of a Kruskal-Wallis test. Where significant, pairwise comparisons between locations were performed via Kolmogorov-Smirnov tests.

Results
Overall, considering the results of all the DISTLM tests, I C , I L , SST AVG and SST VAR all appeared highly correlated between each other (R > 0.78; Supplementary Table S4). Thus, disentangling the relative portions of explained variation of each variable in the study area was not possible, while I S resulted as the most independent predictor variable (Supplementary Table S4). percent cover. The percent cover of A. calycularis decreased by ~38% between the two extremes of the cumulative human impact gradient from Loc 1 (low I C ) to Loc 5 (high I C ) (Figs 1-3).
The DISTLM marginal tests revealed that I S significantly explained ~53% of the observed variation followed by SST AVG and SST VAR (~39% each) and I L (~33%) ( Table 1). Given the above mentioned high correlation among predictor variables, the DISTLM sequential test, fitting I C , SST AVG and SST VAR within the same model, revealed that SST VAR and SST AVG significantly explained ~39% and ~12% of the observed variation respectively, while the test for I C was not significant after the forward selection procedure ( Table 2). The DISTLM sequential test, fitting I L , I S , SST AVG and SST VAR within the same model, revealed that I S significantly explained ~53% of the observed variation followed by SST VAR (~11%), I L (~8%) and SST AVG (~2%) ( Table 3). colony and polyp biometry. Colony area decreased by ~36% between the two extremes along the gradient of I C , from Loc 1 to Loc 5 (Figs 2 and 4).
I S explained ~25% of colony area variation according to the DISTLM marginal test, followed by I C (~19%), I L (~17%), SST AVG (~11%), while the test for SST VAR was not significant ( Table 1). The DISTLM sequential test, fitting I C , SST AVG and SST VAR within the same model, revealed that I C significantly explained ~19% of the observed variation, while the tests for SST VAR and SST AVG were not significant ( Table 2). The DISTLM sequential test, fitting I L , I S , SST AVG and SST VAR within the same model, revealed that I S significantly explained ~25% of the observed variation while the tests for the other variables were not significant (Table 3).
Similarly to colony area, polyp length decreased from Loc 1 to Loc 5 (Fig. 5). According to the DISTLM marginal test, I S ranked first among the predictor variables, explaining ~13% of the variation, followed by I C (~10%), I L and SST AVG (~9%), and SST VAR (~6%). The DISTLM sequential test, fitting I C , SST AVG and SST VAR within the same model, revealed that I C explained ~10% of the observed variation, while the tests for the other variables were not significant ( Table 2). The DISTLM sequential test, fitting I L , I S , SST AVG and SST VAR within the same model, revealed that I S and SST AVG significantly explained ~13%, and ~2% of the observed variation respectively, while the tests for SST VAR and I L were not significant (Table 3).
Polyp size classes ranged from 2 to 13 mm. The minimum number of polyps measured in a colony was 13 (at Loc 4), the maximum was 134 (at Loc 1). The gradient of human impact had an effect on the median of the size distribution of polyps among locations (H 4 = 320.64, p < 0.001, n = 2367). Pair-wise comparisons showed that the size frequency distributions of polyps were not significantly different between Loc 1 and Loc 2, while they differed among the remaining localities, and among these and the former two (see Supplementary Table S1). Following the gradients of I C , I L , and I S , the observed pattern of distribution of polyp size classes appeared coherent with the I S gradient in the study area, with Loc 1 and Loc 2 showing comparable I S values (8.5 and 8.8 respectively), higher than those observed at Loc 3 and Loc 5 (9.4 and 9 respectively), while the highest I S value (10.6) was found at Loc 4 (see Supplementary Table S2, Fig. 6). Average polyp length, standard deviation and number of size classes were higher at Loc 1 and Loc 2, and lowest at Loc 4 and Loc 5 (see Supplementary Table S5). The modal class ranged www.nature.com/scientificreports www.nature.com/scientificreports/ from 6 to 7 mm, except for Loc 2 where the size distribution had two modes at 6 and 9 mm. Skewness was generally positive except for Loc 5, with absolute maximum values at Loc 1 and Loc 5. Kurtosis was generally negative except for Loc 4 where size classes were fewer and frequencies around the mode dropped abruptly, particularly to the right end of their distribution (Fig. 6).

Discussion
Anthropogenic pressure negatively relates to the health status of A. calycularis at all the organization levels here considered, from polyp to population. Our research suggests that both land-and sea-based threats impact the species performance, although local variation in sea surface temperature may have a role in explaining the observed pattern. However, the linearly decreasing trend in I S recorded in the study area had a stronger explanatory power than I L , I C , and SST spatial-temporal variation in describing the patterns observed at all organization levels.
In this study, I S showed a stronger relationship with percent cover, explaining more than 50% of its variability. The aim of this study was not to identify the weight of every single agent of change, as the study was not designed for this purpose. However, it is interesting to note that artisanal (i.e., inshore) fishing effort appeared strongly related with the reduction in A. calycularis cover along the analyzed gradient (data not shown). Inshore fishing is common in the investigated area, and nets thrown in proximity of the rocky substrate might damage the colonies. Professional and recreational fishing activities, as well as lost fishing gear have been shown to negatively affect not www.nature.com/scientificreports www.nature.com/scientificreports/ only the abundance but also the morphological and mechanical characteristics of Mediterranean hard-bottom communities 29,84,85 . Entanglement on habitat-forming macro-invertebrates, such as corals, damage their structures through mechanical friction, thus favouring parasitic colonization of the exposed skeleton ultimately leading to death 86 . However, the most likely explanation for the observed correlation in A. calycularis could be that fishing and I S in general could represent proxies of a generalized overexploitation in densely populated areas, thus increasing the possibility of accidental contacts and subsequent damage or even detachment of the colonies. Recreational, and often illegal harvesting of shallow-water sea urchins, abalones and other traditional seafood, as well as spearfishing, are common activities in the study area. These activities often imply direct human contact with the reefs and sometimes the illegal use of poison (e.g., copper sulphate to catch octopuses), as reported by local fishermen. The reduction in percent cover by ~38% along the cumulative human impact gradient could also influence the orange coral bioconstructions at a functional level. Polyps belonging to neighbouring colonies are known to cooperate to catch large prey, and it has been hypothesized that this behaviour may explain the tendency of this species to aggregate and form biogenic reefs 50 . If aggregation enhances the ability of this species to catch large jellyfish, the reduction of percent cover in impacted areas is expected to limit A. calycularis ability in capturing such prey, affecting its nourishment and fitness and, ultimately, impairing an ecological function of the biogenic reefs. Following this hypothesis, such indirect, potential effect of human impact could be widespread among other cnidarians. In fact predatory cooperation among polyps has been observed in other gregarious/ colonial species [87][88][89] .  Table 2. Non-parametric multiple regression analysis. DISTLM-forward analysis sequential tests fitting cumulative human impact (I C ), average sea surface temperature (SST AVG) and average annual variation of sea surface temperature (SST VAR) in the same model in order to analyze their combined contributions in explaining the variation of Astroides calycularis percent cover, colony area, and polyp length in the study area. www.nature.com/scientificreports www.nature.com/scientificreports/ Besides percent cover, also colony area and polyp length significantly decreased by ~36% and ~13% respectively along the gradient (Figs 4 and 5). I S showed significant relations with these morphometric variables. Additionally, polyp size frequency distributions at the most impacted locations showed a decrease in the number of size classes with a lack of the largest ones, and significant transitions of size, suggesting a reduction in recruitment efficiency 90 . Also, the colony with the largest number of polyps was found at Loc 1 (N = 134), while the smallest ones were found at Loc 4 (N = 13) and Loc 5 (N = 19). Thus, the observed decline in the biotic variables is a likely effect of the overexploitation of the coastal area previously discussed for percent cover. Also, mechanical damage with accidental removal of colonies in highly impacted areas, besides reducing the overall coral cover, may cause the reduction in number of larger colonies in favour of new, smaller ones that would grow in the space made available by the damaged colonies. This would also explain the reduction of polyp length as well as the reduction in the number of size classes along the gradient, since larger polyps are usually found at the centre of larger colonies 43 .
However, the possible contributions of land-based human impact and of spatial temporal variation of SST in the study area should not be overlooked. As far as the SST and its annual variation are concerned, it is well known that A. calycularis is a thermophilic stenotherm species and large temperature fluctuations might impair its performance 16 . In the study area, the largest difference in average annual SST among locations was about 1 °C, while the difference among Loc 1 (Marettimo) and Loc 4 (Capo Gallo), where the species showed the best and the worst performance, respectively, was 0.66 °C. Annual SST variation ranged from 12 °C at Loc 2 to 13.51 °C at Loc 5, but at Loc 4 it was only 0.47 °C higher than Loc 1. Winter minimum SST below 14 °C is known to limit the distribution of this species 16 , but this is not the case of the study area. Moreover, A. calycularis had its best performance at Loc 1 and Loc 2 where the lowest minimum SST occurred. The highest maximum SST were detected at Loc 4 and Loc 5, which were about 1 °C higher than at Loc 1. It should also be noted that the highest summer average SST detected in the study area (about 28 °C in August at Loc 5) is similar to other more southern areas where A. calycularis thrives (e.g. the Pelagie Archipelago) 53 . In summary, although SST and its annual variation may be related to the health status of the species in the study area, it is unlikely that they may have had a leading role in determining the observed variation in A. calycularis health status. The observed high correlation among SST AVG, SST VAR, I c and I L could likely be spurious.
As far as the land-based human impact is concerned, the coastal area between Loc 4 and Loc 5 occurs within the Gulf of Palermo, which is among the most densely populated areas in Sicily (above 4000 citizens per km 2 ). This area is characterized by the presence of rivers that have been proved to discharge xenobiotic contaminants to the sea (e.g., chromium-enriched sediments found in river estuaries) 91 . The Palermo area hosts a variety of industrial, agricultural, commercial and above all harbor-related activities, with the latter considered the main local source of contaminants for the coastal zone 91,92 . Increasing levels of contaminants (e.g., polycyclic aromatic hydrocarbons -PAHs), far exceeding national and international regulatory guidelines have been detected in the area 93 . The population density of coralligenous species, including corals, decreases with increasing pollution 94 and in general pollution due to metals and organic compounds is detrimental for the marine biota 95,96 . For example, orthophosphate ions (used in fertilizers) inhibit calcification 97 , which may be heavily detrimental for stony corals. Thus, toxic metal compounds or fertilizers that reach the seawater could have affected the fitness of A. calycularis at the more impacted locations, contributing in part to explain the observed patterns. Further studies are needed to test these hypotheses.
Whatever the cause, the overall effect of increasing human impact is the reduction of the species performance at all the considered organizations levels. The regeneration potential and adaptability of corals 98,99 are reflected in their ability to display, in some cases, a mixed mode of sexual and asexual reproduction 100    www.nature.com/scientificreports www.nature.com/scientificreports/ reproduction can help organisms to adapt to changing environments 101 and can be seen as a "best-of-both-worlds" scenario. Generally speaking, in a stable environment asexual reproduction is employed to exploit suitable conditions for survival as all the new individuals will be adapted to that environment 102 . However, when environmental conditions are less favourable sexual reproduction is expected, as it ensures a mixing of the gene pool, leading to genetically diverse offspring that can exploit a broader spectrum of resources and adapt to a broader range of environmental conditions 102 . Fecundity of A. calycularis (i.e., the number of mature oocytes produced per female polyp) is known to decrease with decreasing polyp and colony size 43 . Consequently, we can hypothesize that at the less favourable sites A. calycularis could allocate more energy toward sexual rather than asexual reproduction, as expected from a species facing adverse environmental conditions. Moreover, the limited amount of offspring produced by smaller colonies and polyps at disturbed locations could not be sufficient to increase the percent cover to a level comparable to less impacted sites. This might contribute to explain the fragmented pattern of distribution observed in Loc 4 and Loc 5 (Fig. 2), with small and dispersed colonies and low percent cover, as opposed to the larger colonies found in Loc 1 and Loc 2 that colonize the available space by polyp budding, resulting in a higher percent cover. Further studies on the ecology and reproductive biology of A. calycularis, also during its earliest life stages, are needed to test these hypotheses.  study has shown that population abundance and biometry at the colony and polyp levels in this species are severely impacted by anthropogenic pressure. Given the importance of this species as a bioconstructor and habitat former, and as an attraction to recreational divers, who represent an important source of income for coastal tourism, the current vulnerability status of A. calycularis should be re-assessed, and dedicated management plans should be considered. In the Egadi Archipelago diving is an important tourist activity with eight diving centers in the area. Diving is common between May and October with a peak between July and August with about 40 dives per day in Favignana, 10 in Levanzo and 80 in Marettimo, with 85-90% of divers coming from outside Sicily. Images of A. calycularis are iconic symbols of the diving attractions in the websites of those diving centers (source: the Egadi Archipelago diving centers).
Human pressure on coastal areas is increasing, with severe impacts on marine biodiversity 3 . Thus field studies like this one, providing baseline information on the status of endangered species and the effect of anthropogenic impacts on them, are crucial to achieve proper information for their management and conservation. The shallow-water sites hosting A. calycularis bioconstructions are generally small coves and cliffs not larger than a few hundred square metres 32 . Properly managed marine protected areas have strong potential to prevent damage www.nature.com/scientificreports www.nature.com/scientificreports/ from intense use of coastal habitats. However, as already suggested 32 , the creation of microprotected areas allowing specific restrictions and ensuring proper communication to the large public could be the best practice to achieve the conservation of coastal bioconstructions. Non protected small coves and cliffs where the orange coral usually thrives appear to be good candidates to test and implement focused protection schemes.

Data Availability
The datasets generated and/or analysed during the current study are available from the corresponding author on reasonable request.