Cross-Correlation Analysis of the Stability of Heterogeneous Slopes

Abstract

Geologic material properties of hillslopes are inherently heterogeneous, with complex layering structures due to geological deposition processes. Lacking detailed sampling of the properties’ spatial distribution has led to the stochastic representation of the properties to address uncertainty in the hillslope stability evaluation. This study treats the spatial distributions of the shear strength parameters, the cohesion (c), and the internal friction angle (φ), in a synthetic two-dimensional slope as stochastic random fields characterized by their means, variances, and correlation scales. This study then evaluates the cross-correlation between these parameters and the factor of safety (FS) using unconditional Monte Carlo simulation (MCS). Different from classical sensitivity analyses, the cross-correlation analysis of FS and the stochastic parameter fields stresses the importance of the orientation of the large-scale geological layering, the correlation between the geological media’s cohesion, and the internal friction angle at the local scale on the probability and uncertainty of failure of the heterogeneous hillslope. The analysis further guides the field sampling strategy to reduce uncertainty in the slope stability analysis due to unknown heterogeneity. More importantly, it suggests the location of stability reinforcement measures. The results of this study provide cost-effective tools for geoengineers to deal with field slope stability analysis under uncertainty.

1. Introduction

Shear strength parameters (cohesion (c) and internal friction angle (φ)) are the critical parameters for evaluating slope stability. These parameters exhibit spatial variability at different scales due to complex geological processes [1,2,3]. Many previous studies have concluded that the spatial variability of geotechnical properties is vital in slope stability evaluations [4,5,6]. For example, Qi et al. [7] investigated the effect of the spatial variability of shear strength parameters on a two-dimensional slope’s critical slip surface distribution. Griffiths et al. [8] and Jiang et al. [9] suggested that ignoring the spatial variability of the shear strength parameters would lead to an underestimation of the probability of slope failure when the coefficient of variation of the shear strength parameters was significant. The slope stability evaluated using the mean value of the parameter could be misleading. Cho [10] emphasized the importance of the spatial variability of soil mechanics’ parameters in evaluating failure probability. Additionally, the effect of the spatial variability of soil properties is crucial when the slope risk assessment is over large areas due to both the variability and lack of measurements [11,12]. Based on genetic algorithms and machine learning, Miao et al. [13,14] performed displacement prediction and landslide susceptibility mapping over a large area and evaluated its uncertainty.

In recent years, geostatistical random field theory has been used to describe the heterogeneity of parameters, using its mean, variance, and correlation scales to describe the most probable values, variability, and spatial structure of the parameter fields [15]. The heterogeneous parameter field is statistically isotropic when the horizontal correlation scale is equal to the vertical correlation scale. The field is statistically anisotropic when the horizontal correlation scale is different from the vertical correlation scale.

Combining random fields and Monte Carlo simulations (MCS), one can conduct a cross-correlation analysis of heterogeneous shear-strength parameters and slope stability. Cross-correlation analysis is a method that reveals the spatial relationship between parameters and processes, considering the variance and spatial structure of the parameters. Recently, Cai et al. [16] proposed an effective sampling strategy in slope stability evaluation based on the cross-correlation analysis of shear strength parameters and the factor of safety (FS), defined as the ratio of the slope’s absolute strength to the actual applied load. Cross-correlation analysis is widely used in many fields, such as underground engineering and groundwater science. Using cross-correlation analysis, Gao et al. [17] investigated the spatial relationship between the rock parameters of the unlined rock caverns and the displacement at a location of interest. Mao et al. [18] studied the relationship between observed heads and hydraulic properties at different times and locations of unconfined aquifers during pumping tests. Sun et al. [19] proposed a temporal sampling strategy based on the cross-correlation analysis of hydraulic parameters and observed hydraulic heads.

