Two sources of uncertainty in estimating tephra volumes from isopachs: perspectives and quantification

Calculating the tephra volume is important for estimating eruption intensity and magnitude. Traditionally, tephra volumes are estimated by integrating the area under curves fit to the square root of isopach areas. In this work, we study two sources of uncertainty in estimating tephra volumes based on isopachs. The first is model uncertainty. It occurs because no fitted curves perfectly describe the tephra thinning pattern, and the fitting is done based on log-transformed square root of isopach area. The second source of uncertainty occurs because thickness must be extrapolated beyond the available data, which makes it impossible to validate the extrapolated thickness. We demonstrate the importance of the two sources of uncertainty on a theoretical level. We use six isopach datasets with different characteristics to demonstrate their presence and the effect they could have on volume estimation. Measures to better represent the uncertainty are proposed and tested. For the model uncertainty, we propose (i) a better-informed and stricter way to report and evaluate goodness-of-fit, and (ii) that uncertainty estimations be based on the envelope defined by different well-fitted curves, rather than volumes estimated from individual curves. For the second source of uncertainty, we support reporting separately the volume portions that are interpolated and extrapolated, and we propose to test how sensitive the total volume is to variability in the extrapolated volume. The two sources of uncertainty should not be ignored as they could introduce additional bias and uncertainty in the volume estimate.


