Northern Hemisphere Glaciation: Its Tectonic Origin in the Neogene Uplift

: The Earth has cooled since the early Pliocene, which was punctuated by accelerated cooling indicative of thresholds. I posit that the cooling was initiated when the Neogene uplift of the Tibetan highland caused it to ice over, augmenting the albedo. I formulate a minimal warm/cold/Arctic box model to test this hypothesis and prognose the Pliocene climate. In particular, based on model physics, I discern three thermal thresholds as Pliocene cools: (1) when the Arctic temperature falls below the marking temperature of the ice front, the East Greenland ice sheet would descend to the sea level and calve into the Nordic Seas; (2) when the Arctic temperature cools to the freezing point, the ice sheet would form and expand over circum-Arctic lowlands to cause a massive deposition of ice-rafted debris marking Northern Hemisphere glaciation (NHG); (3) when glacial state persists through low eccentricity, it would cause a transition from obliquity-to eccentricity-dominated glacial cycles. Aligning these thresholds with the observed ones around 3.5, 2.7, and 1 million years ago, the model produces a temporal evolution of the Pliocene temperature as well as its driving albedo change. Since the latter can be accommodated by the observed one, it supports the Neogene uplift as the tectonic origin of NHG.


Introduction
The latest cooling of the Earth on the tectonic timescale began in the early Pliocene about 4 million years ago (Ma) [1,2].It was punctuated by glaciation of West Greenland [3,4] before culminating in Northern Hemisphere glaciation (NHG) manifested in a massive increase of ice-rafted debris (IRD) [5].Possible origins of the Pliocene cooling have been widely explored and extensive critiques of some proposed mechanisms have been put forth.Rather than reiterating these criticisms, I shall raise below additional points that are often overlooked.
First, one should be cautious in assigning causative roles to internal properties of the climate system in probing its response.One prominent example is the partial pressure of the atmospheric carbon dioxide (pCO 2 ), which would equilibrate with temperature over millennia [6]; hence, it may not be independently prescribed.As a case in point, the two track each other closely during glacial cycles [7], but since the temperature signal is induced by the orbital forcing via energy balance-with or without pCO 2 -the latter signal is a response rather than a driver of the temperature change.Over the tectonic timescale, pCO 2 has decreased to the preindustrial level since 20 Ma [8] and stayed low through the Pliocene, so it cannot account for the observed cooling of about 10 • C [4] based on the calculated sensitivity [9].
One well-known tectonic event is the closure of the Panama Isthmus, which has warmed the subtropics and increased poleward moisture transport [10], with the latter argued to propel the ice growth despite the warmth [11].However, this causality is opposite to that found in observations and climate models [12,13] hence remains unsupported.Moreover, the seaway closure was largely complete around 3. 6 Ma, which preceded NHG by 1 Myr [14], and, given the instantaneous (decadal) ocean response, this timing mismatch is difficult to bridge.As such, the seaway closure may cause only a temporary disruption of the Pliocene cooling [15].
Another tectonic candidate is the widespread uplift of the northern continents during the Neogene [16], which may enhance weathering to draw down pCO 2 , thus cooling the planet [17], but the effect is likely minor [18] and, in any case, weathering constitutes negative feedback to further muddle causality.An additional effect of the uplift well demonstrated by climate models is the amplification of the planetary wave [19], which would cool the northern lands, but this is predominantly a winter phenomenon [20], while it is the summer temperature that regulates the ice extent [21,22].
Another readily discerned effect of the uplift is the increase in planetary albedo caused by ice-topped highlands.The author of [23] shows that the Tibetan highland was ice-covered during the last ice age, accompanied by about 10 • C cooling, underscoring the efficacy of its albedo effect because of its vast size.This effect, however, is largely overlooked in the literature, and not even mentioned in most reviews of the subject [24,25].It is to quantify this albedo effect due to uplift and assess its possible role in the Pliocene cooling that prompts the present inquiry.
For this purpose, I shall present a minimal warm/cold/Arctic box model, which is nonetheless thermodynamically closed to allow a prognosis of the Pliocene temperature and glacial regimes.Specifically, because of inherently turbulent planetary fluids, a climate state is a macroscopic manifestation of a nonequilibrium thermodynamic (NT) system, and hence subject to maximum entropy production (MEP), a generalized second law.Readers are referred to [26,27] for reviews of MEP, whose physical basis is continually strengthened by its linkage to the fluctuation theorem [28,29], as the latter has been tested in the laboratory [30].In addition, recent eddy-resolving and direct numerical simulations of horizonal convection [31,32] have reproduced its signature feature (a mid-latitude front [33]), which is absent from coarse-grained general circulation models [34], and its application has replicated broad climate features and resolved long-standing glacial puzzles [35,36].With these supports, I am justified to posit MEP as a working hypothesis in formulating our model closure (see further discussion in Section 6).
The organization of the paper is as follows: I first formulate the box model in Section 2 and then link temperature to the radiative forcing in Section 3. In Section 4, I discern the thermal thresholds of changing glacial regimes.Aligning these thresholds with observed ones, I calculate in Section 5 the temporal evolution of the Pliocene climate, which supports the uplift hypothesis.I provide additional discussions in Section 6 and conclude the paper in Section 7.

