RUSSIAN JOURNAL OF EARTH SCIENCES, VOL. 14, ES1003, doi:10.2205/2014ES000537, 2014

Evolution and variability of the pleistocene ice ages: A new view

V. A. Bol'shakov1, Ya. V. Kuzmin2

1Faculty of Geography, M. V. Lomonosov Moscow State University, Moscow, Russia

2Institute of Geology & Mineralogy, Siberian Branch of the Russian Academy of Sciences, Novosibirsk, Russia



The evolution and variability of Pleistocene glaciations is studied by comparison of the benthic $\delta^{18}{\rm O}$ stack LR04 and the Orbital Climatic Diagram (OrCD). Taking into account that the OrCD construction is based on a very simple principal, it corresponds well to the LR04 data within the last 1.25 Myr. It is shown that along with cases of good evidence for climatic influence of joint action by precession and obliquity variations there are some concrete situations in the paleoclimatic record of the last 1 Myr when this evidence is absent. A shift in the main periodicity of glacial cycles (the Middle Pleistocene transition, MPT) took place at around 1.239 Myr ago. Some possible mechanisms of the MPT are considered, and we argue that parametric resonance is the underlying mechanism responsible for the MPT phenomenon.


Almost four decades ago, Hays et al. [1976] showed for the first time that at least 80% of global climatic fluctuations during the last 500 kyr are represented by systematic variations with periodicities of about 100, 41, 23, and 19 kyr. Exactly these time periods characterize the quasi-periodic changes of the Earth's orbital elements, namely eccentricity, obliquity, and precession. Fluctuations of these elements cause changes in insolation, i.e. the solar radiation which comes to the upper boundary of the Earth's atmosphere. At the same time a phase correspondence between insolation and climatic changes was found for separate frequency bands by Hays et al. [1976]. For example, the paleoclimatic component with 40-kyr periodicity was slightly behind in phase of obliquity (also "axial tilt" or "$\varepsilon$ angle") changes as expected, because of inertia in the climate system. At the same time, a minimal value of the $\varepsilon$ angle was found to correspond to times of cooling in climatic records. It is worth noting that these findings came more than one hundred years after Croll [1867] published the first proposed relationship between obliquity and climate, who noted in particular that a decrease in $\varepsilon$ will result in global cooling. Thus, in the mid-1970s the orbital hypothesis of paleoclimate, described originally in the nineteenth century [Adhémar, 1842; Croll, 1864, 1875], was confirmed. According to this hypothesis, global climatic changes (i.e. glaciations and interglacials) are connected with orbital fluctuations of insolation.

