sampling the basal friction

[nicolearetz]

Hello,

I was wondering

1. Is there currently a way in ISSM to treat the basal friction as a random field and draw samples from it? I've written some thoughts on that below to illustrate what I mean.
2. Are there any common or accepted approaches on how to specify its uncertainty? E.g., which distribution to choose (log-normal, truncated-normal, Poisson, ...)

Thank you for your help.

Best regards,
Nicole

Ad 1) These are some options I thought about, but I don't know to which extend these would currently be feasible.

A) Bayes: I think the ideal case would be to model the log-friction-coefficient $\beta$ as random field, take the surface velocity data, apply Bayes theorem, and sample from the posterior distribution. However, the involved non-linearities and the high / infinite dimension of $\beta$ are prohibitive for standard MCMC methods.
In [1], they suggest to use a local Gaussian approximation of the posterior as proposal distribution in Metropolis-Hastings MCMC, with low-rank approximation of the Hessian (see B below).

B) Gassianized posterior: One option would be to approximate the posterior with a Gaussian. The mean could be chosen by solving the deterministic inverse problem (which I know is implemented in ISSM). The question is how to choose the covariance. In [2] they choose a low-rank approximation of the Hessian around the solution of the inverse problem. But that requires evaluating the action of the Hessian or of the prior-preconditioned Hessian. Is there a way of doing that?

C) Partitioning: It's my understanding that in the QMU tutorial, the domain is partitioned into small regions on which the basal friction is assumed constant. Is there a way to choose the mean and covariance for each partition based on the surface velocity data or the inverse problem?

[1] Petra, Noemi, James Martin, Georg Stadler, and Omar Ghattas. "A computational framework for infinite-dimensional Bayesian inverse problems, Part II: Stochastic Newton MCMC with application to ice sheet flow inverse problems." SIAM Journal on Scientific Computing 36, no. 4 (2014): A1525-A1555.
[2] Isaac, Tobin, Noemi Petra, Georg Stadler, and Omar Ghattas. "Scalable and efficient algorithms for the propagation of uncertainty from data through inference to prediction for large-scale problems, with application to flow of the Antarctic ice sheet." Journal of Computational Physics 296 (2015): 348-368.

[MathieuMorlighem]

Hi Nicole

Nice to hear from you! These are some very interesting (and important!) questions that we have not been really investigated so far. To me, these two papers from Noemi and Toby remain "state of the art" and it looks like nobody has continued their work. Dakota (through UQ) can do nonuniform distributions but it will not have a spatial covariance, which is probably what you need here. If you have ideas, don't hesitate to contact us. We are happy to implement new features to make this sort of studies possible.
Mathieu

[nicolearetz]

Hi Mathieu,

Thank you for the encouraging reply. I'm glad I didn't overlook any features.

Being able to do some sort of Gaussian approximation around a chosen mean (e.g., the solution of the inverse problem) would already be very helpful. If I see it correctly, that requires being able to evaluate the Hessian action and specifying a spatial covariance. I imagine the Hessian action is already implemented for solving the inverse problem. The covariance could be specified through a Laplacian-like differential operator (e.g., [3]), that's at least the go-to in my bubble of the world. That requires some knowledge on how to implement PDEs in the ISSM source code but should be straight-forward -- at least I was able to do it very quickly in FEniCS after exporting an ISSM mesh.

I'll think about whether there are any other ingredients. I'd be happy to help as well.

Best regards,
Nicole

[3] Bui-Thanh, Tan, Omar Ghattas, James Martin, and Georg Stadler. "A computational framework for infinite-dimensional Bayesian inverse problems Part I: The linearized case, with application to global seismic inversion." SIAM Journal on Scientific Computing 35, no. 6 (2013): A2494-A2523.