Introduction
Calculating the volume of tephra fall deposits is important for the study of explosive volcanic eruptions (Pyle 1989; Communicated by A. Tadini. This paper constitutes part of a topical collection: Uncertainty Quantification in Volcanology: Observations, Numerical Modelling and Hazard/Risk Assessment. Fierstein and Nathenson 1992). Its value is critical to the eruption VEI or magnitude assessment (Newhall et al. 1982;Pyle 2015), is tied to physical processes of eruptions, and helps to constrain other eruption source parameters (e.g., Legros 2000;Mastin et al. 2009;Yang et al. 2020). Conventionally, tephra volume is estimated based on isopachs using the method proposed by Pyle (1989). New statistical and engineering methods have been proposed to construct isopachs, generate tephra thickness distributions and estimate tephra volumes (e.g., Engwell et al. 2015;Green et al. 2016;Yang et al. 2016;White et al. 2017;Yang et al. 2020Constantinescu et al. 2022;Rougier et al. 2022) with the motivations of quantifying the uncertainty, reducing the subjectivity, and/or increasing the efficiency in the calculation. Regardless of the method used, however, interpolation and extrapolation, which construct a smooth surface or curve based on discrete points, are inevitable in estimating tephra volumes.
Despite the care, strictness, and robustness in studies that attempt to estimate tephra volumes, it is widely accepted that we should still interpret the estimated volume and quan-tified uncertainty with caution, because of the difficulty in objectively quantifying different sources of uncertainty (e.g., assigning appropriate value to quantify the measurement uncertainty). Despite the development of novel methods that do not use isopachs to estimate tephra volume (e.g., Engwell et al. 2013;Yang et al. 2016;Rougier et al. 2022), the isopach-based method is still being widely used to estimate or validate new methods for tephra volume (e.g., Pyle 2016;Yang et al. 2019;Buckland et al. 2020;Williams et al. 2020;Prival et al. 2020;Takarada et al. 2020;.

Motivation
There are different ways to examine the isopachs of tephra deposits on a 2D biplot, such as the log(area)−log(thickness) and log(thickness)− √ isopach area plots (Walker 1980;Pyle 1989). With the latter, the hypothesis that tephra thins exponentially away from the vent is tested, and shown to work well for many tephra isopach datasets by Pyle (1989).
The tephra volume is estimated in the following way with the help of the log(thickness) − √ isopach area plot: the isopach data are plotted on the log(thickness)− √ isopach area plot, and then fitted to different curves. The tephra volume is estimated through the integration of these fitted curves (Pyle 1989;Fierstein and Nathenson 1992; more details introduced below). In addition to the (segmented) exponential curve proposed in Pyle (1989), power-law (Bonadonna and Houghton 2005) and Weibull (Bonadonna and Costa 2012) curves have been proposed and used for tephra volume estimation.
Methods to quantify the tephra volume uncertainty using the isopach-based method have been proposed and widely used (e.g., Daggitt et al. 2014;Biass et al. 2013Biass et al. , 2019, but there is still one source of uncertainty and one question left unresolved (which is also considered generally as a source of uncertainty) with the method. The first source of uncertainty (model uncertainty) occurs due to the fact that there is not a curve that always fits perfectly well to the data, which introduces model uncertainty to the volume estimate. The misfit may seem trivial, but it leads to a paradox: we use the fitted curves that are partially inconsistent with the isopachs to calculate the tephra volume, while the isopachs are constructed with caution and experts' knowledge, which might be (locally) more reliable than the fitted curves in reflecting the tephra thickness. This uncertainty is often omitted or considered compensated or overridden by other sources of uncertainty such as the uncertainty in constructing isopachs, while the uncertainty in constructing isopachs should be categorized as data uncertainty and thus distinguished from the model uncertainty.
The question left unanswered, namely the second source of uncertainty considered in this work, has been proposed in Klawonn et al. (2014a), which suggested that the uncertainties for tephra volumes associated with interpolated thickness (i.e., those between isopachs) and extrapolated thickness (i.e., those within the thickest and outside of the thinnest isopach) differ and should be treated separately. For the extrapolated volume that is within the thickest isopach, Pyle (2016) analyzed the volume distribution of an exponential cone (which can be used to approximate the thickness distribution of tephra deposits), and found that a large portion of the volume of the cone was enclosed within a small distance from the vent. Sulpizio (2005) has analyzed the difficulty and reliability of estimating distal tephra volume (e.g., the extrapolated volume that is outside the thinnest isopach) given different types of known information about a tephra deposit (e.g., proximal isopachs and sparse distal observations and proximal isopachs and one distal isopach). The works of Pyle (2016) and Sulpizio (2005) and the references therein provide critical insights on the theoretical and practical difficulties (and also possibility) of estimating the extrapolated tephra volumes. We should be cautious with the extrapolated volume and its associated uncertainty, but the most appropriate way to view the challenge posed by quantifying the extrapolated volume remains unclear.
Moreover, the importance of the above two sources of uncertainty might be greater than how they appear on the log(thickness) − √ isopach area plot: neither of the y-and xaxes scale linearly with the volume (as they are transformed from the original isopach thickness and area). A small thickness or isopach area difference on the plot may correspond to a much greater difference in the volume. The two sources of uncertainty will be introduced with greater detail later in the text.

Hierarchical propagation of uncertainty in tephra volume estimation
The model uncertainty and uncertainty from extrapolation listed above should not be considered unimportant when other sources of uncertainty (e.g., thickness measurement error) are present and considered in estimating tephra volumes when the isopach-based method is used. This is because (1) different factors and processes contribute to and propagate the uncertainty in tephra volume estimation following a hierarchical order; (2) the model uncertainty and the uncertainty from extrapolation are introduced in the curve fitting process, which is before integration that leads to the final volume estimation. They are thus the last sources of uncertainty in the hierarchy; (3) their inevitable presence (as stated above The volume difference could be larger than it appears here as the x-and y-axes do not scale linearly with the volume Processes affecting tephra thickness b c d e Fig. 1 a Sketch showing the difference between uncertainty being propagated (yellow polygon) following the hierarchy and uncertainty being "summed" (i.e., taking the union of different uncertainties; pink polygon). The model uncertainty and uncertainty from extrapolation could further propagate other sources of uncertainty (yellow polygon), but if neglected, they would lead to additional bias and uncertainty (pink polygon in a); b-e example showing if we neglect the model uncertainty, and only focus on the uncertainty in constructing isopachs (and excluding other sources of uncertainty), the estimated tephra volume range or uncertainty does not cover the true tephra volume (golden area in e).
Note that in d, to account for the uncertainty in constructing the thickest isopach, the three fitted curves are generated, but they all underestimate the tephra thickness at the thickest isopach level and introduced more later) indicates that they could distort the uncertainty inherited from previous steps in the hierarchy if neglected. The factors and processes that contribute to tephra volume uncertainty are listed below (Fig. 1a), and it is important to note that they occur following a hierarchical order, and the hierarchical order cannot be switched: before the thickness was measured in the field, processes such as reworking, compaction, effects of slope and vegetation would modify the primary tephra thickness (e.g., Engwell et al. 2013;Blong et al. 2017;Cutler et al. 2018); The measured thickness during field work is subject to measurement uncertainty (Engwell et al. 2013;Kawabata et al. 2013;Green et al. 2016); To use the isopach-based method to estimate tephra volume, the construction of isopachs also introduces additional uncertainty after collecting the thickness observations (Engwell et al. 2013;Klawonn et al. 2014b); At last, as mentioned earlier, the model uncertainty and the uncertainty from extrapolation are introduced right before the volume estimation.
The above demonstration suggests that even if the impact of one source of uncertainty is much greater than the others, rigorously speaking, this greater uncertainty will still propagate and get augmented by the subsequent uncertainty factors following the hierarchy (Fig. 1a), and it is thus theoretically inaccurate to assume that the final uncertainty is the sum or union of the individual uncertainties. We are aware that it is practically impossible to quantify the uncertainty in this way, but this perspective is important for the two uncertainties studied in this work as they are the last in the hierarchy. Neglecting them could dampen, exaggerate, or distort the uncertainty inherited from previous factors and processes in the hierarchy, and hence introduces additional, irrelevant bias and uncertainty to the final volume estimate. A sketch in Fig. 1b-e shows such a situation in which the estimated volume range does not cover the true tephra volume when the model uncertainty is neglected.
The above arguments demonstrate the importance of the two sources of uncertainty on a theoretical level, but in practice, it is possible that their impacts on the volume estimate are relatively small. Whether it is the case depends on the specific isopach dataset and perspectives to view and interpret the two sources of uncertainty, which will be illustrated in the following text.
In this work, we study and demonstrate the presence of the two sources of uncertainty when the isopach-based method is used for tephra volume estimation and uncertainty quantification, and propose measures and perspectives to address and interpret them. We do not consider the situation when the potential uncertainty is so great (e.g., due to limited or poorly constrained information about a tephra deposit) that it is of no use to quantify the volume uncertainty (e.g., the volume estimate is only used to determine an order of magnitude). We first introduce the log(thickness) − √ isopach area plot method, datasets, and parameters used in this work. Then we demonstrate the presence of the two, and propose measures to address them. The main contributions of this work include: (1) explicitly pointing out the two sources of uncertainty and their importance and (2) the proposed measures to address and interpret the two sources of uncertainty. They present a more accurate and more logically consistent way to capture and interpret the two sources of uncertainty when the isopach-based method is used.
Throughout this work, we assume that a set of isopachs, rather than individual thickness measurements, is available for each tephra deposit to be analyzed. We assume that all isopachs are uncertainty-free unless otherwise specified, namely they represent the true isopach areas of the studied deposits at the corresponding thicknesses. In other words, all misfit between the isopach data and fitted curves belongs to the model uncertainty that is of interest in this work. With these assumptions, we neglect uncertainties such as those from measuring tephra thickness in the field, the effect of post-eruption weathering, movement or compaction, and from constructing isopachs, and assume that they can be analyzed in other stages of the uncertainty quantification. These assumptions are necessary as otherwise we cannot exclude the impact of other uncertainties in our analysis. We also do not attempt to use our knowledge on the physics of tephra transport to constrain the thinning pattern of tephra deposits as is done in some previous works (e.g., Carey and Sparks 1986;Bursik et al. 1992;Bonadonna et al. 1998;Koyaguchi et al. 2001). This is independent from the two sources of uncertainty, and represents another uncertainty source that should be included in other stages of the uncertainty quantification. In this work, for a volume to be estimated (e.g., total tephra volume or the volume of a subset of a tephra deposit), we quantify the uncertainty as the difference between the maximum and minimum estimates, referred to here as volume variability.

Background
The isopach-based method Pyle (1989) proposed to use the one-or multi-segment exponential curve to describe how tephra thins with distance from the vent on the log(thickness)− √ isopach area plot. This was later reviewed, proved to be appropriate, and modified (to a more general formula that is not constrained by the shape of the isopachs) by Fierstein and Nathenson (1992). The volume (V ) of tephra deposits can be calculated as (Pyle (1989); Fierstein and Nathenson (1992)): where T is tephra thickness, and A is the isopach area. The tephra volume can be calculated using Eq. 1 if a relationship exists between T and A 1/2 . The relationship can be determined by fitting certain curves to the isopach data. Three types of curves have been proposed so far for the fitting. The first one is the one-segment or multi-segment exponential functions proposed in Pyle (1989) and examined by Fierstein and Nathenson (1992). The one-segment form is written as: where T 0 is the extrapolated thickness when A = 0, and −k is the slope of the line on the log(thickness) − √ isopach area plot. The multi-segment form, as its name suggests, is composed of two or more connected functions each defined by Eq. 2 with different intercepts and slopes on the log(thickness) − √ isopach area plot. In particular, the two-segment form is defined by two extrapolated thicknesses (T 0 and T 1 ), two slopes on the log(thickness)− √ isopach area plot defining the different thinning rates (k and k 1 ), and a value defining where the two intersect (A 1/2 i p ). Three empirical methods established based on the two-segment exponential curve have been proposed and tested to esti-  Fig. 2 a Tephra volume is discretized as "vertical bars" and calculated with Eq. 5 using the trapezoidal rule; b how V prox , V int , and V dist are defined; c how Regions A, B, and C are defined in Klawonn et al. (2014a). The shaded area is the difference between V prox and Region C and also the difference between V int and Region A; d and e with two fitted curves that both fit well (based on a specified criterion) to the isopach data, instead of using the volumes from the two individual curves (as volumes of the solid and dashed lines in d), we propose to use the union or the envelop thickness defined by the two curves to constrain the volume from the curves (as volumes of the green and yellow lines in e); f sketch showing how V prox , V dist , R prox , and R dist are defined mate the distal volume (similar to the extrapolated volume outside the thinnest isopach) of tephra fall deposits by Sulpizio (2005). The power-law relationship is proposed in Bonadonna and Houghton (2005), and can be written as: where T pl is a constant, and m is the power-law coefficient. Bonadonna and Costa (2012) used Weibull function to describe the relationship between T and A 1/2 : where λ is a characteristic decay length scale denoting deposit thinning, θ corresponds to a thickness scale, and n is a shape parameter. By substituting Eqs. 2-4 and the segmented form of Eq. 2 each to Eq. 1, the total volume of tephra deposit can be integrated analytically. The corresponding parameters of the curves ( Eqs. 2-4) are derived through curve fitting based on hand-drawn isopachs or isopachs from interpolation techniques. For the power-law function (Eq. 3), (both proximal and distal) limits need to be specified during integration to prevent the total volume from going to infinity (Bonadonna and Houghton 2005). It is noted by Fierstein and Nathenson (1992) that we could change the limits in Eq. 1 from 0 to infinity to any other √ isopach area ranges, say isopachs A and B with areas A a and A b and thicknesses T a and T b . If we assume that tephra thins exponentially between √ A a and √ A b , the tephra volume (V a−b ) between the two isopachs is written as (the equation below is from Eq. 13 of Fierstein and Nathenson 1992): where k a−b and log(T 0a−b ) are the slope and intercept of the line defined by the two isopachs on the log(thickness) − √ isopach area plot. With the help of the trapezoidal rule, we could use Eq. 5 to calculate the volume within any pair of isopachs for a certain fitted curve (Fig. 2a).

(N-1)-segment exponential curve
Another form of the multi-segment exponential function needs to be introduced as it will assist the analysis in this work. It simply extends the one or two exponential segments to (N-1) segments on the log(thickness) − √ isopach area plot. Here N refers to the number of isopachs for the deposit. It has been adopted in previous studies before (Fierstein and Nathenson 1992), but is not widely adopted in more recent studies. We refer to it as the (N − 1)-segment exponential curve in this work. Each segment of this curve is defined by the straight line connected by a pair of neighboring isopach data on the log(thickness) − √ isopach area plot (Fierstein and Nathenson 1992). For ranges from zero to the √ thickest isopach area and from the √ thinnest isopach area to infinity, the curve can be simply defined by extending the first and last segments. The (N-1)-segment exponential curve is unique because it is defined by the isopach data, and it does not require additional curve-fitting procedure.

Datasets
We use isopach datasets from six well-studied tephra deposits as examples in this work. The deposits are the tephras from the 1815 Tambora eruption (Kandlbauer and Sparks 2014), Taupo Pumice Fall (Walker 1980), Cotapaxi Layer 5 (Biass and Bonadonna 2011;Biass et al. 2019), Hatepe tephra (Walker 1981;Fierstein and Nathenson 1992), Minoan tephra (Pyle 1990;Daggitt et al. 2014), and tephra from the 1980 Mt. St. Helens eruption (Sarna-Wojcicki et al. 1981;Fierstein and Nathenson 1992). They are chosen as their volumes and extents are different in magnitude, and the number of their isopachs varies from three to twelve. They represent tephras from island settings, where tephra measurements are typically limited where it falls in the sea (e.g., Taupo tephra) and where access to the volcanic island is limited (e.g., tephra from the 1815 Tambora eruption), and in continental settings, where a more complete distribution of isopachs are available (e.g., the 1980 Mt. St. Helens tephra). The selected tephra deposits also differ in the quality, uncertainty, or reliability of their thickness measurements and isopachs. For example, older tephra (e.g., Taupo Plinian deposit) might be incompletely preserved, poorly accessible, while tephra from well-observed or very recent eruption (e.g., 1980 Mt. St. Helens tephra) is likely to be better preserved, the extent better constrained, and with more isopachs "mappable". Under the assumption in this work (i.e., zero uncertainty in the isopachs), this would not affect any analysis in this work, but we stress that in reality such differences in the deposits and the corresponding datasets could lead to different expectations of what might be possible to reduce the uncertainty. We consider the selected datasets covering characteristics of most isopach datasets. More information about these datasets is given in Table 1, and the isopach data are presented in Tables 5 and 6. [1] Kandlbauer and Sparks, 2014;[2] All data from Table 1 of Bonadonna and Costa, 2012. Volume estimates from exponential and power-law curves have slight deviations compared to data listed in Fierstein and Nathenson (1992) and Bonadonna and Houghton (2005). We select data in Table 1 of Bonadonna and Costa (2012) for consistency, and they would not affect any conclusions in this work  St. Helens (12) 0.0126 Bonadonna and Costa (2012).
[1] linear regression based on the log-scaled thickness and square root of isopach area; [2] estimated using Excel spreadsheet from Bonadonna and Costa (2012)

Parameters
We apply the one-, two-, and (N-1)-segment exponential, power-law, and Weibull curves to the six isopach datasets. Parameters of the fitted curves are given in Table 2, and the fitted curves are shown in Fig. 3. Some of them are directly referenced from previous studies, and the others are updated in this work. For the 1815 Tambora tephra which only has three isopachs, the parameters of its fitted curves are from Kandlbauer and Sparks (2014). Parameters of fitted curves for the Cotapaxi Layer 5 are from Biass and Bonadonna (2011) and Bonadonna and Costa (2012). Parameters of the oneand two-segment exponential curves for the rest of the datasets are from Fierstein and Nathenson (1992). The m and T pl in the power-law curve (Eq. 3) are updated for the Hatepe dataset. The power-law curve does not fit well to the Taupo Pumice Fall, Minoan, and 1980 Mt. St. Helens tephra datasets, and is thus not adopted here. Parameters of the Weibull curves for the Taupo Pumice Fall, Hatepe, and Minoan deposits are updated to best fit the data using the Excel Spreadsheet provided in Bonadonna and Costa (2012). Parameters of the Weibull curve for the 1980 Mt. St. Helens tephra are referenced from Bonadonna and Costa (2012). We note here that a high standard is adopted to determine whether a curve fits the isopach data well: the curves need to fit the data visually well on the log(thickness)− √ isopach area plot, and the predicted thicknesses from the curves are characterized by high correlation coefficients with the isopach data (Table 4; all greater than 0.936 in non-logged thickness). This is necessary because we are interested in the presence of the two sources of uncertainty when the fitted curves seem to perform well. In other words, there is no point studying the curves if they obviously fit poorly to the data.

Volume calculation and distribution with respect to isopach area
We calculate the total volume (Table 3) and volume distribution with respect to √ isopach area based on the fitted curves using the trapezoidal rule (Fig. 2a). The fitted curves are discretized to 0.5 or 1 km-length √ isopach area segments, and the volume within each segment is calculated using Eq. 5. The maximum √ isopach area for each deposit in this calculation, i.e., the upper integration limit for the volume calculation, is given in Table 1, which is at least ∼4 times greater than the √ thinnest isopach area of the respective datasets. Different integration limits are specified for the power-law and Weibull (only for the 1815 Tambora tephra as the corresponding Weibull curve leads to significantly great thickness) curves to be consistent with previous works and to avoid unrealistically large thickness values. Different proximal integration limits (e.g., 1, 2, 4 km √ isopach area) are specified for the power-law curves of the 1815 Tambora tephra, Cotopaxi Layer 5, and Hatepe tephra datasets to account for its impact on the volume estimation. Proximal integration limits of the 1815 Tambora tephra are specified based on different thicknesses (thickness below 200, 150, and 120 cm), rather than by √ isopach area, to be consistent with the work of Kandlbauer and Sparks (2014). The same proximal integration limit is applied to the fitted Weibull curve of this deposit. Distal integration limits for the power-law curves are specified from previous works as 300 km for the Cotapaxi Layer 5 (Biass and Bonadonna 2011), 1000 km for the Hatepe tephra (Bonadonna and Houghton 2005), and 1500 km for the 1815 Tambora tephra (Kandlbauer and Sparks 2014).
All volume estimates from this work are consistent with previous estimates calculated from analytical integration (Table 3). The volume distributions with respect to √ isopach area based on the fitted curves are given in Fig. 4 for each deposit.

The model uncertainty
Misfit occurs because deviation between the isopach data and the fitted curves inevitably exists. This is due to uncertainty in both the isopach data and model. In the context of this work, where the former is temporarily assumed to be zero, the latter becomes more apparent (which does not mean that it is not significant otherwise): we cautiously construct isopachs such that they represent the true thicknesses and isopach areas of the deposit at the corresponding levels, but then we use the fitted curves that inevitably deviate or even fit poorly (see paragraph below) to the isopach data to estimate the tephra volume.
Practically, the curve fitting process is done based on log-transformed thickness and square root of isopach area. A curve that seems to fit well on the log(thickness) − √ isopach area plot does not necessarily fit well to the original isopach data (e.g., a predicted √ isopach area that is 120% of the original √ isopach area corresponds to 144% of the original isopach area). This could be further exaggerated by the maximum thickness (at most hundreds of meters) to extent (a few to even more than a million square kilometers) ratio of tephra fall deposits being extremely low. A thickness difference of 0.5 cm between the isopach data and the fitted curve may seem small, but the small difference might span an extensive area, which could greatly affect the volume estimation. For the power-law curve, different proximal integration limits are used to test whether they are sensitive to the total volume. The Weibull curve for the Tambora tephra also adopts a proximal integration limit to avoid unrealistically large thickness and volume prediction [1] Integration limits to be consistent with Kandlbauer and Sparks, 2014 [2] Different from values listed in Bonadonna andCosta, 2012, but consistent with Fierstein andNathenson, 1992 [3] Different from value estimated in Bonadonna and Costa, 2012, but based on a better fitted curve Fig. 3 The log(thickness) − √ isopach area plot showing the original isopach data and the fitted curves with their parameters given in Table 2 One-segment exponential (N-1)-segment exponential

Problem demonstration
We examine the fitting between the fitted curves and the isopach data to demonstrate that the misfit as a result of model uncertainty occurs commonly and could affect the tephra volume estimation. The high correlation coefficients (0.936-> 0.999; Table 4) between the non-logged thickness predictions of the fitted curves and the original data seem to suggest that all fitted curves are consistent with the isopach data.
The original isopach data and predictions from the fitted curves are compared in Tables 5 and 6. Excluding the one-segment exponential curve applied to Cotapaxi Layer 5 (because it can be better described by the two-segment exponential curve), out of the 98 predictions from different curves and for different deposits, there are 35 thickness predictions that are outside the 90-110% range with respect to non-logged thicknesses of the original datasets. Examined based on the isopach area, there are 44 predictions that are outside the 90-110% range. Among these outlier predictions, Ratios of prediction to original data are also given. Predictions that are 80-90% or 110-120% with respect to the original data are marked in blue, and the ones that are below 80% or above 120% with respect to the isopach data are marked in red there are eight thickness and twenty isopach area predictions that are outside the 80-120% range with respect to the original datasets with the maximum deviation of 55% (20-cm isopach for the Mt. St. Helens tephra with the Weibull curve) and 142% (200-cm isopach for the Hatepe tephra with the power-law curve) in terms of thickness, and 16% (20-cm isopach for the Mt. St. Helens tephra with the Weibull curve) and 167% (200-cm isopach for the Minoan tephra with the Weibull curve) in terms of isopach area. These results suggest that the deviation between the isopach data and predictions from the fitted curves occurs commonly even when the fitted curves seem to perform well (e.g., as indicated by the high correlation coefficients), and hence prove the common presence of misfitted curves. The above results also show that for the same datasets, more predictions with greater misfit can be detected when examined based on the isopach area, showing practically the importance of examining the goodness-of-fit based on area rather than thickness.
The tephra volume range defined by different fitted curves for each deposit is shown in Table 7. The ratio of maximum and minimum volume difference to the maximum volume is also given. The volume estimates that require integration limits (i.e., the power-law or Weibull curve for the Tambora tephra) are not included here because their values could greatly affect the volume estimate (Fig. 6), making the corresponding results non-comparable. Except for the Taupo Plinian deposit (0.6%), the ratios range from 10.5 to 20.7%, and the low ratio for the Taupo deposit will be discussed later in the text. Here we focus on the 1815 Tambora and 1980 Mt. St. Helens tephra datasets, as they are the end members with respect to the volume difference ratio (20.7 and 10.5%), and also because the one-segment exponential (for the Tambora tephra) and Weibull (for the Mt. St. Helens tephra) curves fit the deposits relatively poorly (correlation coefficients of 0.967 and 0.936, respectively). For the 1815 Tambora eruption dataset with just three isopachs, the one-segment exponential curve greatly underestimates the thickest isopach and greatly overestimates the second thickest one (curve predictions: 73.9% and 142.4% of the thicknesses, and 53.9% and 141.4% of the isopach areas Table 6 Original isopach data and isopach areas predicted by different fitted curves for the tephra datasets Ratios of prediction to original data are also given. Predictions that are 80-90% or 110-120% with respect to the original data are marked in blue, and the ones that are below 80% or above 120% with respect to the isopach data are marked in red respectively). As a result, the volume estimate based on the one-segment exponential curve, 103.4 km 3 (Table 2), should not be used to characterize the deposit volume.
The 1980 Mt. St. Helens dataset provides an interesting case in which, depending on the perspective, the impact of the model uncertainty could be considered negligible or noteworthy. The 1980 Mt. St. Helens tephra dataset has 12 isopachs. The thickness and area ranges of the 12 isopachs are 20 to 0.05 cm and 200 to 167,000 km 2 , respectively. The two-and (N-1)-segment exponential curves are highly consistent with the isopach data (correlation coefficient greater than 0.999) and with each other, leading to volume estimates of 1.13 and 1.14 km 3 (Fig. 5). The volume estimate from the Weibull curve (correlation coefficient: 0.936) is 1.02 km 3 . The volume difference of −0.11(=1.02−1.13) km 3 seems small, but closer examination in Fig. 5 shows that it occurs For the former, the volume variability divided by the maximum volume is given for reference. Volumes whose calculation requires the specification of integration limit are not included here to avoid additional complexity. For the 1980 Mt. St. Helens tephra dataset, crossed-out calculation is done including the Weibull curve. See text for more details almost entirely due to the misfit between the isopach data and the fitted Weibull curve. Figure 5c shows the volume per √ isopach area difference between the two-and (N-1)segment exponential curves, and between the Weibull and (N-1)-segment exponential curves. The absolute volume difference of the former is a lot smaller than that of the latter. For the latter pair, the largest volume difference arises from the misfit in the √ isopach area ranges of ∼ 0-24 and ∼ 60-220 km, which encloses the 20-and 10-cm, and 3-, 2-, and 1-cm isopachs, respectively. The Weibull curve estimates thicknesses of 11.03 and 8.31 cm for the 20-and 10-cm isopachs, and 2.48, 1.77, and 0.86 cm for the 3-, 2-, and 1-cm isopachs, respectively (see Table 6 for the predicted isopach area data). The misfit in the two √ isopach area ranges contributes to −0.12 km 3 of volume difference between the Weibull and (N-1)-segment exponential curves (this is greater than the total volume difference of −0.11 km 3 because for some √ isopach area, the Weibull curve thickness is greater; e.g., the red shading area in Fig. 5c). For the second range (i.e., √ isopach area from 60 to 220 km), the misfits, namely −0.52, −0.23, and −0.14 cm for the 3-, 2-, and 1-cm isopachs respectively, may seem small, but due to the low thickness to extent ratio, this difference spans a wide area, leading to a non-negligible volume difference of −0.072 km 3 , about 6.4% of the total volume.
Viewing the problem from another perspective, we can also state that the fitting of the Weibull curve is not ideal, but a volume estimate of 1.02 km 3 from it, which is not greatly different from other estimates, shows the robustness of the isopach-based method. In addition, the magnitude of other sources of uncertainty is likely to be greater, and hence overrides the uncertainty from the misfit. However, if not carefully treated and acknowledged, the model uncertainty from the Weibull curve (in this particular case) might introduce bias or augment the uncertainty in the final volume estimate. For example, if we consider the uncertainty in constructing isopachs, and quantify the 1980 Mt. St. Helens tephra volume uncertainty including the Weibull curve estimate, the volume uncertainty would partially "propagate around (as we consider the uncertainty in constructing isopachs here)" the Weibull curve estimate that is already in accurately underestimated, leading to overestimated volume uncertainty. Moreover, omitting the model uncertainty makes the tephra volume uncertainty non-comparable among different deposits, as the impact of the model uncertainty on the total volume uncertainty may vary by deposits.
In contrast to the two examples discussed above, the fitted Weibull curve for the Minoan tephra does not fit well to its isopach data, but its impact on the volume estimate is negligible. Its predictions for the 400-, 300-, and 200-cm isopachs are 84.8%, 90.7% and 119.8% of the respective isopach thicknesses and 62.3%, 75.6%, and 166.9% of the respective areas, but they are extremely small areawise rela- Fitted twoand (N-1)-segment exponential, and Weibull curves compared pairwise. Note that the y-axis shows the non-logged thickness. Selected thickness and isopach area predictions from the curves are labeled; c the volume per √ isopach area difference between the two-and (N-1)segment exponential curves (yellow line) and between the Weibull and (N-1)-segment exponential curves (dark blue line). The volumes that the shaded areas correspond to are marked tive to the other isopachs. For example, the 200-cm isopach has an area of 124 km 2 ( √ isopach area: ∼11.1 km), while the next thinner isopach, the 30-cm isopach, has an area of 21710 km 2 ( √ isopach area: ∼147.3 km). The misfit for the thicker isopachs thus has negligible impact on the total volume (Fig. 4e). This example is briefly presented here to show that the impact of the model uncertainty on tephra volume estimation also depends on properties of the specific isopach dataset (e.g., isopachs thicknesses and areas, their numbers and spacing).

Proposed measures to better quantify the model uncertainty
We propose to use isopach area, rather than correlation coefficient, thickness, or log-thickness, to examine the goodnessof-fit for the fitted curves, and all fitted isopach areas should be reported. Curves that do not fit well to the isopach data should not be used for volume calculation, and the criteria for using or not using a certain curve should be specified. We acknowledge the practical difficulty of determining what is "fitting well to the data", and refrain from drawing a hard line on it. However, one intuitive, reasonable, and bottom line criterion is provided here: if half or more of the predictions from a certain fitted curve are outside the ±40% range with respect to the corresponding isopach areas, the curve should not be used for the volume calculation. Essentially, we advocate for maximizing the clarity of how the fitted results are reported such that poorly fitted curves would not be used for the volume estimation, and that the most-informed interpretation on the volume estimates can be made.
For a set of fitted curves that pass the above or other stricter criterion, we should use the envelope or union of the thickness from these curves to define the volume range and variability, rather than using volumes calculated from individual curves. The (N-1)-segment exponential curve should also be included. The union thickness of the curves (Fig. 2d and e) is defined by the range between the maximum and minimum thicknesses among all well-fitted curves for each √ isopach area value. The bounds of the union thickness for all √ isopach area value are two curves which can be used to define the tephra volume range and variability. This can be realized with the help of Eq. 5 and discretizing the curves as done in this work.
The justification of this proposition is that if a fitted curve can be used to calculate the total tephra volume, any subset of the curve should be qualified to calculate the local tephra volume. In this way, the model uncertainty could be better captured and quantified independently from the uncertainty in constructing isopachs. Indeed, this proposition neglects the individual thinning pattern of different curves. Even though the segmented exponential, power-law, and Weibull curves are proposed with certain justifications, none of them are always better than the others. The complexity of plume dynamics and tephra dispersal and deposition suggests that each of these curves might be a good, but definitely not always perfect, approximation to the true thinning pattern, justifying the proposition.
The proposed measures do not consider other sources of uncertainty, and hence can be easily coupled with methods that quantify the other uncertainty sources in tephra volume estimation with the isopach-based method (e.g., Biass et al. 2013;Daggitt et al. 2014;Biass et al. 2019) without interrupting the uncertainty propagation.

Application
We compare volume ranges of the six deposits defined by individual curves (the (N-1)-segment exponential curve included) and by the proposed measure in Table 7. Curves that require the specification of integration limit are excluded due to their significant impact on tephra volume (Fig. 6). The onesegment exponential curve for the Tambora tephra which has been shown to fit poorly to the isopach data is excluded. This only leaves the two-segment exponential curve for the dataset (which is also the (N-1)-segment exponential curve as it only has three isopachs), and the volume range for the deposit is thus not calculated. For the Mt. St. Helens tephra, we calculate the volume variability (max-min volume) including and excluding the Weibull curve which does not fit well to the isopach data to show the effect of including ill-fitting curves. The volume variability is smaller if the Weibull curve is excluded using the proposed measure (volume variability of 0.05 and 0.15 km 3 , respectively, excluding and including the Weibull estimate), and it is even smaller than the volume variability calculated based on individual curves (0.12 km 3 ; Weibull curve included). This shows that neglecting the model uncertainty might lead to overestimated tephra volume uncertainty. For the other four datasets, their volume variabilities defined based on the proposed measure are all greater and theoretically more robust than those calculated based on the individual curves.

Uncertainty from extrapolation
The uncertainty from extrapolation is well-recognized in time series and spatial data analysis. Tephra volumes can be separated into volumes estimated based on interpolation and extrapolation. In this work, we separate the tephra volume into three portions, namely the portion that is within the thickest isopach, the portion that is within the thickest and thinnest isopachs, and the portion that is outside the thinnest isopach. Their volumes are denoted as V prox , V int , and V dist , respectively (Fig. 2b). V int is the interpolation volume, and the sum of the other two corresponds to the extrapolation volume. The definition of V prox , V int , and V dist is different from how Klawonn et al. (2014a) defined the three regions of tephra volume (sketch shown in Fig. 2c), but this would not affect any conclusions in this work: V dist defined here is equivalent to Region B in Klawonn et al. (2014a), and the difference between V prox and Region C equals to the difference between Region A and V int , which is a constant solely   Fig. 6 V prox , V int , and V dist estimated from different fitted curves for the six tephra datasets shown as histograms with integration limits given determined by the thickest isopach (shaded area in Fig. 2b and c). The interpolated thickness and volume can be examined by leave-one-out validation. However, this cannot be done for those from extrapolation. As the uncertainties associated with interpolation versus extrapolation are non-comparable, Klawonn et al. (2014a, b) suggested that better estimation of the volume would come from strategies that realistically extrapolate deposit thickness and volume for the proximal and distal portions of the deposit, rather than focusing on the best fit to the thickness versus square-root area values.

Problem demonstration
As the presence of extrapolation uncertainty has been pointed out in Klawonn et al. (2014a), we briefly demonstrate it with the six datasets. V prox , V int , and V dist of each deposit based on each fitted curve are plotted as histograms in Fig. 6. V prox calculated based on different integration limits are also marked. It is well-known that the distal integration limits for the power-law are difficult to specify and justify, we simply specify them based on previous works as stated previously. Figure 6 shows that ratios of the interpolation and extrapolation volumes to the total volume vary by the isopach datasets. Similarly, variabilities of V prox , V int , and V dist based on different fitted curves also vary greatly for each dataset. These corroborate arguments from Klawonn et al. (2014a, b), which suggest that the extrapolation volume and its uncertainty could have a big impact on the total volume estimation.
In addition, Fig. 6 shows that the variability of V prox could be significantly affected by the proximal integration limit when it needs to be specified for a certain curve. This is shown in the 1815 Tambora, Cotapaxi Layer 5, and Hatepe tephra datasets (Fig. 6a, c, and d).

Underestimated uncertainty
Not properly addressing the extrapolation uncertainty could also lead to underestimated volume uncertainty. For a set of isopach data, two different fitted curves could provide similar estimates for the extrapolated thickness and volume, but this does not necessarily suggest that the uncertainty on the extrapolated volume is small. This can be illustrated with the Taupo tephra dataset. The fitted one-segment and Weibull curves both fit the isopach data well (Tables 5 and 6), and the predicted thinning patterns from the two curves are highly consistent with each other (Fig. 3b). Ranges of V prox , V int , and V dist defined by the two curves are 0.36-0.52, 5.33-5.59, and 1.84-1.90 km 3 , respectively. These seem to suggest that variabilities of V prox and V dist are small, and the extrapolated volume is subject to limited uncertainty. These statements are not accurate, because the one-segment and Weibull curves both fit well to the isopach data, but it is likely that they happen to be consistent with each other for the extrapolated thickness. The extrapolated thickness from the two is possibly not consistent with the true thinning pattern of the deposit. In such circumstances, the uncertainty of the extrapolated volume is underestimated or the volume estimate is subject to bias.

Proposed measures to better understand uncertainty from extrapolation
We concur the proposition by Klawonn et al. (2014a) that tephra volumes from interpolation and extrapolation should be reported separately. We recommend that measure to address model uncertainty proposed in this work should be applied to report ranges of V prox , V int , and V dist .
We also propose that instead of trying to quantify the uncertainty for V prox and V dist , it is more objective and accurate to describe the uncertainty from extrapolation as uncertainty that cannot be better and robustly quantified based on the given isopach data. Hence, we can only test whether the total tephra volume is sensitive to the potential variability of V prox and V dist , i.e., treat it as a sensitivity test.
To implement the sensitivity test, we propose to first calculate the tephra volume changes by manually setting the maximum thickness (for V prox ) and extrapolated isopach area for a certain thickness (for V dist ) to different values, and calculate ratios of the volume changes to the total tephra volume. We denote the two volume differences as V prox and V dist and the two ratios as R prox and R dist (Fig. 2f). As R prox and R dist are ratios, their values are comparable among different isopach datasets. Larger R prox and R dist indicate that the total volume is more sensitive to the potential variability of V prox and V dist , respectively. We define T 0,N −1 and A 0.5 * Nth,N −1 as the maximum thickness and the isopach area of half of the thinnest isopach thickness inferred based on the (N-1)-segment exponential curve, respectively.
V prox and V dist are defined as (Fig. 2f): • V prox : the difference between the volumes calculated based on the (N − 1)-segment exponential curve assuming T 0,N −1 to be (a) its original value calculated based on the (N-1)-segment exponential curve and (b) five times of its original value; • V dist : the difference between the volumes calculated based on the (N − 1)-segment exponential curve assuming A 0.5 * Nth,N −1 to be (a) its original value calculated based on the (N-1)-segment exponential curve and (b) one and a half times of its original value.
R prox and R dist are defined as the ratios of V prox and V dist to the total tephra volume calculated based on the (N-1)-segment exponential curve. The denominator is chosen such that curve fitting process can be avoided. In this way, the impact of misfit would not affect values of R prox and R dist . T 0,N −1 and A 0.5 * Nth,N −1 are important for calculating R prox and R dist , but they are defined based on the (N-1)-segment exponential curve. That is to say, their values only depend on the two thickest and two thinnest isopachs. This is again a compromise we have to take to avoid curve fitting. We have tried defining V dist based on manually changing the isopach area of the 0.01-cm isopach, and the resultant values of R dist are similar compared to the current way of defining V dist .
It is possible that the true maximum thickness and isopach area of half of the thinnest isopach thickness exceed what are assumed in defining V prox and V dist . This is likely to occur especially for a scenario in which the deposit can be characterized by a two-segment exponential curve, and the existing isopachs only cover the proximal or distal portion of the deposit. If that happens, it can be recognized by experienced users. Moreover, in such cases, values of R prox and R dist could still be alarmingly large because their values depend not only on the assumed thickness and isopach area ranges, but also on the total tephra volume (calculated based on the given isopachs), thicknesses and areas of the thickest and thinnest given isopachs. This will be demonstrated below. The scenario with well-constrained proximal thinning pattern and poorly constrained distal thinning pattern has been carefully studied and analyzed by Sulpizio (2005) based on 23-32 isopach datasets.
As mentioned earlier, large R prox and R dist indicate that the total volume is sensitive to the potential variability of V prox and V dist . Based on our analysis with the six datasets given below, we are confident that the total volume is greatly sensitive to V prox or V dist if its value is above 0.4. This does not suggest that the volume is insensitive if otherwise.

Application
We calculate R prox and R dist for the six tephra datasets plus the Minoan and 1980 Mt. St. Helens tephra datasets with selected isopachs (Table 8). The latter two tephra datasets can be well-characterized by two-segment exponential curves. In addition to R prox and R dist based on all of their isopachs, we calculate R prox and R dist based on their isopachs that only display the proximal (four thickest isopachs for both; which is similar to the scenario analyzed by Sulpizio (2005)) and distal (four and six thinnest isopachs for the two datasets, respectively) thinning patterns to check how R prox and R dist respond to such circumstances. This is necessary because, as mentioned earlier, the true maximum thickness and true isopach area of half of the thinnest isopach thickness might exceed the ranges assumed in defining R prox and R dist given only the proximal or distal isopachs.
For the datasets with the complete isopachs, we highlight datasets with the greatest R prox and R dist . The 1815 Tambora tephra and the Cotopaxi Layer 5 lead to the greatest R prox (0.681 and 0.693). They have large R prox because: the 1815 Tambora tephra dataset has its thickest isopach (20 cm) with a relatively large area of 144,964 km 2 . The uncertainty of its thinning pattern within the thickest isopach could significantly affect the total volume; For Cotopaxi Layer 5, its thickest (100 cm) and second thickest (50 cm) isopachs imply that the deposit may have a very rapid thinning rate within the thickest isopach (Fig. 3c). Its V int is relatively small, and the true V prox value could take up a large portion of the total volume (Fig. 6c). Its potential variability thus would greatly affect the total deposit volume. The Tambora and Cotapaxi Layer 5 tephra datasets have similarly large R prox values, suggesting that the sensitivity of their total volumes to their V prox is great and at the same level. The two deposits are different in volume, thinning pattern, and number of isopachs, but the proposed measure enables directly comparing the total volume sensitivity to V prox between the two. The Taupo Plinian deposit is characterized by the greatest R dist (0.437). Its thinnest isopach corresponds to a thickness of 12.5 cm. How the deposit thins beyond this isopach and the potential variability of V dist are uncertain and could greatly impact the total volume estimate. The high value of R dist for the Taupo deposit suggests that the total volume is sensitive to V dist . As mentioned earlier, this cannot be reflected if we simply examine V dist estimated based on the fitted curves that possibly happen to be consistent with each other, showing the effectiveness of the proposed measure.
The Minoan and 1815 Tambora tephra datasets have the lowest R prox (0.001) and R dist (0.041), respectively. For the former, its thickest isopach is small areawise (600 cm, 9 km 2 ), and the thinnest isopach (5 cm) has an area of 191,710 km 2 . A significant portion of its volume is from V int and V dist . Regardless of the thickness distribution within the thickest isopach, the total volume of the deposit would not be greatly affected by it and thus not sensitive to V prox . For the 1815 Tambora tephra dataset, the thinnest isopach is thin and large areawise (0.1 cm, 4,288,784 km 2 ), which means that its V dist has to be very small relative to the total volume, and hence would not greatly affect the total volume estimate.
Similarly, the 1980 Mt. St. Helens tephra dataset has 12 isopachs. Its thickest (20 cm) and thinnest (0.05 cm) isopachs have areas of 200 and 167,000 km 2 , respectively. A large portion of the tephra volume belongs to V int , leading to relatively low R prox (0.065) and R dist (0.101) for the deposit. Values of R prox and R dist for other deposits that are not mentioned above range from 0.1 to 0.4.
Given just proximal isopachs, how the distal Minoan and 1980 Mt. St. Helens tephras thin is unknown to us. Their proximal isopachs suggest a great thinning rate, which means that in such circumstances, the assumed isopach area range for calculating R dist might be too small compared to the true values. R dist for the two increase from 0.389 and 0.101 with the complete datasets to 0.609 and 0.575 given only the proximal isopachs, respectively. The latter two values are alarmingly large, indicating that the assumed volume variability for V dist could take up more than half of the total volume given just the proximal isopachs. The total volume is sensitive to the potential variability of V dist in such circumstances. R dist could be alarmingly large here because R dist depends on the total volume and V dist , which is a function of the thinnest isopach thickness and area in addition to the assumed isopach area range. Similarly, R prox of the Minoan and 1980 Mt. St. Helens tephras increase greatly from 0.001 and 0.065 with the complete datasets to 0.238 and 0.442 given only the distal isopachs, respectively. This indicates that the total volume would be a lot more sensitive to the potential variability of V prox if the proximal isopachs are unavailable, suggesting the effectiveness and robustness of the proposed measure in such circumstances. The above results indicate that values of R prox and R dist could effec- tively indicate whether the total volume is sensitive to the potential variability of V prox and V dist given isopach datasets of various coverage and quality.

Conclusions
In this work, we study two sources of uncertainty in estimating tephra volumes using the isopach-based method when the given information about a tephra deposit is relatively sufficient such that the goal is to estimate the tephra volume and quantify its uncertainty (i.e., excluding the scenario when the potential uncertainty is so great that it is of no use to quantify it). The two occur because fitting certain curves to the isopach data on the log(thickness) − √ isopach area plot is needed for the method. The first source of uncertainty is the model uncertainty. It occurs because (1) there is no curve that could always fit perfectly well to the isopach data whether the data uncertainty exists or not and (2) the fitting is done based on log-transformed thickness and square root of isopach area, and as a result, curves that fit poorly to the isopach data could be used for the volume estimation without being noticed. If omitted, this source of uncertainty could introduce additional bias, or lead to under-or overestimated uncertainty for the volume estimate. The second source of uncertainty is from extrapolation, as originally proposed in Klawonn et al. (2014a). It occurs because the predicted thickness for each fitted curve is partially from interpolation and partially from extrapolation. The total tephra volume is the sum of volumes from extrapolation and interpolation, but the two are not comparable because the extrapolated thickness or volume cannot be validated based on data.
The two sources of uncertainty may not always greatly affect the volume estimate, especially in well-constrained datasets, but their importance can be proved theoretically. Different sources of uncertainty propagate hierarchically in tephra volume estimation. The two sources of uncertainty occur in the last step, i.e., during the curve fitting process, before the volume calculation (i.e., volume integration). If they are omitted, the sources of uncertainty in previous steps of the hierarchy might not be properly inherited, potentially leaving the estimated uncertainty subject to misrepresentation (Fig. 1).
We use six tephra isopach datasets to demonstrate the presence of the two sources of uncertainty and show their impact on tephra volume estimation. Propositions to address them are given. For the model uncertainty, the goodness-offit should be evaluated based on isopach areas, and curves that do not fit well to the isopach data should not be used to characterize the tephra volume. We recommend a bottom line criterion that if half or more of the predictions from a fitted curve are outside the ±40% range with respect to the corresponding isopach areas, the curve should not be used for the volume estimation. For a set of curves that satisfy the above or stricter criterion, we propose to use the envelope (i.e., union) of the curves to define the volume range, instead of using volumes estimated from individual curves. Thus the model uncertainty is more accurately captured. This proposed measure can be easily incorporated into methods that quantify other sources of uncertainty in estimating tephra volumes with the isopach-based method.
For the uncertainty from extrapolation, we concur to Klawonn et al. (2014a) that volumes from interpolation and extrapolation should be reported separately. We propose that the uncertainty from extrapolation should be addressed as a sensitivity test. We calculate tephra volume changes by assuming different maximum thicknesses and different isopach areas for half of the thinnest isopach thickness, and use the ratios (R prox and R dist ) of the two volume differences to the total tephra volume to show if and how the total volume is sensitive to the extrapolated volumes within the thickest and outside the thinnest isopachs, respectively. We propose that R prox or R dist being greater than 0.4 indicates strong sensitivity of total volume to the volume within the thickest (V prox ) or outside the thinnest (V dist ) isopachs. Proposed measures to address the two sources of uncertainty are tested against the six isopach datasets, and are proved to be effective.
We hope that this work could help quantify tephra volume uncertainty in a more robust and accurate way, and make tephra volume uncertainty comparable across different tephra deposits in future works when the isopach-based method is used.