It was a pleasure to read this well-documented review and synthesis of the Quaternary Milankovitch theory. The author's objective is to provide a historical overview of the understanding of Milankovitch forcing of glacial cycles over the Quaternary, and provide a conceptual model that is based on CLIMBER to outline the main questions associated with the astronomical control of ice ages, provide plausible solutions, and outline  a  future research programme. 

Overall, the text reads particularly well, with a good historical context and as far as I could tell most of the relevant references are cited. This is a landmark paper, which beautifully fits the tradition of Milankovitch medal lectures. 

Perhaps the two main controversial aspects are the title, and a few technical aspects, which I would encourage the author to consider, but from my point of view do not require major revision. 

The approach for the current review is therefore to discuss the title first, then proceed line by line mixing editorial comments with more scientific ones, and consider the more sensitive technical aspects in the end. 

General considerations

----------------------

It is relevant to consider the science reviewed in the present contribution as a development of Milankovitch's effort, because Milankovitch's purpose was to provide a theory of ice ages. Yet, I would argue,  this is not the only possible way of 'generalizing' Milankovitch's theory.  For example, a stream of modern research is focused on understanding orbital control of climate for other periods in the past, down to the Paleozoic. This other form of generalisation would give more importance to other climatic and environmental effects of the orbital forcing than their effects on ice sheets, such as ocean circulation and nutrient supply. Therefore, I would suggest to be a bit more explicit about the Quaternary ice ages in the title (see suggestion below). 

