Responses of tundra plant community carbon flux to experimental warming, dominant species removal and elevation

1. Rising temperatures can influence ecosystem processes both directly and indi -rectly, through effects on plant species and communities. An improved under standing of direct versus indirect effects of warming on ecosystem processes is needed for robust predictions of the impacts of climate change on terrestrial eco system carbon (C) dynamics. 2. To explore potential direct and indirect effects of warming on C dynamics in arctic tundra heath, we established a warming (open top chambers) and dominant plant species ( Empetrum hermaphroditum Hagerup) removal experiment at a high and low elevation site. We measured the individual and interactive effects of warming, dominant species removal and elevation on plant species cover, the normalized difference vegetation index (NDVI), leaf area index (LAI), temperature, soil mois ture and instantaneous net ecosystem CO 2 exchange. 3. We hypothesized that ecosystems would be stronger CO 2 sinks at the low elevation site, and that warming and species removal would weaken the CO 2 sink because warming should increase ecosystem respiration (ER) and species removal should reduce gross primary productivity (GPP). Furthermore, we hypothesized that warm ing and species removal would have the greatest impact on processes at the high el evation where site temperature should be most limiting and dominant species may buffer the overall community to environmental stress more compared to the low elevation site where plants are more likely to compete with the dominant species.


| INTRODUC TI ON
Global environmental changes are shaping ecosystem carbon (C) feedbacks to the atmosphere, especially in high latitude terrestrial ecosystems that currently serve as important C sinks (Tarnocai et al., 2009) but that are undergoing rapid warming (IPCC, 2013). The rate of CO 2 uptake by plants as represented by gross primary productivity (GPP) and ecosystem CO 2 loss through ecosystem respiration (ER; both autotrophic and heterotrophic) are the key fluxes that drive terrestrial net ecosystem exchange of CO 2 (NEE) between the biosphere and the atmosphere (Luo, 2007). While GPP and ER can be temperature sensitive, responses of NEE to rising temperatures are challenging to predict and often vary among locations, plant communities and the duration of warming treatments (e.g. Cahoon, Sullivan, & Post, 2016;Dorrepaal et al., 2009;Leffler, Klein, Oberbauer, & Welker, 2016;Sharp, Sullivan, Steltzer, Csank, & Welker, 2013). Due to this complexity, improved understanding of the direct effects of temperature change versus effects of other environmental factors on C flux would help to refine predictions of C dynamics under global warming.
Plant growth is more constrained by low temperatures than is photosynthesis (e.g. Körner, 2013), and responses of GPP to warmer temperatures are thus often governed by the availability of other plant growth limiting factors such as soil moisture and nutrient availability (Hobbie & Chapin, 1998). As a result, when weak or no responses of GPP and NEE to experimental warming occur, factors other than temperature may be limiting GPP. When positive responses of GPP occur, they are often associated with increased leaf area, plant growth and/or changes in plant community composition (Hobbie & Chapin, 1998;Leffler et al., 2016;Sharp et al., 2013). With regard to ER, positive responses of ER to experimental warming in tundra (Hobbie & Chapin, 1998;Luo, 2007) can occasionally result in ecosystems becoming net CO 2 sources (Welker, Brown, & Fahnestock, 1999). However, a number of factors influence how ER responds to warming, including plant community properties and processes and nutrient and water availability (Cahoon, Sullivan, Shaver, Welker, & Post, 2012;Luo, 2007). Furthermore, short-and long-term experimental warming can increase the contribution of autotrophic respiration to ER (Hicks Pries et al., 2015), while the magnitude of soil respiration response (i.e. heterotrophic and plant roots) to warming often declines over time (Melillo et al., 2017;Rustad et al., 2001). Such a decline in the sensitivity of soil respiration to warming may result from soil microbial acclimation to warmer temperatures, especially when temperatures are already warm (Luo, Wan, Hui, & Wallace, 2001;Melillo et al., 2017). Additionally, the decomposition rate of soil organic matter can have a higher temperature sensitivity when the organic matter is more recalcitrant (Davidson & Janssens, 2006) and when temperatures are lower (Kirschbaum, 1995).

