The question of ecosystem stability is central to ecology and paleoecology and is of particular importance for estuaries, which are environmentally highly variable, considered as geologically short lived, and among the most degraded modern ecosystems of our planet. Understanding their ecological dynamics over geological time scales requires paleontological data in a sequence stratigraphic framework, which allows evaluation of paleocommunity dynamics in an environmental context. A 445-m-thick estuarine succession in a satellite basin of the Vienna Basin (Austria) shows continuous sedimentation over 700 k.y. and can be divided into two transgressive systems tracts and a highstand systems tract. In contrast to expectations, no major physical disturbances of the ecosystem involving abrupt changes in diversity and biofacies composition occurred at flooding surfaces and at the sequence boundary. Accommodation space remained remarkably constant over the depositional history of the basin, and all changes between depositional environments were therefore more or less gradational. Biotic change along the studied succession can be described as a gradual faunal replacement in response to habitat tracking, a process also reported for some normal marine shelf environments. Benthic assemblages in the estuarine succession were strongly dominated by a few taxa and developed along two indirect gradients, water depth and hydrodynamic energy. These gradients show subtle long-term trends, corresponding to the sequence stratigraphic architecture. Tectonics affected the sequence architecture in this particular marginal marine setting: it controlled accommodation space and sedimentary input, and provided stable boundary conditions over hundreds of thousands of years. Our study demonstrates for the first time that estuaries, which are under great environmental pressure today, are resilient to natural environmental perturbations and can persist over geological time scales.

The fossil record is rich in examples of ecosystem stability exceeding millions of years, mostly in Paleozoic and Mesozoic marine level bottom communities, interrupted by brief episodes of strong turnover (DiMichele et al., 2004). Ecosystem stability at the scale of tens to hundreds of thousands of years is a typical feature of Pleistocene coral reefs (Pandolfi and Jackson, 2006), but has only rarely been reported from Phanerozoic soft bottom assemblages (Holterhoff, 1996). High-resolution studies of this kind require a stratigraphic framework to evaluate paleocommunity dynamics in their environmental context (Holland and Patzkowsky, 2004; Scarponi and Kowalewski, 2004). Here we report for the first time on ecosystem stability over 700 k.y. in a Miocene estuarine succession (Fig. 1).

Estuaries are semienclosed coastal water bodies with salinities that differ from the open sea. They exhibit distinct biota (Whitfield and Elliott, 2011), are environmentally highly variable (Elliott and Quintino, 2007), are considered as geologically rather short lived (Wolff, 1983), and belong to the most degraded modern ecosystems of our planet (Lotze et al., 2006). Because of their highly variable physicochemical characteristics, estuaries also count as naturally stressed areas: their biotas have the ability to adapt to various stressors and the ecosystem can compensate for changes in the environment, a feature termed homeostasis (Elliott and Quintino, 2007). The study of turnover and community stability in such a fluctuating ecosystem over hundreds of thousands of years is therefore of major ecological and paleoecological interest.

The Korneuburg Basin, a satellite basin of the Vienna Basin, is located in Lower Austria and is 20 km long and 7 km wide (Fig. 1). It was formed by Burdigalian pull-apart movements within the Alpine-Carpathian thrust belt. The basin sediments mostly belong to the Miocene Korneuburg Formation, which was dated into nannoplankton zone NN4, paleomagnetic chron C5C, and mammal zone MN5 (Harzhauser and Wessely, 2003). Faunal composition and stable isotopes show that the Korneuburg Basin was divided into a southern estuarine part, probably connected to a huge, south-north–trending river system in the adjoining Vienna Basin, and a northern marine part (Harzhauser et al., 2002; Latal et al., 2006). A connection to the large epicontinental Paratethys Sea was most probably established along the northeastern margin of the Korneuburg Basin (Fig. 1) (Harzhauser and Wessely, 2003). The foraminifera in parts of the basin were adapted to brackish-water conditions and indicate a very shallow water environment through basin history, with maximum water depth not exceeding 30 m (Rögl, 1998). The fish fauna (e.g., Gobiidae, Sparidae, Dasyatidae, Myliobatidae) in the southern Korneuburg Basin indicate littoral and shallow sublittoral conditions in a subtropical to tropical environment with freshwater influence (Reichenbacher, 1998).

