Model Uncertainty in Predicting Facing Tensile Forces of Soil Nail Walls Using Bayesian Approach

The model uncertainty in prediction of facing tensile forces using the default Federal Highway Administration (FHWA) simplified equation is assessed in this study based on the Bayesian inference method and a large number of measured lower and upper bound facing tensile force data collected from the literature. Model uncertainty was quantified by model bias which is the ratio of measured to nominal facing tensile force. The Bayesian assessment was carried out assuming normal and lognormal distributions of model bias. Based on the collected facing tensile force data, it is shown that both the on-average accuracy and the spread in prediction accuracy of the default FHWA simplified facing tensile force equation depend largely upon the distribution assumptions. Two regression approaches were used to calibrate the default FHWA simplified facing tensile force equation for accuracy improvement. The Bayesian Information Criterion was adopted to quantitatively compare the rationality between the competing normal and lognormal statistical models that were intended for description of model bias. A case study is provided in the end to demonstrate both the importance of model uncertainty and the influence of distribution assumptions on model bias in reliability-based design of soil nail walls against facing flexural limit state.


Introduction
Soil nails, the excavated ground, and the facing structure are the three main components of a typical soil nail wall system. Due to complex soil-nail-facing interactions, the design of a soil nail wall must consider various limit states, including global, external, internal, and facing stabilities. This study focuses on facing limit states.
The facing limit states of a soil nail wall mainly refer to facing flexural/bending, punching shear, and headed-stud tensile limit states [1,2]. An appropriate estimation of the facing tensile forces under working conditions is crucial to the success of a rigorous while cost-effective reliability-based design of these limit states as the facing could take up to onethird of the cost of the whole wall [3]. In US, the Federal Highway Administration (FHWA) soil nail wall design manuals [1,2] provide a default simplified model for calculation of the nominal facing tensile forces for preliminary deterministic facing design.
In order to extend the FHWA simplified model to facing design within a probabilistic context, its model uncertainty must be quantified. The importance of model uncertainty in reliability-based geotechnical designs has been thoroughly discussed in [4]. Dithinde et al. [5] recently summarized methods for model uncertainty evaluation and model calibration. The quantification of model uncertainty can be done by characterizing the statistics of model bias which is the ratio of measured to nominal value of the variable of interest. Extensive studies on model uncertainty characterization for different geotechnical models have been reported in the literature, e.g., shallow and deep foundations [6][7][8][9][10][11] and retaining structures [12][13][14][15].
This study is focused on characterization of the model bias in predicting facing tensile forces using the default FHWA simplified model. The model bias is computed as the ratio of measured to nominal facing tensile force. Usually, the measured facing tensile force is taken as the nail head tensile load [1,2], which is difficult to measure directly. This 2 Mathematical Problems in Engineering is because a soil nail is a composite cylinder comprising a steel bar as the nail tendon encapsulated by a grout column. When measuring the nail loads, strain gauges are always mounted on the steel bar; as a result, the readings only account for the tensile loads that are carried by the steel bar, whereas the amounts carried by the grout column are not reflected. Wentworth [16] and Banerjee et al. [17,18] developed a method to estimate the nail head loads considering both contributions from the steel bar and the grout column. Their estimated nail head loads were then used by Liu et al. [19] for evaluation of the model uncertainty of the FHWA simplified facing tensile force equation. Expectedly, the evaluation outcomes depend heavily on the accuracy of the method developed by Wentworth [16] and Banerjee et al. [17,18] for interpretation of the strain gauge readings; however, the accuracy was not assessed.
At nail head, the true tensile load should be no less than the tensile load of the steel bar but also no larger than the sum of the nail head steel load plus the grout column load. That is to say, the former is the lower bound whereas the latter is the upper bound. Lower and upper bound nail tensile load data were reported in [16][17][18], the maximum of which were used in [20] for estimation of the model uncertainty of the FHWA simplified nail load model. The objective of this study is to evaluate the model uncertainty of the default FHWA simplified facing tensile force equation using the lower and upper bound nail head load data reported in [16][17][18]. The model uncertainty is quantified by the model bias whose statistics are estimated using Bayesian inference technique. The Bayesian method has been widely demonstrated to be a powerful tool for estimation of statistical model parameters given measured data [21][22][23][24][25][26][27][28][29][30][31][32][33]. The model bias of the default FHWA simplified facing tensile force model is characterized as a normal or a lognormal random variable, and the competence of the two statistical models is compared using the Bayesian Information Criterion (BIC). This study also conducted calibration of the default FHWA simplified model for accuracy improvement using two regression approaches summarized in [5]. A reliability-based design example of facing flexural limit state is provided in the end to demonstrate the importance of model uncertainty and show the influence of distribution assumption of model bias on design outcomes. The present work is a valuable step towards reliabilitybased analysis and design of facing limit states of soil nail walls.