Positive interactions among plants often increase with increas-
ing stress and harsher climatic conditions (Callaway et al., 2002). Some experiments in the Arctic found no, or negative effects of species removals on the growth of neighbouring species (Chapin, McGraw, & Shaver, 1989;Jonasson, 1992), even under experimental warming (Shevtsova, Haukioja, & Ojala, 1997). Under the climatic conditions that characterize high elevations and latitudes, removal of a dominant species may result in net reductions in the growth of the remaining plants, and potentially of GPP and autotrophic respiration. In turn, such reductions would result in no net shift in NEE if heterotrophic respiration is similarly reduced. In partial support of this, Nielsen et al. (2017) found that NEE was unresponsive to shrub removal and its interaction with experimental summer warming in a Greenland fen. Pairing plant removal (Díaz, Symstad, Chapin, Wardle, & Huenneke, 2003) and global change experiments can thus provide important information on how global 4. The instantaneous CO 2 flux, which reflected a weak CO 2 sink, was similar at both elevations. Neither experimental warming nor dominant species removal significantly changed GPP or instantaneous net ecosystem CO 2 exchange even though species removal significantly reduced ER, NDVI and LAI. 5. Our results show that even the loss of dominant plant species may not result in significant landscape-scale responses of net ecosystem CO 2 exchange to warming. They also show that NDVI and LAI may be limited in their ability to predict changes in GPP in these tundra heaths systems. Our study highlights the need for more detailed vegetation analyses and ground-truthed measurements in order to accurately predict direct and indirect impacts of climatic change on ecosystem C dynamics.

K E Y W O R D S
carbon, ecosystem respiration, global warming, gross primary productivity, leaf area index, normalized difference vegetation index, plant-plant interactions change factors may mediate plant-plant interactions to influence high latitude terrestrial ecosystem processes (Aerts, 2010;Nielsen et al., 2017).
Arctic landscapes are spatially heterogeneous and contain pronounced variation in climatic conditions, such as those that occur with changes in elevation. Consequently, plant communities that occur at contrasting elevational sites have experienced long-term differences in climatic conditions (Körner, 2007), and variation in elevation is often associated with changes in plant biomass, plantplant interactions and community composition (e.g. Callaway et al., 2002;Sundqvist, Sanders, & Wardle, 2013). Experiments that disentangle the role of direct versus indirect (e.g. plant community change) effects of temperature on C flux for communities at high and low elevation sites can therefore provide information for C model parametrization about long-and short-term controls over ecosystem C flux dynamics (Ostle et al., 2009;Saleska et al., 2002).
We used a dominant species removal experiment (Díaz et al., 2003) coupled with experimental warming (OTCs; Dorrepaal et al., 2009) at a high and low elevation site in a subarctic tundra in northern Sweden to explore the direct and indirect effects of warming on C flux in an arctic landscape. Specifically, our aim was to quantify standardized instantaneous CO 2 flux responses across the experimental treatments at these sites (Metcalfe & Olofsson, 2015;Wardle, Jonsson, Mayor, & Metcalfe, 2016) and compare treatment responses between elevational sites. Hence, to ensure comparability among treatments and measurements, we conducted all measurements at the period of maximum plant biomass. We related instantaneous NEE standardized at the same light level and temperature, to plant and abiotic properties measured over the same period. Our study is thus intended to provide standardized, comparable data on the role of experimental warming and dominant plant species removal, and their interaction, in regulating ecosystem CO 2 flux in a tundra heath landscape. This design allowed us to test three inter-related hypotheses: 1. Across treatments, higher temperatures, more plant biomass and higher leaf area at the low elevation site relative to the high elevation site will result in a greater net CO 2 sink at the lower elevation.
2. Both GPP and ER will be stimulated by short-term experimental warming, but ER will be more stimulated than GPP and dominant plant species removal will reduce GPP more than ER. Hence, both experimental warming and dominant plant species removal will result in a weaker net CO 2 sink.
3. Decreases in the net CO 2 sink resulting from warming and dominant species removal should be greater at high elevations than at low elevations. We predict this pattern because dominant plant species removal should be most detrimental for neighbouring plant cover at the high elevation site. Furthermore, the response of ER to experimental warming should be greatest at the high elevation due to lower acclimatization of soil respiration and a higher temperature sensitivity of soil organic matter decomposition.
By testing these hypotheses, we aim to advance the understanding of direct and indirect (via plant species and community responses) regulatory effects of temperature on C flux.

