Located in the northernmost part of the central Mediterranean, the semi-enclosed Adriatic Sea is divided into the northern, central and southern Adriatic. The Jabuka Pit depression (280 m depth) separates the northern shallow shelf (depths up to 80 m) from the deeper central Adriatic, while the northwest perimeter of the South Adriatic Pit (1200 m depth) separates the central and the southern Adriatic ( Figure 1 ).

The most important physical parameter affecting lobster abundance, and thus landings and CPUE, is sea bottom temperature (Green et al., 2014; Zhao et al., 2019). Throughout the Adriatic, this is dominantly affected by bathymetry and seasonal changes of heat flux (Artegiani et al., 1997). Over the shallow northern Adriatic, temperatures vary from 6°C during severe winter cooling events to 20°C during periods of vertical mixing in the autumn. Deeper regions of the central Adriatic exhibit less seasonal changes with temperatures ranging from 13 to 17°C (Buljan and Zore-Armanda, 1976; Lipizer et al., 2014), and within deepest regions of southern Adriatic temperatures vary even less, from 12.4 to 13.7°C (Lipizer et al., 2014; Cardin et al., 2020).

Temperature data were retrieved from measurements and the numerical model Regional Ocean Modelling System (ROMS) (Shchepetkin and McWilliams, 2005; Shchepetkin and McWilliams, 2009). Temperature measurements were collected at the stations Stončica (depth: 95 m) in the central Adriatic (43°0’0”N, 16°20’0”E) and SJ107 (depth: 35 m) in the northern Adriatic (45°2′52″N, 13°19′0″E) in the period 2008–2017 ( Figure 1 ). These stations are parts of two oceanographic transects monitored for at least a half of century: Stončica by the Institute of Oceanography and Fisheries (Split) and SJ107 by Institute Ruđer Bošković (Rovinj) (Marić et al., 2012; Vilibić et al., 2012). The stations are highly representative in term of depth for the northern and central Adriatic. Also, the depths at which the measurements were performed correspond to the typical inhabitation and catch depths of the analysed lobster species in the northern (20–35 m) and central (50–100 m) Adriatic. Temperature measurements were mostly carried out once a month or once every two months. At Stončica measurements were performed in 64%, and at SJ107 in 75% of months within the study period. Measurements were taken at Stončica using CTD probes (IDRONAUT 316, SeaBird-25 and 911+) with an accuracy of ±0.003°C, and at SJ107 using protected reversing thermometers (Richter and Wiese, Berlin; precision ±0.01°C) and reversing digital thermometers (SiS RTM 4002; precision ±0.003°C) attached to Niskin bottles. For this research, we used the deepest sampling, which was usually taken 2 m above the seabed.

Daily values of the ROMS were used to reproduce temperature changes in the Adriatic Sea between 2008 and 2017. The ROMS model is a 3D hydrostatic, nonlinear, free surface, sigma coordinate, time splitting finite difference primitive equation model (Shchepetkin and McWilliams, 2005; Shchepetkin and McWilliams, 2009). The lateral boundary conditions were taken from the AREG model of the Adriatic Forecasting System (AFS) (Oddo et al., 2006), while atmospheric forcing was prescribed via bulk formulation (Fairall et al., 1996), using all the required variables from the operational local area model ALADIN/HR (Tudor et al., 2013). Horizontal resolution was 2 km, with 20 sigma vertical layers in the model. River discharges were introduced to the ocean model following climatology by Vilibić et al. (2016) for all rivers except the River Po, for which real daily discharges were used. All other details about the model setup are provided by Janeković et al. (2014) and Vilibić et al. (2016). Temperatures modelled at the lowest sigma coordinate, roughly corresponding to sea bottom temperatures, were analysed in the paper.

Additionally, bottom temperature averages were calculated for each fishing zone and year/season by considering each model point within the fishing zone polygons. Trend significance of sea bottom temperature was tested by the Mann-Kendall nonparametric test.