Performance Functions of Facing Limit States
Facing limit states are classified as structural strength limit states, including flexure limit state (Figure 1(a)), punching shear limit state (Figure 1(b)), and headed-stud tensile limit state (Figure 1(c)). The facing structure is considered as a continuous two-way slab, which endures the flexural moments due to the lateral earth pressures at the soil-facing interface.
The moments cause deformation of the facing and result in tensions. The performance function of the flexure limit state can be written as where F m and T m are measured uncensored ultimate facing flexural capacity and facing tensile force under working conditions, respectively; F n and T n are nominal ultimate flexural capacity and nominal facing tensile force, respectively; and F and T are model biases correcting prediction errors in F n and T n , respectively. Accordingly, F and T are defined as F = F m /F n and T = T m /T n , which are the ratios of measured to nominal values. The FHWA soil nail wall design manuals [1,2] provide simplified equations to compute the nominal values F n and T n , as where C F is the uniformity factor for the soil pressures behind the facing, e.g., C F is the 1 for permanent soil nail walls; h is the thickness of the facing; a n and a m are the reinforcement cross-sectional area per unit width at the nail head and at midspan, respectively, in the vertical or horizontal direction; S h and S v are the horizontal and vertical soil nail spacing, respectively; f y is the reinforcement tensile yield strength; is the empirical spacing factor accounting for the soil arch effect and expressed as = 0.
where is the face batter angle, e.g., is the 0 for vertical facing; is the effective soil friction angle; is the back slope angle, e.g., = 0 for horizontal back slope; and is the interface friction angle between the facing and the soil. As the facing tensile forces increase, a conical failure surface would form locally around the nail head. This local failure mode is called punching shear limit state. Depending on the facing thickness and the type of nail-facing connection, the conical failure surface typically extends behind the bearing plate or headed studs and punches through the facing at an angle of about 45 ∘ [1,2]. The performance function for this limit state is Mathematical Problems in Engineering where P m and P n are measured uncensored and nominal ultimate facing punching capacity, respectively, and P is model bias for P n . Headed studs are routinely used to provide anchorage of the nail into the permanent facing. The headed studs would yield tension once the facing tensile force exceeds their tensile capacity. The performance function for the head stud tensile limit state is where H m and H n are measured uncensored and nominal ultimate tensile capacity of the headed studs, respectively; H is model bias for H n .

Mathematical Problems in Engineering
The load terms, T and T n , are the same for the three facing limit states, as indicated by the performance functions (i.e., (1), (5), and (6)). This study is focused exclusively on the characterization of T for T n where T n is computed using the default and improved FHWA simplified facing tensile force models. The improved model is a calibrated version of the default model for accuracy improvement and will be introduced later in this study in detail. The characterization outcomes for T (i.e., mode, mean, COV, and distribution) will be used to perform a reliability-based example design of facing flexural limit state presented in the end of this paper to demonstrate both the importance of model uncertainty and the influence of model bias distributions on reliability-based design outcomes for soil nail wall facing limit states. The expressions for the nominal punching shear capacity P n and the nominal headed-stud tensile strength H n are not explicitly given here for brevity because they will not be used hereafter. Readers are directed to [1,2] for full formulations of P n and H n .
The analysis of T from a statistical point of view requires a broad database of measured facing tensile forces (T m ) as T is defined as T m /T n . The database of facing tensile forces (T m ) used in this study for characterization of T is described in the next section.