| Study site and experimental set-up
We established this experiment at a low (500 m) and a high ( The bedrock consists of salic igneous rocks and quartic and phyllitic hard schists. All of the study plots (n = 40) at both elevations were placed in heath vegetation, which is a common circumpolar vegetation type (Tybirk et al., 2000). Prior to experimental manipulations, all plots were dominated by Empetrum hermaphroditum Hagerup.
(M ± SE of ground cover = 66.2 ± 1.8%; n = 40), with Betula nana L. Warming, Warming × Removal with the restriction that plots with different treatments were placed ≥3 m apart, and plots with the same treatment were >10 m apart, rendering five blocks at each elevation. While the distance among plots should ensure spatial independence for soil microbial processes among plots (Baldrian, 2014), we further constructed semivariograms for our CO 2 flux data from each elevational site, which verified that fluxes across our plots were not spatially autocorrelated ( Figure S1). For plots designated for removal, the above-ground biomass of the dominant species, E. hermaphroditum (as defined by the species with the greatest per cent cover within each elevation), was physically removed over 18-25 July 2014 across the entire 2 × 2 m area of each plot. Removal treatments were then maintained in late June in subsequent years. The removed biomass was dried at 48°C until constant mass (≥72 hr) and weighed. Removal resulted in roughly 3.3 times more biomass removed at the low than at the high elevation plots; the average total E. hermaphroditum biomass removed (M ± SE; n = 10) was 719.4 ± 37.8 g/m 2 per plot for the low elevation plots, and 218.7 ± 14.0 g/m 2 for the high elevation plots.

| Air temperature, soil temperature, relative air humidity and soil moisture
Between July 1 and August 31, 2016, we measured air and soil temperature and relative air humidity hourly in all plots at 5 cm aboveground (air temperature and humidity) and 5 cm below-ground (soil temperature) at the centre of each plot (Thermochron & Hygrochron Ibuttons, Maxim Integrated Corp.). These measurements were used to calculate mean daily and mean monthly values for each variable. Three measurements of volumetric soil moisture were made at 12 cm soil depth (Campbell HS2 soil water probe, Campbell Scientific) in each plot on each of three occasions during the growing season (22-21 June, 13-14 July, 16-17 August) including at the same time of day as CO 2 flux measurements were made. For each time of measurement, the three measurements per plot were averaged to retain a single mean value of soil moisture per plot. to the ground with plastic skirts and chains (Metcalfe & Olofsson, 2015). We recorded the change in CO 2 (in p.p.m.) within the chamber over a 2-min period with an infra-red gas analyser (LiCor 7500, LICOR Biosciences). Electric fans were used to mix the air inside the chamber during CO 2 measurements. We measured photosynthetically active radiation (PAR, μmol m −2 s −1 ) during each measurement with a light sensor placed above the vegetation orientated directly upwards at ~0.55 m height inside the chamber. Four measurements of ecosystem CO 2 flux, each at a different light level, were taken for each plot using mesh covers and an opaque cover to reduce the light inside the chamber; one measurement under full ambient light, two measurements under different levels of reduced light and one in complete darkness. We calculated NEE for each light level using the following equation (e.g. :

| CO 2 flux measurements
where ρ is the air density (mol/m 3 ), V is the chamber volume (m 3 ), dC/dt is the slope of chamber CO 2 concentration against time (μmol mol −1 s −1 ) and A is the ground surface area enclosed by the chamber (m 2 ).
These measurements represent instantaneous NEE of CO 2 from the ground surface (μmol CO 2 m −2 s −1 ). For each measurement in light, NEE = GPP − ER, where a positive number for NEE represents a net uptake of CO 2 from the atmosphere to the ecosystem at the plot level.
For each measurement taken in darkness GPP = 0 and NEE = ER.
To allow direct comparison of light-dependent CO 2 fluxes across our treatments, we standardized our measurements of NEE and GPP across all plots to a single PAR value (600 μmol m −2 s −1 ). To do so we fitted a rectangular hyperbola to the measured relationship between PAR levels and the measured corresponding NEE values  for each plot, and freely varied A max and k until the root mean square error between predicted and observed data was minimized, using Excel solver: where A max is the light saturated rate of GPP (μmol CO 2 m −2 s −1 ) and k is the half-saturation light constant (μmol m −2 s −1 ). There was a strong linear and proportional relationship between resultant predictions of NEE based on fitted values of A max and k and observed NEE from the study site ( Figure S2). Fitted values of A max and k for each plot were then used to estimate plot-specific NEE at 600 μmol m −2 s −1 (NEE 600 ).
GPP at 600 μmol m −2 s −1 (GPP 600 ) for each plot was then calculated from the formula: with a positive value of NEE 600 indicating a net uptake of CO 2 from the atmosphere to the ecosystem at the plot level at PAR levels of 600 μmol m −2 s −1 , and a negative value of NEE 600 indicating a net release of CO 2 from the ecosystem to the atmosphere. In addition, for each plot we standardized ER at a temperature of 10°C (ER 10 ) by using the relationship between temperature derived from the LiCor 7500 and the NEE measurements made under full darkness across plots. For each plot, we used these ER 10 values to calculate NEE 600at10 by correcting NEE 600 for the difference in ER and ER 10 . As the relationship between temperature and NEE measurements made in darkness was statistically significant (p = 0.023) but had a low R 2 of 0.132 (df = 1,38), we report fluxes both as unstandardized and standardized for temperature at 10°C.

| Plant community measurements
The cover of each plant species was determined in each plot be- Instruments, Llandrindod Wells). We used this NDVI data to calculate Leaf Area Index (LAI, m 2 leaf m −2 ground) for each plot using a relationship previously developed for E. hermaphroditum heath in this region : These measurements were used to explore the effect of our treatments on NDVI and LAI. To further explore differences in GPP 600 due to variation in NDVI and LAI across our plots, we divided GPP 600 for each plot by the corresponding plot measurement of NDVI measurements and by LAI calculated from our NDVI data for August  to derive GPP 600 per unit NDVI (GPP 600NDVI ) and per unit leaf area GPP 600LAI .

| Statistical analyses
To explore how temperature (air and soil), air humidity, soil moisture, NDVI and LAI responded to our treatments across the two elevations and among the 3 months of the growing season (June, July and August), we used linear mixed-effects models (LMMs) with summer month, elevation, Removal and Warming as fixed factors, and block and plot as random factors. We further used LMMs to test for effects of time since experimental establishment (i.e. pre-establishment 2014 vs. 2016), elevation, Removal and Warming on the plant cover data for all vascular species that had an average cover of ≥3% in at least one treatment for each elevation. LMMs were also used to test for the main and interactive effects of elevation, Removal and Warming on NEE 600 , GPP 600 , ER, NEE 600at10 , ER 10 , GPP 600NDVI and GPP 600LAI . Furthermore, we calculated the unstandardized mean difference (D) ± 95% confidence interval between treatment means and control plots. Within each elevation, we used Spearman's rank correlation to examine monotonic relationships between NEE 600 , GPP 600 , ER, NEE 600at10 , ER 10 , and the following variables: NDVI measured from 10 to 13 August 2016, per cent cover measured from 2 to 13 August 2016 for each plant species with an average cover of ≥3% in at least one treatment, air temperature, air humidity and soil temperature averaged from 10:00 a.m. to 5:00 p.m. on the day that CO 2 -flux measurements were taken, and soil moisture from the day that

| Temperature, relative air humidity and soil moisture
In 2016, mean monthly summer air and soil temperatures were highest at the low elevation site compared to the high elevation, and the Warming treatment significantly increased air temperatures (Appendix S1, Tables S1 and S2). Warming and Removal influenced mean monthly summer soil temperature but these effects varied by elevation and among treatments, resulting in no differences among treatments at the high site while soil temperatures were higher in Warming × Removal plots compared to Control plots at the low site in July (Appendix S1, Tables S1 and S2).

Relative air humidity varied most among low elevation plots
where it was reduced in Warming × Removal plots and Removal plots, but not in Warming plots, relative to the Control plots (Tables S1 and S2). In contrast, at the high site relative air humidity was only reduced in Warming plots relative to the Control plots (Tables S1 and S2). Soil moisture varied over time between 31 ± 5% and 25 ± 1% (M ± SE; n = 5; Control plots), and decreased from June to August, at the high elevation site. At the low elevation site, soil moisture varied over time between 18 ± 2% and 12 ± 1% (M ± SE; n = 5) and was lowest in July. Removal consistently reduced soil moisture at the high elevation site but not at the low elevation site (Tables S1 and S2, significant Removal × elevation interaction).

| Plant properties
Removal had a significant and positive effect on V. vitis-idaea cover, which increased most in Removal and Warming × Removal plots (Tables S3 and S4). Furthermore, B. nana cover increased in treatment plots compared to Control plots with the highest increase in Warming × Removal plots (Tables S3 and S4). There was a positive effect of warming on the cover of V. uliginosum which increased in LAI = 0.0259 × e (5.087×NDVI) .
Warming × Removal plots at both elevations, and in Warming plots at the high elevation, relative to Control plots (Tables S3 and S4).

Removal reduced both NDVI and LAI and the effect of Removal on
LAI varied by month (Figure 1; Tables S5 and S6; Figure S3). Warming slightly reduced LAI only at the high but not the low elevation site ( Figure S3; Tables S5 and S6) Table S6 for treatment effects within months).

