Cavitation damage on a mercury target vessel for a pulsed spallation neutron source is induced by a proton beam injection in mercury. Cavitation damage is one of factors affecting the allowable beam power and the life time of a mercury target vessel. The prediction method of the cavitation damage using Monte Carlo simulations was proposed taking into account the uncertainties of the core position of cavitation bubbles and impact pressure distributions. The distribution of impact pressure attributed to individual cavitation bubble collapsing was assumed to be Gaussian distribution and the probability distribution of the maximum value of impact pressures was assumed to be three kinds of distributions: the delta function and Gaussian and Weibull distributions. Two parameters in equations describing the distribution of impact pressure were estimated using Bayesian optimization by comparing the distribution of the cavitation damage obtained from the experiment with the distribution of the accumulated plastic strain obtained from the simulation. Regardless of the distribution type, the estimated maximum impact pressure was 1.2–2.9 GPa and existed in the range of values predicted by the ratio of the diameter and depth of the pit. The estimated dispersion of the impact pressure distribution was 1.0–1.7 μm and corresponded to the diameter of major pits. In the distribution of the pits described by the accumulated plastic strain, which was assumed in three cases, the delta function and Gaussian and Weibull distributions, the Weibull distribution agreed well with the experimental results, particularly including relatively large pit size. Furthermore, the Weibull distribution reproduced the depth profile, i.e., pit shape, better than that using the delta function or Gaussian distribution. It can be said that the cavitation erosion phenomenon is predictable by adopting the Weibull distribution. This prediction method is expected to be applied to predict the cavitation damage in fluid equipment such as pumps and fluid parts.