Most current studies have investigated the effect of the coefficient of variation of parameters on slope stability with statistical isotropic media. However, due to sedimentation processes, geotechnical material often displays a layered structure [7,20,21,22,23]. The failure types of slopes with layered structures are generally translational slides, mudslides, and creep-fatigue [24,25,26]. That is, the correlation scales in different directions of the parameters are different, which leads to the fact that parameter fields with statistical anisotropy are standard. Moreover, Cai et al. [27] developed an analytical approach for the reliability analysis of infinite slope stability in the presence of spatially variable shear-strength parameters. They concluded that the p_f increases when the correlation between c and tanφ increases. Griffiths et al. [8] and Jiang et al. [9] obtained similar conclusions by studying two-dimensional heterogeneous slopes using the random finite element method. The correlation between c and tanφ significantly influences the probability of slope failure. However, the influence of the correlation between c and tanφ on slope stability considering parameter statistical anisotropy has not been comprehensively studied. Likewise, studies have compared the effect of the conditional random field with the unconditional random field on slope stability analysis using MCS with the finite element method. They concluded that conditional random field simulations could address the reduction in uncertainty due to conditioning with sampled parameters in evaluating slope stability [28,29,30,31]. Nonetheless, these studies have not addressed the vital issue of reducing uncertainty by selecting sampling locations. This issue (where to sample and how many sampling locations) to reduce the uncertainty is critical since only a limited number of samples is available in field situations.

This paper investigates the effect of statistical anisotropy of shear strength parameters on slope stability. We first introduce the influence of the correlation scale on the parameter distribution and the method of cross-correlation analysis between the parameters and FS. Then, we use a two-dimensional slope with statistical isotropy to analyze the cross-correlation between parameters and FS. The comparison of the correlation analysis, sensitivity analysis, and traditional limit equilibrium method comes next. The analyses of the cross-correlations between parameters and FS, and the influences of the correlation scales on slope stability in statistically anisotropic media are investigated. Lastly, using conditional random fields with MCS, we compared the uncertainty in evaluating slope stability for four different sampling schemes to validate the results of the cross-correlation analysis.

2. Methodology

2.1. Random Field Modeling of Heterogeneity

This study adopts the random field theory to describe the heterogeneity of shear strength parameters. The impossibility of obtaining the parameter values at every slope location recommends that we consider the c (or φ) of each position as a random variable. A collection of these random variables of c (or φ) in the entire slope becomes a random field characterized by a joint probability density function with mean, variance, and autocorrelation function, which describes the probability of a parameter value at any slope location. The mean represents the most likely value of the parameter, and the variance represents its average deviation from the actual value (i.e., uncertainty) due to variability and lack of measurements. On the other hand, the autocorrelation function quantifies the parameter’s spatial structure (average spatial distribution of the clusters of the geotechnical properties in the slope).

Suppose the slope has n × n random variables, and i and j = 1,2…, n, a two-dimensional autocorrelation function is represented as follows:

$${\rho }_{ij}=\exp \left[-2\left(\frac{\left|x_i-x_j\right|}{{\lambda }_x}+\frac{\left|y_i-y_j\right|}{{\lambda }_y}\right)\right]$$

(1)

where ρ_ij is the autocorrelation coefficient between the parameter at location (x_i, y_i) and location (x_j, y_j), and λ_x and λ_y are the horizontal and vertical correlation scales, respectively. Many other forms of the autocorrelation function are available. They all are ensemble statistics (i.e., general knowledge). We, therefore, chose the most simplistic one. This study considers the correlation scales of c and tanφ as the same.

Physically, the correlation scale represents the average dimensions (e.g., length, width, and thickness) of heterogeneity (e.g., layers or stratifications) within the domain [32]. With specified mean, variance, and correlation scales, numerous realizations of parameter distributions can be generated with different random seeds by a spectral representation method [33,34]. In order to avoid negative values while generating the random fields, the natural log-normal distribution of c and tanφ are employed. Figure 1 shows six realizations of the cohesion parameter field with the same mean and variance (μ_c = 15 kN/m², σ_c = 7.5 kN/m²), but with different correlation scales. These fields are called statistically isotropic when the λ_x is the same as the λ_y (Figure 1a,b), and statistically anisotropic when the λ_x is different from the λ_y (Figure 1c–f). In Figure 1c,d, the correlation scales are 1 m in the vertical direction and 5 m and 10 m in the horizontal direction, respectively. As the λ_x becomes large, the strong and weak zones extend greatly in the horizontal direction, and the slope shows an apparent horizontal layered structure. Figure 1e,f display the cases where λ_x = 1 m, λ_y = 5 m and λ_x = 1 m, λ_y =10 m, respectively. The slope exhibits a vertically layered structure as the λ_y becomes greater than the λ_x. In this study, we consider different cases as follows. For statistical isotropy: λ_x = λ_y = 1 m; for statistical horizontal anisotropy: λ_y = 1 m, λ_x = 5 m, 10 m, 20 m, 40 m, 80 m, respectively; for statistical vertical anisotropy: λ_x = 1 m, λ_y = 5 m, 10 m, 20 m, 40 m, respectively.

