The intense seismic scattering seen in Apollo lunar seismic data is one of the most characteristic features, making the seismic signals much different from those observed on the Earth. The scattering is considered to be attributed to subsurface heterogeneity. While the heterogeneous structure of the Moon reflects the past geological activities and evolution processes from the formation, the detailed description remains an open issue. Here we present a new model of the subsurface heterogeneity within the upper lunar crust derived through a full 3D seismic wave propagation simulation. Our simulation successfully reproduced the Apollo seismic observations, leading to a significant update of the scattering properties of the Moon. The results showed that the scattering intensity of the Moon is about ten times higher than that of the heterogeneous region on the Earth. The quantified scattering parameters could give us a constraint on the surface evolution process on the Moon and enable the comparative study for answering a fundamental question of why the seismological features are different on various planetary bodies.