Bayesian Network Analysis reveals resilience of the jellyfish Aurelia aurita to an Irish Sea regime shift

Robust time-series of direct observations of jellyfish abundance are not available for many ecosystems, leaving it difficult to determine changes in jellyfish abundance, the possible causes (e.g. climate change) or the consequences (e.g. trophic cascades). We sought an indirect ecological route to reconstruct jellyfish abundance in the Irish Sea: since zooplankton are jellyfish prey, historic variability in zooplankton communities may provide proxies for jellyfish abundance. We determined the Bayesian ecological network of jellyfish–zooplankton dependencies using jellyfish- and zooplankton-abundance data obtained using nets during a 2-week cruise to the Irish Sea in 2008. This network revealed that Aurelia aurita abundance was dependent on zooplankton groups Warm Temperate and Temperate Oceanic as defined by previous zooplankton ecology work. We then determined historic zooplankton networks across the Irish Sea from abundance data from Continuous Plankton Recorder surveys conducted between 1970 and 2000. Transposing the 2008 spatial dependencies onto the historic networks revealed that Aurelia abundance was more strongly dependent over time on sea surface temperature than on the zooplankton community. The generalist predatory abilities of Aurelia may have insulated this jellyfish over the 1985 regime shift when zooplankton composition in the Irish Sea changed abruptly, and also help explain its globally widespread distribution.


Scientific Reports
| (2021) 11:3707 | https://doi.org/10.1038/s41598-021-82825-w www.nature.com/scientificreports/ and between B and C, it is likely that there is also an indirect dependency between A and C: the BNIA would only report the direct dependencies (as two edges), not the indirect one. BNIAs can detect direct positive correlations, which for spatial data indicate spatial co-occurrence, that are quantif by a positive Influence Score (IS). Negative correlations have a IS < 0 and represent a negative correlation or spatial segregation such that when the abundance of one species is high, the abundance of the other species will be low. Associations between the two species will be either positive or negative, dependent on the species abundances are coded within BNs as IS = 0. Thus, BNs presents information on both magnitude and direction of interactions between variables and select only those interactions that are direct. These generated BNs can be used to infer how the network is likely to change under different scenarios.
In this study, we first determine jellyfish-zooplankton interacts by using BNIAs to infer contemporary (2008) spatial jellyfish-zooplankton networks in Irish Sea field data. Secondly, we built zooplankton networks (using weekly data grouped into three decades 1971-2000) using historic Irish Sea data from the Continuous Plankton Recorder, and then thirdly we inferred historic jellyfish abundances  by applying the contemporary jellyfish-zooplankton interactions to the historic zooplankton networks. It was our hope that these analyses would enable us to infer how jellyfish abundances may have changed historically in the Irish Sea over a period when a regime shift occurred and the herring fishery collapsed 27, 45 .