The section near the village of Stetten in the Korneuburg Basin can be divided into Stetten West (total thickness of 445 m) and Stetten East (total thickness of 165 m) (Fig. 1; Fig. DR1 in the GSA Data Repository1). Thin coaly deposits, washed-in land-, and freshwater snails indicate marginal marine conditions. Sand packages with trough cross-bedded sets are mostly interpreted as tidal sand waves of the shoreface. Pelitic sediments mostly show even lamination to wavy bedding or thinly alternating sandy and muddy layers, indicative of tidal flat deposits (Reineck and Singh, 1975; Boyd et al., 2006). The succession is dominated by upward-fining parasequences, which are typical for a tidal flat to subtidal environment on a muddy, tide-dominated shoreline (Van Wagoner et al., 1990). The astronomical tuning of the gamma ray record fixed the deposition to the time interval from 17.0 to 16.3 Ma, yielding a total duration of 700 k.y. (our data).

We took a total of 118 bulk samples, weighing ∼1 kg each, from 96 shell beds. The sediment material >1 mm mesh size was quantitatively picked under a binocular microscope for all biogenic components. For sponges, every bioeroded biogenic hard part with distinct traces of the ichnogenus Entobia was counted as a sponge colony. Each balanid plate, decapod claw, echinoid spine, tooth, and otolith, and for molluscs every shell (gastropods, scaphopods) and every isolated valve (bivalves), was counted as an individual. In addition to these quantitative samples, each layer was scoured for large-sized species, which are likely underrepresented in bulk samples, and a species was added as present to the data matrix where appropriate. Our approach biases relative importance toward taxa with disarticulated body parts. The patterns of our study, however, are consistent when corrections for disarticulation are applied. Abundance and species richness are strongly dominated by molluscs (Table DR1 in the Data Repository). Each sample was classified as either pelitic or sandy and assigned to a sequence and systems tract (ST) (Table DR2).

Species contributing <0.1% to the total assemblage, samples with <20 individuals, and one sample consisting solely of washed-in river snails were removed from the final data matrix (108 samples, 45 species, and 97.6% of the original data). To standardize samples, percentages of the abundances were calculated and square root transformed to de-emphasize the influence of the most abundant taxa.

The benthic assemblages were explored for differences between sediment composition and STs using analysis of similarity (ANOSIM). Because of the large sample size, differences between sample groups were significant, although low R-values indicate only poor discrimination. Results of ANOSIM are virtually identical when corrections for disarticulation are applied. Two-way cluster analysis (Ward’s method) was applied to characterize biofacies and to detect hierarchical groupings within the data set; detrended correspondence analysis (DCA) and nonmetric multidimensional scaling were used as complementary ordination methods to detect ecological gradients (Patzkowsky and Holland, 2012, p. 66). Cluster analysis and ANOSIM are based on the Bray-Curtis similarity coefficient. Logarithmic scale rank abundance plots and rarefaction curves were used to compare community organization between the total assemblages of STs and biofacies. Mean rarefied species richness of samples (n = 20, according to smallest sample size in the data set) was used to track diversities within and between STs. Statistical analyses were performed with the software packages PRIMER version 6.1.6 (Clarke and Warwick, 1994) and PAST (Hammer et al., 2001).

Four overlapping biofacies are developed along two major environmental gradients (Fig. 1; Figs. DR2 and DR3). DC1 is interpreted as water depth, because taxa with low scores (<130) typically inhabit subtidal environments, and those with higher scores inhabit the intertidal or very shallow subtidal environments (Appendix DR1 in the Data Repository). The complex indirect gradient water depth is therefore the dominant source of faunal variation in this study and most likely reflects the combined influence of many environmental parameters that are highly correlated (Patzkowsky and Holland, 2012, p. 56). Sediment grain size is a poor predictor for the composition of benthic assemblages (R = 0.12, p = 0.001; Fig. DR3) and DC2 is therefore best interpreted as hydrodynamic energy, where taxa with high scores (>100) prefer higher water energy and firm or hard substrata, and taxa with low scores (<100) prefer lower water energy and soft substrata. Accordingly, low DC1 and DC2 sample scores represent quiet water and mostly pelitic subtidal environments, which are characterized by the Nassarius-Turritella-Corbula biofacies (NCTB). High DC1 and low DC2 values represent quiet water environments of an inner tidal flat, characterized by the Agapilia biofacies (AB). Low DC1 and high DC2 scores, in turn, represent relatively high energy environments of the shoreface, mostly characterized by sandy substrates and the Nassarius-Paphia-Loripes biofacies (NPLB). High DC1 and DC2 sample scores represent somewhat higher energetic conditions of the outer tidal flat, characterized by the Granulolabium-Agapilia biofacies (GAB) (Fig. 1; Fig. DR3).

