We selected regions of interest (ROI) with different spectral signatures for sampling and analysis of C and N contents and for analysis of molecular structure information by solid-state ¹³C-CPMAS-NMR spectroscopy. Regions of interest were sampled based upon score projections of the first three components of a principal component analysis (PCA) of the 186 bands in the VNIR spectra, which generally accounted for more than 95% of spectral variability in the cores. We sampled two ROI per core that received straw as a treatment (N=16). Because of the greater visual heterogeneity of the green manure, which could be roughly identified as wood-like or non-woody particles, we increased the number of ROI in these samples (N=31), sampling one to two ROI of visually identifiable wood-like material and two ROI of visually identifiable non-wood material per core (i.e. bark). Wood samples with a diameter larger than 5 mm were not sampled. In total, we sampled 47 ROIs, corresponding to straw or green manure at different stages of decomposition. All material was milled using an agate mortar and pestle prior to further analysis.

Solid-state ¹³C cross-polarisation magic angle spinning nuclear magnetic resonance spectra of the ROI samples were obtained at a frequency of 50.3 MHz (Bruker DSX 200). The powdered samples were filled into a zirconium dioxide rotor and spun at a magic angle under a magnetic field of 6800 Hz with an acquisition time of 0.01024 s. A ramped ¹H pulse was applied during 1 ms contact time to prevent Hartmann-Hahn mismatches. An average of 5500 scans were accumulated with a delay time of 1 s according to the amount of sample and the carbon content.

