This study proposes a subset simulation (SS)-based approach for efficient slope reliability analysis involving copula-based cross-correlated random fields of cohesion (c) and friction angle (ϕ) of soils. First, the copula theory for modeling the cross-correlation between c and ϕ is briefly introduced. The algorithms for generating the copula-based cross-correlated random fields of c and ϕ are detailed. Then, the SS for efficient slope reliability analysis involving copula-based cross-correlated random fields of c and ϕ is explained. Finally, two slope examples with the same geometry but different sources of probability information are presented to illustrate and demonstrate the proposed approach. The results indicate that the proposed approach has both good accuracy and efficiency in slope reliability analysis involving the copula-based cross-correlated random fields of c and ϕ at low failure probability levels. The copula theory for characterizing the cross-correlated random fields can consider both the Gaussian and non-Gaussian dependence structures between c and ϕ. The copula selection has a significant impact on slope reliability with spatially variable c and ϕ. The probabilities of slope failure produced by different copulas differ considerably. This difference increases with decreasing probability of slope failure. The commonly-used Gaussian copula may lead to a significant underestimate of the probability of slope failure. The reasonable identification of the best-fit copula for characterizing the cross-correlated random fields of c and ϕ based on the test data is highlighted in practical slope reliability analysis.