Database.
Tensile forces that the facing of a soil nail wall sustains are usually measured indirectly through the strain gauges mounted on the nail steel bar near the nail head. Typically, a pair of strain gauges is symmetrically mounted on the top and the bottom of the steel bar, as shown in Figure 2, and the average of the two readings is used to compute the nail tensile load at nail head. This nail head load is then taken as the facing tensile force [34].
However, the measured steel loads do not always correctly reflect the true loads that the facing structure is subjected to, especially when the diameter of the grout column is large. This is because the grout column also carries tensile loads which are not recorded by the strain gauges on the steel bargrout interface. The computed facing tensile forces using the steel bar strain gauge readings can only be taken as the lower bounds of the true facing tensile forces.
On the other hand, the maximum tensile loads that the grout column can sustain should not exceed the maximum tensile capacity of the grout column; otherwise, the grout column cracks and all the tensile loads it carries will transfer to the steel bar. Here, the maximum tensile capacity of the grout column is calculated as the product of the tensile strength of the grout and its cross-sectional area. The actual nail head load (facing tensile force) is smaller than or equal to the sum of the measured steel load plus the calculated maximum tensile capacity of the grout column. In other words, this is the upper bound of the actual facing tensile force.
Wentworth [16] and Banerjee et al. [17,18] reported a large number of measured lower bound and upper bound nail load data for nails from ten well-instrumented soil nail  walls across the United States. All the walls were under working conditions. The wall geometry, soil properties, and nail arrangement were described in detail in the source documents [16][17][18]. The descriptions of the walls are not repeated here for brevity. However, the key information of those walls summarized in [20] is reproduced as Table S1 in the Supplementary Materials of this paper for completeness. Each nail in those walls has several measurement points along the nail, and each point has a measured lower bound nail load and a measured upper bound nail load. All the measured nail loads are summarized in Table S2 in the Supplementary Materials. Nominal facing tensile forces T n using the default FHWA simplified model (3) are also provided in Table S2 for comparison. The maximum lower and upper bound nail load data given in Table S2 were used in [20] for characterization of the model bias of the default and modified FHWA simplified nail load models, while in this study only the load data from the strain gauges that were nearest to the facing are adopted for further analyses. According to the source documents [16][17][18], most of the strain gauges were mounted at a distance of about 0.15 to 0.6 m from the nail heads; a few were at a distance up to 1 m.
A total of n = 112 measured lower and upper bound facing tensile force data points were provided in Table S2. Each dataset has n = 56 data points. Since soil nail walls are rarely built in soft to very soft soils due to the concerns of additional stability and settlement issues [1,2,35], the measured facing tensile force data (n = 11) for walls W6 and W7 (Table S1) are abandoned from further analyses. The source documents also identified questionable data for walls W9 (n = 3) and W10 (n = 6) due to construction problems. Interested readers are directed to [16][17][18] for detailed explanations. These n = 9 questionable data points are also excluded from the analyses to follow. After filtering, the number of facing tensile force data adopted in this study to perform Bayesian inference analyses is n = 72 (bold in Table S2); each dataset has n = 36 data points. These load data were collected from soil nail Mathematical Problems in Engineering  walls built in frictional and cohesive-frictional soils. The soil friction angles are typically from 30 to 40 ∘ , and soil cohesions are typically less than 10 kPa. The data are considered as longterm facing tensile forces because they were collected several months to several years after the end of construction of the soil nail walls.

Initial
Screening. The current FHWA soil nail wall design manual [2] recommends the use of / = 1.0 for computation of K a using (4). This recommendation of / = 1.0 is adopted in the analyses to follow. Nevertheless, later in this study a sensitivity analysis is also carried out to examine the influence of / on the estimation outcomes of the model bias statistics. Figure 3 shows plots of normalized measured facing tensile force ( ) versus normalized nail depths (z/H). The normalized measured facing tensile force values are computed as Nominal facing tensile force, ４ Ｈ (kN)  Figure 5. Bias data are examined in two forms, including raw bias values and logarithm of bias values. Both raw and logarithm of bias values seem to follow normal distributions as they (at least visually) can be adequately described using first order polynomials, regardless of lower or upper bound bias data. This is quantitatively confirmed by applying the Kolmogorov-Smirnov tests to the four datasets in Figure 5. Hence, both lower and upper bound bias values can be taken as normally and lognormally distributed. Based on these findings, the model bias of the default FHWA simplified equation is assumed to be a normal random variable and a lognormal random variable. Quantitative outcomes of accuracy evaluation on the default FHWA facing tensile force equation using the Bayesian method are presented in the following.

