Impacts drive lunar rockfalls over billions of years

On airless bodies like the Moon, erosion was traditionally attributed mainly to space weathering, but high-resolution imagery has revealed widespread mass-wasting features such as granular flows, slides, slumps, creeps, and rockfalls. Rockfalls are identifiable by displaced boulders and associated tracks, which persist on the order of ~1.55–35 Ma, with an upper limit of ~150–300 Ma, indicating geologically recent activity. Prior regional studies proposed drivers including shallow and deep moonquakes, impact-induced shaking, and thermal fatigue. It has also been hypothesized that basin-forming impacts (e.g., Orientale, ~3.8 Ga) caused extensive seismic shaking that smoothed slopes and reduced mass wasting on ancient (pre-)Nectarian terrains, which are thought to have shallower slopes and minimal mass wasting. However, past work was limited to small regions and/or coarse sampling. The authors aim to produce a high-resolution global map of lunar rockfalls to evaluate their spatial distribution, drivers, and implications for lunar erosion activity by analyzing >2 million LRO NAC images with a CNN-based approach.

Previous studies mapped mass-wasting features in selected lunar regions (equatorial, polar highlands, mare, pyroclastics, permanently shadowed areas) and suggested rockfall drivers such as moonquakes, impact-induced shaking, and thermal fatigue. Apollo-era observations and sample analyses constrained rockfall track survivability (~1.55–35 Ma; upper limit ~150–300 Ma), tying observed tracks to recent processes. Basin-forming impacts were proposed to produce global seismic effects that steeply degraded slopes early, implying older terrains would now be shallow and largely free of mass-wasting features. Nonetheless, prior analyses were constrained by limited spatial coverage and coarse data, leaving the global distribution and dominant drivers uncertain.

The study applied a convolutional neural network (CNN) to detect rockfalls in >2 million images from NASA’s Lunar Reconnaissance Orbiter Narrow Angle Camera (LRO NAC). The network (M5) is a RetinaNet single-stage dense detector with a ResNet101 backbone, trained on 809,550 augmented rockfall image chips using open-source tools (Python, Keras, TensorFlow). On the testing set, M5 achieved average precision 0.89, recall 0.44–0.69, and precision 0.98–1.00 at 50–60% confidence; a 60% confidence threshold was selected to balance high precision with acceptable recall. Data ingestion used pre-calibrated, compressed NAC pyramid files (8-bit) to maximize resolution and efficiency. Image selection was automated through a customized JPL Moon Trek API, prioritizing high spatial resolution and solar incidence angles of 30–60°, while minimizing inter-image overlap and underlap using the Klee algorithm. Selected images were tiled, processed by the CNN, and detections were georeferenced using LROC metadata. Rockfall sizes were approximated from CNN bounding boxes calibrated against landmark boulders across varying resolutions (~0.4–1.5 m/pixel). Computational resources included 2 local machines and 9 Google Cloud AI Notebooks (30 CPUs, 15 GPUs: NVIDIA Titan Xp, GTX 1080Ti, Tesla K80, Tesla V100). In total, 240,401 NAC images were processed (~64 Tpixels; ~7.5 TB). Global analysis integrated detections with auxiliary datasets in QGIS: LRO WAC mosaic, a global DTM (GLD100), and a USGS global geologic map (May 2019) refined for consistency and reprojection (equirectangular, polar stereographic). Additional vector datasets included highland/mare boundaries, crater catalogs (1–5 km and ≥20 km), Copernican craters, lobate scarps, tectonic features (graben, wrinkle ridges, faults), volcanic edifices (shield volcanoes, domes), rilles, and vents. Overlapping geomorphic polygons were dissolved to avoid double counting. Rockfalls were tallied across eight geologic epochs (Copernican, Eratosthenian, Eratosthenian–Imbrian, Imbrian, Imbrian–Nectarian, Nectarian, Nectarian–pre-Nectarian, pre-Nectarian). Spatial densities were computed per 1°×1° or 2°×2° quadrangles for visualization and per km² for normalization. An empirical relationship between spatial rockfall density and age was fit with a power law using epoch-binned counts: y_rockfall = 2.6358 a^(-0.813), where a is age in Ma and y_rockfall is rockfalls per km²; a normalized form relative to Copernican density was also provided. Limitations addressed included NAC coverage gaps, polar illumination issues, image quality, CNN misclassifications/missed detections, and auxiliary map uncertainties.