For management purposes and data collection, the Croatian marine fishing area has been administratively divided into smaller units (11 fishing zones; OG, 2011) ( Figure 1 ). All available data on H. gammarus landings in the Adriatic Sea for fishing zones A, B, C, D, E, F, G (coastal zones) and H, I, J, K (offshore zones) for the period 2008–2017 were obtained from the Fisheries Directorate (Croatian Ministry of Agriculture) based on fisher’s logbooks (fishery dependant data). These data correspond with open season for lobsters (May–August) and the MED EU minimum legal size of 105 mm carapace length (CL) (EU Regulation 1967/2006). During the study period, insignificant landings (less than 50 kg per year) were reported in fishing zones H, I, J and K and therefore these zones were excluded from further analysis. For analysis purposes, data were standardised as the catch per unit effort (CPUE), expressed as the biomass of H. gammarus caught per fishing trip of a single fisher. For this study, the available data at the scale of fishing zones were aggregated to three general regions: northern, central and southern following the natural geographical division of the Adriatic Sea. Each fishing zone has boundaries expressed by geographical position. Thus, fishing zones A, B and E corresponds to the northern, zones C, F and G to the central, and zone D to the southern Adriatic. As with temperature, the trend significance in H. gammarus landings, fishing effort and CPUE was tested with the Mann-Kendall nonparametric test. Pearson correlation analysis was used to test the linear relationship between sea bottom temperature and CPUE (kg/day) data of H. gammarus.

A total of 22.83 t of H. gammarus was landed on the eastern Adriatic coast during the study period (data of the Fisheries Directorate for 2008–2017), distributed by region as follows: northern 13.93 t, central 8.65 t and southern Adriatic 0.25 t. Particularly, H. gammarus was mostly captured by pots and gillnets in the fishing zones A and E in the northern Adriatic ( Table 1 ). In that area, as shown in Figure 2A , landings averaged over three fishing zones (A, B, E) fluctuated from 0.21 t (2008; landings = 631.5 kg; CPUE = 2.5 kg/day) to 0.75 t (2016; landings = 2244 kg; CPUE=3.41 kg/day). Landings averaged over the three fishing zones corresponding to the central Adriatic (C, G, F) ranged from 0.23 t (2010; landings = 685.7 kg; CPUE = 1.77 kg/day) to 0.37 t (2014; landings = 1112 kg; CPUE = 2.18 kg/day). In the same period, landings in the southern Adriatic (zone D) were less than 50 kg per year, namely from 14 kg (2008; 10 fishing days corresponding to CPUE of 1.4 kg/day) to 41.1 (2010; 25 fishing days and CPUE of 1.6 kg/day) and therefore the data from that zone, i.e. the southern Adriatic, were not included in the further analysis. Temporal trends (2008–2017) of H. gammarus landings in the Adriatic Sea indicate a statistically insignificant and weak increase in the central Adriatic (4.5 kg/year; p > 0.05; R² = 0.008) and a statistically significant increase in the northern Adriatic (160.9 kg/year; p < 0.05; R² = 0.755).

The temporal trend across the study period (2008–2017) of fishing effort ( Figure 2B ) was weakly negative and statistically insignificant (-3.1 days/year, p > 0.05; R² = 0.021) in the central Adriatic, while fishing effort was positive and statistically significant (40.6 days/year, p < 0.05; R = 0.728) in the northern Adriatic. As a consequence, the decennial temporal trend of H. gammarus of catch per unit effort followed a similar pattern as landings: the CPUE trend over time in the central Adriatic was weakly positive but statistically insignificant (0.03 kg/effort, p > 0.05, R² = 0.307) while CPUE was positive and statistically significant (0.09 kg/effort, p < 0.05, R² = 0.643) in the northern Adriatic ( Figure 2C ). The negative trend of fishing effort in the central Adriatic affects the CPUE calculation in the same area and should be carefully considered.

Time series of the yearly averages of measured bottom temperatures at the SJ107 and Stončica stations and modelled bottom temperatures in the selected fishing zones are shown in Figure 3 . The yearly averages of measured temperatures at SJ107 for 2012 and 2017 were not calculated, since data were lacking at SJ107 for February to June 2012, i.e. during the part of year characterised by the lowest bottom temperatures, while in 2017, measurements are lacking for September to December, i.e., the part of the year characterised by the highest bottom temperatures. Including these values would thus result in a significant overestimation (or underestimation) of the measured yearly averages for these two years.

The bottom temperature trends at both stations (0.15°C/yr, p < 0.05 at SJ107; 0.14^°C/yr, p < 0.05 at Stončica) and over fishing zones were positive and statistically significant. The temperature trend is apparently governed by a pronounced jump of bottom temperatures values, starting in 2013 ( Figure 3 ). This jump is evident in both the measurements and the model, in both the central and northern Adriatic, though it is more pronounced in the northern Adriatic (~1.5-2.5°C) than in the central Adriatic (~1°C). Figure 3 clearly shows that the model is biased when it comes to reproducing absolute values of bottom sea temperature, in particularly in the shallow northern Adriatic ( Figure 3A ). Nonetheless, the reproduction of the variability of yearly ( Figures 3C, D ) and seasonal ( Figure 4 ) temperature was satisfactory for both areas, and thus we chose to use the model for more detailed analysis.