Bayesian Inference Technique for Model
Bias Statistics Bayes' theorem expressed as where is the updated parameter vector and usually referred to as a posteriori; (D | ) is the likelihood of the observations D conditioned on the a priori ; and k is a normalizing constant. The observations D can be uncensored data, censored data, or both. Censored observations include right censored, left censored, and interval censored data. The likelihood is the joint probability densities of all observations conditioned on . The formulation of depends upon the distribution of X and the data types of D. Suppose a number of n interval censoring data are available for X as D = ( 1 < 1 < 1 , 2 < 2 < 2 , . . . , < < ), then the log form of the likelihood can be written as (e.g., [20,[36][37][38][39]) if X is normally distributed, and as ln ( , | ln , ln ) if X is lognormally distributed.
Here, and are log-mean and log-standard deviation of X, respectively; erf( ) is the error function. and can be easily transformed to the corresponding and using (e.g., [40]).
= ln − 0.5 2 (10) Substitute (10) and (11) into (9); then the log-likelihood function can be uniformly written as (D | ) where = ( , ), regardless of normal or lognormal X. For the prior joint distribution ( ), it can be simply taken as 1 if noninformative. As such, the a posteriori ( ) is exclusively related to the likelihood function. Discussion on selection of a priori for Bayesian updating can be found in [41,42]. The likelihood ratio test [43,44] can be employed to determine the confidence intervals of the log-likelihood functions as: where ln is the mode of the log-likelihood function; ln is the log-likelihood value corresponding to ; is the desired level of significance for the confidence interval; 2 ;1 is the 100(1− ) percentile point of a Chi-Square distribution with one degree of freedom. For example, taking the confidence intervals of 50%, 90%, 95%, and 99%, then the 2 ;1 values are equal to 0.455, 2.706, 3.841, and 6.635, respectively.
It is usually difficult to derive analytical solutions for ( ); numerical simulation techniques must be adopted. The most commonly used algorithm is the Markov Chain Monte Carlo (MCMC) approach. An excellent review and the technical details of the MCMC method can be referred to [45].

Comparison between Competing Statistical Models.
Two statistical models, i.e., normal and lognormal distributions, are intended for the random variable X. The relative suitability of the two distribution models can be quantitatively compared based on the observed data i and i using the Bayesian Information Criterion (BIC) [46]. The BIC value of a statistical model is calculated as: where n = number of observed data; m = number of model parameters. In this paper, m = 2 since both the normal and lognormal distribution models have two model parameters, i.e., and . The criterion states the smaller the BIC value of a statistical model, the better the model describes the observed data. To be more intuitive, the BIC value of a statistical model can be further used to compute the probability of a model being the best among all competing models. This probability of being the best model is calculated as (e.g., [47]): Where Δ (BIC) = BIC − min{BIC }, = 1, 2, . . . ; = number of candidate models; min{BIC j } = minimum value among all the BIC values.  (7), respectively. The a priori of the parameter vector is taken as 1, i.e., noninformative. The Metropolis-Hastings algorithm (e.g., [45]) was employed to generate a single Markov chain with 200,000 samples for each model parameter, as shown in Figure 6. The first 50,000 samples were discarded as burn-in, whereas the remaining 150,000 were used for further analyses. Figure 7 shows full distributions (histograms) of the a posteriori and and the likelihood ln with fitted curves using the Kernel density smoothing technique. Based on Figures 7(a) and 7(b), the modes and means were found to be 0.90 and 0.87 for , and 0.400 and 0.425 for , respectively. The modes and means are close, suggesting that both and are approximately normally distributed. Taking the modes of and as the mean and standard deviation of the model bias ( T ), the default FHWA simplified method would on average overpredict the facing tensile forces by about 10%. The spread in prediction accuracy expressed as coefficient of variation (COV) of T is about 44% which falls within typical range of 30% to 50% identified in [48]. The Pearson's correlation coefficient between and was computed as , = 0.197. The database used in this study was built upon facing tensile force data collected from soil nail walls across the United States. As a result, the reported analysis outcomes should be taken as default statistics for T , which do not specifically apply to any soil nail wall projects. However, for design of a specific soil nail wall the default T statistics can be taken as the a priori and be updated using project-specific bias data to reflect the unique project conditions. From this sense, although most studies typically reported the estimated modes or means of parameters, we argue that all the modes, means, COVs, and correlations among parameters be fully reported for a Bayesian analysis. The mode of the log-likelihood was determined as −41.799, as demonstrated in Figure 7(c). Equation (12) was used to compute the 50%, 90%, 95%, and 99% confidence intervals for the log-likelihood, which were −42.027, −43.152, −43.720, and −45.117, respectively. Taking the 95% confidence interval as an example, it means that any pairs of and or and resulting in a log-likelihood value falling between −43.720 and −41.799 cannot be rejected as a reasonable pair of the model parameter estimates. Figure 8 shows the contour plots of versus corresponding to those confidence intervals based on the 150,000 pairs of samples (c.f. Figure 6). Specific values of , , and are summarized in Table 1. The 95% confidence intervals are from about 0.74 to about 1.07 for mean of T and 0.335 to 0.629 for COV of T . From a practical point of view, these 95% confidence intervals are not unreasonably wide, suggesting that the Bayesian analysis is efficient.

