Unfolding Jellyfish Bloom Dynamics along the Mediterranean Basin by Transnational Citizen Science Initiatives

Abstract

Science is addressing global societal challenges, and due to limitations in research financing, scientists are turning to the public at large to jointly tackle specific environmental issues. Citizens are therefore increasingly involved in monitoring programs, appointed as citizen scientists with potential to delivering key data at near to no cost to address environmental challenges, therein fostering scientific knowledge and advising policy- and decision-makers. One of the first and most successful examples of marine citizen science in the Mediterranean is represented by the integrative and collaborative implementation of several jellyfish-spotting campaigns in Italy, Spain, Malta, and Tunisia starting in 2009. Altogether, in terms of time coverage, geographic extent, and number of citizen records, these represent the most effective marine citizen science campaigns thus far implemented in the Mediterranean Sea. Here, we analyzed a collective database merging records over the above four countries, featuring more than 100,000 records containing almost 25,000 observations of jellyfish specimens collected over a period of 3 to 7 years (from 2009 to 2015) by citizen scientists participating in any of the national citizen science programs included in this analysis. Such a wide citizen science exercise demonstrates a valuable and cost-effective tool to understanding ecological drivers of jellyfish proliferation over the Western and Central Mediterranean basins, as well as a powerful contribution to developing tailored adaptation and management strategies; mitigating jellyfish impacts on human activities in coastal zones; and supporting implementation of marine spatial planning, Blue Growth, and conservation strategies.

1. Introduction

Jellyfish have been acknowledged as a natural component of marine ecosystems since antiquity, as the fossil record shows that gelatinous blooms were already occurring several hundred millions of years ago [1,2]. There is evidence to suggest that jellyfish blooms may occur following decadal cycles [3,4], presenting interannual fluctuations over climate-related cycles [5,6,7,8,9,10], and that their blooms may be driven by physical forcing in the sea [11]. In the past years, evidence also suggests that jellyfish blooms are fostered in many coastal marine ecosystems over the world by increasing human-related pressures, such as overfishing, eutrophication, climate change, habitat modification (“ocean sprawl”), and species translocation (reviewed in [7,12,13,14,15]).

At a global scale, there is a paucity in historical datasets covering large temporal and spatial scales, such that it is not possible to conclude whether a global increase of gelatinous zooplankton has occurred or not [4,16,17]. However, some coastal marine ecosystems may have been undergoing a long-term increase in jellyfish biomass proliferation in recent decades (e.g., the Mediterranean Sea [4,16,18,19,20]). This lack of data on gelatinous zooplankton is partly due to the fact that jellyfish were not targeted species in fisheries or oceanographic research, as they were poorly known by scientists and considered ecologically unimportant. Another reason behind such data deficiency is that the fragility of most gelatinous species and their patchy distribution makes conventional plankton nets inadequate for sampling this taxon: gelatinous organisms may be damaged beyond recognition, therefore critically limiting the resulting information on species composition, abundance, and distribution [21,22]. Hence, recent surveys and research campaigns on jellyfish also relied on non-traditional data collection methods, including citizen science, to assess patterns of jellyfish diversity, abundance, seasonality, and distribution (e.g., shoreline surveys and stranding observations, aerial surveys, or interviews to stakeholders, among others [22,23,24,25,26,27]).

In recent years citizen science has expanded the scope of data collection on the presence/absence of jellyfish, either in open waters or stranded along coastal areas [28,29,30,31,32]. Long-term and broad-scale data are necessary to determine jellyfish population dynamics, and citizen science is considered as one of the best methods to contribute beyond the local scale. Citizen science programs are, in fact, generally cost-effective and allow for the establishment of monitoring programs at broader geographic scales and for longer periods [33,34,35,36]. Moreover, citizen science has the potential to increase biodiversity awareness and contribute to understanding of spatial and temporal distribution of species, especially when their biogeographical ranges are broad [37,38].

