Human Health Risk Assessment (HHRA) is a widely applied method to make decisions about the environmental status of sites affected by toxic substances. Its conclusions are affected by the variability and uncertainty of the input variables in the HHRA model. The aim of this work is to apply an algorithm based on 2D Monte Carlo simulations to integrate the variability and uncertainty of exposure factors, concentration, and bioaccessibility, reported by various information sources, to assess and compare their influence on the risk outcome. The method is applied to a specific case study of exposure of children to arsenic from accidental soil ingestion in a residential setting in the city of Madrid (Spain) by combining information from 12 studies. The consideration of the variability and uncertainty of the exposure parameters in the Baseline Risk Assessment (BRA, deterministic) resulted in a greater reduction in the numerical value of risk estimations than that produced by considering only the bioaccessibility factor. The results of the Probabilistic Risk Assessment (PRA) showed that the risk distribution was more sensitive to the variabilities of the accidental soil intake rate and the total arsenic concentration than to other variables such as bioaccessibility. In this case study, the uncertainty introduced by using the "default" reasonable maximum exposure factors in the HHRA model and the variability of the concentration term produce overestimates of risk that are at least in the range of those produced by omitting the bioaccessibility term. Thus, the inclusion of bioaccessibility is, alone, insufficient to improve the HHRA since the selection of the exposure factors can significantly affect the estimates of risk for the soil ingestion pathway. In other sites or for other contaminants, however, the role of the uncertainties associated with the bioaccesible fraction could be more pronounced. The method applied in this work may be useful in updating exposure factors to reduce uncertainties in HHRAs.