Prediction of peat properties from transmission mid-infrared spectra

A better understanding of peatland dynamics requires more data on more peat properties than provided by existing databases. These data needs may be addressed with resource-efficient measurement tools, such as models that predict peat properties from mid-infrared spectra (MIRS). High-quality spectral prediction models are already used for mineral soils, but similar developments for peatland-focused research lag behind. Here, we present transmission-MIRS prediction models for peat that are openly available, easy to use, include quality checks to assess prediction quality, and propagate prediction errors. The models target element contents (C, N, H, O, P, S, K, Ca, Si, Ti), element ratios ($\mathrm{C}/\mathrm{N}$, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$), isotope values (δ¹³C, δ¹⁵N), physical properties (bulk density, loss on ignition (LOI), macroporosity, non-macroporosity, volume fraction of solids, hydraulic conductivity, specific heat capacity, dry thermal conductivity), thermodynamic properties (Gibbs free energy of formation ($\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$)), and nominal oxidation state of carbon (NOSC). They are representative for a more diverse set of peat samples than currently existing peat-only models while having a competitive predictive accuracy. Relatively accurate predictions can be made, for example, for many element contents (C, N, O, S, Si, Ca, $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, $\mathrm{O}/\mathrm{C}$, $\mathrm{H}/\mathrm{C}$, bulk density, and LOI). Many of these properties are not predicted by existing high-quality prediction models focusing on mineral soils. For some of the target variables, high-quality prediction models focused on mineral soils exist. These models may be more accurate, but reported predictive accuracies are not directly comparable because the training data is imbalanced in the number of organic versus mineral soil samples. We suggest that some soil properties are easier to predict for peat, whereas others are easier to predict for mineral soils, emphasizing that we need new approaches to meaningfully compare prediction errors of spectral models computed on datasets with variable amounts of organic soil samples. Our tests also indicate that replacing δ¹³C and δ¹⁵N measurements with MIRS models probably is unlikely to be feasible due to large prediction errors. Future studies should address the lack of open training and validation data for some peat properties (O, H, NOSC, $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, LOI, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$), the lack of mineral-rich peat samples, and improve and standardize model validation and comparison for models trained on data with very different proportions of peat soils. This study is a step to catch up with high quality standards set by models for mineral soils and provides models for several peat properties for which we could not find descriptions of previous models in the literature. By filling data gaps in the Peatland Mid-Infrared Database, we make a step towards providing the data required to better understand peatland dynamics.

1 Introduction

Northern peatlands are projected to both become a larger sink and source for greenhouse gases (Frolking et al., 2011; Qiu et al., 2022). However, these projections have large uncertainties: They range from an estimated net emission of 0.2 to an estimated net carbon sequestration of 0.1 Pg C yr⁻¹ (Qiu et al., 2022). Also peatland C stock estimates are highly uncertain, ranging from ∼ 270 to 1045 Pg (Yu, 2012; Nichols and Peteet, 2019). Some fraction of this uncertainty is caused by a lack of knowledge on processes, their incomplete representation in models (Loisel et al., 2021; Qiu et al., 2022), a high spatial variability of peat properties (Frolking et al., 2011; Loisel et al., 2014, 2017; Normand et al., 2021; Qiu et al., 2022), and a lack of systematic data collections.

