The operations needed to decarbonise our energy systems increasingly involve faulted rocks in the subsurface. To manage the technical challenges presented by these rocks and the justifiable public concern over induced seismicity, we need to assess the risks. Widely used measures for fault stability, including slip and dilation tendency and fracture susceptibility, can be combined with Response Surface Methodology from engineering and Monte Carlo simulations to produce statistically viable ensembles for the analysis of probability. In this paper, we describe the implementation of this approach using custombuilt open source Python code (pfs – probability of fault slip). The technique is then illustrated using two synthetic datasets and two case studies drawn from active or potential sites for geothermal energy in the UK, and discussed in the light of induced seismicity focal mechanisms. The analysis of probability highlights key gaps in our knowledge of the stress field, fluid pressures and rock properties. Scope exists to develop, integrate and exploit citizen science projects to generate more and better data, and simultaneously include the public in the necessary discussions about hazard and risk.