2.2. Slope Stability Analysis

MCS is one of the widely utilized methods to estimate the probability of failure and the reliability of a slope with the general knowledge of the mean, variance, and correlation scales of the geotechnical properties of a hill slope. This study, different from many previous studies, adopts the MCS to investigate the effects of correlation scales on the analysis of the probability of failure and reliability of a slope and guides the sampling scheme to reduce the uncertainty of our estimated probability and reliability.

For this purpose, a two-dimensional synthetic slope model for plane strain analysis (Figure 2) with a slope height H = 10 m and slope inclination α = 26.6° was discretized into 385 elements with 1 m × 1 m in size, and some of them are truncated because of the slope surface. The left and right boundaries of the model are zero horizontal displacements (u_x = 0 m), the bottom boundaries are zero horizontal and vertical displacements (u_x = 0 m, u_y = 0 m), and the slope surface is free displacement. The slope is assumed to be subjected to gravity loads only and consists of elastic–perfectly plastic soils following the Mohr–Coulomb failure criterion. Specifically, the loading stress at each element is the total weight of the element above, and the shear strength of each element follows the Mohr–Coulomb failure criterion. Other complex or advanced constitutive models [35,36,37] or numerical simulation methods [38] could be used with corresponding randomized input parameters.

Table 1. Prior statistics of shear strength and other parameters for the numerical model.

Parameters	Values
Mean of cohesion, μc	15 kN/m2
Coefficient of variation of cohesion, COVc	0.5
Mean of friction angle, μφ	10°
Coefficient of variation of friction angle, COVφ	0.5
Dilation angle, ψ	0°
Young’s modulus, E	1 × 105 kPa
Poisson’s ratio, υ	0.3
Unit weight, γ	20 kN/m3

lists the statistics of soil mechanical parameters, except for their correlation scales. Subsequently, the corresponding FS for the entire slope was evaluated based on the finite element strength reduction method (SRM) [39,40]. The program for calculating FS in this study is mainly based on the program p64 [40], and the main difference lies in the automatic MCS and the search for the critical strength reduction factor. The SRM has been widely used due to its practicality and reliability. The factor of safety (FS) is defined as the proportion by which c and tanφ must be reduced in order to cause slope failure. The strength reduction based on the Mohr–Coulomb criterion is shown in the following equation [39,40]:

$$c_{\mathrm{trial}}=\frac{c}{SRF}$$

(2)

$$\tan {\phi }_{\mathrm{trial}}=\frac{\tan \phi }{SRF}$$

(3)

where SRF is the strength reduction factor; c_trial and tanφ_trial are the trial shear strength parameters, which decrease with an increase in SRF. Several gradually increasing values of the SRF are tested, and the updated c_trial and tanφ_trial are used for elastoplastic analysis. When the algorithm does not converge for 1000 iterations, the slope is considered as having failed in this study. The smallest value of SRF causing failure is then interpreted as the factor of safety FS. For example, the FS of this slope using the mean values (Table 1) as input parameters, considering the slope as homogeneous, is 1.141.

Using this synthetic hill slope model and MCS, we evaluated the probability of slope failure p_f and the reliability index β according to the following formula [16]:

$$p_f=\frac{N_{FS<1}}{N}$$

(4)

$$\beta =\frac{{\mu }_{FS}-1}{{\sigma }_{FS}}$$

(5)

where N is the number of realizations in MCS; in this study, N = 500 realizations. N_FS<1 is the number of realizations whose FS value is less than 1 (i.e., the slope fails). A small value of N_FS<1 implies that the probability of failure of the slope is small. In Equation (5), μ_FS is the mean value and σ_FS is the standard deviation of FS values of N realizations of MCS. (μ_FS − 1) represents the slope stability and σ_FS represents the uncertainty in the evaluating FS. Therefore, the larger β is, the more reliable the estimated FS is and the smaller the probability of slope failure.

2.3. Cross-Correlation Analysis