Tetramethylsilane was equalised at 0 ppm as a reference for the chemical shifts. Phase and baseline correction of the acquired spectra were applied after Fourier transformation. The spectra were then integrated in following chemical shift regions according to ((Beudert et al., 1989) with slight adjustments according to Mueller & Koegel-Knabner (2009): -10 to 45 ppm (alkyl C), 45 to 110 ppm (O/N-alkyl C), 110 to 160 ppm (aromatic C) and 160 to 220 ppm (carboxyl C), spinning sidebands were included. We calculated the ratio of alkyl C to O/N-alkyl C (hereafter noted alkyl ratio) based on these integrated shift regions.

Total carbon and nitrogen contents of the ROI samples were measured by dry combustion (HEKAtech EuroEA 3000). Calibration was made against sulphanilamide (C₆H₈N₂O₂S, 41.8 % C and 16.3 % N) and BBOT (C₂₆H₂₆N₂O₂S, 72.5 % C and 6.5 % N). C/N ratio (hereafter noted C/N) is expressed as the mass ratio of the two elements. All measurements were performed with at least two analytical replicates. An additional replicate was analysed when the difference between analytical replicates exceeded 2 mg C g^-1.

The alkyl ratio and C/N were modelled using the 186 band standard normal variate transformed VNIR reflectance spectra as predictor variables in random forest and artificial neural networks algorithms. The consistency of the spectra of a sample in the profile before and after powdering was checked before modelling. After removing five samples with an Euclidean distance of > 1.5 between the spectra from the core and the powder, four samples corresponding to features with overexposure or shadows leading to low quality spectral information, and four samples from a core replicate with inaccurate spectra, 34 samples were retained as the training dataset. The VNIR spectral data of each ROI were extracted and the mean spectra calculated, which were then used as predictors for POM composition.

Predictive modelling was performed using two algorithms implemented with R software version 4.0.3 (R Core Team, 2021). Firstly, an unconstrained random forest (RF) algorithm was implemented after optimisation for the number of trees in the forest using the ‘party’ package (Hothorn et al., 2005; Strobl et al., 2007, 2008). Secondly, we implemented a random ensemble of fully-connected artificial neural networks (ensemble ANN) and optimized for the number of models in the ensemble. Individual ANN modelling was performed using the ‘neuralnet’ package (Fritsch et al., 2019). Lastly, we evaluated model averaging from the RF and ensemble ANN results, using both weighted and unweighted estimate averaging, with weighting applied proportional to the goodness-of-fit of the individual model estimates.

To implement the random ensemble ANN we independently trained an ensemble of neural network models, in which each model was initialised with random starting weights and an architecture of two hidden layers with a random number of nodes. The number of nodes was randomly selected between 3 and 185 for the first hidden layer and 2 to one less than the number of nodes in the first layer for the second hidden layer. The ensemble approach was applied in order to overcome the instability and high sensitivity of ANN model performance to initialization which is to be expected for our small datasets (Kolen & Pollack, 1991). The random structure of the ANN models was selected to avoid optimisation and over-fitting issues given the small data-set used for model fit. The random ensemble was performed using n-out-of-n bootstrapping and was independently optimized for the number of models for each of the variables using the out-of-bag samples in each bootstrap. Optimization was performed using the MSE and percent variance explained by the model and the optimal number of models was selected according to the stability of model variance with increasing number of models. In the ANN models, the target variables were scaled to a range of 0-1 prior to model fitting and estimates rescaled to the original scale after prediction.

Model evaluation was performed using n out of n bootstrap estimation for both the RF and ANN models, with evaluation done using the out-of-bag model estimates using percentage of variance explained by the models as well as root mean-squared error (RMSE) and mean average error (MAE) to determine the best models. For the ANN, we evaluated both the mean and median of the out-of-bag estimates for predictive performance.

Before applying the predictive model to the images of whole soil cores, pixels were classified into two groups (POM or mineral soil) using a spectral angle mapper classification algorithm in ENVI version 5.2 (Exelis Visual Information Solutions, Boulder, Colorado) with a spectral angle threshold at 0.2 rad. The models were then used to predict the alkyl ratio as well as the C/N for each pixel classified as POM. After predicting the distribution of alkyl ratio or C/N in each core, values were checked for plausibility against the calibration range as well as against published literature values for these soil parameters.

Statistical differences in the chemical characteristics of the organic materials were tested using 3- way ANOVA to evaluate the effects of incubation, type of added OM and depth on organic material chemical composition. Orthogonal contrasts were subsequently computed to estimate the effect of incubation on the different types of OM at the two depths before and after incubation. For the predictions of C/N and alkyl ratio based on hyperspectral imaging, the differences in the distributions between treatments and times were tested based upon the data from all pixels classified as POM using the Kolmogorov-Smirnov non-parametric test.

The alkyl ratio, alkyl, aryl and carboxyl functional group contributions respectively increased by 0.04-0.06, 1.8-2.6%, 3.3-6.0% and 2.2-2.6% as a result of incubation, while O/N alkyl and C/N decreased by 8.2-10.7% and 37-46. The changes in alkyl ratio and alkyl contributions resulting from the incubation were greater for green manure, while the changes in C/N, O/N alkyl, aryl and carboxyl contributions were greater for straw (Table 1). For all molecular indicators, we measured a higher variability in the data for green manure compared to straw, as highlighted by the large standard deviations presented in Table 1.

The ANOVA of the molecular indicators measured on ROI indicated significant differences (p<0.05) between OM composition at T0 and T1, whereas differences between green manure and straw were significant only for O/N alkyl and aryl contributions, with higher O/N alkyl and lower aryl contributions for straw compared to green manure (Table 1). No significant differences in POM composition between depths were found. No significant interactions between incubation time, depths and type of added OM were found so 3-way ANOVA models without interactions were applied for all variables presented in Table 1.

For predicting C/N, the random forest was optimal with 200 trees, whereas the random ensemble ANN was optimal with 100 ANNs in the ensemble. For predicting alkyl ratio, the random forest was optimal with 400 trees, and the random ensemble ANN was optimal with 200 neural networks in the ensemble. Both modelling approaches yielded accurate predictions, with high to very high coefficients of determination (R² > 0.7). The predictions for C/N were better than for alkyl ratio, and in general, the random ensemble ANN outperformed the RF models, and the median estimates of the ANN were better than the mean estimates (Table S1).

For C/N, model averaging improved the predictions, with the best estimates obtained using a weighted average of the RF and median of the random ensemble ANN (R²=0.90, RMSE 14.3, MAE 10.4). For alkyl ratio predictions, the random ensemble ANN outperformed the RF and model averaging yielded no improvement as assessed by the goodness-of-fit parameters. In predicting the measured values close to the maximum and minimum, the RF models tended to underestimate the higher and overestimate the lower, while the ensemble ANN models did not show this tendency (Figure 1 and Fig. S1).