All four biofacies occur in both sections. The AB mostly occurs in the highstand systems tract (HST) and the GAB mostly in transgressive systems tract 1 (TST1) of the lower sequence. A few samples of the GAB were also present in the HST of the lower sequence, but it was not present in TST2 of the upper sequence. This is probably because the lower (i.e., shallower) part of this ST, where this biofacies would be expected to occur, is missing. The NCTB is almost exclusively present in the two TSTs. Similarly, the NPLB is almost exclusively restricted to the two TSTs (Fig. 1; Table DR2).

Benthic assemblages are barely separable between STs (R = 0.189, p = 0.001; Fig. DR3), and DCA axes scores and diversity of samples show only weak trends across the studied sections. DC1 scores decrease from bottom to top of TST1 in both sections. Apart from some distinct outliers, they are uniformly high in the HST and decrease again in TST2. Along the transect DC2 scores are similarly high in the two TSTs, but typically lower in the HST. Rarefied species richness of samples is more variable in the two TSTs and often above 8; it is less variable and, with few exceptions, below 7 in the HST (Fig. 1; Fig. DR1). STs and biofacies all show strong dominance of a few taxa, especially in the GAB and the AB and in TST1 and the HST (Fig. 2). Evenness of the total assemblage of the two TSTs and of the NCTB and the NPLB is much higher than that of the HST and the GAB and the AB.

Biofacies changes are usually gradational within parasequences, but abrupt across flooding surfaces, typically either within TSTs or when a sequence boundary is associated with a major shift in facies (Holland, 2000). The ordination in this study (Fig. 1), however, indicates no disjunct change in biofacies composition within the TSTs or across the SB. Our interpretation is that accommodation space remained remarkably constant over the depositional history of the basin. Water depth varied only within the narrow limits of few tens of meters, and all changes between depositional environments are therefore more or less gradational. Accordingly, biodiversity also changed only little over the studied succession, and most changes can be related to the presence or absence of environments within STs. The biotic change along the studied succession is termed faunal replacement, a process that can only be observed in continuous sections. Such replacement typically occurs over thousands to hundreds of thousands of years and most likely reflects habitat tracking, where species migrate to preferred habitats along environmental gradients (Brett et al., 2007). Ecosystem stability here is exemplified by recurrence of biofacies in STs. The strong dominance of a few taxa in all STs and biofacies indicates that the benthic habitats covered in our study are physically controlled. This is typical for estuaries, where the environment fluctuates constantly (Elliott and Quintino, 2007). As can be expected in such settings, evenness is higher in the two subtidal biofacies and in the two TSTs.

During the studied interval the strong tectonic reorganization of the central Paratethys resulted in rapid subsidence of the Korneuburg Basin (Wessely, 1998). Tectonics therefore affected sequence architecture in this particular marginal marine setting by controlling subsidence and sedimentary input, and by providing stable boundary conditions over 700 k.y.

Ecosystem stability over geological time scales is shown here for the first time for environmentally highly variable estuaries. We demonstrate that under appropriate tectonic boundary conditions, these ecosystems, which are strongly affected by anthropogenic impact today and considered as geologically short lived, can persist over hundreds of thousands of years, a pattern previously only reported for some normal marine level bottom assemblages of the shelf and coral reefs. The stability of such a low-diversity ecosystem over such a long time and across flooding surfaces and a sequence boundary may support the hypothesis of homeostasis formulated for modern estuaries.

We thank S. Ćorić, K. Kleemann, A. Kroh, K. Rauscher, O. Schultz, and W. Sovis for help with field work and species identifications; P. Strauss for help with the sequence stratigraphic interpretation; S. Dominici, M. Stachowitsch, and A. Tomasovych for discussions; and C.E. Brett, S. Holland, and M. Kowalewski for stimulating reviews. This study was supported by project P19013-B17 of the Austrian Science Fund (FWF).

1GSA Data Repository item 2014002, Table DR1 (abundance and species richness of higher taxa), Table DR2 (environmental and stratigraphic affiliation of samples in biofacies), Figure DR1 (section Stetten East), Figures DR2-DR4 (results of two-way cluster analysis and ordination [nMDS]), and Appendix DR1 (actualistic comparison of abundant taxa), is available online at, or on request from or Documents Secretary, GSA, P.O. Box 9140, Boulder, CO 80301, USA.