More data on peat properties would help reduce errors in estimating peat stocks and greenhouse gas emission from peatlands. Such data could be used to interpret and compare individual peat cores, quantify global peat properties at high spatial resolution, constrain parameter estimates in peatland models, and test peatland models and hypotheses on peatland processes (Loisel et al., 2014). Existing databases on global and regional peat properties partly address these needs (Zoltai et al., 2000; Charman et al., 2013; Loisel et al., 2014; Treat et al., 2016; Gallego-Sala et al., 2018; Leifeld et al., 2020). However, they still have large gaps in terms of spatial coverage (e.g. Loisel et al., 2017) and the peat properties considered. For example, one of the largest database on peat properties (Loisel et al., 2014) contains peat ages, bulk densities, C and N contents, but other peat properties, such as nutrient contents, hydraulic properties, and availability of labile organic matter also play a pivotal role in peatland models (e.g., Yu et al., 2001; Bauer, 2004; Frolking et al., 2010; Baird et al., 2012; Mahdiyasa et al., 2022; Qiu et al., 2022) and are not part of this database. Even though similar databases exist for such data, they are much smaller, containing not only fewer samples, but also covering smaller geographical areas (e.g., Liu and Lennartz, 2019; O'Connor et al., 2020).

A major obstacle towards addressing these limitations is that data collection is difficult in practice: Collection and analysis of peat cores are costly and resource limitations constrain both the number of peat samples collected and the set of peat properties measured. Examples for such resource limitations are restricted access to measurement devices, lack of funding, lack of time, or simply the fact that the amount of peat in a given sample is limited; this makes it difficult to measure all peat properties at high spatio-temporal resolution (e.g., Zaccone et al., 2018).

An approach to address this problem has been developed for mineral soils: Soil spectroscopy combined with spectral prediction models can reduce the time and cost of the measurements of diverse soil properties at an accuracy partly competitive to traditional measurement methods (Viscarra Rossel et al., 2008; Nocita et al., 2015; Hartmann and Nopmanee, 2019; Sanderman et al., 2020). This approach requires comparatively few sample preprocessing steps and only one infrared spectrum needs to be measured to predict soil properties as diverse as carbon content, bulk density, and pH value (e.g., Nocita et al., 2015; Dangal et al., 2019). In addition, spectral prediction models can be combined with pedotransfer functions to predict soil properties which are only rarely measured together with spectra (e.g., McBratney et al., 2006). In case of peat, this would be in particular useful to predict peat hydraulic and thermal properties since these are rarely measured together with spectra (e.g., Teickner et al., 2025a). Within the last years, databases, modeling approaches, and user interfaces have been improved to an extent that many traditional measurements can be replaced or supplemented at a trade-off for accuracy that is acceptable for some applications (e.g., Shepherd et al., 2022; Hengl et al., 2021).

Even though aforementioned databases and models also comprise peat, peat soils are still underrepresented in prediction models (even the largest databases contain no more than several hundred to 1000 samples which may be classified as peat based on their C content, see for example Sanderman et al., 2020; Helfenstein et al., 2021; Hengl et al., 2021). Moreover, existing mineral soil databases target soil properties linked to agriculture and forestry. Peatland research often targets stocks, indicators for peat degradation, and modeling ecosystem processes and for this requires additional soil properties, for example loss on ignition, total sulfur contents, stable isotope abundances, and porosity. Finally, many existing studies use different spectroscopic measurement techniques and spectral ranges (transmission mid-infrared spectroscopy, visible and near-infrared spectroscopy, or mid-infrared diffuse reflectance spectroscopy). In summary, existing spectral libraries and prediction services can be very useful to predict some peat properties (e.g., C, N content, bulk density) with some spectral data, but there still is a need for complementary spectral libraries and prediction models dedicated to the specific questions and approaches of peatland research.

There are many spectral prediction models exclusively for peat, but these are still in an early phase of development in comparison to those for mineral soils (Downey and Byrne, 1986; Bergner and Albano, 1993; McTiernan et al., 1998; Chapman et al., 2001; Ludwig et al., 2006; Rinnan and Rinnan, 2007; Artz et al., 2008; Terhoeven-Urselmans et al., 2008; Lang et al., 2009; Laiho et al., 2014; Hayes et al., 2015; Pérez-Rodríguez et al., 2016; Hodgkins et al., 2018; Straková et al., 2020; Helfenstein et al., 2021; Teickner et al., 2022). Major limitations in comparison to state of the art studies for mineral soils are the following: First, many of the computed models are probably unrepresentative for many new samples and likely have low reliability and predictive accuracy because they were computed with samples from only few peat cores (e.g., McTiernan et al., 1998; Chapman et al., 2001). Second, even though there exist many prediction models, there are some often used peat properties for which no prediction models for peat have been developed yet (e.g., some nutrient contents, C and N isotope values, peat hydraulic, thermal and thermodynamic properties). Third, several newer and important concepts in state of the art mineral soil studies have not yet been applied in the development of prediction models for peat: prediction domains (Wadoux et al., 2021; see also Sanderman et al., 2020 and Shepherd et al., 2022 for some recent examples), propagation of prediction errors (e.g., Viscarra Rossel et al., 2008; Shepherd et al., 2022; Padarian et al., 2022), open accessibility of model data and code (Hengl et al., 2021; Shepherd et al., 2022), and software that allows non-experts to generate predictions and check their reliability (Hengl et al., 2021; Shepherd et al., 2022). We will discuss each of these concepts in turn.

The prediction domain (sometimes referred to as calibration space, predictor space, or feature space) is the range of the spectra for which a model has been computed and validated (Wadoux et al., 2021). If a new spectrum is outside this range, a prediction model may still produce accurate predictions, but there is no guaranty for this because the model was not trained and tested on such extreme samples. Checking whether a new spectrum is within the prediction domain is therefore a first rule of thumb to decide whether to trust predictions and their uncertainties (e.g., Shepherd et al., 2022).

Many studies which use spectral prediction models to analyze peat ignore prediction errors. This can lead to overly confident inferences (e.g., Hodgkins et al., 2018; Moore et al., 2019; Baysinger et al., 2022). Prediction models should propagate errors, such that they can be appropriately included in subsequent analyses. Computing prediction errors is both a matter of developing statistical approaches which allow to compute reliable estimates of prediction errors and of software which allows scientists to easily access prediction errors.

Making raw spectra and code for spectral preprocessing and computation of models accessible is one of the largest current efforts in the development of spectral prediction models for mineral soils (Hengl et al., 2021; Shepherd et al., 2022). This allows others to use, to criticize, and to improve prediction models. The advantages therefore are improved reliability, more engagement between developers and users of models, and faster, more efficient, improvement of models. In contrast, only few of the peatland models we are aware of are openly accessible, including the underlying data (Hodgkins et al., 2018; Teickner et al., 2021, 2022; Teickner and Knorr, 2022a).

Finally, recent projects target the development of estimation services to allow non-experts to easily generate predictions (Shepherd et al., 2022). These estimation services are web applications which estimate soil properties based on spectra you have uploaded, while taking care of spectral preprocessing and reliability checks in the background (e.g., whether spectra are within the prediction domain) (Shepherd et al., 2022). Ultimately, this makes spectral prediction models widely accessible to scientists who can then quickly fill existing data gaps.

To summarize, spectral prediction models for mineral soils currently have higher quality standards than existing models for peat, and are openly accessible and easier to use. Our study aims to develop spectral prediction models for many peat properties relevant in peatland research, which use transmission mid-infrared spectra and meet quality standards set by models for mineral soils. More specifically, we want to:

1. Develop openly available spectral prediction models for many peat properties relevant in peatland research.

2. Provide prediction domains for all these models as reliability check for predictions.

3. Propagate errors from data (only for $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$) and parameters (all target variables) to predictions.

4. Fill data gaps in the pmird database, one of the largest open databases focusing on peat and mid-infrared spectra (Teickner et al., 2025a), to support peatland research.

To this end, we computed spectral prediction models for element contents (C, H, N, O, P, S, K, Ca, Si, Ti), element ratios ($\mathrm{C}/\mathrm{N}$, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$), isotope values (δ¹³C, δ¹⁵N), physical properties (bulk density (BD), loss on ignition (LOI) macroporosity, non-macroporosity, volume fraction of solids, saturated hydraulic conductivity (K_s), specific heat capacity (c_p), dry thermal conductivity (K_T)), standard Gibbs free energy of formation ($\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$), and nominal oxidation state of C (NOSC) (Masiello et al., 2008). To our knowledge, our study is the first test to compute prediction models for some of these variables (O, K, Ca, Si, Ti contents, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$, δ¹³C, δ¹⁵N, $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, NOSC, LOI). The models were computed and validated with data from one of the largest open accessible peat spectral libraries (Teickner et al., 2025a) which makes them representative for more diverse peat samples than most existing models that focus on peat. While some of our models may only generate approximate predictions (especially δ¹³C and δ¹⁵N), many have prediction errors small enough to be useful in diverse applications.

Bayesian data analysis was used to propagate errors from model parameters and data (where measurement errors are known) to predictions. We made the models openly accessible via the R packages irpeatmodels (Teickner, 2025b) and irpeat (Teickner and Hodgkins, 2025). For all computed models and those already implemented in both packages, we also computed prediction domains. Finally, we extended irpeat to allow making predictions with a few lines of code, where spectral preprocessing and basic checks for the reliability of predictions are performed in the background.

Figure 1 Workflow how to use the models developed in this study to estimate peat properties.

The models can be used to estimate peat properties in a workflow as shown in Fig. 1. We demonstrate the usefulness of the models and the workflow to fill gaps in databases by predicting these peat properties for samples in the pmird database. This makes pmird to one of the largest open access collections of, for example, peat element contents (C, H, O, N, P, S, K, Ca, Si, Ti) and thermodynamic properties ($\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$).

2 Methods

2.1 Target variables

Table 1 summarizes the peat properties (target variables) for which we computed models. The variables were selected because they can be used to understand and quantify important processes in peatlands and because many of them cannot be predicted with existing models (e.g., Sanderman et al., 2020; Shepherd et al., 2022).

(e.g., Koerselman and Meuleman, 1996; Olde Venterink et al., 2003; Wang and Moore, 2014; Wieder, 2022)(e.g., Wieder, 2022)(e.g., Moore et al., 2005; Wieder et al., 2016)(e.g., Koerselman and Meuleman, 1996; Olde Venterink et al., 2003; Wang and Moore, 2014; Wieder, 2022)(e.g., Olde Venterink et al., 2003; Wang and Moore, 2014; Wieder, 2022)(e.g., Waughman, 1980; Wang et al., 2015b)(Shotyk, 1996; Hölzer and Hölzer, 1998; Gałka et al., 2022a)(Nadelhoffer and Fry, 1988; Lerch et al., 2011; Biester et al., 2014)(Williams and Flanagan, 1996; Schmidt, 2014; Xia et al., 2020)(Larmola et al., 2010,  )Schmidt (2014)(e.g. Williams and Flanagan, 1996; Asada et al., 2005a; Bragazza and Iacumin, 2009; Drollinger et al., 2019)Drollinger et al., 2019Asada et al., 2005bBiester et al., 2014(Asada et al., 2005b)Masiello et al. (2008)(Masiello et al., 2008)(e.g., Worrall et al., 2018)(e.g., Malmer and Holm, 1984; Kuhry and Vitt, 1996; Biester et al., 2014; Leifeld et al., 2020)(e.g., Bauer, 2004)(e.g., Bader et al., 2018; Zaccone et al., 2018; Moore et al., 2018; Leifeld et al., 2020)Bader et al., 2018; Kim et al., 2003(e.g., Bader et al., 2018; Zaccone et al., 2018; Moore et al., 2018; Leifeld et al., 2020)Bader et al., 2018; Kim et al., 2003(Liu and Lennartz, 2019)(e.g., Chambers et al., 2011)Liu and Lennartz, 2019(Liu and Lennartz, 2019)(e.g., Weiss et al., 2006)Liu and Lennartz (2019)(Liu and Lennartz, 2019)(e.g., Weiss et al., 2006)(e.g., Weiss et al., 2006)Liu and Lennartz (2019)(Liu and Lennartz, 2019)(e.g., Granberg et al., 1999; Weiss et al., 2006)O'Connor et al. (2020)(e.g., Granberg et al., 1999; Weiss et al., 2006)Gnatowski et al. (2022)(e.g., Granberg et al., 1999; Weiss et al., 2006)

Table 1 Overview on the peat properties for which we computed prediction models. “Has model” indicates whether we computed a spectral prediction model for the target variable (“yes”) or a pedotransfer function which uses bulk density or N content to predict the variable (“no”). “Likelihood” is the likelihood used in the prediction model for each target variable (see Sect. 2.4).

2.2 Data sources

The prediction models were computed with a subset of the pmird database (Teickner et al., 2025a; Münchberger, 2019; Münchberger et al., 2019; Schuster et al., 2022; Drollinger et al., 2019, 2020; Agethen and Knorr, 2018; Kendall, 2020; Harris et al., 2023; Harris and Olefeldt, 2023; Pelletier et al., 2017; Teickner et al., 2021, 2022; Broder et al., 2012; Hömberg, 2014; Diaconu et al., 2020; Gałka et al., 2022b, a; Boothroyd et al., 2021; Worrall, 2021; Reuter et al., 2019a, b, 2020; Moore et al., 2005; Turunen et al., 2004). The pmird database is a heterogeneous collection of infrared spectra and other chemical and physical properties of peat, peat forming vegetation, and dissolved organic matter (DOM) compiled from previous studies. For many peat variables, the database has data from several different sites and covers global gradients of conditions under which peat is formed. This makes the models more representative for many peat properties than existing models (Bergner and Albano, 1993; Chapman et al., 2001; Artz et al., 2008; Helfenstein et al., 2021). The spectra that we used are mid-infrared spectra measured in transmission mode on dried and ground peat pressed to pellets together with potassium bromide. Details on the measurement devices and other measurement metadata are summarized, for each model, in Table S1 in the Supplement. C and N, O, H, and some of the S contents were measured with combustion elemental analyzers. Contents of all other chemical elements and of the remaining S measurements were measured using wavelength-dispersive X-ray spectrometry on peat pellets pressed from dried and milled material. Isotope values were measured with isotope ratio mass spectrometry coupled to the elemental analyzers used in the respective studies. Bulk density was measured from known dimensions of peat layers cut from cores and dry masses. Loss on ignition was measured by combustion. Details on analytical methods for individual studies may be found in the metadata stored in the pmird database.

We did not use spectra that were classified as already baseline corrected in the pmird database (Teickner et al., 2025b), except for those from dataset 13 where we checked that the corrected spectra are similar to the result of our procedure here, to avoid that differences in preprocessing would decrease the predictive accuracy of the models. The pmird database contains a comparatively small number of vegetation and dissolved organic matter (DOM) samples and we therefore do not consider our models applicable to DOM or vegetation in general, except for Sphagnum which forms the bulk undecomposed peat material in bogs. Since DOM can have spectral properties different from peat, we did not include DOM samples in our models. Except for some cores, the peat is from ombrotrophic bogs, and peat with larger mineral contents and fen peat is thus underrepresented. An overview on the geographical distribution of the peat cores used to train and test each model is given in Table S2.

The pmird database does not contain data for saturated hydraulic conductivity, total porosity, macroporosity, volume fraction of solids, specific heat capacity, and dry thermal conductivity for samples with MIRS (Teickner et al., 2025a). For these variables, we used modified versions of models from Liu and Lennartz (2019), Liu et al. (2020), Gnatowski et al. (2022), and O'Connor et al. (2020) to predict these variables from bulk density or N content. These models are described in Sect. S1. As suggested in previous studies (e.g., McBratney et al., 2006), these models (also known as pedotransfer functions) can be used with bulk density and N contents predicted from MIRS to predict these physical peat properties, whereby all relevant errors are propagated.

The pmird database also does not contain $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$ estimates. We predicted $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$ from element contents (at least C, H, N, O) using modified versions of the models from Thornton (1917), Patel and Erickson (1981), and Battley (1999) for the enthalpy of combustion and the entropy of formation, as described in Popovic (2019). These models are described in Sect. S2.

2.3 Spectral preprocessing

All computations for this manuscript were made in R 4.3.0 (R Core Team, 2022). For each of the target variables, except the physical peat properties mentioned in the previous section, we computed three spectral prediction models. The three models use differently preprocessed spectra (no derivative, first derivative, second derivative spectra), but otherwise were computed in the same way. All spectral preprocessing was done with the ir package (Teickner, 2022). To harmonize the spectra, we interpolated them to unit wavenumber resolution and clipped them to the range 650 to 4000 cm⁻¹. Next, we conducted an atmospheric correction of water vapor and CO₂ artifacts using the approach suggested in Perez-Guaita et al. (2013). First, we subtracted a baseline created from a Savitzky-Golay smoothed version of the spectra where regions with strong CO₂ peaks (645 to 695 cm⁻¹ and 2230 to 2410 cm⁻¹) (Wallace and NIST Mass Spectrometry Data Center, 1997) were linearly interpolated and then we used CO₂ and water vapor spectra from the pmird R package (Teickner, 2025a) (see also Teickner et al., 2025b) to perform the atmospheric correction as described in Perez-Guaita et al. (2013). Due to differences in devices and measurement conditions, this procedure attenuated CO₂ artifacts, but did not remove them completely. Thereafter, the corrected spectra were baseline corrected using a convex hull procedure (Beleites and Sergo, 2021), normalized using the signal normal variate (SNV), the three versions of derivative spectra were computed, and all spectra were binned with a bin width of 10 cm⁻¹ to reduce the number of redundant predictor variables and reduce possible wavenumber shifts between measurements from different devices. Finally, we excluded intensities from 2250 to 2400 cm⁻¹, to avoid that remaining CO₂ peaks confound predictions.

2.4 Prediction models

We used normal, gamma, and beta distributions as likelihoods (Table 1) and used Bayesian statistics to compute all prediction models. All models were computed with brms (Bürkner, 2018), using a logit (beta regression), log (gamma regression) or identity link function (normal regression), assuming a constant shape parameter (beta, gamma) or standard deviation (normal), using a normal prior for the intercept, gamma priors for the shape parameter or standard deviation, and regularized horseshoe priors (Piironen and Vehtari, 2017a, b) for the slopes (for each predictor variable). The regularized horseshoe prior shrinks coefficients to zero except where they are strongly related to the response variable, conditional on other predictors. To reduce overfitting, we defined a large amount of shrinkage, by assuming that 5 of the 321 predictor variables have non-zero coefficients (Piironen and Vehtari, 2017b). The regularized horseshoe prior can lead to a complex posterior geometry that is difficult to sample from even with efficient sampling algorithms and to run most of the models without divergent transitions, we had to increase the degrees of freedom of the student-t distribution of the horseshoe prior from 1 to 3 or 4. This reduced the number of divergent transitions, but also leads to less regularization (deviation from the horseshoe shape) (Piironen and Vehtari, 2017b) which may lead to overfitting and less interpretable model coefficients. In our case, less regularization was not a critical limitation because our aim was not to interpret model coefficients, but to optimize predictive accuracy, and because our model validation did not indicate overfitting (Table 2).

The posterior distributions were estimated with Markov Chain Monte Carlo (MCMC) sampling with Stan (Stan Development Team, 2021), using 4 chains, 3000 warmup iterations and 2000 sampling iterations per chain. Chains were initialized with pathfinder (Zhang et al., 2022). Maximum Monte Carlo standard errors (Vehtari et al., 2021) for predictions of the target variables by the best models (see the next section) for each target variable are shown in Table 2. As mentioned above, some models had divergent transitions, but at least one model per target variable did not and we only evaluate and interpret models without divergent transitions. The largest rank-normalized $\widehat{R}$ for model parameters was 1.01, indicating convergence of the chains (Vehtari et al., 2021).

Table 2 Overview on the predictive properties of the best models for each target variable. “Unit” is the measurement unit for all numeric variables. “Derivative” indicates which derivative of the spectra was used in the model. “n_train” and “n_test” are training and testing sample sizes. “RMSE_test” is the root mean square error computed for the test data. “bias_test” is the average difference between predictions and measurements. “RMSE_train−RMSE_test” is the difference between the training RMSE and the testing RMSE (RMSE_test). “MCSE_mean”, “MCSE_sd”, “MCSE_lower”, and “MCSE_upper” are maximum estimated Monte Carlo standard errors for the average, standard deviation, lower, and upper 95 % prediction interval for predicted values across the training and testing data. “BD” is bulk density and “LOI” is loss on ignition.

2.5 Model validation and data filtering

We used the Kennard–Stone algorithm as implemented in the prospectr package (Stevens and Ramirez-Lopez, 2013) to split the observations for each target variable into a training and a testing dataset, using the euclidean distance between the underived preprocessed spectra. The number of observations assigned to the training dataset was defined as $min(\mathrm{0.8}n,n_{\text{max}})$, where n is the number of available observations for a target variable and n_max=200. Observations that were not part of the training data were used for model testing. All models for the same target variables use the same observations for training and testing such that the models are comparable.

This procedure was chosen because our aim is to develop prediction models that are applicable to as diverse peat samples as possible, that is, to compute one prediction model with maximum prediction domain (the value ranges covered by all predictor variables) (Wadoux et al., 2021) and smallest possible prediction error across this prediction domain. The Kennard–Stone algorithm maximizes the distance between spectra covered by the training data and therefore selects a diverse training data set.

An ideal test of the prediction models would use test data that covers the whole spectral range of the training data and is independent of the training data, which in the case of peat samples are samples from different peatland sites. We had to deviate from this ideal because of the heterogeneity of the pmird data. In particular, because there are only few samples from few peatland sites with large carbonate or silicate contents, it would have been possible with independent observations either only to test the models over a much smaller range of spectral variation if only independent test data would have been used, or to test the predictive accuracy for spectral conditions the model was not trained on. In the first case, we would risk overfitting in the untested spectral range, and in the second case, the predictive accuracy would be underestimated due to extrapolation. Therefore, as an alternative, we did not separate observations from the same cores or sites when defining training and test data. This allowed us to test the models across a much larger spectral range within their prediction domains. For many target variables, the overfitting risk should be small because both the training and the test data have samples from many different sites.

To compare models, we used the expected log predictive density (ELPD) (e.g., Vehtari et al., 2017) computed on the test data. Model evaluation was performed with the loo package (Vehtari et al., 2019). Following rules of thumb (Sivula et al., 2022), we assumed models to have equivalent predictive performance (according to the capability of our evaluation) when the difference of their ELPD (ΔELPD) is smaller than 4, and otherwise when ΔELPD is larger than two times its standard error (using normal approximation for ΔELPD). Models with divergent transitions were not considered during model evaluation. To give an easier to interpret performance metric, we also computed the root mean square error (RMSE).

We do not interpret model coefficients and how this may reflect causal links between molecular structures and target variables (1) because our model coefficients are not intended to estimate causal effects, (2) because it is very likely that they do not represent causal effects, and (3) because specific wavenumbers cannot be assigned unambiguously to molecular structures (e.g., Stuart, 2004). For those interested in model coefficients, we show a plot of the model coefficients with the best model for each target variable and a table listing possible assignments to molecular structures for coefficients with a posterior probability of being larger than 0 of at least 90 % or a posterior probability of being smaller than 0 of at least 90 % in Sect. S5.

2.6 Prediction domains

A regression model interpolates a target variable within the range of predictor values – the prediction domain (Wadoux et al., 2021). If such a model is used for prediction with new data that are outside the prediction domain, it is unclear how large prediction errors are, particularly for models with high dimensional prediction domain, such as spectral prediction models. Consequently, it should be checked that new data are within the prediction domain of the model (e.g., Roberts et al., 2017), even though this is no guaranty for accurate predictions.

For this reason, we computed the training prediction domain (Wadoux et al., 2021) for each model as the range of the predictor variable values across all training samples (training prediction domain), and a prediction domain for the test samples (test prediction domain) for each model as the range of predictor variable values across all testing samples. When samples are outside the prediction domain, predictions may be less reliable than estimated by the model validation.

The difference between testing and training prediction domain shows where the models need further testing. The difference between training (or testing) prediction domain and the prediction domain formed by all relevant spectra in the pmird database indicates whether the model covers the spectral variability in the pmird database, as approximation of the spectral variability of peat in general, and therefore indicates where additional data can improve the models.

When making predictions with the models, irpeat checks whether the input data are within the testing or training prediction domain. This is a safety device to avoid a misuse of models and it provides information for those who want to improve our models.

This test of MIRS against prediction domains is only a first test because even spectra within the prediction domain may have spectral properties and values for the predicted peat property that are different from the training or test data. With a large enough training and test dataset, such edge cases become more and more unlikely. Additional checks to be provided by future studies are based on a list of error sources identified through targeted tests of the models against such possible edge cases.

2.7 Testing the influence of confounding factors

Based on previous experience in the interpretation of peat MIRS and on the peaks caused by silicates, carbonates, amides, carbohydrates, aromatics, and lipids (Stuart, 2004; Tatzber et al., 2007; Parikh et al., 2014), we suggest that the main gradients in peat chemistry that control spectral variation are (1) the content of silicates, (2) the content of carbonates, (3) the initial vegetation composition that controls differences in the initial content of amides, carbohydrates, and aromatics, and (4) the degree of decomposition, which increases the relative contents of amides, lipids, and aromatics, and decreases the overall content of carbohydrates (e.g., Cocozza et al., 2003). Previous studies suggest that differences in amide contents and silicates can bias predictions (Broder et al., 2012; Teickner and Knorr, 2022b) and similar effects are likely for carbonates, because carbonates cause dominant peaks that overlap with peaks caused by aromatics and amides (Tatzber et al., 2007) and because large carbonate contents usually indicate higher pH values and therefore shifts in carboxyl peaks due to deprotonation (Ellerbrock and Gerke, 2021). To test for such confounding factors, we plotted model residuals versus Ca, Si, and N measured for the same samples (residuals were not plotted for samples where Ca, Si, and N, respectively, were not measured).

2.8 Filling data gaps in the pmird database

To fill data gaps in the pmird database, we used the best models for each target variable (Table 2) to predict missing values, including prediction errors as samples from the posterior predictive distribution, for the target variables where samples have MIRS. These predictions are stored in a published data table. Moreover, we created two additional data tables that indicate, for each prediction, whether the MIRS is in the training or testing prediction domain for the respective model. We restricted gap filling to peat and litter samples with absorbance-FT-MIR spectra. In contrast to the model development, we included spectra that may have already been baseline corrected, since the prediction domains can be used to screen spectra that are not similar to the data used to train and test the models. For variables that can be predicted without MIRS if other data are available (e.g., C, H, O, N, bulk density) with the additional models developed here ($\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, saturated hydraulic conductivity, total porosity, macroporosity, volume fraction of solids, specific heat capacity, and dry thermal conductivity) or previously published models (NOSC, $\mathrm{C}/\mathrm{N}$, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$), we created an additional data table with predictions without MIRS. For $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, we required C, H, O, and N contents to be measured for this; contents of other elements were included if available and otherwise the contents were set to 0 g g⁻¹ when computing $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$.

3 Results and discussion

We developed spectral prediction models for many peat properties relevant in peatland research. In the next subsections, we evaluate these models in terms of their prediction errors estimated on test data, in terms of confounding factors, and in terms of how much of the spectral variation in the pmird database is covered by the models (prediction domains). Where other models for the same target variable are available, we compare prediction errors of our models to that from other studies. Our models are openly available via the irpeatmodels package, functions to use them are available from the irpeat package, and we show how to use these packages. Finally, to provide data that can be useful for future research, we summarize the gap-filling of the pmird database.

3.1 Predictive accuracy

The predictive accuracy for the best models for each target variable is summarized in Table 2 and plots of measured versus predicted values are shown in Fig. 2. Estimates for the predictive accuracy are both worse and better than that of previously published models using spectra in the visible, near infrared or mid-infrared range (Fig. 3), but these estimates are not directly comparable because of different modeling approaches and differences in the variability of data used to train and test the models. Studies that use large databases can use modeling approaches that require more training data but may outperform linear models, such as cubist (Hengl et al., 2021; Sanderman et al., 2020). Another reason for differences in the predictive accuracy is the chemical diversity of the training and testing data. For example, one reason for the better predictive accuracy for C contents, N contents, and $\mathrm{C}/\mathrm{N}$ reported in Chapman et al. (2001) certainly is that samples are from one site only which leads to less confounding between predictors and C content. Similarly, data from Terhoeven-Urselmans et al. (2008), Bergner and Albano (1993), and Artz et al. (2008) also cover smaller gradients in peat properties than the pmird database, in particular no peat with large mineral fractions (Fig. 3).

Not only differences in the range of chemical properties, but also the distribution of observations along chemical gradients can lead to differences in estimated predictive accuracies. For example, models from Sanderman et al. (2020) and Hengl et al. (2021) for C and N contents have a better predictive accuracy and were computed with many peat samples, but the majority of observations is from mineral soils with small C and N contents. It has been repeatedly observed that prediction errors are larger for larger C contents than for smaller C contents (Dangal et al., 2019; Helfenstein et al., 2021). These heterogeneous prediction errors are probably caused by two factors: Firstly, spectra of mineral soil have prominent mineral peaks which allow a more accurate estimation of small C contents, whereas at large OM contents there is a much more complex and diverse pattern of peaks caused by organic matter molecular structures. Secondly, whenever a variable is positive and the majority of values is small, prediction errors are smaller due to the positivity constraint (e.g., Bellon-Maurel et al., 2010). This does not only explain worse estimated predictive accuracy of our models for C and N, but can also explain why our model for bulk density has a better estimated predictive accuracy than models from Sanderman et al. (2020) and Hengl et al. (2021) because the majority of the peat samples have a small bulk density and the positivity constraint therefore implies smaller prediction errors. A last reason for a better predictive accuracy is outlier removal in previous studies not based on specific theoretical considerations (Terhoeven-Urselmans et al., 2008; Sanderman et al., 2020). We did not remove outliers here because we wanted to develop prediction models that are applicable to a diverse range of peat samples, while outlier removal may lead to better predictive accuracy for peat samples with specific chemical properties. Overall, the modeling approach and data properties of our studies are most directly comparable to Helfenstein et al. (2021), who focused on C contents, and here, our model performs similarly well.

In summary, our models have a roughly similar or better predictive performance for some variables as have previous studies focusing exclusively on peat samples. Albeit direct comparison of prediction errors to high-quality models computed with large spectral libraries is not possible in terms of peat, it is very likely that the predictive accuracy for peat properties could be improved with more flexible modeling approaches, such as localization (e.g., Viscarra Rossel et al., 2024; Helfenstein et al., 2021) or different machine learning algorithms, which would however, in many cases, require more balanced data and, for some variables, more data in general.

Figure 2 Measured versus fitted and predicted values for all target variables. For $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, standard errors of prediction estimated by the auxiliary model (Sect. S2) are shown as error bars. Diamond-shaped points with red borders are outliers (discussed in the text).

Figure 3 RMSE for test data from the best models compared to RMSE for test data from other studies (a) and ranges of target variable values in the training data for the models (b). In (a), points are average values and error bars 95 % confidence intervals for individual models (our models) or the ranges of average RMSE for different strata of the data, where different analyses were performed. When studies report only one point estimate, this point estimate is shown was vertical tic. In (b), error bars are ranges for target variable values in the training data.

3.2 Outliers and confounding factors

There are some interesting patterns in the plots of measured versus predicted values (Fig. 2): For C, there are two outliers, one with a measured C content >0.6 g g⁻¹ and one with a measured C content <0.1 g g⁻¹ (Fig. S20). A peat sample with >0.6 g_C g⁻¹ should be decomposed because litter initially has large O and H contents that cause small relative C contents and because preferential decomposition of organic matter fractions with large O and H contents (carbohydrates, phenols) leads to a relative accumulation of C (Schellekens et al., 2015; Moore et al., 2018; Leifeld et al., 2020), but the spectrum does not have pronounced aromatic or lipid peaks one would expect for a peat sample with such high C content (Fig. S20). A peat sample with <0.1 g_C g⁻¹ must have a comparatively large mineral content because undecomposed peat forming litter has much larger C contents, yet the spectrum does not have typical silicate peaks (Fig. S20). We therefore assume that either the C measurements in the pmird database are not correct for these two measurements or that spectra were incorrectly assigned to these samples. For H, the plot indicates overestimation for smaller H contents and underestimation for larger H contents. Samples with smaller H contents have larger Si contents (Fig. S18) which suggests that predictions are confounded by silicate peaks.

For P and K, several samples with the largest measured P and K contents have relatively large prediction errors. These samples are Juncus effusus samples from a short-term NPK fertilization experiment with high nutrient loads (Agethen and Knorr, 2018). While it appears that the models can predict P and K contents reasonably well for these samples, more samples would be required to evaluate whether the model overfits to spectral characteristics of J. effusus litter (Fig. S20) rather than spectral characteristics for high P contents in vegetation.

For Si, there are several samples with larger Si content for which the model overestimates Si contents, besides the four samples with maximum Si content for which the model underestimates Si contents. The overestimated samples have a large Ca content and the underestimated do not (Fig. S16). Presence of carbonates therefore biases predictions of Si contents, even though this bias seems to be small for the training data. The one observation marked as outlier in Fig. 2 may be a measurement error of Si contents or erroneously assigned spectrum to this sample because the sample does not have typical characteristics for silicate rich peat, such as a pronounced peak around 1100 cm⁻¹ and peaks around 1900 cm⁻¹ (Fig. S20) (Parikh et al., 2014).

For S, larger Ca contents do not bias predictions, but they increase the residual variance indicating, similarly to Si, that peaks caused by carbonates confound predictions (Figs. 2, S16). One reason for the difficulty to predict S contents here is that samples in the data with large S contents also contain large Ca contents, but only some samples with large Ca content have carbonate peaks. Samples with large carbonate peaks probably are influenced by groundwater with relative high pH value under which calcite can precipitate with sufficiently high Ca²⁺ concentrations. In contrast, samples with large Ca contents, but without carbonate peaks probably have Ca²⁺ bound to carboxyl groups, which leads to a small carboxyl peak (around 1730 cm⁻¹) and a more pronounced peak around 1640 cm⁻¹ which has contributions by carboxylates (Ellerbrock and Gerke, 2021) (Fig. S19).

It is interesting that prediction of other variables (e.g., C, N, K, P, Ti, bulk density) is possible without such bias. Since the majority of samples with large Ca contents (ca. >15 000 µg g⁻¹) are from cores from one permafrost peatland site, this may be due to overfitting, indicating that more peat MIRS from carbonate rich samples need to be published to improve development and testing of transmission-MIR prediction models for carbonate-rich peat.

Ti is the only target variable for which the training RMSE is significantly larger than the testing RMSE (Table 2). This indicates that Ti contents were more variable in the training data than in the testing data. Since Ti is unlikely to cause detectable peaks in peat MIRS, prediction of Ti concentrations probably relies mainly on a similar atmospheric deposition across the analyzed peatlands and residual enrichment as peat is lost due to decomposition or fires. That even small Ti concentrations can be estimated from MIRS therefore supports application of Ti concentrations as decomposition indicator, even though the model also suggests that there are other sources of variation. In fact, the most useful application of the models may not be the accurate prediction of Ti contents, but the detection of conditions where Ti concentrations are controlled by other factors than residual enrichment from decomposition. Such conditions can be detected by comparing predictions of the model to Ti measurements.

The models for δ¹³C and δ¹⁵N have prediction errors much too large for most current applications of such isotope measurements. C and N isotope values are controlled by many different processes that can easily lead to a large variation in δ¹³C and δ¹⁵N values despite similar spectral properties. For δ¹³C, these factors are differences in δ¹³C signatures of assimilated CO₂ due to isotope fractionation (Williams and Flanagan, 1996; Hobbie and Werner, 2004), the Suess effect (Keeling, 1979), and differences in the fraction of CO₂ assimilated from methanotrophy (Larmola et al., 2010; Schmidt, 2014). For δ¹⁵N, these factors are different pathways via which N is assimilated by peat forming vegetation (Asada et al., 2005b) and the opposite effects of isotope fractionation (Nadelhoffer and Fry, 1988; Lerch et al., 2011; Asada et al., 2005a; Alewell et al., 2011) and enrichment of OM fractions with negative δ¹³C values during aerobic decomposition (Bowling et al., 2008; Xia et al., 2020), which agrees with weak correlations of δ¹⁵N values with peat decompostion indicators (Broder et al., 2012; Biester et al., 2014; Mathijssen et al., 2019; Serk et al., 2022). It may therefore be the case that δ¹³C and δ¹⁵N measurements cannot be estimated accurately with MIRS prediction models. However, many of the δ¹³C and δ¹⁵N measurements in the pmird database do not correct for blank effects (Teickner et al., 2025b), which biases δ¹³C and δ¹⁵N values proportionally to the C and N mass (Langel and Dyckmans, 2017). According to blank corrections for other projects, this bias has a magnitude of ca. 0.2 ‰ to 0.5 ‰ (depending on the sample C mass) for δ¹³C and a magnitude of ca. 0.2 ‰ to 1 ‰ (depending on the sample N mass) for δ¹⁵N. This is only a fraction of the estimated RMSE (Table 2) and we therefore currently assume that unbiased measurements would still result in large prediction errors relative to measurements of δ¹³C and δ¹⁵N. For δ¹⁵N, the contribution of this bias is certainly larger than for δ¹³C. A factor that contributes to the large prediction errors therefore are biased measurements due to uncorrected blank effects.

Figure 4 Prediction domains of the best model for each target variable. The training and testing prediction domains are created from the maximum and minimum intensity values at each wavenumber of the preprocessed spectra across the training or test data. The shaded region “All spectra” is formed in the same way, but with all data that were used for model development. For this plot, to allow easier comparison of prediction domains between target variables, the spectra were all preprocessed and scaled in the same way (zero derivative spectra, scaling of intensities as for C).

For $\mathrm{C}/\mathrm{N}$, prediction errors are larger for samples with larger $\mathrm{C}/\mathrm{N}$. Two factors probably contribute to this pattern: Firstly, large $\mathrm{C}/\mathrm{N}$ values imply small N contents and at large $\mathrm{C}/\mathrm{N}$ values, very small changes in N contents cause large changes in $\mathrm{C}/\mathrm{N}$ values. Such slight changes in N contents probably cause only small changes in peak intensities that are not easy to detect in MIRS and hence not easy to predict. Secondly, some of the samples with large $\mathrm{C}/\mathrm{N}$ values have large silicate contents (Fig. S18). The large differences in spectra between undecomposed Sphagnum peat and mineral-rich peat very likely confounds linear relations present in peat without large mineral contents.

To summarize, for many variables, accurate predictions are possible with the models developed here, whereby the accuracy needed will of course depend on the specific purpose of the analysis. The models for δ¹³C and δ¹⁵N are probably not accurate enough for any analysis of isotope values, even if blank effects are corrected. Besides silicates and nitrogen, calcium – either in the form of carbonates or in the form of Ca²⁺ bound to carboxylates – is an important contributor to spectral variability in peat samples and makes it difficult to predict S, Si, and Ca contents, at least with the modeling approach used here. Since the pmird database contains Ca-rich samples only from few sites, future tests of the models with additional Ca-rich peat samples would be useful.

Table 3 Summary of gap filling the pmird database. For each variable, the table shows the number of measured values (no gap filling needed) (“Measured”), the number of gaps (“Gaps”), the number of gaps that can be filled with the pedotransfer functions (with bulk density or N measurements), with the auxiliary model for $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$ (with element contents), or with formulas for element ratios and NOSC (with element contents) (“Auxiliary models”), the number of gaps that can be filled with MIRS predictions, where spectra are in the training or in the testing prediction domain (“Training domain” and “Testing domain”; some spectra are in both prediction domains and then are counted twice), the number of gaps that could be filled in total (“Filled”), and the fraction of gaps that could be filled (“Filled (%)”). In total, there were 3379 samples with spectra for each target variable.

3.3 Prediction domains

A comparison of training and testing prediction domains for our models shows that the testing prediction domains covers only a small fraction of the training prediction domain for H, O, NOSC, $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$, $\mathrm{C}/\mathrm{N}$, and LOI (Fig. 4). For these peat properties, only few samples were left for model testing and we could therefore not test the predictive accuracy of the models for some parts of the training prediction domain. Future studies should test these blind spots. The difference between testing and training prediction domain can be used to identify samples useful to test the models in the future. Similarly, the difference between the training prediction domain and the domain formed from all spectra identifies samples that would be useful additions to the training data if the target variable would be measured for these samples. For example, for H, O, NOSC, $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$, and $\mathrm{C}/\mathrm{N}$ there is a lack of mineral-rich samples in both the training and testing data, whereas for LOI there are mineral-rich samples in the training data, but not enough mineral-rich samples were left for model testing (Fig. 4). The prediction domains are available from the irpeatmodels package.

3.4 Filling gaps in the pmird database

Table 3 summarizes the results of our gap filling. Even though we filled all gaps with MIRS predictions, predictions that are outside the training and testing prediction domain may be unreliable and therefore we consider all observations for which the spectra is neither inside the training nor the testing prediction domain as unfilled gaps for our evaluation here.

Since the pmird database contains many bulk density and N measurements, but no measurements for porosity, hydraulic conductivity, specific heat capacity, and thermal conductivity, many missing values could be estimated with the pedotransfer functions and the remaining gaps filled with MIRS-predicted bulk density and N. Similarly, a large fraction of fillable gaps for element contents could be imputed, except for H and O for which training and testing prediction domains cover a smaller fraction of the spectral variability. Much fewer gaps could be filled for O, NOSC and $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$ than for H because of differences in the preprocessing for these models: the models for H uses underived spectra, whereas the models for O, NOSC and $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$ use first or second derivative spectra, where high frequency features are more emphasized and therefore many spectra are not within the prediction domains (Table 2). For observations, where spectra are outside the prediction domain only because a few of the variables exceed the prediction domain boundaries by small values are still reliable; the amount of useful predictions is therefore probably underestimated. Compared with existing databases (Zoltai et al., 2000; Charman et al., 2013; Loisel et al., 2014; Treat et al., 2016; Gallego-Sala et al., 2018; Leifeld et al., 2020), this makes the gap-filled pmird database one of the largest available data sources for contents of many elements, hydraulic and thermal properties, and peat chemistry, in particular NOSC, and $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$. The gap filling predictions are available from Zenodo (Teickner and Knorr, 2025a).

3.5 Implementation of the models in the irpeatmodels and irpeat R packages

The models can be used with R by installing the irpeat (Teickner, 2025a) and irpeatmodels (Teickner, 2025b) packages. The irpeatmodels package contains the models itself and irpeat contains functions to interact with the models. This design was chosen to account for faster development cycles for code to interact with the models and code for other functions of the irpeat package and also to account for size limitations for software packages in online repositories.

4 Conclusions

Our aim was to develop transmission-MIRS prediction models for peat with the quality standards of state of the art prediction models for mineral soils. We have developed models for element contents (C, N, H, O, P, S, K, Ti), element ratios ($\mathrm{C}/\mathrm{N}$, $\mathrm{H}/\mathrm{C}$, $\mathrm{O}/\mathrm{C}$), isotope values (δ¹³C, δ¹⁵N), physical properties (bulk density, macroporosity, non-macroporosity, hydraulic conductivity, specific heat capacity, dry thermal conductivity), Gibbs free energy of formation ($\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$), and nominal oxidation state of carbon (NOSC) of bog peat with data from the pmird database, one of the largest spectral databases available for peat soils. The models have predictive accuracies comparable to or better than existing models for peat, while covering a more diverse range of peat spectral chemistry. All models are available via the R packages irpeat and irpeatmodels. irpeat automatically propagates prediction errors and flags predictions where spectra are outside the prediction domain of the models.

A limitation to be addressed in future studies is that more diverse test data and, for some variables, more training data are required to reduce prediction errors and cover more representative chemical gradients of peat, in particular for O and H contents, NOSC and $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$, and peat samples with larger silicate, carbonate, and, in general, mineral contents. This would also make it possible to use more flexible modeling approaches, such as localization or different machine learning algorithms, which could reduce prediction errors. Moreover, it is recommendable to replace pedotransfer functions for peat hydrological and thermal properties by spectral prediction models to reduce prediction errors.

Our model evaluation lets us draw the following conclusions: First, transmission-MIRS probably are not suitable to predict peat δ¹³C and δ¹⁵N very precisely and we attribute this mainly to the diversity of processes that control δ¹³C and δ¹⁵N and their contrasting effects. Second, a comparison to existing models for mineral soils suggests that we need more targeted approaches to meaningfully compare the predictive accuracy of models evaluated on databases with contrasting fractions of organic versus mineral soil samples.

Finally, we used our models to fill gaps in the pmird database, making it one of the largest databases for many peat properties such as hydraulic properties, Ti, O, H, and P contents, NOSC, and $\mathrm{\Delta }{\text{G}}_{\text{f}}^{\mathrm{0}}$. This illustrates that the models represent a first step to provide the data required to quantify diverse peat properties at high spatial resolutions which will in particular be useful to estimate and test process models.