Modelling the habitat preference of two key Sphagnum species in a poor fen as controlled by capitulum water content

Current peatland models generally treat vegetation as static, although plant community structure is known to alter as a response to environmental change. Because the vegetation structure and ecosystem functioning are tightly linked, realistic projections of peatland response to climate change require the inclusion of vegetation dynamics in ecosystem models. In peatlands, Sphagnum mosses are key engineers. Moss community composition primarily follows habitat moisture conditions. The known species habitat preference along the prevailing moisture gradient might not directly serve as a reliable predictor for future species compositions, as water table fluctuation is likely to increase. Hence, modelling the mechanisms that control the habitat preference of Sphagna is a good first step for modelling community dynamics in peatlands. In this study, we developed the Peatland Moss Simulator (PMS), which simulates the community dynamics of the peatland moss layer. PMS is a processbased model that employs a stochastic, individual-based approach for simulating competition within the peatland moss layer based on species differences in functional traits. At the shoot-level, growth and competition were driven by net photosynthesis, which was regulated by hydrological processes via the capitulum water content. The model was tested by predicting the habitat preferences of Sphagnum magellanicum and Sphagnum fallax – two key species representing dry (hummock) and wet (lawn) habitats in a poor fen peatland (Lakkasuo, Finland). PMS successfully captured the habitat preferences of the two Sphagnum species based on observed variations in trait properties. Our model simulation further showed that the validity of PMS depended on the interspecific differences in the capitulum water content being correctly specified. Neglecting the water content differences led to the failure of PMS to predict the habitat preferences of the species in stochastic simulations. Our work highlights the importance of the capitulum water content with respect to the dynamics and carbon functioning of Sphagnum communities in peatland ecosystems. Thus, studies of peatland responses to changing environmental conditions need to include capitulum water processes as a control on moss community dynamics. Our PMS model could be used as an elemental design for the future development of dynamic vegetation models for peatland ecosystems.


Fig. 1 Framework of Peatland Moss Simulator (PMS).
Calculation of Sphagnum growth 161 To model grid cell biomass production and height increment, we assumed that capitula were the 162 main parts of shoots responsible for photosynthesis and production of new tissues, instead of the 163 stem sections underneath. We employed a hyperbolic light-saturation function (Larcher, 2003) to 164 calculate the net photosynthesis, which was parameterized based on empirical measurements (1) 168 where subscript 20 denotes the variable value measured at 20 °C; Rs is the mass-based 169 respiration rate (µmol g -1 s -1 ); Pm is the mass-based rate of maximal gross photosynthesis (µmol 170 g -1 s -1 ); PPFD is the photosynthetic photon flux density (µmol m -2 s -1 ); Bcap is the capitulum 171 biomass; and αPPFD is the half-saturation point (µmol m -2 s -1 ) for photosynthesis.

172
By adding multipliers for capitula water content (fW) and temperature (fT) to Eq. (1), the net The increase in shoot biomass drove the shoot elongation: where Hc is the shoot height (cm); Hspc is the biomass density of Sphagnum stems (g m -2 cm -1 ) 205 and Sc is the area of a cell (m 2 ).  We implemented Tests 5-6 to test the importance of parameters that directly control the species 348 ability to overgrow another species with more rapid height increment (i.e. Pm20, Rs20, αPPFD and 349 Hspec) in lawn and hummock conditions. We eliminated the species differences in the parameter 350 values to be same as those in S. magellanicum (Test 5) and same as those in S. fallax (Test 6).

351
The effects of the manipulation were compared against those from Tests 3-4. For each of Tests 352 3-6, 80 Monte-Carlo simulations were run using the setups described in Test 1.

353
Test 7-8 were implemented to separate the effects of photosynthetic water-response 354 parameters from the effects of the water retention of capitula. We set the photosynthetic water-355 response parameters to be the same as those in S. magellanicum (Test 7) and same as those in S.

