Geotechnical engineering design has relied upon deterministic methods of analysis whereby values for analysis parameters and conditions are selected subjectively based on judgment with the intent of providing acceptable margins of safety. The objective of this research was to improve the use of probabilistic slope stability analysis in practice so that the design of slopes can be made on a consistent and probabilistic basis. The current research involved the development of a methodology for the measurement and modeling of the anisotropic autocorrelation distance of cohesive soils, which was demonstrated at Dyke 17 West of the McArthur Falls Generating Station. In-situ testing using the piezocone and laboratory testing was conducted to characterize the spatial variability of the effective-shear strength envelope. Vertical (down-hole) and horizontal (cross-hole) geostatistical analysis was conducted to assess the anisotropy of the semivariogram. The investigation identified that heterogeneous inclusions had significant impacts on the results, but that simplistic (visual) identification and filtering procedures were adequate. The effective-stress shear strength envelope was statistically characterized as a random field, which was simulated as a first-order Markov process using customized add-in functions in a limit-equilibrium slope stability analysis. The analysis accounts for the spatial variability of shear strength and is capable of simulating both isotropic and anisotropic autocorrelation functions. The study showed that the critical slip surface geometry and the probability of failure can be significantly different when the anisotropy of spatial correlation is accounted for. The study also showed that neglecting spatial correlation may over-estimate the probability of failure, however this finding was noted to be likely case-specific. The primary conclusion of the study was that appropriate representation of spatial correlation is essential to calculating the probability of failure. Finally, convergence of the probabilistic simulation was evaluated using bootstrapping of the simulated factor of safety distribution to assess the standard error in the mean factor of safety, standard deviation of factor of safety and the probability of failure. A convergence criterion based on the percentage standard error in the probability of failure was proposed and used to define the number of Monte-Carlo iterations required.