In the next step, we investigate the sensitivity of a slope’s FS to the parameters’ heterogeneity at every part of the slope. Cross-correlation analysis is the sensitivity analysis of system response cast in a stochastic framework with the consideration of variability (variance) and spatial structure (correlation scales) of the heterogeneous slope (e.g., [18,19]). Specifically, a cross-correlation map represents the most likely relationship between system responses and spatial variability in system properties. This study applies the cross-correlation analysis to examine the relationship between shear strength parameters and slope stability:

$${\rho }_{FSc}\left(x_i\right)=\frac{1}{N}{\displaystyle \sum }_k^N\frac{\left(F{S}_k-{\mu }_{FS}\right)\times \left(c\left(x_i,k\right)-{\mu }_c\left(x_i\right)\right)}{{\sigma }_{FS}{\sigma }_c\left(x_i\right)}$$

(6)

$${\rho }_{FS\tan \phi }\left(x_i\right)=\frac{1}{N}{\displaystyle \sum }_k^N\frac{\left(F{S}_k-{\mu }_{FS}\right)\times \left(\tan \phi \left(x_i,k\right)-{\mu }_{\tan \phi }\left(x_i\right)\right)}{{\sigma }_{FS}{\sigma }_{\tan \phi }\left(x_i\right)}$$

(7)

where $x_i$ is the position vector of the parameter considered; ρ_FSc($x_i$) and ρ_FStanφ($x_i$) are the cross-correlations between the cohesion and the internal friction angle at the location $x_i$ and FS, respectively. FS_k is the FS at the k_th realization; c($x_i$,k) is the cohesion at the location $x_i$ and k_th realization; tanφ($x_i$,k) is the tangent of the internal friction angle at the location $x_i$ and k_th realization. Lastly, σ_c($x_i$) is the standard deviation of the cohesion; σ_tanφ($x_i$) is the standard deviation of tanφ.

3. Results of Cross-Correlation Analysis for Statistical Isotropy

This section examines the cross-correlation between FS and two independent variables (cohesion c and the internal friction angle tanφ) with correlation scales of 1 m in both horizontal and vertical directions (statistically isotropic cases). The cross-correlation maps between the shear strength parameters and FS (Figure 3) reveal that ρ_FSc is positive at the toe and the top of the slope (Figure 3a), suggesting a large c value at these areas leads to the greater FS value of the slope. On the other hand, the interior of the slope areas has a positive correlation between tanφ and FS (Figure 3b). The remaining areas have correlation values close to zero, meaning that shear strength in these areas (i.e., most of the slope surface and the back of the slope.) has little effect on the slope stability.

Comparing Figure 3a,b, we observe that the locations of the positive correlation regions of the two parameters are different. Such a difference implies that the mechanisms of c and tanφ on slope stability are different within the slope. Physically, c is the maximum shear stress a rupture surface can carry when the normal stress is absent. On the other hand, tanφ is the coefficient that converts the normal stress at the rupture surface to frictional force. According to the Mohr–Coulomb yield criterion, the shear strength of a geotechnical material is equal to the sum of c and the normal stress multiplied by tanφ. In the interior of the slope, due to the gravity of the overlying material of the slope, the material in the interior is subjected to high normal stresses, and the coefficient tanφ becomes the key factor transforming the normal stress into the anti-slip force. The larger tanφ is, the greater the anti-slip force under the same normal stress. Therefore, the value of tanφ in the slope’s interior greatly influences the slope stability.

On the other hand, at the toe and top of the slope, the geotechnical mass experiences low normal stress, and the material’s cohesion c dictates the shear strength and influences slope stability. The above discussions should explain the differences in the cross-correlation patterns in Figure 3a,b. For this reason, one must recognize that the impact areas of these two parameters are different when evaluating slope stability.

Comparing the potential sliding surface (the red line in Figure 3) calculated by the limit equilibrium method [41] using the mean values of parameters, we observe that the region of positive correlation partially overlaps the surface. This result demonstrates that the cross-correlation analysis of FS and parameters pinpoint the slope’s critical areas. Moreover, the cross-correlation analysis identifies potential sliding zones in contrast to a single sliding surface obtained by the limit equilibrium method.

