Introduction

Dating of findings is the basis for the chronological coordination of archaeological work. Historians also depend on dating analyses to prove the information of archive documents. Besides stratigraphic methods for relative dating, there are several approaches for chronometric dating. Radiometric dating methods are based on the decay of alpha and beta emitters within the sample. They comprise thermo—and radioluminescence used for dating inorganic materials (Aitken 1994; Huntley et al. 1985; Roberts and Lian 2015). For organic matter, in particular, the decay of 14C—referred to as radiocarbon dating—has become the most important standard dating method (Arnold and Libby 1949; Becker 1992; Libby 1961). Dendrochronology is based on another approach. The relative width of tree rings depends on climatic growth conditions. Therefore, sequences of tree rings can be linked together in chronologies. The term was first introduced by A. Douglass. The first relevant chronologies were built at the beginning of the 20th century in different parts of the world (Douglass 1941; Speer 2010; Stokes and Smiley 2008; Webb 1983).

The processes used for the present dating tool can be seen as the first steps within this transition during sub-fossilization. The molecular decay of wood in the context of subfossil preservation conditions has been well recorded. To reduce the influencing factors of the decay, the present research focused on conservation conditions that are suboptimal for microbial activity (in comparison to biodegradation in litter and soil). These preservation conditions can be divided into cold, waterlogged (Björdal et al. 1999, 2000; Björdal and Nilsson 2002), dry (Blanquette 2003), and toxic (Rowell and Barbour 1990). Often, combinations of preservation conditions lead to synergistic effects and even better reduction in microbial activity (e.g., cold and waterlogged or waterlogged and toxic). Aging processes differ according to the preservation conditions, but dry and cold conditions, in particular, seem to preserve the original cell wall principle components (but also secondary metabolites) quite well over a long time. A rough prediction model for age detected by changes in wood density has been presented by Guyette and Stambaugh (2003).

Fourier transform infrared (FTIR) spectroscopy has been proven to be a powerful tool for detecting chemical changes over time (Fink 2017). The basic principle is the absorption of IR radiation by the chemical bounds at their specific resonant frequency. The method provides the chemical fingerprint of a material. Trojanowicz (2008) concluded that infrared spectroscopy is one of the most promising methods to be employed in archaeometry for dating purposes. The method was applied to establish a dating tool for the aged paper of historical documents (Trafela et al. 2007). It has been widely used in geosciences to characterize not only organic remnants, but also inorganic matter of biological origin in fossils (Bobroff et al. 2016). A prediction model for age has been presented by Inagaki et al. (2018) including two buried logs of wood and covering the time period from 220 to 630 AD.

For statistical analysis of spectral data, the random forest method was applied. In the field of machine learning, tree-based methods for regression are well established (Gareth et al. 2013; Hastie et al. 2017). Although tree-based methods are simple and useful for interpretation, they lack prediction accuracy, i.e., they produce good predictions on the training set but are likely to overfit the data, leading to poor test set performance. The random forest method was introduced by Breiman (2001) to overcome these drawbacks. Like bagging, this approach generates multiple trees, which are then combined to yield a single consensus prediction. Combining a large number of trees reduces the variance and will often result in dramatic improvements in prediction accuracy.

The objective of this work is the establishment of a dating tool for wooden artefacts under specified storage conditions based on the relation between chemical characteristics of the decay and time. The combined biotic and abiotic chemical decay is revealed by means of infrared spectroscopy—a rapid and cheap analytical method. Infrared spectroscopic investigations lead to huge data pools. PLS (partial least squares) regression models are state of the art to analyze such data pools (Brereton 2007). An alternative modeling approach for data of this kind is the random forest method. Dendrochronology serves as the reference method.

Materials and methods

Materials

Samples were taken from existing sample pools out of the Tree Ring Laboratory at the Institute of Wood Technology and Renewable Materials, University of Natural Resources and Life Sciences, Austria, representing tree species/genera Picea abies, Larix decidua, Quercus sp., and Abies alba and the region of Austria (Table 1, Online Resource Figures S-2 to S-5).

Table 1 Origin, preservation conditions and number of samples