Lognormal Model
Bias. Similar analyses were conducted for the lognormal case based on (9). The analysis outcomes are shown in Figures 9-11. The modes and means were determined to be 1.03 and 1.13 for , and 0.946 and 1.251 for , respectively. The deviations of modes from means suggest that both and are not normally distributed, which can be visually confirmed based on the long tails shown in Figures 10(a) and 10(b). Based on the modes of and , the COV of T was computed as = 0.918. These estimation results suggest that the default FHWA facing tensile force calculation model is satisfactorily accurate  on average; however, the scatter in prediction accuracy is extremely large. These outcomes are very different from those obtained assuming normal model bias. The large bias COV in this case can be explained by Figure 5 where the lower tail of the cumulative distribution of the logarithm of T [ln(T m /T n )] deviates significantly from the overall trend. The Pearson's correlation coefficient was computed as , = 0.895, indicating high correlation between and . The confidence intervals are also investigated and the outcomes are shown in Figure 11 and Table 1.

Model
Competence. Earlier in this study it was demonstrated through the Kolmogorov-Smirnov tests that T for the default FHWA simplified method (3) can be treated as a normal or a lognormal random variable. To compare the relative feasibility between these two statistical models, (13) is used to compute the BIC values for the normal and lognormal models and then the probabilities of being the best model are computed using (14). The BIC values are computed as 90.76 for the normal case and 106.20 for the lognormal case. The corresponding P best value is 99.96% for the normal assumption whereas it is only 0.04% for the lognormal assumption. This means that if only normal and lognormal statistical models are intended for the model bias T of (3), then based on the collected interval censoring bias data it is almost certain that the normal model is the better one. In other words, the T of (3) should be treated as a normal random variable. However, the general concern of assuming a normal model bias is that reliability analysis using Monte Carlo simulation technique could generate negative bias values which violate the physical definition of T , i.e., ratio of measured to nominal facing tensile force which must be positive. FHWA recommendation of / = 1. The influence of / is demonstrated insignificant on estimated modes of and , regardless of the normal or lognormal case, as shown in Figure 12. Hence, the adoption of / = 1 in the previous analyses is justified and will be also used in the analyses to follow.