The sensitivity analysis results of FS to the cohesion (the sensitivity of FS to c, J_c, Figure 4a) and internal friction angle (the sensitivity of FS to tanφ, J_tanφ, Figure 4b) display similar patterns to those in Figure 3a,b, derived from MCS. The sensitivity analysis of FS to the parameter takes the following steps. First, the parameter at each location is set as the mean value, and the FS of the homogeneous slope is evaluated by the finite element strength-reduction method. Then, we applied a perturbation of the parameter at a spatial location x_i, keeping the parameter at the other locations as the mean value. We subsequently evaluated the FS corresponding to this perturbation. The ratio of the change of FS to the perturbation is the sensitivity of FS to the parameter at this position. After calculating the sensitivity at each location, we derived a sensitivity map. The map shows that the sensitivity of FS to c is greater at the foot and top of the slope compared to other regions, and the sensitivity of FS to tanφ is more significant in the interior of the slope. Notice that the sensitivity analysis, based on the perturbation method, aims at the change of FS per change in the given mean value of the parameter, ignoring the variability (variance) and spatial structure (correlation scale) of the parameters [42]. Specifically, the cross-correlation analysis considers many possible slopes with heterogeneous parameter fields with the same mean parameter value but different perturbations and spatial structure patterns. It then summarizes the results statistically. Consequently, the sensitivity analysis results (Figure 4) differ from the cross-correlation analysis (Figure 3), and the cross-correlation analysis is most appropriate for cases where spatial parameter values are unknown (i.e., realistic field situations).

4. Results of Cross-Correlation Analysis for Statistical Anisotropy

4.1. Statistical Horizontal Anisotropy (Horizontal Correlation Scale > Vertical Correlation Scale)

This section investigates the results of the cross-correlation analysis between FS and shear strength parameters (c and tanφ) with different λ_x values (Figure 5) while λ_y = 1 m. Figure 5a,c show the cross-correlation maps between FS and c, for λ_x of c equal to 5, 10, and 20 m, respectively. The cross-correlation maps between FS and tanφ, for λ_x of tanφ equal to 5, 10, and 20 m are presented in Figure 5d,f, respectively. We observe that the positive areas of ρ_FSc and ρ_FStanφ expand as λ_x increases, but the areas are confined to the areas at the slope toe. Comparing the results to the cross-correlation map of the statistically isotropic parameters (Figure 3), we notice that ρ_FSc develops from the toe. In contrast, ρ_FStanφ develops from the inside of the slope.

Figure 6a,d illustrate the probability of failure (p_f), reliability β index, the mean of FS (μ_FS), and the standard deviation of FS (σ_FS) as a function of the normalized horizontal correlation scale, respectively. The normalized horizontal correlation scale is λ_x/(H/tanα), the ratio of the horizontal correlation scale to the horizontal projection of the slope length, H/tanα. Notice that α is the slope inclination angle. Since the correlation between c and tanφ (φ is the friction angle) is generally unclear [8,9,16], this study also examines the effect of perfectly positively, zero, and negatively correlated c and tanφ perturbation relationships and they are indicated by the red, green, and blue lines in these figures, respectively.

According to these figures, if c and tanφ are positively correlated or uncorrelated (i.e., the red and green lines), the probability of failure in Figure 6a and the standard deviation of FS (σ_FS) in Figure 6d rapidly increase as λ_xtanα/H approaches one and stabilize afterward. On the other hand, the reliability (β) (Figure 6b) and the mean of FS (μ_FS) (Figure 6c) decrease exponentially. These results stem from the fact that for horizontally layered slopes, the layer with the lowest parameters (such as the weak interlayer or the stratum with highly developed joints and fractures) controls the stability of the slope. A longer correlation scale means that the layer with the weakest strength covers most of the slope, and the slope is less stable. On the other hand, from the physical meaning of the correlation between c and tanφ, the higher the correlation is, the lower c is at a location, and the lower the tanφ at the same location. Therefore, the slope is less stable (i.e., the red line is higher than the green line in Figure 6a).

The blue lines in all figures depict the behaviors of these quantities for the case where c and tanφ are negatively correlated. The blue lines show that as λ_xtanα/H increases, the p_f value increases slightly but remains very low (about 0.05), while the reliability β index increases and remains high at about 2.4. The value of μ_FS increases from 1.075 to 1.085 and remains constant over the rest of λ_x(tanα/H). The value of σ_FS decreases first at λ_x(tanα/H) = 0.5 and remains almost constant at a small value. The trends are distinctly different from when ρ_ctanφ is 1 or 0. The negative correlation means that a large c is at one location, and a small tanφ is at the same location or vice versa. As a result, the slope stability no longer decreases significantly or even increases slightly.