The seasonal changes of near-bottom temperature at SJ107 (34 m) and Stončica (95 m) are shown in Figure 4 . These changes were more pronounced in the northern Adriatic, where their range reaches ~7°C, than in the deeper central Adriatic (~3°C in range). The northern Adriatic is much colder during most of the year, particularly during the winter period (February-March). October is the only month in which the sea bottom temperature is higher in the deep northern Adriatic than in the deep central Adriatic.

The modelled values were overestimated at SJ107 by 0.6-0.9°C during the winter months and even more during the period of developed thermocline (May through October; up to 4.0°C). In contrast, the model slightly underestimated the bottom temperature at the Stončica station (model-to-measurements bias is up to -0.8°C in June). The offset of model values might be due to systematic offsets of values of atmosphere-ocean heat fluxes, or of lateral boundary conditions propagating from Otranto towards the central and northern Adriatic. Also, it might be due to an inadequate reproduction of vertical mixing processes. Nonetheless, it should be noted that the phase of the seasonal signal and its variability and range are well reproduced. Standard deviations of bottom temperatures are consistent between model and measurements, and are higher in the northern Adriatic (~1°C) than in the central Adriatic (~0.5°C). The model best reproduced the JFM (January-February-March) and AMJ (April-May-June) seasonal bottom temperatures. This is particularly important for this study, as JFM and AMJ were found to be the two most relevant periods in which changes can affect lobster landings and CPUE. JFM and AMJ also correspond to periods preceding and during the open fishing season.

Correlations between bottom sea temperature and CPUE (kg/day) for the fishing areas (Northern Adriatic - zones A, B and E; central Adriatic - zones C, G and F) for the period from 2008 to 2017 were found to be statistically significant both for winter (JFM) and spring (AMJ) ( Table 2 ). The simultaneous time series of JFM and AMJ temperatures and CPUE (kg/day) of European lobster are shown in Figure 5 for the northern and in Figure 6 for the central Adriatic. A high correspondence between the two variables is evident: in both areas, the CPUE values significantly increases in 2013/2014, at the same time as the bottom temperatures. Most individual peaks in CPUE also correspond with individual peaks of the temperature time series (2014 and 2016 for the northern Adriatic; 2011, 2014 and 2016 for the central Adriatic).

The observed temperature changes might be the result of ongoing climate changes or just a local or regional phenomenon. Analyses of satellite-derived sea surface temperature documented the largest trend in June, with a rate of 4.3°C over 100 years over the whole Mediterranean and some 30% higher rates over the northern Adriatic (Pastor et al., 2018). This has been attributed quasi-equally to multi-decadal oscillations (like Atlantic Multidecadal Oscillation, Knight et al., 2006) and to real warming trends (Iona et al., 2018). Vilibić et al. (2019) found a substantial increase in temperature in the northern Adriatic between 1979 and 2017, particularly high on the surface and during the summer season. Modelled bottom temperatures in the present study show a substantial increase in winter (JFM) and spring (AMJ) values over the whole Adriatic between the first (2008–2012) and the last (2013–2017) five years of the simulation. These changes are clearly visible on the corresponding difference plots for JFM and AMJ ( Figure 7 ). The largest increase in both seasons was found in the shallow northern Adriatic, particularly in zone A (also the zone where most lobster is caught), thereby confirming the findings of Pastor et al. (2018). The increase is less pronounced in the deeper central and southern parts of the Adriatic, as the seasonal thermocline is far shallower than the ocean depth there (Buljan and Zore-Armanda, 1976; Lipizer et al., 2014). The JFM temperature increase in zone A was between 1 and 2°C, and between 0 and 1°C over most of the central Adriatic. The AMJ temperature increases were even higher in both areas: up to 2.5°C in zone A of the northern Adriatic, and up to 1.5°C over the central Adriatic zones. Despite the observed changes, JFM temperatures remained lower in the northern than in the central Adriatic. On the other hand, the northern AMJ temperatures remained higher than central ones, with the temperature difference between the two areas increasing over time, due to the more efficient heat transfer towards the bottom in shallow waters (Artegiani et al., 1997).