Box Model
I consider a box configuration (Figure 1; all symbols are listed in Appendix A), in which the coupled ocean/atmosphere are divided into warm/cold boxes separated by the subtropical front (30 • N) and an Arctic Ocean box appended at the Fram Strait (80 • N).The sea-surface temperature (SST; thick solid line) represents a first-order description of the observed one, which may also be justified as an asymptotic state of infinite eddy mixing when subject to differential solar forcing.Such mixing would homogenize the temperature gradient within each box, which would be zero (hence isothermal) in warm/Arctic boxes because of the vanishing gradient at the equator/pole of a sphere, but retain a finite value in the cold box due to its anchoring at a possibly different Arctic temperature, a profile that is nonetheless fully specified by the derived cold-box average (thick dashed line).Since there is no blockage of the Arctic air, the troposphere comprises only two isothermal boxes.
Being a minimalistic model, I retain only dominant energy balances: the system is heated differentially by the incoming shortwave (SW) flux (q * i ), which after atmospheric absorption and reflection by the planetary albedo is absorbed by the ocean (q i , referred as 'radiative forcing') to differentiate the SST.The latter, in turn, differentiates the surfaceair temperature (SAT) by the convective flux (q c ) neglecting by comparison latitudinal variation of the surface-and outgoing-longwave (LW) fluxes.Because of the thermal inertia of the ocean, I neglect its seasonality.Air, however, is strongly warmed in summer by atmospheric absorption, which remains bounded by the underlying SST in maintaining radiative-convective equilibrium of the troposphere.This 'convective bound' is a distinct feature in observation [37] (their Figure 10.7, lower panel), so I shall take SST as a proxy of the co-zonal summer SAT that controls the ablation, and hence glacial margin.Being a minimalistic model, I retain only dominant energy balances: the system is heated differentially by the incoming shortwave (SW) flux ( * ), which after atmospheric absorption and reflection by the planetary albedo is absorbed by the ocean ( , referred as 'radiative forcing') to differentiate the SST.The latter, in turn, differentiates the surface-air temperature (SAT) by the convective flux ( ) neglecting by comparison latitudinal variation of the surface-and outgoing-longwave (LW) fluxes.Because of the thermal inertia of the ocean, I neglect its seasonality.Air, however, is strongly warmed in summer by atmospheric absorption, which remains bounded by the underlying SST in maintaining radiative-convective equilibrium of the troposphere.This 'convective bound' is a distinct feature in observation [37] (their Figure 10.7, lower panel), so I shall take SST as a proxy of the co-zonal summer SAT that controls the ablation, and hence glacial margin.
In our convention, global means are overbarred, deviations (referred to as 'differential') are primed and subscripted 1/2/3 for warm/cold/Arctic boxes.The latitudinal distance is defined as  = sin (latitude), so equal intervals span equal surface area on a sphere, and the ocean boundaries at 30 and 80° N have  = 0.5 and 0.98, respectively.Given the smallness of the Arctic box (∆ = 0.02), I shall neglect its effect on the other box-averaged temperature, which thus is simply the two-box solution derived in [29].Readers are referred to [29] for mathematical details, but the abbreviated derivation given below should aid the self-containment of the paper.

Differential Temperature
In the following, I shall first link the differential SAT to SST via atmospheric MEP.With the dominant heat balance stated earlier, the irreversible entropy production of the atmosphere equals the entropy flux through its lower boundary, which is of the form ( −  )/ , where  is the air-sea exchange coefficient set to 10 Wm K [29].Given the small temperature deviation from its global mean (in Kelvins), the inverse temperature in the above expression can be linearized, so the differential entropy production (that is, removing the global mean) varies as (the symbol "~") which is expressed in the cold-box variables as the warm box has the same value.Maximizing this entropy production against  , yields The differential incoming SW flux (q * i ), after atmospheric absorption and reflection by planetary albedo, is absorbed by the ocean (q i ) to differentiate the SST, which in turn differentiates the SAT by the convective flux (q c ).The latitudinal coordinate is defined as y = sin(latitude), so ocean-box boundaries lie at y = 0.5 and 0.98.
In our convention, global means are overbarred, deviations (referred to as 'differential') are primed and subscripted 1/2/3 for warm/cold/Arctic boxes.The latitudinal distance is defined as y = sin (latitude), so equal intervals span equal surface area on a sphere, and the ocean boundaries at 30 and 80 • N have y = 0.5 and 0.98, respectively.Given the smallness of the Arctic box (∆y = 0.02), I shall neglect its effect on the other box-averaged temperature, which thus is simply the two-box solution derived in [29].Readers are referred to [29] for mathematical details, but the abbreviated derivation given below should aid the self-containment of the paper.