Overall, an underestimation of λ_x can result in an overestimated slope stability when evaluating the stability of horizontally layered slopes—the importance of identifying the spatial structure of the slope is clear.

4.2. Statistical Vertical Anisotropy (Horizontal Correlation Scale < Vertical Correlation Scale)

Figure 7 shows the cross-correlation map between FS and shear strength parameters (c and tanφ) with different λ_y values when the slope has vertically stratified formations. Figure 7a–c illustrate the cross-correlation maps between FS and c, and λ_y of c with 5, 10, and 20 m, respectively, while their horizontal correlation scales are 1 m. These figures indicate a distinctly positive correlation between FS and c at the toe of the slope. Furthermore, as λ_y increases, the high correlation area becomes more concentrated and vertical, and the cross-correlation value weakens slightly.

Figure 7d–f demonstrate the cross-correlation maps between FS and tanφ, and λ_y of tanφ equal to 5, 10, and 20 m, respectively. We observe that FS and tanφ are positively correlated at the area x = 20 to 25 m and y = 0 to 5 m, and the cross-correlation decreases slightly as λ_y increases.

The behaviors of p_f, β, μ_FS, and σ_FS in the slopes with longer vertical correlation scales than the horizontal one as a function of λ_y/H (the vertical correlation scale normalized by the height of the slope, H) are displayed in Figure 8a–c, and d, respectively. First, we notice that the value of probability failure (p_f) in this case is much less than that in the horizontal layering slope (i.e., Figure 8a vs. Figure 6a), regardless of the effects of various factors as in Figure 6. In other words, vertical stratification (the orientation of the large-scale structures) plays a more dominant role than the others do in slope stability.

Nevertheless, Figure 8 shows that when ρ_ctanφ is 1 (the red line), with the increase in λ_y, p_f and σ_FS decrease, and μ_FS and β increase, indicative of the fact that as λ_y increases, the stability of the slope increases and the uncertainty of the evaluation decreases. This result stems from the fact that when a slope is vertically layered, the high-strength layer controls the stability, similar to anti-slip piles. The longer λ_y means that the layer with high strength is extensive, and when ρ_ctanφ is 1, the c is large and so is φ large, and the anti-slip pile can be effective, leading to high slope stability.

When ρ_ctanφ is 0 or −1 (the green and blue line in Figure 8), the increase in λ_y leads to increases in p_f and σ_FS, and it decreases μ_FS and β values, suggesting that increasing λ_y worsens the stability of the slope and increases the uncertainty of the evaluation. These trends are the opposite of when ρ_ctanφ is 1, likely because c and tanφ are uncorrelated or perfectly negatively correlated, and the effect of anti-slip piles weakens.

4.3. Effects of the Number of Realizations in MCS

All the results above are from 500 MCS. To ensure that the number of simulations is sufficient to obtain representative results, we plot Figure 9 to show the effect of the number of realizations on p_f and μ_FS at λ_x = 10 m and 20 m, respectively, and λ_y = 1 m, in these conditions, σ_FS are the maximums. As shown in the figure, the mean values of p_f and μ_FS fluctuate widely within 150 realizations but stabilize after more than 300 realizations, certifying the adequacies of the number of realizations used in the MCS and the results’ representativeness.

5. Effects of Conditional Random Fields

To demonstrate that sampling in highly correlated zones (conditioning on the stochastic fields) reduces the uncertainty in the FS assessment, we conducted numerical experiments considering two cases where only the cohesion field is a stochastic process. Case 1 used the distribution of c in Figure 1a as the reference field, and Case 2 used the c distribution in Figure 1c. The mean, variance, and λ_y of the random fields in the two cases were identical, except that λ_x was different.

We considered four sampling schemes (Figure 10). In Scheme 1, samples were taken at the toe of the slope at three 5 m deep vertical boreholes (Figure 10a). In Scheme 2, the three boreholes were near the slope top (Figure 10b). Scheme 3 took samples over a 15 m borehole horizontally into the slope at the corner of the slope (Figure 10c). The fourth scheme sampled a 15 m borehole horizontally near the top of the slope (Figure 10d). These four schemes took a soil sample every meter to obtain 15 samples. As such, we could evaluate the effects of the same number of parameter values in different correlation areas for conditioning on the estimates of FS of the slope.