Methods and materials
Contemporary (2008) Irish Sea spatial jellyfish and zooplankton data. The data for the contemporary jellyfish-zooplankton network were collected over a 2-week cruise in the Irish Sea on the RV Prince Madog in May 2008 (Fig. 2).
Net sampling. To determine the mesozooplankton community composition, and abundance of jellyfish (pelagic medusae of the Cnidarian Classes Hydrozoa and Scyphozoa-these are technically zooplankton, but we call them jellyfish here for the sake of clarity), net sampling was undertaken using a horizontally-towed midwater MIK net, a vertically hauled 1 m diameter ring net, and an otter trawl (jellyfish only). A total of 144 stations were sampled (Fig. 2). Zooplankton were identified (we encountered a total of 111 species) and counts per sample were standardised to counts per unit volume with reference to flowmeter data that enabled volumes of water filtered during each haul to be determined. Two MIK net samples contained Aurelia aurita and Cyanea capillata, as did 65 ring net samples (see supplementary data). Of the 23 otter trawl deployments, (see Visual, acoustic and aerial sampling for jellyfish. In addition to net sampling, efforts were made to obtain abundance data for jellyfish by visual counts from the foredeck cf. 4 , scientific echosounding (18,38 and 120 kHz; Fig. 3 cf. 46 ) and aerial survey 32,47 . Visual observations of jellyfish were made from the bow during daytime only. To derive categorical estimates of abundance for each station the visual observations that occurred 750 m before the station and 750 m after the beginning of each station were extracted, then the counts per distance were placed into Low, Medium and High categories.
Visual observation data were available in daylight only, and were heavily influenced by sea state and glare, so were not used in network analysis. Acoustic data were also ultimately discarded because the combination of transducer depth (transducers were in the vessel hull) and acoustic near-field effects 48 left the upper c. 10 m   49 . Further details of CPR methods are given in Refs. 23,50,51 . Each CPR sample records the number of zooplankton present in approximately 3 m 3 of water. 537 species of zooplankton are routinely identified by the analysts, and counts are standardised to take account of differing volumes of water filtered per sample 52 . For this study we considered abundant 24 morpho-groups routinely detected by the CPR in the Irish Sea. They consisted of: 9 species, 7 larvae groups, Isopods, Chaetognatha, Copepod nauplli, Calanus i-iv, Pseudocalanus adults, Para-pseudocalanus spp. and PCI (Chlorophyl) as a proxy for phytoplankton concentration (see supplementary data). For simplicity these groups (apart from PCI) will heron be referred to as zooplankton. We collapsed them into the five biogeographic groups identified previously 28,53 : Warm Temperate, Warm Temperate Oceanic, Temperate Oceanic, Shelf Sea and Sub Arctic (Fig. 4).
Data considerations for Bayesian Network Inference Algorithms. The Bayesian Network Inference Algorithm (BNIA) ingests discrete data to ensure that noise is masked and only the relative abundances of each taxon are important 34,35,39,40,54 . Calculation of the best Bayesian network given a data set is computationally intractable, so BNIA uses a search technique which explores different networks and returns the highest-scoring network encountered as the solution 55 . Zooplankton with no strong biogeographical affinities were removed from the analysis, namely Copepod nauplii, Calanus i-iv, Isopoda, Larvacea, Lamellibranchia, Echinoderm and Polychaete larvae, Chaetognath and Limacina retroversa. We combined the remaining 14 zooplankton species into five groups as described below (Table 1) 54 . In order to avoid Type I errors, often associated with high-levels of zero count data, we performed contingency-test filtering. This filtering disallowed an edge between two variables whose joint distribution showed no evidence of deviation from the distribution expected from their combined marginal distributions (chi-squared tests, p > 0.25) 35,39,40,54 i.e. those which were correlated due to high levels of mutual zeros. The data were grouped as follows: www.nature.com/scientificreports/ Step 1: Zooplankton Grouping. The zooplankton species were grouped by their biogeographical association (cf. 28,29,45 ; see Table 4) to enable us to determine the effect (if any) of the North Atlantic regime shift on the networks.
Step 2: Weekly Averaging. Weekly averages were taken of the CPR count data to ensure consistent networks.
The number of CPR samples varied between zero and fifteen per week.
Step 3: Temporal Grouping. Previous research indicates that the regime shift occurred in the mid-1980s, with the exact year of varying between species examined and analysis type 28 . We split the data into three temporal groups, 1971-1980, 1981-1990 and 1991-2000 in order to ensure that the entire potential regime shift window was contained within one temporal group (Fig. 4).
Incorporation of environmental data in BNIA. We included sea-surface temperature (SST) and an index of the North Atlantic Oscillation (NAOI) as environmental variables in the zooplankton networks. SST data were obtained as monthly medians from the UK Meteorological Office 56 , and weekly values were calculated by interpolation from those. The NAOI was defined as the difference between the normalised sea level pressure over Gibraltar and the normalised sea level pressure over the southwest UK, and monthly indexes were obtained from The Climate Research Group, University of East Anglia (http://www.cru.uea.ac.uk/ftpda ta/nao.dat). The data we used were as calculated by Ref. 57 .
Bayesian network inference. Bayesian network inference (BNI) was performed in Banjo 41 . Data preparation for Banjo (grouping and discretisation) and statistical analyses were performed in R v 3.3.2 cf. 58 . Further analysis of Banjo outputs, when required, used the functional language Haskell 59 . The scripts are available on Github: https ://githu b.com/egmit chell /boots trap. R code is available a thttps ://githu b.com/egmit chell /jelly fish.
To determine each Bayesian network, the four steps below were followed as per Ref. 35 : 1. Discretisation. We split the data into three intervals; zero counts, low counts and high counts. Low counts consisted of counts below the median for the species group and high counts were counts over the median. Medians were used rather than means because for some groups the high counts were very high and use of the mean would have resulted in a very small number of samples grouped in the highest interval. Employing a large number of bins preserves resolution in the dataset, while fewer bins provide more statistical power, and greater noise masking: three bins has been found to be optimal for ecological data 35,39,54 . Zeros were treated as a separate entity because the absence of individuals (i.e. zero presence) is ecologically very different to the presence of just one: this difference is in contrast to zero versus one in measures of gene expression, which varies on a continuous scale between zero and one 54 ( Table 2). 2. Contingency-test filtering. The CPR dataset (6010 trawls obtained from between 1 and 4 times per week, with 23 zooplankton and PCI (Chlorophyll Concentration) contained a large number of zeros (56%), which could potentially have created problems with false positives (Type I errors) due to autocorrelation between these zeros. We combated this problem in three ways. 1) By grouping zooplankton into biogeographic groups, instead of just using the individual species data. 2) By taking weekly averages and having zero counts as a separate bin we reduced the number of zeros to 39%: the BNIA we used requires that the data are evenly distributed across all the bins, so 39% is adequate for a 3 bin system. 3) We used Chi-squared filtering to remove edge pairs which showed no evidence of deviation from the distribution expected from their combined marginal distributions (chi-squared tests, p > 0.25 cf. 35 ). 3. BNIA. Banjo was run on each data set to search 10 million possible networks using a greedy search, whereby at each search point the optimal edge addition to the network was used. The number of possible parents for each node (the edges that feed into a node) was limited to three to help eliminate artefacts 54 . For each edge the influence score (IS) was calculated. 4. Model Averaging. If a network suggested by Banjo represented the underlying ecological network, we would expect the same network to persist even if a few sample points (under 1%) were removed from the input data. Note that this bootstrapping removes the weekly samples, rather than individual species or groups. We found, however, that there was variation within the suggested networks when a small number of sample points were removed, so we applied a model averaging to the networks to overcome this stochasticity. For each time period, we took 100 random samples, consisting of 90% of the total available data points. For each edge that occurred, the probability of occurrence was calculated as the number of times the edge appeared over the total number of 100-random-sample bootstraps conducted. This probability distribution of node frequencies was bimodal for each of the network (2008, 1970s, 1980s, 1990s). These bimodal distributions www.nature.com/scientificreports/ suggest that there are two Gaussian distributions of edges, the rare/low occurrence ones, and the highly probable edges. The mean IS for each edge over all bootstraps was also calculated.
For the BNI analyses, jellyfish count data from net samples were split into two groups: Aurelia aurita, which was the most abundant species (97.8% by number), and 'other taxa' which consisted of Cyanea spp. (1.9%) and Cosmetira pilosella (0.3%). Rhizostoma were also caught, but due to very low numbers they were not included in the analyses. The zooplankton were grouped by biogeographic affinity cf. 28 (as per Table 1) to ensure consistent networks. This biogeographic grouping was required because not all zooplankton species were present in all of the 144 samples, and zero-inflated data are not valid with the Bayesian priors used in this BNIA.
Inference of historical jellyfish abundance. One of the most powerful aspects of BNs is the ability to make inferences of how the probability of one node (taxa or physical variable) being in each state (zero, low or high for a taxon) is likely to change given that another node is in a given state (zero, low or high for a taxon) 38 . This inference is made by calculating the probability of node A being in a given state given node B is in a given state. For example, the probability of a given species being in a zero, low or high abundance state can be inferred for SST being in either a low or high state. This inference can be used to calculate the likely consequence for species abundance of changing variables such as SST.
All nodes that have dependencies between A and B are included in the calculation: The N are the total number of nodes in the network and n and m are the indices for the chain of nodes of length N connecting the first and last nodes. S are the number of discrete states for each node, which are indexed s. In order to infer the likely change of one taxon's (A) abundance on another (B), for example if A changes from a high abundance state to a zero abundance state, the probabilities of B existing in all states are calculated given a high abundance state for A, and the probabilities for B existing in all states are calculated given a zero abundance state for A. The inferred change in one taxon's abundance is then the difference between these probabilities.
The abundance probabilities for the biogeographical zooplankton groups (Warm and Cold Temperate) were calculated (e.g. P(High abundance|Low Temperature)) given the SST state for each of the three decades. The www.nature.com/scientificreports/ mean temperature in the Low SST state was approximately 4 °C lower than the High SST state. An increase in the probability of a species/group given the probability of another species/group having a zero count represents the increase or decrease in abundance in that area, for that temperature change.
In order to calculate historic jellyfish abundance probabilities, the dependencies of jellyfish on zooplankton were assumed to remain static over the three decades. Note the implicit assumption here that Aurelia continue to feed on the same biogeographic zooplankton groups, although the relative rate of feeding may change. For each decade the abundance probabilities of the biogeographic zooplankton groups that the jellyfish nodes depended on were calculated. The jellyfish node probability abundances (zero, low, high) given either Low SST or High SST were calculated for each decade given the dependencies found in the spatial network inferred from the 2008 field data.
Using Naïve Bayes Classifiers to infer the relative importance of different zooplankton groups to herring. In order to determine the relationships between zooplankton groups and herring, we used naïve Bayes Classifiers to determine the relative importance of each zooplankton biogeographic group and environmental group to herring. The herring data were recorded as annual catch data for the period 1970-2000 27 , which meant there were only 10 samples for each decade: this was insufficient for direct incorporation directly into the CPR BNs, but naïve Bayes Classifiers can determine and rank the importance of each group relative to each other 60 . Herring have a multi-year life-cycle, and there are lags between spawning and recruitment and capture. Impact by jellyfish on herring is possible as direct predation on eggs and larvae, and by competition for zooplankton prey 1 . To calculate a calendar years' worth of zooplankton data we took yearly averages of zooplankton starting in September (herring spawning time) 61 , discretized into low and high categories based on median values for each decade, and compared them to the zooplankton data the following year (Table 2). To rank the groups, the herring catch data were included in the zooplankton abundance network, linked via an edge to a zooplankton (and environmental) group, and the network score calculated. This scoring was repeated for all groups across the three decades. For each decade the rank of the biogeographical zooplankton group was given by the change in network score between the CPR networks without herring and the networks that included a forced link between the herring and the group. Note that these analyses make no assumptions or suggestions to whether the herring did depend on any zooplankton, instead they rank the relative importance of each group, so provide a different approach to that used with the jellyfish data.