Preservation conditions at the sampling sites

The oldest available living trees of Central and Eastern Austria were sampled. Spruce and larch samples came from “Dachstein” (47° 34′ N, 13° 38′ E), oak from “Wienerwald” (48° 11′ N, 16° 14′ E) and fir from “Kalkalpen National Park” (47° 44′ N 14° 27′ E). Construction wood was sampled in buildings mainly from roof or ceiling constructions within historical building research. Samples from different buildings in the cities of Salzburg (47° 48′ N, 13° 02′ E), Hallein (47° 41′ N, 13° 06′ E), St. Stephen’s Cathedral in Vienna (48° 12′ N, 16° 22′ E), and the castle in Heidenreichstein (48° 52′ N, 15° 07′ E) were included. The small crater-like lake “Schwarzer See” (47° 31′ N, 13° 49′ E; 1450 m a.s.l.) is surrounded by steep rocky slopes and cliffs. Dead and broken trees fell into the cold water body; preservation conditions minimized microbial attack immediately (Kłusek et al. 2015). The Roman harbor near Weyregg at Lake Atter (47° 54′ N, 13° 34′ E; 482 m a.s.l.) is situated close to the borderline in shallow water. Mean surface water temperature from 1976 to 2015 was 10.5 °C (± 6 °C SD, ehyd.gv.at, last accessed: 2018-10-22). Samples from Pottenbrunn (48° 14′ N, 15° 42′ E; 245 m a.s.l.) were buried in gravel below the groundwater level. Mean groundwater temperature from 1993 to 2015 was 10.5 °C (± 1 °C SD, ehyd.gv.at, last accessed: 2018-10-22). Samples from Hallstatt came from the prehistoric salt mine. Details are given in Grabner et al. (2015) and Tintner et al. (2016).

Methods

Dendrochronological reference

All samples were collected and analyzed within the last 25 years and dated according to standard procedures (Grabner et al. 2018; Speer 2010; Stokes and Smiley 2008). From living trees, usually two cores per tree taken at breast height were available. Roof constructions were sampled using a dry wood corer—one core per piece of timber. The same procedure was used for the prehistoric mining timber from Hallstatt. Partially, disks were also available. Disks from the tree trunks found in the lake (Schwarzer See) were carefully dried and measured under dry conditions. The same procedure was applied to archeological samples (Pottenbrunn, Weyregg). All samples were stored in the laboratory under dry conditions without any conservation agents or preservatives. No extraction step was made prior to the measurements.

Fourier transform infrared (FTIR) spectroscopy

FTIR spectra were recorded in the ATR (attenuated total reflection) mode in the mid infrared area (4000–400 cm−1) with an optical crystal of a Bruker ® Helios FTIR micro-sampler (Tensor 27). This device allows spot measurements with a spatial resolution of 250 Âμm. Cores and cross sections of the disks were measured tree ring-wise on randomly chosen tree rings, so the age detected by means of dendrochronology could be assigned to a certain spectrum. It was attempted to cover the widest span of tree rings available against the background of limitation discussed in the section “When does the MD-Dating model fail?”. The crystal was placed directly on the polished surface of a certain tree ring. Thirty-two scans were collected at a spectral resolution of 4 cm−1. Spectra were vector-normalized using the OPUS © (version 7.2) software. Smoothing and second derivative spectra were obtained using The Unscrambler × 10.1 (© Camo) by applying the Savitzky–Golay algorithm (Savitzky and Golay 1964) with a polynomial order of 2 and 13 smoothing points.

Statistical methods

All statistical analyses were done using the statistical computer software language (R Development Core Team 2011). Partial least squares (PLS) regression was done using the R package pls (Mevik et al. 2016). Further, the R package randomForest (Liaw and Wiener 2002) was used to fit a random forest model to the data. A supplementary method description is given in the Online Resource document.

Results

Development of the “MD-Dating” tool

