Half a century of monitoring macrobenthic animals on tidal flats in the Dutch Wadden Sea

Macrobenthic animals living in a tidalflat area in the westernmost part of the Wadden Sea were monitored for 50 yr (1970−2019) using consistent methods. About 100 papers were published on this project. We review a number of results and conclusions on observed changes and their possible underlying causal processes. The most significant changes in population sizes and growth rates of several species could be attributed to climate warming (by about 2°C), along with increasing trends in species richness and total late-winter zoobenthic biomass. In the initial years, eutrophication (doubling of nutrients and chlorophyll concentrations) resulted in a doubling of zoobenthic biomass. The subsequent de-eutrophication after the mid-1980s was reflected only in the biomass values observed in late summer. A long-term trend in food supply for birds was not observed. Disturbances from fisheries were intermittent and modest. In several bivalve species, magnitudes of production and biomass were determined primarily by recruitment variation, which was mainly caused by spring abundance of epibenthic predators (shore crabs and shrimps). Their abundance in creased with temperatures in the preceding winter. In contrast to this top-down regulation, bottom-up processes apparently played only a minor role in the determination of bivalve biomass. Rarely occurring extremely high bivalve numbers resulted in reduced rates of growth and production. We conclude that the uniquely long monitoring of the tidal-flat macrozoobenthos yielded data series which not only indicated several longterm trends, but also contributed to our insight in processes underlying the observed trends. Most of the observed trends were related to climate change and eutrophication followed by de-eutrophication. OPEN ACCESS


INTRODUCTION
For a long time, investigators generally regarded monitoring of ecosystems as a rather dull job of limited scientific value. However, monitoring is not routine and trivial, as it (1) characterizes ecosystem shifts by providing evidence for biological changes as a consequence of large-scale and long-term environmental changes (such as caused by climate ABSTRACT: Macrobenthic animals living in a tidalflat area in the westernmost part of the Wadden Sea were monitored for 50 yr (1970−2019) using consistent methods. About 100 papers were published on this project. We review a number of results and conclusions on observed changes and their possible underlying causal processes. The most significant changes in population sizes and growth rates of several species could be attributed to climate warming (by about 2°C), along with increasing trends in species richness and total late-winter zoobenthic biomass. In the initial years, eutrophication (doubling of nutrients and chlorophyll concentrations) resulted in a doubling of zoobenthic biomass. The subsequent de-eutrophication after the mid-1980s was reflected only in the biomass values observed in late summer. A long-term trend in food supply for birds was not observed. Disturbances from fisheries were intermittent and modest. In several bivalve species, magnitudes of production and biomass were determined primarily by recruitment variation, which was mainly caused by spring abundance of epibenthic predators (shore crabs and shrimps). Their abundance in creased with temperatures in the preceding winter. In contrast to this top-down regulation, bottom-up processes apparently played only a minor role in the determination of bivalve biomass. Rarely occurring extremely high bivalve numbers resulted in reduced rates of growth and production. We conclude that the uniquely long monitoring of the tidal-flat macrozoobenthos yielded data series which not only indicated several longterm trends, but also contributed to our insight in processes underlying the observed trends. Most of the observed trends were related to climate change and eutrophication followed by de-eutrophication. warming, eutrophication, or fisheries), (2) generates data on complex multiannual processes taking place in ecosystems, (3) can be used to test various hypotheses and models (e.g. on productivity and on stress), and (4) offers a means of distinguishing between natural and anthropologically induced changes in ecosystems. Thus, coupling monitoring results to proper analyses of mechanisms behind changes in the studied ecosystem will lead to a better understanding of ecosystem functioning and alterations therein, fully warranting monitoring efforts. For the evaluation of some ecological problems caused by ecosystem changes, monitoring is the only way to obtain meaningful results.
The aim of the present paper is to show, based on a 50 yr data set of the numbers and growth rates of some 50 species of benthic animals sampled twice per year at 15 fixed stations in a 50 km 2 tidal-flat area, that data obtained by monitoring can yield much more evidence than simply recording fluctuation patterns and trends. These data allow a study of the relationships among numbers, growth rates, production, and biomass of zoobenthos with available long-term data on various environmental factors such as (1) air and water temperatures providing evidence for effects of climate warming, (2) chlorophyll concentrations providing evidence for effects of (de)-eutrophication, and (3) fisheries providing evidence for effects of disturbance and subsequent recovery. Further, the data offer a basis for studying the effects of changes in sea level. Data on diversity, population dynamics, species interactions, and productivity allow us to study the underlying mechanisms of change. The obtained data are also relevant for developing and testing various hypotheses and ecosystem models. Over the course of this 50 yr monitoring effort, some 100 papers on the project have been published. Here we summarize and review the results of these papers and provide conclusions based on these results.

Study area
Since the early 1970s, sampling has been carried out twice annually (in late winter and in late summer) at 15 fixed stations on Balgzand, a 50 km 2 tidal-flat area in the northwest of the Netherlands, at about 53°N, 6°E. The position of the permanent stations (distributed throughout Balgzand), their intertidal level (from 0.5 m above to 0.9 m below mean tide level), and their sediment composition (from 0.1 to 30% fine particles of < 60 µm) were described by Beukema & Cadée (1997a). The main environmental gradient runs from high (sheltered and muddy) coastal areas via intermediate tidal flats (with moderate elevation and mixed sediments) to low (exposed and sandy) areas far off shore. The permanent sampling sites cover the entire gradient.

Sampling procedure
The bottom fauna was sampled by taking 50 cores of 0.019 m 2 (winter) or 0.009 m 2 (summer), each to a depth of about 30 cm at 20 m intervals along 12 transects of 1 km length. Moreover, 9 samples of 0.105 m 2 each were taken from 3 plots of 900 m 2 . The bottom samples were sieved in the field through 1 mm mesh screen. All samples were sorted in the laboratory while the animals were still alive.
Abundance was expressed both in density (n m −2 ) and in ash-free dry mass (AFDM, g m −2 ). Bivalves were aged by counting external year marks on their shells. This sorting by age class allowed calculations of rates of annual growth, survival, and production.
Data on surface water temperatures were obtained from a Royal Netherlands Institute for Sea Research (NIOZ) data base on daily records at a jetty in the Marsdiep tidal inlet, ca. 5 km from Balgzand. Data on chlorophyll a concentrations (mostly 1 or 2 samples per month) were available from Marsdiep data sets of NIOZ and of the Directorate-General for Public Works and Water Management (Rijkswaterstaat) .

