In normally consolidated, shallow (depth \< 3 km) siliciclastic sequences, faults develop clay smears. Existing models include the dependence of permeability on the clay fraction, but improved predictions of fault permeability should account for uncertainty and anisotropy. We introduce PREDICT, a methodology that computes probability distributions for the directional components (dip-normal, strike-parallel, and dip-parallel) of the fault permeability tensor from statistical samples for a set of geological variables. These variables, which include geometrical, compositional, and mechanical properties, allow multiple discretizations of the fault core to be populated with sand and clay smears, which can be used to upscale the permeability to a coarser scale (e.g., suitable for reservoir modeling). We validated our implementation with experimental data and applied PREDICT to several stratigraphic sequences. We show that fault permeability is controlled by the clay smear configuration and, crucially, that it typically exhibits multimodal probability distributions due to the existence of holes. The latter is a unique feature of our algorithm, which can be used to build fault permeability scenarios to manage and mitigate risk in subsurface applications.