| CO 2 flux
The plots overall represented a net sink for CO 2 at the time of our measurements. There were no main or interactive effects of elevation, Warming, or Removal on NEE 600 or NEE 600at10 but a significant interactive effect of Warming and Removal on NEE 600at10 (Figure 2; Figure S4; Table S7); NEE 600at10 was higher in Warming plots and Removal plots, and lower in Warming × Removal plots, than in Control plots (Table 1). Removal reduced ER by 0.85 ± 0.62 (μmol CO 2 m −2 s −1 ; D ± 95% CI; n 1 = n 2 = 20) and ER 10 by 0.55 ± 0.55 (μmol CO 2 m −2 s −1 ; D ± 95% CI; n 1 = n 2 = 20; Figure 2; Figure S4; Table S7). Furthermore, there were significant interactive effects F I G U R E 1 Normalized difference vegetation index (NDVI) values in response to warming and removals at high (a-c) and low (d-f) elevation in arctic tundra heath vegetation in June, July and August, 2 years after treatments were imposed. The boundaries of the boxes represent the 25 and 75 percentiles and the error bars indicate the 5 and 95 percentiles; filled boxes denote intact vegetation, open boxes denote vegetation where the dominant vascular plant species Empetrum hermaphroditum was removed. *Significance at p < 0.05; **significance at p < 0.01; ***significance at <0.001. Linear mixed effects model results are given in Table S5 F I G U R E 2 Net ecosystem CO 2 exchange (NEE), gross primary productivity (GPP) and ecosystem respiration (ER) measured at the plot-level in response to warming and removals at a high (a-c) and a low (d-f) elevation in arctic tundra heath vegetation in August 2016. NEE and GPP are standardized to 600 PAR. Positive NEE 600 and GPP 600 signifies net CO 2 uptake from the atmosphere into the ecosystem. The boundaries of the boxes represent the 25 and 75 percentiles and the error bars indicate the 5 and 95 percentiles; blue boxes denote ambient temperature, red boxes denote warming treatment by OTCs, filled boxes denote intact vegetation, open boxes denote vegetation where the dominant plant species Empetrum hermaphroditum has been removed. *Significance at p < 0.05; ***significance at <0.001. D = unstandardized mean difference between treatment (W = Warming, R = Removal, W × R = Warming × Removal, n = 10) and control plot means, n = 10, ±95% confidence interval (CI) across the study system. Linear mixed effects model results of main and interactive treatment effects are given in Table S7 of Warming and Removal on GPP 600 , GPP 600NDVI and ER (Table S7).
Both GPP 600 , GPP 600NDVI were reduced in Warming × Removal plots, but not in Removal plots and Warming plots, relative to the Control plots, while the strongest reduction in ER in response to Removal occurred in Warming × Removal plots (Table 1, Figure 2; Figure S5).
For the high elevation plots, NEE 600 was negatively correlated with soil moisture and positively correlated with B. nana cover (Table 2). Furthermore, GPP 600 and ER were negatively correlated with soil temperature, while ER and ER 10 were positively correlated with NDVI. For the low elevation plots, NEE 600 was positively correlated with V. vitis-idaea cover, NEE 600at10 was positively correlated to B. nana cover and GPP 600 , ER and ER 10 were positively correlated to NDVI and B. nana cover (Table 2).

