Marc Kennedy Clive Anderson, Stefano Conti, Tony O’Hagan

1. Applications of Bayesian sensitivity and uncertainty analysis to the statistical analysis of computer simulators for carbon dynamics Marc Kennedy Clive Anderson, Stefano Conti, Tony O’Hagan Probability & Statistics, University of Sheffield

2. Outline • Uncertainties in computer simulators • Bayesian inference about simulator outputs • Creating an emulator for the simulator • Deriving uncertainty and sensitivity measures • Example application • Some recent extensions

3. Uncertainties in computer simulators • Consider a complex deterministic code with a vector of inputs and single output • Use of the code is subject to: • Input uncertainty • Code uncertainty

4. Input uncertainty • The inputs to the simulator are unknown for a given real world scenario • Therefore the true value of the output is uncertain • A Monte Carlo approach is often used to take this uncertainty into account • Sample from the probability distribution of X • Run the simulator for each point in the sample to give a sample from the distribution of Y • Very inefficient…not practical for complex codes

5. Code uncertainty • The code output at a given input point is unknown until we run it at that point • In practice codes can take hours or days to run, so we have a limited number of runs • We have some prior beliefs about the output • Smooth function of the inputs

6. Bayesian inference about simulator outputs • Bayesian solution involves building an emulator • Highly efficient • Makes maximum use of all available information • A single set of simulator runs is required to train the emulator. All sensitivity and uncertainty information is derived directly from this • The inputs for these runs can be chosen to give good information about the simulator output • A natural way to treat the different uncertainties within a coherent framework

7. Inference about functions using Gaussian processes • We model as an unknown function having a Gaussian process prior distribution • h(.) is a vector of regression functions and are unknown coefficients Prior expectation of the model output as a function of the inputs

8. Inference about functions using Gaussian processes • We model as an unknown function having a Gaussian process prior distribution • c(.,.) is a correlation function, which defines our beliefs about smoothness of the output and is the GP variance Prior beliefs about covariance between model outputs

9. Choice of correlation function • We use the product of univariate Gaussian functions: • Where is a measure of the roughness of the function in the kth input

10. roughness = 0.5

11. roughness = 0.2

12. roughness = 0.1

13. roughness = 0.01

14. Conditioning on code runs • Conditional on the observed set of training runs, is still a Gaussian process, with simple analytical forms for the posterior mean and covariance functions

15. 2 code runs

16. 2 code runs

17. 2 code runs Large b Small b

18. 3 code runs

19. 5 code runs

20. More about the emulator • The emulator mean is an estimate of the model output and can be used as a surrogate • The emulator is much more… • It is a probability distribution for the whole function • This allows us to derive inferences for many output related quantities, particularly integrals

21. Inference for integrals • For particular forms of input distribution (Gaussian or uniform), analytical forms have been derived for integration-based sensitivity measures • Main effects of individual inputs • Joint effects of pairs of inputs • Sensitivity indices

22. Example Application

23. Sheffield Dynamic Global Vegetation Model (SDGVM) • Developed within the Centre for Terrestrial Carbon Dynamics • Our job with SDGVM is to: • Apply sensitivity analysis for model testing • Identify the greatest sources of uncertainty • Correctly reflect the uncertainty in predictions

24. Photosynthesis Loss • Terrestrial carbon source if NEP is negative • Terrestrial carbon sink if NEP is positive Plantrespiration Loss Soilrespiration Net Ecosystem Production (CARBON FLUX)

25. Leaf life span Leaf area Budburst temperature Senescence temperature Wood density Maximum carbon storage Xylem conductivity Soil clay % Soil sand % Soil depth Soil bulk density Some Inputs Parameters

26. Main Effect: Leaf life span If leaves die young, NEP is predicted to be higher, on average. Why?

27. Main Effect: Leaf life span (updated) If leaves die young, SDGVM allowed a second growing season, resulting in increased carbon uptake. This problem was fixed by the modellers

28. Main Effect: Senescence Temperature Large values mean the leaves drop earlier, so reduce the growing season Small values mean the leaves stay until the temperature is very low

29. When soil bulk density was added to the active parameter set, the Gaussian Process model did not fit the training data properly

30. NEP 80 60 40 20 0 -20 0 500000 1000000 1500000 Bulk density Bulk density Before… After… Error discovered in the soil module

31. Our GP model depends on the output being a smooth function of the inputs. The problem was again fixed by the modellers

32. SDGVM: new sensitivity analysis • Extended sensitivity analysis to 14 input parameters (using a more stable version) • Assumed uniform probability distributions for each of the parameters • The aim here is to identify the greatest potential sources of uncertainty

33. NEP (g/m2/y) NEP (g/m2/y)

34. Leaf life span 69.1% by investing effort to learn more about this parameter, output uncertainty could be significantly reduced Water potential 3.4% Maximum age 1.0% Minimum growth rate 14.2% Percentage of total output variance

35. Extensions to the theory

36. Multiple outputs • So far we have created independent emulators for each output • Ignores information about the correlation between outputs • We are experimenting with simple models linking the outputs together • This is an important first step in treating dynamic emulators and in aggregating code outputs

37. Dynamic emulators • Physical systems typically evolve over time • Their behaviour is modelled via dynamic codes • wherexare tuning constants andztare context-specific drivers • Recursive emulation ofytover the appropriate time span shows promising results

38. CENTURY output ( ) and dynamic emulator ( )

39. Aggregating outputs • Motivated by the UK carbon budget problem • The total UK carbon absorbed by vegetation is a sum of individual pixels/sites • Each site has a different set of input parameters (e.g. vegetation/soil properties), but some of these are correlated • This is a multiple output code • Each site represents a different output • Bayesian uncertainty analysis is being extended, to make inference about the sum

40. References • For Bayesian analysis of computer models: • Kennedy, M. C. and O’Hagan, A. (2001). Bayesian calibration of computer models (with discussion) J. Roy. Statist. Soc. B, 63: 425-464 • For Bayesian Sensitivity analysis: • Oakley, J. E. and O’Hagan, A. (2004). Probabilistic sensitivity analysis of complex models: A Bayesian approach. J. Roy. Statist. Soc. B, 66: 751-769