Results
Irish Sea contemporary spatial network from cruise data. The Bayesian Network for the field data collected in 2008 connected all nodes (Fig. 5). Aurelia was dependent on Cold Temperate zooplankton, Oceanic and Warm Temperate zooplankton, and on the Cyanae spp. and Cosmetira pilosella jellyfish group. Note that because these two zooplankton groups are connected to other zooplankton groups, the implication is that Aurelia will be effected indirectly by all these other zooplankton groups too 38 . The dependency of Aurelia with Warm Temperate Zooplankton was negative: high Warm Temperate Zooplankton corresponded to low Aurelia abundances whereas dependence on the Cold Temperate Oceanic Zooplankton group was non-monotonic, i.e. had different effects dependent on zooplankton abundances. Three of the edges between the zooplankton, Chlorophyll Concentration (indicating phytoplankton) and jellyfish groups showed positive dependences, three negative and two non-monotonic (Fig. 5).
Historic zooplankton networks from CPR data. The Bayesian Networks for CPR-sampled zooplankton across the three decades have different dependencies (Fig. 6). There are two sub-networks for the 1970s and 1980s networks, but these join into a single network for the 1990s. Three dependencies are present throughout the three networks: Shelf Sea Zooplankton with Cold Temperate zooplankton; Cold Temperate Zooplankton with SST, and Warm Temperate Oceanic zooplankton with Sub Arctic Zooplankton. The merge of the two subnetworks from the 1970s and 1980s in the 1990s occurs via the inclusion of the Warm Temperature Zooplankton in the 1980s network, and with the connection of Warm Temperate Oceanic Zooplankton to Shelf Sea zooplankton in the 1990s.
In order to investigate whether the changes in network structure were statistically significant, tests were conducted on the occurrence of edges for the three different time periods. There was no significant difference between the 1980s and 1990s (p = 0.3076), but the 1970s were significantly different to the 1980s (p = 0.0211) and 1990s (p = 0.0156). The network for 1981-1990 represents the time period over which zooplankton community composition has been recorded as changing elsewhere in the North Atlantic 49,53,62 . The 1981-1990 network is the only network where the Chlorophyll Concentration and NAOI variables are connected to other variables: this connection is via the Warm Temperate Oceanic Zooplankton group.
Changes to the ecological networks over time were also evident as changes in the strengths of the dependencies, as measured by the Interaction Strength (IS). The Shelf Sea Zooplankton and Cold Temperate zooplankton IS did not change significantly (as assessed by Chi-squared test) between 1970-1980 and 1981-1990, however it did change significantly (p < 0.0001) between the last two decades. Warm Temperature Oceanic Zooplankton and Sub Arctic Zooplankton changed significantly (p < 0.0001) between all three periods, with mean ISs of 0.44, 0.17 and 0.53 for the 1970s, 1980s and 1990s respectively. The IS for Cold Temperate Zooplankton and SST increased significantly over the three time periods, increasing sequentially from 0.24 to 0.47 and then to 0.64 (p < 0.0001).