The total annual European landings of H. gammarus over the last 10 years has been approximately 5000 tonnes (FAO, 2021). However, northern European countries have reported considerably higher landings than those in the Mediterranean basin, and a pattern of low H. gammarus abundance is evident throughout the Mediterranean when compared with Atlantic stocks. However, statistics of H. gammarus landings obtained by small-scale fisheries in general should be considered with caution across the Mediterranean, since they are difficult to evaluate and are often underestimated (Lloret et al., 2018; Pere et al., 2019). The Fisheries Department of the Croatian Ministry of Agriculture have been collecting data in a uniform manner since 2008 due to the pre-requisites of entering the EU. Thus, the ten years of data selected for the present study can be considered to have a higher degree of reliability. However, the reliability of catch statistics may be insufficient, since fishery-independent data were not available for validation of the landings and CPUE data with reference to population abundance, as suggested by Salas et al. (2007). As an example of how drastically different official statistics can be from actual ones, shows the study of Kleiven et al. (2012) where total estimated catch of European lobster was 14 times higher than officially reported. Also, CPUE data could be biased due to improvement of fishing technology, a phenomenon known as technological creep (Kleiven et al., 2022). For sure, on board surveys representing fishery-independent data together with scientific surveys will increase the reliability of catch data issue. Fluctuations in lobster abundance may occur as a consequence of the combination of environmental and fishery-related processes (Lizárraga-Cubedo et al., 2015). Since an increase in landings in the northern Adriatic was also associated with an increase in CPUE, this suggests higher abundance and/or increased catchability (Bueno-Pardo et al., 2020).

In this study, a positive and strong correlation was observed between winter (JFM) and spring (AMJ) sea bottom temperatures and CPUE. Similarly, McCleese and Wildner (1958) detected a strong correlation between long-term catch rates of the American lobster and sea surface temperature (SST) at the largest spatial scales, with lags of 0–3 years. More recently, Zhao et al. (2019) also reported that a temperature rise in Gulf of Maine led to increased catchability of American lobster over many years, with an expanded juvenile habitat in the north (Steneck and Wahle, 2013; Tanaka and Chen, 2016). These environmental changes have also been accompanied by the decline of large predators (Le Bris et al., 2018). On the contrary, warming waters have been associated with declined landings related to decreased juvenile habitat availability (Tanaka and Chen, 2015; Wahle et al., 2015) and increased prevalence of epizootic shell disease in the southern Gulf of Maine (Glenn and Pugh, 2006). Although there is no reference about the temperature range limits for European lobster in the published literature, Caputi et al. (2013) and Green et al. (2014) warn that H. gammarus may be sensitive to climate change, as rising temperature is the most important factor driving shifts in its distribution range, with increased abundance found at higher latitudes and decreased abundance at lower latitudes (Le Bris et al., 2018 and references therein). Sea temperature is one of the most important environmental factors affecting the fluctuations of lobster abundance, but how it contributes to total variation in the catch rate, and whether are there any temporal and spatial differences, remains unclear.

In homarid lobsters, sea temperature influences behaviour, which in turn affects their availability to fisheries (Lizárraga-Cubedo et al., 2015). Catchability is related to lobster movement and affected by numerous factors, including feeding behaviour and moulting status (Wahle et al., 2013), both of which are closely related to temperature (Green et al., 2014). Higher water temperatures affect lobster movements and catchability by encouraging increased movement, thereby increasing catchability (Smith et al., 1999; Moland et al., 2011). Rising seawater temperatures at the beginning of fishing season causes higher mobility and lobster feeding activity (Bennett and Lovewell, 1977; Lizárraga-Cubedo et al., 2015). This is followed by reduced mobility attributable to moulting activity (Miller, 1990; Sheehy et al., 1999), and then increased mobility due to higher activity observed during the post-moult period in summer. Moland et al. (2011) indicated limited movement of adult European lobsters, while American lobster exhibits higher movement and migration patterns (Estrella and Morrissey, 1997). Also, Moland et al. (2011) reported that seasonal variation in H. gammarus activity was correlated to water temperature, where lobster activity declined during the winter, with a minimum during February and March, and resumed again in April. This is in line with the winter (JFM) and spring (AMJ) temperatures chosen in this study, preceding and during the fishing season. Temperature is also positively correlated with growth rate among crustaceans, due to the within species thermal tolerance (Hartnoll, 2001). Higher temperature can positively stimulate growth rate by decreasing the time of the intermoult period or by increasing the moult increment (Green et al., 2014). This of course reduces the period in which lobsters are most susceptible to predatory mortality, while also increasing their fitness and the potential to affect growth rate and size at maturity.

