Main

On 1 December 2020, Chang’e-5 successfully landed in Northern Oceanus Procellarum at (51.92° W, 43.06° N)15. The landing site is located in the geologic unit Em4 mapped by Qian et al.16 (see their fig. 1b), which belongs to the Eratosthenian-aged mare. The landing area is characterized as a young phase of volcanism16 with the basalt thickness varying from 39 to 63 m overlying the old Imbrian-aged mare units17. The thickness of the lunar regolith in Em4 is ~3–15 m (ref. 18). Before landing and sampling, our team identified and mapped craters >200 m in diameter in Em4 (Fig. 1b) and other units in the landing area19.

Fig. 1: Lunar sampling sites and mapped craters in the Chang’e-5 landing area.
figure 1

a, Distribution of all of the lunar sampling sites. b, Counted craters (marked in red) in the Chang’e-5 landing area. North is up.

Automated sampling was the core task of the Chang’e-5 mission. Two sampling methods were employed: surficial sample collection by a robotic arm and underground sample collection by a drill. The Chang’e-5 samples that have been radiometrically measured were scooped from the lunar regolith surface20. The samples were sealed in a storage device to ensure that they were unaffected by the external environment when returning to Earth. The samples were carried by the ascender and transferred to the returner in lunar orbit. The returner was subsequently separated from the orbiter and finally carried the samples back to Earth on 17 December 2020. The total weight of the returned samples was 1,731 g. The samples collected by the robotic arm have been radiometrically measured at the Institute of Geology and Geophysics, Chinese Academy of Sciences14,20.

Mineral and radiometric analyses were conducted at the Institute of Geology and Geophysics, Chinese Academy of Sciences, using 3 g of the surficial sample. Basalt clasts used for radiometric measurements have similar major mineral constituents and formed during the same final crystallization stage of the magma14. Radiometric measurements using the Pb–Pb dating technique were taken on both Zr-rich minerals and phosphates, and reveal the same age of 2.030 ± 0.004 Gyr ago (Ga), which represents the absolute age of the basalt within the landing area of Chang’e-520. This age has been used here to update the lunar cratering chronology function that is widely used at present.

We updated the chronology function established by Neukum13 using the the radiometric age and corresponding terrain N(1) value; that is, the total number of craters with diameters ≥1 km per square kilometre. The crater size–frequency distribution in Em4 has been studied previously and different N(1) values have been derived17,19,21,22. Here we use the result from Jia et al.19 (that is, N(1) = (1.74 ± 0.32) × 10−3 km−2), which was derived from crater measurements in a high-resolution orthophoto map of the Chang’e-5 landing area generated from more than 700 Narrow Angle Camera images from the Lunar Reconnaissance Orbiter Camera. In studying the size–frequency distribution with high-resolution images, Hiesinger et al.23 updated several N(1) values of the young units sampled. Based on the data used in Neukum13, the new data from Chang’e-5 and the updated data for Copernicus crater in Hiesinger et al.23 were used to fit the chronology function model. As shown in Fig. 2, the data point from Chang’e-5 does not change the chronology function model greatly from that of Neukum13; thus a function with an exponential decline and a linear rate is again used to fit the data presently available. Our new lunar cratering chronology function is as follows, where t is the model age:

$$N\left( {1,t} \right) = 1.089 \times 10^{ - 13}\left( {{\mathrm{e}}^{6.757t} - 1} \right) + 7.660 \times 10^{ - 4}t$$
(1)
Fig. 2: New lunar chronology model and comparison with the widely used Neukum model.
figure 2

a, The new lunar chronology model (solid red line) and the model from Neukum13 (dashed black line). b, The model age difference with respect to N(1) (our model − Neukum model). c, The model age difference with respect to the age from Neukum13 (our model − Neukum model).

Figure 2 also shows the model age difference (that is, new model − Neukum model) with respect to N(1). The new chronology model gives older ages for all N(1) values except 0.055 km−2 to 0.0065 km−2 (with corresponding model ages from 3.98 Ga to 3.59 Ga), for which the new chronology model gives slightly younger ages. The maximum difference appears at 2.55 Ga, with a corresponding N(1) of 0.0021 km−2, and the maximum model age difference is 0.24 Gyr. It is worth noting that the norm of the residuals decreases quickly with the iteration times in the least-squares fitting (see details in the Methods), demonstrating that the residuals have been stabilized with the iteration. A detailed comparison with other N(1) values and other chronologies is presented in the Methods.