For band localization, it was started with the complete wavenumber range of the vector-normalized spectra. The PLS regression model was used to identify relevant bands for the final prediction models. Main areas are the band regions between 1771 and 1690 cm−1, and 1270 and 800 cm−1. Additionally, the band region of methylene groups between 2970 and 2800 cm−1 was taken into consideration. All regions have already been used to describe molecular decay processes in wood (Tintner et al. 2016; Traoré et al. 2016). All three band regions were further used to build up the dating tools for Picea, Abies and Quercus. To visualize changes in the spectra and especially in the region between 1780 and 1700 cm−1, spectra of randomly selected measurements on randomly selected cores are displayed in Online Resource Figure S-1. The region can be assigned to C=O valence vibrations of acetyl or carboxyl groups (Schwanninger et al. 2004). The chemical decay of these functional groups is paralleled by a decrease in the band intensity. This decay is dominated by the factor “time”; microbial degradation plays a subordinate and comparable role among the preservation conditions included in the model. According to the loadings plot of the PLS regression (not shown), an additional region between 1690 and 1610 cm−1 was included for the Larix model. The band region is affected by the carbonyl-stretch vibration of resin acids (Tintner and Smidt 2018). As the root-mean-squared errors (RMSE) of the random forest models were always smaller than the ones of PLS regression, only the results of the random forest approach are presented here.

The combination of ATR-FTIR spectroscopy and the adequate statistical evaluation tool, as applied to a series of ancient wood, leads to a new dating tool. The MD-Dating tool (abbreviation for Molecular Decay Dating) is applicable to the following general conditions: It covers approximately 3000 years. Chemical changes start directly after the formation of wood. This was demonstrated by the oldest available living trees (labeled in the Figures with “living”). The prediction of these living trees cannot be separated from the predictions of construction wood measurements. Furthermore, the conditions of construction wood (“constr”) in buildings and of waterlogged wood (“water”) at mean temperatures below 12 °C are covered. The three conditions cannot be separated within the final models (Fig. 1). It demonstrates that they lead to comparable degradation processes. Chemical differences between heartwood and sapwood, earlywood and latewood, as well as different amounts of extractives reasoned by different growing situations do not affect the molecular changes over time in a way that age prediction by the models becomes meaningless. All these effects contribute to the RMSE (Table 2). Tintner et al. (2016) already displayed that molecular decay over time and chemical differences between only earlywood measurements and narrow rings can be separated in a principal component analysis as they load on different components.

Fig. 1
figure 1

Calibrated models for age based on molecular decay versus dendrochronological age: Dashed lines indicate the 95% confidence bands, and dotted lines mark the 95% forecast bands. Different symbols indicate the different preservation conditions. aPicea abies, n = 1324, bQuercus spp., n = 767, cLarix decidua, n = 895, dAbies alba, n = 496

Table 2 Cross-validated RMSE of prediction and calibration (in years)

Out of the spectral data and dendrochronology as a reference, random forest models for the prediction of age were created. The random forest model underestimates the age of wood dated by dendrochronology before the beginning of Common Era (see Fig. 4, although the figure is displayed for another purpose discussed later). This effect is true not only for the tree species Picea abies but also for all other species. Hence, the final models were subjected to a further calibration step to improve the prediction quality [described in the Supplementary Information, especially Eqs. (1 and 2)]. This procedure is referred to as MD-Dating model.

The model quality is represented by the RMSE. Both approaches, the MD-Dating model and the random forest model alone, were compared by tenfold cross-validation. The averaged RMSE is presented in Table 2.

Further, it was assessed how applicable the models are to different species/genera. In Table 3, the RMSE of each fitted random forest model applied to the other three tree species/genera is stated. The columns in Table 3 refer to the different models, while the rows indicate the different species/genera. Again, the RMSEs for both approaches, the MD-Dating model and the random forest model alone, are given. It became evident that models of different species/genera lead to significantly weaker results.

Table 3 RMSE of prediction and calibration (in years)