The best example of marine citizen science in the Mediterranean is represented by the integrative and complementary implementation of the CIESM JellyWatch Programme (http://www.ciesm.org/marine/programs/jellywatch.htm [accessed on 16 June 2021]), the “Occhio alla Medusa–Meteomeduse” campaigns in Italy, the “iMedjelly” project in Spain, the “Spot the Jellyfish” campaign in Malta, and the “Decouvrons Les Meduses” citizen science campaign in Tunisia from 2008 to 2020. Altogether, in terms of time coverage, geographic extent, and number of citizen records, these represent the most effective marine citizen science campaigns thus far implemented in the Mediterranean Sea.

On the other hand, jellyfish blooms have been causing detrimental impacts to many coastal ecosystems [12,13,39,40,41], as well as to the conduction of several human activities and to several socially important sectors such as fisheries, aquaculture, tourism, public health, and coastal industries [7,12,13,14,42,43,44]. All these impacts translate into negative effects on marine ecosystem services and important economic losses within the affected human activities/enterprises [14,42,45,46,47,48,49].

As a result, the development of adaptation and management strategies to address jellyfish blooms is imperative in order to mitigate the costs and effects of a massive jellyfish presence within coastal areas [49,50]. Management of jellyfish blooms need to address various disparate aspects, including information and education campaigns, citizen science-based programs, systems for the removal/recovery of jellyfish in extreme cases of infestation, short-term beach closures, usage of temporal anti-jellyfish nets, training of lifeguard personnel and medical aid provision in case of stings, monitoring of stranded jellyfish, and different methodological alternatives for forecasting the occurrence of future blooms [50,51,52,53]. In fact, a forecasting capacity can become a great management tool in the form of an “early warning” system [54,55,56]. This allows coastal managers to identify high-risk days and to be aware of the possible arrival and occurrence of jellyfish blooms in coastal areas, enabling them to take action and mitigate the potential negative impact of jellyfish [53]. This capacity is especially important in coastal areas afflicted by life-threatening species where reliable forecasting capabilities can counsel the closure of beaches in order to prevent stings during high-risk days [47].

Within this context, the present work focused on compiling and analyzing the jellyfish abundance data originating from integrated, collaborative, citizen science-based programs conducted in four Mediterranean countries (Italy, Malta, Spain, and Tunisia), in order to understand the dynamics of jellyfish blooms at a sub-basin-scale and to find tools to mitigate the impact of jellyfish blooms in coastal areas. To achieve such an objective, we conducted a thorough analysis of the four corresponding databases in order to standardize the data and effort. In addition, the spatial and temporal dynamics and prevalence of the targeted jellyfish species were characterized through a habitat-modeling approach, which was applied to the coastlines of the monitored countries. This was done in order to understand the operational environmental forcing and to predict high probability areas for massive beach stranding of jellyfish. The jellyfish bloom forecasting system output is presented as a preventive and mitigation tool for citizens and coastal stakeholders, aiming to reduce the jellyfish blooms socio-economic impact in coastal areas through a feasible and powerful management strategy.

2. Materials and Methods

2.1. Study Area and Jellyfish Citizen Science-Based Programs

The study was conducted along the coastlines of the central and western Mediterranean within the framework of the ENPICBC MED MEDJELLYRISK Project, where participatory (citizen science-based) data on the occurrence of jellyfish species was gathered in four countries: Italy, Malta, Spain, and Tunisia (Figure 1).

To ensure a high degree of accuracy and quality of the data being collected, from the early stages, we mass-distributed species recognition tools, such as identification guides and posters (Figure S1), within the MED-JELLYRISK project to all salient coastal stakeholders.

Moreover, training sessions in Spain and Tunisia were held for coastal lifeguard services, and the information was also included within different digital platforms that facilitated access in the field, including mobile apps (Focus Meteo Meduse App and iMedJelly App) and social networks such as Facebook (https://www.facebook.com/Spot-the-Jellyfish-258981024173161 [accessed on 16 June 2021], https://m.facebook.com/Meduses.Tunisie/ [accessed on 16 June 2021], @meteomeduse [accessed on 16 June 2021] in Italy). All relevant information about the species most likely to be encountered along the targeted Mediterranean coastlines (at least 29 species) was included, i.e., a photo or illustration, most recognizable characteristics, frequency and season of appearance, and its stinging potential. These materials were updated whenever necessary, such as the inclusion of additional gelatinous species and the updating of existing information as a result of the release of additional knowledge on the subject. A compilation of all this information was included in an identification guide (available online at: http://193.188.45.233/jellyfish/docs/Englishguide.pdf [accessed on 16 June 2021]) and in a jellyfish sting treatment booklet (available online at: http://193.188.45.233/jellyfish/docs/firstaid.pdf [accessed on 16 June 2021]), specifically developed within the framework of the MEDJELLYRISK project. The devised data collection protocols were simple and consistent within the four targeted coastal regions, with participants being asked to submit jellyfish records containing basic attributes through the on-line (mail, Facebook, app) citizen science report submission form: species observed, abundance (quantitative and categorical scales), date, time, and location of the observation.

Once the citizen science reports were gathered, experts from the different institutions involved in the project validated all data. Data were filtered for eliminating error and bias (e.g., erroneous reports, duplicates, and records with missing information were excluded). The emerging databases were then treated to standardize the sampling effort, the abundance data (by adopting common categories ranging from 0 to 3), the days of observation, and the sampling sites in order to be able to compare the data between the 4 citizen science programs.

In general, citizen science-based data followed two different monitoring approaches. The first approach, followed by Italy and Malta, reported only presence records of jellyfish. Data from Italy were gathered from a citizen science-based program called “MeteoMeduse” coordinated by the University of Salento in Lecce. It covers the shallow and deep seas of the entire Italian coast: the Adriatic, Ionian, Tyrrhenian, and Ligurian seas (Boero 2013). Data analyzed covered the period from 2009 to 2014 and was validated by experts from the University of Salento. Although the dataset was collected from the Italian citizen science-based campaign, some observations were located in waters near other countries (e.g., Slovenia and Croatia), but they were considered in the Italian analysis. Data from Malta originated from a citizen science-based program launched in 2010 called “Spot the Jellyfish”, which was run by the International Ocean Institute (IOI) and the University of Malta and included coastal sightings made between 2011 and 2014 from the entire Maltese archipelago. Experts from the University of Malta validated the data. For these two countries, the total reports equal the number of positive (presence) reports. The 0 (zero) values that represented the absence of jellyfish within coastal areas were not considered for Italy and Malta’s citizen science report databases.

Conversely, the second citizen science-based monitoring approach, followed by Spain and Tunisia, gathered the information for both presence and absence of jellyfish, and thus the 0 (zero) values representing the absence of jellyfish were also considered in these two databases. This information was available due to a higher effort invested by the public in surveying the respective coastal areas. For the evaluation of the medium-term arrivals of jellyfish along the Spanish Catalan coast, people (mainly from rescue services and from coastal municipalities) were trained to implement a citizen science-based monitoring program that has been carried out since 2007. The main objective was to evaluate the presence of jellyfish along the Catalan coast from May to September from 2007 until 2013 through a daily sampling scheme involving a total of 243 beaches. Data from Tunisia for the 2013–2015 period originated from a citizen science campaign established during the project. Experts from the Faculty of Science of Bizerte validated the corresponding data.

Intensity index = (NPR · TR ⁻¹ · NBS ⁻¹) × 100

where NPR represents the number of positive (presence) reports; TR the total reports, including presence and absence (when available); and NBS the total number of beaches surveyed. Higher intensity index values represent conditions of high prevalence of jellyfish in particular coastal areas, since index values do not only consider the positive records over the total records but also incorporate the monitoring effort invested by including the number of beaches surveyed.

2.2. Predictor (Environmental) Data

The oceanographic variables (i.e., SST, salinity, chlorophyll, nitrate, phosphate, oxygen, and temperature) were obtained through the European “Copernicus Marine Service” (CMEMS-http://marine.copernicus.eu/ [accessed on 16 June 2021]). All the oceanographic variables were downloaded in separate datasets, with a resolution of six degrees of latitude, and these were uploaded to the server on a weekly basis. The motu-client framework (developed by Copernicus as a tool for tapping into their database) was used to recall the data when necessary. For each oceanographic variable, the nearest point to a reported stranded jellyfish coordinate was extracted from the dataset to obtain the value inputted within the machine-learning algorithm (see next section). Following Benedetti-Cecchi et al. [57], the distance (in kilometers) from the nearest marine canyon was included as a predictor in the models as a “slope index”, calculated as the distance to the nearest 1000 m depth isobath (except for Malta, where 800 m was the deepest value) for each site using the marmap package in R [58] (Figure S2). The summary statistics for all the predictor (environmental) data is reported in the (Supplementary Materials Table S1).

2.3. Data Analysis

The oceanographic variables used as potential explanatory variables inside machine learning algorithms were analyzed a priori using the Pearson correlation “rho” to avoid collinearity (more important for GAM and ANN models [59]), where pairs of correlations with a rho value higher than 0.7 were discarded from the analysis (Figure S3). These “dropped” variables were those variables with less-known effects or indirect ecologically related effects on the jellyfish ecophysiology [60,61]. The response variable was defined as the number of jellyfish blooms (the sum of categories 2 and 3, as these categories represent incidences of “blooming”) per each sampled beach. The machine learning workflow (averaging protocol and evaluation) closely followed the routine proposed by Thuiller et al. [62]. Different machine learning methods were applied to the (participatory) jellyfish data, so as to model the environmental forcing over selected species. In addition, in order to identify some potential areas of massive jellyfish occurrence and/or stranding and as a way to create a management tool to assist in the mitigation of jellyfish bloom negative impacts (i.e., people get stung), an ensemble potential distribution “predictive” map was also generated and made available through the iMedjelly App. The algorithms used in the machine learning approach were generalized additive models (GAM), artificial neural networks (ANN), and random forest (RF). All these models were fitted using the ensemble platform for species distribution modeling “biomod2” package [63].

The generalized additive model (GAM) works in a similar fashion to a linear model; however, in GAMs, the linear predictor involves a sum of smooth functions of covariates (explanatory variables), allowing this method to model (fit) more complex relationships between the response (jellyfish data) variable and the selected covariates [64]. This flexibility is achieved since the relationships are defined as “smooth functions”, rather than detailed parametric relationships, reducing the over-fitting problem. Moreover, the “generalized” part of the model allows for a family distribution of the errors, which is different from the Gaussian distribution, making it ideal for non-negative and/or continuous data, including, as in this case, count-like data. A second approach considered here was the use of a different family of models based on the decision trees, in particular the random forest (RF) approach, which is an alternative to cope with the non-linear relationships commonly found in natural systems. This model is based on the generation of a large number of regression trees as a (random) subset of the original data by means of permutation [65]. For the prediction of the response value (abundance category), the covariates were included in each tree, producing a set of response values where the final value is taken as the average value among the trees considered in the forest [66]. Finally, we used an artificial neural network (ANN) statistical modeling approach. This model is based on the analogy to the human brain’s way to learn [67] from the relationships between an input layer, hidden layers, and output layers. At each layer, the inputs values were processed independently in nodes. The output is reached by iteratively connecting two pairs of nodes, with the back-propagation algorithm being one of the most used in the training process. In general, ANN has proved to be robust for noisy data given their ability to handle linear and nonlinear relationships [68]. In order to compare modeling output, we retained only significant variables, and the importance of each environmental variable was calculated through a permutation process, where the model was fitted by excluding (in a stepwise fashion) one environmental variable at a time and calculating how decrease in accuracy when the variable is excluded. Thus, when an “important” variable is excluded, the model’s accuracy decreases to a higher degree than when a “non-important” variable is removed. The estimation of variable importance is given for a zero-to-one (0–1) range, where the zero represents a variable, which has no influence, and the value 1 represents a variable, which is essential.

3. Results

3.1. Jellyfish Dynamics

Jellyfish bloom data varied in the way it was collected and in its temporal coverage. In the case of Spain and Tunisia, all the data collection was conducted by trained people, and reports were submitted both when jellyfish were present (positive reports) but also when they were absent (negative reports). On the other hand, reports from Italy and Malta came from on-line citizen science platforms where people submitted only images of jellyfish species encountered within the coastal area, equivalent to positive records only (without considering the negative = absent records). In order to use this data, experts from the different institutions involved in the project validated the sightings. Temporal coverage of the jellyfish bloom data varied across countries, where Spain, Italy, Malta, and Tunisia conducted seven, six, four, and three years of data collection, respectively. Similarly, the average (per year) number of sampling days and surveyed beaches for Italy were 194.5 and 965, and 111 and 193 for Spain, respectively. Malta and Tunisia surveyed, to a much lesser extent, for an average of 89.8 and 48 days, and on 32 and 33 beaches, respectively (

Table 1. Citizen science data summary for all the countries. NDS: number of days surveyed, PPR: percentage of positive (presence) reports, NBS: number of beaches surveyed, NPR: number of positive (presence) reports, TR: total reports (presence + absence). TOTAL represents the sum of variables, except for PPR and intensity index, which represents the average (and standard error) values (* --: due the small sampling effort of this year, this data were not considered in the analysis).

Country	Year	NDS	PPR	Intensity Index	NBS	NPR	TR
Italy	2009	117	100	0.51	197	303	303
	2010	201	100	0.12	811	1818	1818
	2011	191	100	0.14	720	3738	3738
	2012	217	100	0.04	2387	2671	2671
	2013	214	100	0.09	1130	5422	5422
	2014	227	100	0.18	542	1779	1779
Malta	2011	72	100	3.13	32	124	124
	2012	128	100	2.94	34	242	242
	2013	84	100	3.57	28	122	122
	2014	75	100	3.03	33	134	134
Spain	2007	101	1.3	0.01	243	204	15,762
	2008	101	2.9	0.01	243	456	15,729
	2009	107	4.9	0.02	242	773	15,797
	2010	108	4.2	0.02	243	709	16,907
	2011	123	30.2	0.38	79	1509	4995
	2012	103	43.5	0.43	101	2982	6863
	2013	134	5.1	0.03	197	1135	22,320
Tunisia	2013	50	16.5	0.50	33	66	400
	2014	46	25.5	0.80	32	61	239
	2015	1	-- *	-- *	2	2	2
TOTAL	--	2400	7329	115,367	7329	24,250	115,367

). A linear model showed the positive effect of the increase in the number of days surveyed over the spatial coverage of the reports, visualized as the number of beaches surveyed (r-squared = 0.6, t-value = 4.73, p-value < 0.01), but this was not consistent across countries (Figure S4).

3.1.1. Jellyfish Temporal Dynamics

In the case of Italy, data from MeteoMeduse included reports from 2009 until 2014, comprising a total of 15,731 jellyfish presence records. No absence data were registered, since data originated from citizens that are accustomed to report only the presence of jellyfish; hence, the percentage of positive reports was always 100% (Table 1). During the 2011–2013 period, Italy reported the highest number of jellyfish sightings, with a maximum of 5422 reports for 2013. Data from Malta also refers only to the occurrence/presence of jellyfish individuals, originating from the 2011–2014 period from the entire Maltese archipelago for a total of 622 records, with a maximum of 242 observations of jellyfish reported in 2012 (Table 1). Data from Spain considered reports covering the 2007–2013 period, collected exclusively from the Catalonian coast. All the reports (98,373) were made during the summer season, of which only 7.9% were positive reports (Table 1). The year with the highest number of jellyfish presence records was 2012 (2982 positive reports). For Tunisia, most citizen science records came from fishers. Tunisian data covered the years 2013 and 2014, with a total of 639 records collected, from which 21% corresponded to jellyfish presence and 80% to jellyfish absence reports (Table 1).

The relationship between the positive reports and the total reports received, considering the monitoring effort as the total surveyed beaches, was analyzed through the intensity index. This index showed that among the four targeted countries, Malta had the highest values (3.57), reflecting the higher frequency of jellyfish reporting on its beaches during the study period, especially during 2013. Even though its dataset was the shortest in duration, Tunisia showed the second level of intensity (0.80), followed by Italy (0.51) and Spain (0.43) (Table 1). For Italy, 2009 was the most important year in terms of intensity of jellyfish sighting reports, with a clear decline in the next year until a minimum was registered in 2012, after which a positive (but weak) trend is evident. In the case of Spanish coastal areas, years 2011 and 2012 were the two most important years in terms of the intensity of jellyfish reports.

The jellyfish species reported from the four target countries included the scyphozoans Rhizostoma pulmo, Pelagia noctiluca, Cotylorhiza tuberculata, Aurelia spp., and Chrysaora hysoscella, as well as the hydrozoans Velella velella, Aequorea forskalea, and Physalia physalis. As far as Aurelia spp. is concerned, it should be noticed that, following Scorrano et al. [69], the only open-water Aurelia species recorded in the Western and Central Mediterranean (including the Adriatic Sea) was A. solida Browne 1905, whereas the Aurelia moon jellyfish occurring in sheltered coastal lagoons and harbors (such as Varano lagoon in Italy, or Empuriabrava in Cataluña) should be referred to as Aurelia coerulea von Lendenfeld 1884.

Within the target countries, different species dominated the report cohort; in general, for those countries in the northern half of the Mediterranean (Italy and Spain), R. pulmo was the most common species reported, followed by P. noctiluca. Conversely, reports from central Mediterranean countries (Malta and Tunisia) were dominated by P. noctiluca, with few (or no) sightings of R. pulmo.

An analysis of the temporal dynamics of jellyfish reports per country and per species showed that, for Italy, frequency of reports for the species Rhizostoma pulmo increased consistently from 2009 to 2013 when they decreased to minimum levels, even lower than those reported in 2009. Reports of the species for abundance category 1 were the most abundant, followed by categories 2 and 3. This observed hierarchy in the prevalence of different abundance categories was evident for different species, with the exception of Velella velella, for which abundance category 3 was the most commonly reported abundance category since 2010, showing a progressively increasing trend in abundance until the end of the dataset (Figure 2). For P. noctiluca, reports showed two periods of high abundance, in 2010 and in 2013.

In the case of Malta, reports were broadly dominated by P. noctiluca, for which a maximum number of records were reported during 2012, with abundance category 2 being as almost as commonly reported as abundance category 1 (Figure 3). A similar increase in the number of reports for V. velella was evident for Italian coastal areas (Figure 2).

Along the Spanish (Catalan) coast, the reports were dominated by R. pulmo, whose reports increased from 2007 to a maximum in 2012 when the abundance category 3 reached a maximum level, with more than 100 reports (Figure 4). However, numbers seemed to decrease during 2013, when no report was submitted for the abundance category 3 (Figure 4). The second most commonly reported species (P. noctiluca) showed a maximum during 2008, with a subsequent decrease in the number of reports until a minimum being registered during 2012. Afterwards, during 2013, the reports of this species seemed to recover again, at least for abundance category 1, which reached the maximum number of reports registered during the entire 7 years (Figure 4). For 2009 and 2010, the reports were broadly dominated by V. velella and A. forskalea, two species that were almost absent within other years (Figure 4).

For the Tunisian coastal areas, the reports were dominated by P. noctiluca, which showed an increase in the number of records for abundance categories 2 and 3 from 2013 to 2014; however, the species was not reported during 2015. Aurelia (cf. A. solida; see also Gueroun et al. [70]) was a commonly reported species for Tunisian waters, with a similar abundance pattern being reported for 2013 and for 2014, with no record of the species being made during 2015. The few reports submitted for R. pulmo showed a marked decrease in category 1 since 2013, with a peak of abundance being registered during 2014, the only year when the abundance category 2 was registered (Figure 5).

3.1.2. Jellyfish Spatial Dynamics

From a spatial (horizontal) perspective, the distribution of the “low” (1) abundance category, which represents an abundance of less than 10 individuals of jellyfish per beach, showed no clear patterns. The “medium” and “high” (categories 2 and 3) abundance categories showed similar horizontal patterns, and thus their values were combined within analysis and were used to represent jellyfish massive sightings and strandings. For the two most commonly reported jellyfish species (i.e., Rhizostoma pulmo and Pelagia noctiluca), the horizontal distribution of the highest abundance categories showed that, for Italy, P. noctiluca blooms were persistent throughout the study period and that the most affected areas were those in the southern and northern Tyrrhenian Sea, the Ligurian Sea, and also the central Adriatic, which contributed several reports of P. noctiluca blooms (Figure 6). R. pulmo blooms were mainly reported from the Ligurian, the north Tyrrhenian Sea, and the north Adriatic Sea, and even occasionally from inside the Gulf of Taranto (north Ionian Sea) (Figure 7).

For Malta’s coastal areas, P. noctiluca blooms were reported from all around the island, with a higher number of reports being submitted for NW and SE coastal areas (Figure 8). R. pulmo was largely absent from these waters, with the exception of a few records made in 2013.

In the case of Spain, the greatest abundance of P. noctiluca and R. pulmo showed contrasting patterns of distribution. P. noctiluca showed its greatest abundance along northern sections of the coast, showing a clear association with the presence of submarine canyons. Along southern parts of the Catalan coast, high numbers of P. noctiluca were recorded to the north of the Ebro River delta, where circulation eddies are common (Figure 9). On the other hand, R. pulmo showed its greatest abundance along the central area of the Catalan coast, near the city of Barcelona, especially for 2011 and 2012, when records for this species were dominated by the high abundance category (Figure 10).

In the case of the Tunisian coastal areas, the northern part (near Bizerte) was the area from where most reports originated. Thus, blooms of P. noctiluca (Figure 11) and R. pulmo (Figure 12) were mostly recorded in this area.

3.2. Analysis for Jellyfish Potential Spatial Distribution

Environmental variables affect jellyfish distribution by altering their reproductive cycle, changing their growth and survival time (affecting the amount of time being drift with currents), and by causing “aggregations” of jellyfish [11] where hundreds of individuals are concentrated in some places by dispersal restriction and currents.

The modeling approach (e.g., machine learning algorithms and ensembles), considering the relationships (linear, non-linear, and interacting effects) of environmental variables over the observed distribution of stranded jellyfish, revealed a high degree of inter-specific and inter-country variation. In order to analyze which variables have a higher important role affecting the spatio-temporal dynamics of jellyfish, we focused on the average variable importance among the three machine learning algorithms. Thus, for the species Pelagia noctiluca, the most important (average value > 0.4) variables in explaining the observed distribution of the species were sea surface temperature (SST), salinity, wind direction, and slope index (

Table 2. Variable importance estimates for each explanatory variable for Pelagia noctiluca distribution models. GAM: generalized additive model, ANN: artificial neural network, RF: random forest. Those variables with an importance estimate > 0.4 are highlighted in bold.

Country	Environmental Variable	GAM	ANN	RF	All Models (Average)
Italy	Sea surface temperature (°C)	0.08	0.09	0.35	0.17
	Salinity	0.62	0.24	0.45	0.44
	Current direction (East, m s−1)	0.08	0.00	0.08	0.05
	Current direction (North, m s−1)	0.06	0.00	0.10	0.05
	Current speed (m s−1)	0.08	0.05	0.12	0.08
	Slope index	0.29	0.93	0.66	0.62
Malta	Sea surface temperature °C)	0.13	0.93	0.23	0.43
	Salinity	0.09	0.01	0.15	0.08
	Phosphate (mmol m−3)	0.50	0.00	0.02	0.17
	Nitrate (mmol m−3)	0.07	0.00	0.11	0.06
	Current direction (East, m s−1)	0.02	0.01	0.11	0.05
	Current direction (North, m s−1)	0.09	0.00	0.08	0.06
	Slope index	0.02	0.36	0.09	0.15
Spain	Sea surface temperature (°C)	0.39	0.78	0.63	0.60
	Salinity	0.17	0.07	0.60	0.28
	Chlorophyll a (mg m−3)	0.04	0.07	0.34	0.15
	Wind direction (East, m s−1)	0.14	0.51	0.41	0.35
	Wind direction (North, m s−1)	0.19	0.62	0.59	0.47
	Current direction (East, m s−1)	0.36	0.01	0.35	0.24
	Current direction (North, m s−1)	0.26	0.01	0.40	0.22
	Current speed (m s−1)	0.51	0.00	0.36	0.29
	Slope index	0.18	0.55	0.47	0.40
Tunisia	Sea surface temperature (°C)	0.54	0.28	0.42	0.41
	Chlorophyll a (mg m−3)	0.23	0.00	0.26	0.16
	Phosphate (mmol m−3)	0.25	0.00	0.02	0.09
	Nitrate (mmol m−3)	0.48	0.15	0.13	0.25
	Wind direction (East, m s−1)	0.25	0.43	0.05	0.24
	Wind direction (North, m s−1)	0.25	0.26	0.11	0.21
	Wind speed (m s−1)	0.23	0.20	0.06	0.17
	Current direction (East, m s−1)	0.10	0.00	0.02	0.04
	Current direction (North, m s−1)	0.27	0.00	0.04	0.10
	Current speed (m s−1)	0.07	0.00	0.03	0.03
	Slope Index	0.97	0.93	0.45	0.78

). In Italy, salinity and slope index were the two most explanatory variables (0.44 and 0.62, respectively), with all the other environmental variables scoring 0.1 or less in terms of importance estimates. The predicted potential areas for high numbers of jellyfish, considering average environmental conditions, showed important areas near the Gulf of Genoa (Ligurian Sea), the south Tyrrhenian Sea (near the Aeolian Islands), around the island of Sardinia, in the Strait of Sicily, and one area in the central Adriatic Sea (Figure 13A). In the case of Malta, the observed distribution of the species could be mostly explained in terms of SST values, with an average significance estimate of 0.43. The predicted important areas blooming of the species within Maltese coastal waters, in concordance with its intensity index, showed high values around the entire archipelago, with the north and southeast areas scoring the highest values (Figure 13B). In Spain, the variables with the highest average importance values were SST, wind direction (“north wind”), and slope index (0.6, 0.47, and 0.4, respectively). Predicted important areas for future blooming of the species consisted of two main areas: (i) those associated with the presence of submarine canyons (which reflects the influence of high slope index values) along the northern Catalan coast near the Cap de Creus National Park, and (ii) the coastal area contiguous to the northern Ebro river delta (Figure 13C). In Tunisia, the two main explanatory variables were SST and slope index, with values of 0.41 and 0.78, respectively. The predicted areas with high potential for P. noctiluca blooming were all circumscribed to the northwest coastal area near Bizerte (Figure 13D).

For the species Rhizostoma pulmo, the relative importance of the environmental variables in explaining the observed distribution was similar to that observed for the previous species (P. noctiluca), albeit in a different variable combination (

Table 3. Variable importance estimates for each explanatory variable for Rhizostoma pulmo distribution models. GAM: generalized additive model, ANN: artificial neural network, RF: random forest. Those variables with an importance estimate > 0.4 are highlighted in bold.

Country	Environmental Variable	GAM	ANN	RF	All Models (Average)
Italy	Sea surface temperature (°C)	0.03	0.10	0.22	0.12
	Salinity	0.29	0.38	0.41	0.40
	Current direction (East, m s−1)	0.12	0.00	0.07	0.06
	Current direction (North, m s−1)	0.10	0.01	0.10	0.07
	Current speed (m s−1)	0.12	0.05	0.12	0.09
	Slope index	0.42	0.8	0.64	0.62
Spain	Sea surface temperature (°C)	0.40	0.65	0.32	0.46
	Salinity	0.36	0.03	0.45	0.28
	Chlorophyll a (mg m−3)	0.01	0.07	0.20	0.09
	Wind direction (East, m s−1)	0.02	0.08	0.40	0.17
	Wind direction (North, m s−1)	0.07	0.25	0.35	0.22
	Current direction (East, m s−1)	0.06	0.03	0.29	0.13
	Current direction (North, m s−1)	0.05	0.02	0.19	0.08
	Current speed (m s−1)	0.02	0.01	0.31	0.11
	Slope index	0.29	0.68	0.39	0.45
Tunisia	Sea surface temperature (°C)	0.85	0.33	0.24	0.48
	Chlorophyll a (mg m−3)	0.70	0.00	0.16	0.29
	Phosphate (mmol m−3)	0.80	0.00	0.00	0.27
	Nitrate (mmol m−3)	0.80	0.25	0.17	0.40
	Wind direction (East, m s−1)	0.84	0.40	0.06	0.43
	Wind direction (North, m s−1)	0.84	0.48	0.11	0.48
	Wind speed (m s−1)	0.83	0.48	0.55	0.62
	Current direction (East, m s−1)	0.49	0.00	0.16	0.22
	Current direction (North, m s−1)	0.52	0.00	0.12	0.21
	Current speed (m s−1)	0.24	0.00	0.12	0.12
	Slope Index	0.91	0.84	0.37	0.71

). Thus, for Italy, the salinity and slope index were the two most important variables (0.4 and 0.62, respectively). The predicted important areas for blooming of the species showed a different spatial distribution than that for P. noctiluca. In this case, these areas were restricted to the Levante Riviera and to the Tuscan Archipelago (east Ligurian Sea), the central area of the Tyrrhenian Sea, north of Naples, the southern coastal areas of the island of Sicily, the coastal areas of the Gulf of Taranto, the southern part of the Adriatic Sea, and the coastal areas of the Gulf of Venice (Figure 14A). In the case of Spain, SST and slope index were the two most important environmental variables (0.46 and 0.45, respectively) that significantly explained the observed distribution of blooms. The spatial prediction of the future occurrence of blooms of R. pulmo in these waters contrasted with that for P. noctiluca, showing high probability values in the central area of the Catalan coast, mainly south of the city of Barcelona (Figure 14B). In Tunisia, the distribution of blooms of R. pulmo was explained mostly by SST, nitrate, wind (speed and direction), and slope index, variables that featured within all the model algorithms. The predicted blooming high probability areas for R. pulmo in these waters showed a similar pattern to those for P. noctiluca, where the coastal area near Bizerte showed the highest values (Figure 14C).

4. Discussion

At a regional, sub-basin scale, our results demonstrate the potential of citizen science data for the unraveling of jellyfish temporal dynamics, as well as the identification of potential blooming areas of the most abundant jellyfish species, through habitat-modeling. To our knowledge, this is the first time that the spatial and temporal dynamics, the role of the environmental variables, and the characterization of the most probable presence areas have been described for the species P. noctiluca and R. pulmo on a broad scale across the Western and Central Mediterranean, from Italy to Tunisia to Spain.

Citizen science has been widely recognized as a great tool for producing large amounts of data, especially data covering broad spatio-temporal extents, without whose application the corresponding data acquisition and analysis might be constrained or even impossible [33,35,37,71,72,73,74,75]. The effort invested in the current study in collating together the jellyfish sighting databases from four different citizen science programs implemented at different levels in different Mediterranean countries allowed us to generate a standardized, large-scale database. This, in turn, has contributed to addressing the lack of knowledge about jellyfish and their spatio-temporal distribution, at least in the central-western Mediterranean.

One of the big challenges in utilizing data from citizen science programs is the generation of a good-quality dataset [76,77], and the current study was not exempt from such a challenge. However, in the present study, training sessions and/or teaching tools (e.g., species identification guides and posters, mobile apps, and social media platforms) were implemented in order to facilitate participation as data collection, as well as to increase the capacities and experience of volunteers. This early stage strategy helped minimize error and to increase accuracy [34,37,75,77,78]. A thorough data management and validation of every single jellyfish sighting report by experts also contributed to database quality [34,78]. Although the most common jellyfish species in the Mediterranean are easily recognizable, a poor status of the specimen observed (missing parts such as tentacles, only umbrella, or very degraded mesoglea) or a poor photographic quality may generate confusion in the reports, making the validation and the data treatment processes essential to reduce bias.

Citizen science monitoring campaigns can provide data about the presence and/or absence of jellyfish and their relative abundance, which can be a preventive tool in themselves by acting as an early warning system. However, the applicability of the same data can be extended, as in the delivery of predictive models that would act as a management tool [50,55,79]. The analyzed database, featuring more than 100,000 records containing almost 25,000 observations of jellyfish individuals, collected over a period of 3 to 7 years by citizen scientists participating in any of the citizen science programs included in this analysis, is a valuable contribution to the understanding of jellyfish ecology. Using different machine learning methods and the available environmental data, we were able to identify key environmental variables and to develop the species distribution modeling. In so doing, we inferred considerable inter-specific variability in occurrence and a well-defined inter-annual variability in abundance for each species. Our data covers a short period of time (3 years in Tunisia and 7 years in Spain); hence, the observed trends may not be representative of larger (decadal) oscillations, although these results are definitely a key aspect for understanding jellyfish dynamics and a starting point for future related work. Jellyfish are known for their highly variable spatio-temporal dynamics and for the periodicity of their blooms [3,28,80]; only large temporal databases can aspire to identify the correct temporal trends of jellyfish blooms [4,16,81]. Regarding the spatial distribution, the low abundance category showed no clear species-specific pattern, reflecting the observation that the presence of small numbers of jellyfish along the coastal areas is a common phenomenon. The medium and high abundance categories (considered together as “blooms”) were analyzed for P. noctiluca and R. pulmo in all the study areas due to their public hazard importance, with P. noctiluca being considered as the most important Mediterranean jellyfish from the toxic point of view [82].

Pelagia noctiluca is a holoplanktonic species that presents an extensive range of distribution in all warm and temperate waters [83]. Within the Mediterranean, the species inhabits coastal and oceanic waters [84,85] and prominent occurrence all year-round, albeit with interannual differences in intensity, emerging in the current study as a result of the collected citizen science data. Rhizostoma pulmo is a meroplanktonic species widely distributed in the whole Mediterranean basin (reviewed in [86]) and one of the second most abundant species along western Mediterranean coasts [87].

In the Adriatic, P. noctiluca was considered a rare species [88] before its anomalously high blooms recorded in 1977 [89,90,91]. In 2010, P. noctiluca dominated the western basin but it was much less observed from the Adriatic [15]. Over the 4-year period of 2010–2013, P. noctiluca had the greatest abundance along the western coast of Italy, but it was increasingly recorded at bloom levels from the Adriatic Sea. Pelagia noctiluca was still not as abundantly recorded as R. pulmo, especially in the northern Adriatic, where the latter species is considered to be native and is thus observed regularly [18] and where its blooms were reported during the entire duration of the current six-year study. Along the rest of the Italian coast, records of P. noctiluca were notably more abundant than in the Adriatic, especially those belong to the high abundance category. The pattern was similar throughout the six-year study, except for 2009, when the fewest observations were reported. In 2009, P. noctiluca abundance peaked early (January–March), according to studies carried out in the Strait of Messina [92]. Given that citizens’ participation is lower during winter months than during the bathing season, this might explain the lower number of reports for the species during 2009, at least for the Tyrrhenian and Ionian Seas. In the Ligurian Sea, P. noctiluca has been reported as present along most of the coastline but exhibiting alternating years of presence and absence [3,28], in full agreement with data collected in the current study. Regarding R. pulmo, the general distribution pattern shown in the present study is consistent with that reported by Boero et al. [15] and Leoni et al. [93]. For P. noctiluca, one of the environmental variables with highest importance was salinity, which is in agreement with Piccinetti Manfrin and Piccinetti [94] and Canepa et al. [30]. Although previous studies in different areas of the Adriatic suggested a direct correlation between sea-water temperature and the abundance of P. noctiluca [89,95,96], with temperature being identified as a major environmental factor affecting population densities [30,85], this factor was identified the third most important in this study for the same species. Another important explanatory factor was the slope index, corroborating the association between jellyfish stranding events and their proximity to marine canyons reported in previous studies, including the one for the Strait of Messina [92]. In the case of R. pulmo, the same factors as for P. noctiluca were relevant, with SST also being identified as the third one of importance in this study, being correlated with high abundance of Rhizostoma spp. in Europe [97].

One of the scyphomedusae known for considerable coastal aggregations in Malta is P. noctiluca, and since its massive blooms in the early 1980s [98], this species is considered to be one of the most common species around the island. As previous studies have shown (on the basis of the same citizen science data), the species is present throughout the year, along both flanks of the island, disappearing only sporadically in autumn [99]. Over the years 2011–2015, Gatt et al. [99] highlighted the year 2012 as the highest in terms of the number of reports submitted for the most commonly reported jellyfish species, with P. noctiluca having the highest abundance values. Conversely, other species such as Rhizostoma pulmo, cf. Aurelia solida, and C. hysoscella were very rarely reported during the same period. The results of this study showed that the outstanding forcing environmental factor for P. noctiluca was sea surface temperature, which is consistent with the results of previous studies [30,99]. In addition, for this species, hydrodynamic forcing (e.g., winds, current, and tidal effects) has been implicated in the formation of coastal and offshore aggregations in the Adriatic [90,100,101] and in Maltese waters [98]. Although in this study these variables did not emerge as important ones, except for the currents, they may explain the local aggregations in coastal areas of Malta due to the island effect [98,101]. Malta’s island status and the relative importance of sea currents within this context could also be the reason why in Malta the intensity index for P. noctiluca was so high when compared to the rest of the countries included in the study.

The spatial distribution of P. noctiluca along the Catalan coast in Spain covered the entire coastline but was highly concentrated at the extremes, especially in the northern area, as was also previously described in Canepa et al. [30], whose analysis was based on the same data series collected up to 2010. This distribution can be explained in terms of the proximity to the canyons in the north that can enhance circulation, making the seasonal occurrence of this species more likely, as hypothesized in previous works [30]. Moreover, in the current study, the slope index was identified as an important environmental factor, consistent with the findings of Benedetti-Cecchi et al. [57], who reported fewer jellyfish outbreaks with increasing distance from canyons. According to the data collected by Fuentes et al. [102], for the years 2007–2009, based on citizen science reports and supplemented with monthly coastal surveys, R. pulmo was the second most observed jellyfish after P. noctiluca, results that agree with our findings. In the same study, the authors indicated that, by 2010, important impacts caused by this species had been reported by fishermen as a result of more frequent and intense blooms. Within our temporal dynamics analysis, a notable increase in R. pulmo records was observed from the year 2010 to 2013, with a considerable peak in 2012, agreeing also with the results shown in Leoni et al. [93]. The spatial distribution observed in this study is also in agreement with Fuentes et al. [102], showing the occurrence of this species along the entire Catalan coast, with the highest concentration in the central area around the Barcelona province. This is supported by the findings of previous coastal surveys that recorded the presence of ephyrae, mainly in the central area of the Catalan coast [102]. To our knowledge, this is the only work describing the spatio-temporal distribution of R. pulmo together with the role of the environmental variables along Spanish coasts and integrating a seven-year-long database. For both P. noctiluca and R. pulmo, the most important explanatory environmental factors were SST and the slope index. For P. noctiluca, the wind direction was also identified as an important variable, in agreement with the findings of Canepa et al. [30], although in the latter study, the explanatory variable behind the stranding of this jellyfish species was southeastern winds, whilst northwesterly winds were similarly identified in the current study.

Pelagia noctiluca in Tunisian coastal areas was the most sighted species (blooms of high abundance category) and mostly in the northern part of the country, in agreement with the findings of previous studies [103,104,105]. Some of the most important explanatory environmental variables for blooms of this species within this region were SST, slope index, and nitrate, in agreement with the findings of Touzri et al. [105]. Another environmental variable positively associated with the presence of P. noctiluca in the bay of Bizerte identified by Touzri et al. [106] was salinity, a variable that was highly correlated (rho = 0.85) with the slope index. R. pulmo is a characteristic species of Tunisian coasts [106], and blooms are described under conditions of high temperature and coastal eutrophication [107]. However, there is no ample information about the spatio-temporal distribution of this species in Tunisian waters, except from the recent work by Leoni et al. [93] that documents a year-round presence and large blooms for this species, neither of the key environmental variables that determine its presence, elevating the current study to the first one of its kind on R. pulmo within these waters.

Several studies along the Mediterranean have shown that high abundances of gelatinous zooplankton are generally related to variations in water mass hydrodynamic variables, particularly in salinity and sea surface temperature [3,85,103,105,108]. The different Mediterranean seas and coastal areas have their own hydrodynamic and physico-chemical peculiarities, as do different jellyfish species in terms of preferred conditions. Hence, it is fundamental to know those peculiarities in order to plan and implement adequate management measures that are region- and species-specific. In this study, we explored those peculiarities and applied the identified key variables within further analysis, testing a predictive model as a potential prevention and mitigation tool, considering that forecasting systems may be useful for anticipating events and facilitating management [50,55]. Our forecasting models appear to be consistent with the data of P. noctiluca (Figure 13) and R. pulmo (Figure 14) collected through a participatory approach, since coastal areas identified as carrying a high probability of bloom occurrence largely coincide with areas reporting the highest frequency of abundance category 2 and 3 citizen science records.

Jellyfish blooms can have a serious socio-economic impact, especially on the sectors of tourism, aquaculture, and fisheries. Blooming species, such as P. noctiluca and R. pulmo, are among the most important species, as a result of their sheer impact, within the western Mediterranean. As a result, investing in citizen science may be an effective approach to provide useful information about jellyfish blooms and to explore potential mitigation tools for jellyfish coastal management purposes. Similarly, but in light of a Blue Growth vision, the potential to predict areas with dense jellyfish proliferation may turn into a positive perspective for the potential exploitation of commercially valuable jelly biomasses [109,110].

Herein, we demonstrate the importance and usefulness of citizen science programs in order to determine jellyfish dynamics at broad spatial scales. In this case, four countries of the Mediterranean covering the Adriatic, Ionian, Tyrrhenian, Ligurian, and Balearic Seas, were able to characterize the spatio-temporal dynamics of the most commonly occurring jellyfish species, and to propose a mitigation and prevention tool for coastal jellyfish bloom management through the data coming from their national citizen science initiatives. Well-designed, -implemented, and -evaluated citizen science programs, conducted with committed participants, can efficiently generate quality-controlled data, contribute to effective management strategies, and help mitigate negative impacts [33,72,111]. In this framework, the implementation of well-structured citizen science programs should be recommended across large marine regions to increase knowledge on jellyfish proliferation mechanisms and distribution, thus contributing to the development of adaptive management strategies in coastal areas.