With the sampled data, 500 realizations of the conditional random field corresponding to each sampling scheme were generated using the Kriging Superposition Approach (KSA) [43,44]. The resulting conditional realizations honored the sampled values at the sampled locations and retained the specified spatial statistics of the random fields. Afterward, we used them to conduct MCS as the previous unconditional MCS, and the results are summarized in

Table 2. The results of conditional MCS.

	Scheme 1	Scheme 2	Scheme 3	Scheme 4
μFS of Case 1	1.033	1.041	1.031	1.040
σFS of Case 1	0.066	0.072	0.058	0.072
μFS of Case 2	0.935	0.939	0.924	0.935
σFS of Case 2	0.094	0.105	0.082	0.112

.

The reference FS of Case 1 is 1.031. The μ_FS of sampling Scheme 3 is the closest to the FS of the reference since the sampling area is primarily in the high correlation region. Sampling Scheme 1, which samples a smaller portion of the highly correlated region, yields μFS value that is the second closest to the reference value. Since sampling Schemes 2 and 4 cover the minimal correlation between the parameters and FS, they yield a μ_FS that differs significantly from the reference value and is close to the μ_FS of the unconditional random simulation, indicative of their ineffectiveness for defining the actual factor of safety.

For Case 2, the reference FS is 0.922. The results of the conditional MCS are similar to Case 1, and the mean and standard deviation of FS are the best for sampling Scheme 3. Moreover, due to the apparent layered structure of Case 2 (λx = 5 m), Scheme 3 with horizontal borehole sampling yields significantly better results than the other three schemes. As expected, the uncertainty in the FS evaluated for sampling Scheme 3, which samples the most highly correlated areas, is the smallest in both Case 1 and Case 2. That is, sampling the high correlation areas’ parameters is critical for evaluating slope stability because it reduces uncertainty significantly more than sampling in other regions. As demonstrated in this example, cross-correlation analysis guides optimal borehole placement locations. Furthermore, we could suggest that focusing on the high correlation areas for slope reinforcement measures is essential since the same reinforcement measures will be more effective in the high correlation areas than in the low correlation areas.

6. Conclusions

Cross-correlation analysis of the shear strength parameters and FS of statistically anisotropic heterogeneous slopes shows that the large-scale heterogeneity structures dominate the stability of a slope. Specifically, with the given mean and variance of the shear strength parameters, the probability of failure of slopes with a long vertical correlation scale (vertical structure) is much smaller than that with a long horizontal correlation scale (horizontal structure).

In statistically isotropic parameter fields, the distributions of the high correlation areas of ρ_FSc and ρ_FStanφ are different: the former distributes at the toe and top of the slope, and the latter in the interior of the slope.

For the slope with the long horizontal correlation scale, the high correlation area of shear strength parameters with FS is located at the toe of the slope and distributed horizontally. The larger λ_x is, the longer the extension of the area in the horizontal direction is. In addition, when c and tanφ are positively correlated or uncorrelated, the larger λ_x is, and the less stable the slope is. When c and tanφ are perfectly negatively correlated, the effect of λ_x on the stability of the slope decreases.

In the slope with a long vertical correlation scale, the high correlation area of the shear strength parameter with FS is located at the toe of the slope and distributed vertically. With the increase in λ_y, the correlation in the high correlation region slightly decreases, and the region’s distribution becomes more vertical than others. Moreover, when c and tanφ are positively correlated, the larger λ_y is, and the more stable the slope is. When c and tanφ are negatively correlated or uncorrelated, the larger λ_y is, and the less stable the slope is.

This study further demonstrates that sampling in high cross-correlation regions can reduce the uncertainty in slope stability analysis. In addition, when the slope has an apparent layered structure (statistical anisotropy), the sampling direction consistent with its structure is recommended.

In summary, to evaluate the stability of a slope, first, one must detect the orientation of the large-scale structures (i.e., long correlation scales). The knowledge of the correlation between c and tanφ becomes essential, which dictates the probability of slope failure in the slopes with horizontal layering. The cross-correlation analysis presented in this study yields the location of the critical areas where shear strength parameters affect slope stability. Furthermore, these high cross-correlation regions guide the optimal borehole placement locations to reduce uncertainty in slope stability analysis and even for the selection of slope stability reinforcement locations.