The method to treat randomly distributed spherical fuels in continuous energy Monte Carlo calculations has been established. In this method, the location of a spherical fuel is sampled probabilistically along the particle flight path from the spatial probability distribution of spherical fuels, called the nearest neighbor distribution. The necessary probability distribution was evaluated by a newly developed Monte Carlo hard sphere packing simulation code, which employs a random vector synthesis method to reduce overlaps of spherical fuels. The obtained probability distribution was validated by comparing a cross-section photograph of a real fuel compact and an X-ray diffraction experimental result. This method was installed in a Monte Carlo particle transport code and validated by an inventory check of spherical fuels and criticality calculations of ordered packing models. Also, an analysis of a critical assembly experiment was performed with the new code. As a result, it was confirmed that the method was applicable to practical reactor analysis. The method established is quite unique in the respect of probabilistically modeling the geometry of a great number of spherical fuels distributed randomly without any loss of the advantage of the continuous energy method.