Ranking of relative importance of zooplankton groups to herring. The naïve Bayes classifiers
found that the ranking of zooplankton groups with respect to influence on herring changed across the regime shift period (    were 4th most influential on herring, but did not change the network score (ΔS = 0.0000), suggesting a loss of influence, while the NAOI remained ranked second and Warm Temperate Oceanic Zooplankton were ranked first (ΔS = − 4.0144). The NAOI influence score changed from positive pre regime-shift to negative post regimeshift. This change coincides with a shift in the NAO from a low phase to a high phase between the 1960s and 1990s 63 .

Inference and probabilities. The Warm Temperate and Cold Temperate Oceanic Zooplankton groups
have higher abundances at high SST ( Table 4). The probability of being in a high abundance state given either high or low SST increases through time for the Warm Temperate Zooplankton group. In contrast, the abundance of the Cold Temperate Oceanic Zooplankton group remained steady throughout the three time periods given low SST, but decrease given high SST. Across the three decades, low SST corresponded to high inferred jellyfish probability, and high SST to low inferred jellyfish probability (Table 4). For both high and low inferred Aurelia abundance the probability of being in a high abundance state remained relatively constant, in contrast to the abundances of the zooplankton groups which changed significantly over the regime shift (Fig. 7).

Discussion
Jellyfish-abundance time-series globally are rare, with less than 40 datasets over 10 years in length as of 2013 2 . Most of these datasets are from the northern hemisphere (87%), and in particular the Atlantic Ocean (17%) and the Mediterranean Sea (17%) 2 , so our understanding of how jellyfish abundances are changing globally through time is extremely limited. The debate on changing frequency of jellyfish bloom occurrences is however lively and ongoing, with a perception that the frequency of blooms is increasing e.g. Ref. 64 due to anthropogenic influences such a overfishing and climate change 6,65 . While increases in jellyfish abundance have been found in some areas 2 , including in the Irish Sea 1, the misreporting of increased jellyfish blooms has led to over-stated generalisations on the occurrence of bloom events 66,67 . Global analyses of the available jellyfish time-series in 2013 showed a negligible overall increase with time, with a weak increasing linear trend that was possibly just an up-phase on a longer-period global oscillation 2 : more data and analyses are needed in order to tease apart the potentially complex interplay of jellyfish with other marine life and environmental variability.
In order to extend time-series of jellyfish abundance to a broader temporal and geographical scale, we propose an eco-palaeontological method for inferring historic jellyfish abundance from historic zooplankton data. We demonstrated this approach here by first inferring ecological zooplankton networks for three decades using historic CPR data from the Irish Sea from 1971 to 2000. Our BNI showed marked changes in zooplankton networks across these three decades, but three dependencies remained throughout: Cold Temperate Oceanic Table 3. The naïve Bayes classifiers for the time-lagged herring catch data with each group. The more negative the change in score, the better the network fit to the data, which reflects an underlying dependency between the herring and the group.  (Fig. 6). In the 1980s, environmental factors became significantly more important-it was the only decade for which significant dependencies between a zooplankton group (Warm Temperate Oceanic Zooplankton) and Chlorophyll Concentration and an index of the NAO were found, possibly due to the regime shift perturbing ecosystem function. The regime shift resulted in fundamental changes to the historic zooplankton networks. Prior to 1990, the Warm Temperate Oceanic and Sub Arctic plankton did not depend on the other zooplankton groups, nor SST. As such, prior to 1990, there are two distinct ecological networks of zooplankton so that changes abundances within one network would not have affected the other network. Then, after the regime shift of the 1980s, in the 1990s Shelf Sea Zooplankton formed a dependency with Warm Temperate Oceanic Zooplankton, after the regime shift of the 1980s, so that changing abundances of any zooplankton group began to have knock-on effects across all the other groups. This increased connectivity likely increased ecosystem stability 68,69 by potentially providing a broader spectrum of prey to predators of zooplankton. We inferred the likely jellyfish-zooplankton networks for each of the three study decades (1970s, 80s and 90s) by adding the jellyfish-zooplankton interactions determined from the 2008 survey data networks to the historic CPR networks. Finally, we inferred Aurelia abundances for each of the three decades from the historic networks (Fig. 7). We found that inferred Aurelia abundances were not markedly impacted by changing zooplankton networks: T the probability of Aurelia existing in a high-abundance state being significantly higher for colder sea surface temperatures. This higher abundance state occurring despite significant changes over time in the abundance probabilities and Bayesian network structures for multiple zooplankton groups (Figs. 5, 7). It thus appears as though Aurelia aurita has been resilient to change in the composition of the zooplankton community, which is its prey. Our finding here that SST in the Irish Sea has a stronger influence on abundance of A. aurita than does the composition of the prey fields is consistent with a previous study that showed, using CPR data, that temperature was a more significant driver of the abundance of shelf jellyfish in the northern Benguela upwelling ecosystem than was food 16 .
From our results alone it is not possible to determine the underlying basis of the apparent resilience of Aurelia aurita to the mid 1980's regime shift. One possibility is whilst the species composition of the Irish Sea zooplankton community changed, the size spectrum of that community and the nutritional content did not change, so Aurelia′s preyfield remained functionally similar. Although the regime shift in the 1980s saw changes in species composition and ecological network structure, ecosystem stability and resilience from the perspective of Aurelia was maintained via various processes including functional and/or trophic redundancy 70,71 , competition trade-offs 72 , stabilizing feedback loops 36 and/or functional complexity 73 . Despite this putative functional similarity in the zooplankton cimmunity from the perspective of Aurelia, there were though major differnces  (Fig. 7). The increased connectivity within the 1990s zooplankton network may be due to the decreasing month-specificity of zooplankton group peaks, with a broadening of the abundance peaks for most zooplankton groups, resulting in an increased temporal overlap (Fig. 4). In the context of Cushing's match-mismatch hypothesis, this may give Aurelia a longer time window each year to surf the wave of zooplankton production. We suggest that the increased connectivity between different biogeographical zooplankton groups created redundancy in prey species, and so enabled prey-switching by Aurelia 74 . Aurelia had dependencies with the Warm Temperate and Cold Temperate Oceanic Zooplankton groups, in which cyphyonautes larvae dominate the Warm Temperate group and Calanus helgolandicus and barnacle larvae dominate the Cold Temperate Oceanic zooplankton group. Therefore, it is likely that these taxa are the key drivers for this network, with Aurelia switching prey post regime shift. In our data, particularly abundant members of the Warm Temperate Zooplankton group included cyphonaute larvae of bryozoans and numerous species of mesozooplanktic copepods. Aurelia spp. are known to be able to thrive on very small zooplankton 75 and indeed feeding experiments on Aurelia aurita from the Black Sea suggest that food items overlooked by researchers appraising prey visually may provide ten times the energy for Aurelia than that obtained from mesozooplankton 76 : it is possible that Aurelia medusa could be sustained by cyphonautes. Another possibility is that environmental conditions that favour strobilation (asexual reproduction) of Aurelia benthic phases also favour bryozoan reproduction. Our own work, for example, has shown that cold winter conditions lead to elevated rates of ephyrae production by strobilation 77 and that abundance of medusae in the North Sea is high in years characterised with low NAO index. A similar association is evident between NAO phase and abundance of cyphonaute larvae in Loch Hyne 78 . There was no clear association between Aurelia abundance probabilities (Table 4) and herring catch (Fig. 1), which suggests that Aurelia abundances and herring recruitment do not directly interact. This lack of interaction can be explained by our 2008 field-data network, which did not include any direct dependencies between Aurelia (or any other jellyfish) and the Sub Arctic Zooplankton group (Fig. 6) which is dominated by Calanus finmarchicus (Table 1), historically the primary food source of herring 79 . The ranked zooplankton groups (Table 3) show that the Cold Temperate Zooplankton (Acartia spp., Podon spp. and Evadne spp.) was the most important group prior to the regime shift, and afterwards changed to Cold Temperate Oceanic (dominated by cyphyonaute larvae) and Sub-Artic Zooplankton group (dominated by Calanus finmarchicus). It appears that Aurelia do not compete directly with herring for larger stages of Calanus finmarchicus, and while small herring may compete to some extent with Aurelia for the cyphyonaute larvae, these larvae historically formed only a relatively small proportion (5.4%) 80 of the herring diet so are unlikely to have a strong influence on herring success.
Our conclusions on jellyfish-herring interactions in the Irish Sea differ from those drawn by previous studies of the North Sea, which found that over the period 1958-2007 jellyfish were negatively impacting herring recruitment 4,7 . The biogeographical distribution of Calanus finmarchicus changed significantly in the North Sea Figure 7. The network used to infer historic Aurelia abundances across the three decades (70s: purple edges, 80s: blue edges and 90s: cyan edges). The red edges and nodes were inferred from the cruise data, and are assumed to be constant over the historic study period  www.nature.com/scientificreports/ during the mid-1980s 24 , suggesting that Irish Sea herring would have experienced a less substantial change in their prey field during this time. In any case our networks suggest minimal interaction between Calanus finmarchicus in the Irish Sea. It would be instructive to run a network analysis of North Sea zooplankton and jellyfish, but that is beyond the scoipe of the present study (Fig. 6). Plankton communities in the Irish Sea, like as in the North Sea, are primarily driven by hydroclimatic conditions 24 , with some areas now showing increased autumnal diatom production due to anthropogenic nutrient enrichment 81 . Future anthropogenic impact and climatic changes are likely to lead to further changes for Irish Sea plankton communities. These changes are likely to include further zooplankton biogeographic shifts, with cold-water species being replaced by the warm-water species as the cold-water species retreat north. It is not possible to predict from our analyses here how these changes will affect jellyfish abundances, and as jellyfish responses would depend on whether further changes to network structure occur. However, jellyfish abundances may not necessarily increase: from our inferred historic jellyfish abundances there are noticeably higher abundances with cold sea surface temperature, so any increased warming may have a negative impact on jellyfish abundances, despite their resilience over the 30 year study period here 5 .
Our BNI technique has provided a window into the past and revealed a new perspective on jellyfish ecosystem interactions in the Irish Sea. Using a combination of historic zooplankton data with contemporary jellyfish data we have been able to demonstrate that, despite a change in the ecological networks of zooplankton between the 1970s and 1990s, the abundance of Aurelia aurita most likely remained stable: we suggest that is a consequence of Aurelia′s generalist ability to take a wide range of prey. Our work suggests a remarkable resilience of Aurelia aurita to ecosystem reorganisation experienced in the Irish Sea over the regime shift of the 1980s. The generalist ability of Aurelia as exhibited in the Irish Sea may help explain the pan-global distribution of the Moon jellyfish species complex and its success in a diversity of locations 82,83 .