| D ISCUSS I ON
Despite differences in temperature between the low and high elevation sites, NEE 600 was similar at both elevations, in contrast to our first hypothesis which predicted that the low elevation would be a stronger net CO 2 sink than the high elevation site. Our second hypothesis, which predicted that short-term experimental warming and dominant species removal would result in a weaker net CO 2 sink, was also unsupported.
Instead, we found that temperature and soil moisture were correlated with CO 2 flux variables only at the high elevation site, potentially indicating a greater sensitivity to changing abiotic conditions for high elevation communities. However, NEE 600at10 , GPP 600 and ER responded to the interactive effect of short-term (3 years) experimental warming and removal. As these responses were similar at both the high and low elevation sites, our third hypothesis, which predicted interactive effects of removal, warming and elevation on NEE, was also unsupported. Hence, our results are more in line with previous studies finding NEE to be unresponsive to warming (Cahoon et al., 2016;Hobbie & Chapin, 1998). However, they reveal that responses of instantaneous CO 2 flux variables to short-term experimental warming can interact with dominant plant species in a similar manner for plant communities that are adapted to different temperature regimes.
Global changes alter the distribution of species, how those species interact with one another and the ecosystem functions that those TA B L E 1 Unstandardized mean difference (D) ± 95% confidence interval (CI) for instantaneous CO 2 fluxes (μmol CO 2 m −2 s −1 ) between treatment and Control plots across two elevational sites  Note: Independent data points are individual plots at each elevation; n = 20 for all except RH and AirT where it is 15 and 17 at the high and low elevation, respectively, and for SoilT where n = 18. * , **Correlation coefficients significantly different to 0 at p ≤ 0.05 and p ≤ 0.01 respectively (in bold).
Abbreviations: AirT, air temperature (°C); ER, ecosystem respiration; ER 10 , ecosystem respiration standardized at 10°C; GPP 600 , gross primary productivity standardized at 600 PAR; NEE 600 , net ecosystem CO 2 exchange standardized at 600 PAR; NEE 600at10 , NEE 600 standardized at 10°C; RH, relative air humidity (%); SoilT, soil temperature (°C). species maintain. We found similar NEE 600 and NEE 600at10 across arctic tundra heath vegetation at a low and a high elevation site at maximum plant biomass during the growing season, where the low elevation have ~2°C higher air temperatures and higher total plant biomass compared to the high elevation (Blüme-Werry et al., 2018). GPP 600 and ER did not differ significantly between elevations across this tundra heath landscape. Furthermore, average annual air temperature has increased ~1°C over the last 100 years in this region (Kivinen, Rasmus, Jylhä, & Laapas, 2017)-a similar increase in air temperature to that simulated in our, and other, warming experiments (e.g. Dorrepaal et al., 2009). However, neither ER, ER 10 or GPP 600 responded to short-term warming, which is likely related to the relatively small (<1°C) change in soil temperature in response to warming. Hence, the lack of any effect of short-term experimental warming on NEE 600 and NEE 600at10 may reflect an overall weak response of GPP 600 to the increase in air temperature imposed by the warming treatment ( Figure 2; Table S1), as well as insulation of the soil by the above-ground vegetation from the increased air temperature caused by the OTCs.
The dominant species in this system, E. hermaphroditum, is a relatively unproductive but widespread species in high latitude ecosystems (Tybirk et al., 2000). Dominant species usually exert large influences on ecosystem functioning (Grime, 1998 (Ylänne et al., 2015). Our results may suggest that other species contributed more to net ecosystem CO 2 exchange at the time of measurement in our system, in line with results from a long-term removal experiment in boreal forest that involved experimental removals of E. hermaphroditum (Wardle & Zackrisson, 2005).
In that study, removal of E. hermaphroditum had only weak effects on soil biota and soil processes compared to faster growing Vaccinium species (Wardle & Zackrisson, 2005). In our study, the cover of the faster growing deciduous shrub B. nana was also related to CO 2 fluxes, reinforcing evidence that B. nana plays a central role in regulating NEE in arctic tundra (Cahoon et al., 2016;Metcalfe & Olofsson, 2015).
Clearly, NDVI and LAI are useful for predicting GPP and modelling NEE and GPP across large environmental gradients in Arctic tundra (Shaver, Street, Rastetter, van Wijk, & Williams, 2007;Street et al., 2007 (Aerts, 2010). Furthermore, E. hermaphroditum is an ericaceous shrub that can reduce seedling establishment, germination and growth of other plants (González et al., 2015) and in line with previous findings (Aerts, 2010) we did not find strong support for rapid increases in cover of neighbouring plant species following E. hermaphroditum removal. Additionally, these tundra heaths are slow growing, and some bare soil patches were still evident 3 years after removal. Therefore, that neighbouring plant species contributed more to GPP 600 following E. hermaphroditum removal was not likely a primary mechanism underlying our findings and our results rather point to an overall low GPP 600 of E. hermaphroditum at peak biomass in this study system.

| CON CLUS IONS
Despite the difference in temperatures caused by elevation and experimental warming, and dominant species removal significantly reducing NDVI and LAI, we found no main effect of any of these factors on NEE 600 or NEE 600at10 . Our results instead revealed that, while warming may interact with dominant species removal to affect several CO 2 flux variables, the dominant species E. hermaphroditum appeared to play a limited role in governing NEE at peak biomass in this ecosystem. These results suggest that warming can influence GPP and ER following removal, or loss, of a dominant species in a similar manner across tundra heath communities adapted to different climatic conditions. They also suggest that vegetation indices that are useful for modelling NEE and GPP across tundra communities and vegetation types spanning a wide range in NDVI Street et al., 2007) are not always suitable for assessing C balances (Valentini et al., 2000) at the scale of our study. Hence, even when substantial changes in NDVI and LAI occur as a consequence of changes in climate, biotic interactions or extreme events in these tundra heath ecosystems (Bokhorst et al., 2015;Callaghan et al., 2013) such changes might not always result in landscape-scale responses of NEE even under global warming.

ACK N OWLED G EM ENTS
The