Appendix A. Uncertainty assessment for annual mean embayment area, volume, and depth.
The propagation of uncertainty to annual mean estimates of embayment area (A, m^{2}), volume (V, m^{3}), and depth (Z, m) demonstrates the use of statistical models, bootstrapping approaches, and classical error propagation. As shown in Fig. 2, the mean area and volume of the embayment were subject to errors from the 3D digital model used to represent the embayment and from the WSE estimate. To address these errors we used a twostep process. First, a universal kriging algorithm was fit to the elevation coverage of 8055 points to generate a 3D continuous surface of the study site using a commercial GIS software (ArcGIS, Geostatistical Analyst extension). The final form of the kriging algorithm was determined iteratively by varying the model to produce the smallest root mean square error (RMSE = 0.06) and also reproduce, as close as possible, a 1:1 relationship between observed and predicted elevations (slope = 0.998). The kriging algorithm generated 8055 predicted elevations and a SE of each prediction. The predicted elevations were next used to create a triangulated irregular network (TIN) hull that could be used in the GIS software to determine embayment area and volume. The uncertainty in the predicted elevations could not be propagated into the estimates of area and volume using classical methods since the analytical expression of the algorithm used by the proprietary commercial software (ArcGIS) to derive the TIN hull from the predicted elevations was not known. Instead a Monte Carlo approach was implemented to propagate the error. Each predicted elevation and its SE was used to generate 20 normally distributed elevation values within the 99.99% CI (t value = 3.891) of the predicted elevation. The BoxMueller equation was used to randomly generate the elevations within the specified CI:
(A.1) 
for j = 1 to 20, where z_{ij} is an i (8055) x j (20) matrix of elevations, predicted z_{i} is the ith elevation predicted by the kriging algorithm, t is the t value, SE_{predictedZi} is the SE of predicted z_{i}, and RAND_{ij} is a random number from 0 to 1. Subsequently, each j set of z_{i} was used to create a TIN hull which resulted in the development of 20 TIN hulls.
The second step involved using the TIN hulls to calculate the area and volume of the embayment (ArcGIS 3D Analyst) based on the mean WSE. Twenty WSEs were randomly selected from the large population of measured and modeled WSEs (see Appendix F). Each of these randomly drawn WSEs also had a SE. The SE for measured values resulted from adjustment to the NAVD88 vertical datum, and the SE for modeled values resulted from prediction when the Gongora gauge was inoperable and adjustment to the NAVD88 vertical datum (see Appendix F). To fully incorporate this error a 99.99% upper and lower limit for each of the 20 WSEs was calculated, resulting in 40 WSEs.
The decision to use only 20 TIN hulls and the initial draw of 20 WSEs was based on the consideration that these had to be manually entered into a graphical user interface for the calculation of area and volume in the GIS software. Thus, with the 40 WSEs there were 20 × 40 (800) data entries with corresponding area and volume estimates to tabulate. The mean area and volume were calculated from the 800 area and volume estimates with SEs calculated by equation (2) and df = 799.
Finally, the embayment average depth (Z) was calculated as Z = V/A with SE_{Z} and df_{Z} calculated by Eqs. 3 and 7, respectively.
(A.2) 

(A.3) 
Mean area, volume, and depth are shown in Table A1.
TABLE A1. Scalar values of mean embayment area (A), volume (V), and depth (Z). SE and df were the standard error and degrees of freedom, respectively.