Please see the attached PDF for the review.

(note my quality/impact/... ratings on the review form are based on the current version)

The Pollard et al submission is an attempt of applying history

matching to the Penultimate Glacial Maximum (PGM) for Eurasia. The

choice of history matching is appropriate, however the implementation

is limited and currently inadequate to match the claims. Critical, the

submissions claims via the title to quantify "the Uncertainty in the

Eurasian Ice-Sheet Geometry at the Penultimate Glacial Maximum", and

via the abstract to "robustly quantify uncertainties" and yet it only

partially does so. It fails to account for the structural uncertainty

of their static perfectly plastic ice sheet model. It also fails to

account for uncertainties in the MIS 6 ice margin. A key point is that

MIS 6 maximum ice extent has poor age control. It is unclear to what

extent parts of the margin represent short term surge events (which

are not going to be well represented by a perfectly plastic ice sheet

model), nor is it clear to what extent the maximum ice margin extent

was synchronous. Nor is the potentially large uncertainty associated

with the assumption of an equilibrium ice sheet addressed.

The work is of potential value, but it first needs to make claims that

are defensible. This includes clarity and accuracy on the extent to

which uncertainties are addressed and that this is a very approximate

history matching as the parameter sampling is far from

complete. History matching typically relies on emulators to adequately

sample the parameter space. 200 samples for 7 parameters is far from

adequate unless the response is very linear with minimal interaction

between parameters (which would have to be shown).

Secondly, the chosen model uncertainty estimate has no

justification. Futhermore, it is clear that there is a bias error that

also needs to be addressed in the implausibility function given the

mean NROY misfits to the GLAC1D and ICE6G ice sheet chronologies.

However, this can be rectified. Select at least 20 of your simulations

that have the least RMSE error for GLAC1D (and separately for ICE6G if

that is fully from a glaciological model).  Use the variance of the

residuals as your minimum structural model variance error estimate and

use the mean bias as your minimum model bias error estimate in the

implausibility. Note, these values will clearly vary around the ice

sheet. You can either make the error estimates a field (ie depending

on relative location in the ice sheet), or you can choose the maximum

value across the ice sheet (easier but at cost of wider

uncertainties).  As both of these ice chronologies have their own

limitations, expanding the resultant variance and bias estimates by

some fudge factor, would still be needed. To account for errors

growing with a larger PGM ice sheet, scale the uncertainties for ice

thickness by the ratio of mean ice height.

Thirdly, the authors are conflating NROY with plausible. This is a

problematic stretch. Showing something as being NROY, ie not

implausible, doesn't necessarily make it plausible.

Finally, the authors do not adequately address the limitations of

the perfectly plastic approximation and make a number of inaccurate claims,

some of which are detailed below (running late on this review, I figured

it's better to give you something to work with now than a completely detailed

evaluation).

# some specific comments

Quantifying the Uncertainty in the Eurasian Ice-Sheet Geometry at

the Penultimate Glacial Maximum (Marine Isotope Stage 6)

# The title is misleading, as uncertainties are inadequately quantified.

51.16±6.13 m sea-level equivalent

# The significant digits given for results are meaningless, correct this to

# an appropriate amount.

We perform Bayesian uncertainty quantification

# History matching is not Bayesian (where do you invoke Bayes Rule?)

# cf https://www.physics.mun.ca/~lev/revCalG.pdf

# for an explanation of what Bayesian is and entails.

Finally, the simple ice-sheet model approach is designed to generate

ice geometries based on simple, steady-state ice-sheet physics for a

prescribed margin

# Misleading as stated. The perfectly plastic ice sheet model is

# is derived from "steady-state ice-sheet physics" but it doesn't

# reflect it, only provides a limited approximation. 

simple ice-sheet model whose minimal input requirements

# The specification of a 2D basal shear stress map for each timeslice

# is far from  minimal. 

However, reliance on poorly constrained rebound data required for GIA

inversion modelling (Lambeck et al., 2006) or assumptions of highly

uncertain climate data used in dynamic ice-sheet simulations

(Abe-Ouchi et al., 2007; Peyaud, 2006) make these approaches

challenging to constrain for the Penultimate Deglaciation and give

only a very limited view of possible pasts with no grasp on the vast

range of plausibility

# The last claim is incorrect. I'm using a glaciological (hybrid

# shallow shelf/shallow ice) model with large degrees of freedom in

