Spatio-temporal statistical models for glaciology

The purpose of this thesis is to develop spatio-temporal statistical models for glaciology, using the Bayesian hierarchical framework. Specifically, the process level is modeled as a time series of computer simulator outputs (i.e., from a numerical partial differential equation solver or an emulator) added to an error-correcting statistical process, closely related to the concept of model discrepancy. This error-correcting process accounts for spatial variability in simulator inaccuracies, as well as the accumulation of simulator inaccuracies forward in time. For computational efficiency, linear algebra for bandwidth-limited matrices is used for evaluating the likelihood of the model, and first-order emulator inference allows for the fast approximation of numerical solvers. Additionally, a computationally efficient approximation for the likelihood is derived. Analytical solutions to the shallow ice approximation (SIA) of the full Stokes equation system for stress balance of ice are used to examine the speed and accuracy of the computational methods used, in addition to the validity of modeling assumptions. Moreover, the modeling and methodology within this thesis are tested on data sets collected by the University of Iceland Institute of Earth Science (UI-IES) glaciology team, including bi-yearly mass balance measurements at 25 fixed sites at Langjökull (a glacier) over 19 years, in addition to 100 meter resolution digital elevation maps. As a byproduct of the construction of the Bayesian hierarchical model, a novel finite difference method is derived for solving the SIA partial differential equation (PDE). Although the application domain of this work is glaciology, the model and methods developed in this thesis can be applied to other geophysical domains. The thesis is structured around three papers. The first of these papers reviews dynamical modeling of glacial flow, introduces a second-order finite difference method for solving the SIA PDE, presents a Bayesian hierarchical model involving this numerical solver, and validates the model with analytical solutions to the SIA PDE. The second of these papers generalizes the statistical model of the first paper, probes higher- order random walks for representing model discrepancy, incorporates first-order emulators, and analyzes methods for efficient log-likelihood evaluation. The third of these papers applies the model framework of the first two papers to mass balance and surface elevation data at Langjökull. The major contributions of the thesis are the derivation of a new numerical method for solving the SIA PDE in two spatial dimensions and time, the use of a random walk to represent model discrepancy (i.e., an error-correcting process), efficient methods for log-likelihood evaluation, and the application of spatio-temporal statistical modeling to Langjökull, one of Iceland’s main glaciers.Markmið þessarar ritgerðar er að þróa tíma- og rýmisháð tölfræðileg líkön fyrir jökla með því að nota stigskipt Bayesísk líkön. Hér er sá hluti stigskipta Bayesíska líkansins sem snýr að undirliggjandi líkani fyrir ferlið sem er verið að skoða, útfærður þannig að tímaraðirnar sem koma frá tölulegum hermi (þ.e. frá tölulegri lausn á hlutafleiðujöfnu eða nálgun á slíkri tölulegri lausn) eru lagðar saman við líkindafræðilegt ferli sem hefur það hlutverk að leiðrétta fyrir mismuninn á milli tölulega hermisins og raunverulega ferlisins, og er nátengd hugmyndinni um misræmi líkana. Þetta líkindafræðilega ferli leiðréttir fyrir rýmisháð frávik og tekur tillit til að frávikin safnist upp yfir tíma. Til að flýta fyrir útreikningum þá er línuleg algebra fyrir rýr fylki notuð til að reikna sennileikafall líkansins og fyrstu gráðu hermar notaðir til að flýta fyrir útreikningum á tölulegum lausnum hlutafleiðujafna eða öðrum kerfum. Að auki er ný reiknisparandi nálgun fundin fyrir sennileikafallið. Fræðilegar lausnir á þunnjökla nálguninni sem byggir á jöfnum Stokes fyrir spennur í jöklum, eru notaðar til að meta reiknihraða og nákvæmni tölulegra nálgana, og hversu vel líkanið fellur að gögnunum. Að auki er líkönunum og aðferðafræðinni í þessari ritgerð beitt á raunveruleg gagnasöfn sem Jarðvísindastofnun Háskóla Íslands hefur sett saman, þar með taldar afkomumælingar á 22-25 föstum stöðum sem teknar eru tvisvar á ári á Langjökli yfir 19 ára tímabil, auk hæðarkorts sem hefur 100 metra upplausn. Aukaafurð sem kom til við smíði á stigskipta Bayesíska líkaninu, er ný mismunaaðferð til að leysa tölulega hlutafleiðujöfnuna fyrir þunnjökla nálgunina. Þó svo að aðferðirnar sem hér eru settar fram séu fyrir jöklafræði þá má útfæra þær fyrir önnur jarðeðlisfræðileg gögn og tilsvarandi líkön. Ritgerðin byggir á þremur vísindagreinum. Í fyrstu greininni er farið yfir þau líkön sem hafa verið þróuð til að lýsa hreyfingu jökla, annarar gráðu mismunaaðferð til að leysa tölulegu hlutafleiðujöfnuna fyrir þunnjökla nálgunina er kynnt sem og stigskipt Bayesískt líkan sem notar tölulegu lausnina, og mat byggt á stigskipta Bayesíska líkaninu er borið saman við fræðilega lausn þunnjökla nálgunarinnar. Í grein tvö er fjallað um nánari útfærslu á tölfræðinni og útreikningunum fyrir stigskipta Bayesíska líkanið, slembigangur af stigi hærra en einn er skoðaður sem líkindafræðilegt líkan fyrir mismuninn á milli tölulega hermisins og raunveruleikans, sýnd notkun á fyrstu gráðu hermum, og ný reiknisparandi nálgun fyrir sennileikafallið er kynnt. Í þriðju greininni er aðferðafræði fyrstu tveggja greinanna beitt á afkomugögn og hæðargögn frá Langjökli. Framlag ritgerðarinnar felst í: (i) nýrri annarrar gráðu mismunaaðferð til að leysa tölulegu hlutafleiðujöfnuna fyrir þunnjökla nálgunina í tveimur rúmvíddum og tíma, (ii) notkun slembigangs til að lýsa mismuninum á milli tölulega hermisins og raunverulega ferlisins, (iii) nýrri reiknisparandi nálgun fyrir sennileikafallið, (iv) að beita nýju tíma- og rýmisháðu tölfræðilegu líkani fyrir jökla við greiningu gagna frá Langjökli, einum af stærstu jöklum Íslands.Icelandic Research Fun

A Bayesian hierarchical model for glacial dynamics based on the shallow ice approximation and its evaluation using analytical solutions

Bayesian hierarchical modeling can assist the study of glacial dynamics and ice flow properties. This approach will allow glaciologists to make fully probabilistic predictions for the thickness of a glacier at unobserved spatiotemporal coordinates, and it will also allow for the derivation of posterior probability distributions for key physical parameters such as ice viscosity and basal sliding. The goal of this paper is to develop a proof of concept for a Bayesian hierarchical model constructed, which uses exact analytical solutions for the shallow ice approximation (SIA) introduced by Bueler et al. (2005). A suite of test simulations utilizing these exact solutions suggests that this approach is able to adequately model numerical errors and produce useful physical parameter posterior distributions and predictions. A byproduct of the development of the Bayesian hierarchical model is the derivation of a novel finite difference method for solving the SIA partial differential equation (PDE). An additional novelty of this work is the correction of numerical errors induced through a numerical solution using a statistical model. This error-correcting process models numerical errors that accumulate forward in time and spatial variation of numerical errors between the dome, interior, and margin of a glacier.The Icelandic Research Fund (RANNIS) is thanked for funding this research.Peer Reviewe