However, Hays et al. [1976] and other authors [see Berger, 1999; Bol'shakov, 2001, 2003a, 2003b; Imbrie and Imbrie, 1986; Imbrie et al., 1993; Muller and MacDonald, 1995; Raymo and Nisancioglu, 2003; Ruddiman et al., 1986] have revealed significant deviations between empirical data and the most widely accepted version of the orbital theory of paleoclimate created by Milankovitch [1930, 1941]. Among the best-known disagreements are "the problem of 100-kyr period" and "the problem of the Middle Pleistocene transition". There are numerous studies related to these issues but a consensus solution has not yet been reached.

A correct paleoclimatic theory is important not only as a means of understanding the past, but also for the prediction of future climatic change. A unique situation exists in the interpretation of Pleistocene climate, especially for the last 1–1.5 Myr. The ingoing signal of variations in insolation, which affects the Earth's climatic system, is well-known. The same is true for the outgoing signal which is represented by a plethora of empirical data derived from both deep-sea cores and long terrestrial records. These data allow us to reconstruct changes in global ice volume, temperature, sea-level changes, and atmospheric circulation, among other parameters. Consequently, when we have this information and take into account different feedbacks, it is possible to look for solution of task about the work of the Earth's climatic machine which transforms the insolation signal into global climatic changes. Obviously, such knowledge could help us in the creation of more precise forecast of global climate changes when considering anthropogenic impacts on the natural environment.

In this paper, the Orbital-Climatic Diagram (hereafter – OrCD) will be used for investigation of the above mentioned problems. The OrCD was proposed as a template of the global climatic changes during the last 1 Myr [Bol'shakov, 2001, 2003b].

The OrCD as an Alternative to Discrete Insolation Curves of M. Milankovitch and his Followers

Controversies in the Use of Discrete Insolation Curves for Paleoclimatic Interpretations and Modeling

It is well-known that the famous Milankovitch diagram [see Milankovitch, 1930, p. 154, Diagram I], widely used in first half of the twentieth century for interpretation of glaciations in the last 600 kyr, represents changes in caloric summer half-year insolation at 65° N. Such changes in insolation we call "discrete" because they are obtained not for the whole year and not for all latitudes on Earth. The curves created by followers of Milankovitch are even more discrete, as they represent changes of the monthly (usually for June or July) or even daily (for the summer solstice) insolation at 65° N.

It is important to keep in mind that variations in both obliquity and precession contribute to discrete insolation changes, and that direct input of the eccentricity variations in insolation is usually neglected in climate models due to its small value. It is also worth noting that the effect of precession on monthly or daily insolation values is much higher than influence from obliquity. This constitutes one of the most significant controversies in paleoclimate theory among the Milankovitch followers, because the largest by amplitude precession harmonics in the orbital insolation signal is subdued in paleoclimatic records [Bol'shakov, 2003a, 2008; Cleaveland and Herbert, 2007; Hays et al., 1976; Imbrie et al., 1993; Karner et al., 2002]. Discrete insolation curves, nevertheless, are widely employed for interpretation of paleoclimatic data and for modeling.

However, already in the nineteenth century A. von Humboldt and J. Herschel pointed out that in order to interpret paleoclimatic data and create a theory of climate changes, full annual variations of insolation for the entire Earth should be used and not discrete changes [see Bol'shakov, 2003b, 2008; Croll, 1875; Imbrie and Imbrie, 1986]. Generally speaking, the same approach is claimed by the followers of Milankovitch. For example, Imbrie [1982, p. 413] wrote: "There has also been a tendency for investigators to believe they could model the response of the system from a radiation curve representing the input at a single latitude and season [e.g., Broeker and van Donk, 1970; Kukla, 1968; Milankovitch, 1941]. Since no one could be sure which insolation curve, if any, was the crucial one, investigators had great freedom to choose a curve that resembled a particular set of data. Understandably the resulting ambiguity did much undermine confidence in the validity of time domain predictions. Starting in 1976, with the advent of numerical models that integrated the effect of insolation changes over all latitudes and seasons, this situation was much improved..." (italics are ours – V. B., Y. K.). Nevertheless, J. Imbrie himself [see Imbrie and Imbrie, 1980; Imbrie et al., 1993] used for modeling and paleoclimatic interpretations monthly (June) insolation at 65° N. A similar case can be found in Berger et al. [1998, p. 616]: "Such time-depended climate models must therefore be forced only by the astronomical variations of insolation for each latitude and day...", but on the next page of the same article they added: "June insolation at 65° N is very often used as a guideline for the analysis of climatic changes and, in particular, for ice volume changes" [Berger et al., 1998, p. 617]. Six years later a similar argument was made by Loutre et al. [2004, p. 2]: "A more general version of the astronomical theory is now widely used, especially in climate modeling, where changes in insolation at all latitudes and times of the year are taken into account. Nevertheless, it is often supposed that insolation at 65° N in June can be used for comparison with most proxy records." It remains unclear why instead of using "the astronomical variations of insolation for each latitude and day" for comparison with proxy records the discreet insolation at 65° N in June has been repeatedly employed. After all, the global whole year insolation variations are fundamentally different from the monthly insolation variations under the separate latitude, because the latter do not show the seasonal and latitudinal contrasts of insolation which forcing the climate really.

Fig 1
Figure 1

Based on the above cited papers, it follows that the interpretation of empirical paleoclimatic data in light of discreet insolation curves is a flawed approach. This observation is supported by the well-known disparities between the Milankovitch insolation diagram and empirical data. The same shortcomings also apply to monthly and daily insolation curves. In order to support this conclusion, one can turn to comparison between changes in July insolation at 65° N [Berger and Loutre, 1991] and global climatic fluctuations as reflected in the benthic $\delta^{18}{\rm O}$ stack LR04 [Lisiecki and Raymo, 2005] (Figure 1). It should be reminded that changes in benthic $\delta^{18}{\rm O}$ values reflect changes in global ice volume and, to a smaller degree, temperature fluctuations. The higher $\delta^{18}$O values reflect an increase in ice volume and/or a drop in temperature, and vice versa. In Figure 1, because ordinate axis $\delta^{18}{\rm O}$ is multiplied by $(-1)$ for convenience in comparing with insolation curve, $\delta^{18}{\rm O}$ graphic minima correspond to relatively cool periods (the largest of which are glaciations, with even Marine Isotope Stage [hereafter – MIS] numbers). The $\delta^{18}{\rm O}$ graphic maxima reflect relatively warm periods, and interglacials are assigned with odd MIS numbers.

By comparison of insolation and oxygen isotope (OI) curve, it is possible to see that significant insolation minima correspond to both glaciations and interglacials, and also to transitions between them; the same can be said about insolation maxima. The July insolation values comparable to the present time coincide with different global climatic conditions, from glaciation to interglacial. Changes in periodicity and amplitude of the $\delta^{18}{\rm O}$ record also reflect the phenomenon of Middle Pleistocene Transition (MPT), dated to ca. 1.24 Myr ago [Bol'shakov, 2013]. The MPT occurred in the absence of any significant changes in character of orbital insolation variations. The correspondence between changes of July insolation and $\delta^{18}{\rm O}$ variations is also absent before 1.24 Myr ago (Figure 1). In other words, quantitative changes of discreet insolation do not correspond to paleoclimatic changes. It should be added that the correlation coefficient between July insolation at 65° N and the $\delta^{18}{\rm O}$ stack LR04 for the last 1.25 Ma is 0.236. It grows to a maximum of 0.34 if insolation curve would be shifted in time by 3.5 ka ahead, i.e. if we accept the delay in climate response to insolation changes is 3.5 ka. However, this last correlation coefficient 0.34 is low too.

Thus, according to Berger et al. [1998, p. 616], and Loutre et al. [2004, p. 2] (see above), it is possible to conclude the following. For paleoclimatic interpretation of empirical data, modeling, and the creation of a diagram that qualitatively accounts for global paleoclimatic changes, it is essential to consider annual (i.e. for all seasons or months) and global (i.e. for all latitudes) insolation changes. It is also necessary to develop for each orbital element individual feedback mechanisms which transform insolation signals to climatic changes; it is shown, for example, by Bol'shakov [2003a, 2003b, 2008]. A comprehensive solution accounting for these essential ingredients has not yet been done completely due to following circumstances: 1) an incomplete knowledge of feedback mechanisms, where the most complicated issues are the influence of clouds and aerosols; and 2) an enormous amount of computing time necessary for calculations of paleoclimatic models in which all of the above mentioned factors are taken into account.

The Construction of OrCD

The difficulty in using a global annual insolation for all latitudes (in order to interpret paleoclimatic data) results from the following. First, the integral annual insolation value is equal to zero at any latitude for precession changes because changes in summer insolation are compensated by counter-phased fluctuations during the winter. This follows directly from the equations of Milankovitch [1930] for the calculation of changes in insolation of summer ($\Delta Qs$) and winter ($\Delta Qw$) caloric halves of year in the Northern Hemisphere, as expressed in equations (1)–(2):

\begin{equation} \tag*{(1)} \Delta Qs = \Delta Ws \times \Delta \varepsilon - m \times \Delta (e \times \sin w) \end{equation}

\begin{equation} \tag*{(2)} \Delta Qw = \Delta Ww \times \Delta \varepsilon + m \times \Delta (e \times \sin w) \end{equation}

Here "$\Delta Ws$" and "$\Delta Ww$" are changes in summer ($\Delta Ws$) and winter ($\Delta Ww$) solar radiation at the variation of the tilt of Earth's axis $\Delta \varepsilon$ at 1°, and they depend on latitude. The "m" is a factor that depends on the latitude; the "w" is the longitude of the perihelion of the Earth relative to the vernal equinox. Thus, the first terms on the right-hand sides of equations (1)–(2) represent the contribution of changes in insolation due to variations of the Earth's tilt, and the second terms are contributions related to precession. In order to obtain the annual insolation change $\Delta Q$, it is necessary to sum up the equations (1)–(2). In this case, the second terms on the right-hand sides will be mutually annihilated, i.e. the contribution of precession to the change in annual insolation is equal to zero. Changes in the annual insolation at certain latitude are determined only by variations of the obliquity $\Delta \varepsilon$:

\begin{eqnarray*} \Delta Q = (\Delta Ws + \Delta Ww) \times \Delta \varepsilon \end{eqnarray*}

Furthermore, precession-induced insolation changes are reversed in phase for different hemispheres. For example, when one half of the globe experiences precession-induced insolation effects characterized by short hot summers and long cold winters, the other hemisphere will experience long cool summers and short mild winters. Thus, the annual global precession signal is equal to zero. Secondly, if we consider the Earth as a sphere, the global insolation signal responding to obliquity variations will also be zero. From this perspective, we have no reason to expect direct global climate responses to global annual insolation signal equal to zero, and generally speaking, variations in precession and obliquity should not cause global climatic changes.

The correct way of solving this problem was described by J. Croll. He introduced for the first time the concept of positive feedbacks [see Croll, 1875, p. 13] which amplify weak insolation signals into global climatic changes. For example, he established a mechanism for the global climatic influence of obliquity variations, by taking into account the positive feedback related to albedo. Changes in the $\varepsilon$ angle result in counter-phased fluctuations of insolation in high and low latitudes of both hemispheres. When the $\varepsilon$ angle decreases, the annual insolation in high latitudes becomes smaller, while in low latitudes the insolation increases. When the $\varepsilon$ angle increases, the high latitude insolation rises, and the low latitude one declines. Because in high latitudes the albedo feedback is much more effective than in low latitudes due to the presence of large fields of snow and ice, global climatic changes connected to obliquity variations are determined by fluctuations in insolation in the high latitudes. This is why the decrease in $\varepsilon$ angle should be accompanied by a temperature drop and an increase of ice volume in the high latitudes of both hemispheres; in other words, lead to global cooling [Bol'shakov et al., 2012; Croll, 1867]. Taking into account the above-mentioned circumstances, Bol'shakov [2001, 2003a] proposed a simplified method for the construction of a diagram to represent global Pleistocene climatic fluctuations, referred to as the Orbital-Climatic Diagram (hereafter – OrCD). Implicit in the name of the diagram is the relationship between orbital insolation variations and Pleistocene climatic changes. In constructing the OrCD, variations of orbital elements and not of insolation were used. As this takes place, climatic influence of full annual and global insolation changes, connected to all three orbital elements, was determined based upon the most general and simple assumptions about the possible impact of full insolation variations on climate. This is important because it allows to interpret the influence of orbital variations on climate change in the most simple and straightforward way.

For example, the decrease of the $\varepsilon$ angle is expected to cause increase in ice volume and global cooling, according to Croll's mechanism (see above). It is assumed that a reduction in eccentricity of the Earth's orbit, which leads to decrease in annual insolation of the entire planet, also causes global cooling and expansion of ice volume. As for precession, a rise of ice volume and cooling in the mostly continental Northern Hemisphere should favor conditions of small insolation contrasts with long cool summers and short mild winters. This corresponds to a maritime type of climate favorable for growth of continental glaciers under decreasing temperature.

Because the influence of precession in different halves of the globe occurs with opposite phase, concurrent Southern Hemisphere seasonal insolation contrasts will be larger than those of the Northern Hemisphere, and long cold winters and short hot summers should result. There are reasons (see details in [Bol'shakov, 2003b, 2010]) to suggest that these particular conditions will facilitate an increase in ice volume and subsequently, cooling in the oceanic Southern Hemisphere. Thus, counter-phased changes of precession insolation in both halves of the globe could lead to similar climatic consequences because of physical-geographic differences between the two hemispheres. This fact can serve as one of the ways to eliminate the long-standing controversy over counter-phased changes of precession insolation in different hemispheres and sinphased climatic fluctuations in both of them. The last circumstance allows us to consider the OrCD as a standard for global climatic changes in the last 1 Myr; the OrCD was initially created to describe climate during this period.

Fig 2
Figure 2

The OrCD is constructed as the sum of variations of orbital elements calculated by Berger and Loutre [1991], normalized in relation to the mean value, and afterwards multiplied by "coefficients of climatic importance" (CCI). These coefficients were found by selection, given the condition that the resulting diagram should correspond well to existing paleoclimatic (OI) curves. This yielded the most suitable CCI values of 1.0, 0.7, and $- 0.55$ for variations in eccentricity, obliquity, and precession, respectively. The relationship of the selected CCI absolute values indicates that the largest relative climatic impact over the last 1 Myr comes from variations in eccentricity, and the smallest from precession, and this corresponds to empirical findings. The OrCD can be considered in terms of conventional relative probability ($\Delta P$) for glaciations (negative $\Delta P$ values) and interglacials (positive $\Delta P$ values) for the last 1 Myr (Figure 2a).

A diagram similar to the OrCD was compiled earlier by Imbrie et al. [1984] for the last 800 ka. It is called the ETP diagram because it represents the sum of normalized changes in eccentricity (E), tilt (T), and precession (P). Imbrie et al. [1984, p. 297] also wrote: "... we have reversed the sign of the precession index so that positive excursions in this core have the same climatic direction in the Northern Hemisphere as positive excursions in eccentricity and obliquity." Due to the same reason, the precession CCI was taken with negative value for construction of OrCD. It should be noted that in spite of using CCI there is no large visual difference between the OrCD and ETP diagrams for the time interval of 0–800 ka.

It is necessary to emphasize that the ETP curve was not used for paleoclimatic purposes; its creation and employment were determined by convenience for conducting a comparative spectral analysis [Imbrie et al., 1984, p. 297]: "Our object in performing cross-spectral analysis against the ETP curve, rather than against the individual orbital curves, is purely one of convenience, namely, to obtain a compact summary of orbital-isotopic relationships across the entire visible spectrum in a single diagram. Within each frequency band of interest, the results using ETP could be duplicated exactly by calculating a cross spectrum against the appropriate individual orbital curve". Thus, the formal difference between the OrCD and the ETP curves is the use of CCI for the former. However, this distinction stems from a principal divergence in interpretation. The ETP curve, unlike the OrCD curve, could not be used (and in fact was not used) for paleoclimatic interpretations [see Imbrie et al., 1984]. This is impossible, because according to Milankovitch theory and Imbrie et al. [1984] who follow it, the immediate climatic influence of variations in eccentricity is negligibly small compared to the orbital influences imposed by obliquity and precession. On the contrary, the OrCD is constructed with an assumption that changes in eccentricity produce the strongest and direct impact on global climate within the last 1 Ma.

The possibility of appreciable and direct global climatic influence from variations in eccentricity, along with the impact of obliquity and precessional changes, has been explained by structural differences of consequent insolation signals [Bol'shakov, 2003a, 2003b]. As noted above, the global annual changes of the insolation caused by variations in obliquity and precession are equal to zero. More than that, the precession variations are in counter-phase in different hemispheres. Only eccentricity variations result in annual insolation changes of the whole Earth, although these are relatively small.

In order to reconcile these facts, one should connect a mechanism of "non-linear amplification" just with global climatic influences from variations in precession and obliquity, because they do not change global annual insolation of the Earth (the incoming insolation signal of these elements is zero). On the contrary, the climatic impact of global eccentricity variations, which is larger than the zero, would be related to a mechanism of "linear amplification", with a relatively large coefficient of it. So, we consider mechanisms of amplification for each orbital element to be "individually non-linear" [Bol'shakov, 2003b].

Comparison of the OrCD and the Benthic ${\delta}^{18}{\rm O}$ stack LR04

A comparison of OrCD and changes in eccentricity and the OI curve LR04 for the last 1.8 Myr is shown on Figure 2. The coincidence for glaciations on the LR04 curve with minima in eccentricity is observed at least for the last 900 kyr. This was first demonstrated by Broeker and van Donk [1970] for the first ten MIS, and later by Hays et al. [1976] for the last 500 kyr. The only exception is the MIS 2–4 which constitute a single glaciation (Weichselian in Western Europe, or Valdai in Eastern Europe), and containing the MIS 3 interstadial warming. This knowledge allows us to distinguish glacial and interglacial stages on the OrCD as well. The Figure 2a shows primarily interglacial OrCD stages (hereafter – OrCDS) that correspond to the eccentricity maxima.

Within the last 900 kyr the OrCDS correlate well with traditional MIS (Figure 2b). The best match is observed for MIS 1–5, 7 (with deep MIS 7.4 minimum in the middle of MIS 7), 13–15, and 17–21. It should be noted that the maximum for OrCDS 19 corresponds to MIS 19, and has an age of about 770–790 kyr. This fits well to the age of Matuyama-Brunhes inversion as 780 kyr, determined with the Ar–Ar method [Cande and Kent, 1995; Spell and McDougal, 1992]. This paleomagnetic inversion occurs in the middle of the MIS 19 [Bassinot et al., 1994; Bol'shakov, 1999; de Menocal et al., 1991; Schneider et al., 1992; Shackleton and Opdyke, 1973; Tauxe et al., 1996]. The OrCD precedes in time the OI curve of global ice volume and temperature, because of inertia in the climatic system response to orbital forcing. This inertial delay has been estimated for well-pronounced and dated marks of the Late Pleistocene – MIS 2.2, 5.5, and 6.0. Radiometric dating methods allowed the determination of their ages as: MIS 2.2 (22–23 kyr); MIS 5.5 (122.6 kyr); and MIS 6.0 (129.3 kyr) [Imbrie et al., 1993; Thompson and Goldstein, 2006]. The ages estimated from the OrCD are: 27 kyr, 128 kyr, and 133.5 kyr, respectively. Thus, the inertial climatic delay can be estimated as about 5 kyr. This allows us to pinpoint chronologically the OI records by temporal correlation of typical extrema in the OI and OrCD curves. It can be easily seen that the OrCD corresponds much better to the LR04 curve than the discreet insolation curve (see Figure 1 and Figure 2). The correlation coefficient between OrCD and the LR04 for the last 1.25 Ma is 0.42. It grows to a maximum of 0.633 if the insolation curve is shifted 5.5 ka forward (i.e., into the future), which is close to the time lag of climate response to insolation changes as defined above. Nevertheless, there are some disagreements, and this should be expected if one takes into account the simplified method of construction the diagram of paleoclimatic changes which is used by us.

The most significant temporal discrepancies between the OrCD and LR04 curves for the last 1 Myr are observed for the deepest $\delta^{18}{\rm O}$ minimum related to MIS 8 (corresponds to the OrCDS 8), and for the $\delta^{18}{\rm O}$ minimum of MIS 22 and the related low OrCDS 22 value (see Figure 2b). The orbital signal forestalls the climatic response for 24 and 20 kyr, respectively. Therefore, by taking into account the 5 kyr temporal lag, it is necessary to assign minima for MIS 8 and 22 in the LR04 curve to 19 and 15 kyr, respectfully. This disagreement reflects differences in the approaches used for dating the OI curve by Bol'shakov [2003b] and Bol'shakov et al. [2005], and by Lisiecki and Raymo [2005].

Both studies employed in fact "orbital tuning" to assign a chronology; however, Bol'shakov's [2003b] method is considered here to be more secure because it takes into consideration not only variations of obliquity and precession, but also 100-kyr oscillations in eccentricity which are well-represented in the OI records of the last 1 Myr. The use of 100-kyr variations in eccentricity was done for the first time by Johnson (1982) and proved to be quite successful; it dated the Matuyama-Brunhes boundary to ca. 790 kyr ago eight years before the famous study by Shackleton et al. [1990]. Previously, the age of this event was determined as ca. 730 kyr ago on the basis of K–Ar dating.

Some differences in shape between the OI and OrCD curves are visible in Figure 2b. The most notable are significant maxima of the OrCD (indicated by stars) which are not detected in the LR04 record. Also, there are discrepancies in correlation of the largest OI signal of MIS 11 and the small maximum of the OrCDS 11 (this reflects well-known "MIS 11 problem"), and of the small intensity signal of OrCDS 6 and deep minimum of the MIS 6. As it is shown by Bol'shakov [2010], these disagreements are partly due to the fact that in the OrCD curve the 400-kyr eccentricity cycle is well-represented while it is absent in the OI records. The OrCD maxima indicated by stars (Figure 2b). corresponds to lower values or even minima in eccentricity. Using this observation, it is suggested that these maxima are connected with coherent changes in the two other orbital elements, precession and obliquity, which should lead to warming. However, the orbital changes are not reflected on the OI curve as significant warmings or decreases in global ice volume. So, it is evident that orbital impact does not influence climatic change in the way it would be expected in these particular cases.

This situation needs to be mentioned in particular, because it touches upon the important problem of establishment of the mechanism of climatic system response to orbital insolation forcing. Possibly, this is connected to peculiarities of processes in the climatic system itself during the examined time intervals, because there are cases when similar signals are reflected in the OI record as warmings. For example, there are OI maxima during MIS 3 and 23, and a peak inside MIS 18. All these three warmings have the same patterns: a) they correspond to minima in eccentricity; b) they reflect interstadial warmer conditions within glacial periods (unlike the generally accepted opinion about MIS 23 as a separate interglacial within the 41-kyr cycle, we consider it to represent a "warm" sub-stage of glacial MIS 22 belonging to the 100-kyr cycle; see below); and c) they are logically connected to joint and the unidirectional influences of two orbital parameters (precession and obliquity) which are expected to cause warming.

It should also be noted that there are cases of opposite insolation-related climatic impacts from precession and obliquity, during maximum values in eccentricity. For example, sub-stages MIS 5.4 and 7.4 are the most significant cooling periods within interglacial MIS 5 and 7 (Figure 2). These sub-stages correspond to maxima in eccentricity, but are determined by changes in precession and obliquity which lead to cooling. These changes are also reflected on the insolation curve as two distinct minima (see Figure 1).

So, returning to an OrCD maxima indicated by stars in Figure 2b, we can state that for the last 1 Ma the joint climatic impact of the variations in insolation, determined by the precession and obliquity, was not always reflected in paleoclimatic records adequately.

Overall, the 100-kyr eccentricity cycles of the OrCD correspond well to the analogous LR04 cycles for the last 1 Ma (see Figure 2). This circumstance allows us to stretch the OrCD further into the past, in order to compare it with the OI curve and to study the problem of the Middle Pleistocene Transition (MPT). The latter is related to the change in periodicity of global glacial cycles at about 1 Myr ago from the dominated 41-kyr period to the 100-kyr one [Clark et al., 2006; Ruddiman et al., 1986]. This phenomenon was accompanied by the increase of global ice volume, drop in temperature, and larger amplitude fluctuations in sea-level, global ice volume, and temperature over glacial-interglacial cycles. Despite the tremendous transformation of the natural environment related to the MPT, there is still no consensus explanation as to its cause(s).

Study of the Middle Pleistocene Transition

The MPT problem was recognized in the early 1980s [Pisias and Moore, 1981]. The essence of the problem is that change in dominant periodicity of global glacial cycles took place in the absence of substantial changes in orbital forcing. As it was pointed out, the adequate explanation of it is still absent. Besides it, the parameters used to define the MPT vary from study to study. Examples include: the timing of the beginning of the MPT; its duration; and changes in the characteristics which are used as indicators of global climate. For example, Clark et al. [2006, p. 3152] write: "... there is wide disagreement in defining when the MPT occurred, with descriptions ranging from an abrupt versus gradual transitions that began as early as 1500 ka and as late as 600 ka."

Clark et al. [2006] explained the above mentioned uncertainties as a consequence of using various deep-sea records. This is a reasonable suggestion, but even the use of the same record does not guarantee that the results will be the same. The use of the most widely cited OI record – the LR04 stack – may serve as a good example. Clark et al. [2006] used this archive to determine changes in global ice volume in the LR04 curve. They write: "... the MPT began 1250 ka with a gradual increase (decrease) in average ice volume (deep-water temperature), accompanied by an increase in the amplitude of the variability, with this transition reaching completion by 700 ka. Analysis of the frequency domain, however, suggests that while the emergence of the low-frequency signal also began 1250 ka, that signal disappeared for $\sim 100$-kyr before reemerging as a persistent signal since 900 ka." [Clark et al., 2006, p. 3180].

However, other researchers who used the LR04 record give different parameters of the MPT. For example, Lisiecki and Raymo [2007] determined the temporal interval for MPT as 1400–800 kyr ago, not 1250–700 kyr ago as suggested by Clark et al. [2006]. Bintanja and van de Wal [2008] give two intervals for the MPT, 1250–750 and 1400–800 kyr ago, centered around 1000 kyr ago. Both Lisiecki and Raymo [2007], and Bintanja and van de Wal [2008] consider the MPT as a transition from the 41-kyr climatic cycle to the 100-kyr one, and they do not talk about the disappearance of the 100-kyr cycle as discussed by Clark et al. [2006]. All these differences, especially in the determination of the MPT chronology, can be explained by dissimilarities in defining the meaning of MPT.

Within the framework of current discussion, let us compare the OrCD and LR04 data for the time interval of 1800–800 kyr ago. If we put aside the difference between the MIS 28–31 in LR04 record and OrCD (see below), it is possible to identify eccentricity 100-kyr cycles on the LR04 curve up to the maximum of MIS 37. For the earlier time frame, the OI record displays comparatively uniform and relatively even amplitude fluctuations with shorter periods, and they cannot be presented as a 100-kyr cycles. According to spectral analysis, the prevailing periodicity of these oscillations is about 41 kyr. Thus, it is reasonable to consider the maximum of MIS 37 as the moment of change in glacial rhythm from a 41-kyr period (determined by obliquity variations) to a 100-kyr one (controlled by changes in eccentricity). According to the LR04 chronology, it happens at 1.24 Myr ago. According to the OrCD data, it can be estimated at 1.239 Myr ago if we take into account the 5 kyr delay discussed above. Since this time, the 100-kyr glacial cycle has not been interrupted (Figure 2b).

This conclusion is drawn despite Clark et al.'s [2006] conclusion that the disappearance of 100-kyr cycle occurred at about 1 Myr. Our inference is drawn directly from analysis of the OI data within time interval of MIS 28–31. It is noteworthy that the phase of MIS 29 is counter to the OrCDS 24; and the same is true for MIS 30 and the related time span on the OrCD (Figure 2b). One should keep in mind that the OrCDS 24 Minimum should correspond to the deepest minimum in the LR04 record, according to the accepted mechanism of climatic impact from the orbital elements. However, the MIS 30 Minimum precedes the OrCDS 24 Minimum significantly, and this should not happen because the effect cannot precede the cause. It is possible to propose two explanations for this feature: 1) the OrCD do not reflect global climatic changes at that time; and 2) the LR04 record of global climate for this time span is distorted.

The first explanation does not look reasonable because the OrCD consistently follows global $\delta^{18}{\rm O}$ variations within the last 1.24 Myr. Also, for this time frame there are no controversial cases when OI extrema, which reflect maxima in cooling or warming trends, forestall corresponding OrCD extreme values. This leads us to favor the second explanation.

There are a number of reasons to suspect the LR04 record. First, the initial apportionment of MIS 27–29 is doubtful. The "cold" MIS 28 has almost the same maximal $\delta^{18}{\rm O}$ graphical value as the "warm" MIS 27. What is important in our case is that the existence of such $\delta^{18}{\rm O}$ maximum extends the temporal interval between MIS 27 and 29, and introduces a shift in the timing of MIS 30 toward the earlier period. It could be possible that for this time span there is a superfluous precession cycle in the LR04 record. Such inconsistencies in orbital tuning, with gap of two precession cycles or one obliquity cycle, have been detected previously [Bassinot et al., 1994; Shackleton et al., 1990]. According to Shackleton et al. [1990] and Bassinot et al. [1994], the inaccuracies distorted the correct age determination of the Matuyama-Brunhes boundary, especially in the SPECMAP curve chronology [Imbrie et al., 1984].

Second, Shackleton's [1995] OI curve was built on the basis of three cores, and the MIS 28 is represented by a minimum in $\delta^{18}{\rm O}$ value located between the maxima of MIS 29 and 27. A similar situation exists for the $\delta^{18}{\rm O}$ record of core MD 97–2140 from the Pacific Ocean [de Garidel-Thoron et al., 2005], and for the ODP 1143 core from the South China Sea [Ao et al., 2011]. In these cases, the maxima of MIS 27 and 29 contain an additional smaller $\delta^{18}{\rm O}$ maximum of definitely subordinate importance. This does not, however, allow us to infer this relatively small $\delta^{18}{\rm O}$ maximum as a separate unit of the 23-kyr precession cycle in the paleoclimatic record. Thus, we may assume that the clear maximum of MIS 28 as reflecting the corresponding precession extreme in the LR04 curve is the false one.

If this is the case, in order to fit the OrCD and LR04 records it is necessary to shift the timing of the MIS 29 Maximum to the time of MIS 28 Maximum. The MIS 30 Minimum will move concurrently, and appear in a position corresponding to the minimum of OrCDS 24. The latter will precede the minimum of MIS 30 as it should be. Also, the counter-phase position of the MIS 29 Maximum in relation to the OrCDS 24 Minimum will disappear. As a result, two 100-kyr eccentricity cycles can be clearly distinguished on the LR04 curve for the time interval of 1120–900 kyr ago, and they correspond to OrCDS 23–24 and 25–26. Further into the past, an additional full eccentricity cycle can be established in the LR04 record corresponding to OrCDS 27–28.

Therefore, after comparison of the OrCD and LR04 records it became evident that the 100-kyr cycle has not been interrupted for the past 1.24 Myr.


Before debating the results obtained, it is necessary to highlight again that the OrCD curve used for comparison with the OI data is compiled on the basis of the most general and simple notions of the mechanism of climatic impact of the insolation variations related to the changes in Earth's orbital elements. This circumstance provides clarity in the interpretation of the OrCD as it pertains to global climatic fluctuations. At the same time, these assumptions place distinct limitations on the application of OrCD for paleoclimatic interpretations [Bol'shakov, 2003b, 2010].

The Response of Climatic System to Orbital Forcing

One of the first conclusions reached when comparing the OrCD and LR04 records is that there are cases of inadequate climatic response to orbital forcing during the last 1 Myr. They are related to an absence in particular times of significant warming peaks which should correspond to OrCD maxima, reflecting the simultaneous contribution of variations in both obliquity and precession. These OrCD maxima occur during recessions or even minima in eccentricity when global ice volume increases.

Our analysis, including cases of logical reflection of orbital variations in paleoclimatic records, allows us to draw the conclusion that the causes of absence of climatic response to insolation are related to the nature of the climatic system, and its conditions at particular times. More detailed study of this phenomenon is needed to understand how the climate machine works, and how it transforms insolation changes into global climate fluctuations. At this point we do not know of similar cases documenting the presence or absence of correspondence between concrete orbital forcing and paleoclimatic records during the last 1 Myr, except for partial temporal agreement between paleoclimatic data and variations of discreet insolation in the Late Pleistocene.

A practical aspect of the above mentioned discrepancies between orbital forcing and paleoclimatic records is that we can suggest that uncertainties in use of the orbital tuning method for chronological purposes can occur. These ambiguities are related to inadequate reflection of variations of the orbital elements; in this case, obliquity and precession. This may lead to errors in the determination of the number of cycles for separate orbital elements, with either larger or smaller quantity. The case of the underestimated age for the Matuyama-Brunhes boundary is a good illustration. The possible existence of an additional precession cycle during MIS 27–29 in the LR04 record can be considered another example (see above).

Determination of Empirical Characteristics and Main Index of the MPT

As noted above, differences in definition of the main characteristics of the MPT, especially in its timing, explain the alternate conclusions drawn for the MPT in different studies. From the view of orbital theory, the main index of the MPT is the change of glaciation periodicity from an obliquity-driven 41-kyr period to an eccentricity-driven 100-kyr period. A common element of all MPT studies is that it was accompanied by an increase in global ice volume and a consequent drop in global temperature. Also, most scholars consider the decrease in temperature as the leading cause for the increase in global ice volume and the onset of the MPT [e.g. Berger et al., 1999; Bol'shakov, 2003b; Clark et al., 2006; Imbrie et al., 1993]. However, de Garidel-Thoron et al. [2005] and Elderfield et al. [2012] came to different conclusions. They studied changes in sea surface temperature (SST) of the eastern equatorial Pacific and deep-water temperature of the Pacific Ocean near New Zealand using the Ca/Mg ratios measured in planktonic and benthic foraminifers, respectfully. It was shown that minimal temperatures corresponding to glaciations did not change during the interval of 1500–600 kyr ago, although the MPT took place at this time. More than that, from data presented by Elderfield et al. [2012, p. 705, Figure 1] one can draw the conclusion that during the MPT the average deep-water temperature in a given region of the Pacific Ocean increased, as well as global ice volume. These findings conflict with the view that decreasing temperature initiated the MPT. The increase in both deep-water temperature and global ice volume seems to be contradictory also.

We assume that global ice volume and temperature (including deep-water one) are interdependent because they are part of a unified climate system. This is why the increase of global ice volume should lead to a decrease in global (and deep-water) temperature. The latter value, which corresponds to glacial epochs, cannot decrease below a minimal number (freezing point), so a drop in average temperature implies a decrease in interglacial temperatures within the glacial-interglacial cycle. But, according to Elderfield et al. [2012], a drop in interglacial temperatures did not occur. This is especially surprising as within each 100-kyr cycle changes of ice volume are observed to follow temperature fluctuations: a rise of temperature is always accompanied by a decrease in ice volume and vice versa.

General state of controversy in terms of ocean temperature changes at the time of MPT should be also noted. According to de Garidel-Thoron et al. [2005], the mean SST in the western equatorial Pacific did not change for the last 1.75 Myr. At the same time, the SST in the eastern equatorial Pacific [see Lawrence et al., 2006; Liu and Herbert, 2004] show a decrease within the interval 1400–150 kyr ago. The deep-sea temperature in the Pacific Ocean, according to Siddall et al. [2010], declined from 3.5 Myr to about 900 kyr ago and afterwards fluctuated near a constant mean value, without increase as it takes place according to Elderfield et al. [2012]. In the North Atlantic, according to Lawrence et al. [2010], a decline of both deep-sea and surface water temperatures occurred over the last 3 Myr. Based on these findings, we conclude that at present there is no clear indication that temperature change is a reliable index for the MPT.

Next we consider changes in global ice volume whose increase during the MPT is unambiguous, except in terms of timing. Lisiecki and Raymo [2005], and Clark et al. [2006] assumed that the change in ice volume took place gradually within a time span of about 1400–700 kyr ago. Elderfield et al. [2012], on the other hand, concluded that a sharp increase in ice volume has occurred around 900 ka. In both cases, the following changes in global ice volume remain not uniform: at some points after 700 kyr ago increases of ice volume were more significant (MIS 16 and 12), and at others were less remarkable (MIS 14 and 8). In this light, the date of 700 ka cannot be accepted as a distinct value typical of the end of the MPT process in terms of global ice volume. A similar conclusion can be drawn about changes in the amplitude of $\delta^{18}{\rm O}$ fluctuations during this period.

Elderfield et al. [2012] recognized a pivotal time of ice volume change that occurred about 900 ka, corresponding to MIS 22. But this point was determined on the basis of a single core ODP 1123, and a sole archive is not sufficient to document global ice volume changes. The globally stacked LR04 record is a better reflection of worldwide changes in ice volume and temperature. As a result of this discussion, one can conclude that the timing of the MPT determined by using the global ice volume and amplitude of $\delta^{18}{\rm O}$ variations as indicators are not without problems, and this is true for both the beginning and end of the MPT.

So, the most robust index of the MPT is the change in duration of the glacial-interglacial cycle. Based on our results, we interpret the MPT as a simultaneous change in the periodicity of glacial cycles that took place around 1.24 Myr ago. We cannot, however, fully exclude alternative view on the timing of the MPT, which is the following. It began about 1.24 Myr ago through two glacial-interglacial cycles that lasted about 100-kyr each, and at around 1070–970 kyr ago a perturbation in the cyclicity occurred. The MPT periodicity was restored about 950 kyr ago. In order to either prove or reject this possible scenario, it is necessary to get independent and precise age estimates for MIS 30, 29, and 27.

Differences in Records of the 100-kyr and 41-kyr Cycles and the MIS Numbering

Next we consider the difference between cycles in eccentricity and obliquity, as reflected in the OI curves. The main dissimilarity is that they belong to distinct scales of climatic fluctuations of the glacial-interglacial cycles. Both periods are easily distinguished with spectral-temporal analyses within the last 1 Myr. But the establishment of the 41-kyr cycle is more difficult because of the presence in the paleoclimate record of all three orbital inputs. That explains, for example, why the LR04 curve for the last 1 Myr is so complicated, with peaks with variable amplitude. The 100-kyr cycles can nevertheless be easily detected even with the naked eye in LR04 record, and these are closely tied to eccentricity extrema in the OrCD (Figure 2). On the contrary, during the time span of 1.8–1.25 Myr ago the 41-kyr component is much more homogenous, and usually consists of one maximum and one minimum. The amplitude of 41-kyr cycles is essentially smaller than the 100-kyr ones (see Figure 2b).

These differences between the 100-kyr and 41-kyr cycles make it impossible to distinguish MIS 23 as an interglacial stage of the 41-kyr cycle in the LR04 curve. This is already clear from a visual analysis of the LR04 record (Figure 1 and Figure 2). The MIS 23 maximum is similar to the MIS 18.3 one, which reflects interstadial warming within a glacial stage MIS 18 controlled by the 100-kyr cycle. It is well-represented on the OrCD as OrCDS 22. Thus, the MIS numbering scheme should be similar to the OrCD scheme up to MIS 37, and Number 29 (instead for No. 37) should be used to mark the beginning of the 100-kyr cycle for Pleistocene climatic events. This conclusion is in good accord with the inference that the Pleistocene change in periodicity from 41-kyr to 100-kyr occurred around 1.24 Myr ago. Although Clark et al. [2006] also dated this transition to about 1.25 Myr ago, they did not contest the incorrect numbering of MIS events around 1250–900 kyr ago. It is not necessary that this should lead to immediate re-numbering of the OI stages for the time span earlier than 900 kyr. Instead, this should serve as a reminder that current numbering of MIS 23–37 does not reflect adequately the process of 100-kyr eccentricity oscillations at corresponding time.

On the Mechanism of the MPT

According to our conclusions (see above), the MPT occurred at 1.239 Myr ago by change in the periodicity of glaciations. Thus, it is necessary to develop a cause-and-effect mechanism to explain the transition from 41-kyr cycles to 100-kyr ones. An attempt was made by Berger et al. [1999] who proposed a model of global ice changes for the last 3 Myr. They considered two external factors that influenced climate: 1) a linear decrease of carbon dioxide (CO$_2$) gas within this time frame; and 2) variations in the 65° N July insolation at the top of the atmosphere.

According to Berger et al. [1999], the decrease in CO$_2$ concentration causes both a gradual reduction of global temperature and an extension of ice volume. As the latter increases, the climatic system reacts in different ways to orbital insolation impact. Up to around 1 Myr ago when the ice volume was not large the Milankovitch mechanism would be the controlling factor, and glaciation may be expected to occur when summer (July sensu [Berger et al., 1999]) insolation reaches its minimum. Berger et al. [1999] assumed that this mechanism allows for the manifestation of 41-kyr obliquity periodic variations.

Another mechanism of insolation impact on climate existed after 1 Myr ago. Global ice volume increased at that time, in association with a reduction in CO$_2$ concentration and a drop in temperature. Berger et al. [1999, p. 8] write: "In this glacial world, the occurrence of interglacials needs very large insolation in summer at high latitudes, i.e. the conjunction of high eccentricity, high obliquity and the Northern Hemisphere summer at perihelion. ... As eccentricity reaches a maximum roughly every 100 ka, the interglacial occur every 100 ka. This is particularly clear over the recent 0.8 Ma: in our simulations, the complete or near-complete melting of the Northern Hemisphere ice sheets occur only when the 65° N July insolation exceeds $\sim 460$ Wm$^{-2}$, which are also the times of high eccentricity."

Besides the absence of justification for a linear decrease in CO$_2$ concentration, we can see three significant defects in the establishment of the MPT mechanism proposed by Berger et al. [1999]. The first defect is that they explain the emergence of interglacials and not glaciations. This distorts the essence of the Cenozoic era (and the Pleistocene epoch in particular) where climate change is personified by global cooling and a gradual increase in ice volume. At around 1 Myr ago and during the MPT, the cooling continued; however, the world of these times cannot be called "glacial". We consider the basic state of the Earth to be an interglacial one. This is why it is necessary to explain the appearance of just glaciations, as was done by founders of orbital theory at different stages of its development, including M. Milankovitch.

In our opinion, the fact that the basic state of the Earth is the interglacial one during the Pleistocene is reflected in the non-symmetric saw-tooth shape of the OI curves for the last 1 Myr. According to it, transition to the glaciations took place more slowly than return to interglacial conditions. Such dynamics are logically related to the fact that the exit from the more sustainable situation to a less stable one requires more energy and time for transition than the return to a relatively steady climate.

A second defect in Berger et al.'s [1999] model is the following. During the times of eccentricity's maxima, the values of monthly insolation are both maximal and minimal; the latter occurred around 11 kyr after the former (Figure 1). This results from a decisive contribution to monthly insolation of precession component which is modulated by the eccentricity. In these cases, it is necessary to explain why an increase in insolation brings about shrinking of the ice sheets but a subsequent insolation decrease does not cause a return to glacial conditions. Berger et al. [1999] do not consider this possibility.

The third defect in Berger et al.'s [1999, see Figure 6] model is that it uses the incoming signal of July insolation for 65° N, and this is not a complete insolation forcing but only a part of it. If the incoming signal does not fully represent the true signal, it is unreasonable to expect correct results from a model built on it. We note that the same authors (see above) emphasized the need to consider a complete picture of insolation to model paleoclimate.

Clark et al. [2006] explained the MPT phenomenon within the framework of a "regolith hypothesis" introduced by Clark and Pollard [1998]. The essence of this idea it that periodic glacial advances in the Early Pleistocene eroded the loose sedimentary rocks (calls "regolith") beneath the glaciers. Because of this, during subsequent glaciation at the beginning of the MPT glaciers moved over the hard bedrock. Periodic advances and retreats of continental glaciers changed the resilient pressure on crystalline bedrock, and according to Clark and Pollard [1998] this favored the increase in amplitude and periodicity of glacial fluctuations. The MPT process is described by Clark and Pollard [1998, p. 7] as this: "... at roughly 100 kyr intervals the ice sheet must become very large and thick with a deep bedrock depression, and thus susceptible to complete and rapid deglaciation in the next warm orbital interval. A deforming sediment layer prevents these conditions before $\sim 1$ Ma and allows them after it has been eroded...."

Several key points in the "regolith hypothesis" remain unclear. First, why does the transition just to 100-kyr eccentricity cycles take place? If Clark and Pollard [1998] assume that thawing happened "... in the next warm-orbit interval", subsequent glacial cycles should be around 21 kyr long. Our conclusion is based on the fact that changes in summer insolation [see Clark and Pollard, 1998, p. 2, Figure 1a] are determined mainly by precession. This is why warm orbital intervals that follow the cold ones should occur around 10–11 kyr later. Second, the conclusion about the exposure of crystalline bedrock during glacial movements after the Early Pleistocene is not fully justified because at least in Central/Eastern Europe, Middle Pleistocene moraine (glacial) deposits are divided by soft interglacial sediments of various origins (eolian, alluvial, lacustrine, etc.).

It is also worth mentioning two mechanisms of MPT proposed by de Garidel-Thoron et al. [2005] and Elderfield et al. [2012]. The former authors concluded that the MPT was caused by an increase in the equatorial temperature gradient in the Pacific Ocean. They assumed that this increase favors the strengthening of meridional transport of moisture to the high latitudes with consequent increases of global ice volume. However, this mechanism does not explain why the glacial cycles became just 100-kyr long. More than that, the mechanism of ice volume increase as proposed by de Garidel-Thoron et al. [2005] is doubtful. The supply of moisture would favor not only increases in volume, but also shrinking of glaciers, because vapor condensation in high latitudes will release latent heat. This conclusion follows from the fact that the Earth's solar climate should have a much higher meridional temperature gradient than it really exists today. Smoothing of the meridional gradient is caused by efficient meridional circulation of the atmosphere and hydrosphere.

Elderfield et al. [2012] assumed that the MPT was initiated by sharp increase in ice volume in the Antarctic at 900 kyr ago. However, this again does not explain the transition just to a 100-kyr glacial cyclicity. Elderfield et al. [2012, p. 707] postulate: "This period is associated with an anomalously low Southern Hemisphere summer insolation across the minor melting MIS 23." But in accord with the above mentioned critical remarks, several points need to be raised: 1) there are no grounds to connect this increase in ice volume, which is well-represented in the benthic $\delta^{18}{\rm O}$ records of ODP 1123 core at 42° S, with changes in discreet summer insolation at 75° S; 2) "the minor melting MIS 23" is not extraordinary event, because MIS 23 is not an interglacial but an ordinary interstadial warming well-represented on the OrCD; and 3) it is well-known that 100-kyr glacial cycles are determined first of all by changes in ice volume of the Northern Hemisphere, whereas fluctuations of ice volume in the Southern Hemisphere were of essentially smaller scale. This is why one should look for the cause of glacial rhythmic change among environmental conditions in the Northern Hemisphere or in the whole Globe.

Parametric Resonance Mechanism of the MPT

In the context of current discussion, the MPT can be related to the mechanism of so-called "parametric resonance" [Bol'shakov, 2001, 2003a, 2003b]. The potential importance of resonance in climatic systems has been suggested by several authors before [e.g. Sergin and Sergin, 1969; Hagelberg et al., 1991]. For our "parametric resonance" mechanism, changes of external conditions cause alterations of a system's parameters and, consequently, their resonance frequencies (periods). External conditions produced directional cooling during the time interval of 2–1 Myr ago; a trend which has begun much earlier, in the Eocene epoch. This cooling should have led to an increase in size of high latitude glaciers. Before about 1.2 Myr ago, however, the ice volume was relatively small, and the Earth's temperature drop was not enough to provide for ice sheet growth on a scale comparable with the Middle Pleistocene. In this case, changes in high latitude ice sheet volume were controlled by variations in $\varepsilon$ values with a relatively short period, in accordance with the empirically determined 41-kyr cyclicity. We can assume that the same duration was typical for periods of internal (resonance) ice oscillations at that time. The impact of eccentricity-driven insolation variations with a longer period, and of global temperature changes that result from them, was too weak to cause the expansion of the ice caps. However, with continued global cooling the volume of ice gradually increased.

Probably, at around 1.2 Myr ago the Earth's temperature and size of high latitude glaciers reached a critical value in relation to eccentricity-driven insolation changes. At this point, the eccentricity determined drop in temperature was large enough to overcome thawing of ice sheets, which were moving from the higher latitudes to the lower ones. Along with the increase of the mass and volume of glaciers, positive feedbacks intensified as well in response to changes in albedo and greenhouse gas concentration in the atmosphere; these factors favor the growth of glaciers [Bol'shakov, 2003a, 2003b]. The time constant (or resonance period) for the growth of glaciation also became longer. Hence, the combined effect of these three factors established a new rhythm for glacial cycles for the last 1 Myr.

Thus, the evolution and dynamics of global glaciations in the last 2 Myr was ruled mainly by joint action of variations in eccentricity and obliquity, and strengthened by positive feedbacks on the background of decrease in global temperature. We understand the hypothetical nature in asserting parametric resonance as a driver of climate. However, there are a number of arguments in favor of it. First of all, it is the similarity of the OrCD variations, which are determined mainly by changes in eccentricity, with 100-kyr cycles of the LR04 curve (Figure 2). The observation that OrCD extrema from the 100-kyr cycle correspond well to those of the LR04 curve (with some delay) is fully consistent with the idea that eccentricity variations directly affected climate over the last 1.24 Myr. Second, the resonance mechanism assumes selectivity. Resonance amplification is possible only for signals which coincide with characteristic periods of a system's internal fluctuations. In the Pleistocene, resonance amplification with 100-kyr periodicity is observed, but not the 400-kyr one. Third, resonance can explain how prolonged increases in ice sheet volume will bring not only expansion of their size but also make the characteristic period of glacial oscillations longer. This is why under the further growth of glaciation one should expect climatic impact of long eccentricity cycle of 400-kyr which was not displayed for the last 2 Myr [see Imbrie et al., 1993; Lisiecki, 2010].

The possibility of this scenario can be confirmed by published data [Heckel, 1986; Veevers and Powell, 1987]. According to them, the 400-kyr cyclicity was revealed in the maximal phase of the Permian-Carboniferous glaciation, and it was reflected in the sea-level changes. It is well-known that the Paleozoic glaciation of Gondwana was much bigger than the Pleistocene one, and glaciers have reached latitude of 30° S [Veevers and Powell, 1987]. Hence, the above mentioned data can confirm the proposed mechanism of links between insolation changes (determined by variations in eccentricity and obliquity), fluctuations of ice sheets, and rhythms of their expansion and shrinking.

One can compare this mechanism with changes in global ice volume and temperature recorded in the LR04 dataset for 1.3–1.1 Myr ago, which do not have evidence of sharp fluctuations in ice volume or temperature corresponding to the glacial phases of MIS 40, 38, and 36. What occurred at that time was the gradual increase of ice volume and a temperature decrease, as this continued until around 700 kyr ago. Therefore, the threshold value of ice volume which allowed for 100-kyr resonance fluctuations was reached about 1.2  Myr ago. One more conclusion drawn from the spectral analysis of the LR04 record [Lisiecri and Raymo, 2005] is that a continuous increase of ice volume did not prevent the response of the climatic system to insolation impact of obliquity variations. More than that, it favored climatic responses to shorter and less powerful (judging from the view of global full annual changes in insolation) precession insolation forcing.


This study of evolution and variability of Pleistocene glacial cycles is based on a correlation between the OrCD and benthic $\delta^{18}{\rm O}$ stack LR04 because we consider the latter as the most reliable record of global ice volume and temperature changes during the Quaternary. The development of the LR04 dataset is an important and necessary step forward in the investigation of global climatic change in the Pliocene–Pleistocene. Nevertheless, it is possible that there are uncertainties in the LR04 record, and it can be refined as suggested by its authors [see Lisiecki and Raymo, 2005, p. 6]. Below are some conclusions of our research.

  1. The OrCD is a better tool than discrete insolation curves for correlation with empirical data which reflect global climatic changes in the last 1.25 Myr. Thus, the OrCD is a reasonable alternative to continuous employment of the discrete insolation changes for paleoclimatic interpretation of empirical data during this time frame. Comparison of the OrCD and LR04 records suggest some possible errors in the LR04 chronology, particularly for MIS 27–29, where there is evidence of a superfluous precession cycle.
  2. Change in periodicity of dominant glacial cycles (i.e. the MPT phenomenon) took place at around 1.239 Myr ago, coinciding with the extreme value of MIS 37. After that, eccentricity cycles were not interrupted. Therefore, the numbering of OI stages beginning with the MIS 23 and ending with the MIS 37 is incorrect, because they reflect the 41-kyr cyclicity instead of the 100-kyr one which really existed at that time.
  3. Along with cases of logical reflection of the orbital forcing in paleoclimatic records (for example, for MIS 2–5, 7.4, 13–15, 18, and 23), the cases of inadequate reproduction of orbital forcing are revealed. The latter are generally characterized by joint one-directional impacts of both obliquity and precession variations, which one would expect to cause interstadial warmings during the glaciations. Evidently, the reason for these discrepancies can be found in the properties of the global climatic system itself at certain times.
  4. The most genuine and adequate indicator of the MPT is the change in periodicity of dominant glacial oscillations. Other characteristics, such as sharp increase of ice volume and oscillation amplitudes, are not as definite as the change in periodicity, and they took place both before and after changes in cycle duration in a chaotic manner. This suggests that these characteristics should not be incorporated into the definition of the MPT. They should be considered separately, possibly as factors which could lead to the emergence of the MPT phenomenon or accompany it.
  5. We believe that the cause of MPT, as well as of creation the ice sheets first in Antarctic and afterwards in the Northern Hemisphere, is a decrease of the Earth's temperature which began in the Eocene. The ultimate reason for this progressive cooling is not yet exactly known. We assume that the most possible mechanism of the MPT is parametric resonance as proposed earlier [see Bol'shakov, 2001, 2003a, 2003b].

In the process of this study, we came across a number of issues which need to be resolved for successful investigation the problem of the evolution and dynamics of the Pleistocene glacial cycles. One such question is the absence of precise and independent dating methods of sediments throughout the entire Pleistocene, i.e. the last 2.6 Ma. This would greatly help also to analyze the age discrepancies between the OrCD and LR04 records.

It is obvious that adequacy of knowledge about events in geological past depends on the reliability of empirical information. The same data are the basis for understanding nature's mechanisms of climate change. In this light, we should focus our efforts on resolving contradictory empirical evidence in the deep-sea core record at the time of MPT, especially with regard to temperature changes. Analysis of reliable information derived from long terrestrial archives would be particularly valuable. This choice is a logical consequence of the fact that glacial sheets originated on the continents, and the conditions for their development were determined largely by landmasses' climate. The most perspective in this respect would be the synthesis of empirical information for long lacustrine records.


This research was supported by the Russian Foundation for Basic Research (RFFI), grant 11-05-00147; Siberian Branch of the Russian Academy of Sciences (Program IV.31.2), and the Russian Academy of Sciences (Program 4.6, Project "Problems of Desertification of Central Asia"). We are grateful to Dr G. S. Burr (NSF-Arizona AMS Laboratory, University of Arizona, Tucson, AZ, USA) for help in translating this paper.


Adhémar, J. A. (1842), Revolutions de la Mer: D luges P riodiques, Carilian-Goeury et V. Dalmont, Paris.

Ao, H., M. J. Dekkers, L. Qin, G. Xiao (2011), An updated astronomical timescale for the Plio-Pleistocene deposits from South China Sea and new insights into Asian monsoon evolution, Quaternary Science Reviews, 30, 1560–1575, doi:10.1016/j.quascirev.2011.04.009.

Bassinot, F. C., L. D. Labeyrie, E. Vincent, X. Quidelleure, N. J. Shackleton, Y. Lancelot (1994), The astronomical theory of climate and the age of the Brunhes-Matuyama magnetic reversal, Earth and Planetary Science Letters, 126, 91–108, doi:10.1016/0012-821X(94)90244-5.

Berger, W. H. (1999), The 100-kyr ice-age cycle: internal oscillation or inclinational forcing? International Journal of Earth Science, 88, 305–316, doi:10.1007/s005310050266.

Berger, A., X. Li, M. Loutre (1999), Modelling northern hemisphere ice volume over the last 3 Ma, Quaternary Science Reviews, 18, 1–11, doi:10.1016/S0277-3791(98)00033-X.

Berger, A. L., M. F. Loutre (1991), Insolation values for the climate of the last 10 million years, Quaternary Science Reviews, 10, 297–317, doi:10.1016/0277-3791(91)90033-Q.

Berger, A. L., M. F. Loutre, H. Gallee (1998), Sensitivity of the LLN climate model to the astronomical and CO$_2$ forcings over the last 200 ky, Climate Dynamics, 14, 615–629, doi:10.1007/s003820050245.

Bintanja, R., R. S. W. van de Wal (2008), North American ice-sheet dynamics and the onset of 100,000-years glacial cycles, Nature, 454, 869–872, doi:10.1038/nature07158.

Bol'shakov, V. A. (1999), On the depth of paleomagnetic record acquisition by deep-sea sediments with reference to the climatic-stratigraphic position of the Matuyama-Brunhes reversal, Izvestiya, Physics of the Solid Earth, 35, 522–526.

Bol'shakov, V. A. (2001), A new concept of the astronomical theory of paleoclimate: One step forwards, after two steps backwards, Izvestiya, Physics of the Solid Earth, 37, 906–916.

Bol'shakov, V. A. (2003a), Modern climatic data for the Pleistocene: Implications for a new concept of the orbital theory of paleoclimate, Russian Journal of Earth Sciences, 5, 2, 125–143, doi:10.2205/2003ES000116.

Bol'shakov, V. A. (2003b), The New Concept of the Orbital Theory of Paleoclimate, (in Russian with English Abstract), Moscow State University Press, Moscow.

Bol'shakov, V. A. (2008), How long will the "precession epoch" last in terms of Pleistocene glacial cycles? Russian Journal of Earth Sciences, 10, ES3004, doi:10.2205/2008ES000299.

Bol'shakov, V. A. (2010), The problem of the 11th Marine Isotope Stage from the viewpoint of the new conception of the orbital theory of the paleoclimate, Oceanology, 50, 215–225, doi:10.1134/S0001437010020074.

Bol'shakov, V. A. (2013), Study of parameters of the Middle Pleistocene transition by comparison of the oxygen-isotope record LR04 with the orbital-climatic diagram, Doklady Earth Sciences, 449, 358–361, doi:10.1134/S1028334X13030197.

Bol'shakov, V. A., E. V. Ivanova, A. G. Prudkovskii (2005), Applying a new method for timing paleoclimatic deep-sea sedimentary records, Oceanology, 45, 867–876.

Bol'shakov, V. A., A. P. Kapitsa, W. G. Rees (2012), James Croll: a scientist ahead of his time, Polar Record, 48, 201–205, doi:10.1017/S0032247411000301.

Broecker, W., J. van Donk (1970), Insolation changes, ice volumes, and the $\delta^{18}{\rm O}$ record in deep-sea cores, Reviews of Geophysics and Space Physics, 8, 169–198, doi:10.1029/RG008i001p00169.

Cande, S. C., D. V. Kent (1995), Revised calibration of the geomagnetic polarity timescale for the Late Cretaceous and Cenozoic, Journal of Geophysical Research, 100, B4, 6093–6095, doi:10.1029/94JB03098.

Clark, P., D. Pollard (1998), Origin of the middle Pleistocene transition by ice sheet erosion of regolith, Paleoceanography, 13, 1–9, doi:10.1029/97PA02660.

Clark, P. U., D. Archer, D. Pollard, J. D. Blum, J. A. Rial, V. Brovkin, A. C. Mix, N. G. Pisias, M. Roy (2006), The middle Pleistocene transition: characteristics, mechanisms, and implications for long-term changes in atmospheric pCO$_2$, Quaternary Science Reviews, 25, 3150–3184, doi:10.1016/j.quascirev.2006.07.008.

Cleaveland, L. C., T. D. Herbert (2007), Coherent obliquity band and heterogeneous precession band responses in Early Pleistocene tropical sea surface temperatures, Paleoceanography, 22, PA2216, doi:10.1029/2006PA001370.

Croll, J. (1864), On the physical cause of the change of climate during geological epochs, Philosophical Magazine, 28, 121–137.

Croll, J. (1867), On the change in the obliquity of the ecliptic, its influence on the climate of the polar regions and on the level of the sea, Philosophical Magazine, 33, 426–445.

Croll, J. (1875), Climate and Time in their Geological Relations: A Theory of Secular Changes of the Earth's Climate, Edward Stanford, London.

de Garidel-Thoron, T., Y. Rosenthal, F. Bassinot, L. Beaufort (2005), Stable sea surface temperatures in the western Pacific warm pool over the past 1.75 million years, Nature, 433, 295–298, doi:10.1038/nature03189.

de Menocal, P. B., W. F. Ruddiman, D. V. Kent (1990), Depth of post-depositional remanence acquisition in deep-sea sediments: a case study of the Brunhes-Matuyama reversal and oxygen isotopic Stage 19.1, Earth and Planetary Science Letters, 99, 1–13.

Elderfield, H., P. Ferretti, M. Greaves, S. Crowhurst, I. N. McCave, D. Hodell, A. M. Piotrowski (2012), Evolution of ocean temperature and ice volume through the mid-Pleistocene climate transition, Science, 337, 704–709, doi:10.1126/science.1221294.

Hagelberg, T., N. Pisias, S. Elgar (1991), Linear and nonlinear couplings between orbital forcing and the marine $\delta^{18}{\rm O}$ record during the Late Neogene, Paleoceanography, 6, 729–746, doi:10.1029/91PA02281.

Hays, J. D., J. Imbrie, N. Shackleton (1976), Variation in the Earth's orbit: Pacemaker of the ice ages, Science, 194, 1121–1132, doi:10.1126/science.194.4270.1121.

Heckel, P. H. (1986), Sea-level curve for Pensylvanian eustatic marine transgressive-regressive depositional cycles along midcontinent outcrop belt, North America, Geology, 14, 330–334, doi:10.1130/0091-7613(1986)14<330:SCFPEM>2.0.CO;2.

Imbrie, J. (1982), Astronomical theory of the Pleistocene Ice Ages: A brief historical review, Icarus, 50, 408–422, doi:10.1016/0019-1035(82)90132-4.

Imbrie, J., J. Z. Imbrie (1980), Modeling the climatic response to orbital variations, Science, 207, 943–953, doi:10.1126/science.207.4434.943.

Imbrie, J., K. P. Imbrie (1986), Ice Ages: Solving the Mystery, Harvard University Press, Cambridge, MA.

Imbrie, J., et al. (1993), On the structure and origin of major glaciation cycles. 2. The 100,000-year cycle, Paleoceanography, 8, 699–735, doi:10.1029/93PA02751.

Imbrie, J., J. Hays, D. Martinson, A. McIntyre, A. C. Mix, J. J. Morley, N. G. Pisias, W. L. Prell, N. J. Shackleton (1984), The orbital theory of Pleistocene climate: support from a revised chronology of the marine $\delta^{18}{\rm O}$ record, A. L. Berger et al. (Eds.), Milankovitch and Climate, Reidel, Dordrecht, 269–305.

Johnson, R. G. (1982), Brunhes-Matuyama magnetic reversal at 790,000 yr. B.P. by marine-astronomical correlations, Quaternary Research, 17, 135–147, doi:10.1016/0033-5894(82)90055-2.

Karner, D. B., J. Levine, B. P. Medeiros, R. A. Muller (2002), Constructing a stacked benthic $\delta^{18}{\rm O}$ record, Paleoceanography, 17, 1030, doi:10.1029/2001PA000667.

Kukla, G. J. (1968), Comment to: Pleistocene epoch and the evolution of man, Current Anthropology, 9, 37.

Lawrence, K. T., Z. Liu, T. D. Herbert (2006), Evolution of the eastern Tropical Pacific through Plio-Pleistocene glaciation, Science, 312, 79–83, doi:10.1126/science.1120395.

Lawrence, K. T., S. Sosdian, H. E. White, Y. Rosenthal (2010), North Atlantic climate evolution through the Plio-Pleistocene climate transitions, Earth and Planetary Science Letters, 300, 329–342, doi:10.1016/j.epsl.2010.10.013.

Lisiecki, L. E. (2010), Links between eccentricity forcing and the 100,000-year glacial cycle, Nature Geosciences, 3, 349–352, doi:10.1038/ngeo828.

Lisiecki, L. E., M. E. Raymo (2005), A Pliocene-Pleistocene stack of 57 globally distributed benthic $\delta^{18}{\rm O}$ records, Paleoceanography, 20, PA 1003, doi:10.1029/2004PA001071.

Lisiecki, L., M. Raymo (2007), Plio-Pleistocene evolution: trends and transitions in glacial cycle dynamics, Quaternary Science Reviews 26, 56–69, doi:10.1016/j.quascirev.2006.09.005.

Liu, Z., T. D. Herbert (2004), High-latitude influence on the eastern equatorial Pacific climate in the early Pleistocene epoch, Nature, 427, 720–723, doi:10.1038/nature02338.

Loutre, M. F., D. Paillard, F. Vimeux, E. Cortijo (2004), Does mean annual insolation have the potential to change the climate? Earth and Planetary Science Letters, 221, 1–14, doi:10.1016/S0012-821X(04)00108-6.

Milankovitch, M. (1930), Mathematische Klimalehre und Astronomische Theorie der Klimaschwankungen, Gebruder Borntraeger, Berlin.

Milankovitch, M. (1941), Kanon der Erdbestrahlung und seine Andwendung auf das Eiszeitenproblem, Royal Serbian Academy, Beograd.

Muller, R., G. MacDonald (1995), Glacial cycles and orbital inclination, Nature, 377, 107–108, doi:10.1038/377107b0.

Pisias, N. G., T. C. Moore (Jr.) (1981), The evolution of the Pleistocene climate: a time series approach, Earth and Planetary Science Letters, 52, 450–458, doi:10.1016/0012-821X(81)90197-7.

Raymo, M. E., K. Nisancioglu (2003), The 41-kyr world: Milankovitch's other unsolved mystery, Paleoceanography, 18, 1011, doi:10.1029/2002PA000791.

Ruddiman, W. F., M. Raymo, A. McIntyre (1986), Matuyama 41,000-year cycles: North Atlantic Ocean and Northern Hemisphere ice sheets, Earth and Planetary Science Letters, 80, 117–129, doi:10.1016/0012-821X(86)90024-5.

Schneider D. A., D. V. Kent, G. A. Mello (1992), A detailed chronology of the Australasian impact event, the Brunhes-Matuyama geomagnetic polarity reversal, and global climate change, Earth and Planetary Science Letters, 111, 395–405, doi:10.1016/0012-821X(92)90192-X.

Sergin, V., S. Sergin (1969), How the Earth glaciations originated, Priroda, 9, (in Russian), 10–17.

Shackleton, N. (1995), New data on the evolution of Pliocene climatic stability, E. S.Vrba, G. H. Denton, T. C. Partridge, L. H. Burckle (Eds.), Paleoclimate and Evolution, with Emphasis on Human Origins, Yale University Press, New Haven, CT, 242–248.

Shackleton, N. J., A. Berger, W. Peltier (1990), An alternative astronomical calibration of the lower Pleistocene time scale based on ODP Site 677, Transactions of the Royal Society of Edinburgh, 81, 251–261, doi:10.1017/S0263593300020782.

Shackleton, N., N. Opdyke (1973), Oxygen isotope and palaeomagnetic stratigraphy of equatorial Pacific core V28-238: Oxygen isotope temperatures and ice volumes on a 105 year and 106 year scale, Quaternary Research, 3, 39–55, doi:10.1016/0033-5894(73)90052-5.

Siddall, M., B. Hönisch, C. Waelbroeck, P. Huybers (2010), Changes in deep Pacific temperature during the mid-Pleistocene transition and Quaternary, Quaternary Science Reviews, 29, 170–181, doi:10.1016/j.quascirev.2009.05.011.

Spell, T., I. McDougall (1992), Revisions to the age of the Brunhes-Matuyama boundary and the Pleistocene geomagnetic polarity timescale, Geophysical Research Letters, 19, 1181–1184, doi:10.1029/92GL01125.

Tauxe, L., T. Herbert, N. Shackleton, Y. Kok (1996), Astronomical calibration of the Matuyama-Brunhes boundary: consequences for magnetic reminisce acquisition in marine carbonates and Asian loess sequences, Earth and Planetary Science Letters, 140, 133–146, doi:10.1016/0012-821X(96)00030-1.

Thompson, W. G., S. L. Goldstein (2006), A radiometric calibration of the SPECMAP timescale, Quaternary Science Reviews, 25, 3207–3215, doi:10.1016/j.quascirev.2006.02.007.

Veevers, J. J., C. M. A. Powell (1987), Late Paleozoic glacial episodes in Gondwanaland reflected in transgressive-regressive depositional sequence in Euramerica, Geological Society of America Bulletin, 94, 475–487, doi:10.1130/0016-7606(1987)98<475:LPGEIG>2.0.CO;2.

Received 28 March 2014; accepted 10 April 2014; published 14 April 2014.

      Powered by MathJax

Citation: Bol'shakov V. A., Ya. V. Kuzmin (2014), Evolution and variability of the pleistocene ice ages: A new view, Russ. J. Earth Sci., 14, ES1003, doi:10.2205/2014ES000537.

Generated from LaTeX source by SemTeXML, v.1.2) software package.