Climate warming
Over the last 50 yr, water temperatures in the western Wadden Sea have risen at a rate of on average 0.4−0.5°C per decade (Beukema & Dekker 2020;Fig. 1). Severe winters (with mean water temperatures in January/February of < 2°C) were recorded 7 times between 1970 and 2000, but have not occurred in the 20 yr since. On the other hand, hot summers (with mean water temperatures of > 20°C in July/ August) were observed only after 2000 (3 times).
In contrast to hot summers with few observations of (minor) influences, serious effects of severe winters on the tidal-flat ecosystem of the Wadden Sea have been observed . Out of the ~50 common zoobenthic species inhabiting the tidal flats, some 10−12 species showed increased and high (up to 100%) mortality rates during the coldest winters (Beukema 1989, Beukema & Dekker 2011b. As ca. 10% of the species completely disappeared in the coldest winters, species numbers were seriously reduced after a cold winter (Beukema & Dekker 2011b). With rising temperatures, species numbers on Balgzand have shown significantly increasing long-term trends (Table 1), both when expressed as number per 14 m 2 (all sites together) or per 0.95 m 2 (sampling site). For a substantial part, this trend was due to winter-sensitive species such as Cerastoderma edule, Nephtys hombergii, Lanice conchilega, Crangon crangon, and Carcinus maenas. Therefore, we think that the significant correlation of species richness with temperatures was a causal relationship. The long-term increase in species richness was partly caused by the declining negative influence of low winter temperatures and partly by the arrival of new species, whereas no species disappeared.
The fauna at sites located relatively high in the intertidal zone was hardly affected by low winter temperatures. The winter-sensitive species tended to be found more frequently at relatively low levels in the intertidal, where total biomass was roughly halved after a severe winter (Beukema , 1985. Winter-sensitive species that did live at higher sites showed decreased survival there than lower in the intertidal zone, although the survival of sensitive species in the coldest winters was very low everywhere (Beukema 1985). At the highest sites, temperatures are more extreme than at lower sites, which are covered most days at high tide by relatively warm flood water (mostly from the North Sea).
In severe winters, cold-sensitive species partly survive in the deeper tidal streams and the nearby North Sea. As a consequence, repopulation of the tidal flats after a cold winter always occurred rapidly, and almost complete recovery to pre-winter distribution patterns and species richness took place within half a year (Beukema & Dekker 2011b), both by migrations from deeper water as well as from augmented recruitment after cold winters .
Total biomass values were also reduced after cold winters. The positive relationships between winter water temperatures and total zoobenthic biomass as well as bivalve biomass at the end of winter were highly significant (Fig. 2). The causation of this relationship was complex. On the one hand, and most importantly, survival of several species was better in mild winters (see above). Higher summer temperatures in recent years have promoted bivalve growth, which would also enhance biomass (Beukema et al. 2017a). On the other hand, winter weight losses in bivalves increased with increasing winter temperatures (Honkoop & Beukema 1997). As to be expected from the long-term temperature trends, total zoobenthic biomass at the end of winter showed positive long-term trends with year ( Fig. 3A, Table 2). Apart from the positive influence of rising temperatures, a number of recently introduced species also contributed to this trend (see Section 3.5).
More detailed investigations of the temperature effects on numbers and growth were done in several bivalve species. Cockles C. edule of all ages almost completely died out in winters with mean Year Water temperature (°C)  vival in adult cockles only, but thus far this has not had a strong effect (Beukema & Dekker 2020). Cold winters promoted recruitment in the subsequent summer. The difference was substantial: recruitment was about 10 times more successful after the coldest than after the mildest winters (Beu kema & Dekker 2020). In the same paper, we showed that summer temperatures did not affect cockle recruitment. The negative effect of higher winter temperatures on recruitment was probably not a direct effect of temperature, but rather an indirect effect via the positive influence of high winter temperatures on the abundance in subsequent spring months of shrimps C. crangon and shore crabs C. maenas, which are important predators of bivalve spat . Recruitment success was negatively related to densities of both shrimp as well as shore crabs in spring. With the rising winter and summer temperatures in the course of the 50 yr observation period, a consistent change in cockle dynamics was observed: recruitment significantly declined and winter survival increased. The net effect of the 2 opposing forces was an absence of any long-term trend in the abundance of adult cockles (Beu kema & Dekker 2020).
As with cockles, the annual growth rates of blue mussels Mytilus edulis and gaper clams Mya arenaria also increased with increasing temperatures in the growing season (Beukema et al. 2017a), and winter weight losses increased with winter temperatures (Honkoop & Beukema 1997). In M. edulis, recruitment declined with increasing winter temperatures (Beukema & Dekker 2014), whereas in M. arenaria, this negative relationship was non-significant (Table 1). In both species, negative relationships were observed between their recruitment success and the abundance of shrimp and/or shore crabs in spring . In the course of the 50 yr period of observation, the decline in annual recruitment was nearly significant in M. edulis, but non-significant in M. arenaria. A positive long-term trend in biomass was observed in M. arenaria (Beukema & Dekker 2019b), probably resulting mostly from increased survival by faster growth, which allows individuals to reach predator-safe sizes in a shorter time. In the Baltic tellin Limecola (Macoma) balthica, growth as well as survival showed relationships with winter and summer temperatures that were very different from those in the above 3 bivalve species. Growth rates were reduced at high summer temperatures . Survival rates were reduced at high temperatures in winter (Beukema et al. 2017c) as well as in summer (Beukema et al. 2009). As in other bivalves, recruitment was more successful after cold than after mild winters (Beuke ma et al. 2009(Beuke ma et al. , 2017c. The same causes  likely played a dominant role, al though Philippart et al. (2003) suggested that a reduced reproductive output and advancement of spawning were also involved. However, Honkoop et al. (1998) proved that the magnitude of reproductive output was non-decisive for recruitment success, and  found that L. balthica phenology did not differ from that of the predators of their post-larvae. Thus, the observed reduced re cruitment after higher winter temperatures could have resulted completely from increased predation pressure by epibenthic predators, namely shrimps and shore crabs . Spat-sized L. balthica migrate in winter from high to low tidal flats and to subtidal ar-  (Beukema & de Vlas 1989), where conditions for survival and growth are better (Beukema 1993a). In colder winters, a higher proportion of the animals migrated ). Long-term changes found in total numerical densities of L. balthica included negative trends in recruitment and adult abundance, but these changes became prominent only after the year 2000 and were probably related more strongly to an epidemic disease than to rising temperatures (Beukema et al. 2017c). Dairain et al. (2020) presented evidence that disseminated neoplasia played a key role in the dynamics of L. balthica and might explain its recent population decline. As a consequence of the strongly negative long-term trend in adult survival, the biomass values in this species showed a significantly negative relationship with year (Beukema & Dekker 2019b). Compared to the other important bivalve species, annual growth rates in L. balthica declined with increasing temperatures during the growing season (Beukema et al. 2009(Beukema et al. , 2017a. Moreover, higher winter temperatures enhanced winter weight losses (Beukema et al. 2009), as they did in the other important bivalve species (Honkoop & Beukema 1997).
Relationships of abundance and growth with temperatures were studied in 2 other tellinid species: Abra tenuis (Dekker & Beukema 1993 and Macomangulus (Tellina) tenuis . Both species were found to be winter-sensitive, showing almost complete population collapses in winters with January/February water temperatures of <1.5°C . In A. tenuis, recruitment was higher when temperatures in summer (their reproductive period) were higher (Dekker & Beukema 1993). Apparently, temperatures in the preceding winter did not affect their recruitment success. Both species are present at spat size only in summer, much later than the other bivalves that produce peak numbers of spat in spring. As opposed to spring densities, summer numbers of predators of bivalve spat were only weakly related to the conditions of the preceding winter. In both M. tenuis as well as in A. tenuis, higher annual growth rates were observed in warm than in cool summers (Dekker & Beukema 1993. In summary, the strongest and most consistent effects of warming climate appear to be long-term increases in species richness (Beukema & Dekker 2011b). Late-winter biomass increased as well (Fig. 3), but as this was mostly due to increasing amounts of large M. arenaria and Ensis leei and not to smaller animals, food availability for birds did not show a positive long-term trend (Beukema et al. 2010, Beukema & Dekker 2019b. Another significant effect of climate warming was the increase in the abundance and earlier spring appearance on tidal flats of epibenthic predators of bivalve spat, i.e. shore crabs (Beukema 1991b) and shrimps (Beukema 1992b), resulting in declining recruitment success with increasing winter temperatures in most bivalve species (Table 1). As recruitment success is an important determinant of subsequent biomass (Beukema et al. 2010, Beukema & Dekker 2019a) and production (Beukema & Dekker 2019a), a long-term decline in bivalve biomass and production might be expected from recruitment declines (Beukema 1992a). For instance, a succession of a number of years with failing recruitment (as in 1988, 1989, and 1990) indeed resulted in very low bivalve biomass values (Beukema & Cadée 1996). However, to date, a general long-term decline in total bivalve biomass has failed to occur (Table 2), probably as a consequence of declining mortality rates (Beukema & Dekker 2019b).
Another consistent effect of climate warming was the increased weight loss of bivalve individuals at higher winter temperatures (Table 1), caused by enhanced metabolic levels in the absence of sufficient food intake. As warmer winters are followed by on average warmer spring/summers periods, these losses were compensated in most species by stronger growth in such years. However, this was not the case in L. balthica, which appears to be the species that is most seriously suffering from climate warming (Beukema et al. 2009). Warmer winters resulted in lower values of condition index (weight at length), lower survival, and lower egg production (Beukema et al. 2001b(Beukema et al. , 2009. After mild winters, the adults were found to reduce their output of eggs to prevent their body condition from lowering (Beukema et al. 2001b). In contrast to this one species that is negatively affected by climate warming, there are several species that benefit from climate change through increased survival and growth, including bivalves, polychaetes, and crustaceans (see above). Philippart et al. (2007) described how nutrient loads to the western Wadden Sea sharply increased around 1980 and declined after 1987. Chlorophyll concentrations and primary production closely followed these changes (Beukema et al. 2002, Philippart et al. 2007. Chlorophyll concentrations in the nearby Marsdiep showed sharp increases in the late 1970s and a subsequent gradual decline starting in the mid-1980s (Fig. 4). In this section, we discuss the possible response of the benthic fauna to these changes in primary food supply.

Eutrophication
The Balgzand tidal flats are rich in bottom animals: in most periods, late-winter biomass amounts of around 30 g AFDM m −2 were present (Table 2). Latesummer biomass values were about 2 times higher (Fig. 3). Late winter is the season with the annual minimal biomass values (Beukema 1974). Molluscs (nearly all of their biomass being bivalves) were the most important group, with 65−70% of the total biomass; the share of polychaetes amounted to 30−35%. These percentages did not show any long-term trends. In this sense, the Balgzand benthic ecosystem appears stable. Results of surveys of the bottom fauna over more extensive parts of the Dutch Wadden Sea (Beukema 1976, Beukema et al. 1978 showed similar characteristics in biomass as well as in species composition of the bottom fauna. Therefore, and based on its broad range of sediment types and intertidal levels, the Balgzand tidal-flat fauna may be considered representative of the much larger tidal-flat areas of almost the entire Dutch Wadden Sea. Late-winter biomass values showed consistently increasing trends for the 1970s to 2000s, both in molluscs and polychaetes ( Table 2). The most spectacular increases (60% in molluscs as well as in poly-chaetes) were observed from the 1970s to the 1980s. In fact, all important (with regard to biomass contribution) species showed such increases (Table 2). Increases after the 1980s were limited to 25% per decade or (mostly) less. The strong increase around 1980 was discussed extensively by Beukema et al. (2002). Almost simultaneously, but with a time lag of 2 yr, the near-doubling of phytoplankton concentrations as well as primary production in the western Wadden Sea resulted in an almost doubling of the zoobenthic biomass (Beukema & Cadée 1986, Beukema 1991a, Beukema et al. 2002, Philippart et al. 2007). The increase in zoobenthic biomass was not found in parts of Balgzand with extremely harsh conditions, where food supply is probably not a limiting condition (Beukema & Cadée 1997a). Therefore, we concluded that the increase in zoobenthic biomass was caused by less intense competition for food resulting from the eutrophication of the area (Beukema et al. 2002).
Zoobenthic biomass as estimated from late-winter sampling apparently responded only to the initial increases that occurred around 1980, but not to the declines in phytoplankton characteristics that occurred after 1987 (Table 2). On the contrary, total late-winter biomass of zoobenthos further increased after the 1980s. This further increase was not only due to better survival of winter-sensitive species, but also to high biomass values in the last 25 yr of some invasive species (E. leei, Magallana gigas, Marenzelleria viridis). Apparently, at the arrival of these invasive species, the tidal-flat ecosystem was not yet occupied to carrying capacity for zoobenthos.
Carrying capacity for bivalves of a Wadden Sea tidal flat ecosystem was extensively studied by Beukema & Dekker (2019a). In that paper, we proposed a practicable definition for the ill-defined concept of carrying capacity: the numerical density (n m −2 ) at which production (g m −2 yr −1 ) is maximal. At low densities, production is limited by numbers, and it increases with numerical density. However, above an optimal density, production declines with further increases in density as a consequence of reduced growth rates at higher densities , 2019a. In practice, such high bivalve densities of over ~400 ind. m −2 during the growing season (resulting in seriously reduced growth) were rarely reached on the Balgzand tidal flats: only twice in the 39 years with available production estimates (Beukema & Dekker 2019a). Thus, in only about 5% of the observation years, bivalve abundance went beyond the carrying capacity of the ecosystem. In these few years, competition for food among bivalves would have been strong. Nevertheless, this severe competition with intensive grazing may have resulted in only 1 year (out of 37) with reduced chlorophyll concentration in the main tidal stream (Beukema & Dekker 2019a). One should realize, however, that phytoplankton is only part of the food supply for bivalves on Wadden Sea tidal flats. Both Kamermans (1994) and Hummel (1985) found that bivalves also consume microphytobenthos, and Jung et al. (2019) found high proportions of freshwater algae in their diet. This diversity of food sources may explain why growth rates in most bivalve species did not increase with increasing chlorophyll concentrations (Beukema et al. 2017a). The exception was L. balthica, in which growth rates increased with increasing chlorophyll concentrations (Beukema et al. , 2017a and in particular with increasing diatom concentrations (Beukema & Cadée 1991).
As contrasted to the late-winter biomass values, those found in late summer did not further increase in the 1990s and the following years (Fig. 3A). If anything, they showed a declining trend after 1992. Trends in summer values can be studied best from 1992 onwards, as biomass was seriously reduced in 1990 (and not yet completely restored in 1991) by intensive fishery, removing all mussel beds and nearly all cockles from the tidal flats. Unfortunately, no reliable data on total summer biomass are available for the initial years (before 1988), because samples were not dug sufficiently deep in the summers of these early years. As a consequence of the opposite trends in the winter and summer estimates of total biomass values in the last 28 yr, the differences between summer and winter values showed a declining trend (Fig. 3B). This decline of the seasonal biomass increment was statistically significant. It might reflect a declining trend in chlorophyll concentrations (Fig. 4), thus being a manifestation of the deeutrophication during this period.
In conclusion, clear eutrophication effects were observed for a short period (late 1970s to early 1980s) with simultaneous increases in nutrient concentrations, phytoplankton and chlorophyll concentrations (Fig. 4), primary production, and (several years later) late-winter biomass of bottom animals (Fig. 3A). The de-eutrophication after the mid-1980s was visible in nutrient and phytoplankton characteristics, but for zoobenthos biomass only in the summer values (Fig. 3A) and in the seasonal increments (Fig. 3B). Bottom-up regulation of this tidal-flat ecosystem does not appear to be strong, whereas top-down regulation was strong in several species of bivalves. As  showed, recruitment success in these important species was governed by the abundance of epibenthic predators (shrimps, shore crabs) which are affected by winter temperatures, and recruitment strength largely determined bivalve abundance in subsequent years and thereby their biomass (Beukema et al. 2010, Beukema & Dekker 2019a) and production (Beukema & Dekker 2019a). As a consequence of the non-response of the late-winter biomass of bottom animals to the recent de-eutrophication, food supply for shellfish-eating birds on the tidal flats in the most critical season (with minimal biomass and maximal need for food) did not show any long-term decline (Beukema & Dekker 2019b), and total zoobenthic biomass values in late winter even showed an increasing trend over the 50 yr period (Fig. 3A, Table 2).

Sea level rise
Sea level shows a rising trend. Along the Wadden Sea coast, this rise is at present limited to about 2 mm yr −1 . As yet, there appear to be no clear signs of an acceleration. In the Wadden Sea, sedimentation is currently sufficient to keep up with such low rates of sea level rise. Some additional concerns arose from gas extraction, but the resulting bottom subsidence (and thus relative sea level rise) proved to be very limited (Beukema 2002). Future rates of sea level rise are uncertain, but are generally expected to accelerate. Van der Spek (2018) predicted that at a certain rate (higher than present, at > 5 mm yr −1 ), sedimentation can no longer keep up. As a consequence, low tidal flats will get drowned and high flats will turn to low ones. Consequences of insufficient compensation of sea level rise (by extra sedimentation) for the benthic animals and their biomass were studied (Beukema 2002). On tidal flats, zoobenthic biomass is maximal at intermediate intertidal levels (immersed for 40 to 80% of the time), and declines coastward (higher levels) as well as farther off shore (lower levels). Net (uncompensated) sea level rise will cause a shift of this biomass maximum to a coastal direction: formerly low flats will become even poorer in benthos, and formerly high flats will become enriched. At substantial net sea level rise, intertidal surface areas will diminish by drowning of lower flats. Intermediate flats will turn to lower and poorer flats, while the area of high flats will diminish as the coastward shift is arrested by the presence of dikes along most of the coast. Therefore, total zoobenthic biomass is bound to decline if sea level rise seriously accelerates; however, as long as sedimentation can keep up with sea level rise, little will happen with total biomass (Beukema 2002).
Species that live predominantly in high tidal flats are most endangered by substantial uncompensated sea level rise, as ultimately, all of these flats will turn into lower flats. No habitat will then be left for species that are limited to higher flats. Examples of such species are A. tenuis (Dekker & Beukema 1993) and Corophium volutator (Beukema & Flach 1995). For the time being, A. tenuis populations are doing increasingly well as a result of the elevated summer temperatures, promoting their growth as well as recruitment (Dekker & Beukema 1993. High flats further act as nursery areas for juveniles of some species that migrate to lower areas later in life, such as L. balthica (Beukema & de Vlas 1989, Beukema 1993a) and lugworms Arenicola marina (Beukema & de Vlas 1979, Flach & Beukema 1994.
In conclusion, to date, no significant effect of sea level rise has been found, but in the future, serious problems might arise if sea level rise accelerates further.

Disturbance and recovery
Fishery appears to be the most marked disturbance on the Balgzand tidal flats. For a long time, digging for lugworms has been, and still is, a licensed activity on restricted parts of Balgzand, executed year-round by a few ships, although at present only 1 ship is active . Fishery for cockles C. edule occurred in about half of the years in late-summer/ autumn until 1993, when it was prohibited (Beukema & Dekker 2018). Fishing for mussels M. edulis was not frequent, but in the summer of 1990 all mussel beds were completely eliminated (Beukema 1993b, Beukema & Cadée 1996. After that disastrous year, this fishery was prohibited on tidal flats. Mechanical lugworm harvesting involved digging to a depth of 40 cm in a small part of Balgzand. This activity not only halved lugworm density in the fished area, but also total benthic biomass . Biomass reduction was caused in particular by the death of large M. arenaria, which accounted for almost half of the total biomass at the start of the fishing period. After the end of the fishing period, it took some 3 yr before total biomass and lugworm abundance reached pre-fishing levels and 5 yr before this was also the case in M. arenaria . Piersma et al. (2001) suggested long-lasting (at least 8 yr) negative effects from mechanical cockle dredging, including reduction of recruitment in C. edule as well as in L. balthica. However, more recent studies (Beukema & Dekker 2005, 2018 did not corroborate these conclusions. On Balgzand, we observed no difference in cockle recruitment after fishing between fished and unfished areas (Beukema & Dekker 2005). Further, survival of bivalve recruits that had settled a few months before the late-summer fishing period did not differ between years with and without fishing, nor between fished and unfished areas (Beukema & Dekker 2018). Thus, we did not find any evidence that the fishery negatively affected recruitment in C. edule, L. balthica, and M. arenaria. The reason for this contradiction with results reported by Piersma et al. (2001) may be that the main conclusion of that paper was based on only 1 comparison between 2 groups of years (viz. the periods 1992−1994 and 1996−1998). Data for 1991 were omitted for no stated reason, and this year was characterized by an above-average spat fall. The data for 1995 were omitted because data from just 1 station were lacking. Within the period 1990−1998, several comparisons between groups of years would have been possible and should have been included by Piersma et al. (2001). The authors did not give any reason for the choice of the single comparisons they reported. This raises the question whether the choice was purely by chance or a case of 'cherry picking'. An analysis of the available data to compare other periods did not yield any result that supported the conclusion from the one comparison chosen by Piersma et al. (2001) to be reported (J. van der Meer pers. comm.). The questionable conclusion by Piersma et al. (2001) on a long-lasting negative influence of mechanical cockle fishing is far from trivial, as it played a decisive role in the termination by the government of mechanical cockle fishery in the Dutch Wadden Sea. We were co-authors of the paper by Piersma et al. (2001), and we regret and apologize that we did not notice at that time that this publication was flawed.
Large-scale blue mussel M. edulis fishing occurred only once on Balgzand, with drastic consequences (Beukema 1993b, Beukema & Cadée 1996. After 3 successive years (1988,1989,1990) of recruitment failure in bivalves, stocks of adult bivalves had become low by 1990 across the Dutch Wadden Sea, and mussel farmers had serious need of mussels to seed their culture plots. With the help of cockle fishers, mussels were (almost) completely removed from tidal flats in the summer of 1990. As a consequence, zoobenthic biomass (and particularly bivalve stocks) were unprecedentedly low in early 1991. At that time, the few bivalves left showed above-average individual weights, probably as a result of superfluous algal food supply in the 1990/91 winter (Beukema & Cadée 1996). Food shortage in shellfish-eating birds (oystercatchers and eiders) resulted in increased mortality in alternative bivalve prey species and depletion of these food sources (Beukema 1993b). As a consequence of food shortage, eiders and oystercatchers left the intertidal area earlier in late winter/early spring than in other years (Beukema & Cadée 1996). Most of the bottom fauna rapidly recovered in the summer of 1991, accelerated by an extraordinarily high recruitment in several species (Beukema & Cadée 1996). However, across the Dutch Wadden Sea, mussel stocks remained low for a long period (9 yr), primarily as a consequence of infrequently occurring successful recruitment (van der Meer et al. 2019).
For a more detailed study of rates of recovery of macrozoobenthos on tidal flats, we killed all bottom animals by suffocation in a number of 120 m 2 quadrats and subsequently followed the recovery of the bottom fauna in these plots for several years ). Within half a year, nearly all species had settled again within the plots. After a summer had elapsed, total numbers of bottom animals equaled those in the surrounding areas. As nearly all of these newly settled species were of a very young age, recovery of biomass took a much longer time, i.e. several years. Immigration of adults into the plots was a slow process. In several species, recruitment of young stages proved to be higher within the formerly de-faunated plots than in the surroundings . Apparently, recruitment inhibition in well populated areas is a frequently occurring process. A similar observation of exceptionally high recruitment all over Balgzand was made in 1991 after a very sparse bottom fauna in early 1991 (Beukema & Cadée 1996).
In conclusion, except for 1 year (1990/91), disturbance by fishery on the tidal flats was of a limited nature. Recovery from prolonged lugworm digging and mussel-bank removal took several years.

Diversity
During our long observation period on the Balgzand tidal flats, the number of macrozoobenthic species found increased by at least 30% (Beukema & Dekker 2011b). During the 40 yr period 1970−2009, the number of species encountered annually rose from an average of 28 to a mean number of 38, this is at a rate of 0.25 yr −1 . The number of species encoun-tered per sampled m 2 increased from about 13 to about 18. Increases were caused partly by invasions of exotic species and partly by movements of species from subtidal locations to tidal flats. A further cause of the increase in species numbers per m 2 was the increasing trend in the abundance of several individual species (increasing their chance to be included in a sample of limited size), as particularly observed in species that are sensitive to cold winters. Thus, 3 groups of species contributed, to a roughly similar extent, to the overall increase in species richness: aliens, local invaders, and species with increasing abundance.
Invasive immigrants included: E. leei, M. viridis, M. gigas, and Hemigrapsus takanoi (Beukema & Dekker 2011b). On Balgzand, the first specimen of E. leei (E. directus, E. americanus) was observed in 1982 (Beukema & Dekker 1995). Numbers of E. leei increased after 1991, but high abundance was found only after 2004 . Spat of this bivalve was found over most of Balgzand, but it survived only near the low water mark (in the transition zone from intertidal to subtidal). In this zone, adult densities of around 100 m −2 and biomass values of around 50 g AFDM m −2 were reached ). In these low, sandy, and exposed parts of Balgzand, the species strongly enriched the benthic fauna (increasing biomass by an order of magnitude) in 2006−2009. Growth was usually rapid and was negatively related to numerical density . M. viridis has been found on Balgzand since 1989 and rapidly spread in the early 1990s, to occur nearly everywhere by 1996 (Beukema & Dekker 2011b). It reached maximal numbers in the late 1990s. Biomass values were high (about 3 g AFDM m −2 in winter and 7 g AFDM m −2 in summer) during the short period 1999−2003, to decline to low values after 2003 (Dekker 2012; Table 2). The species thus showed a typical pattern of population development known in invasive species: a lag-phase following the first observation, an explosive increase, a stabilization, and a phase of decline, just as observed in this species some 5 to 10 yr earlier in the Dollard estuary by Essink & Dekker (2002). M. gigas arrived around 2000 on Balgzand and established on mussel beds. As yet, it is mostly restricted to a few sampling sites with mussel beds.

H. takanoi has been observed on Balgzand since 2005 (Beu kema & Dekker 2011b).
In earlier papers (Beukema & Dekker 2011a), we criticized a study by Compton et al. (2008), who used numbers of species in very small samples (of 0.02 m 2 ) to characterize local species richness. We showed that numbers in small samples do not predict species numbers in larger areas. On the Balgzand tidal flats, species numbers in small samples were relatively high at near-shore sites and low at offshore sites. For larger samples, the reverse was true. This reversal was caused by the very high densities (n m −2 ) occurring in near-shore areas in a few species, such as C. volutator and Peringia (Hydrobia) ulvae. High densities increase the chances of a species to be captured in a small sample. Such high-density species are generally lacking in off-shore areas.
Here, the number of species in larger samples was significantly higher than in near-shore areas. In conclusion, more species live in off-shore areas, but nearly all of them are present at such low densities that only large samples can fully represent them.
To summarize, the species richness of the tidal-flat fauna has increased, in particular by introduction of alien species and by the declining frequency of cold winters, promoting both success of migration from subtidal areas as well as abundance of species already present.

Productivity
Bottom animals play a significant role in the Wadden Sea ecosystem, both as consumers and as producers. Therefore, it is well worth to attempt a quantification of their role in the functioning of the ecosystem. At the very start of the work of J.J.B. on marine benthos, in the late 1960s, his research task was to estimate the benthic secondary production of the Wadden Sea, within the framework of the then developing International Biological Program (IBP). It soon became clear that production estimates in populations of bottom animals showed high year-to-year variability. To obtain representative values, J.J.B. launched a long-term monitoring programme for macrozoobenthic species in a limited part of the Wadden Sea: the 50 km 2 Balgzand tidal-flat area. This programme yielded the necessary data for estimates of annual production by several bivalve species. The advantage of the choice for bivalves for these studies is the ease of assessing age in individual specimens by counting the year marks of their shells. For each sampling, the samples thus yielded for each year class present separate data on their n m −2 and on their mean individual weights. Multiplying mean density by mean weight change over a certain period gives of an estimate of production for this period. The usual way of calculating secondary production, from mean numerical density and individual weight changes, in fact underestimates true production if parts of the animals are lost to predators and subsequently regenerated, at the expense of body growth. De Vlas (1979Vlas ( , 1985 showed that such cropping amounted to substantial parts of total production of the cropped species, e.g. lugworm A. marina tails and L. balthica siphons. We reported estimates of annual soft-part production for long periods (several decades) for some important bivalve species: C. edule (Beukema & Dekker 2006), L. balthica , M. edulis , and E. leei . The observed long-term means of annual production of these species on Balgzand amounted to about 7, 3, 5, and 1 g AFDM m −2 yr −1 , respectively. The summed bivalve production of 17 g AFDM m −2 yr −1 represents 4 or 5% of the annual primary production of the area (Beukema & Dekker 2019a). In all of these species, positive estimates were limited to the growing season (March/ April to July/September). In the autumn/winter period, the animals usually lost weight (Honkoop & Beukema 1997), and production estimates turned to negative values. Bivalve production was maximal at intertidal levels between mid-tide and the low-water level, where (in most species) growth rates and also adult densities were maximal. Over almost the entire range of numerical densities, production estimates increased with numbers of bivalves and with magnitude of recruitment in preceding years (van der Meer et al. 2001, Beukema & Dekker 2019a). Only at (rare) extremely high densities, production declined as a consequence of growth reduction (Beu kema & Dekker 2019a). As this reduction was ob served in all bivalve species studied, it likely had a common cause, namely lack of food as a result of severe competition.
The most proper unit to express production rates to be used in studies of energy budgets and energy transfer is not g m −2 , but cal m −2 or J m −2 . Therefore, we studied calorific values of dried soft parts of the bivalve L. balthica. Based on chemical composition (Beukema & De Bruin 1977) as well as direct calorimetry, calorific values of this material were determined and possible causes of the consistent differences between the outcomes of the 2 methods were discussed (Beukema & De Bruin 1979). The best estimates for the annual mean calorific value were 5.5 and 5.6 cal mg −1 for the direct (combustion in a micro-calorimeter) and indirect (calculation from chemical composition) method, respectively. A review paper on publications of calorific values of marine bivalves and other marine invertebrates (Beu -kema 1997) revealed that a high proportion of the reported records were improbably low. Based on acceptable records, we concluded that 5.5 cal (or 23 J) mg −1 may be the most representative value.
Estimates of the production of shell material on the Balgzand tidal flats by L. balthica (Beukema 1980) and C. edule (Beukema 1982) amounted to annual means of 13 and 118 g m −2 yr −1 , respectively. Mean total shell production on tidal flats of the entire Dutch Wadden Sea for these 2 species was estimated at 11 and 156 million kg yr −1 , respectively. As part of these shells are either fragmented or too small to be exploitable, an estimated average of 88 million kg of cockle shells was added annually to the exploitable stock of the Dutch Wadden Sea (Beukema & Cadée 1999). This amount did not differ much from the amounts annually extracted from the Wadden Sea stocks by shell fishing. We advised Rijkswaterstaat to continue their present policy of issuing permits for shell fishing based on mean annual shell production (Beukema & Cadée 1997b).

Population dynamics
Population sizes fluctuated from year to year and from place to place, but in some species more than in others (Beukema et al. 1983). Suspension feeders (such as most bivalves) fluctuated more than deposit feeders (such as most polychaetes). Generally, species fluctuating strongly from year to year also fluctuated strongly from place to place (Beukema et al. 1983). A possible explanation is that deposit feeders rely on a more stable and evenly distributed food source (competing for microphytobenthos and detritus in the immediate surrounding area) than suspension feeders do, being dependent on intermittent superfluous pulses of phytoplankton irregularly supplied by tidal currents. Thus, competition with neighbours will be more important in deposit feeders than in suspension feeders. More even distributions indeed appear to be realized in deposit feeders. Moreover, in suspension feeders, recruitment may generally be more variable and less gregarious than in deposit feeders.
If the same environmental factor governs the abundance of 2 or more species, the abundances of these species will fluctuate in parallel over the area where this factor is exercising its influence. Such an environmental factor may be the winter conditions of a given year: severe winters not only kill high proportions of sensitive species , 1989, Beukema & Dekker 2011b), but also promote bivalve recruitment in the following spring ; Table 1). Winter conditions in any year are similar over extensive ranges (up to 2000 km) along the European coasts . As a consequence of their similar responses to winter temperatures, a high proportion of benthic species (including some bivalves as well as some polychaetes) showed synchronized population sizes throughout the Dutch and part of the German Wadden Sea . If a certain local population does not share the general large-scale fluctuation pattern, this aberration may point to a local (non-natural) disturbance . The study of synchrony in fluctuation patterns may then offer a way to distinguish between natural and manmade changes.
In the 4 most important bivalve species on Wadden Sea tidal flats (C. edule, M. edulis, M. arenaria, and L. balthica), annual recruitment was synchronized over a vast area of the Wadden Sea as a consequence of their similar response to winter conditions and subsequent abundance of epibenthic predators . Synchrony was observed not only between distant (100−200 km or even longer distances) populations of the same species (Beukema et al. 2001a(Beukema et al. , 2017b, but also between the 4 studied species living together in the same area (Beukema et al. 2001a). For several reasons, this recruitment synchrony did not manifest in a possible synchrony of annual bivalve biomass values of the various species. Reasons for this included: absence of between-species synchronization in annual survival rates (Beukema & Dekker 2019b), multi-age composition of species biomass (Beukema et al. 2001a), and between-species differences in ages at which maximal biomass values are reached (Beukema et al. 2001a). Consequently, available food sources for birds showed fluctuations that were not as heavy as they would have been if all recruitment synchronization was passed fully to biomass (Beukema & Dekker 2019b). Biomass synchronization was present throughout the Wadden Sea tidal flats in a number of wintersensitive species (including the bivalve C. edule and the polychaetes N. hombergii and L. conchilega), resulting in relatively low food supply for some bird species at the end of severe winters (Beukema et al. 1993). After the coldest winters, the average zoobenthic biomass on the Balgzand tidal flats amounted to around 20 g AFDM m −2 , compared to around 40 g AFDM m −2 after the mildest winters (Fig. 2).
An important mechanism explaining the dynamics of populations is the relationship between the size of the adult stock and annual recruitment strength.
In the cockle C. edule, recruit densities showed a negative relationship with sizes of adult stock (Beukema & Dekker 2018). The abundance of cockle adults also showed negative relationships with recruit numbers of 2 other bivalves: L. balthica and M. arenaria (Beukema & Dekker 2018). In contrast, populations of L. balthica showed a positive relationship between the size of the adult stock and the density of recruits (Beukema et al. 2017c), but this relationship was limited to a small range of very low adult densities.

Species interactions
In a study on relationships between fluctuation patterns in the biomass values of 14 important benthic species on Balgzand, Williams et al. (2004) revealed several interactions, most of them being negative. The strongest effects found were between the abundance of the predatory polychaete N. hombergii and its main prey species, i.e. smaller polychaetes such as Scoloplos armiger and Heteromastus filiformis. As these negative interactions had already been described in earlier papers (Beukema 1987, the above negative relationships were expected. N. hombergii is winter-sensitive, in contrast to its main prey species. As a consequence, the abundance of prey was particularly high in summers after cold winters . Abundance of the predator proved to be governed primarily by winter temperatures and not by food supply, whereas prey abundance was governed primarily by predator abundance , van der Meer et al. 2000. Despite the good correlation found between winter temperatures and N. hombergii abundance for the period before 2000 (van der Meer et al. 2000), the predictive value of high winter temperatures for N. hombergii abundance after 2000 proved to be poor (van der Meer et al. 2013). Apparently, some (unknown) factor prevented the expected further population growth in this species. Williams et al. (2004) reported only 1 highly significant positive interaction: between A. marina and L. balthica. This might have been caused by a common behavioural trait: juveniles of the 2 species migrate in winter from high to low tidal flats (Beukema & de Vlas 1979, 1989, Beukema 1993a, Flach & Beukema 1994 and do this for higher proportions in cold than in mild winters (Flach & Beukema 1994. Surprisingly, the exercise described by Williams et al. (2004) did not reveal strong negative interactions of A. marina or C. edule with several other species. Such interactions were observed in field experiments by Flach (1992Flach ( , 1996.

Model testing
Distribution patterns of macrozoobenthic species on Balgzand were used (Beukema 1988) to test the abundance/biomass comparison (ABC) model developed by Warwick (1986) to detect pollution effects on marine zoobenthic communities. As a consequence of the frequent occurrence of very high numerical densities of 2 small-sized species (P. ulvae and C. volutator), particularly at high intertidal levels, the ABC method would indicate such areas as polluted. There are no other indications of pollution on Balgzand and, moreover, the observed high densities of these species in coastal areas have existed for a very long time. We concluded that the ABC method is useless for assessing the pollution state of Wadden Sea tidal-flat communities.
In a recent study , we tested a model developed by Brey (1990) to estimate secondary production from biomass values and simultaneous mean individual weights, using a 33 yr series of direct bi-annual estimates of production, biomass, and mean individual weights of the Balgzand population of L. balthica. The model satisfactorily predicted production on the basis of estimates obtained by late-winter samplings, but the late-summer samplings resulted in overestimates (roughly a doubling) compared to the true production. The base of the model, namely the derivation of estimates of instantaneous rates of mortality (and thus values for the production:biomass ratio) from the mean individual weights appeared to be sound. However, this derivation did not account for year-to-year variation in mortality rates at the same individual mean weight. Thus, for a sufficiently precise estimate of production, we consider real estimates of annual mortality to be indispensable. As the model was designed for communities rather than individual species, a further study based on estimates of groups of species would be necessary to justify the use of the Brey model.

MAIN CONCLUSIONS
The Balgzand benthic ecosystem developed favourably during the 50 yr monitoring period: species numbers (Beukema & Dekker 2011b) as well as latewinter biomass values Fig. 3A) showed rising trends by 0.25 yr −1 and 0.5 g AFDW yr −1 , respectively. Some environmental factors showed long-term trends as well: water temperatures rose by about 2°C (Beukema & Dekker 2020), and chlorophyll concentrations in the main tidal stream steadily declined after an increase until the mid-1980s (Philippart et al. 2007).
Several effects of climate change were detected: (1) It contributed to the rise in species numbers (Beukema & Dekker 2011b) and biomass (Fig. 2) encountered on the tidal flats.
(2) It promoted the abundance of shrimps and shore crabs in spring (Beukema 1991b, 1992b, Beuke ma & Dekker 2014. (3) These increases in epibenthic predator abundance reduced recruitment success in several bivalve species ), but possible effects on bivalve biomass were compensated by better survival and/or growth (Beukema & Dekker 2019b; Table 1).
(4) It advanced the periods of spawning and settlement in spring of both predators and their prey to roughly the same extent , so no mismatches arose.
Effects of eutrophication were particularly clear around 1980, when primary production and chlorophyll concentrations roughly doubled within a few years (Philippart et al. 2007), and doubling of zoobenthic biomass followed with a time lag of 2 yr (Beukema & Cadée 1986, Beukema et al. 2002. This biomass increase was observed only in areas with already high biomass (Beukema & Cadée 1997a), pointing to relief of food limitation (Beukema et al. 2002). The de-eutrophication occurring after the mid-1980s (with chlorophyll concentrations declining gradually to a halving in recent years) was reflected only in zoobenthic biomass values obtained in summer and in the annual seasonal biomass increments from winter to summer (Fig. 3). High biomass values of introduced species will have added to the lack of effect of de-eutrophication.
Sea level rise was small during the period 1970− 2019, amounting to a total of about 10 cm, and was compensated by extra sedimentation in most parts of Balgzand. Thus far, none of the changes in zoobenthic populations could be attributed to sea level rise. Disturbance by the lugworm fishery was clear, but it was local. The biomass reduced by lugworm digging required several years for complete recovery (Beu kema 1995). Disturbance by mussel fishery occurred only once (in 1990) and had severe effects (Beukema 1993b, Beukema & Cadée 1996. Disturbances by the frequent, but meanwhile terminated, mechanical cockle fishery were restricted to local depletions of adult stocks, but had no demonstrable long-term effects (Beukema & Dekker 2005, 2018. Immigration of species (either exotic or from nearby deeper waters or more southern areas) contributed significantly to the richness of tidal flats (Beukema & Dekker 2011b). Bivalve productivity was high (Beukema & Dekker 2019a). It was mainly dependent on recruitment success in a few preceding years (van der Meer et al. 2001, Beukema & Dekker 2019a), but hardly on availability of algal food (Beukema & Dekker 2019a). Apparently, topdown processes dominated bottom-up ones as causative mechanisms for production variability on a species level (bivalves). Population dynamics of several species were strongly affected by the magnitude of their annual recruitment and subsequent survival. In most bivalves, the same environmental factor (winter conditions) governed re cruitment success , resulting in synchronized year-to-year fluctuations over extensive areas in numeric abundance of young stages (Beukema et al. 2001a). Such synchrony, however, did not extend to adult populations. Therefore, it did not affect fluctuations in benthic biomass or in food available for birds (Beukema & Dekker 2019b).
To date, the Balgzand monitoring programme has seldom been criticized. Only Kuipers et al. (1981) made some remarks. They emphasized that (1) smaller organisms are wrongly ignored in the Balgzand monitoring, (2) macrozoobenthos uses only a relatively small part of the available primary energy, and (3) the non-studied mortality of young benthic stages of macrobenthos (rather than physical factors and food competition) might control macrobenthic production. Although most of these points are basically valid, it would have been unfeasible to include micro-and meiobenthos (and the young stages of macrobenthos belonging temporarily to these groups) in the already extensive programme.

EPILOGUE
We are now looking back on our life's work of 50 yr of monitoring efforts. The resulting unique, uninterrupted, and consistent 50 yr data series recorded in particular what happened to an intertidal soft-bottom community during a period of climate change. Rising temperatures explained a major part of the observed ecosystem changes. This major conclusion is similar to the one for the long-term research in the English Channel, which also indicated the overwhelming importance of climate change (manifested as temperature change) for a variety of ecosystem changes (Southward et al. 2005).
Prospects for the future of the Wadden Sea tidalflat ecosystem appear favourable. Faunal diversity will keep increasing by introduction of new species and continuation of climate warming, reducing overwinter mortality in a high proportion of the species. Biomass will probably stabilize at a high level as a result of opposing forces: higher temperatures enhance survival and growth, but reduce bivalve recruitment, and de-eutrophication lowers algal concentrations. The main thread may be the accelerated rise in sea level, resulting in drowning of the tidal flats. We recommend a minimum of interventions in the management of this World Heritage region and suggest allowing it to develop along natural lines. realization of our programme. Large institutes such as NIOZ bear a great responsibility for continuous funding and execution of long-term monitoring programmes. The fate of too many monitoring programmes has been an untimely premature end caused by lack of continued funding (Duarte et al. 1992). Unfortunately, we have to mention that the present director decided to terminate our uniquely long data series in 2020 with the near-retirement of R.D. (20 years after that of J.J.B.). We are grateful to Pauline Kamermans (WUR), P. Jacobs (NIOZ) and 3 reviewers (Karsten Reise, Erik Bonsdorff, Karin Troost) for their constructive comments. We thank numerous assistants and students for their dedicated help in field and laboratory work. The data of the project are stored in files at NIOZ and are accessible. They will prove to be a rich treasury that is far from being exhausted by the papers we have published to date.