Although the collection of generated trees is much more difficult to interpret than a single tree, one can obtain an overall summary of the importance of each predictor using prediction accuracy. The random forest algorithm estimates the importance of a variable by looking at the degree to which the prediction error increases when data for that variable are permuted while all others are left unchanged. The necessary calculations are carried out tree-by-tree as the random forest is constructed. The decrease in accuracy as a result of this permuting is averaged over all trees and is used as a measure of the importance of that variable in the random forest. For each tree species, the numbers of the 30 most important variables, i.e. wavenumbers, used in the random forest models divided according to the different spectral regions are given in Table 4. The results explain why the use of specific models for different species/genera is necessary. While the molecular changes in Picea are reflected in the region between 1771 and 1690 cm−1, the Quercus model takes most variables out of the region between 1271 and 800 cm−1 corresponding to numerous bonds originating from all principal components of the cell wall. The detailed lists of the most important wavenumbers are given in the Online Resource Figure S-6.

Table 4 Numbers of the 30 most important wavenumbers used in random forest models

Additionally, approximate confidence and forecast bands were calculated. The final models are presented in Fig. 1. The calibrated years forecasted by the MD-Dating model versus the dendrochronologically measured age for all four tree species together with the corresponding confidence and forecast bands are shown. The dashed lines that are only visible in Fig. 1d indicate the 95% confidence bands while the dotted lines mark the 95% forecast bands. The different plotting characters indicate the different preservation conditions. In particular, the measurements of the living trees are strongly straightened on the line of best fit. Therefore, the symbols can hardly be distinguished.

When does the MD-Dating model fail?

For the practical application before including samples into a prediction model, it has to be checked whether the frame conditions are fulfilled. Three instances can lead to a failure of prediction. Two effects were observed on samples that disturb the decay kinetics in the significant band regions. Therefore, spectra influenced by these effects were excluded from the data set as outliers. The third one refers to a preservation condition different from the ones included in the MD-Dating models presented here.

Brittleness

The first effect is fungal deterioration leading to brittleness. Very brittle parts of the sample lead to strongly different—partly lower and partly increased—band intensities. Gray spectra in Fig. 2 represent spectra from a very obvious brittle part of the core. The black spectra were recorded from the parts directly adjacent to them, but not or less brittle. Typical bands of celluloses and hemicelluloses are diminished and typical lignin bands are increased.

Fig. 2
figure 2

Exemplary spectra of a sample with brittle parts (gray spectra) and not or less brittle parts (black), a whole spectral range and b zoom into the important spectral region around 1735 cm−1

Brittleness in this context is understood as dried-out waterlogged wood that contains quite heavy decay by either soft rot and/or erosion bacteria, which creates a very fragile network of the wood structure. Upon air-drying the wood collapses because the cells are too fragile to withstand the capillary forces of the evaporation water. Between this fragile material that collapses and appears brittle after drying, there is very often a broad zone of intermediate decayed wood. Some cells are decayed others are solid (or partly solid). This partly decayed wood will only shrink, not collapse, and does not appear brittle. If such wood is used for dating, the changed chemical composition will most likely disturb the dating result. Therefore, the measurements were currently restricted to the macroscopically intact wood that did not display any deformation under the pressure of the ATR crystal. Further research will have to clarify how to detect microbial attack that changes the molecular decay to an extent that the measurements cannot be used in the present models.

Compartment air

The second effect occurs at tree rings of construction wood near to the wood surface that is exposed to the compartment air (e.g., the room under the roof in case of wooden roof structure). Wood is affected up to a depth of about 0.5 cm. The accelerated decrease in the significant band at 1735 cm−1 can be explained by the increased access of oxygen (Fig. 3).

Fig. 3
figure 3

Exemplary spectra of a sample with an intensified aging of the tree ring nearest to compartment air (spectrum in black)

Salt environment

Preservation of wood in living trees and under waterlogged and dry conditions leads to comparable chemical changes. Other conditions can lead to different speeds or even different behavior over time. Samples were included from a prehistoric salt mine in Hallstatt, Upper Austria, where the preservation conditions are characterized by a highly salty clay environment. Figure 4 clearly reveals a different behavior of the samples stored under these conditions. Looking at the gray data points, the model shows an underestimation of samples older than the beginning of Common Era, visible by an upwards bending of the data cloud. The bending itself is less problematic, as it can be considered by a calibration step (as done for the final prediction model). However, samples from salt environment seem to follow a systematically biased degradation process. Even if we do not expect fundamentally different parameter estimates for the resulting MD-Dating model, it is suggested that a particular model should be built for this preservation condition. Figure 4 shows the estimated years predicted by the random forest model versus the years measured by dendrochronology for the tree species Picea abies including 51 measurements from salt environment.