Testing model robustness 394
Test 2 addressed the model robustness to the uncertainties in several parameters that closely 395 linked to hydrology and growth calculations. Modifying most of the parameter values by +10% 396 or -10% yielded marginal changes in the mean cover of species in either hummock or hollow 397 habitat (Table 4). Reducing the moss carpet and peat hydraulic parameter n had stronger impacts 398 on S. fallax cover in hummocks than in lawns. Nevertheless, changes in simulated cover that 399 were caused by parameter manipulations were generally smaller than the standard deviations of 400 the means i.e. fitting into the random variation. In Tests 5 and 6, the species differences in the growth-related parameters were eliminated.  In Tests 9, the model failed to simulate the preference of S. magellanicum to hummocks (  could be considered as an elemental design for future DVM development. 590 We see PMS as an elemental design for the future development of dynamic vegetation models 591 for peatland ecosystems, yet there are certain uncertainties and features that should be developed further. We used a gas-exchange-based method to quantify the simultaneous changes in capitula 593 water potential, water content and carbon uptake of Sphagnum moss capitula. It should be noted 594 that, the measurements mainly covered the changes from RWCopt towards RWCcmp (Table 1  conditions. In our modelling, we used the volumetric water content of moss carpet to estimate 600 RWC as an approximation for wet conditions (Eq. 17). The accuracy of such approximation for 601 high RWC conditions remains ambiguous and more information is still required.

602
We assumed that tissue structure did not change during the measurement process, and that the Consequently, PMS could underestimate capitula water potential towards the drying end for 608 those species, if a constant ra is derived from the maximal evaporation rate (Em, Eq. 5; Fig 3C).

609
The water-retention relationship in PCM may not sufficiently capture water potential changes at 610 wet and dry extremes (e.g., S. magellanicum in Fig. 4C). Water retention functions developed

Wcmp
Capitulum water content at the compensation point g g -1

Wmax
Maximum water content of capitula g g -1

Wopt
Optimal capitulum water content for photosynthesis g g -1 Wcf field-water contents of Sphagnum capitulum g g -1 Wsf field-water contents of Sphagnum stem g g -1 WTm Measured multi-year mean of weekly water     grid cell level. In Test 9-10, the photosynthetic water-response parameters (i.e. aW0, aW1 and aW2. 946 See Table 2) were set to be the same as those in S. magellanicum (Test 9) and same as those in S. where θs is the saturated water content; θr is the permanent wilting point water content; α is a 976 scale parameter inversely proportional to mean pore diameter; n is a shape parameter; and m=1-977 1/n. where Ts is the temperature of the moss carpet; z is the thickness of the moss layer (z = 5 cm).

1011
The latent heat flux was calculated by the "interactive scheme" (Daamen and McNaughton, where Δ is the slope of the saturated vapor pressure curve against air temperature; λ is the latent 1015 heat of vaporization; E is the evaporation rate; VPDb is the vapor pressure deficit at zb; rb,x is the 1016 total resistance to water vapor flow, the sum of boundary layer resistance (rsa,x) and surface 1017 resistance (rss); and A is the available energy for evapotranspiration and Ax = Rnx -Gx. 1018 The calculations of γ, λ and VPDb require the temperature (Tb) and vapor pressure (eb) at the 1019 mean source height (zb). These variables were related to the total of latent heat (∑λEx) and 1020 sensible heat (∑Hx) from all surfaces using the Penman-type equations: where ρaCp is the volumetric specific heat of air; raero is the aerodynamic resistance between zb 1025 and the reference height za, and was a function of Tb accounting for the atmospheric stability and was calibrated in Section 2.5.

1043
The value of Ei was calculated as: where subscript 20 denotes the variable value measured at 20 °C; Rs is the mass-based dark 1158 respiration rate (µmol g -1 s -1 ); Pm is the mass-based rate of maximal gross photosynthesis (µmol 1159 g -1 s -1 ); and αPPFD is the half-saturation point (µmol m -2 s -1 ), i.e., PPFD level where half of Pm is 1160 reached. The measured morphological and photosynthetic traits are listed in Table 2. when A dropped to less than 10% of its maximum. Each measurement lasted between 120 1175 and180 minutes. Each sample was weighed before and after the gas-exchange measurement, then dried at 40°C for 48 h to determine the biomass of capitula (Bcap).

1177
The GFS-3000 records the vapor pressure (ea, kPa) and the evaporation rate (E, g s -1 ) 1178 simultaneously with A at every second (Heinz Walz GmbH, 2012). The changes in Wcap with 1179 time (t) was calculated as following: We assumed that the vapor pressure at the surface of water-filled cells equaled the saturation where λ is the latent heat of vaporization; γ is the slope of the saturation vapor pressure -1187 temperature relationship; ra is the aerodynamic resistance (m s -1 ) for vapor transport from inter-1188 leaf volume to headspace; rs is the surface resistance of vapor transport from wet leaf surface to 1189 inter-leaf volume. The bulk resistance for evaporation (rbulk) was thus calculated as ra+rs. 1190 We assumed that the structures of tissues and pores did not change during the drying process 1191 and assumed ra to be constant during each measurement. A tended to increase with time t until it 1192 peaked (Am) and then decreased (Fig. 2B). The point A=Am implied the water content where 1193 further evaporative loss would start to drain the cytoplasmic water, leading to the decrease in A.

1194
The response of A to Wcap was fitted as a second-order polynomial function (Robroek et al.,  Similarly, the evaporation rate (E) increased from the start of measurement until maximum 1227 evaporation Em, and then decreased (Fig. B2C). The point E=Em implied the time when the wet 1228 capitulum tissues were maximally exposed to the air flow. Therefore, ra was estimated as the 1229 minimum of bulk resistance using Eq. (B5), by assuming ei(t)≈es when E(t) = Em: Based on the calculated ei(t), we were able to derive the capitulum water potential (h) where R is the universal gas constant (8.314 J mol -1 K -1 ); M the molar mass of water (0.018 kg 1235 mol -1 ); g is the gravitational acceleration (9.8 N kg -1 ); ei/es is the relative humidity; h0 is the 1236 water potential due to the emptying of free-moving water before measurement (set to 10 kPa 1237 according to Hayward and Clymo, 1982).