Most crustaceans have synchronised spawning and time their reproduction based primarily on temperature (Lawrence and Soame, 2004). Harmonising hatching with food abundance increases the larval survival rate (Cushing, 1972). In this study, the reported positive trends in seawater temperature could possibly be responsible for more successful spawning and increased recruitment in the northern Adriatic, particularly due to the fact that positive trends were observed even before 2008 (Vilibić et al., 2019). The egg-bearing females of H. gammarus spawn in late summer (Tully et al., 2001), and egg hatching occurs in late spring to early summer (Phillips, 2013). The duration of egg development is largely influenced by temperature, with increasing temperature shortening the egg incubation process (Green et al., 2014). In the context of climate change in the North Sea, regimes with elevated temperatures (mild winters) resulted in a strong seasonal forward shift of larval hatching, while experiments showed that larval duration decreased and survival increased significantly at higher temperatures (Schmalenbach and Franke, 2010). The first few weeks post-hatching are characterised by a pelagic phase and the duration of this phase is temperature-dependent and reported to last for 14–35 days (Jørstad et al., 2001; Browne et al., 2009). Although specific observations of benthic post-hatch larvae of European lobsters in the wild are still lacking (Linnane et al., 2001), it is assumed they settle and remain cryptic in shelter-providing rocky substrata and emerge from their shelters only once they reach capapace lengths (CL) between 25 and 40 mm (Linnane et al., 2000a; Linnane et al., 2000b; Ball et al., 2001). Therefore, high temperatures in shallow coastal waters during summer in the Mediterranean could play a major role in juvenile survival and the recruitment success (Pere et al., 2019). Moreover, a recent study in the northern Adriatic reported a high number of pre-adult lobster (Pavičić et al., 2021).

Regulatory measures of lobster management in Croatia include a minimum landing size, closed season and prohibition of catching berried females (OG, 2016). Currently, the fishing season is open from 5 May until 1 September. Given the possible earlier increased movement due to elevated temperature, a trend of illegal lobster catching was observed in Croatia before the open season, mainly in March and April (pers. comm). Modifications of the existing time frame of the fishing ban should be considered in the future. Also in recent years, V-notching schemes for berried females (Tully, 2001) and protected areas were implemented (Moland et al., 2021) in Europe, and the implementation of these measures should also be considered in Croatia. Recently established MPA network for European lobster in Norway provides good example of management measures (Knutsen et al 2022). After establishment of MPAs protection effects started to manifest, including effects on density, growth, demography, behaviour, and phenotypic diversity. Climate change will surely provoke changes in all fishery sectors, professional and recreational. Artisanal and industrial professional fishers may adapt to these changes mainly through the expansion of fishing grounds following the distribution of target species, which will consequently increase operation costs. The effect on the fishing community is highest when the socio-ecological system is already under pressure, such as with overfishing in the Mediterranean (Miller et al., 2010; Colloca et al., 2011; Pranovi et al., 2013) related to decreased demographic structure, limitations of geographic distribution and diversity loss (Rijnsdorp et al., 2009; Perry et al., 2010; Planque et al., 2010). It is well known that differences in lobster size at maturity, fecundity and population size structure between different areas may occur as a response to the local environmental conditions and fishing strategies (Lizárraga-Cubedo et al., 2015). Since lobster fisheries are most often regulated only through the minimum landing size (MLS) (Pere et al., 2019), fishery managers need to be sure that MLS regulation currently in force corresponds to real size at maturity, throughout the distribution range. Thus, the EU Directive recognises and prescribes different minimum landing size for H. gammarus for northern European countries and the Mediterranean (European Union, 2006). It is possible that this should also be considered at the subregional level, since different genetic populations have been documented in the distribution range (Triantafyllidis et al., 2005; Ellis et al., 2017), while Adriatic populations are panmictic (Pavičić et al., 2020b). It can be expected that if these warming trends continue, the differences in the main biological points between Atlantic and Mediterranean stocks will become even more pronounced. In this study, we confirm that the northern Adriatic as a particularly vulnerable area to climate change. However, it is configured as a cul-de-sac (Ben Rais Lasram et al., 2010; Pranovi et al., 2016), preventing further northward migration of temperate and boreal affinity species like H. gammarus. Considering its shallowness, it is questionable how these species will behave in the future. Further on, recently, the American lobster (Homarus americanus) was reported in Adriatic Sea (Pavičić et al., 2020a) and Aegean Sea (Kampouris et al., 2021). If American lobster establish population in the Mediterranean Sea it can negatively affect European lobster populations, since these two species compete for the same habitat and hybridization is possible (Jørstad et al., 2007). For certain, the observed changes require a deeper and more complex analysis of the lobster stocks and updating of current management measures in the Mediterranean region.