The word 'theory' is used in the present contribution to designate the ensemble of theoretical and empirical considerations that leads to our current understanding of the astronomical control on ice ages. The "theory" is illustrated by the "Model 3" conceptual model. The latter has obvious similarities with previously published models (especially from the Paillard's school) but the author explains that it is formulated such as to be compatible with CLIMBER-2 output, along with  general considerations about the bi-stability of ice sheets at the global level. The key point, here, is that the "Model 3" does not explicit introduces a 80ka time scale as Paillard 99 did, which in effect, puts more responsibility on the succession of precession cycles, modulated on eccentricity, to generate the 100ka signal. This definitely is an interesting point, and in that sense I can see why the author  presents this work as a 'generalisation' of Milankovitch's theory. Model 3 features a dynamics that wasn't in the original Milankovitch proposal. Yet, let us admit,  many authors of models and theories over the last 3 decades have introduced a dynamical component, and in that sense generalize Milankovitch's theory. Furthermore, as  I will discuss below, there are still some uncomfortable theoretical hurdles. For example, the linear relationship the derivative of ice volume (in m3/s) and the orbital forcing (in W/m2) is not straightforward to explain from the theory of ice sheets. 

In summary, I would concede that it is adequate to use the word "theory" in the title but  I would be more explicit about the focus on Quaternary glacial cycles and on the fact that this contribution is part of a long-term development. Here is my proposal: 

A proposal: "steps towards generalizing the Milankovitch theory of Quaternary ice ages"

I leave it with the editor, based on the comments of the other reviewers. 

Line by line comments

---------------------

p. 2: 50: "The comprehensive theory cannot be simple". It is a nice punchy sentence but perhaps a bit tautological. In fact, can we ever have 'a' general comprehensive theory of ice ages ? Perhaps this would be similar to ask for a general theory of the blood circulation. One can convincingly sketch the general idea (this is one objective of "Model 3"), and then go into endless details about the dynamics of fluids, chemistry, function and shape of heart etc. There will never be "one" theory because we speak about a complex system, which can be attacked from many angles. 

My suggestion (first draft): "Despite significant progress in understanding climate dynamics and [...] glacial cycles, questions have remained about the mechanism of importance of glacial dynamics, carbon cycle, and scaling relationships  between different variables, and how they concurred to synchronise glacial cycles on short eccentricity cycles. These are the targets for progressing towards better theories of glacial cycles. "

Somewhere in the text I would also see adequate to pay some tribute to the visionary paper by MacAyael 1979. 

l. 125, typo (starting with ".")

Section 3.2 : Perhaps this is a place to note a difference of status between different forms of low-order models. Paillard 1999 is clearly "inductive" in the sense that the existence of three states, plus the associated time scales, emanate from inspection of the data. Paillard 1999 summarises what he sees in the form of a small model, and then draws the consequences of it. Verbitsky et al. 2018 start from theoretical considerations about glacial scaling laws and deduces the ice volume trajectory, with obviously a bit of (reasonable) fine tuning. Saltzman's work from the late eighties /early nineties are somewhere in between, depending on how one looks at it. These different models have different functions, all useful,  in the construction of our understanding of the astronomical control of ice ages. 

l. 296: the word "massive" may be unnecessary. 

ll. 320 - 322 : CLIMBER-2 itself was tuned; it implicitly includes observations, in that sense it is not quite certain that CLIMBER-derived constraints should be considered as fully independent. 

l. 404: "erasing the memory". :this is correct, but the phenomenon is already implicitly there in MacAyael 1979. 

l. 449: "Model 3 represents" -> "Model 3 is"

l. 481: The author may also consider a reference to the Huybers - Tziperman 2008 paper. 

l. 605: I must admit having been unconvinced about the "Quantum tunnel" analogy (in which case the potential barrier is crossed by a form of delocalisation). For the modelling of ice ages, the basic limitation of the potential barrier image is that it is 1-dimensional, while the dynamics for the relaxation imply at least another degree of freedom. This is not what happens in quantum dynamics. I leave the author with these considerations without any intention to fight on this point. 

ll. 650 - 655 : Basal sliding related to thermal balance (as encoded by Verbitsky et al. 2018) is another potentially important mechanism not mentioned here. From informal conversations with glaciologists, I understand that this is a very plausible, major ingredient for the deglaciation catastrophe. 

p. 24 before section 5.8 : The author may consider adequate to mention already at this point the Paillard-Bouttes theory about sequestration and release of carbon due to change in the formation of AABW (it is alluded to later in the text, admittedly)

p. 26 ll. 820 - 824 : regolith over CO2. Fair point, but pre-800 ka CO2 estimates remain uncertain, especially that one needs a trend over the mean state (set by the balance between outgassing and weathering), which may differ from what individual glacials or interglacials estimate. 

l. 866 : The long interglacials indeed bring an interesting constrain in the decision between 'self-sustained' vs 'driven' constraints and this section is relevant. However, what model simulations give in this respect (l. 867) may bring a tautological argument, especially when it comes for future climate simulations because the anthropogenic CO2 emissions would have broken the self-sustained oscillation, if there was one, anyway. 

sect. 5.11 : The discussion in this paragraph is reasonable, but I would like to use this opportunity to clarify one point. Several models with self-sustained oscillations have indeed sensitive dependence on parameters and/or additive stochastic noise, in the sense that certain terminations may be triggered one precession cycle in advance, or delayed, delaying the whole sequence. Mathematically, this indeed occurs as a manifestation of 'non-chaotic strange attractors' (Mitsui and Aihara, 2014; Crucifix, 2013). But in all cases, these models do not display sensitive dependence to _initial conditions_. I agree that the surprising efficiency of 'simple rules' suggest that the timing of deglaciations is less sensitive to details or stochastic elements than these models may suggest, but on the other hand, I observe that the catastrophic character of some deglaciations (especially termination V) is hard to capture by those models which are the most well-behaved (like CLIMBER). So we have to explain a paradox here: On the one hand terminations would be highly catastrophic (which indeed suggest a domino effect, in essence very sensitive to an initial trigger), and on the other hand their timing would be very robustly set by the astronomical forcing. 

l. 960: it is quite clear what the author means here, but perhaps the semantics could be polished ( "simply by a non-linear response" , "an arbitrary non-linear response", "a very special type of non-linearity"). Perhaps there is a way to nail it a bit more explicitly. The key point is that a linear response (or weakly non-linear) response to a non-linear transformation of the orbital forcing will generate all the eccentricity spectrum, including its 100ka, 400ka and even longer components, simply because it will merely rectify the signal. What we need is a strongly non-linear internal dynamics, that is, internal feedbacks (or, to put it otherwise: dynamical mechanisms) which are triggered by the internal state of the system. Perhaps one could contrast 'non-linear dynamics' to 'non-linear response' but I concede that this is not quite satisfactory either. 

About Model 3

-------------

There is one fundamental difficulty with many conceptual Quaternary ice sheet models, and, admittedly, many low-order ice sheet models (include from this reviewer, see Martínez Montero et al. 2022): why would the ice volume derivative be proportional to its volume ? 

Indeed, if we admit that mass balance is area times net accumulation rate, it would perhaps be natural to except the net accumulation rate to be proportional to climate factors, which themselves are excepted to be proportional to the _area_ (not volume) of the ice sheets. 

These considerations concur to the scale relationship of Verbitsky et al. 2018 (dS/dt proportional to S^(3/4)*accumulation, S the area). Yet, I agree, the linear relationship between bulk accumulation and volume works well in simple models, and this is perhaps why, besides its simplicity, it is so popular. The problem is briefly alluded to in Verbitsky and Crucifix, CPast 2023 (original version in https://doi.org/10.5194/cp-2023-30). At least there is a little theoretical challenge, here. 

Interglacial metrics

--------------------

The difference between eq. A1 and the original Huybers metric (l. 1034) is n*I_0, where n is the number of days with I>I0. We indeed expect this latter term to be a function of obliquity, which then cancels part of the obliquity component in the original Huybers metric. Both the Huybers original formulation, and the Milankovitch caloric insolation, work in similar ways. Huybers and Tziperman 2008 are quite explicit about why they do not clip the insolation (just after their equation 2); I would like to bring the additional arguments: 

 - when we calibrate one of the good old Saltzman's models we obtained posterior distributions with about equal weights of precession and obliquity (Carson et al., 2013, Figs. 3 and 5, compare gamma_P and gamma_E)

 - in LOVECLIM, with interglacial conditions, precession and obliquity have more or same-order-of-magnitude effects on the "GDD" (but actually equivalent to PDD) at high latitudes (Bouncer et al., 2015, Fig. 7) with, I would concede, more precession in North America. 

So the discussion clearly has merit and the author is undoubtedly right in tackling this issue, but perhaps the conclusion lines 1046 deserve some caveat. 

The Domino effect

-----------------

Such a domino effect is indeed what is encoded, more or less explicitly, in many conceptual models, but how 'irreversible' or 'catastrophic' the domino needs to be is not straightforward. Some deglaciations are deeper than others, or sometimes stalled (those leading to 7e, or the strange 15c /15 a duet) which conceptual models tend to overdo. 

References

----------

MacAyeal D. (1979), A catastrophe model of the paleoclimate, Journal of Glaciology, (24) 245-257 doi:10.3189/S0022143000014775

Martínez Montero M., M. Crucifix, V. Couplet, N. Brede and N. Botta (2022), SURFER v2.0: a flexible and simple model linking anthropogenic CO2 emissions and solar radiation modification to ocean acidification and sea level rise, Geoscientific Model Development, (15) 8059–8084 doi:10.5194/gmd-15-8059-2022

Carson J., M. Crucifix, S. Preston and R. D. Wilkinson (2018), Bayesian model selection for the glacial-interglacial cycle, Journal of the Royal Statistical Society: Series C (Applied Statistics), (67) 25-54 doi:10.1111/rssc.12222

Bounceur N., M. Crucifix and R. D. Wilkinson (2015), Global sensitivity analysis of the climate-vegetation system to astronomical forcing: an emulator-based approach, Earth System Dynamics, (6) 205--224 doi:10.5194/esd-6-205-2015

This manuscript presents the author’s current understanding of Quaternary ice age cycles based on a discussion of several conceptual models, on more elaborated models (in particular Climber2) together with a thorough discussion of the most important associated physics. It is a very interesting and well-written paper and I liked it very much. I have therefore no hesitation to recommend publication in Climate of the Past. I only have a few minor comments that the author may use in his final version before publication.

1 - I do not agree with all the views presented here by the author. I certainly do agree on the main mechanisms behind the 100-kyr-cyclicity. But I don’t like the so-called “regolith hypothesis”: now, even its first promoter (Peter Clark) explains why it was not a good idea after all. In the paper, it is clearly presented as a hypothesis, which fills our lack of knowledge on the origin of long-term trends: may be, for completeness, the author could state that this remains a controversial hypothesis.

2 - Figure 14 is not referenced in the text. Besides, in my opinion, it is entirely useless and does not help the reader to understand the paper. I would suggest removing it.

3 - When discussing the speed of the forcing behind the MPT, ie. the rate of change of the critical ice volume Vc, around lines 815-820: an interesting paper was written on this specific point by Legrain, Parrenin and Capron (Nature communications Earth & Environment, 4, 2023) using a rather similar “threshold based” model; with the conclusion that a gradual change appears more likely when using random parameters.

4 - Equation (2) line 300. A minus sign is missing for k=2 (v should decrease during terminations).

5 - line 574. I am not convinced at all that the 100-kyr-cycle is a “peculiar regime”. Such a 100-kyr periodicity appears in many different contexts and in many pre-Quaternary Earth’s paleoclimatic records. For instance, in Pälike et al. (Science 2006) there is a clear 100-kyr cycle in the ¹⁸O that might be linked to Antarctic ice-sheet variations. Of course, interpretations are more difficult for these earlier periods, but the Quaternary is certainly much too short a time span to talk about “peculiar” or “ordinary” regimes.

6 - line 560: « phase locking  has a different meaning … ». Of course, “model MiM” and “model3” are different, but in both cases the locking is directly linked to threshold crossing. It is not clear to me why « phase locking… » should have a different meaning.

7 - Figure 8e: there is a single very long “cycle” of about 150-kyr in the histogram. I am wondering which one it is… and how is it possible with this model: a few details on this particular cycle could be helpful.

8 - Line 350: “Climber-2 has a problem simulating timing of TV while Model 3 does not”… I find this quite interesting! In the author’s opinion, is it pure chance? Or could conceptual models be “in some way” more robust than physical models?

9 - some typos:

Figure 12: Qusi-linear ->Quasi-linear

Legend Fig 7: artifitial -> artificial