Improved FHWA Simplified Model for Facing Tensile Forces
Dithinde et al. [5] summarized two regression approaches that have been widely used in the literature for geotechnical model calibrations. The first approach is to regress the model bias against the nominal facing tensile forces. Although simple and straightforward, this approach does not provide insights into the sources of statistical correlations and thus it may not be the best from the perspective of accuracy improvement. The second approach is to conduct a multiple regression of the model bias against the input parameters of the default FHWA facing tensile force equation. This approach usually gives better accuracy improvement; however, the regression expression could be too complicated to be practical if the input parameters are too many. Both approaches are used to calibrate the default FHWA facing tensile force model in this study. Figure 13 shows the plots of measured lower and upper bound model bias data against nominal facing tensile forces computed using the default FHWA method (  the datasets in Figure 13. The Spearman's values are negative and the p-values are less than 0.05 for both cases, implying the negative correlations of the datasets at a level of significance of 5%. As such, the current default FHWA facing simplified model equation can be modified as

Approach 1: Regression to Nominal Facing Tensile Force.
where M is a dimensionless empirical modification factor. When regression to the nominal facing tensile force (T n ) is conducted, the modification factor can be generally written as M = f(T n /A t /P a ), where A t = 1.5 × 1.5 m 2 = 2.25 m 2 is defined as typical tributary area [20,35]; P a = 101 kPa is the atmospheric pressure. The introduction of A t and P a makes M nondimensional.
In order for (15) to be practical, the expression of M should not be too complicated. In this study, four simple expressions for M are examined, including linear, exponential, logarithm, and power formulations. It was found that the logarithm formulation for M yields the minimal model bias COV in this case, which is written as where a and b are empirical constants to be determined. The optimal values of a and b are determined as the pair that meets two requirements: (1) the mode of model bias is equal to 1.00; (2) the mode of model bias standard deviation is as small as possible. Rounding to two decimal places, the optimal values for the two constants are determined as a =   Figure 11: Confidence intervals of 50%, 90%, 95%, and 99% for estimated mean and COV for model bias (assuming lognormal model bias). and lognormal cases while the corresponding COV values are 0.431 and 0.896. Hence, the improved FHWA model is accurate on average and the spread in prediction accuracy is about 2% less compared to the default FHWA model, i.e., 0.431 versus 0.444 for the normal case and 0.896 versus 0.918 for the lognormal case.

Approach 2: Regression to Input Parameters.
Theoretically, to allow the best accuracy improvement, the model bias of the default FHWA facing tensile force calculation equation should be regressed against all input parameters and the modification factor can be written as The penalty of such a full multiple regression is that the improved FHWA equation could be too complicated to be practical. To strike the balance point, this study takes that M can be representatively written as a function of only one input parameter of (3), i.e., M = f(x) where x is , , K a , H+q s , or S h S v . Table 2 summarizes the results of Spearman's rank tests that were applied to the measured lower and upper bound model bias (T m /T n ) against input parameters , , K a , H+q s , and S h S v . The results suggested that M is highly unlikely to be a function of at a level of significance of 5% because the p-values using both lower and upper bound bias data are higher than 0.05. Input parameters , K a , H+q s , and S h S v were found to be candidate regression variables for x in M = f(x) because for these parameters at least one of the two pvalues is less than 0.05.
Again, the four simple expressions used earlier for M are examined for each of the four candidate regression variables, resulting in a total of 16 cases. It was found that taking M as a linear function of ( H+q s ) gave the lowest model bias COV. The formulation of M is as follows.
Here, a and b are regression constants to be determined. The optimal values for b in (17) Table 3.

Case Study
This section presents a reliability-based design of the facing flexural limit state of a real case provided in the FHWA soil nail wall design manual [1]. The design considers seven different cases of the model bias, and the design outcomes are compared to that by [1]. Through this case study, the importance of model uncertainty and the influence of assumption of model bias distributions in reliability-based design of soil nail wall facing limit states are demonstrated. The example wall taken from [1] was H = 10 m in height with vertical facing and horizontal back slope. The nails were spaced at S h = S v = 1. = 0 kPa. The facing was constructed using CIP concrete with a thickness of h = 200 mm, which gave C F = 1 in (2). Inside the concrete slab was the steel bars mesh spaced at 300 mm each way. The nominal yield strength of the steel bars was f y = 420 MPa. The performance function for the facing flexural limit state, i.e., (1), has four random variables, including F , T , , and . In this case study, the mean values of both and are taken as their design values while the corresponding COV values are taken as 0.10 and 0.05, respectively, which are typical [48,49]. The model uncertainty of the nominal ultimate flexural capacity, F n , is assumed to primarily stem from the uncertainty in f y . If so, then the mean and COV of the model bias F for F n in (1) can be commonly assumed to be 1.10 and 0.10, respectively [50,51]. The random variables F , , and are assumed to be lognormally distributed. For T , the statistics and distributions were summarized in Table 3.
The main design parameter for the facing flexural limit state is the total reinforcement cross-sectional area per unit width, (a n + a m ). Controlled by the empirical depth factor , both the default and modified nominal facing tensile forces computed using (3) and (15) are trapezoidal curves with depth (e.g., Figure 3). With uniform facing thickness (i.e., h = 200 mm in this case study), the most critical depth would be from z/H = 0.2 to z/H = 0.7 where reaches its maximum value of 0.75. As a result, the analyses to follow were conducted using z/H = 0.5.
In the FHWA manual [1], the total reinforcement crosssectional area per unit width was determined as (a n + a m ) = 1326 mm 2 /m, resulting in a minimum factor of safety (FS) over 4.3 which far exceeds the minimum requirement of FS = 1.5. In this case study, the (a n + a m ) was redetermined on a target reliability index ( T ) basis. Three target reliability indices were selected as T = 3.0, 3.5, and 4.0. The outcomes of reliability-based redetermination of (a n + a m ) satisfying different T values are summarized in Table 4 and shown in Figure 14. The design was based on the performance function (1) and the model uncertainty ( T ) was considered for seven cases.
The first observation is that the facing flexural limit state would be underdesigned (unsafe) if model uncertainty is not taken into account. For example, by using the default FHWA facing tensile force equation, the (a n + a m ) value to satisfy T = 3.5 is 324 mm 2 /m without considering model uncertainty (i.e., Case 1 in Table 4). This value would be increased to 536 mm 2 /m or 2240 mm 2 /m if taking the model bias as a normal or lognormal random variable, respectively.
The second observation is that the assumption of distribution of the model bias has a huge impact on the design outcomes. For example, to achieve the same T values the required (a n +a m ) would be much less based on normal model bias than those based on lognormal model bias, i.e., Cases 2, 4, and 5 against Cases 3, 6, and 7 in Table 4. The design by the FHWA manual [1] where (a n + a m ) was determined as 1326 mm 2 /m would be judged to be very conservative based on the normal model bias assumption; however, opposite judgement would be made if based on the lognormal model bias assumption. This highlights the importance of model bias in reliability-based design of facing limit states.     15) and (16) Modified FHWA, Eqs. (15) and (17) Lognormal bias Default FHWA, Eq. (3) Modified FHWA, Eqs. (15) and (16) Modified FHWA, Eqs. (15) and (17) Figure 14: Design outcomes of (a n +a m ) with respect to different target reliability indices using default and modified FHWA simplified facing tensile force models.

