## 1 Introduction  
Slope failure modes include planar, wedge, circular, and toppling failures [1]. The methods used to explore these failure modes include field investigations [2; 3; 4], the traditional limit equilibrium method back-analysis of stability [5; 6; 7; 8], and the numerical simulation method [9; 10; 11; 12]. However, the characteristic information is hard to obtain effectively in a forward analysis because of no failure occurrence. It is necessary to explore potential failure modes, and the current understanding of this is insufficient.  
In recent years, many researchers have reported that different failure modes may occur in a slope [13; 14; 15]. They have studied the influence of weak interlayers (such as joints, bedding, foliation, and other discontinuities) on slope failure modes using numerical simulations and physical experiments [16; 17; 18; 19; 20]. On the other hand, experiment studies of structural plane materials [21; 22; 23] have explored the slope failure mechanism from a microscopic point of view. Additionally, some novel numerical methods for modeling the structure plane were proposed to assist in analysis [24; 25]. However, few articles have mentioned the effects of weak interlayers at different locations. Experiments or simulations based on a single weak interlayer may not be able to characterize the failure mechanism comprehensively. It is necessary to investigate slope failure modes based on multiple potential sliding surfaces. Additionally, in forward analysis, the numerical method may make it difficult to identify all the potential sliding surfaces since the failure has not yet occurred. The traditional limit equilibrium method (LEM) performs well in this regard due to its clear concept and convenience.  
Identifying the triggering factors of different failure modes has a guiding significance for practical engineering. Contact between water and rock reduces the rock's strength; this is called softening and is summarized by [PERSON] [26]. This mechanical property of rock is influenced by saturation [27; 28; 29]. Anchor prestress loss is one of the main factors affecting the stability of anchored rock masses and has been widely investigated [30; 31; 32; 33]. Earthquakes are factors that cannot be ignored for the stability of slopes [34; 35; 36]. All these studies are focused on a specific factor. However, slope failure may occur under the influence of multiple factors. It is worthwhile to explore whether the trigger factors of different failure modes are the same.  
The factor of safety (FOS) is used to quantify the probability of slope failure. Because the FOS has a strong linear correlation, it is hard to identify the triggering factors of each failure mode from the FOS directly. The present study attempts to use statistical methods to analyze the influence of the above-mentioned factors on the FOS and to identify the triggering factors of different failure modes. In multivariate analysis, exploratory factor analysis can effectively reduce the dimension of factors and eliminate the correlation between factors [37; 38]. The factor load matrix can explain the foundations of classification [39]. The present study regards FOS as a sample. This sample is divided into two categories, according to the different failure modes. EFA was used to analyze these two samples to obtain different factor load matrices. These factor load matrices were then used to explain and identify the triggering factors of each failure mode.  
The present study has two research objectives. The first one is to study the influence of potential sliding surfaces at different locations on the slope failure modes. The second one is to capture the main triggering factors of each failure mode. This work contains three parts. The first is the detection of potential sliding surfaces (Section 3.1). The second is the numerical evaluation of failure characteristics (Section 3.2). The last is a multivariate statistical analysis of the FOS values (Section 3.3).