Late-Holocene variability in chironomid functional assemblages and carbon utilization in a tundra lake food web

High latitude freshwater systems are facing changes in catchment-mediated allochthonous input, as well as physical and chemical controls triggered by on-going climate change, which may alter their carbon processing and ecological characteristics. To explore changes in chironomid functional responses and carbon utilization in relation to long-term environmental change, we studied a sediment core covering ca. 2000 years from a tundra lake in northern Finland, which was analysed for sediment geochemistry, isotopic composition of chironomid remains and their functional assemblages. We aimed to relate changes in chironomid functional feeding assemblages and resource utilization, based on Bayesian stable isotope modelling, and determined that the long-term resource utilization was more controlled by sediment geochemistry (resource availability) and climatic variables, reflecting changes in habitat and lake ontogeny, rather than the functional feeding assemblage composition. Change horizons were observed for both sediment geochemistry and functional assemblage composition. However, different timing of these changes suggests different drivers affecting the dynamics of primary production and chironomid community functionality. We also compared the recent warming period to Medieval Climate Anomaly (MCA), observing divergent patterns, which suggests that MCA may not be a good analogue for changes induced by on-going climate warming.


Introduction
The ecological, physical and chemical characteristic of northern small and shallow lakes are especially susceptible to environmental change, including Handling editor: Jasmine Saros and Koen Martens Electronic supplementary material The online version of this article (https://doi.org/10.1007/s10750-019-04151-7) contains supplementary material, which is available to authorized users. changes in climate conditions, catchment processes and human impact (Smol, 2016). These particularly common freshwater systems often feature low-nutrient conditions, short food webs and dominance of benthic production (Christoffersen et al., 2008;Rautio et al., 2011). However, climate warming and associated increased (allochthonous) catchment input may cause functional changes in high latitude aquatic ecosystems, such as boosting of primary production and reorganizing of community structure (Smol et al., 2005;McGowan et al., 2018), as well as alteration in carbon processing (Tranvik et al., 2009;Mariash et al., 2018;Wauthy et al., 2018). Understanding these changes in energy flow and ecosystem dynamics is essential for assessing the future of northern freshwaters. Here, we explore the theme of carbon processing and functional changes with paleolimnological methods to gain temporal insight. While paleolimnology is a well-established tool in research of long-term environmental change (Smol & Douglas, 2007), aspects of ecological functioning have only recently emerged in paleoecological studies (Gregory-Eaves & Beisner, 2011;Luoto & Nevalainen, 2015;Luoto & Ojala, 2018). Defining ecological structures according to similar biological functions, rather than taxonomic relationships, and exploring them in a spatio-temporal context may open new insight into understanding ecosystem dynamics, such as ecosystem's resilience, functional responses to environmental perturbations and biodiversity changes (Gregory-Eaves & Beisner, 2011;Oliver et al., 2015).
Chironomid larvae (Diptera: Chironomidae) are key players in benthic carbon processing of high latitude lakes: they participate in recycling of nutrients and decomposition of organic matter (Hansen et al., 1997;Vanni, 2002), and where present, provide carbon transfer (nutrition) to higher trophic levels (Sierszen et al., 2003). At the same time, our understanding of feeding behaviour of chironomids is limited, especially considering the selectivity and plasticity expressed in their feeding, and the driving factors triggering adaptations in their feeding habits (Berg, 1995). For secondary producers, habitat attributes (such as oxygen levels, substrate organic matter quantity and quality) and physiological factors (body size, developmental stage) have been noted important in defining the feeding mode and resource utilization (Merritt et al., 2008;Craig et al., 2015). In addition, Kivilä et al. (2019) recently showed that occurrence of chironomid functional feeding groups (FFGs) was connected to habitat-related factors, such as DOC and nutrient concentrations. In this study, we use FFGs, which are defined by feeding mode but also have their general preferred type of nutrition (Mandaville, 2002;Merritt et al., 2008), to capture temporal patterns in functional assemblage composition and feeding modes, and assess their connections to environmental changes and macrobenthic resource utilization.
Resource utilization is intimately linked with resource availability, however, the nutritional quality (Crenier et al., 2017) and habitat-limited aspects (Craig et al., 2015) may alter the proportional consumption of each nutrition pool. Exploring longer time scales over different climatic periods may improve the understanding of environmental controls on resource utilization. Resource utilization in contemporary ecology is commonly studied by consumer stable isotope (SI) signatures, which have been increasingly applied to chitinous invertebrate remains for disentangling various paleoecological questions . Paleoecological work on chironomid SI signatures has demonstrated timeand space-sensitive consumption of both autochthonous and allochthonous resources (Belle et al., 2017(Belle et al., , 2018Kivilä et al., 2019), as well as feeding on isotopically distinct methane-oxidizing bacteria (van Hardenbroek et al., , 2013). Yet, questions like what in the end controls the resource utilization, what might be its cascading effects for the food web functions and can we find such insight into aquatic carbon processing from paleolimnology, remain largely unanswered requiring further research.
In this study, we utilize functional paleoecology to explore SI signatures, carbon utilization and changes in functional feeding assemblages of benthic chironomid larvae in a subarctic lake for the past ca. 2000 years. We hypothesize that resource availability controls the long-term carbon utilization, but characteristics of the functional composition may alter it, either by habitat-mediated differences or by selective feeding behaviour. We also assume that the long-term environmental, climatic and the lake's ontogenic changes are expressed as changes in the ecosystem's functional composition. We aim to identify chironomid functional behaviour in relation to environmental perturbations, and decipher the character of ecological responses in relation to environmental change. We also investigate the relationship between functional assemblage structure and resource utilization, to assess the role of FFG structure for the type of nutrition consumed under varied conditions.

Study site, sampling and chronology
The study site lake Loažžejávri (69°53 0 N, 26°55 0 E) is located in Utsjoki, northern Finnish Lapland (Fig. 1). The shallow (max. depth 1.2 m), 3.0 ha lake is situated at an elevation of 255 m a.s.l. and features a ca. 2 km 2 catchment comprised of subarctic shrub tundra and small minerotrophic fen areas. At the sampling time, lake water was neutral, oligotrophic and transparent (Secchi depth [ 1.2 m = bottom) (Table 1). Hydrologically, the lake is connected to one small outlet stream and receives water from surrounding springs and a smaller spring fed lake. Bedrock in the area is composed of Proterozoic quartz feldspar paragneiss and amphibolite. The present-day climate on the region is subarctic, with mean annual temperature at -1.3°C, mean July temperature at 12.3°C and mean annual precipitation at 450 mm (Finnish Meteorological Institute).
A 38 cm sediment core was retrieved in late summer 2014 with a Limnos gravity corer (Kansanen et al., 1991) from mid-basin (1.2 m). Samples were sectioned at 1-cm intervals and kept refrigerated. In addition to sediment sampling, limnological data were A B Fig. 1 a Location of the study site Lake Loažžejávri in northern Finnish Lapland and b a detailed catchment map for the lake  (Rantala et al., 2016), which are reviewed in Table 1. The sediment core chronology was originally published in Luoto et al. (2017) and can be viewed in Supplementary file 1. The age-depth model is based on two AMS radiocarbon dates from plant macrofossils and was built with Clam 2.2 (Blaauw, 2010) in R (R Core Team, 2018) based on a linear interpolation function and the IntCal13 calibration curve (Reimer et al., 2013). 210 Pb activities were unfortunately too low to aid the age-depth model (see Luoto et al., 2017). As the chronological markers used for the age-depth model are few, the model should be interpreted with caution. Further on, we refer to the age model with a centennial resolution.

Chemical and biological analyses
In this study, the core was analysed for geochemical variables of sediment organic matter (OM), isotopic and functional assemblage composition of chironomid remains and the ratio of planktonic to littoral cladocerans (Crustacea: Cladocera). The amount of OM was determined with loss-on-ignition (LOI), where sediment samples were dried at 105°C overnight to determine water content, then burned in a muffle furnace for 2 h at 550°C for OM content, and consequently for 2 h at 950°C to determine the amount of carbonates (Dean, 1974;Heiri et al., 2001). Elemental (C org %, N tot %, C/N mass ratio) and isotopic (d 13 C OM , d 15 N OM ) composition of OM were analysed on freeze-dried sediment. Before analysis, carbonates (\ 4%) were removed with acid fumigation, where moistened subsamples (ca. 45 mg) were placed in a desiccator with a 100 ml beaker of 12 M HCl for 18 h and after fumigation dried in 60°C until dry (Ramnarine et al., 2011;Wotherspoon et al., 2015). Fumigated sediments were homogenized and transferred into tin capsules for elemental and SI analysis. The results are an average of duplicate runs and are corrected for the weight difference caused by fumigation process (Ramnarine et al., 2011;Wotherspoon et al., 2015). For SI analysis of chironomid head capsules (d 13 C HC , d 15 N HC ) standard procedures (van Hardenbroek et al., 2010;Heiri et al., 2012) were applied. Sediment was sieved (100 lm) and rinsed with ultrapure water, and head capsules were picked from the residue under a stereomicroscope (940). Head capsules (0.2-0.5 mg) were transferred to pre-weighed tin cups, dried in 40°C and weighed. Closed tin cups were stored in a desiccator until analysis. The sample material represents bulk assemblage composition of chironomids, due to scarcity of remains as well as high taxonomic diversity and evenness (see Luoto et al., 2017). SI analyses of both sediment OM and chironomid remains were performed at the University of Jyväskylä with a FlashEA 1112 elemental analyser coupled with a Thermo Finnigan DELTA plus Advantage mass spectrometer (Thermo Electron Corporation, Waltham, MA, USA) using internal standards (fish for chironomids, birch leaf for OM). The SI values are expressed as delta notations d = (R sample / R standard -1) 9 1000, where R = 13 C/ 12 C or 15 N/ 14 N and standardised to VPDB for carbon and atmospheric concentration for nitrogen. Discrimination factor was calculated as an offset between chironomid and OM isotope signature (e.g. Dd 13 C ¼ d 13 C HC À d 13 C OM ) for both carbon and nitrogen. These offsets are used as selectivity indices (Belle et al., 2017).
Results for the chironomid taxonomic assemblages in the sediment profile are published in Luoto et al. (2017). In this study, every taxon presented in Luoto et al. (2017) was assigned a FFG following Mandaville (2002) and Merritt et al. (2008). Relative proportions of FFGs were calculated from the full assemblage. For cladoceran remains, sediment subsamples were heated in 10% KOH, sieved (51 lm mesh) and centrifuged (4000 rpm) for 10 min (Szeroczyńska & Sarmaja-Korjonen, 2007;Rautio & Nevalainen, 2013). Samples were mounted on microscope slides, and for this study, the cladoceran remains were identified into planktonic or littoral forms and their relationship was calculated as planktonic to littoral (P/L) ratio (Alhonen, 1970;Sarmaja-Korjonen, 2001).
For temperature reference, tree-ring-based summer temperature (June-July-August, JJA) data with annual temporal resolution from Matskovsky and Helama (2014) was utilized. A temperature range for each sample was calculated as a 25-year average for the year indicated by our age-depth model. In addition, a chironomid-based summer temperature reconstruction from the same core  is presented for reference and to assist geochronological comparisons between the data and the independent (tree-ring) temperature reconstruction.

Statistical methods
Based on the SI measurements from the chironomid remains, a three-source Bayesian mixing model was built with simmr (Parnell & Inger, 2016) in R (R Core Team, 2018) to examine the relative proportions of different nutritional sources incorporated into the chironomid assemblage. The model with the same attributes has been previously applied on a surface sediment data set by Kivilä et al. (2019). The endmembers were based on published measurements from close-by region or similar environments (for further details see Kivilä et al., 2019, Appendix 2). We used measurements from bacterial mats and biofilm to estimate benthic signal (mean ± SD: d 13-C = -21.1 ± 3.7, d 15 N = 0.7 ± 1.3), zooplankton for the pelagic signal (d 13 C = -32.1 ± 1.6, d 15-N = 3.8 ± 1.3), and terrestrial leaves, litter and soil for the terrestrial signal (d 13 C = -28.3 ± 1.1, d 15-N = -1.7 ± 2.9). Methane-derived carbon, although known to substantially incorporate into chironomids (e.g. Jones & Grey, 2011), was not included in the model, as the d 13 C HC did not show values depleted enough in 13 C to indicate methane influence. The used literature-based end-member estimates allow for sensible spatio-temporal variation beneficial for a core study, excluding the problem of using site-specific one-time measurements as end-members. In addition, Bayesian models are relatively insensitive to variation in source values, as long as the sources are clearly separated from each other (Tanentzap et al., 2017). Chironomid-specific fractionation coefficients of 0.50 ± 0.56% and 1.5 ± 0.50% for carbon and nitrogen, respectively (Goedkoop et al., 2006), were applied for terrestrial and benthic sources. Pelagic end-member composing of zooplankton measurements did not require trophic correction. All observations fall within the area defined by the end-members in the isospace plot, supporting the suitability of the model (Phillips et al., 2014). The model was applied for each sample to observe variation within the core and then separately to three depth intervals where the modelled source contributions indicated differences.
Statistical relationships were examined with Spearman's rank correlations (q) and multivariate analyses. Prior to multivariate analysis the data were log 10 transformed and a constant was added to variables with negative values. Due to linear nature of the data, redundancy analysis (RDA) was performed to explore relationships of environmental variables with the chironomid functional assemblage, their isotopic signatures and modelled average source contributions. As explanatory parameters we used sediment geochemistry (LOI, C/N, d 13 C OM , d 15 N OM ) and the treering-based mean summer temperature (T JJA ), which all yielded low inflation factors (\ 10). Forward selection with 999 permutations was used for identification of significant connections. To identify major changes over time in the composition of functional assemblages, d 13 C HC and d 15 N HC , a segmented regression analysis was performed for defining breakpoints. Correlations were calculated with PAST 3.0 (Hammer et al., 2001), RDA were performed with Canoco 5 (Lepš & Š milauer, 2003) and breakpoint analysis was performed with SegReg software (Oosterbaan, 2005). Further editing of images was performed with CorelDraw.

Results
The results show that both sediment chemistry and biological variables have changed considerably during the investigated interval. In sediment chemistry, no major changes occurred before ca. 1600 CE, with OM quantity around 30%, C/N ranging between 10 and 11, d 15 N OM approximately 3% and d 13 C OM ca. -23% with slight variations (Fig. 2). At ca. 1600 C.E. proportion of OM (from 27 to 51%), C org % and N tot % started to rapidly increase. Concurrently, a clear enrichment in d 13 C OM and a depleting trend in d 15 N OM were observed, as well as a gradual lithological colour change from red-brown to bright green gyttja.
The chironomid functional assemblages featured five different FFGs: collector-gatherers (C-G), collector-filterers (C-F), shredders (SHR), scrapers (SCR) and predators (PRD) (Fig. 3). C-G and C-F were the most common groups throughout the core, their relative proportions ranging between 30-68 and 19-52%, respectively. At * 1400 C.E. a functional shift identified with the breakpoint analysis occurred for C-F and C-G. There C-F increased rapidly (breakpoint significance 99.9%) from an average of 26% (pre-1400 CE) to 41% (post 1400 CE) and concurrently C-G experienced a decrease (breakpoint significance 99.0%) that persisted to present day. SHRs featured a maximum abundance of 29% at the bottom of the core and gradually decreased towards the present day, while the trend in PRD abundances was fluctuating throughout the study interval ranging between 0 and 18%. Occurrence of SCR was sporadic, occurring in only three samples in very low abundances (\ 2%), and the guild composed of only one morphotype (Cricotopus bicinctus-type). Cladoceran P/L ratio showed mainly low values, but featured a distinct peak of planktonic fauna at 16-20 cm depth, with decreasing and increasing trends before and after the peak, respectively.
Carbon SI signatures of the chironomid remains (d 13 C HC ) ranged between -25.3 and -22.7% (average -24.1%) and d 15 N HC between 0.9 and 3.3% (average 1.7%) (Fig. 4). They followed crudely the same trend as the OM isotopes, d 13 C HC featuring first stable and then 13 C-enriching trend, while d 15 N HC depicted a mildly 15 N-enriching trend followed by a depleting trend. The chironomid-OM isotopic offsets (Dd 13 C, Dd 15 N) ranged between 0.8-2.3% and -2.4 to 3.3% for carbon and nitrogen, respectively. Dd 13 C fluctuated throughout the core but has been above average ([ 1.6%) since ca. 1800 CE, while Dd 15 N featured an overall decreasing trend.
The modelled source contributions of benthic, pelagic and terrestrial components for the chironomid diet throughout the core are shown in Fig. 4 as mean values, 50% probability range and 97.5% probability range. Benthic source dominated throughout the core with average contribution of 62 ± 3%, while proportions of pelagic and terrestrial components were considerably smaller, 31 ± 3% and 8 ± 3% on average, respectively. While the source contributions do not show large changes within the stratigraphy, three sections where the sources differ can be identified (Fig. 5). In the early core (before ca. 800 CE), the terrestrial component was slightly elevated compared to the rest of the core (probability 68%) (Fig. 5d). Between ca. 800 and 1600 C.E. the pelagic component showed elevated values (Fig. 5c), and highest pelagic contribution within the core occurred during this time with probability of 81%. Importance of benthic source has been increasing especially since ca. 1600 C.E. (probability 76%, Fig. 5b), and consequently the pelagic and terrestrial components have experienced a further decrease towards the present day. T JJA (tree-ring) T july (chironomid) Fig. 3 Relative percentage of functional feeding groups (Collector-gatherers: C-G, collector-filterers: C-F, shredders: SHR, predators: PRD and scrapers: SCR, last one with 10 times exaggeration line) and the cladoceran P/L ratio, indicative of water level fluctuations. In addition, chironomid-based mean July temperature reconstruction (T July , open circles) from the same core ) and a regional tree-ring based air temperature reconstruction (T JJA , black circles, from Matskovsky and Helama, 2014) are shown. Major climate phases and the functional shift are marked for reference as in Fig. 2  Statistical bivariate relationships between the geochemical variables and chironomid functionality, isotope signatures and resource utilization are shown in Table 2. The RDA (gradient length 0.8 SDU) revealed that environmental connections with the chironomid FFGs explained 24.2% (adjusted 12.4%) of the variance and the only significant variable picked out with forward selection was LOI, which contributed for 57.7% (P = 0.001) (Fig. 6). For the chironomid SI signatures explored with RDA, environmental connections explained 53.4% (adjusted 46.1%) of the variation and the most influential variables were d 13 C OM (64.4%, P = 0.001), d 15 N OM (20.6%, P = 0.002) and C/N (9.8%, P = 0.02). For source contributions, the explanatory variables accounted for 62.5% (adjusted 56.3%), whereas d 13 C OM (91.8%, P = 0.001) and C/N (5.6%, P = 0.05) were considered significant with forward section. In all of the RDA ordinations ( Fig. 6) uppermost 12 cm of the sediment core differed from the rest of the core as they were grouped separately in the ordination spaces.

Discussion
Mixed organic matter sources in the early core contributions within the whole core material, and contribution probabilities for b benthic, c pelagic and d terrestrial sources during three differing intervals in the core century in the Arctic), which in northern Scandinavia was characterized by gradually cooling conditions and southward retreat of the tree line (Kaufman et al., 2009;Luoto et al., 2017). Hydroclimate at the time is not fully resolved, but Linderholm et al. (2018) suggest that it may have been a drier period in northern Scandinavia. In our sediments, this was reflected as a geochemically fairly stable period with OM content of ca. 28%, and OM isotopic signature (d 13 C OM ca. -22%, d 15 N OM ca. 3%) suggesting a combination of allochthonous and autochthonous input. A decreasing trend in d 13 C OM and an increasing trend in C/N, which persists until ca. 1000 CE, may indicate relative increase of terrestrial input into the lake system, or a respective decrease in autochthonous production. In addition, the relative lake level was gradually declining during the time, as suggested by the cladoceran P/L ratio, which gives support to limited pelagic production and may indicate drier climate. The cladoceran signal may relate to actual water table fluctuation or alternatively to gradual filling in of the lake basin, but regardless is indicative of relative reduction in pelagic niche space. Considering the modelled resource utilization, benthic source was dominant in the early parts of the sediment profile, being so throughout the whole core (Figs. 4, 5a). The terrestrial nutrition component was also slightly elevated in the lower part of the core compared to overall average and modern conditions (Fig. 5d), which supports the presumed larger relative allochthonous influence. As input of terrestrial carbon into Arctic lakes is typically enhanced by increased temperature and precipitation (Tranvik et al., 2009), neither of which occurred at the time interval, we presume that suppression of the in-lake production leads to a wider utilization of terrestrial nutrition components. It has been suggested that terrestrial components are more likely to support food webs when high-quality food components are limited (Solomon et al., 2011;Tanentzap et al., 2017). However, terrestrial input was rather low throughout the studied interval, and therefore, terrestrial sources do not seem to be of high importance in supporting the food web.
The distribution between FFGs was most even and diversity highest at the bottom of the core (Fig. 3). Occurrence of SCR and highest abundance of SHR may indicate suitable habitats differing from present day, and the overall conditions must have been more favourable for diverse functional communities. SHR are generally associated with coarse particulate organic matter (CPOM), meaning living or decaying vascular plants, macrophytes and woody debris, while SCR utilize mostly periphyton as their food source (Merritt et al., 2008). This may indicate more pervasive macrophyte flora in the shallower lake, as many taxa belonging to SHR or SCR groups are generally associated with macrophyte habitats (Brooks et al., 2007;Merritt et al., 2008). Due to the relatively omnivorous character of chironomids, availability of specific food items may not be a major controller of presence and abundance of the different FFGs. Instead, habitat-related factors, such as water depth or DOC and nutrient concentrations, have been formerly shown to guide occurrence of FFGs (Luoto & Ojala, 2014;Kivilä et al., 2019).

Increasing pelagic influence
During the Medieval Climate Anomaly (MCA, ca. 900-1150 CE), which in northern Scandinavia was characterized by warm and dry summers (Linderholm et al., 2018), the T JJA showed clearly elevated temperatures for the region especially between ca. 600-1000 CE. In the sediment geochemistry, this period is characterized by more depleted d 13 C OM and increasing C/N until ca. 1000 C.E. (Fig. 2), which is most likely indicative of combined effects of increased terrestrial input and increased planktonic production.  Terrestrial input increases during warm periods due to increased bio-activity of catchment soils (Tranvik et al., 2009), contributing to both of the biogeochemical signals. Further on, planktonic production in northern lakes is enhanced during warm periods, benefitting from warmer water temperatures, longer open water period and increased catchment nutrient inputs (Rautio et al., 2011;Seekell et al., 2015). After ca. 1000 CE, as T JJA featured increased variability and a cooling trend towards the Little Ice Age (LIA), C/N began to decrease and d 13 C OM experienced a slight enrichment, which supports an increase in the relative proportion of autochthonous production. The signal may also be indicative of decreasing terrestrial input instead of actual increase of in-lake production during the rather cool climate conditions. However, despite the shallowness of the lake and rather harsh climatic conditions, the sedimentary evidence shows that the lake has retained enough pelagic primary production to host a pelagic food web throughout its studied history .
Our resource modelling shows a similar pattern: the pelagic component rose slightly soon after the temperature increased and featured elevated contribution for the rest of the MCA, extending some time to the LIA (ca. 800-1600 CE), while concurrently the benthic component experienced its lowest contribution (Figs. 3, 5b, c). It has been suggested that autochthonous, high-quality algal nutrition sources are preferred by at least collector-gathering chironomid Sergentia coracina when available (Belle et al., 2017), and that even small amounts of high-quality food significantly improves the growth of detritivores in laboratory (Crenier et al., 2017). If there is a true preference for high-quality nutrition, it could be expected that the relative availability of such would be reflected into the resource utilization, even if already small amounts contribute to improved growth and survival of the larvae.
During the MCA, given our chronology is adequate, the C-F group featured low abundance (Fig. 3), in support with former work showing that C-F group strongly prefers cold and oligotrophic conditions (Luoto and Nevalainen, 2015). The dominant FFG at the time was C-G, which appears as the most cosmopolitan of the FFGs; they are common and have a wide tolerance for many environmental variables, for instance higher DOC and nutrients (Luoto & Ojala, 2014, 2018Kivilä et al., 2019). Also PRD, which favour warmer conditions (Luoto & Nevalainen, 2015), showed an abundance peak. The distinct FFG composition during the MCA is also discernible in the RDA, where samples from 29 to 25 cm depth, corresponding to the warmer T JJA temperatures, form a separate cluster (Fig. 6a). Accordingly, in addition to dominance of C-G, T JJA is associated with the MCA cluster, which supports that climatic factors like temperature may have been important in guiding the FFGs (Luoto and Nevalainen, 2015). In addition, the same samples are connected to higher C/N in the isotope and source RDAs, and are aligned with the pelagic source (Fig. 6b, c), further supporting the effects of temperature-mediated, increased pelagic production supporting the food web.
According to the age-depth model, the MCA is located in the lower half of the core where the chronological control is weaker. In addition, the warmer temperatures are not pronounced in the within-core chironomid-based temperature record at the same time as the T JJA suggest a warm period. This may be due to the spatial heterogeneity of the MCA or differences between the terrestrial and aquatic temperature responses, as discussed in Luoto et al. (2017), but while some chronological mismatch cannot be completely ruled out, one of our age markers at 25 cm (997 ± 30 cal BP) suggests that timing of the MCA should not be far off.

Functional changes and boosting benthic production
The most drastic changes in our data occurred during the Little Ice Age (LIA) period. The LIA (ca. 1400-1850 CE) was a spatially and temporally heterogeneous cold period, which, according to most paleoenvironmental reconstructions, is considered to feature dry winters and wet summers in northern Scandinavia, especially at its early phase (Linderholm et al., 2018). A major shift in the composition of functional assemblages occurred at ca. 1400 CE, most explicitly expressed as the sudden increase of C-F, accompanied by a relative decrease in C-G (Fig. 3). Concurrently occurring rapid increase in relative water depth inferred from cladoceran P/L ratio suggests hydrological changes in the lake system, for instance opening of an inlet stream or a temporary drainage channel. Such hydrological changes may well have occurred during the hydroclimatic fluctuations of the early LIA. Accordingly, northern Finland has been identified as an area where precipitation was above the European average during this time (Cook et al., 2015). The relative abundance of C-F remained high since the rapid increase until present, despite the cladoceran P/L ratio indicating that the high water table was a temporary event of a little over a 100 years. However, as several external, habitat-related factors guide the occurrence of the FFGs , for instance the prevailing cool temperatures (T JJA ) until ca. 1900 C.E. may partly explain the C-F dominance, given their colder temperature optimum (Luoto and Nevalainen, 2015), while other conditions, such as more oxygenated water , may also influence. Overall, it seems likely that the environmental changes (e.g. water level fluctuations) triggered habitat change or diversification in the system, which better accommodated the C-F group.
Major changes in the sediment geochemistry occurred at ca. 1600 CE, approximately two centuries later than the functional shift (Fig. 2). Rapid increase in organic content (LOI, % C org ) at the same time likely indicated boosting in-lake production, and a concurrent 13 C enrichment of d 13 C OM (from ca. -23 to ca. -21%) showed a strong benthic signal (Hecky & Hesslein, 1995). These changes occur after the coldest peak of the LIA in the record, and may thus relate to the slowly increasing temperatures, or other more favourable habitat conditions following the environmental perturbations of the early LIA. Increasing productivity, associated with natural and anthropogenic post-LIA warming, is a common feature to many lakes in the Arctic zone . Decreasing trend in d 15 N OM , originating already at 1400 CE, suggests increasing fixation of atmospheric nitrogen into the biomass and may thus be a signal of cyanobacterial production. Cyanobacteria are common members of autotrophic benthic communities in shallow northern lakes (Rautio et al., 2011), becoming increasingly abundant since ca. 1800 C.E. (Taranu et al., 2015). The deposition of reactive nitrogen may also contribute to the near-surface-negative trend (Holtgrieve et al., 2011). Post-depositional alteration of organic matter may play some role in the biochemistry of the uppermost sediments, especially if these sediments are younger than implied by the age-depth model. This may cause slight overestimation of the OM content due to fresh and non-degraded material and contribute to the lower C/N ratio, however, the post-depositional alteration of isotope signals likely is small (Meyers, 1994).
The modelled source contributions seemed to follow changes in the sediment geochemistry rather than those in the FFG composition, as the major shift in resource utilization occurred at ca. 1600 C.E. along with the sediment SI signatures. This implies that the availability of nutrition has more control over the chironomid resource utilization than the food acquiring method (FFG). This may not be surprising given the generally omnivorous feeding habits of chironomids, alike to most aquatic insects (Marcarelli et al., 2011), which allows them to utilize various types of nutrition (Serra et al., 2016). External controls may thus be most important in guiding composition of nutrition consumed . As the geochemical changes at ca. 1600 C.E. seem to link to changes in amount and composition of benthic production, it further suggests that benthic primary production and macroinvertebrate functionality are sensitive to different environmental conditions.
In all of the RDA ordinations, the top part of the sediment formed a separate cluster from the rest of the core along RD1 axis values (Fig. 6). In the isotope and source RDAs (Fig. 6b, c) this corresponded to the top 12 cm (since ca. 1650 CE), while only the top 5 or 6 cm (since ca. 1800 CE) were separately clustered in the RDA for FFGs (Fig. 6a). The cluster is aligned with LOI and d 13 C OM , suggesting a positive association with increased benthic primary productivity. In addition, the uppermost 12 cm correspond roughly to the increasing temperatures after the coldest LIA, and the top 5 cm roughly to the recent rapid increase in temperature according to T JJA (Fig. 3), which is in line with the increasing productivity. Temperature-mediated changes (e.g. reduced ice-cover) have been formerly identified as main drivers of ecological change for primary producers in the Arctic (Griffiths et al., 2017) and in Finnish Lapland (Sorvari et al., 2002). This supports the climate-mediated effects on productivity and consequent higher organic conditions, even though T JJA showed no particular gradient, nor association with the top 12 cm cluster in our RDAs (Fig. 6). In addition to temperature, deposition of atmospheric nitrogen may have contributed to the increasing productivity during the industrial era (Bergström & Jansson, 2006). Overall, the RDAs confirm that the younger and older parts of the sediment core differ from each other considerably. In addition, it is notable that the recent samples and those representing the MCA were located in nearly opposite ends of the RD1 axis (Fig. 6). The difference was clearest with the FFGs (Fig. 6a), where the more recent (ca. 1800 C.E. onward) samples were associated with LOI whereas the cluster of MCA samples was associated with T JJA . This implies that the benthic functional assemblage was different by composition and environmental controls between the MCA and recent warming, suggesting that the MCA may not be a good analogue for functional changes caused by the on-going climate warming. The reason for this discrepancy remains uncertain, however, we suggest it is related to changes in productivity and habitat differences. Compared to the MCA the modern increasing benthic production and microhabitat diversification, resulting from climate-mediated changes (e.g. reduced ice-cover) and possibly amplified by atmospheric N-deposition, reflect different proportional availability of nutrition sources, which is further affecting the isotope signatures. The FFGs on the other hand tend to be sensitive to several habitat-related factors (Luoto & Ojala, 2014;Kivilä et al., 2019), and their occurrence may be related to for instance differences in substrate quality.
In this upper part of the core, the chronology can be considered more reliable; however, there might be some temporal mismatch in the more recent sediments, due to lesser compaction. The functional change horizon at 21 cm is temporally well constrained (age marker at 20 cm, 1415 ± 25 cal BP), in addition, the two temperature reconstructions agree fairly well with the onset and duration of the LIA (Fig. 3)  . The sediment at 21-20 cm depth was particularly rich in mossy plant macrofossils, and it is assumed that the rejected dates also originate from this layer.

Implications for feeding behaviour and ecosystem functions
Generally, the observed changes in stable isotope composition of chironomids were relatively small, reflecting the availability of food items with a rather narrow variability in isotopic signatures throughout the core. To some extent, homogeneity of isotope measurements might be due to us using bulk assemblage measurements for the chironomid SIs and thus averaging out taxon-specific feeding. For example Frossard et al. (2014) noted different d 13 C values between taxa. However, they observed considerable differences between different time periods, and particularly between the littoral and deep zone of the lake, whereas the variability of actual d 13 C values of different littoral taxa from the same time window was small. In addition, they observed rather similar temporal patterns between taxa in the littoral zone. This supports that the resource availability (changing between time and habitats) is a larger control over chironomid d 13 C than individual feeding habits, as also noted by Kivilä et al. (2019). Most of the taxa in our sequence compose of detritus feeders (C-G and C-F) according to Merritt et al. (2008), and therefore, large isotopic differences in nutrition between the taxa are unlikely in the shallow study lake where the sources of organic matter seem to remain rather stable. In the end, the bulk assemblage measurements provide a well-homogenized community level measurement for the carbon utilization of chironomids in this setting, however, the discrepancy in resolution when comparing with FFG or taxonomic composition should be noted.
We hypothesised that characteristics of the functional assemblage such as selective feeding behaviour may alter the chironomid resource utilization. Looking at the selectivity level of chironomids (Dd 13 C), the offsets are elevated at the bottom (ca. 100-500 CE) of the core and again at the top (since ca. 1800 CE) (Fig. 4), while the FFG assemblage differs between these two time intervals. At the bottom, specialized feeding groups (SHR, PRD and SCR) were more common, which may explain the larger offset from detritus. At the top, however, the assemblage is clearly dominated by collectors (C-F and C-G), and it seems likely that the offset is created by selective feeding behaviour. The lowering of sediment C/N occurred at ca. 1800 C.E. concurrently with the beginning of increasing trend in Dd 13 C, indicating that the chironomids preferred to assimilate something other than the bulk of benthic microbiota. The C/N change may indicate increase of bacterial activity in the benthos, as their stoichiometry commonly reflects a lower C/N ratio (Goldman et al., 1987). This is in accordance with Kivilä et al. (2019), who suggested that the selectivity may be higher in environments where many food components of differing nutritional value are easily available, such as benthic mats. It has also been noted that chironomids select for higher quality nutrition in excess of poor quality food (Belle et al., 2017), which suggest that the lower quality components (e.g. terrestrial carbon) are incorporated into the food web most when the higher quality components are not present in large amounts. This is understandable given the requirements for growth and reproduction that require high-quality nutrition, for example essential fatty acids (Brett and Müller-Navarra, 1997;Crenier et al., 2017). It seems that the selectivity (as Dd 13 C) is elevated when high-quality nutrition is not a major component in the bulk OM, regardless if the habitat was productive or unproductive. However, the d 13 C HC follows the trends of the dominating component of OM, as the generally omnivorously feeding chironomids utilize less valuable nutrition sources for maintaining their life functions. Overall, the chironomid isotope signatures appear to be more directed by sediment geochemistry (food source), and to a much lesser extent mediated by selectivity implied by different FFGs.
Resilience of ecosystems is a key question in assessment of future responses to the on-going climate change. Generally, many factors improve the system resilience either by being resistant to change or allowing a quick recovery after a perturbation (Peterson et al., 1998). The shift in functional composition in our data is not evident by looking at the taxonomy it is based on (see Luoto et al., 2017), which highlights that the stability characteristics of species composition are not always relevant for the resilience of ecosystem functions (Oliver et al., 2015). The environmental perturbations experienced in the early LIA (ca. 1400 CE) seemed to be enough to push the FFG assemblage into an alternative state, which has persisted similar through the changes that caused the primary production to rearrange (ca. 1600 CE), as well as under the current warming. In support of our hypothesis, that the functional assemblages would reflect changing environmental conditions and perturbations, our data suggests that examining the FFG assemblages from down-core sites provides different ecological information than the taxonomic data and may improve our understanding of past ecosystem dynamics in relation to environmental perturbations.

Conclusions
Our study featured ca. 2000 years of variability in geochemistry and functional assemblages of chironomids, where major reorganization horizons occurred at different times. This suggests that changes in functional composition of chironomids and sediment geochemistry, at most reflecting changes in autochthonous production, are guided by different external variables. Resource modelling based on chironomid stable isotopes revealed that the utilization of different carbon pools followed primarily sediment geochemistry, and was thus controlled more by resource availability than the functional composition of the assemblage. The effects of Late-Holocene climate fluctuations were observed for functional assemblages, and interestingly present and past warm periods shared little similarity. As taxonomic and functional assemblage dynamics differed throughout the core, it seems that functional paleoecology may improve our knowledge of past ecosystem dynamics, however, more research is required to unravel the full potential of functional approaches in down-core studies.