The impact rate resulting in craters of 1 km in diameter and larger is given by the derivative of the impact chronology function with respect to t; that is, \(\phi = \frac{{\partial N}}{{\partial t}}\). From the new lunar cratering chronology model (equation (1)), the derived impact rate for diameter D ≥ 1 km is as follows:

$$\phi \left( {1,t} \right) = 7.358 \times 10^{ - 13}{\mathrm{e}}^{6.757t} + 7.660 \times 10^{ - 4}$$
(2)

Similar to the rate from Neukum13, the impact rate resulting in craters of 1 km in diameter and larger from the new model also exponentially decreases before ~3.0 Ga, and then an almost constant value is obtained. The difference in the impact rate resulting in 1 km and larger craters between the two models is also shown in Fig. 3. Before 3.866 Ga, the cumulative cratering rate from the new model is always smaller than that from Neukum13, with the difference becoming larger for older ages. From 3.866 Ga to ~3.156 Ga, the new model gives a slightly higher cumulative impact rate. From that time until today, the two models are almost same.

Fig. 3: Cratering rate and the difference between the new model and the Neukum model.
figure 3

a, Cumulative cratering rates with respect to model age for our model and the Neukum model. b, The difference in the cumulative cratering rates between the two models.

The differences in ages and cratering rates between the new and previous models are shown in Fig. 4. The models from Neukum13, Wagner et al.24, Le Feuvre and Wieczorek25 and our new model are generally close, with maximum differences of ~0.2 Gyr. These four models have significant differences from those of Marchi et al.26 and Robbins27, with maximum differences of up to 1.2 Gyr. So far, the model from Neukum13 has been the most widely used in lunar and planetary research. Our new model, which was updated from Neukum13 with the new radiometric age of the Chang’e-5 sample, should be more reliable and applicable in future research.

Fig. 4: Comparison of our lunar chronology model and previous models.
figure 4

a,b, Age differences between our model and previous models13,24,25,26,27 relative to N(1) (a) and age (b), respectively.

In summary, a new chronology model based on the radiometric age of Chang’e-5 samples is presented and compared with previous models. The radiometric age from Chang’e-5 samples falls within the gap between 3.0 and 1.0 Ga and serves as a ‘golden spike’ for lunar chronology. The result confirms that a summation of an exponential equation and a linear equation is applicable to describe the chronology model. The new model refines the previous model and will enhance the reliability of ages estimated with the crater size–frequency distribution method. We propose that the new chronology model is used in lunar surface dating, and that the chronology models for inner Solar System bodies (for example, Mars, Mercury and asteroids) are revised on the basis of this new lunar model while accounting for differences in their impact rate, impact velocity, gravity and target properties compared with the Moon28,29,30.

Methods

Crater counting has long been used to determine the model age of units of the lunar surface31. The rationale of the approach is to fit the observed crater size–frequency distribution of a given surface unit to a known crater production function1,32, which is further used to derive the absolute model ages along with a chronology function chronology function calibrated to radiometric dating from lunar samples13,33.

The production function refers to the standard crater size–frequency distribution of the lunar surface, which was established by combining a number of crater size–frequency distributions of measured craters on geologic units of different ages and in overlapping crater diameter ranges10,13,33,34. The production function widely used at present for the Moon was derived by Neukum13, Hartmann34 and Neukum et al.10 based on the hypothesis that the craters are randomly distributed on the lunar surface. Possible asymmetries in the distribution of impact craters on the Moon were proposed and their effects on the production function were modelled in several studies24,25,35,36,37. However, these models are mutually inconsistent, given the asymmetries in crater distribution, and they are not necessarily widely accepted at present23.