• Detected 136,610 rockfalls globally between 80°N and 80°S; mean density ~2 rockfalls per 1°×1° quadrangle. 91% lie in highland terrains. Nearside/farside and hemispheric distributions are broadly balanced.
• Pronounced clusters occur across terrains, with asymmetric intra-crater distributions in Tsiolkovsky, Atlas, and Plato. Peak density: 272 rockfalls per 1°×1° quadrangle in Bürg crater. Other dense clusters include Mare Orientale basin and Montes Rook/Cordillera, and craters Tsiolkovsky, Pasteur D, Atlas, Aristoteles, Milne N/L, and Crookes. Clusters may indicate regions of recent seismic activity.
• Boulder diameters range ~3–25 m, most commonly ~7–10 m. Older regions host boulders with systematically larger volumes, potentially due to wider spacing of impact-induced radial fractures as craters erode.
• Impacts directly trigger rockfalls (impact-ejected boulders with radial tracks; small impacts on crater slopes displacing boulders) and, importantly, create persistent fracture networks that precondition slopes for rockfalls over geologic timescales. Other processes (tectonic, volcanic) may act as additional triggers but are secondary.
• Spatial context: ~84% of rockfalls occur in impact craters or basins; ~0.8% in graben/faults/scarps; ~0.4% on volcanic edifices/vents/rilles; ~14.8% on regular/unclassified slopes. Excluding 1–5 km craters, ~72.7% are in craters/basins and ~26.1% on regular slopes; the crater fraction is likely underestimated due to lack of a global <1 km crater map.
• Rockfall density in craters <25 km diameter is on average ~3× higher than in other geomorphic regions (~13E-3 vs ~4E-3 rockfalls/km²). No general increase near tectonic features was observed, suggesting either limited endogenic seismic triggering or insufficient fracture preconditioning in those contexts; localized exceptions include Vallis Alpes, Aristarchus plateau, and a pyroclastic vent west of Baily.
• Temporal/geologic distribution: ~46% of rockfalls are in Imbrian terranes; a substantial ~25% occur in ancient Nectarian and pre-Nectarian terrains (>3.8 Ga), contradicting the assumption that such terrains lack mass wasting. Rockfall spatial density peaks in young (Copernican) terrains and declines with age, with a ~50% reduction by ~1.5 Ga relative to Copernican levels.
• Track survival times (~1.55–35 Ma; maximum ~150–300 Ma) imply that observed rockfalls in (pre-)Nectarian terrains reflect active erosion within the last 300–~1.55 Ma, likely ongoing today. This suggests 3.8 Ga is insufficient to fully erase relief on airless bodies and has implications for other bodies (e.g., Mercury, Vesta).
• An empirical power-law fit relates rockfall spatial density to age: y_rockfall = 2.6358 a^(−0.813) rockfalls/km²; normalized to Copernican density gives relative densities ~1 at 500 Ma, ~0.5 at 1.5 Ga, and ~0.25 at 3.2 Ga.

The global map and statistical analysis demonstrate that impacts are the primary driver of lunar rockfalls. Impact events not only directly eject boulders and trigger rockfalls but also generate extensive fracture networks (radial, concentric, spall, shear) in cratered terrains that precondition slopes to fail over billions of years. The asymmetric rockfall distributions within certain craters (e.g., Tsiolkovsky) match inferred oblique impact trajectories, supporting the role of impact-induced discontinuity patterns. Fracture morphology and density also influence rockfall magnitudes, which are larger on average than on Earth. The detection of significant rockfall activity in ancient (pre-)Nectarian terrains indicates a long-lived legacy of early impacts, potentially amplified by a higher impactor flux during the Late Heavy Bombardment. These findings revise assumptions about the erosional maturity of ancient lunar surfaces, show that mass wasting persists, and provide a quantitative link between rockfall density and topographic age that can complement other age-dating methods. Additionally, rockfall clusters may pinpoint regions of recent seismicity and help localize ongoing tectonic or volcanic activity.

This study presents the first global, high-resolution map of lunar rockfalls, revealing 136,610 events and establishing that impacts and their long-lived fracture networks predominantly drive lunar mass wasting. Rockfalls occur across all terrains, including a substantial fraction in Nectarian and pre-Nectarian units, indicating that erosion has remained active from the late Copernican to the present. The empirical relationship between rockfall spatial density and age offers a complementary tool for assessing erosional state and relative ages of lunar topography. The results have broader implications for the evolution of airless bodies with similar impact histories. Future work should refine auxiliary datasets (e.g., global crater catalogs <1 km), improve polar coverage and illumination handling, enhance CNN recall while maintaining precision, integrate mechanical modeling of fracture networks, and apply these methods to other planetary bodies. Rockfall clustering patterns can be leveraged to investigate localized recent seismic or volcanic activity.

Key limitations include incomplete NAC coverage and polar illumination/shadowing issues leading to potential underdetection near the poles; variable image quality; CNN performance constraints (recall ~0.44–0.69 at 60% confidence, precision ~0.98), implying missed detections and rare false positives; potential spatial gaps from image selection and tiling; uncertainties and inconsistencies in auxiliary datasets (global geologic map classifications, incomplete tectonic/volcanic feature catalogs, and crater maps—particularly <1 km diameters); difficulty distinguishing endogenic seismic triggers from lack of impact-induced fractures near tectonic features; and uncertainties in rockfall track survival times that affect timing interpretations. The empirical age-density fit depends on the geologic timescale used, map accuracies, and the specific CNN (M5) performance and should be applied cautiously.