-
The coupled Yarkovsky–YORP evolution model is not specified with sufficient mathematical and algorithmic detail to be assessed or reproduced (Sec. 2.6–2.7; Sec. 3.3). In particular, the manuscript does not provide the explicit forms of ${\rm d}a/{\rm d}t$, ${\rm d}\omega/{\rm d}t$, and ${\rm d}\epsilon/{\rm d}t$ (or equivalent), whether diurnal/seasonal components are included, how heliocentric-distance scaling enters, how diameter/density/albedo/thermal inertia enter, and how numerical integration is performed (time step, scheme, random seeds, number of realizations). Closely related: $C_{\rm yark}$ and $C_{\rm yorp}$ are not defined mathematically (units/dimensions, baseline normalization), so reported best-fit magnitudes cannot be audited or compared to prior work.
Recommendation: In Sec. 2.6, add an explicit “Model equations and parameters” subsection that (i) writes the governing evolution equations used for $a(t)$, $\omega(t)$, and $\epsilon(t)$, with citations; (ii) defines $C_{\rm yark}$ and $C_{\rm yorp}$ precisely (including units or explicit statement that they are dimensionless scalings relative to a stated baseline model), and lists all symbols with dimensions; (iii) states whether diurnal and seasonal Yarkovsky are both modeled and how thermal parameters are represented (explicitly or absorbed into $C_{\rm yark}$); and (iv) documents the integration method, timestep, duration (= family age), random seeding, and number of Monte Carlo realizations per grid point. Provide the full ($C_{\rm yark}$, $C_{\rm yorp}$) grid ranges/spacing in a table and ensure Sec. 3.3 and figure captions reference these exact definitions.
-
Definition of the family center $a_c$ and the derived variable $\Delta a$ may be biased and is not aligned with standard family analyses unless carefully justified (Sec. 2.3; Sec. 3.1). The manuscript sets $a_c$ as the mean of member semimajor axes and uses $\Delta a_i = a_i - a_c$, but the mean can be shifted by Yarkovsky-evolved tails, truncation by resonances, and interlopers; and it is unclear whether proper or osculating elements are used. Since $\Delta a$ is central to every correlation and to the simulation matching, an inconsistent anchor can systematically distort both the observed slopes and inferred efficiencies.
Recommendation: In Sec. 2.3, explicitly state whether $a$ is proper or osculating and justify the choice. Replace or benchmark the mean-based $a_c$ against at least one robust alternative (median/sigma-clipped mean; using only large $D$ objects less affected by drift; or adopting published family centers from a standard family catalog). Perform a sensitivity test for a subset of key families (including one with apparent anomalies, e.g., Gefion) showing how (i) $\epsilon$–$\Delta a$ correlations, (ii) $D$–$|\Delta a|$ relations, and (iii) best-fit ($C_{\rm yark}$, $C_{\rm yorp}$) regions change under alternative $a_c$ definitions and/or outlier clipping. Report the outcome explicitly in Sec. 3.1/3.3.
-
Signed $\Delta a$ and $|\Delta a|$ are used inconsistently in interpreting size-dependent drift (Sec. 3.1; related figures across Sec. 3). The physical expectation from Yarkovsky scaling is primarily a relationship between size and drift magnitude ($|\Delta a|$ or V-shape envelope), while the sign of $\Delta a$ depends on spin-axis orientation (e.g., $\cos \epsilon$). Several places motivate a diameter–$\Delta a$ correlation with arguments that actually correspond to diameter–$|\Delta a|$, which risks confusion and can generate misleading null/positive correlations if both drift directions are present.
Recommendation: Audit Sec. 3.1–3.2 and all relevant plots to state explicitly, for each size–drift test, whether $\Delta a$ is signed or absolute. For the core “size dependence of drift” claim, prioritize $D$ vs $|\Delta a|$ (and/or a V-shape boundary method) and report those results consistently. If any signed $D$–$\Delta a$ analysis is retained, provide a clear physical rationale (e.g., asymmetric obliquity distribution, resonance removal on one side, selection effects) and demonstrate that the inference does not hinge on this choice.
-
The statistical inference and goodness-of-fit framework used to select best-fit ($C_{\rm yark}$, $C_{\rm yorp}$) is not defined in a unified, reproducible way (Sec. 2.7; Sec. 3.3). It is unclear exactly which 1D distributions are compared with which tests, how $\chi^2$ is computed (binning/normalization), whether any multi-dimensional comparison is attempted, how regression-slope comparisons enter, and what scalar objective function determines “best”. Use of KS $p$-values as an optimization target (especially a min-$p$ rule) is ad hoc and sensitive to sample size; additionally, near-flat contour maps in some families could indicate weak identifiability, plotting-range issues, or an implementation bug.
Recommendation: In Sec. 2.7, define a single objective function used to rank grid points (or a clearly specified multi-criteria decision rule), and justify it. Specify: (i) which observables are compared ($\Delta a$, spin period, obliquity; possibly transformed variables like $\cos \epsilon$); (ii) the exact test statistics used, including binning/smoothing choices for any $\chi^2$; and (iii) how statistics are combined into one score. Strongly consider adding or substituting a distance metric better suited for optimization (e.g., Wasserstein distance for 1D; energy distance/MMD for 2D like $(\Delta a, \cos \epsilon)$). For each family, quantify identifiability by showing example simulated distributions at several grid points and report whether the objective varies meaningfully; if contours are flat, state that constraints are weak and explain why (data sparsity, degeneracy, or model insensitivity).
-
Measurement uncertainties and observational selection effects—especially for obliquity—are not incorporated, risking biased correlations and biased ($C_{\rm yark}$, $C_{\rm yorp}$) inference (Sec. 2.1–2.5; Sec. 3.1–3.3). The requirement of complete property sets likely selects larger/brighter objects with dense photometric coverage, and the availability of spin-axis solutions is strongly non-random. Family ages also have non-negligible uncertainties that directly trade off against inferred drift efficiencies.
Recommendation: In Sec. 2.1–2.2, provide data provenance and typical uncertainties for $D$, $P$, $\epsilon$, and ages, with citations. For each family, report: total catalog members vs those used after requiring ($a$, $D$, $P$, $\epsilon$, taxonomy), and summarize how the used subset differs in $D$ (or $H$), inclination, and $a$. In Sec. 2.4–2.5 and Sec. 2.7, propagate uncertainties via bootstrap/Monte Carlo resampling (including perturbing ages within quoted uncertainties) to produce confidence intervals on correlations, regression slopes, and best-fit parameter regions. At minimum, repeat the key $\epsilon$–$\Delta a$ and model-fitting steps on a diameter-limited subset where spin-axis completeness is closer to uniform, and discuss any changes in Sec. 4.
-
Key quantitative results are not comprehensively reported across all 19 families, limiting robustness assessment and reproducibility (Sec. 3.1–3.4). The text highlights examples, but there is no single table listing, per family, $N$ used per analysis, correlation coefficients (with uncertainty), regression parameters, and the inferred best-fit ($C_{\rm yark}$, $C_{\rm yorp}$) with uncertainty/acceptable ranges and fit-quality metrics.
Recommendation: Add summary tables (Sec. 3.1–3.3): (i) per-family sample sizes for each observable ($N$ with $D$; with $P$; with $\epsilon$), and Spearman/Pearson (or preferred) correlations for $\epsilon$–$\Delta a$ (or $\cos \epsilon$–$\Delta a$), $D$–$|\Delta a|$, and any spin-related relations, with $p$-values and bootstrap confidence intervals; (ii) per-family best-fit ($C_{\rm yark}$, $C_{\rm yorp}$) with uncertainty regions (e.g., $\Delta$objective threshold), alongside the associated goodness-of-fit values. Ensure Sec. 3.3 and Sec. 4 base comparative claims (C-type vs S-type; young vs old) on these tables rather than selective examples.
-
The physical interpretation “constraints on thermal properties” is not quantitatively supported because the mapping from ($C_{\rm yark}$, $C_{\rm yorp}$) to thermal inertia/conductivity/albedo/roughness/density is left qualitative and is likely degenerate (Sec. 3.3–3.4; Sec. 4). Without an explicit thermophysical model and stated priors/assumptions, a constrained efficiency coefficient does not uniquely imply thermal inertia; similarly, $C_{\rm yorp}$ is strongly shape/roughness dependent.
Recommendation: Either (a) explicitly map inferred ($C_{\rm yark}$, $C_{\rm yorp}$) ranges to thermal inertia (and/or conductivity) ranges using stated formulae and assumed parameters (density, albedo, beaming/roughness, shape distribution), including a degeneracy discussion; or (b) narrow claims to “relative efficiencies” and avoid direct thermal-property language. For a few benchmark families (e.g., Themis, Koronis, Vesta or other well-studied cases), provide an order-of-magnitude conversion with uncertainty bands and compare to published thermophysical constraints. Clearly flag where inference is model-dependent (Sec. 4).
-
Several key modeling/analysis choices that strongly affect long-term spin-state distributions are missing or un-justified: collisional/stochastic reorientation, YORP self-limitation/cycles, spin-barrier/fission and any reset rule, and the handling of resonances or dynamical depletion that can reshape $\Delta a$ distributions (Sec. 2.6–2.7; Sec. 3.3–3.4). Without these, apparent fits may be non-unique, and outliers (e.g., anomalous trends such as Gefion’s reported size–drift behavior) cannot be interpreted confidently.
Recommendation: In Sec. 2.6–2.7, explicitly list included vs neglected physical processes and justify omissions with timescale arguments (size- and heliocentric-distance dependent). If collisional reorientation and/or YORP cycling is omitted, discuss how this might bias obliquity clustering and spin-period distributions and thus inferred $C_{\rm yorp}$. In Sec. 3.4 or Sec. 4, add a limitations/alternatives subsection addressing resonance truncation, initial ejection-velocity fields, and dynamical diffusion; for families with anomalous behavior (e.g., Gefion noted in Sec. 3.1), add a short quantitative check (resonance proximity, asymmetry in $\Delta a$ distribution, interloper contamination, or age uncertainty sensitivity).