The chronology function is usually expressed as the relationship between the radiometric ages of the returned samples and the cumulative frequency of craters larger than some particular diameter, such as 1 km, that is, N(1) in the corresponding sampled area. However, the values of N(1) were often calculated through extrapolation using the production function instead of direct measurements, thus the production function also influences the chronology function. The current widely used chronology function is based on the values of N(1) from Neukum13, and corresponding radiometric ages. Hiesinger et al.23 updated N(1) for young ages with the highest-resolution lunar orbital images so far and indicated that the updated values (especially for Copernicus crater) are more consistent with the chronology function from Neukum13. Therefore, N(1) values from Neukum13 and the updated N(1) for Copernicus crater are used in this study. However, the updated N(1) from Hiesinger et al.23 was calculated using the production function from Neukum et al.10. It is better to derive the values of N(1) using the same production function. The production functions of Neukum13 and Neukum et al.10 have a similar relationship with different coefficients. If the diameter of the craters used to derive the N(1) is D1, the ratio of N(1) values from the two production functions can be expressed as follows:

$$\frac{{N\left( 1 \right)^ \ast }}{{N\left( 1 \right)}} = 10^{\mathop {\sum}\nolimits_{i = 1}^{11} {\left( {a_i - a_i^ \ast } \right)\left( {{{{\mathrm{log}}}}_{{{{\mathrm{10}}}}}D_1} \right)} ^i}$$
(3)

where N(1)* is from Neukum13 and N(1) is from Neukum et al.10; the coefficients of ai and \(a_i^ \ast\) are for Neukum et al.10 and Neukum13, respectively. Supplementary Fig. 1 shows the above equation for crater diameters ranging from 10 m to 300 km, which is also the recommended diameter interval for the two chronology functions. It is clear that for crater diameter less than 1.0 km the values of N(1) from the two models are very similar. For Copernicus crater in Hiesinger et al.23 (See Fig. 10 therein), the above ratio is between ~1.01 and 1.02. Thus N(1) is also very similar for the two production functions and was directly used here, as the accurate ratio is unknown and the difference plays a small role in generating the new chronology function.

Through the distribution of all of the radiometric ages, including that from the Chang’e-5 samples, the chronology model can be expressed as the summation of an exponential equation and a linear equation. Thus a function similar to Neukum13 was adopted. The trust-region-reflective algorithm38,39 of the nonlinear least-squares fitting was used to solve the coefficients of the chronology function by minimizing the sum of squared residuals through iteration. The uncertainties of the ages and N(1) values were not incorporated in the model fitting for simplicity and easy comparison with previous models. The coefficients of the Neukum model13 were used as the initial values, and the iteration was executed 1,250 times, when the norm of the residuals stabilized at the level of 1.42242 × 10−5.

In addition to the N(1) value provided by Jia et al.27 of (1.74 ± 0.32)×10−3 km−2, Wu et al.22 and Qian et al.17 gave N(1) values for Em4 of 1.24 × 10−3 km−2 and 1.28 × 10−3 km−2, respectively. The production function of Neukum et al.10 was used for both studies. If the above values are directly used to fit the new chronology function, as we have done here, the results are as follows:

$$N\left( {1,t} \right) = 1.223 \times 10^{ - 13}\left( {{\mathrm{e}}^{6.729t} - 1} \right) + 7.363 \times 10^{ - 4}t$$
(4)
$$N\left( {1,t} \right) = 1.218 \times 10^{ - 13}\left( {{\mathrm{e}}^{6.730t} - 1} \right) + 7.378 \times 10^{ - 4}t$$
(5)

Equations (4) and (5) were derived with the N(1) of Em4 unit from Wu et al.22 and Qian et al.17, respectively. Comparing the two equations with equation (1), the differences are not large. Supplementary Fig. 2 shows the differences between equations (1), (4) and (5) and the chronology function from Neukum13. These curves have similar shapes; that is, the maximum difference appears at a model age of about 2.52 Ga, and the modifications are small at young and old ages for all three models. The differences between the three models result from the different N(1) values in the Chang’e-5 landing area. This indicates the importance of utilizing high-resolution images and the subjectivity of identifying craters, including the careful identification and elimination of secondary craters when applying the crater size–frequency distribution method to estimate the age of a lunar surface unit.

Besides the model from Neukum13, there are other chronology function models24,25,28,29 and their differences from the new model in this research are shown in Fig. 4. In summary, the largest differences in the ages and cumulative cratering rates are between the new model and that from Robbins27; the new model is most similar to that from Neukum13.