Fig. 4
figure 4

Uncalibrated prediction model for Picea abies including also 51 measurements from the prehistoric salt mine in Hallstatt, Upper Austria. The solid line indicates perfect fit, different symbols indicate the different preservation conditions

Discussion

A big strength of the method is the strictly monotonous decay over different (but not all) preservation conditions. This distinguishes it from other methods, especially from dendrochronology, where it might happen that different results far away from one other appear equally probable. With regard to 14C dating, there are time spans at which calibration curves bend into a plateau as well (Aitken 1997). In comparison with 14C dating, the comparably far lower costs here are another main advantage. As FTIR technology is a low-cost technology and the measurements can be taken within some minutes, costs per measurement can be kept low, or a higher number of measurements can be done for the same price. Costs per sample might be less than 10% of 14C dating costs. Dendrochronology might fail, as the tree ring pattern does not fit into the chronology with a sufficient probability. Reasons can be different factors of growth apart from climate (frequent harvests in the forests nearby the tree, very special microclimatic effect, pedogenic reasons, not enough tree rings, etc.). Such material can be dated via molecular decay anyhow.

It is advisable to measure on more than just a single spot. RMSE values of Table 1 refer to a single measurement. Multiple measurements can be taken with almost the same effort. The exact calculation of confidence intervals of the averaged value is difficult, as taking additional measurements on one artefact only reduces errors due to variability within a piece of wood and does not take into consideration the variability among different artefacts. It is evident, however, that an averaged value of more measurements will have a smaller error. Besides the described problem with near surface layers, there is no need for measuring in different tree rings as required with dendrochronology. Artefacts with broad extension but only few tree rings can also be dated.

Anyhow, plenty of waterlogged wood suffers from severe fungal deterioration leading to very brittle sections. Such material absconds from the present models. The problem of compartment air excludes most artefacts of artistic-historical interest like small statues, panel paintings and lots of furniture. The chemical wood composition must not have undergone strong alterations due to microbiotic agents (bacteria or fungi = brittle wood) or insects (there must be no pest damages in the analyzed wood), but neither due to superficial oxidative events related to storage conditions. Currently, all models are referenced via dendrochronology. It will become an interesting task to establish models for species, where no chronologies exist. Therefore, sample sets referenced by other methods will have to be compiled.

A careful focus has to be placed on different preservation conditions. The presented models include only construction wood and wood from cold, sweet waterlogged conditions. In particular, marine wood or wood from bog sites cannot be dated by the presented models. Future work will have to focus on such material, but also different preservation situations have to be investigated in more detail. Pedersen (2015) and Pedersen et al. (2015) present FTIR spectra of spruce wood from an early 17th century moat that was filled up with soil after few years of service with a strongly diminished band between 1770 and 1700 cm−1. Such situations have to be particularly assessed from a geochemical and microbial point of view.

Anyhow, there is a considerable amount of samples that can be assessed by the present models. Pizzo et al. (2015) demonstrated a good behavior of the holocellulose/lignin-ratio independent of different storage conditions over a time span corresponding to the present research. FTIR measurements were taken in a very comparable way. It proves that time-dependent processes can be assessed monotonously, even if processes of molecular changes are influenced by a large number of factors.

Conclusion

This article presents a novel method to date four tree species common in Europe with the perspective of further prediction models for other species and other spatial validity. Future work will focus on extension in a spatial sense and toward further preservation conditions, especially wooden artefacts frozen in glaciers and buried in peat bogs. Models for further species/genera will also be established. As molecular decay processes of lignocellulosic materials might be comparable, the method is thought to be applicable to other plant materials. The combination of the spectral method and random forests statistics is a promising approach to stimulate and support the work of building historians, archaeologists, and even environmental scientists.