# the climate forcing. Though not perfect, I suspect I've generated a

# larger range of history matched simulations for last glacial cycle Eurasian ice

# sheet evolution (in early process of write-up) than is offered by

# the methodology of this submission, as evidenced by my current 50

# ensemble parameters (versus 7 in this submission).

a simple ice-sheet model whose minimal input requirements enables the

production of large ensemble simulations with controlled sources of

uncertainty (Gowan et al., 2016a).

# If the sources of uncertainty are "controlled", then they should be

# fully assessed.

generate physically plausible ice-sheet reconstructions

# Given the approximations involved, I don't see how you can call

# these physically plausible. 

since our simulations

505 are process based

# The perfectly plastic approximation is not a process based model but

# a diagnostic tool. What processes are you actually modelling?

The model has been successfully applied where large uncertainty in

inputs required for dynamic ice-sheet models, such as climate, have

reduced the confidence in using the outputs of such models as inputs

to sea-level models due to misfits against ice extent and volume

distributions that impact GIA, and where large numbers of runs are

required making computation efficiency 105 paramount, such as in the

exploration of variable global ice-sheet configurations (Gowan et

al., 2021)

# What do you mean by "successfully"? Instead of vague descriptors,

# be precise or do you want to be stuck with my definition of "successful"?

Limited constraints on climatic conditions, the requirement for large

ensemble simulations to explore the range of plausible scenarios, and

a need for well-defined sources of uncertainty make ICESHEET an ideal

choice for exploring uncertainty in ice sheet configurations during

the PGM.

# Again, given the statement above, a full uncertainty assessment

# should be provided to match the claim. Without seeing that assessment,

# I see no basis for calling ICESHEET an "ideal" choice, or even a defensible choice.

In order to account for GIA, we assume that the Eurasian ice sheet at

the PGM had been at its maximum extent sufficiently long for the solid

Earth underneath to be at (or close 180 to) an isostatic equilibrium

with the ice load, an assumption we consider reasonable given the lack

of constraints during this time

# this is a large assumption, contributing another source of unassessed uncertainty

# One main reason for your lack of constraints, is your missing of full glaciological

# physics of the ice and earth system which provide some pretty strong constraints

# on the system. 

We run 200 simulations for each reconstruction and each of the 4

selected time periods (22, 20, 18, 16 ka), totalling 1600 simulations

(Figure 5 and Figure A1).

# so you are assuming minimal interactions between your parameters,

# but do not provide evidence to support this. If this were not the

# case, then even just a 3 value min/median/max grid search over 7

# parameters would entail 3^7= 6561 simulations for 1 timeslice.

eq 1

# The denominator for implausibility should include variance for model

# structural error. The numerator should include model structural bias

# error. You are choosing to specify model structure error as simply

# some fraction of ensemble variance. On what basis do you justify

# such a choice?  To understand why this can be problematic, cf the

# simple example in subsection 2.5 of

# https://doi.org/10.5194/egusphere-2022-1410

# The above reference also provides guidance some guidance on how model

# structural uncertainty can be defensibly assessed.

Bayesian History Matching

# Looking at the first page of googled hits for "Bayesian History Matching", the actual papers

# describe Bayesian emulation used in history matching. As you are not using emulators, your

# description is inaccurate. This needs to be corrected throughout.

We restrict our NROY space to parameter values that correspond to

models runs with implausibility I(ˆp) higher than 3, following the

Pukelsheim (2012) three-sigma rule typically used in Bayesian History

Matching

# Should be less than 3 not "greater than 3" for NROY. You should also

# state clearly what assumptions you are making about the residual distribution

# to justify the choice 3 sigma (the modelling community relies too much on "this

# is what others do"). cf the above https link for the assumptions made.

Figure 5

# please add a few contours to each frame as its hard to discern the

# colour map values within even 300 m Also, I can't make sense of the

# colour scale for the first two columns. How can you have negative

# ice thickness or ice surface elevation (not clear what is being

# plotted) a km below sealevel?

Our work has expanded this methodology to include the cold-based ice

and active ice streaming basal processes which have had a strong

impact on the implausibility metric, improving the simulation fit

during history matching when applied to the Last Deglaciation, with

the exception of the British ice sheet (Figure 4) where simulation

mismatch is likely due to discrepancies in ice-margin extraction.

# You need to be honest about exactly what the perfectly plastic

# approximation entails and how far what you implement is from

# the physics of ice streams.