The PDF report

The PDF report offers a visualization of the numerical benchmarking results. It is structured by the different metrics used in the settings. Each metric has its own section where the result from all outputs using this metric are displayed.

The front page shows some meta information on the benchmarking. Of particular importance is the simulation table, that maps short simulation identifiers to the runs considered in the benchmarking. Usually these identifier are abbreviations of the simulation names, as defined in the lpjml config. If these names are not unique, abbreviations of the simulation paths are used. The identifiers are employed because the simulation names and paths are often quite long, making it more convenient to use short abbreviations to refer to the different runs in the report, especially when comparing the baseline with multiple under test runs.

The benchmark result object

The benchmarking function also outputs its numerical results as a benchmark result object. Accessing this data is a bit technical, so skip this section if you are only interested in the PDF report.

The benchmark result behaves similarly to a nested R list. The top layer is again structured by the different metrics. The data that is generated by a metric is stored in its var_grp_list attribute. It contains the results for all the variables that the metric evaluated as a list of so called variable groups.
Each variable group consists of the processed baseline and (possibly multiple) under test outputs of that variable, as well as the comparisons of all under test outputs against the baseline. See ?Metric for details on how the metrics organize their processed data.

Say we are interested in how the time series of global vegetation carbon has changed in the newer version. The data from the baseline run can be accessed via

vegc_old_ts <- benchmark_result$GlobSumAnnAvgTimeseries$var_grp_list$vegc$baseline
vegc_old_ts

Note that this is an LPJmLData object of lpjmlkit. Evaluating an LPJmLData object in the R Console will provide you with a lot of useful meta data, about the origin and processing history of the object. Note that subset is TRUE and space aggregation is weighted_sum, which reflects the processing steps done by the benchmarking function.

You will encounter several elements of the lpjmlkit package in lpjmlstats - in fact lpjmlstats completly builds upon lpjmlkit.

This also means that we can employ the functionality of this latter package to work with the benchmark results. For example, to output the first 5 years of vegc_old_ts we use

subset(vegc_old_ts, time = 1:5)$data

For under test results the simulation identifiers are used to access the data. This is needed to differentiate the under test outputs if multiple runs are benchmarked against a baseline run.

The simulation table that relates the identifiers to simulation names, as well as other meta data can be accessed via:

get_benchmark_meta_data(benchmark_result)

Note that the simulation with the new version has the identifier 5.8.. Thus, to access the new vegc time series of the new LPJmL version we use

vegc_new_ts <- benchmark_result$GlobSumAnnAvgTimeseries$var_grp_list$vegc$under_test$`5.8.`

To check the difference between the two time series we use

subset(vegc_new_ts, time = 1:5)$data - subset(vegc_old_ts, time = 1:5)$data

Observe that the newer version produced less vegetation carbon in the first five years.

The data from the benchmarking object in combination with the functionality of lpjmlkit is also useful to make comparisons to literature data. Say we have a reference dataset of global soil nitrogen given as a terra object. Again the tools from lpjmlkit come in handy to convert the benchmark data to terra:

soiln_lpjmldata <- benchmark_result$TimeAvgMap$var_grp_list$soiln$under_test$`5.8.`

soiln_terra <- lpjmlkit::as_terra(soiln_lpjmldata)

soiln_terra

# read in validation data terra

# ...

Another convenient feature of the benchmark data object is that all plot functions used for the report generation can also be called individually from here. For example the table of the GlobSumTimeAvgTable metric can be generated via

benchmark_result$GlobSumTimeAvgTable$plot()

This enables a terminal-only benchmarking without generating plots or reports. See ?benchmark for how to skip the PDF report generation in the benchmarking.

Sometimes individual plots of the benchmarking may be needed in another context. The following chunk generates the vegc plot.

benchmark_result$TimeAvgMap$plot()[[1]] # vegc is the first output defined in the settings

With the exception of the table, all plots generated in the benchmarking are ggplot objects, which enables modifications of all visual plot elements. For example a different title can be added to the litter carbon plot with the following code.

vegc_map <- benchmark_result$TimeAvgMap$plot()[[3]] # litc is the third output of the settings

vegc_map + ggplot2::ggtitle("Difference in litter carbon; lpjml 5.8.6 - lpjml 5.7.1")

Maybe the preceding discussion convinced you that working with the benchmark object is fun and useful. It provides a lot of interfaces to further processing steps and functionality. It may thus be valuable to save the object for later analysis. This last step of the tutorial can be done for example with

saveRDS(benchmark_result, file = "benchmark_result_vignette.rds")

This file can also be used to recreate the benchmark report. See ?create_pdf_report for details.