Summary and Conclusions
This study presents a Bayes theory-based assessment on the model uncertainty of the default FHWA simplified facing tensile force equation using measured lower and upper bound data of facing tensile forces collected from the literature. The model uncertainty was quantified by model bias which is the ratio of measured to nominal facing tensile force. The computed model bias values are divided into lower and upper bound datasets, which were later demonstrated to be both normally and lognormally distributed through Kolmogorov-Smirnov tests. The Bayesian assessment of the model uncertainty was carried out for two cases: (1) model bias is a normal random variable; (2) model bias is a lognormal random variable.
Based on the collected data, facing tensile forces by the default FHWA simplified facing tensile force equation on average would be overestimated by about 10% based on the normal bias assumption but underestimated by about 3% based on the lognormal bias assumption. Correspondingly, the spreads in prediction accuracy expressed as the COV of bias were about 45% and 92%, respectively. The assessment outcomes are very different, depending on the assumption of the model bias distributions. Two regression approaches widely adopted for model calibration were used to modify the default FHWA facing tensile force equation for accuracy improvement. The modified FHWA simplified equations were proved to have bias mean values equal to one and smaller bias COV values for both normal and lognormal cases.
The Bayesian Information Criterion (BIC) was used to further quantify the rationality between the two competing statistical models (normal and lognormal) that were intended for description of model bias. The BIC values of the two candidate models were also used to compute their probabilities of being the best. The results showed that model bias for the default facing tensile force prediction equation is much better described as a normal random variable.
A case study of reliability-based design of facing flexural limit state was presented to highlight the importance of model uncertainty and the influence of distribution assumptions for model bias on design outcomes. The work presented in this study is an important step towards both direct rigorous reliability-based design and calibration of resistance factors for load and resistance factor design for facing limit states of soil nail walls.
Finally, the statistics of model biases quantified in this study should be understood as default values and be used where project-specific measured facing tensile force data are not available. For cases where project-specific data are available, the model bias statistics should be updated to reflect the features of the projects. The updating can be performed using Bayesian inference technique in which the priors are the default model bias statistics determined in this study.

Data Availability
All the facing load data used to support the findings of this study were collected from the literature cited in the article. All the data are included within the article or within the supplementary information file. There are not any restrictions for the reader on data access. Table S1. Summary of wall geometry, soil properties, and nail arrangement for soil nail walls reported in the source documents (after [20,35] and source data from [16][17][18]). Table S2. Summary of measured lower and upper bounds of nail load reported in the source documents [16][17][18] at each measurement point and the corresponding nominal T n using the default FHWA simplified facing tensile force model (Supplementary Materials)