Differential Temperature
In the following, I shall first link the differential SAT to SST via atmospheric MEP.With the dominant heat balance stated earlier, the irreversible entropy production of the atmosphere equals the entropy flux through its lower boundary, which is of the form α(T − T a )/T a , where α is the air-sea exchange coefficient set to 10 Wm −2 K −1 [29].Given the small temperature deviation from its global mean (in Kelvins), the inverse temperature in the above expression can be linearized, so the differential entropy production (that is, removing the global mean) varies as (the symbol "~") which is expressed in the cold-box variables as the warm box has the same value.Maximizing this entropy production against T ′ a,2 yields a well-known result of MEP [26] because of mutual compensation of the thermodynamic flux (heat transport, the bracketed term) and thermodynamic force (differential temperature).

22
Having linked SAT to SST, I next apply the oceanic MEP to link SST to the radiative forcing q ′ i .Since the differential heating of the cold ocean box is the radiative forcing minus the convective flux, the oceanic counterpart to (1) is Applying (2) and maximizing (3) against T ′ 2 , we obtain the differential SST is thus linear in the differential forcing.As a cursory check, a forcing range of 200 Wm −2 [37] would produce an SST range of 20 • C, not unlike that observed [38].
In addition, the model-deduced meridional overturning circulation depends only on α to yield a transport of 11 Sv for the 6000 km wide North Atlantic [36], which is of the same order as the observed one [39].These observational comparisons constitute a veritable test of MEP as they contain no free parameters-unlike previous comparisons involving tunable parameters [40].
Given the narrowness of the Fram Strait, which is the sole gateway to the Arctic Ocean since the Miocene [41], I shall neglect the heat transport into the Arctic Ocean, so that Substituting from ( 2) and ( 4), we derive As seen in Figure 1, the Arctic Ocean could be 50% colder than the averaged subpolar temperature (thick dashed line) with respect to the global mean, not unlike the present case [37] (their Figure 8.2a), and since the Arctic temperature anchors the cold-box SST (thick solid line), the latter has a finite gradient in halting the glacial advance.
Having linked temperature to radiative forcing, I shall next examine how upliftinduced albedo change affects radiative forcing, and hence temperature.The radiative forcing of the cold box is where q * 2 is the incoming SW flux, b the (constant) atmospheric absorption, and a 2 the planetary albedo; only the last is altered by uplift.Let i be the area rendered ice-covered by uplift (in fraction of the cold-box area), then the cold-box forcing varies as where ∆r, set to 0.7, is the excess albedo of ice over land.Since the Neogene uplift is dominated by mid-latitude lands within the northern cold box, they are assumed to be free of low clouds concentrated in high latitudes [42] (high clouds have little albedo), and with the northern cold box comprising one quarter of the global surface, the global forcing thus varies as δq = δq 2 /4, ( and the differential forcing varies as That is, because of the concentration of uplift in the northern cold box, its effect on the differential forcing is three times that on the global forcing.The latter induces global cooling via δT = sδq, where the 'global sensitivity' s has a range of [0.5, 1.5] K Wm −2 −1 based on climate models [9].Substituting from ( 9) and (8), the global cooling is thus linked to the ice cover via δT = −(sq * 2 ∆r)δi/4, and from ( 4), ( 6), (10), and (11), the global cooling would cool the cold and Arctic boxes via It is seen that the cold box is cooled more strongly than the globe temperature due to contribution from the differential forcing (10), and the cooling is further amplified in the Arctic Ocean because of blockage of the incoming ocean heat.
With the variations derived above, one may calculate the temperature evolution given their initial values and observed uplift, a prognostic scheme that is, however, thwarted by the uncertain uplift history.On the other hand, there is a proxy time series of the subpolar SST, which exhibits abrupt cooling episodes indicative of thresholds, so if we can identify these thresholds from our model, we may then fix their timing to diagnose the uplift history.Taking this approach, I shall next consider the thermal thresholds from our model.

Thermal Thresholds
To illustrate the thermal thresholds as the Pliocene cools, I draw in Figure 2 three cold/Arctic-box SST profiles (same as the summer SAT) representing the initial condition of the early Pliocene (EP), the expansion of the Greenland ice sheet (GIS), and NHG.For a constant lapse rate set to γ = 6 • C, the temperature profiles also represent the snowline height, both are tick-marked on the right ordinate.The highland above the snowline would be ice-covered, and the net accumulation would propel the ice down the slope until the ice flux is depleted by ablation, an ice edge defining the 'glacial line' and its temperature T g referred to as the 'glacial marking temperature' (GMT).This temperature is drived in Appendix B, which turns out to be primarily an intrinsic property of the ice sheet, and its deduced value of about 5 • C is consistent with the present observation [43], which is indicated in a dashed line.It is recognized that the GMT derived in Appendix B over a slope would be altered for a sea-level ice sheet [36], but since we are not prognosing the ice extent on the sea level, it is drawn as level for illustration.
Dark-shaded columns represent Greenland and the Tibetan plateau, which were uplifted during the Neogene.Greenland is perhaps elevated in two steps [44], and being located in high latitudes with lower snowline, its eastern part was already ice-covered by the late Miocene [45], and hence had no impact on the subsequent albedo change.For the Tibetan highland, I subscribe to its uplift history synthesized by [46] based on vegetation data, which shows a notable lag of its northern from southern provinces (see also [47]), but despite the timing spread, they converge and reach beyond the snowline in the early Pliocene (EP) to initiate the cooling.It is noted that [48] proposed an earlier uplift (~8 Ma) based on assumed crustal/mantle thickness, which the author acknowledged as speculative, and it conflicts moreover with the early Pliocene uplift (~4.5 Ma) based on sedimentary evidence [49].Also of relevance is the fact that the Tibetan plateau was ice-covered during the last ice age when the snowline fell below the plateau height to cause about 10 • C cooling, underscoring the potency of its albedo effect because of its vast spatial coverage [23].Given the phased uplift of the Tibetan highland, one expects a gradual cooling on the tectonic timescale and, as the cooling proceeds, the lowering of the snowline would expose lesser highlands, such as the American West [24] and margins of the North Atlantic Ocean [50], to glaciation [51], so it is this combined uplift area that determines its albedo effect and whether it is sufficient to cause NHG.Dark-shaded columns represent Greenland and the Tibetan plateau, which were uplifted during the Neogene.Greenland is perhaps elevated in two steps [44], and being located in high latitudes with lower snowline, its eastern part was already ice-covered by the late Miocene [45], and hence had no impact on the subsequent albedo change.For the Tibetan highland, I subscribe to its uplift history synthesized by [46] based on vegetation data, which shows a notable lag of its northern from southern provinces (see also [47]), but despite the timing spread, they converge and reach beyond the snowline in the early Pliocene (EP) to initiate the cooling.It is noted that [48] proposed an earlier uplift (~8 Ma) based on assumed crustal/mantle thickness, which the author acknowledged as speculative, and it conflicts moreover with the early Pliocene uplift (~4.5 Ma) based on sedimentary evidence [49].Also of relevance is the fact that the Tibetan plateau was icecovered during the last ice age when the snowline fell below the plateau height to cause about 10 °C cooling, underscoring the potency of its albedo effect because of its vast spatial coverage [23].Given the phased uplift of the Tibetan highland, one expects a gradual cooling on the tectonic timescale and, as the cooling proceeds, the lowering of the snowline would expose lesser highlands, such as the American West [24] and margins of the North Atlantic Ocean [50], to glaciation [51], so it is this combined uplift area that determines its albedo effect and whether it is sufficient to cause NHG.
The first threshold (labeled GIS for Greenland ice sheet) is when the Arctic temperature cools to the GMT (dashed line), so it is defined by When this threshold is crossed, the ice sheet over East Greenland would descend to the sea level and advance southward toward GMT, which, however, is cut short by the shoreline to calve into the Nordic Seas, resulting in minor IRD events.It is important to note that this threshold leading to full glaciation of Greenland is predicated on the existing ice sheet atop East Greenland, which provides the seeding for its expansion to sea level.Without such seeding, there can be no sea-level ice sheet so long as the surface The first threshold (labeled GIS for Greenland ice sheet) is when the Arctic temperature cools to the GMT (dashed line), so it is defined by When this threshold is crossed, the ice sheet over East Greenland would descend to the sea level and advance southward toward GMT, which, however, is cut short by the shoreline to calve into the Nordic Seas, resulting in minor IRD events.It is important to note that this threshold leading to full glaciation of Greenland is predicated on the existing ice sheet atop East Greenland, which provides the seeding for its expansion to sea level.Without such seeding, there can be no sea-level ice sheet so long as the surface temperature remains above the freezing point, leading to the second threshold delineated next.
The second threshold (labeled NHG) is when the Arctic Ocean is cooled to the freezing point or approximately NHG : At this juncture, the convective flux has exhausted its capacity in countering the weakening radiative forcing by cooling the surface, leading to the formation of perennial sea ice, which would expand to the whole Arctic Ocean and expel the snowline to the surrounding lowlands.Subsequently, even a slight cooling perturbation would allow a foothold of positive ice mass balance, which would then expand the polar cap southward until its leading temperature reaches GMT (that is, the intersect of solid and dashed lines), as indicated by the light shade.The polar cap would calve into the surrounding ocean to produce massive deposition of IRD-the defining marker of NHG.
The generation of sea/land ice during NHG would augment the albedo to continue the Pliocene cooling, but also significantly, with the polar cap now encroaching on the contiguous northern lands, it would activate the ice-albedo feedback to implement the precession forcing, which otherwise would be nil on account of Kepler's second law [22]: a lower summer insolation, for example, would be countered by a longer summer to yield zero annual mean in the radiative forcing.As posited by this author [36], this is the reason that the 41 ky obliquity signal prevails prior to NHG, but the 21 ky precession signal would emerge during NHG to be followed by distinctive mid-Pleistocene transition (MPT) to 100 ky ice-age cycles paced by eccentricity [52].Ref. [36] noted that MPT remains unresolved for lack of evidence of putative pCO 2 decrease and regolith removal [52][53][54] and provided an alternative mechanism of a bistable ocean, to which readers are referred for details.A brief derivation is nonetheless provided below to further quantify MPT and address causality between NHG and glacial cycles [55].
[36] postulated that the augmented differential forcing beyond the convective bound (Section 2) would vault the cold box into the freezing point by MEP, rendering a bistable temperate/freezing-point cold box, which would translate to a bistable polarcap/Laurentide ice sheet (LIS).The bistable forcing interval is set by the global convective flux q c , which has a centerline lying about 1.65 q c below the global forcing.The global surface heat balance is of the form where q and q LW are the global radiative forcing and the surface LW flux, respectively.Defining a 'greenhouse' parameter and setting it to s g ≈ −2 [42] (their Figure 1 over temperature range [283 K, 293 K] in accordance with Figure 3), we derive The global convective flux thus decreases during the Pliocene cooling due to both the weaker radiative forcing and drier air.As shown in [36], the threshold of MPT is when the 'excess' convective flux defined below vanishes, so for the third threshold, MPT : q e ≡ 1.65 When this threshold is crossed, the glacial state would persist through low eccentricity, thus allowing ice growth to LIS that defines the ice age, and the glacial cycles would transition from obliquity-to eccentricity-dominated periodicity.Substituting from ( 20), ( 19), (10), and (11), the excess convective flux varies as follows: Therefore, the MPT threshold is specified by its initial value and global cooling.Having identified these thresholds, we shall next discern them from the observed temperature evolution.

Pliocene Climate
Ref. [4] has provided a high-resolution subpolar SST from 4 Ma to the present.The thick line shown in Figure 3 represents my interpretation of its broad features, which can be characterized as gradual cooling punctuated by three more abrupt episodes of cooling, around 3.5, 2.7, and 1 Ma.Although [56] has inferred a mid-Miocene warmth (~3 Ma) from their data spanning 3-1.7 Ma, it is better interpreted, based on the longer time series of [4], as slowing of the cooling when sandwiched between two rapid cooling episodes.Since the temperature is inferred from proxy data of the eastern North Atlantic with summer bias, it is expected to be quite warmer (about 5 • C) than the averaged cold-box temperature (T 2 ) [57], but the abrupt changes dovetail the modeled thresholds to allow fixation of their timing (dashed vertical lines).
Ma) from their data spanning 3-1.7 Ma, it is better interpreted, based on the longer time series of [4], as slowing of the cooling when sandwiched between two rapid cooling episodes.Since the temperature is inferred from proxy data of the eastern North Atlantic with summer bias, it is expected to be quite warmer (about 5 °C) than the averaged coldbox temperature ( ) [57], but the abrupt changes dovetail the modeled thresholds to allow fixation of their timing (dashed vertical lines).For the initial condition, it is noted that the early Pliocene (EP) is slightly warmer than the Holocene [58], so I set the following values based on the present observation [37] (their Figures 6.  Together with the variations derived in sections 3 and 4, I calculate the temperature and ice cover at the thresholds.Connecting them by smooth curves and incorporating more rapid cooling at the thresholds, I produce the modeled time series shown in Figure Figure 3.Time evolution of the modeled global (T), averaged cold-box (T 2 ), and Arctic (T 3 ) temperatures and uplift-induced ice cover (i) when thermal thresholds (dashed vertical lines) are aligned with abrupt cooling of the observed subpolar SST (thick solid line, adapted from [4]).The shaded portions are solution spans for global sensitivity ranging from 0.5 (thin solid line) to 1.5 (thin dashed line).The solid square is the estimated ice cover by uplift, which is seen sufficient to cause NHG.
For the initial condition, it is noted that the early Pliocene (EP) is slightly warmer than the Holocene [58], so I set the following values based on the present observation [37] (their Figures 6.3 Together with the variations derived in Sections 3 and 4, I calculate the temperature and ice cover at the thresholds.Connecting them by smooth curves and incorporating more rapid cooling at the thresholds, I produce the modeled time series shown in Figure 3 (thin lines).I have varied global sensitivity (s) over its outer range, between 0.5 (thin solid line) and 1.5 (thin dashed line) • C Wm −2 −1 [9], so solution spans (shaded portions) underscore their considerable uncertainty, but the general features remain robust.As these time series represent the prognosed cooling driven by the uplift-induced ice cover, I offer the following interpretation of the Pliocene climate.
Since the uplift of the Tibetan and American West highlands is largely confined to the subpolar zone (that is, north of 30 • N), the resulting albedo increase has a greater effect on differential than global forcing-by three-fold (10)-leading to stronger subpolar than global cooling (13).The cooling is further augmented for the Arctic Ocean due to blockage of the incoming ocean heat, resulting in Arctic cooling that could double the rate of global cooling (14).It is significant to note that this Arctic amplification is comparable to that caused by anthropogenic warming [60] even without ice-covered Arctic Ocean and the accompanying ice-albedo feedback.
When the Arctic temperature falls below the GMT of 5 • C around 3.5 Ma (the threshold labeled GIS), the pre-existing ice sheet over East Greenland would expand to West Greenland, descend the slope to the sea level, and advance toward the GMT.The glacial advance, however, is cut short by the shoreline to calve into the Nordic Seas.The rapid ice expansion on the millennial timescale may account for the observed abrupt cooling, which would be accompanied by minor IRD events, whose source has been traced to Greenland [3].
The slow uplift and attendant ice-albedo feedback continues to cool the climate, and when the Arctic Ocean is cooled to the freezing point (the threshold labeled NHG), perennial ice would grow to cover the Arctic Ocean with only a small open area maintained by nonzero ocean heat entering through the Fram Strait, as is the current case.The icecovered Arctic Ocean would expel the snowline to the surrounding lowlands to allow a foothold of positive ice mass balance when perturbed, which would then grow the polar cap toward the GMT.The difference from the GIS threshold is that without seeding by the mountain glaciers, the polar cap can be initiated only when the sea level has cooled to the freezing point (not just the GMT), which represents a more stringent requirement.This is the reason that there is no widespread glaciation over the circum-Arctic lands until NHG [61], which must thus be distinguished generically from the GIS threshold.Again, the rapidity of the ice expansion causes the more abrupt cooling around 2.7 Ma, and the polar cap expanding southward over the contiguous lands would calve into both the Arctic Ocean and Nordic Seas to dispense massive IRD, the defining marker of NHG [5].Contrary to earlier analysis [62], our model suggests that the Arctic Ocean should be largely free of perennial ice until NHG, although its observational evidence remains debated [63,64].
As expected, smaller global sensitivity (thin solid lines) would reduce the global and subpolar cooling relative to Arctic cooling, so a more expansive ice cover is needed for attaining NHG, but even the low end of this sensitivity requires only 17% of the cold-box area to be uplifted.If we assume that the Tibetan and American West highlands comprise 30% of the subpolar land [19], which occupies 60% of the cold box, the uplifted area would be 18% of the cold-box area (solid square), which hence is sufficient to cause NHG and even more so for greater global sensitivity.In addition, as we noted earlier, the cooling would lower the snowline to expose lesser highlands to glaciation, including the West Greenland ice sheet formed during the first threshold, further augmenting the uplift effect.With the above, I reckon that the Neogene uplift and its boost of the planetary albedo are sufficient to cool the Pliocene to NHG.
As discussed in Section 4, the generation of a polar cap over contiguous northern lands during NHG would activate the ice-albedo feedback to implement precession forcing and amplify obliquity forcing.Together with the augmented ice volume response, the glacial cycles would exhibit a stronger δ 18 O signal, as observed [65].Supported by the physics delineated above, I reckon, therefore, that the observed amplification of glacial cycles is the outcome rather than the driver of NHG.Although the reverse causality has been suggested [56], no physical mechanism has been advanced.Nor can it be attributed to changing orbital parameters, which simply do not vary on such timescales.
With continuing cooling by the uplift of circum-Arctic lands [66] and lowering of the global convective flux due to drier air, a threshold (labeled MPT) would be crossed when the glacial state would last through low eccentricity to enable LIS that defines the ice age.Again, ice expansion occurs over short mass balance timescale, so the threshold may account for the more abrupt cooling around 1 Ma, which represents MPT from obliquity-dominated glacial cycles to 100 ky ice-age cycles paced by eccentricity [52].Ref.
[36] provides a detailed discussion of MPT, and Figure 3 offers a tangible constraint on its timing, namely, when the global temperature cools to below about 10 • C.

Discussion
Although climate models are widely used to simulate past climate, their limitations should not be overlooked.In particular, the ocean is teeming with eddies, whose effect is incorporated through sub-grid parameterization finely tuned to yield the current climate [67], which necessarily downgrades their prognostic utility for a vastly different paleoclimate.The sensitivity can be lessened if climate models have resolved eddies, which is computationally challenging and remains to be attempted in paleoclimate studies.Our box model, while extremely crude, has nonetheless accounted for eddies by a closure based on probability laws of the nonequilibrium thermodynamics, namely the MEP-a generalized second law-and its deduced differential temperature and meridional overturning circulation, which involve no free parameters, are commensurate with those observed in support of the MEP.Without resolving eddies to entail MEP, climate models may not be used to arbitrate our deductions, which can only be tested against observation.Given the crudeness of our box model and the uncertainty in observations, this comparison is limited to broad qualitative features, as depicted in Figure 3.To the degree that the modelled time series resemble the observed ones, including three abrupt cooling episodes, it is suggested that the model has captured the governing physics of the Pliocene cooling.
A critical parameter of our model is the glacial marking temperature (GMT) when net accumulation above the snowline is depleted by ablation down the slope.I have provided a derivation of this temperature in Appendix B, which turns out to be primarily an intrinsic property of the ice sheet, and the deduced value of ~5 • C is consistent with the present observation of the Greenland ice sheet [43] (their Figure 8).A fixed GMT sets the vertical drop (about 830 m) of glacial line from snowline, so the two would move in tandem by a changing climate, which moreover is through temperature, not moisture transport.The short-circuiting of the latter is because only the one crossing the snowline of freezing point is of relevance, which is thus insulated from the SST.

Conclusions
Given the many questions raised about previous proposals, the tectonic origin of the Pliocene cooling remains unresolved.Taking note of the Neogene uplift of the Tibetan highland, I posit that its icing above the snowline, hence augmented albedo, may initiate the Pliocene cooling.To test this hypothesis and more generally to explain the time evolution of the Pliocene climate, I formulate a minimal warm/cold/Arctic box model that is thermodynamically closed to prognose the Pliocene SST (same as the summer SAT on account of the well-validated convective bound).After the initial trigger, the ice-albedo feedback would continue the cooling, and three thermal thresholds are discerned: (1) when the Arctic SST falls below the GMT, the ice sheet over the East Greenland would descend to the sea level, accompanied by episodic IRD events; (2) when the Arctic SST cools further to the freezing point, the ice cap would form and expand on circum-Arctic lowlands, as manifested in massive IRD deposition that defines NHG; (3) when global convective flux is sufficiently lowered by the global cooling, the glacial state may last through low eccentricity to cause transition from obliquity-to eccentricity-dominated glacial cycles.
Given the short timescale governing the ice mass balance, the above thermal thresholds should appear as more rapid cooling punctuating the gradual Pliocene cooling, as indeed seen in the observed subpolar SST.Fixing timing of the modelled thresholds to the observed ones around 3.5, 2.7, and 1 Ma, I produce a temporal evolution of the Pliocene climate as well as its driving albedo change.Since the latter can be accommodated by the observed uplift and growing ice cover, the model provides a causal explanation of the Pliocene cooling via the interplay of uplift-induced albedo, temperature, and glacial regimes.This supports our hypothesis that the Neogene uplift of the Tibetan highland is the tectonic origin of the Pliocene cooling and NHG.Absorbed SW flux of box i q Global-mean absorbed SW flux q c

Appendix
Convective flux q r Net radiative flux through tropopause q OLR Outgoing LW radiation q SW Incoming SW flux ∆r Excess reflectance of ice over land s Slope ( = 10 I provide here a derivation of the marking temperature of the glacial line, a key parameter of the model.As sketched in Figure A1 (see for example [43], their Figure 2a, July panel), I consider the north-south section of an idealized plateau when the snowline (thick solid line) is lower than the plateau height, so the moisture transport crossing the snowline would all manifest in the accumulation (no need to distinguish snowline and equilibrium line for our purpose).This moisture transport is seen later to be dominated by the summer season, so moisture and energy transports considered below pertain to averages of the summer-half of the year (referred to as 'summer' for short), which moreover are assumed to be in balance-neglecting seasonal storage of air and ice on account of their small thermal inertia.Assuming north-south symmetry of the ice sheet over the plateau, only half of the accumulation would descend the southern slope until it is depleted by ablation, which defines the glacial line (thick dashed line).The aim is to determine its temperature, referred to as the GMT.
averages of the summer-half of the year (referred to as 'summer' for short), which moreover are assumed to be in balance-neglecting seasonal storage of air and ice on account of their small thermal inertia.Assuming north-south symmetry of the ice sheet over the plateau, only half of the accumulation would descend the southern slope until it is depleted by ablation, which defines the glacial line (thick dashed line).The aim is to determine its temperature, referred to as the GMT.I first link the moisture transport to the energy transport, which is assumed to be dominated by eddies as katabatic-type overturning circulation is appreciable only in winter and confined to the slope [68] and since quasi-geostrophic eddies do not transport potential energy [37] (their Section 13.3.3.2), the energy transport consists mainly of sensible-and latent-heat, and with their ratio denoted by B (referred to as 'lateral Bowen ratio'), the moisture transport (same as the accumulation  ) can be linked to the latent-heat ( ) hence the atmospheric energy ( ) transports via Figure A1.North-south section of a plateau (of meridional span l) when the snowline (SL, thick solid) is lower than the plateau, so its overhead moisture transport would cause accumulation (with annual rate p), and the southward ice flux would be depleted by ablation at the glacial line (GL, thick dashed).The moisture transport crossing the snowline is linked to energy transport (F a ) by the Clausius-Clapeyron relation, which in turn balances the radiative cooling through the tropopause (q r ).The aim is to derive the marking temperature of GL.
I first link the moisture transport to the energy transport, which is assumed to be dominated by eddies as katabatic-type overturning circulation is appreciable only in winter and confined to the slope [68] and since quasi-geostrophic eddies do not transport potential energy [37] (their Section 13.3.3.2), the energy transport consists mainly of sensible-and latent-heat, and with their ratio denoted by B (referred to as 'lateral Bowen ratio'), the moisture transport (same as the accumulation A c ) can be linked to the latent-heat (F L ) hence the atmospheric energy (F a ) transports via is a 'moisture-content' parameter with ρ w being the water density and L s , the latent heat of sublimation.Since the vertical vapor profile tracks the saturated one, [69] has derived an analytical expression of the lateral Bowen ratio as a function only of the surface temperature T and its saturation vapor pressure e s (his Equation ( 23)).Applying the Clausius-Clapeyron (C-C) relation approximated by [70] (their Equation ( 25)) e s = 6.1 e 0.074T mb, (A3) we derive B = (3.86+ 0.014T)e −0.074T (A4) Over the snowline defined by the freezing point, (A4) yields B = 3.86, which agrees with its observed value [69] (his Figure 4), yielding µ = 2.3 × 10 −3 m 3 W −1 y −1 .I next link the atmospheric energy transport to the radiative forcing, which states This deduced GMT is consistent with the observed Greenland ice sheet [43] (their Figure 8).Although the parameters of (A12) are uncertain, the sensitivity is dampened by the ½ power dependence.More significantly, it underscores the primary control of the GMT by ice-sheet properties, so changing climate would displace snow-and glacial lines in tandem, which, however, is through temperature, not moisture transport.The short-circuiting of the moisture transport occurs because only the one crossing the snowline (defined by the freezing point) causes the accumulation, which is thus insulated from the SST.

Glacies 2024, 1 ,Figure 1 .
Figure 1.A model of warm/cold/Arctic oceanic boxes designated by numerals 1/2/3 separated by the subtropical front (30° N) and the Fram Strait (80° N).The SST (thick solid line) is anchored by the Arctic temperature to retain a linear profile in the cold box.The differential incoming SW flux ( * ), after atmospheric absorption and reflection by planetary albedo, is absorbed by the ocean ( ) to differentiate the SST, which in turn differentiates the SAT by the convective flux ( ).The latitudinal coordinate is defined as  = sin(latitude), so ocean-box boundaries lie at  = 0.5 and 0.98.

Figure 1 .
Figure 1.A model of warm/cold/Arctic oceanic boxes designated by numerals 1/2/3 separated by the subtropical front (30 • N) and the Fram Strait (80 • N).The SST (thick solid line) is anchored by the Arctic temperature to retain a linear profile in the cold box.The differential incoming SW flux (q *i ), after atmospheric absorption and reflection by planetary albedo, is absorbed by the ocean (q i ) to differentiate the SST, which in turn differentiates the SAT by the convective flux (q c ).The latitudinal coordinate is defined as y = sin(latitude), so ocean-box boundaries lie at y = 0.5 and 0.98.

Figure 2 .
Figure 2. Profiles of cold/Arctic-box summer SAT (doubling as the snowline, both tick-marked on the ordinate) to illustrate thermal thresholds as the Pliocene cools.The initial condition (labeled EP for the early Pliocene) is when the Tibet plateau rises above the snowline to initiate the Pliocene cooling.The first threshold (labeled GIS for Greenland ice sheet) is when the Arctic temperature falls below the glacial marking temperature ( , dashed) to fully glaciate Greenland.The second threshold (labeled NHG for Northern Hemisphere glaciation) is when the Arctic Ocean is cooled to the freezing point, so the polar cap (light-shaded) would form on the circum-Arctic lowlands and advance toward  .

Figure 2 .
Figure 2. Profiles of cold/Arctic-box summer SAT (doubling as the snowline, both tick-marked on the ordinate) to illustrate thermal thresholds as the Pliocene cools.The initial condition (labeled EP for the early Pliocene) is when the Tibet plateau rises above the snowline to initiate the Pliocene cooling.The first threshold (labeled GIS for Greenland ice sheet) is when the Arctic temperature falls below the glacial marking temperature (T g , dashed) to fully glaciate Greenland.The second threshold (labeled NHG for Northern Hemisphere glaciation) is when the Arctic Ocean is cooled to the freezing point, so the polar cap (light-shaded) would form on the circum-Arctic lowlands and advance toward T g .

Figure 3 .
Figure 3.Time evolution of the modeled global ( ), averaged cold-box ( ), and Arctic ( ) temperatures and uplift-induced ice cover (i) when thermal thresholds (dashed vertical lines) are aligned with abrupt cooling of the observed subpolar SST (thick solid line, adapted from [4]).The shaded portions are solution spans for global sensitivity ranging from 0.5 (thin solid line) to 1.5 (thin dashed line).The solid square is the estimated ice cover by uplift, which is seen sufficient to cause NHG.

Figure A1 .
Figure A1.North-south section of a plateau (of meridional span l) when the snowline (SL, thick solid) is lower than the plateau, so its overhead moisture transport would cause accumulation (with annual rate ), and the southward ice flux would be depleted by ablation at the glacial line (GL, thick dashed).The moisture transport crossing the snowline is linked to energy transport ( ) by the Clausius-Clapeyron relation, which in turn balances the radiative cooling through the tropopause ( ).The aim is to derive the marking temperature of GL.