Northern Hemisphere sea ice variability in a transient CGCM simulation of the past 2.6 Ma

Domain: earth_science
Published: 2025-01-02
Editorial outcome: accept

Abstract
The recent sea ice changes in the Northern Hemisphere (NH), necessitate elucidating the sea ice variability over the past 2.6 million years (Ma), when the Earth’s glacial cycles transitioned from ∼41 to ∼100 kyr periodicity, following the Mid-Pleistocene Transition (MPT) period (0.7–1.2 Ma). Here, we analyze a coupled general circulation model (CGCM) simulation to understand how the NH sea ice responds to changes in the transient orbital, greenhouse gas (GHG), and ice-sheet forcings. We find that the Earth’s axial tilt (obliquity) and axial wobble (precession) strongly influence the variability in high-latitude (>70° N) perennial sea ice and mid-latitude (35° N–70° N) seasonal sea ice, respectively, by modifying the net surface shortwave radiation. Meanwhile, the GHG forcing affects the glacial–interglacial sea ice predominantly in the Labrador Sea, Irminger–Iceland basin sector, and Central North Pacific regions during the MPT and post-MPT (0.0–0.7 Ma) periods by modulating the downwelling longwave radiation. Additionally, we confirm that variability with longer periodicity (∼100 kyr) from GHG and ice-sheet forcings is most pronounced in NH sea ice during the MPT and post-MPT periods.

Introduction
Sea ice is sensitive to climate variations and is critical for understanding past and future climate changes. Sea ice amplifies climate change in high–latitude regions through its high albedo and minimal thermal conductivity, together with other processes (e.g., temperature feedbacks, changes in water vapor, cloud cover, and atmospheric and ocean heat transport)1,2 that alter the radiative budget in the Arctic. Additionally, sea ice dynamics influence the mid-latitude jet stream3, oceanic salt and freshwater fluxes4, global intermediate/deep water circulation5, and North Atlantic Deep Water formation6. Over the past several decades, the Northern Hemisphere (NH) sea ice has been significantly reduced, in response to the increasing greenhouse gas (GHG) concentrations and global warming7. However, the responses of the NH sea ice to internal processes and external climate forcings in the past remain largely unknown. To comprehend historical and recent changes in NH sea ice8–11, understanding the responses of sea ice to external climate forcings and their feedbacks during the Quaternary (past 2.6 million years (Ma)) is crucial.

The Quaternary period is characterized by glacial–interglacial cycles. The Earth’s orbital parameters, namely (1) orbital eccentricity (∼100 thousand years (kyr)), (2) axial tilt (obliquity; ∼41 kyr), and (3) axial wobble or longitude of perihelion (precession; ∼23/19 kyr), significantly vary the latitudinal and seasonal distribution of incoming solar radiation12, which is the primary external forcing of the glacial–interglacial cycles13–15. Additionally, GHG concentration changes and variations of the extent and topography of the NH ice-sheets, together with their associated climate feedback mechanisms, amplify changes initiated by orbital variations16 especially by altering the global mean surface temperature (GMST)17 and influencing the NH mid-latitude westerlies18. The NH ice-sheets experienced cyclic waxing and waning, indicating glacial–interglacial cycles in response to GHG concentration changes19. These glacial–interglacial cycles transitioned from a relatively low amplitude and shorter periodicity (∼41 kyr cycle) during the Early Quaternary (1.2–2.6 Ma; pre-Mid-Pleistocene Transition (MPT)) to a high amplitude and longer periodicity (∼100 kyr cycle) in the Late Quaternary (0.0–0.7 Ma; post-MPT)19–22, following the MPT (0.7–1.2 Ma)23–27. Therefore, in this study, we divided the Quaternary into three distinct periods—post-MPT, MPT, and pre-MPT—to comprehend the evolution of the NH sea ice.

Paleoclimate modeling has shown that incoming solar radiation has driven the NH sea ice variations during interglacials over the past 0.8 Ma17. Additionally, various proxy records and modeling studies that are primarily confined to the post-MPT period identified ∼23/19 kyr (precession) and ∼41 kyr (obliquity) signals in the NH sea ice variability, with some emphasizing the role of precession28–30, and others highlighting the role of obliquity17,31–33. Recently, a transient paleoclimate simulation covering the past 3 Ma presented a general overview of the Northern and Southern Hemisphere sea ice variability, revealing precession, obliquity, and CO2/ice-sheet (∼60–100 kyr) signals and suggesting the combined effect of the GHG and ice-sheet forcings on sea ice variability34. Despite these findings, uncertainties remain regarding the spatiotemporal variability and mechanisms underlying the response of the NH sea ice to the orbital, GHG, and ice-sheet forcings. Furthermore, whether precession or obliquity strongly influences the sea ice remains unclear.

Here, we analyze the 3Ma experiment and explain the spatiotemporal variability and mechanisms of NH sea ice change during the Quaternary. The 3Ma experiment is a coupled atmosphere–ocean–sea ice simulation performed by Yun, et al.34 employing the Community Earth System Model version 1.2 (CESM1.2), forced by transient prescribed orbital35, GHG19 (including CO2, CH4 and N2O concentrations), and NH ice-sheet forcings19. Although the GHG concentrations and ice-sheets are not external forcings in the natural system, they are treated as such in the 3Ma experiment. Therefore, in this study, we describe them from an external forcing perspective. The 3Ma experiment provides an unprecedented view of the global Quaternary climate variability, as it is able to reproduce consistent glacial–interglacial variability36 and key paleo–proxy–based estimates34,36. Nevertheless, the 3Ma experiment design has a few caveats, such as the absence of prescribed ice–sheet melted freshwater forcing, a constant Last Glacial Maximum (∼21 ka) bathymetry, and an acceleration technique with a factor of five applied to the external (orbital, GHG, and ice-sheet) forcings.

Results
Evolution of Northern Hemisphere sea ice during the Quaternary
The prescribed orbital forcing35 (Fig. 1a–c) influence the intensity of summer (June–July–August; JJA) incoming solar radiation at 65°N37, which is the latitude at which the NH ice-sheets most sensitive to orbital forcing15,38. The incoming solar radiation at 65°N consistently varies throughout the Quaternary period (Fig. 1d) with predominant ∼23/19 kyr (precession) and ∼41 kyr (obliquity) signals (Supplementary Fig. 1). The Quaternary climate shows gradual cooling towards the post-MPT; however the trend of JJA incoming solar radiation at 65°N does not change considerably during the past 2.6 Ma (0.0001 Wm−2/kyr) (Fig. 1d). Meanwhile, the prescribed CO2 concentration19 (GHG forcing) (Fig. 1e) (Supplementary Fig. 2a–c) and prescribed global ice volume in sea level equivalent19 (ice-sheet forcing) (Fig. 1f) (Supplementary Fig. 2j–l) show a dominant ∼41 kyr periodicity during the pre-MPT period, ∼60–100 kyr (CO2/ice-sheet signal) periodicity during the MPT period, and a longer ∼100 kyr periodicity during the post-MPT period, consistent with the deep sea proxy record estimates21 (Supplementary Fig. 3).

We investigated the NH (35°N–90°N) sea ice extent (ocean surface covered with more than 15% sea ice concentration) during the post-MPT (0.0–0.7 Ma), MPT (0.7–1.2 Ma), and pre-MPT (1.2–2.6 Ma) periods. Evidently, the simulated sea ice extent increases throughout the MPT and post-MPT periods (Fig. 1g) and expands equatorward, particularly in the West Pacific and Atlantic Ocean sectors (Fig. 2a, b). Here, we define >70°N as high-latitude region and 35°N–70°N as mid-latitude region. Furthermore, perennial sea ice, i.e., sea ice that survives during summer and persists for over a year, is defined by the region north of the 90% contour during the MPT and post-MPT (purple line in Fig. 2a, b) and also corresponds to the region north of the 80% contour (purple line in Fig. 2c) during pre-MPT (cf. Supplementary Fig. 4). Evidently, the high-latitude region is primarily covered with perennial sea ice (north of the purple line in Fig. 2) throughout the Quaternary, whereas, the mid-latitude region is covered with seasonal sea ice (south of the purple line in Fig. 2), i.e., sea ice that forms in winter and melts in summer.

Evidently, the dominant periodicity of sea ice extent changes from ∼41 kyr during the pre-MPT period to a longer ∼100 kyr during the post–MPT period (Fig. 3a–c), which is aligned with climate shifts observed across the MPT21,39 (Supplementary Fig. 3). Spectral analysis of the sea ice extent indicates a latitudinal dependency of sea ice variability throughout the Quaternary (Fig. 3a–c). Meanwhile, the simulated perennial (Fig. 1i) and seasonal sea ice exhibit a strong significant correlation during the Quaternary (r=0.86); however, variations are evident in their periodicity modes. To understand the sea ice variability, remembering that the ∼41 kyr (obliquity) and ∼23/19 kyr (precession) signals do not necessarily represent direct responses to orbital forcing is essential. This is because the prescribed CO2 (GHG) (Supplementary Fig. 2a–c) and ice-sheet (Supplementary Fig. 2j–l) forcings also have a major ∼41 kyr signal during the pre-MPT period and a minor ∼23/19 kyr signal throughout the Quaternary. However, longer periodicity of ∼100 kyr and ∼60–100 kyr (CO2/ice-sheet signal) would indicate a response to the GHG/ice-sheet forcings, as obliquity and precession have distinct ∼41 kyr and ∼23/19 kyr signals, respectively.

Our analysis shows that the pre-MPT high-latitude region covered with perennial sea ice (north of purple lines in Figs. 2c and 3c) has a dominant ∼41 kyr (obliquity) signal (Fig. 3c, f). Meanwhile, the pre-MPT mid-latitude region covered with seasonal sea ice (south of purple lines in Figs. 2c and 3c) has a dominant ∼23 kyr (precession) signal (Fig. 3c). These pre-MPT sea ice signals indicate strong influences of the orbital and GHG forcings, as ice-sheet forcing during the pre-MPT period is weaker than that during the post-MPT period. During the MPT period, perennial sea ice has a dominant ∼41 kyr (obliquity) signal, indicating the influence of the orbital forcing, and a major ∼63 kyr (CO2/ice–sheet) signal (Fig. 3b, e), indicating the influences of the GHG and ice–sheet forcings. The MPT seasonal sea ice shows varying ∼63–56 kyr (CO2/ice–sheet), ∼41 kyr (obliquity), and ∼23/19 kyr (precession) signals (Fig. 3b), indicating the influences of the orbital, GHG, and ice-sheet forcings. Further, the post-MPT period has a primary ∼100 kyr (eccentricity) signal throughout the NH latitudes (Fig. 3a), suggesting strong influences of the GHG and ice-sheet forcings. Meanwhile, the post-MPT perennial sea ice (Fig. 3a, d) also shows a secondary ∼41 kyr (obliquity) signal and seasonal sea ice shows a secondary ∼23 kyr (precession) signal, indicating the influence of the orbital forcing. Thus, our results indicate that NH sea ice variability with longer periodicity (∼100 kyr) during the MPT and post-MPT comes from the GHG and ice-sheet forcings which amplified the eccentricity forcing. Additionally, obliquity plays a stronger role in high-latitude perennial sea ice variability, and precession plays a stronger role in mid-latitude seasonal sea ice regions throughout the Quaternary.

Furthermore, the sea ice extent variability over the Atlantic and Pacific Ocean sectors is similar to the NH sea ice signals, particularly that of perennial sea ice (Fig. 4). However, significant differences are evident in the Atlantic sea ice extent. For example, the Atlantic sea ice extent has a dominant ∼23/19 kyr (precession) signal during the pre-MPT and MPT periods (Fig. 4b, c) and a ∼23/19 kyr (precession) signal (∼40% variance) during the post-MPT period (Fig. 4a) with respect to the NH sea ice variability. Consequently, compared to that of the Pacific Ocean sector, the Atlantic Ocean sector shows a stronger influence of the precession forcing.

One indicator of the spatiotemporal variability of sea ice is the first Empirical Orthogonal Function (EOF) of sea ice concentration during the post-MPT, MPT, and pre-MPT periods. The first EOF of sea ice shows the major glacial–interglacial sea ice changes (Fig. 5a–c). During the pre-MPT period, the major leading pattern shows sea ice variability in the Labrador Sea, especially along the coast of Newfoundland (Fig. 5c). During the MPT period, an additional major leading pattern is evident in the Irminger–Iceland basin sector, especially along the coast of Norway (Fig. 5b). These two regions have open-ocean deep convection sites that are also evident in this model simulation. Following the MPT (post-MPT) period, three major leading patterns are identified, with an additional sea ice variability emerging over the Central North Pacific (CNP) (Fig. 5a). Such a sequential increase in the number of major leading patterns of sea ice variability is related to the southward expansion of sea ice extent (black solid line in Figs. 2 and 5a–c) associated with a gradual decrease in CO2 concentration (GHG forcing) (Fig. 1e).

The EOF results reveal three major regions with strong variability in sea ice throughout the Quaternary. The Labrador Sea, Irminger–Iceland basin, and Central North Pacific regions, which are covered with seasonal sea ice, show a dominant ∼23/19 kyr (precession) sea ice signal during the pre-MPT and MPT periods (Fig. 5e, f, h, i, respectively), whereas the Labrador Sea and Irminger–Iceland basin show a ∼23 kyr (precession) signal (∼40% variance) during the post-MPT period (Fig. 5d), which closely resembles the post-MPT Atlantic sea ice variability (Fig. 4a). Our results show that perennial sea ice is sensitive to obliquity (Fig. 3). Meanwhile, as sea ice expands southward (Figs. 2 and 5a–c), the mid-latitude seasonal sea ice shows greater sensitivity to precession. Proxy studies based on marine sediment core MD01–2414, from the mid-latitude region in the central Okhotsk Sea (53°N, 149°E; Fig. 5a), reveal a dominant precession sea ice variability during the last 130 ka30. Meanwhile, sediment core HLY–6, from the high-latitude region in the western Arctic Ocean (78°18´N, 177°W; Fig. 5a), shows a significant obliquity sea ice variability and episodes of perennial sea ice coverage28. Additionally, the presence of detrital iron oxide mineral grains in the Arctic Coring Expedition (ACEX) composite core from the Lomonosov Ridge (87.5°N and 138.3°W) indicates that the Arctic has been covered with perennial sea ice over the past 14 Ma40. Therefore, the proxy records agree with our findings that the high-latitude Arctic is covered with perennial sea ice and is strongly influenced by obliquity. However, a secondary effect of precession is also evident in the model results during the MPT and post-MPT periods. In contrast, the mid-latitude region is covered with seasonal sea ice and is more sensitive to precession.

Atmospheric and oceanic temperature variability during the Quaternary
Present-day observations, reanalysis data, and the global climate model have shown a tight coupling (negative correlation) between annual mean surface air and ocean temperatures and sea ice41. The simulated surface air temperature and sea ice extent are significantly (negatively) correlated during all three time periods and the surface air temperature leads sea ice by 1 kyr during the pre-MPT and MPT periods and shows no lag during the post-MPT period (Fig. 6d–f). Notably, the simulated pre-MPT zonal mean surface air temperature shows a dominant ∼41 kyr (obliquity) signal in the high-latitude region, while the mid- and low-latitude regions (0° N–30° N) show both ∼41 kyr (obliquity) and ∼23/19 kyr (precession) signals (Fig. 6c). Also, the cospectrum between surface air temperature and sea ice during the pre-MPT period shows strong ∼41 kyr and ∼23 kyr signals (Fig. 6i), indicating strong influences of the orbital and GHG forcings, as the pre-MPT ice-sheet forcing is weaker than the post-MPT ice-sheet forcing.

The MPT surface air temperature shows a major ∼41 kyr (obliquity) signal in the high- and low-latitude regions and a dominant ∼100–56 kyr (CO2/ice-sheet) signal in the mid- and high-latitude regions (Fig. 6b). Meanwhile, the cospectrum shows a strong CO2/ice-sheet (∼100–56 kyr) signal alongside ∼41 kyr and ∼23/19 kyr signals (Fig. 6h). The presence of a dominant CO2/ice-sheet signal during the MPT period suggests a more robust response of the surface air temperature to the GHG and ice-sheet forcings than that of sea ice. Moreover, similar to the case of sea ice variability, the post-MPT period has a dominant ∼100 kyr (eccentricity) signal throughout the NH latitudes (Fig. 6a), a ∼41 kyr signal in the high-latitude region, and a ∼23 kyr signal in the mid- and low-latitude regions (Fig. 6a), indicating the influences of the orbital, GHG and ice-sheet forcings.

Furthermore, the simulated surface air and ocean temperatures have similar latitudinal periodicity variations (Supplementary Fig. 5a–c) and are significantly (positively) correlated with no phase lag (Supplementary Fig. 5d–f) throughout the Quaternary. Our results suggest that the surface air and ocean surface temperatures are tightly coupled at orbital and millennial timescales. Consequently, oceanic heat released to the atmosphere sustains a strong link between sea ice and surface air temperature variability.

Response to the orbital forcing during the Quaternary
The net surface shortwave radiation over the ocean (positive downward) in the mid-latitude region shows a predominant precession (∼23/19 kyr) signal (Fig. 7a–c), whereas the high-latitude region shows a distinct obliquity (∼41 kyr) signal throughout the Quaternary (Fig. 7a–c). Consequently, these changes in the net surface shortwave radiation signals influence the surface air temperature and sea ice variability. Additionally, the response to the obliquity forcing in the high-latitude region is further amplified by the ice–albedo feedback17. Next, we applied obliquity and precession bandpass filters (30–70 and 5–30 kyr, respectively) on the simulated summer incoming solar radiation, surface air temperature, and sea ice extent to investigate the responses of sea ice and surface air temperature variability to the orbital forcing (obliquity and precession). The bandpass-filtered summer incoming solar radiation at 65°N is significantly correlated with the obliquity (positively) and precession (negatively) indices. Furthermore, summer incoming solar radiation at 65°N does not lag the obliquity index (purple line in Fig. 7d) during the Quaternary; however, it lags the precession index by 4 kyr (purple line in Fig. 7e). Therefore, an increase in the obliquity index or a decrease in the precession index increases the NH incoming summer solar radiation at 65°N. Specifically, perennial and seasonal sea ice are significantly correlated with the obliquity (negatively) and precession (positively) indices, respectively. Perennial sea ice lags the obliquity index by ∼5 kyr (green line in Fig. 7d), and seasonal sea ice lags the precession index by ∼10 kyr (green line in Fig. 7e). Hence, the NH mid-latitude region covered with seasonal sea ice is strongly influenced by the precession forcing, whereas the high-latitude region covered with perennial sea ice is more strongly influenced by the obliquity forcing throughout the Quaternary.

Furthermore, the Earth’s orbital forcing has a ∼100 kyr (eccentricity) signal which modulates the strength of the precessional cycle. Naturally or in a more realistic model setup with dynamical carbon cycle and ice-sheets, the ∼100 kyr signal in post-MPT sea ice and surface air temperature would have been explained by eccentricity, amplified by feedbacks involving GHG and ice-sheets. However, without the amplification of these physical quantities through feedbacks within the Earth’s climate system, the variability of sea ice due to the orbital forcing alone is likely to be relatively weaker. Meanwhile, this experiment treats GHG and ice-sheet as external forcing. Therefore, a more reasonable interpretation of the model results is that the signal of the orbital forcing was amplified through feedbacks within the climate system, and these consequent GHG and ice-sheet changes more effectively influenced sea ice variability.

Response to the GHG forcing during the Quaternary
The Earth’s surface and lower atmosphere are warming in response to increased GHG concentrations. As the temperatures increase, more water vapor is released into the atmosphere, intensifying the greenhouse effect by causing increased emissions of longwave radiation from the atmosphere to the surface, resulting in additional warming43 (i.e., the water vapor feedback44). Cloud cover also affects downwelling longwave radiation43. Here, we show that the changes in prescribed CO2 concentration can explain the major regions of glacial–interglacial sea ice variability (Fig. 8a–c), as shown in Fig. 5a–c. Furthermore, changes in downwelling longwave radiation (clear-sky conditions) caused by GHG and low-level water vapor changes explain the glacial–interglacial major sea ice regional changes in the (a) Labrador Sea, Irminger–Iceland basin and Central North Pacific during the post-MPT period (Fig. 8d), (b) Labrador Sea and Irminger–Iceland basin during the MPT period (Fig. 8e), and Labrador Sea during the pre-MPT period (Fig. 8f) as shown in Fig. 5a–c. Additionally, changes in downwelling longwave radiation explain the intensification of sea ice extent during the MPT and post-MPT periods (Fig. 8d, e).

In contrast to the orbital forcing, the GHG forcing has a uniform zonal distribution. The prescribed CO2 (GHG) concentration (Supplementary Fig. 2a–c) and the simulated GMST (Supplementary Fig. 6) have similar periodicity signals throughout the Quaternary. Meanwhile, the pre-MPT prescribed CO2 concentration (Fig. 1e) and the simulated GMST (Fig. 1h) show gradual (significant) downward trends, which intensify across the MPT by approximately two and three times, respectively. In contrast, the simulated sea ice extent shows an upward (significant) trend during the pre-MPT (1.5 thousand km2/kyr), which also intensifies by a factor of approximately 2.2 (3.4 thousand km2/kyr) during the MPT period (Fig. 1g), indicating the influence of the GHG forcing on sea ice. Proxy sea ice records from sediment core U1343 (57° 33.4′N, 176°49.0′W; 1,950 m) located in the Bering Sea shelf region are in agreement with our result, highlighting a significant increase in sea ice extent throughout the MPT period27,45.

The prescribed CO2 concentration (Fig. 1e) and simulated GMST (Fig. 1h) are significantly (positively) correlated with no lag during the Quaternary (Fig. 8g). Consequently, this indicates that as GHGs and GMST increase simultaneously, warmer summers and enhanced sea ice melting occur. Furthermore, the NH sea ice extent with GMST and high-latitude/perennial sea ice with high-latitude surface air temperature are significantly (negatively) correlated with no lag during the Quaternary (Fig. 8g). Therefore, our results show that compared to that during the post-MPT period, the higher (26%) prescribed GHG forcing during the pre-MPT increased GMST by 2.2°C (21%), causing warmer seasons. Consequently, the pre-MPT period would have enhanced summer time sea ice melting and reduced winter time sea ice formation, causing a comparably smaller extent of sea ice (Fig. 2c). In contrast, during the post-MPT period, the relatively lower GHG forcing and GMST facilitate a comparably greater extent of sea ice (Fig. 2a) because of somewhat colder seasons.

Furthermore, our analysis shows that the GHG forcing can largely explain the simulated NH sea ice variability, while a secondary influence of precession forcing is also evident during the Quaternary (Table 1). Notably, the obliquity forcing strongly influences high-latitude/perennial sea ice, and a secondary influence of the precession forcing cannot be ignored. Additionally, the orbital indices lead to sea ice change, while GHG changes lag behind sea ice change, suggesting that changes in the orbital forcing initiate sea ice change, which is then amplified by the GHG forcing.

Methods
Model and simulation
The 3Ma experiment34 is a transient paleoclimate model simulation covering the past three million years, employing the coupled CESM version 1.2 (CESM1.2). The experiment comprises fully coupled atmosphere (Community Atmosphere Model version 4; CAM4), ocean (Parallel Ocean Program version 2; POP2), sea ice (Los Alamos Sea Ice Model; CICE4), and land (Community Land Model version 4; CLM4.0; prescribed vegetation) components. The atmospheric model has latitude–longitude horizontal resolution (T31) of approximately 3.75°×3.75° with 26 vertical levels. The ocean and sea ice models have a 100×116 horizontal grid (gx3v5), with the ocean model comprising 25 vertical levels. The ocean and sea ice output data are regridded and analyzed at latitude–longitude horizontal resolution of 3.75°×3.75°, identical to that of the atmospheric outputs. NCAR Command Language (NCL)66 version 6.6.2 was used to analyze data and generate all figures.

The 3Ma experiment is forced by transient prescribed orbital35 forcing, which changes the astronomical insolation, ice-sheet19 forcing, which changes the NH ice-sheet extent and topography, while the Southern Hemisphere ice-sheet is kept constant, and GHG (including CO2, CH4 and N2O concentrations)19,21 forcing, which changes the radiative forcing. The GHG forcing includes CO2 concentration estimates from various datasets. From 0.8 to 3.0 Ma, the estimated GHGs consist of simulated CO2 concentration from a previous CLIMBER-2 transient simulation19 containing an interactive carbon cycle. Meanwhile, from 0.0 to 0.8 Ma, the GHGs estimate includes CO2 concentration measurements from air bubbles from the European Project for Ice Coring in Antarctica (EPICA) Dome C ice core67. The GHG forcing also includes N2O and CH4 concentration estimates as a CESM forcing obtained by regressing the EPICA Dome C ice core N2O and CH4 concentrations on the LR04 stack21 deep sea proxy record post 0.8 Ma and extending the linear statistical model to 3Ma using the LR04 stack deep sea record. The LR04 stack21 deep sea proxy record is constructed by the graphic correlation of an average of globally distributed 57 benthic foraminiferal calcite δO18 records for the past 5.3 million years (unit is in permil), and it measures the global ice volume and deep ocean temperature. The ice-sheet topography changes are obtained from the simulated CLIMBER-2 orography data from the previous transient simulation19. The 3Ma experiment also uses constant ocean bathymetry from the Last Glacial Maximum (∼21 ka) and a prescribed vegetation component, where the plant functional types (i.e., vegetation) are changed only in the NH land ice-covered areas.

The 3Ma simulation uses an acceleration factor of five for the external (orbital, ice-sheet, and GHG) forcings, which reduces the 3Ma orbital year forcing to 0.6 Ma model years. The orbital, ice-sheet, and GHG forcings are updated every 100 orbital years, corresponding to 20 model years. The CESM1.2 model in the 3Ma experiment resolution exhibits a comparatively low standard equilibrium climate sensitivity of 2.4°C warming for a two-fold CO2 increase with respect to the pre-industrial CO2 concentration. This is approximately 1.5 times smaller than the likely range of estimates (3.7±1.2°C) obtained within 37 other climate model simulations conducted as part of the Coupled Model Intercomparison Project Phase 6 (CMIP6)68. Therefore, to achieve a more realistic response to past longwave radiative forcings and indirectly account for the radiative impacts of other CO2-correlated forcings from vegetation, dust, N2O, or CH4, the CO2 forcing19 anomalies relative to preindustrial conditions are amplified by a factor of 1.5 in the 3Ma experiment. The resulting effective equilibrium climate sensitivity is ~3.8°C, which is quite consistent with the CMIP6 estimate. Additionally, it is in reasonable agreement with prior paleoclimate estimates (i.e., 3.2°C)69,70, which were derived from reconstructions of GMST and radiative forcings during the past 0.784 Ma. The re-scaling of the CO2 forcing in CESM1.2 resulted in a more realistic representation of the magnitude of tropical, global mean, and Antarctic temperature changes34,36.

Furthermore, the 3Ma experiment does not explicitly simulate the climate–carbon cycle feedback. In reality, the GHG forcing is a response and byproduct of the feedback from the carbon cycle. Therefore, an interactive carbon cycle model would be more realistic. However, under the current model’s capability, including a carbon cycle model, may produce other model biases. Therefore, by prescribing the GHG forcing, which involves CO2, CH4 and N2O concentrations, we can assess its impact on sea ice without invoking a complex causality in the Earth’s carbon cycle system. Additionally, because the GHG and ice-sheet forcings are treated as external climate drivers, the interactive feedback processes between the GHGs and ice-sheets were not incorporated. Incorporating these interactions would provide a more physically realistic approach. However, achieving a balanced representation of climate change associated with the natural carbon cycle and the climate feedback resulting from ice-sheet dynamics presents significant challenges at this stage. Therefore, the initial approach was to consider these factors as climate forcing to minimize potential biases in the 3Ma experiment.

The 3Ma simulation data were compared against large-scale climate reconstructions and various paleo–proxy datasets. The 3Ma simulation has been shown to reproduce the magnitude and timing of glacial–interglacial variability34,36 consistent with paleoclimate data constraints and estimates from the intermediate complexity Earth system model. The 3Ma simulation agrees well34,36 with key paleoclimate proxy records, such as the magnitude and timing of GMST variations, annual mean ocean surface temperatures over the tropics and elsewhere (western South and North Pacific, southern Indian, and North Atlantic), Greenland and Antarctic temperatures, the East Asian summer monsoon, and the East Africa hydroclimate.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Dear Gagan Mandal et al., please find comments to your manuscript in the enclosed review document.\n(Remarks on code availability)\nAs far as I see from the submitted files no codes are included in the submission.", "reviewer_2": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature\nCommunications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career\nResearchers who co-review manuscripts.", "reviewer_3": "In this paper, Mandal et al. have analyzed the Northern Hemisphere sea ice variability in a global climate model experiment\nthat covers the last 3 million years. To my knowledge, this experiment is the first in its sort and provides an unprecedented\nview of global Quaternary climate variability. In addition, sea ice is a crucial component of the climate system and it is\nimportant to improve our understanding of its variability at all relevant time scales. The topic of this paper is therefore timely\nand certainly of potential interest to a broad readership. Hower, as discussed in detail below, there are a few concerns that\nkeep me from advising to publish this paper in Nature Communications.\nMain comments\nFirst, the main results of this experiment have already been published in extensive paper by Yun et al. (2023), including the\nNorthern Hemisphere sea ice variability and an analysis of the impact of different orbital parameters on sea ice, temperature\nand CO2 through a wavelet power spectra. Yun et al. have also already looked at differences in spectral variance before and\nafter the Mid-Pleistocene Transition, for instance in Greenland temperature.\nSecond, the design of the experiment makes it less suitable for a detailed analysis as performed by Mandal et al. This has to\ndo with three issues:\n1) In the experiment, an acceleration technique with a factor 5 is applied to the external forcing (orbital, greenhouse gas\nlevels and ice-sheets). This implies that in each simulated model year the external forcings are refreshed with a time-step of\n5 years, so that in total 600 000 model years have been run instead of 3000 000 years. As mentioned by Yun et al. (2023), previous\nstudies have shown that this acceleration technique causes biases in deep ocean variables and some surface\nocean variables at high latitudes that are linked to deep water formation and the Atlantic Meridional Overturning Circulation\n(AMOC). Sea ice cover is one of these surface variables.\n2) In addition, the experiment does not include the impact of variations in meltwater fluxes from the North American,\nEurasian and Greenland ice sheets. These meltwater fluxes have varied considerably on the time scale considered here,\nwith a direct effect on the sea surface salinity and on the deepwater formation, and thus the AMOC. Yun et al. (2023) note that\nignoring the meltwater fluxes likely caused an underestimation of the AMOC variability. This underestimation has\nimportant consequences for the variability in NH sea ice cover.\n3) Moreover, in the experiment the land-sea distribution is fixed at Last Glacial Maximum (LGM) conditions, implying that\nareas like the Barents Sea and the Bering Strait are permanently defined as land in the model. One consequence of fixing\nthe land-sea mask at the LGM is therefore that the area available for sea ice formation is underestimated during warmer\ninterglacial phases with a higher sea level. Yun et al. (2023) also mention that using the LGM land-sea mask can in principle\ngenerate large-scale climate biases in their experiment.\nThese three issues and their consequences for the analysis are not yet discussed by Mandal et al. In my opinion, this is a\nserious weakness of the paper, limiting the usefulness of the presented sea ice variability analysis. For instance, in my view\nit does not make much sense to perform a detailed analysis of sea ice variability in the Labrador Sea and Irminger-Iceland\nbasin, as these are key regions for deepwater formation and it can be expected that sea ice variability is affected here by the\ndiscussed issues.\nOther comments\nPage 3, first paragraph. “The sea ice amplifies climate change in high–latitude areas, known as Arctic amplification, while its\nhigh albedo and minimal thermal conductivity influence energy budgets at low and high latitudes.” This suggests that Arctic\namplification is solely caused by sea ice effects. However, this a too limited view, as there are also processes involved, for\ninstance changes in atmospheric heat transport and variations in cloud cover and water vapor that alter the radiative budget\nin the Arctic (e.g., Serreze and Barry, 2011).\nPage 2. I suggest explaining the meaning of MPT here.\nPage 3, line 68. I suggest clarifying that greenhouse gas (GHG) forcing also involves other important gases besides CO2,\nmost notably CH4 and N2O. The latter two greenhouse gases are also included in the GHG forcing in the analyzed\nexperiment.\nPage 3/4, starting line 70. “This period has seen the cyclical waxing and waning of the NH\ncontinental ice sheets accompanying the variations in global mean surface temperature that follow the aforementioned orbital periodicities (∼100\nkyr, ∼41 kyr, and ∼23/19 kyr). Despite these changes, the global mean radiation balance changed by less than ∼0.5 Wm−2,\nwith no noticeable variation in the key orbital parameters.” I find the last sentence confusing, as the orbital parameters clearly\nvaried in time, as also mentioned in the previous sentence. I suggest rephrasing.\nPage 5, line 96. I suggest making clearer here that this experiment has been performed by others and that many results for\nthis experiment were published previously in Yun et al. (2023).\nPage 8, 9. What is causing the “transition periodicity” at 63-56 kyr? I suggest to clarify.\nPage 13, line 245. Figures 3j-l are missing.\nPage 16, line 295. “The prevalence of the predominant transition periodicity during the MPT period suggests a more robust\nresponse of SAT to GHG and continental ice sheet forcings than sea ice.” Please explain why a clear response to the 63-56\nkyr transition periodicity suggest a more dominant response of SAT to GHG and ice sheet forcing.\nPage 25. I propose mentioning here that the GHG forcing is amplified by a factor 1.5 to compensate for a relatively weak\nclimate sensitivity in the appliec model.\nPage 26, line 507. The authors conclude from Supplementary Figure 8 that the modelled sea ice cover in the Arctic over that\npast 120 ka is in good agreement with a proxy for Atlantic sea ice derived from the NEEM Greenland ice core. In my opinion\nthe agreement is rather poor and cannot be used to conclude that the model is performing well in simulating Arctic sea ice\nvariations for the last 3 Ma.\nAdditional references.\nSerreze, M.C., Barry, R.G. (2011) Processes and impacts of Arctic amplification: A research synthesis, Global and Planetary\nChange 77, 85-96.", "reviewer_4": "Review\nNorthern Hemisphere sea ice variability in a transient CGCM simulation of the past 2.6 Ma\nMandal et al\n--------------------------------------------------------------------------------------------------------------------\nKey results\nThis work uses results from a fully coupled general circulation model (CGCM) to investigate the influence of orbital (orbital\neccentricity, axis tilt or obliquity, and axis wobble – or precision), atmospheric carbon dioxide (CO2), and changes in\ncontinental ice sheets on Northern Hemisphere (NH) sea ice from just before the Mid-Pleistocene Transition (MPT) from\n(~2.6 Million years ago (Ma) through the MPT and afterwards, to the present. There are 4 key findings in this work: 1) NH sea\nice is largely determined by global mean surface air temperature; 2) variability in high latitude (>70°N) Arctic perennial sea\nice (or “quasi-permanent sea ice” as defined by the authors) is strongly influenced by the tilt in the earth’s axis; 3) lower-\nlatitude seasonal NH sea ice is more heavily influenced by the earth’s wobble (precision); and 4) low-frequency variability\n(~100 thousand years, ka) from CO2 and continental ice shelves is most pronounced in NH sea ice post-MPT and at high\nlatitudes in the perennial sea ice area during the MPT.\nValidity\nThe data and analysis w.r.t. power spectra, EOF analysis are valid and robust. The interpretation and conclusions from the\nanalysis are challenging for me to evaluate (see additional comments below).\nSignificance\nI am not intimately familiar with current paleoclimate research (I am not a paleo scientist – see below re. my expertise), so I’m\nnot entirely sure if these results are “new” or “significant”. It seems upon my reading of this submission along with a reading\nof Yun et al., 2023 (Y2023) that some of the findings were already presented in Y2023 (from the same model simulations).\nFor example, Y2023 states that sea ice variability (both NH and SH) show a shift in power spectrum in the post-MPT time\nperiod with lower frequency variabilities of CO2 and continental ice shelves playing a larger role. Additionally, that sea ice is\ninfluenced by orbital forcings is not new (e.g. Howell et al., 2016 and references therein).\nA primary difference between this submission and Y2023 is the presentation of NH sea ice in different latitudinal/seasonal\nregimes: a “quasi-permanent” sea ice area (where sea ice concentrations lie above 80% pre-MPT and above 90% MPT and\nlater) and a seasonal ice pack with lower concentrations and generally south of 70°N. That sea ice variability the perennial\nand seasonal ice regions may be influenced by different orbital parameters is new to me.\nData and methodology\nThe approach of calculating and analyzing power spectra, EOFs and also looking at NH sea ice from the perspective of\nperennial (or “quasi-permanent”) sea ice vs. seasonal sea ice is valid. The results, as presented, read as a narrative that\ngives quantitative results (e.g. summarizing the values given in Figures, etc.). There are two things I find missing from this,\nhowever: 1) a narrative/storyline that helps tie everything together using Figures as supporting material\nand data presentation for the narrative, and 2) mechanisms and context. For example, the primary mechanism (as I understand it)\nthrough which orbital parameters can influence sea ice, global air temperatures, etc. is through incoming radiation (e.g.\neccentricity will change amount of solar radiation earth receives), geographical location of orientation (e.g. tilt in the earth’s\naxis will directly impact seasonal changes), and heat transport (atmospheric, oceanic) of the incoming radiation, yet there is\nno discussion of these mechanisms. It is intriguing that different orbital parameters may influence sea ice in perennial/high\nlatitude and seasonal/lower latitudes differently. Why?\nSuggested improvements\nThere is very little description of the climate model or forcings, other than to refer the reader to Y2023. There should be\nat least a basic summary of the model and the forcings for readers not intimately familiar with paleoclimate forcings and climate\nmodels. Details can of course be left to Y2023, but context for the results presented here would help (I had to read Y2023\nmore than usual for a reference material in order to understand how the climate simulation works).\nSea ice in the climate system plays important roles in, and responds to, transfer of energy from the sun to the earth. It serves\nto reflect incident radiation back (albedo), influences conductive heat fluxes at the surface, etc. As part of a global energy\nbudget, sea ice is important in volume and not just surface area and yet no mention of sea ice thickness or volume is made\nin this work. The focus can remain on sea ice concentrations at the surface – however it would be best to include sea ice\nvolume/thickness at least in discussion if not also within analysis.\nSmaller editorial suggestions:\n1. Define “MPT” before using it (e.g. abstract)\n2. Give some non-disciplinary labels to “precission, obliquity and eccentricity” early on to help non-paleoscientist readers\n(e.g. orbital eccentricity, axis tilt or obliquity, and axis wobble – or precision)\n3. Line 50 Sea ice loss is indeed an important mechanism for Arctic amplification of current global warming. It is not the only\nmechanism (particularly if one is considering AA at the top of the atmosphere rather than the earth’s surface; eg. Pithan &\nMauritsen 2014), although this sentence implies sea ice is THE mechanism. Reword.\n4. Define “LR04 stack” and give context/description? This applies to Supplementary Figure 3 as well)\n5. “Quasi permanent” sea ice is confusing. The analysis investigates sea ice variability – and “permanent” suggests no\nvariability. It is confusing to read about “variability” of “quasi-permanent” sea ice. Suggest “seasonal sea ice” for sea ice that\nforms in winter, melts in summer, and “perennial sea ice” rather than “quasi-permanent”…..permanent is a contradiction to\ninvestigating variability? I suggest changing to”perennial sea ice” and “seasonal sea ice” (such as in Vernal & Hillaire-\nMarcel) although that’s just a suggestions and authors may choose other descriptors.\n6. Figure 1. Sea ice concentration (SIC) contours in j-l suggest that SIC goes below 0! Suggest changing contour\nvalues/scale such that scale goes to 0 not -20%\n7. Line 167 “Considering that monthly one million km2 is widely accepted as the threshold of ice-free Arctic38, we found that\nthe NH high-latitudes (∼70°N–90°N) retained annual mean sea ice cover as it was persistently covered by quasi-permanent\nsea ice throughout the Quaternary (Fig. 1h)” This sentence does not make sense to me….and not sure how the definition of\nan “ice free Arctic” is relevant to this work?\n8. Line 243 only here does the reader learn that Bering strait is closed. This can significantly impact ocean heat transport\n(e.g. Dörr et al, 2024 and references therein) I suggest showing a map of coastline or land or ocean so the reader can\nassess how this is different than present day, and that the Bering strait is closed (perhaps in Figure 1 j-lwhere the Bering\nStrait appears open?)\n9. Line 255 Typically sea surface temperatures directly below sea ice are at freezing temperature (or close) and have little to\nno variability. How then does one assess variabilities in sea ice covered areas?\n10. Line 390 Sea ice extent leads GHG forcing? This is interesting….mechanism/why/how? This actually contradicts one\nconclusion that sea ice is determined in large degree by GHG?\n11. Line 440 The mid-latitude radiation and relative heat capacities of land/ocean influence mid-latitude response as well?\nWouldn’t the geographic locations due to land/ocean in midlatitude in NH play a part too?\n12. Line 471 “high-latitude sea ice varied coherently with periods of intensified NH glaciation during early pre–MP and MPT ,\nindicating the influence of NH continental ice-sheet forcings“ . Correlation does not imply causation. Have continental ice\nsheet changes caused changes in sea ice? Or are they both influence by radiative heating and NH heat transfer? Or? Why?\n13. Figure 5. Define dashed lines in correlation plots (d-f)?\n14. Reference 30 (Yu2023) has a link that is not valid. Correct.\n15. Supplementary Figure 4 Average SIC at 72°N missing labels a, b, c – and two dashed lines are in the figure? There is no\n70°N latitude line? Or is everything north of 70°N considering “high latitude” and south “mid-latitude” (for sea ice)?\nClarity and context\nThe acronyms and abbreviations used in this submission are extensive – so much so that I had to write many of them down\nand make a chart to refer to as I read. I strongly suggest reducing some of the discipline-specific abbreviations and/or some\nplain-English wording in the manuscript to help non paleoscientist readers.\nReferences\nDörr, J. S., Årthun, M., Eldevik, T., & Sandø, A. B. (2024). Expanding influence of Atlantic and Pacific Ocean heat transport\non winter sea-ice variability in a warming Arctic. Journal of Geophysical Research: Oceans, 129(2), e2023JC019900.\nhttps://doi.org/10.1029/2023JC019900\nF.W. Howell, A.M. Haywood, H.J. Dowsett, S.J. Pickering. Sensitivity of Pliocene Arctic climate to orbital forcing,\natmospheric CO2 and sea ice albedo parameterization\nEarth Planet. Sci. Lett., 441 (2016), pp. 133-142, 10.1016/j.epsl.2016.02.036\nPithan, F., & Mauritsen, T. Arctic amplification dominated by temperature feedbacks in contemporary climate models.\nNature Geosci 7, 181–184 (2014). https://doi.org/10.1038/ngeo2071\nVernal, A. & C. Hillaire-Marcel. Quaternary Arctic sea-ice ocver: Mostly perennial with seasonal openings during\ninterglacials. Past Global Changes Magazine, 30 (2), 94-95 (2022), https://doi.org/10.22498/pages.30.2.94\nMy expertise\nMy area of expertise lies in sea ice, sea ice in the climate system, and predominantly through the historical period and in\nclimate model future simulations. I do not consider myself a paleo-scientist. My understanding of paleoclimates, forcing and\nvariability comes from a relatively broad educational background as a climate scientist. I am not, for example, familiar with\nspecific proxy records of paleoclimate, or sources of prescribed climate forcings used in climate model simulations of\npaleoclimates. As such I found I had to devote a significant amount of time and effort into skimming background material\nabout the model and forcings used here, particularly from Yun et al., 2023, in order to put the results presented in this work in\ncontext.", "reviewer_1_attached": "Response to Reviewer #1 and #2 comments\n1\n2 We sincerely thank you for your valuable comments and suggestions on our manuscript. The\n3 comments and suggestions helped revise and improve our paper and the significance of our\n4 research. We studied the comments carefully and revised the manuscript as suggested. We have\n5 also improved our writing and grammar, making the manuscript clearer and scientific. We hope\n6 it meets your expectations with approval.\n7 The primary modifications in the paper and the responses to your comments are as follows:\n8 (a) Original reviewer comments are in italics,\n9 (b) Our response to the comments/suggestions are highlighted in blue color,\n10 (c) References from the revised manuscript are in red.\n11 Reviewer 1 and 2 comments:\n12 I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This\n13 is part of the Nature Communications initiative to facilitate training in peer review and to\n14 provide appropriate recognition for Early Career Researchers who co-review manuscripts.\n15 As far as I see from the submitted files no codes are included in the submission.\n16 Response: Thank you for all the suggestions. We have mentioned in the revised manuscript\n17 that all the code and data are available upon request. Meanwhile, some of the 3 Ma experiment\n18 data is available on the OPeNDAP ICCP climate data server https://climatedata.ibs.re.kr (3Ma\n19 transient simulation with CESM1.2; Dataset doi: 10.22741/iccp.20230001).\n20 Line #556-558: \"The code used in this study is available from the authors upon request. The 3\n21 Ma experiment data is available on the OPeNDAP ICCP climate data server\n22 https://climatedata.ibs.re.kr (3Ma transient simulation with CESM1.2; Dataset\n23 doi: 10.22741/iccp.20230001).\"\n24\n25 General comments\n26 In general, please carefully check statements that you cite from previous work, also in those\n27 cases where comments below do not explicitly raise that point. A rough scan over some\nPage 1 of 48\n28 references suggests that statements attributed by authors to previous publications are not\n29 always accurate.\n30 Response: We appreciate the concern. We have extensively edited the manuscript to improve\n31 clarity. We also made sure that references are properly cited.\n32 In particular the last sentence of the abstract sounds particularly odd if formulated in past\n33 tense. Is there a good reason not to use present tense? This statement applies also to other parts\n34 of the manuscript.\n35 Response: We are grateful for the insight. We have removed redundant texts from the revised\n36 manuscript and have formulated the sentences in the present tense.\n37 The employed model geography could be presented in more detail and more prominently. The\n38 reader can find some hints here and there, but a comprehensive description that is easy to\n39 locate is still missing. At the moment it is hidden in lines 240ff as an addendum to the results\n40 and via some remarks in the summary and discussion. This is probably not appropriate towards\n41 a proper understanding of the model simulation and of the model's capabilities and limitations.\n42 Furthermore, major simplifications appear to have been made with regard to paleogeography\n43 employed as a boundary condition. This reviewer understands that throughout the simulation\n44 the land sea mask has been as for the LGM. This decision should certainly be\n45 motivated/justified in the manuscript and the impacts on interpretation of modelling results,\n46 that will certainly be substantial, must be critically discussed. Keeping the Bering Strait closed\n47 (see line 243) and the Barents Sea shelf subaerial (Figure 1) throughout the whole simulation\n48 will at times lead to sea ice being influenced by the various forcings and dynamics in a different\n49 way than it has been the case in reality, where land sea masks continuously changed. No model\n50 setup will be perfect, but the feeling is that more discussion should explore potential\n51 ramnifications.\n52 Response: We are grateful for the insight. The 3 Ma model experiment uses a constant ocean\n53 bathymetry and land-sea boundary from the Last Glacial Maximum (∼21ka), which causes\n54 some unrealistic land-sea boundaries and has been discussed in our revised manuscript.\n55 Firstly, Figure 1j-l is redrawn as Figure 2a-c in the updated manuscript, which shows subaerial\n56 land regions in grey dots (especially subaerial Barents Shelf and closed Bering Strait) over\nPage 2 of 48\n57 modern land-sea boundaries.\n58 Line #133-140: \"\n59 Fig. 2: The climatology of sea ice during the post–MPT (0–0.7 Ma), MPT (0.7–1.2 Ma), and\n60 pre–MPT (1.2–2.6 Ma) periods. The shaded blue color shows the sea ice coverage over modern\n61 land-sea boundaries. The units are in percentage sea ice concentration. The grey dots show\n62 subaerial land regions (especially subaerial Barents Shelf and closed Bering Strait) as the\n63 experiment uses a constant Last Glacial Maximum land-sea boundary. The black and purple\n64 lines demarcate the sea ice extent and perennial sea ice boundary, respectively.\"\n65 Secondly, we have introduced the model limitations in the introduction section and have\n66 discussed the potential ramifications of using a constant Last Glacial Maximum (∼21ka) land-\n67 sea boundary in the discussion section.\n68 Line #85-89: \" Meanwhile, the 3 Ma experiment design has a few caveats, such as the absence\n69 of prescribed ice-sheet melt freshwater forcing, a constant Last Glacial Maximum (∼21 ka)\n70 bathymetry, and an acceleration technique with a factor of five applied to the external forcings\n71 (orbital, GHG, and ice-sheet). The impact of the experiment's design on the simulated sea ice\n72 is covered in the last section.\"\n73 We show that the high-latitude region in the 3Ma experiment has some unrealistic land-sea\n74 distribution, such as a closed Bering Strait and subaerial Barents Shelf (Figure 2). We note that\n75 a blocked Bering Sea could cause the Atlantic Meridional Overturning Circulation to exhibit\n76 increased hysteresis behavior during glacial periods (Hu et al., 2012) and impair seawater\n77 transfer via the Bering Strait between the Arctic and the Pacific. However, we find that the\n78 simulated Pacific sea ice extent south of the Bering Strait (Supplementary Fig. 7) and the\n79 complete Pacific Ocean (Fig. 4d–f) exhibit similar seasonal and perennial sea ice variability\n80 similar to NH sea ice variability (Fig. 3). This is possible because the simulated NH meridional\nPage 3 of 48\n81 heat is transported equatorward mainly through the atmosphere for high obliquity and high\n82 atmospheric CO 2 levels (Yun et al., 2023) and it undermines the influence of a closed Bering\n83 Strait on sea ice in the 3 Ma simulation.\n84 Additionally, the presence of Barents Shelf as a subaerial land mass would reduce sea ice\n85 formation during glacial-interglacial periods, influence surface albedo changes(Cao et al.,\n86 2019), and could warm the high-latitude(Cao et al., 2019). However, Yun et al. (2023) state that\n87 employing a volume-conserving ocean model with transient sea level and ocean volume would\n88 result in higher computational costs and perhaps unexplained biases and ambiguities that have\n89 not been thoroughly investigated in studies.\n90 Line #471-486: Further, the 3 Ma experiment uses a constant ocean bathymetry and land-sea\n91 boundary from the Last Glacial Maximum (∼21ka), resulting in ∼135 meters lower sea level\n92 than the present(Willeit et al., 2019). Consequently, we note that high-latitude has some\n93 unrealistic land-sea mask, such as a closed Bering Strait and subaerial Barents Shelf (Figure\n94 2). A closed Bering Sea would disrupt seawater transport between the Arctic and the Pacific\n95 via the Bering Strait, and it could lead to the emergence of stronger hysteresis behavior in the\n96 Atlantic Meridional Overturning Circulation during glacial periods(Hu et al., 2012). However,\n97 we find that the simulated Pacific sea ice extent south of the Bering Strait (Supplementary Fig.\n98 7) and the complete Pacific Ocean (Fig. 4d–f) show similar seasonal and perennial sea ice\n99 variability akin to NH sea ice variability (Fig. 3) possibly because the simulated NH meridional\n100 heat is transported equatorward mainly through the atmosphere for high obliquity and high\n101 atmospheric CO 2 levels(Yun et al., 2023). Additionally, the presence of a subaerial Barents\n102 Shelf would reduce sea ice formation during glacial-interglacial periods, influence surface\n103 albedo changes(Cao et al., 2019), and could warm the high-latitude(Cao et al., 2019). Then\n104 again, using a volume-conserving ocean model with time-evolving ocean volume and sea level\n105 would increase computational costs and may result in unidentified biases and ambiguities that\n106 studies have not fully explored(Yun et al., 2023).\n107\n108 Regarding the interpretation of periodicities of greenhouse gases (Supplementary Figure 2):\n109 here it must be taken into account that GHG is not only a forcing, but also subject to other\n110 changes in climate, and it indeed seems to be derived from output from a different simulation\nPage 4 of 48\n111 rather than from proxy-records. It is therefore (at in the opinion of this reviewer) not very\n112 surprising that one finds a change in the periodicity of GHG variability that follows the\n113 periodicity of glacial-interglacial cycles (41 kyr pre-MPT, 100 kyr post-MPT). Hence it would\n114 be useful to point out what we actually learn from such analyses. My apologies if this point has\n115 been missed by this reviewer.\n116 Response: We appreciate the concern. We have discussed the 3 Ma experiment design in detail\n117 in the Methods section and the potential influence on sea ice in the discussion section.\n118 The 3 Ma experiment is forced with transient prescribed greenhouse gases (GHG)(Willeit et\n119 al., 2019) forcing along with transient prescribed orbital (A. L. Berger, 1978) and ice-\n120 sheet(Willeit et al., 2019) forcings and a constant ocean bathymetry and land-sea boundary\n121 from the Last Glacial Maximum (∼21ka).\n122 Line #520-522: \"The 3 Ma experiment uses bathymetry during the Last Glacial Maximum\n123 (∼21 ka), and it is forced with transient prescribed orbital(A. L. Berger, 1978), ice-sheet(Willeit\n124 et al., 2019), and greenhouse gases (including CO 2 , CH 4 , and N 2 O concentration)\n125 (GHG)(Willeit et al., 2019) forcings.\"\n126 We have clarified that the prescribed GHG forcing involves CO 2 , CH 4 , and N 2 O concentrations.\n127 The GHG forcing includes estimates of CO 2 concentration from different datasets, as the 3 Ma\n128 experiment does not include an interactive carbon cycle. From 0.8 to 3 Ma, the estimated\n129 greenhouse gases consist of simulated CO 2 concentration from a previous CLIMBER-2\n130 transient simulation(Willeit et al., 2019), which contains an interactive carbon cycle.\n131 Meanwhile, from 0 to 0.8 Ma, the GHG estimate includes air-bubble measurement of CO 2\n132 concentration from the Dome C ice core(Jouzel et al., 2007), European Project for Ice Coring\n133 in Antarctica (EPICA). The GHG forcing also includes estimates of N 2 O and CH 4\n134 concentrations as a CESM forcing obtained by regressing the EPICA Dome C ice core N 2 O\n135 and CH 4 concentrations on the LR04 global Benthic Stack proxy record(Lisiecki & Raymo,\n136 2005) (δ18O of benthic foraminiferal calcite, which measures global ice volume and deep\n137 ocean temperature) post 0.8 Ma and extending the linear statistical model to 3 Ma using the\n138 LR04 δ18O stack record.\n139 Line #80-83: “The 3 Ma experiment is an atmosphere-ocean-sea ice coupled simulation\n140 performed by Yun et al. (2023) employing the Community Earth System Model version 1.2\nPage 5 of 48\n141 (CESM1.2), which is forced with transient prescribed orbital(A. Berger, 1978), GHG(Willeit\n142 et al., 2019) (including CO 2 , CH 4, and N 2 O concentration), and ice-sheet forcings(Willeit et al.,\n143 2019).”\n144 Line #523-531: \"The GHG forcing includes estimates of CO 2 concentration from different\n145 datasets. From 0.8 to 3 Ma, the estimated greenhouse gases consist of simulated CO 2\n146 concentration from a previous CLIMBER-2 transient simulation(Willeit et al., 2019), which\n147 contains an interactive carbon cycle. Meanwhile, from 0 to 0.8 Ma, the GHG estimate includes\n148 air-bubble measurement of CO 2 concentration from the Dome C ice core(Jouzel et al., 2007),\n149 European Project for Ice Coring in Antarctica (EPICA). The GHG forcing also includes\n150 estimates of N 2 O and CH 4 concentrations as a CESM forcing obtained by regressing the EPICA\n151 Dome C ice core N 2 O and CH 4 concentrations on the LR04 global Benthic Stack proxy\n152 record(Lisiecki & Raymo, 2005) (δ18O of benthic foraminiferal calcite, which measures\n153 global ice volume and deep ocean temperature) post 0.8 Ma and extending the linear statistical\n154 model to 3 Ma using the LR04 δ18O stack record.\n155 We also discuss that the strong association between sea ice and CO 2 concentrations in our\n156 results could be because the GHG forcing anomalies relative to preindustrial conditions are\n157 amplified by a factor of 1.5 to compensate for relatively weak climate sensitivity in the applied\n158 model.\n159 Line #459-462: “This strong association between sea ice and CO 2 concentration could be\n160 because the GHG forcing anomalies relative to preindustrial conditions are amplified by a\n161 factor of 1.5 to compensate for relatively weak climate sensitivity in the applied model.”\n162 Line #535-539: \"Further, in order to compensate for the impact of unresolved and highly\n163 uncertain dust forcing and to partially offset the relatively weak climate sensitivity of CESM1.2\n164 compared to multi-model ensemble-based estimation from the Coupled Model\n165 Intercomparison Project Phase 6 (CMIP6), the greenhouse gas forcing anomalies relative to\n166 preindustrial conditions have been amplified by a factor of 1.5.\"\n167 Meanwhile, in our study, we performed spectral analysis on the forcing data and found that the\n168 orbital forcing has a precession and obliquity signal throughout the Quaternary. In contrast, the\n169 prescribed CO 2 concentration (GHG forcing) and prescribed global ice volume in sea level\n170 equivalent (ice-sheet forcing) have a dominant ∼41 kyr periodicity during the pre–MPT period,\nPage 6 of 48\n171 ∼60–100 kyr (CO 2 /ice-sheet signal) periodicity during the MPT period and a longer ∼100 kyr\n172 periodicity during the post–MPT period. The significance of our analysis is that the presence\n173 of CO 2 /ice-sheet signal (∼60–100 kyr) during MPT and ∼100 kyr signal during post–MPT in\n174 the simulated sea ice/surface air temperature/sea surface temperature directly reflects the\n175 influence of GHG/ice-sheet forcings. Additionally, we agree with the reviewer that during the\n176 Quaternary period, the GHG is a forcing, which is a response and a byproduct of the feedback\n177 from the carbon cycle. Therefore, an interactive carbon cycle model would be more realistic.\n178 However, under the current model's capability, including a carbon cycle model may produce\n179 other model biases. Therefore, by prescribing GHG forcing, which involves CO 2 , CH 4 , and\n180 N 2 O concentrations, we could assess the impact of GHG forcing on sea ice without invoking a\n181 complex causality residing in the Earth's carbon cycle system.\n182 Line #94-96: \"The incoming solar radiation at 65°N consistently varies throughout the\n183 Quaternary period (Fig. 1d) with a predominant ∼23/19 kyr (precession) and ∼41 kyr\n184 (obliquity) signal (Supplementary Fig. 1).\"\n185 Line #98-104: \"The prescribed CO 2 concentration(Willeit et al., 2019) (GHG forcing) (Fig. 1e)\n186 (Supplementary Fig. 2a–c) and prescribed global ice volume in sea level equivalent(Willeit et\n187 al., 2019) (ice-sheet forcing) (Fig. 1f) (Supplementary Fig. 2j–l) show a dominant ∼41 kyr\n188 periodicity during the pre–MPT period, ∼60–100 kyr (CO 2 /ice-sheet signal) periodicity during\n189 the MPT period and a longer ∼100 kyr periodicity during the post–MPT period …..\"\n190 Line #147-154: \" Before understanding the sea ice variability signals, it is important to keep\n191 in mind that although the prescribed orbital forcing has a characteristic ∼41 kyr (obliquity) and\n192 ∼23/19 kyr (precession) signal throughout the Quaternary, it does not necessarily indicate a\n193 direct response to orbital forcing. This is because the prescribed CO 2 (GHG) (Supplementary\n194 Fig. 2a–c) and prescribed ice-sheet (Supplementary Fig. 2j–l) forcings also have a major ∼41\n195 kyr signal during the pre–MPT and a minor ∼23/19 signal across the Quaternary. Furthermore,\n196 the CO 2 /ice-sheet signal (∼60–100 kyr) during MPT and ∼100 kyr signal during post–MPT\n197 indicate a direct response to GHG/ice-sheet forcings.\n198\n199\n200 Overall: Avoid using \"~ \"for stating the latitudes considered by the authors. If the analyses are\nPage 7 of 48\n201 done for northward of 70N, it is precisely northward of that latitude as the authors chose this\n202 as a reference.\n203 Response: Thank you. As per suggestion, \"~ \"is avoided.\n204\n205 The authors could start the results section clarifying their domain of analyses (35N to 90N),\n206 their definitions of sea ice extent and other variables, and what is considered NH high latitude.\n207 This would decrease redundancy and the usage of the NH acronym and other explanations.\n208 Response: Thank you for the suggestion. As advised, we have clarified the domain and\n209 definitions earlier in the result section.\n210 Line #119-122: \"We investigate sea ice extent (ocean surface covered with more than 15% sea\n211 ice concentration) in the Northern Hemisphere (NH; 35°N–90°N) during the post–MPT (0–0.7\n212 Ma), MPT (0.7–1.2 Ma) and pre–MPT (1.2–2.6 Ma) periods, as our study domain. The high-\n213 latitude is defined as the region north of 70°N, and the mid-latitude is between 35°N and 70°N\n214 latitudes.\"\n215 Line #127-130: \"Here, we define sea ice concentration exceeding 80–90% are similar to\n216 perennial sea ice which refers to sea ice that survives during summer and persists for over a\n217 year. Meanwhile, less than 80–90% of sea ice concentration is similar to seasonal sea ice that\n218 forms in winter and melts in summer.\"\n219 The authors also need to determine what is variance and what is variability, sea ice area, extent,\n220 concentration, and use it consistently. There is also a general inconsistency in describing what\n221 drives certain processes. To say that the variability of a specific variable is more pronounced\n222 under precession-like periodicity or timescales is not the same as to affirm that the variable is\n223 precession-driven.\n224 Response: Thank you for the insight. We have extensively edited the manuscript, removed\n225 inconsistencies, and made it clearer and more scientific, using the Figures as supporting\n226 material.\n227\n228 At summary and discussion the authors do not show a specific link between the sentences in\nPage 8 of 48\n229 the first paragraph. Previous studies reconstruct or simulate sea ice extent/concentration (not\n230 variability – it varies over time) with limited timescale, and therefore this manuscript provides\n231 novel and important results? Past modeling studies also cover SAT and SST. How is the link to\n232 authors' findings on linkages between SAT, SST, and sea ice in the Arctic?\n233 Response: Thank you for the insight. We have extensively edited the manuscript, removed\n234 inconsistencies, and made it clearer and more scientific. Subsequently, in the discussion section,\n235 our revised manuscript notes that previous studies are mainly limited to the post-MPT period,\n236 while our study examines in detail the spatiotemporal variability and the mechanism of NH sea\n237 ice variability during the post–MPT, MPT, and pre–MPT periods. Our study also emphasizes\n238 that our analysis of the extremely long 3 Ma simulation data provides an excellent test bed to\n239 explore the effects of orbital, GHG, and ice-sheet forcings on sea ice in more detail and with\n240 more precision compared to previous studies limited to the post-MPT period.\n241 Line #416-422: \"Previous proxy(Cronin et al., 2017; Cronin et al., 2013; Lo et al., 2018; Polyak\n242 et al., 2013; Sha et al., 2016; Stein et al., 2017; Turney et al., 2015) and modeling studies(Lo\n243 et al., 2018; Mantsis et al., 2011; Tuenter et al., 2005; Wu et al., 2020, 2022; Yin & Berger,\n244 2010) have mostly recorded or simulated sea ice concentration limited to the post-MPT. A\n245 recent study (3 Ma experiment) has presented a broad outline of the Northern and Southern\n246 Hemisphere sea ice variability, which shows a ∼21 kyr (precession), ∼41 kyr (obliquity), and\n247 ∼60–120 kyr (CO 2 /ice-sheet) periodicity, suggesting the combined effect of GHG and ice-sheet\n248 forcing(Yun et al., 2023) on sea ice variability. Here, we analyze the 3 Ma experiment to\n249 investigate the spatiotemporal and latitudinal variability and the mechanism of NH sea ice\n250 variability during the post–MPT, MPT, and pre–MPT periods..\"\n251 How can the authors affirm that their study agrees with almost all geological SST records? The\n252 references indicated in the manuscript are surely neither descriptive of all datasets that include\n253 SST reconstructions during the Quaternary nor consider most boundary conditions employed\n254 in previous experiments.\n255 Response: Thank you. We actually meant that most of the geological sea surface temperature\n256 records show a 100 kyr cycle during the post-MPT period (McClymont et al., 2013), which\n257 agrees with our result. Meanwhile, we have extensively edited the manuscript and removed\n258 this redundant information to make our manuscript clearer and more scientific.\nPage 9 of 48\n259\n260 What is the 'it' referred to in line 431? Also, all of the statistical descriptions could be translated\n261 into the effect of the forcing onto sea ice in a concise manner.\n262 Response: Thank you for the suggestion. As suggested, we have moved the statistical\n263 descriptions in the results section and have clearly discussed it using Figure 8 in the revised\n264 manuscript.\n265 Line #353-358: \"We find that the perennial sea ice covering the high–latitude (Fig. 2) has a\n266 strong Pearson correlation of −0.79 and standardized regression coefficient of −0.76 with the\n267 obliquity index (lags obliquity index by ∼5 kyr) (green line in Fig. 8a). Meanwhile, the seasonal\n268 sea ice covering the mid-latitude (Fig. 2) has a strong Pearson correlation of 0.70 and\n269 standardized regression coefficient of 0.72 with precession index (lags precession index by\n270 ∼10 kyr) (green line in Fig. 8b). Therefore, we identify a strong coupling between simulated\n271 sea ice with orbital forcings throughout the Quaternary.\"\n272 Sea ice distribution? coverage? The authors need to define the sea ice measurement they are\n273 exploring in the manuscript and keep it consistent. Also, there can be more description of what\n274 it entails to do seasonal studies (computational costs, resolution, boundary conditions,\n275 calendar corrections) to provide a clear but broad view for future studies.\n276 Response: Thank you for the suggestion. As suggested, we have removed the inconsistencies\n277 regarding sea ice and have clarified the definitions of sea ice extent and seasonal sea ice earlier\n278 in the result section. We have also replaced quasi-permanent sea ice with perennial sea ice.\n279 Additionally, we have not discussed possible future studies involving seasonal sea ice to\n280 maintain conciseness.\n281 Line #119-121: \"We investigate sea ice extent (ocean surface covered with more than 15% sea\n282 ice concentration) in the Northern Hemisphere (NH; 35°N–90°N) ……..\"\n283 Line #127-130: \"Here, we define sea ice concentration exceeding 80–90% are similar to\n284 perennial sea ice which refers to sea ice that survives during summer and persists for over a\n285 year. Meanwhile, less than 80–90% of sea ice concentration is similar to seasonal sea ice that\n286 forms in winter and melts in summer.\"\n287\nPage 10 of 48\n288 In the methods section the description of data could clearly be improved, it appears very\n289 selective at the moment. What about the other data, than model data, that was mentioned\n290 during the text?\n291 Response: We appreciate the suggestion. We have only analyzed the 3 Ma simulation data and\n292 have elaborated the 3 Ma climate model experiment design and forcings in the Methods section.\n293 Line #510–554: “The 3 Ma experiment(Yun et al., 2023) is a transient paleoclimate model\n294 simulation for the past 3 million years employing the coupled CESM version 1.2 (CESM 1.2).\n295 The experiment comprises a fully coupled atmosphere (Community Atmosphere Model version\n296 4; CAM4), ocean (Parallel Ocean Program version 2; POP2), sea ice (Los Alamos Sea Ice\n297 Model; CICE4), and land (Community Land Model version 4; CLM4.0; prescribed vegetation)\n298 components. The atmospheric model has about 3.75° × 3.75° latitude-longitude horizontal\n299 resolution (T31) with 26 vertical levels. The ocean and sea ice model has a 100 × 116\n300 horizontal grid (gx3v5), and the ocean model has 25 vertical levels. The ocean and sea ice\n301 output are gridded to a 3.75° × 3.75° latitude-longitude horizontal resolution identical to\n302 atmospheric outputs.\n303 The 3 Ma experiment uses bathymetry during the Last Glacial Maximum (∼21 ka),\n304 and it is forced with transient prescribed orbital(A. L. Berger, 1978), ice-sheet(Willeit et al.,\n305 2019), and greenhouse gases (including CO 2 , CH 4 , and N 2 O concentration) (GHG)(Willeit et\n306 al., 2019) forcings.\n307 The GHG forcing includes estimates of CO 2 concentration from different datasets.\n308 From 0.8 to 3 Ma, the estimated greenhouse gases consist of simulated CO 2 concentration from\n309 a previous CLIMBER-2 transient simulation(Willeit et al., 2019), which contains an interactive\n310 carbon cycle. Meanwhile, from 0 to 0.8 Ma, the GHG estimate includes air-bubble\n311 measurement of CO 2 concentration from the Dome C ice core(Jouzel et al., 2007), European\n312 Project for Ice Coring in Antarctica (EPICA). The GHG forcing also includes estimates of N 2 O\n313 and CH 4 concentration as a CESM forcing obtained by regressing the EPICA Dome C ice core\n314 N 2 O and CH 4 concentration on the LR04 global Benthic Stack proxy record(Lisiecki & Raymo,\n315 2005) (δ18O of benthic foraminiferal calcite, which measures global ice volume and deep\n316 ocean temperature) post 0.8 Ma and extending the linear statistical model to 3 Ma using the\n317 LR04 δ18O stack record. The ice-sheet topography changes are taken from the simulated\nPage 11 of 48\n318 CLIMBER-2 orography data from the previous transient simulation(Willeit et al., 2019).\n319\n320 Further, in order to compensate for the impact of unresolved and highly uncertain dust\n321 forcing and to partially offset the relatively weak climate sensitivity of CESM1.2 compared to\n322 multi-model ensemble-based estimation from the Coupled Model\n323 Intercomparison Project Phase 6 (CMIP6), the greenhouse gas forcing anomalies relative to\n324 preindustrial conditions have been amplified by a factor of 1.5. The 3 Ma simulation also uses an\n325 acceleration factor of 5 to the external forcings (orbital, ice-sheet, and GHG forcings), which\n326 reduces the 3 Ma orbital years forcing to 0.6 Ma model years. The orbital, ice-sheet, and GHG\n327 forcings are updated every 100 orbital years, which means 20 model years.\n328 The 3 Ma simulation data has been compared against large-scale climate\n329 reconstructions and numerous paleo-proxy data. The 3 Ma simulation has been shown to\n330 reproduce the magnitude and timing of glacial–interglacial variability(Timmermann et al.,\n331 2022; Yun et al., 2023) consistent with paleoclimate data constraints and estimates from the\n332 intermediate complexity Earth system model. The 3 Ma simulation agrees well(Timmermann\n333 et al., 2022; Yun et al., 2023) with key paleoclimate proxy records such as the magnitude and\n334 timing of global mean surface temperature variations, annual mean sea surface temperatures\n335 over the tropics and elsewhere (western south and north Pacific, southern Indian and north\n336 Atlantic), Greenland and Antarctic temperatures, East Asian summer monsoon and eastern\n337 African hydroclimate.\n338 The 3 Ma simulation data are available in the annual average millennial mean. Yun et\n339 al. (2023) provide the complete details on the 3 Ma model set-up and experimental design.”\n340\n341 Comments on Figures:\n342 Figure 1: The discussion fostered by the authors based on figure 1 would be more clear if the\n343 time-series used different colors for the parameters examined instead of a solid black line for\n344 some and a dashed red line for others. In this way, key results such as the surface air\n345 temperature are hardly visible due to clutter from the overlying time-series. The authors can\n346 also describe in the figure caption the direction of the time axis instead of using the arrows as\nPage 12 of 48\n347 it increases the clutter in the image and decreases understanding of the processes. The slope\n348 of the curves between the post-MPT, MPT, and pre-MPT periods could be described in a table\n349 or just mentioned in the text.\n350 Furthermore: you do not provide any reference to any of the data sources. Providing them\n351 would enable readers to better understand your plots without going back and forth between\n352 manuscript text and figure caption.\n353 Response: We are grateful for the suggestion. Firstly, as suggested, the slope is mentioned in\n354 the text.\n355 Line #96–98: \"The Quaternary climate shows a gradual cooling towards the post–MPT but the\n356 trend of the incoming solar radiation did not appreciably change during the past 2.6 Ma (0.0001\n357 Wm−2/kyr) (Fig. 1d).\"\n358 Line #390–398: \"We find that the prescribed CO 2 concentration (GHG) (Supplementary Fig.\n359 2a–c) and simulated global mean surface temperature (Supplementary Fig. 6) have similar\n360 periodicity signals throughout the Quaternary. During the pre–MPT, the prescribed CO 2\n361 concentration (Fig. 1e) and simulated global mean surface temperature (Fig. 1h) have a gradual\n362 downward trend, intensifying across the MPT by about 2 and 3 times, respectively. On the\n363 contrary, the simulated sea ice extent accompanies an upward trend (significant at a 95%\n364 confidence level) during the pre–MPT (1.5 thousand km2/kyr), which intensifies by a factor\n365 of about 2.2 (3.4 thousand km2/kyr) during the MPT (Fig. 1g) indicating the influence of GHG\n366 and/or ice-sheet forcing on sea ice.\"\n367 Secondly, Figure 1 shows time series for different parameters in colors. We have provided\n368 references to data sources in the figure caption and main manuscript.\n369 Line #80–83: “The 3 Ma experiment is an atmosphere-ocean-sea ice coupled simulation\n370 performed by Yun et al. (2023) employing the Community Earth System Model version 1.2\nPage 13 of 48\n371 (CESM1.2), which is forced with transient prescribed orbital(A. Berger, 1978), GHG(Willeit\n372 et al., 2019) (including CO 2 , CH 4, and N 2 O concentration), and ice-sheet forcings(Willeit et al.,\n373 2019).”\n374 Line #93–94: “The prescribed orbital forcing(A. L. Berger, 1978) (eccentricity, precession, and\n375 obliquity) (Fig. 1a–c) influences the intensity of summer (June–July–August; JJA) ……..”\n376 Line #98–101: “The prescribed CO 2 concentration(Willeit et al., 2019) (GHG forcing) (Fig.\n377 1e) (Supplementary Fig. 2a–c) and prescribed global ice volume in sea level equivalent(Willeit\n378 et al., 2019) (ice-sheet forcing) (Fig. 1f) (Supplementary Fig. 2j–l) show a dominant ∼41 kyr\n379 periodicity ……………………………..”\n380 Line #105–113: \"\n381 Fig. 1: The 3 Ma experiment forcings and simulations using the CESM1.2 model. The time\n382 series shows the prescribed orbital forcing: (a) eccentricity index(A. L. Berger, 1978),\n383 (b) precession index(A. L. Berger, 1978) (𝑒sinω, where 𝑒 is the eccentricity and ω is the longitude\n384 of the perihelion), (c) obliquity index(A. L. Berger, 1978) (units are in degrees); (d) summer\n385 (JJA) incoming solar radiation at 65°N (units are in watts per square meter), (e) prescribed CO 2\n386 concentration(Willeit et al., 2019) (units are in ppmv), (f) prescribed global ice volume in sea\n387 level equivalent(Willeit et al., 2019) (units are in meters), (g)……..\"\nPage 14 of 48\n388\n389 Figure 1j-l should be made a separate Figure (2a-c) to improve understanding and discussion\n390 of the separate representations (one focuses on variability, the other on the mean-state).\n391 Regarding Fig. 1 j-l: a colorbar of sea ice concentration that starts with -20% is very confusing\n392 - unless negative values should illustrate model land sea mask or similar (see my comment\n393 below regarding land sea masks). The meaning of values below zero should then be clearly\n394 documented, and ideally using a color different than white (which currently just means \"no sea\n395 ice \").\n396 Next, why do you apply different definitions of the permanent sea ice boundary pre- and post-\n397 MPT? Doing so makes the plots difficult to compare. Is the fact, that the isolines look very\n398 similar between all three time periods, a sideeffect of different definitions?\n399 Furthermore, there is one problem with the visualization: I hypothesize the permanent ice free\n400 area in the Barents Sea stems from a model land sea mask with subaerial Barents Shelf that is\n401 not shown on the plot? In particular, I assume the model considers a subaerial Barents Shelf,\n402 i.e. the white color shown there does not imply that there are ice free ocean regions, am I right?\n403 I think it would be very important to illustrate the land sea mask as it is in the model rather\n404 than overlaying a (high-resolution) modern land sea mask - otherwise readers will be confused\n405 regarding sea ice dynamics. In this respect, more information on the paleogeography used in\n406 the model, and potential differences between the studied time periods, should be presented\n407 somewhere in the manuscript. This is particularly relevant for continental ice sheet forcing as\n408 you argue that this forcing would control high latitude sea ice in the Northern Hemisphere.\n409 Response: Thank you for the suggestion. Figure 1j-l is drawn separately as Figure 2a-c in the\n410 updated manuscript. As suggested, the contour value/scale is from 0 to 100% sea ice\n411 concentration. Figure 2 in the updated manuscript also shows subaerial land regions (especially\n412 subaerial Barents Shelf and closed Bering Strait) as the experiment uses a constant Last Glacial\n413 Maximum land-sea boundary. The discussion section also discussed the potential ramifications\n414 of using a constant Last Glacial Maximum (∼21ka) land-sea boundary.\n415 Line #133-140:\"\nPage 15 of 48\n416 Fig. 2: The climatology of sea ice during the post–MPT (0–0.7 Ma), MPT (0.7–1.2 Ma), and\n417 pre–MPT (1.2–2.6 Ma) periods. The shaded blue color shows the sea ice coverage over modern\n418 land-sea boundaries. The units are in percentage sea ice concentration. The grey dots show\n419 subaerial land regions (especially subaerial Barents Shelf and closed Bering Strait) as the\n420 experiment uses a constant Last Glacial Maximum land-sea boundary. The black and purple\n421 lines demarcate the sea ice extent and perennial sea ice boundary, respectively.\"\n422 Line #471-486: \"Further, the 3 Ma experiment uses a constant ocean bathymetry and land-sea\n423 boundary from the Last Glacial Maximum (∼21ka), resulting in ∼135 meters lower sea level\n424 than the present(Willeit et al., 2019). Consequently, we note that high-latitude has some\n425 unrealistic land-sea distributions, such as a closed Bering Strait and subaerial Barents Shelf\n426 (Figure 2). A closed Bering Sea would disrupt seawater transport between the Arctic and the\n427 Pacific via the Bering Strait, and it could lead to the emergence of stronger hysteresis behavior\n428 in the Atlantic Meridional Overturning Circulation during glacial periods(Hu et al., 2012).\n429 However, we find that the simulated Pacific sea ice extent south of the Bering Strait\n430 (Supplementary Fig. 7) and the complete Pacific Ocean (Fig. 4d–f) show similar seasonal and\n431 perennial sea ice variability akin to NH sea ice variability (Fig. 3) possibly because the\n432 simulated NH meridional heat is transported equatorward mainly through the atmosphere for\n433 high obliquity and high atmospheric CO 2 levels(Yun et al., 2023). Additionally, the presence\n434 of a subaerial Barents Shelf would reduce sea ice formation during glacial-interglacial periods,\n435 influence surface albedo changes(Cao et al., 2019), and could warm the high-latitude(Cao et\n436 al., 2019). Then again, using a volume-conserving ocean model with time-evolving ocean\n437 volume and sea level would increase computational costs and may result in unidentified biases\n438 and ambiguities that studies have not fully explored(Yun et al., 2023).\"\n439 Next, we applied a slight difference in defining perennial sea ice boundaries during pre-MPT\n440 (above 80% sea ice concentration) and MPT/post-MPT (above 90% sea ice concentration). We\nPage 16 of 48\n441 reason that during the pre-MPT period, the sea ice is less expanded than in the post-MPT period\n442 because of relatively higher pre-MPT greenhouse gasses concentration and global mean surface\n443 temperature. Therefore, we used slightly different definitions and argued that, in general,\n444 perennial sea ice covered the high-latitude (north of 70°N latitude) regions throughout the\n445 Quaternary. Furthermore, there is not much difference between the terms \"above 80%\" and\n446 \"above 90%\" in the definition, as shown in Supplementary Figure 4 in the revised manuscript.\n447 Supplementary Line #52-59: \"\n448 Supplementary Fig. 4: The time series of zonal mean sea ice concentration at 72°N and 76°N\n449 for the past 2.6 Ma. The purple dashed line shows the average sea ice concentration at 72°N\n450 during the post–MPT, MPT, and pre–MPT periods. The units are in percentage sea ice\n451 concentration. We notice that sea ice concentration in the high-latitude region (north of 70°N\n452 latitude), such as those at 72°N and 76°N, is covered with more than 80% sea ice concentration\n453 during the pre-MPT and more than 90% sea ice concentration during the MPT and post-MPT\n454 periods.\"\n455\n456 Figure 2: Define the units of the colorbar.\n457 Response: Thank you for the suggestion. The units show standardized variance and are defined\n458 in the caption for Figure 3.\n459 Line #157-159: \"Fig. 3: The periodicity-variance spectrum of NH sea ice and perennial sea ice.\n460 The standardized periodicity-variance spectrum of sea ice is shown in colors (units show\nPage 17 of 48\n461 standardized variance) during…..\"\n462 Figure 3: The domain analysed in this manuscript is from 35N to 90N. The repetition is\n463 redundant.\n464 Response: Thank you for the suggestion. As suggested, we have removed the redundancy and\n465 have clarified the domain analyzed in our manuscript earlier in the result section.\n466 Line #119-121: \"We investigate sea ice extent (ocean surface covered with more than 15% sea\n467 ice concentration) in the Northern Hemisphere (NH; 35°N–90°N) ……..\"\n468\n469 Figure 4: Split the figures (similar to the suggestion above). Consider showing more latitude\n470 circles in Fig. 4a-c to make it easier for the reader to identify the progressive shift in sea ice\n471 extent.\n472 Response: Thank you for the suggestion. Figure 4 is drawn separately and clearly as Figure 5\n473 in the revised manuscript. As suggested, we have shown more latitude circles in Figure 5a-c\n474 and presented Figure 5d-i more clearly.\n475 Figure 5a-c shows the major glacial-interglacial sea ice variability in the Labrador Sea,\n476 Irminger-Iceland basin, and central Pacific regions. Meanwhile, Figure 5d-i shows the\n477 periodicity-variance spectrum of sea ice in those major sea ice variability regions, as shown in\n478 Figure 5a-c. Therefore, we have shown them together in Figure 5 to make it easier for the reader\n479 to identify the progressive shift in sea ice extent.\n480 Figure 6: Since these two graphs represent the same calculation but for different indices, the\n481 authors can just refer to the indices in the title and keep only one legend relative to the region\n482 (using solid lines only instead of both solid and dashed lines).\n483 Response: Thank you. We have followed the suggestions, and Figure 6 is clearly drawn as\n484 Figure 8 in the revised manuscript.\nPage 18 of 48\n485 Line #322-322:\n486 \"Fig. 8: Cross-correlation of summer incoming solar radiation at 65°N, surface air\n487 temperature, and sea ice with obliquity and precession indices (orbital parameters) from 0–\n488 2.6 Ma. (a) ….\"\n489\n490 Figure 7: This could be one figure with two y-axes. The negative and positive correlation\n491 coefficients are misleading.\n492 Response: Thank you. We have followed the suggestions, and Figure 7 is clearly drawn as\n493 Figure 9g in the revised manuscript.\n494 Line #378-380: \"\n495 Fig. 9: Relationship between CO 2 (GHG) forcing and sea ice during the Quaternary…...\"\nPage 19 of 48\n496\n497 Specific comments:\n498 line 31-34: to me the sentence seems difficult to understand; maybe: \"By altering the latitudinal\n499 distribution of solar insolation, orbital modulation of NH surface air temperature and sea ice\n500 extent (ocean surface covered with more than 15% sea ice concentration) is dominated by the\n501 precession cycle at mid latitudes and by the obliquity cycle at high latitudes. \"\n502 Response: Thank you. As suggested, we have added the mechanism and rephrased it clearly.\n503 Line #30-33: \"The orbital forcing modulates incoming net surface shortwave radiation,\n504 influencing surface air temperature and sea ice. For instance, the mid-latitude region, which is\n505 covered with seasonal sea ice, is primarily driven by the Earth's axis wobble (precession).\n506 In contrast, the high-latitude (>70°N), which is covered with perennial sea ice, is mainly driven\n507 by the tilt in the Earth's axis (obliquity).\"\n508\n509 line 35-38: the sentence is quite complicated due to multiple use of brackets that are used in\n510 different context; furthermore, the statement regarding intensifying GMST and sea ice extent is\n511 unclear at least to me; what does it physically mean if GMST intensifies by a factor of ~3? It\n512 does not get three times warmer, right? And also: \"followed by no appreciable change during\n513 the post–MPT (0–0.7 Ma) period. \"This reads a bit as if there were no Pleistocene glacial-\n514 interglacial climate changes, right? I assume that this is not the message that authors would\n515 like to convey here.\n516 Response: Thank you. As suggested, we have added the mechanism, removed redundant\n517 information, and rephrased it clearly.\n518 Line #34-35: \"Further, the GHG forcing changes downwelling longwave radiation, which\n519 influences global mean surface temperature and explains the spatiotemporal sea ice variability.\"\n520\n521 line 50-51: Is this really the definition of Arctic amplification? Arctic amplification does not\n522 imply that sea ice must be the driver. In fact, other drivers may be more important (e.g. Pithan\n523 and Mauritsen, 2014).\nPage 20 of 48\n524 Response: Thank you for the insight. As suggested, the sentence is reworded, and other\n525 processes such as temperature feedback, changes in water vapor, cloud cover, and atmospheric\n526 and ocean heat transport with citations have been added.\n527 Line #39-41: \"The sea ice amplifies climate change in high–latitude region, through its high\n528 albedo and minimal thermal conductivity, together with other processes (temperature feedbacks,\n529 changes in water vapor, cloud cover, and atmospheric and ocean heat transport)(Pithan &\n530 Mauritsen, 2014; Serreze & Barry, 2011), alter the radiative budget in the Arctic.\n531\n532 line 51-53: please clarify the link between sea ice and low latitude energy budget\n533 Response: Thank you. The sentence is reworded, and other processes influencing the radiative\n534 budget in the Arctic have been added, such as temperature feedback, changes in water vapor,\n535 cloud cover, and atmospheric and ocean heat transport with citations.\n536 Line #39-41: \"The sea ice amplifies climate change in high–latitude region, through its high\n537 albedo and minimal thermal conductivity, together with other processes (temperature feedbacks,\n538 changes in water vapor, cloud cover, and atmospheric and ocean heat transport)(Pithan &\n539 Mauritsen, 2014; Serreze & Barry, 2011), alter the radiative budget in the Arctic.\n540\n541 line 53: Consider changing from \"sea ice melting stratifies the near-surface...\" to \"sea ice\n542 melting increases stratification in the near-surface…\"\n543 line 55-56: This sentence needs a better connection with the previous one.\n544 Response: Thank you. As suggested, the sentence is rewritten clearly, and concisely with a\n545 better connection\n546 Line #42-45: \"Moreover, sea ice dynamics influence the location of the mid-latitude jet\n547 stream(Kidston et al., 2011), oceanic salinity structure(Abernathey et al., 2016), surface heat,\n548 and freshwater fluxes anomalies, which cause buoyancy anomalies affecting global\n549 intermediate/deep water circulation(Goosse et al., 2013) and North Atlantic Deep Water\n550 formation(Bretones et al., 2022) and, therefore, the Atlantic Meridional Overturning\n551 Circulation(Li et al., 2021; Sévellec et al., 2017).\"\nPage 21 of 48\n552\n553 line 59: Increased atmospheric CO2 concentration\n554 Response: Thank you. As suggested, corrections have been made.\n555 Line #45-47: “ …. coinciding with increased greenhouse gas (GHG) concentration and ……”\n556\n557 line 60: Less is understood in respect to what? Is it \"less\" or \"little\n558 Response: Thank you. As suggested, corrections have been made.\n559 Line #47: \"…. However, little is understood about the response ……\"\n560\n561 line 60: \"internal or external climate forcings \"- it sounds quite strange to speak about an\n562 internal climate forcing; normally, forcings are considered to be external factors impacting the\n563 climate system. Therefore, only external climate forcings make sense from a terminological\n564 point of view, and carbon dioxide must be one of these (unless predicted by the model, but then\n565 it is not a forcing). Please be clear in your terminology. Do you mean to write here \"internal\n566 processes/variability or external climate forcings \"?\n567 Response: Thank you. As suggested, corrections have been made.\n568 Line #47-48: \"…. response of NH sea ice to internal processes and external climate forcings.\"\n569\n570 line 61: Consider changing \"thus\" to \"additionally\" or to elaborate further on the importance\n571 of the Quaternary period. What is your motivation to focus your work on the Quaternary?\n572 You could have chosen another period to derive transient simulations. Consider also changing from\n573 \"era of modern humans\" to \"current climate\" or any other expression that is more related to\n574 climate science (and more precise - I am not sure which time period the \"ear of the modern\n575 humans \"actually refers to).\n576 Response: Thank you. As suggested, \"thus\" is replaced by \"additionally\" and the sentence is\n577 reworded concisely.\nPage 22 of 48\n578 Line #48-50: \"Additionally, to comprehend the historical………\"\n579\n580 line 66-68: The periodicity of the orbital forcing should be detailed in an additional sentence.\n581 Here it is put under parenthesis, but it is not clear by reading it the first time that the authors\n582 refer to periodicity.\n583 Response: Thank you. As suggested, the sentence is rewritten clearly and concisely.\n584 Additionally, the periodicity of orbital forcing is clearly mentioned in the results section.\n585 Line #52-56: \"The Earth's orbital parameters (orbital eccentricity (∼100 thousand years (kyr)),\n586 axial tilt or obliquity (∼41 kyr), and axis wobble or longitude of perihelion, known as\n587 precession (∼23/19 kyr)) significantly vary the latitudinal and seasonal distribution of the\n588 amount of incoming solar radiation(Masson-Delmotte et al., 2013) which is the dominant\n589 external forcing ……\"\n590 Line #94-96: \"The incoming solar radiation at 65°N consistently varies throughout the\n591 Quaternary period (Fig. 1d) with a predominant ∼23/19 kyr (precession) and ∼41 kyr\n592 (obliquity) signal (Supplementary Fig. 1).\"\n593\n594 line 68: Is the GHG forcing only CO2 concentration? If yes, you are analysing sensitivity to\n595 CO2 concentration, not GHG forcing overall.\n596 Response: Thank you for the insight. In the 3 Ma experiment, GHG includes CO 2 , CH 4 , and\n597 N 2 O concentration. We have mentioned it clearly in the introduction and methods sections. In\n598 the methods section, we have also elaborated on how the 3 Ma experiment estimates the CO 2 ,\n599 CH 4 , and N 2 O concentrations from different datasets.\n600 Line #80-83: “………….which is forced with transient prescribed orbital(A. Berger, 1978),\n601 GHG(Willeit et al., 2019) (including CO 2 , CH 4, and N 2 O concentration),……………….”\n602 Line #520-532: \" The 3 Ma experiment uses bathymetry during the Last Glacial Maximum\n603 (∼21 ka), and it is forced with transient prescribed orbital(A. L. Berger, 1978), ice-sheet(Willeit\n604 et al., 2019), and greenhouse gases (including CO 2 , CH 4 , and N 2 O concentration)\n605 (GHG)(Willeit et al., 2019) forcings.\nPage 23 of 48\n606 The GHG forcing includes estimates of CO 2 concentration from different datasets. From 0.8 to\n607 3 Ma, the estimated greenhouse gases consist of simulated CO 2 concentration from a previous\n608 CLIMBER-2 transient simulation(Willeit et al., 2019), which contains an interactive carbon\n609 cycle. Meanwhile, from 0 to 0.8 Ma, the GHG estimate includes air-bubble measurement of\n610 CO 2 concentration from the Dome C ice core(Jouzel et al., 2007), European Project for Ice\n611 Coring in Antarctica (EPICA). The GHG forcing also includes estimates of N 2 O and CH 4\n612 concentration as a CESM forcing obtained by regressing the EPICA Dome C ice core N 2 O and\n613 CH 4 concentration on the LR04 global Benthic Stack proxy record ..\"\n614 line 69: What in \"continental ice sheets\"? and \"played significant roles\" in what?\n615 Response: Thank you. The NH ice-sheet topography influences NH mid-latitude transient and\n616 stationary eddy. As suggested, the sentence is rewritten clearly and concisely.\n617 Line #56-59: \"…………and NH ice-sheet, which influences NH mid-latitude transient and\n618 stationary eddy, together with their associated internal climate feedback mechanisms amplify\n619 changes initiated by the orbital…...\"\n620\n621 line 70: Perhaps consider rewriting the sentence. Suggestion: \"During this period, NH\n622 continental ice sheets experienced cyclic waxing and waning in response to variations in global\n623 mean surface temperature.\"\n624 Thank you. As suggested, the sentence is rewritten clearly and concisely.\n625 Line #59-61: \"The NH ice-sheet experienced cyclic waxing and waning indicating glacial-\n626 interglacial cycles in response to variations in GHG concentration……\"\n627\n628 line 71-72: This is redundant as the authors already mentioned the periodicity shortly before\n629 the sentence.\n630 Thank you. As suggested, the sentence is rewritten clearly and concisely.\n631\n632 line 73 -74: \"with no noticeable variation in the key orbital parameters \"; this is a\nPage 24 of 48\n633 misformulation and potentially a wrong citation. The statement in the publication cited by you\n634 (Clark et al., 2006) refers in fact to a statement made in earlier work and it would likely be\n635 more appropriate to cite that work here. Furthermore, orbital parameters of course varied, that\n636 is the whole idea behind Milankovitch cycles. Clark et al. (2006) state, in citation of previous\n637 work, that \"variations in the Earth's orbital parameters did not change over the last 2000 ka \".\n638 There is an important difference between this formulation and your rephrasing - the cited\n639 literature states that, while orbital parameters changed with time, *the nature of their\n640 variation* did not change.\n641 Thank you. As suggested, the sentence is rewritten clearly and concisely.\n642 Line #52-56: “The Quaternary period is characterized by glacial-interglacial cycles. The\n643 Earth's orbital parameters (orbital eccentricity (∼100 thousand years (kyr)), axial tilt or\n644 obliquity (∼41 kyr), and axis wobble or longitude of perihelion, known as precession (∼23/19\n645 kyr)) significantly vary the latitudinal and seasonal distribution of the amount of incoming solar\n646 radiation(Masson-Delmotte et al., 2013) which is the dominant external forcing for the glacial-\n647 interglacial cycles (Milankovitch hypothesis)(Hays et al., 1976; Imbrie et al., 1992;\n648 Milankovitch, 1941).”\n649\n650 line 74: What is the difference between astronomical and orbital parameters? The authors\n651 should keep consistency. Also, either refer to \"incoming solar radiation\" or simply \"insolation.\"\n652 Thank you. As suggested, consistency is maintained throughout the revised manuscript.\n653\n654 line 76-80: There is no previous mention of what MPT refers to before it is mentioned. It is\n655 described only after mentioning it for the first time. Perhaps a good paragraph introducing the\n656 MPT and why to study it, or on studies on its importance, would give a better view of why the\n657 authors chose to analyse this period.\n658 Response: Thank you for the suggestion. As suggested, we have defined the Mid-Pleistocene\n659 Transition (MPT) period (0.7–1.2 Ma). and the shift in the glacial cycles from 41 kyr to 100\n660 kyr across the MPT in the abstract.\nPage 25 of 48\n661 Line #24-26: \"Given the recent changes in the Northern Hemisphere (NH) sea ice, it is essential\n662 to understand sea ice variability over the past 2.6 million years (Ma), when the Earth's glacial\n663 cycles transitioned from ∼41 kyr to ∼100 kyr periodicity, following the Mid-Pleistocene\n664 Transition (MPT) period (0.7–1.2 Ma).\"\n665\n666 line 82-92: Very confusing paragraph, including concepts like periodicity, variation, cycles, all\n667 in terms of orbital parameters, not in terms of processes varying at x, y, z timescales. Some\n668 grammar issues also reduce clarity.\n669 line 85-87: \"Past studies have also shown [...] the past 800 ka 31 \"You speak of studies in\n670 plural but cite only one.\n671 line 87-90: The second part of the sentence could be reformulated for clarity: \"Additionally,\n672 various proxy records and modeling studies that are confined to the post–MPT periods have\n673 pointed to the significance of precession and obliquity cycles in NH sea ice variability, wherein\n674 some studies emphasized the role of precession, while others highlighted the role of obliquity\n675 \".\n676 Response: Thank you. As suggested, the paragraph is rewritten clearly and concisely.\n677 Line #67-74: \"Paleoclimate modeling study have shown that incoming solar radiation (orbital\n678 forcing) dominated the interglacial variations in NH sea ice during the past 0.8 Ma(Yin &\n679 Berger, 2012). Additionally, various proxy records and modeling studies that are primarily\n680 confined to the post–MPT period found ∼23/19 kyr (precession) and ∼41 kyr (obliquity)\n681 signals in NH sea ice variability, wherein some emphasized the role of precession(Cronin et al.,\n682 2013; Lo et al., 2018; Wu et al., 2020), while others did the role of obliquity(Mantsis et al.,\n683 2011; Tuenter et al., 2005; Turney et al., 2015; Yin & Berger, 2012). Recently, a transient\n684 paleoclimate simulation for the past 3 Ma presented a general overview of the Northern and\n685 Southern Hemisphere sea"}, "author_response": "Response to Reviewer #1 and #2 comments\n1\n3 comments and suggestions helped revise and improve our paper and the significance of our\n4 research. We studied the comments carefully and revised the manuscript as suggested. We have\n5 also improved our writing and grammar, making the manuscript clearer and scientific. We hope\n6 it meets your expectations with approval.\n8 (a) Original reviewer comments are in italics,\n9 (b) Our response to the comments/suggestions are highlighted in blue color,\n10 (c) References from the revised manuscript are in red.\n13 is part of the Nature Communications initiative to facilitate training in peer review and to\n14 provide appropriate recognition for Early Career Researchers who co-review manuscripts.\n17 that all the code and data are available upon request. Meanwhile, some of the 3 Ma experiment\n18 data is available on the OPeNDAP ICCP climate data server https://climatedata.ibs.re.kr (3Ma\n19 transient simulation with CESM1.2; Dataset doi: 10.22741/iccp.20230001).\n22 https://climatedata.ibs.re.kr (3Ma transient simulation with CESM1.2; Dataset\n23 doi: 10.22741/iccp.20230001).\"\n24\n27 cases where comments below do not explicitly raise that point. A rough scan over some\nPage 1 of 48\n28 references suggests that statements attributed by authors to previous publications are not\n29 always accurate.\n31 clarity. We also made sure that references are properly cited.\n33 tense. Is there a good reason not to use present tense? This statement applies also to other parts\n34 of the manuscript.\n36 manuscript and have formulated the sentences in the present tense.\n38 reader can find some hints here and there, but a comprehensive description that is easy to\n39 locate is still missing. At the moment it is hidden in lines 240ff as an addendum to the results\n40 and via some remarks in the summary and discussion. This is probably not appropriate towards\n41 a proper understanding of the model simulation and of the model's capabilities and limitations.\n43 employed as a boundary condition. This reviewer understands that throughout the simulation\n44 the land sea mask has been as for the LGM. This decision should certainly be\n45 motivated/justified in the manuscript and the impacts on interpretation of modelling results,\n46 that will certainly be substantial, must be critically discussed. Keeping the Bering Strait closed\n47 (see line 243) and the Barents Sea shelf subaerial (Figure 1) throughout the whole simulation\n48 will at times lead to sea ice being influenced by the various forcings and dynamics in a different\n49 way than it has been the case in reality, where land sea masks continuously changed. No model\n50 setup will be perfect, but the feeling is that more discussion should explore potential\n51 ramnifications.\n53 bathymetry and land-sea boundary from the Last Glacial Maximum (∼21ka), which causes\n54 some unrealistic land-sea boundaries and has been discussed in our revised manuscript.\n56 land regions in grey dots (especially subaerial Barents Shelf and closed Bering Strait) over\nPage 2 of 48\n57 modern land-sea boundaries.\n60 pre–MPT (1.2–2.6 Ma) periods. The shaded blue color shows the sea ice coverage over modern\n61 land-sea boundaries. The units are in percentage sea ice concentration. The grey dots show\n62 subaerial land regions (especially subaerial Barents Shelf and closed Bering Strait) as the\n63 experiment uses a constant Last Glacial Maximum land-sea boundary. The black and purple\n64 lines demarcate the sea ice extent and perennial sea ice boundary, respectively.\"\n66 discussed the potential ramifications of using a constant Last Glacial Maximum (∼21ka) land-\n67 sea boundary in the discussion section.\n69 of prescribed ice-sheet melt freshwater forcing, a constant Last Glacial Maximum (∼21 ka)\n70 bathymetry, and an acceleration technique with a factor of five applied to the external forcings\n71 (orbital, GHG, and ice-sheet). The impact of the experiment's design on the simulated sea ice\n72 is covered in the last section.\"\n74 distribution, such as a closed Bering Strait and subaerial Barents Shelf (Figure 2). We note that\n75 a blocked Bering Sea could cause the Atlantic Meridional Overturning Circulation to exhibit\n76 increased hysteresis behavior during glacial periods (Hu et al., 2012) and impair seawater\n77 transfer via the Bering Strait between the Arctic and the Pacific. However, we find that the\n78 simulated Pacific sea ice extent south of the Bering Strait (Supplementary Fig. 7) and the\n79 complete Pacific Ocean (Fig. 4d–f) exhibit similar seasonal and perennial sea ice variability\n80 similar to NH sea ice variability (Fig. 3). This is possible because the simulated NH meridional\nPage 3 of 48\n81 heat is transported equatorward mainly through the atmosphere for high obliquity and high\n82 atmospheric CO 2 levels (Yun et al., 2023) and it undermines the influence of a closed Bering\n85 formation during glacial-interglacial periods, influence surface albedo changes(Cao et al.,\n86 2019), and could warm the high-latitude(Cao et al., 2019). However, Yun et al. (2023) state that\n87 employing a volume-conserving ocean model with transient sea level and ocean volume would\n88 result in higher computational costs and perhaps unexplained biases and ambiguities that have\n89 not been thoroughly investigated in studies.\n91 boundary from the Last Glacial Maximum (∼21ka), resulting in ∼135 meters lower sea level\n92 than the present(Willeit et al., 2019). Consequently, we note that high-latitude has some\n93 unrealistic land-sea mask, such as a closed Bering Strait and subaerial Barents Shelf (Figure\n94 2). A closed Bering Sea would disrupt seawater transport between the Arctic and the Pacific\n95 via the Bering Strait, and it could lead to the emergence of stronger hysteresis behavior in the\n97 we find that the simulated Pacific sea ice extent south of the Bering Strait (Supplementary Fig.\n98 7) and the complete Pacific Ocean (Fig. 4d–f) show similar seasonal and perennial sea ice\n99 variability akin to NH sea ice variability (Fig. 3) possibly because the simulated NH meridional\n100 heat is transported equatorward mainly through the atmosphere for high obliquity and high\n101 atmospheric CO 2 levels(Yun et al., 2023). Additionally, the presence of a subaerial Barents\n103 albedo changes(Cao et al., 2019), and could warm the high-latitude(Cao et al., 2019). Then\n104 again, using a volume-conserving ocean model with time-evolving ocean volume and sea level\n105 would increase computational costs and may result in unidentified biases and ambiguities that\n106 studies have not fully explored(Yun et al., 2023).\n107\n109 here it must be taken into account that GHG is not only a forcing, but also subject to other\n110 changes in climate, and it indeed seems to be derived from output from a different simulation\nPage 4 of 48\n111 rather than from proxy-records. It is therefore (at in the opinion of this reviewer) not very\n112 surprising that one finds a change in the periodicity of GHG variability that follows the\n113 periodicity of glacial-interglacial cycles (41 kyr pre-MPT, 100 kyr post-MPT). Hence it would\n114 be useful to point out what we actually learn from such analyses. My apologies if this point has\n115 been missed by this reviewer.\n117 in the Methods section and the potential influence on sea ice in the discussion section.\n119 al., 2019) forcing along with transient prescribed orbital (A. L. Berger, 1978) and ice-\n120 sheet(Willeit et al., 2019) forcings and a constant ocean bathymetry and land-sea boundary\n121 from the Last Glacial Maximum (∼21ka).\n123 (∼21 ka), and it is forced with transient prescribed orbital(A. L. Berger, 1978), ice-sheet(Willeit\n124 et al., 2019), and greenhouse gases (including CO 2 , CH 4 , and N 2 O concentration)\n125 (GHG)(Willeit et al., 2019) forcings.\"\n128 experiment does not include an interactive carbon cycle. From 0.8 to 3 Ma, the estimated\n129 greenhouse gases consist of simulated CO 2 concentration from a previous CLIMBER-2\n130 transient simulation(Willeit et al., 2019), which contains an interactive carbon cycle.\n132 concentration from the Dome C ice core(Jouzel et al., 2007), European Project for Ice Coring\n133 in Antarctica (EPICA). The GHG forcing also includes estimates of N 2 O and CH 4\n134 concentrations as a CESM forcing obtained by regressing the EPICA Dome C ice core N 2 O\n135 and CH 4 concentrations on the LR04 global Benthic Stack proxy record(Lisiecki & Raymo,\n136 2005) (δ18O of benthic foraminiferal calcite, which measures global ice volume and deep\n137 ocean temperature) post 0.8 Ma and extending the linear statistical model to 3 Ma using the\n140 performed by Yun et al. (2023) employing the Community Earth System Model version 1.2\nPage 5 of 48\n141 (CESM1.2), which is forced with transient prescribed orbital(A. Berger, 1978), GHG(Willeit\n142 et al., 2019) (including CO 2 , CH 4, and N 2 O concentration), and ice-sheet forcings(Willeit et al.,\n143 2019).”\n145 datasets. From 0.8 to 3 Ma, the estimated greenhouse gases consist of simulated CO 2\n146 concentration from a previous CLIMBER-2 transient simulation(Willeit et al., 2019), which\n147 contains an interactive carbon cycle. Meanwhile, from 0 to 0.8 Ma, the GHG estimate includes\n148 air-bubble measurement of CO 2 concentration from the Dome C ice core(Jouzel et al., 2007),\n150 estimates of N 2 O and CH 4 concentrations as a CESM forcing obtained by regressing the EPICA\n152 record(Lisiecki & Raymo, 2005) (δ18O of benthic foraminiferal calcite, which measures\n153 global ice volume and deep ocean temperature) post 0.8 Ma and extending the linear statistical\n154 model to 3 Ma using the LR04 δ18O stack record.\n156 results could be because the GHG forcing anomalies relative to preindustrial conditions are\n157 amplified by a factor of 1.5 to compensate for relatively weak climate sensitivity in the applied\n158 model.\n160 because the GHG forcing anomalies relative to preindustrial conditions are amplified by a\n161 factor of 1.5 to compensate for relatively weak climate sensitivity in the applied model.”\n163 uncertain dust forcing and to partially offset the relatively weak climate sensitivity of CESM1.2\n164 compared to multi-model ensemble-based estimation from the Coupled Model\n166 preindustrial conditions have been amplified by a factor of 1.5.\"\n168 orbital forcing has a precession and obliquity signal throughout the Quaternary. In contrast, the\n169 prescribed CO 2 concentration (GHG forcing) and prescribed global ice volume in sea level\n170 equivalent (ice-sheet forcing) have a dominant ∼41 kyr periodicity during the pre–MPT period,\nPage 6 of 48\n171 ∼60–100 kyr (CO 2 /ice-sheet signal) periodicity during the MPT period and a longer ∼100 kyr\n172 periodicity during the post–MPT period. The significance of our analysis is that the presence\n173 of CO 2 /ice-sheet signal (∼60–100 kyr) during MPT and ∼100 kyr signal during post–MPT in\n174 the simulated sea ice/surface air temperature/sea surface temperature directly reflects the\n175 influence of GHG/ice-sheet forcings. Additionally, we agree with the reviewer that during the\n177 from the carbon cycle. Therefore, an interactive carbon cycle model would be more realistic.\n179 other model biases. Therefore, by prescribing GHG forcing, which involves CO 2 , CH 4 , and\n181 complex causality residing in the Earth's carbon cycle system.\n184 (obliquity) signal (Supplementary Fig. 1).\"\n186 (Supplementary Fig. 2a–c) and prescribed global ice volume in sea level equivalent(Willeit et\n187 al., 2019) (ice-sheet forcing) (Fig. 1f) (Supplementary Fig. 2j–l) show a dominant ∼41 kyr\n188 periodicity during the pre–MPT period, ∼60–100 kyr (CO 2 /ice-sheet signal) periodicity during\n189 the MPT period and a longer ∼100 kyr periodicity during the post–MPT period …..\"\n191 in mind that although the prescribed orbital forcing has a characteristic ∼41 kyr (obliquity) and\n192 ∼23/19 kyr (precession) signal throughout the Quaternary, it does not necessarily indicate a\n193 direct response to orbital forcing. This is because the prescribed CO 2 (GHG) (Supplementary\n195 kyr signal during the pre–MPT and a minor ∼23/19 signal across the Quaternary. Furthermore,\n196 the CO 2 /ice-sheet signal (∼60–100 kyr) during MPT and ∼100 kyr signal during post–MPT\n197 indicate a direct response to GHG/ice-sheet forcings.\n198\n199\nPage 7 of 48\n201 done for northward of 70N, it is precisely northward of that latitude as the authors chose this\n202 as a reference.\n204\n206 their definitions of sea ice extent and other variables, and what is considered NH high latitude.\n209 definitions earlier in the result section.\n211 ice concentration) in the Northern Hemisphere (NH; 35°N–90°N) during the post–MPT (0–0.7\n213 latitude is defined as the region north of 70°N, and the mid-latitude is between 35°N and 70°N\n214 latitudes.\"\n216 perennial sea ice which refers to sea ice that survives during summer and persists for over a\n217 year. Meanwhile, less than 80–90% of sea ice concentration is similar to seasonal sea ice that\n218 forms in winter and melts in summer.\"\n220 concentration, and use it consistently. There is also a general inconsistency in describing what\n221 drives certain processes. To say that the variability of a specific variable is more pronounced\n222 under precession-like periodicity or timescales is not the same as to affirm that the variable is\n223 precession-driven.\n225 inconsistencies, and made it clearer and more scientific, using the Figures as supporting\n226 material.\n227\nPage 8 of 48\n229 the first paragraph. Previous studies reconstruct or simulate sea ice extent/concentration (not\n230 variability – it varies over time) with limited timescale, and therefore this manuscript provides\n231 novel and important results? Past modeling studies also cover SAT and SST. How is the link to\n232 authors' findings on linkages between SAT, SST, and sea ice in the Arctic?\n234 inconsistencies, and made it clearer and more scientific. Subsequently, in the discussion section,\n235 our revised manuscript notes that previous studies are mainly limited to the post-MPT period,\n236 while our study examines in detail the spatiotemporal variability and the mechanism of NH sea\n237 ice variability during the post–MPT, MPT, and pre–MPT periods. Our study also emphasizes\n238 that our analysis of the extremely long 3 Ma simulation data provides an excellent test bed to\n239 explore the effects of orbital, GHG, and ice-sheet forcings on sea ice in more detail and with\n240 more precision compared to previous studies limited to the post-MPT period.\n242 et al., 2013; Sha et al., 2016; Stein et al., 2017; Turney et al., 2015) and modeling studies(Lo\n243 et al., 2018; Mantsis et al., 2011; Tuenter et al., 2005; Wu et al., 2020, 2022; Yin & Berger,\n244 2010) have mostly recorded or simulated sea ice concentration limited to the post-MPT. A\n245 recent study (3 Ma experiment) has presented a broad outline of the Northern and Southern\n247 ∼60–120 kyr (CO 2 /ice-sheet) periodicity, suggesting the combined effect of GHG and ice-sheet\n248 forcing(Yun et al., 2023) on sea ice variability. Here, we analyze the 3 Ma experiment to\n249 investigate the spatiotemporal and latitudinal variability and the mechanism of NH sea ice\n250 variability during the post–MPT, MPT, and pre–MPT periods..\"\n252 references indicated in the manuscript are surely neither descriptive of all datasets that include\n254 in previous experiments.\n256 records show a 100 kyr cycle during the post-MPT period (McClymont et al., 2013), which\n257 agrees with our result. Meanwhile, we have extensively edited the manuscript and removed\n258 this redundant information to make our manuscript clearer and more scientific.\nPage 9 of 48\n259\n261 into the effect of the forcing onto sea ice in a concise manner.\n263 descriptions in the results section and have clearly discussed it using Figure 8 in the revised\n264 manuscript.\n266 strong Pearson correlation of −0.79 and standardized regression coefficient of −0.76 with the\n267 obliquity index (lags obliquity index by ∼5 kyr) (green line in Fig. 8a). Meanwhile, the seasonal\n268 sea ice covering the mid-latitude (Fig. 2) has a strong Pearson correlation of 0.70 and\n269 standardized regression coefficient of 0.72 with precession index (lags precession index by\n270 ∼10 kyr) (green line in Fig. 8b). Therefore, we identify a strong coupling between simulated\n271 sea ice with orbital forcings throughout the Quaternary.\"\n273 exploring in the manuscript and keep it consistent. Also, there can be more description of what\n274 it entails to do seasonal studies (computational costs, resolution, boundary conditions,\n275 calendar corrections) to provide a clear but broad view for future studies.\n277 regarding sea ice and have clarified the definitions of sea ice extent and seasonal sea ice earlier\n278 in the result section. We have also replaced quasi-permanent sea ice with perennial sea ice.\n280 maintain conciseness.\n282 ice concentration) in the Northern Hemisphere (NH; 35°N–90°N) ……..\"\n284 perennial sea ice which refers to sea ice that survives during summer and persists for over a\n285 year. Meanwhile, less than 80–90% of sea ice concentration is similar to seasonal sea ice that\n286 forms in winter and melts in summer.\"\n289 selective at the moment. What about the other data, than model data, that was mentioned\n290 during the text?\n292 the 3 Ma climate model experiment design and forcings in the Methods section.\n294 simulation for the past 3 million years employing the coupled CESM version 1.2 (CESM 1.2).\n296 4; CAM4), ocean (Parallel Ocean Program version 2; POP2), sea ice (Los Alamos Sea Ice\n298 components. The atmospheric model has about 3.75° × 3.75° latitude-longitude horizontal\n299 resolution (T31) with 26 vertical levels. The ocean and sea ice model has a 100 × 116\n300 horizontal grid (gx3v5), and the ocean model has 25 vertical levels. The ocean and sea ice\n301 output are gridded to a 3.75° × 3.75° latitude-longitude horizontal resolution identical to\n302 atmospheric outputs.\n304 and it is forced with transient prescribed orbital(A. L. Berger, 1978), ice-sheet(Willeit et al.,\n305 2019), and greenhouse gases (including CO 2 , CH 4 , and N 2 O concentration) (GHG)(Willeit et\n306 al., 2019) forcings.\n309 a previous CLIMBER-2 transient simulation(Willeit et al., 2019), which contains an interactive\n310 carbon cycle. Meanwhile, from 0 to 0.8 Ma, the GHG estimate includes air-bubble\n311 measurement of CO 2 concentration from the Dome C ice core(Jouzel et al., 2007), European\n313 and CH 4 concentration as a CESM forcing obtained by regressing the EPICA Dome C ice core\n315 2005) (δ18O of benthic foraminiferal calcite, which measures global ice volume and deep\n316 ocean temperature) post 0.8 Ma and extending the linear statistical model to 3 Ma using the\nPage 11 of 48\n319\n321 forcing and to partially offset the relatively weak climate sensitivity of CESM1.2 compared to\n322 multi-model ensemble-based estimation from the Coupled Model Intercomparison Project\n324 have been amplified by a factor of 1.5. The 3 Ma simulation also uses an acceleration factor of\n325 5 to the external forcings (orbital, ice-sheet, and GHG forcings), which reduces the 3 Ma orbital\n326 years forcing to 0.6 Ma model years. The orbital, ice-sheet, and GHG forcings are updated\n327 every 100 orbital years, which means 20 model years.\n329 reconstructions and numerous paleo-proxy data. The 3 Ma simulation has been shown to\n330 reproduce the magnitude and timing of glacial–interglacial variability(Timmermann et al.,\n331 2022; Yun et al., 2023) consistent with paleoclimate data constraints and estimates from the\n332 intermediate complexity Earth system model. The 3 Ma simulation agrees well(Timmermann\n333 et al., 2022; Yun et al., 2023) with key paleoclimate proxy records such as the magnitude and\n334 timing of global mean surface temperature variations, annual mean sea surface temperatures\n335 over the tropics and elsewhere (western south and north Pacific, southern Indian and north\n339 al. (2023) provide the complete details on the 3 Ma model set-up and experimental design.”\n340\n343 time-series used different colors for the parameters examined instead of a solid black line for\n344 some and a dashed red line for others. In this way, key results such as the surface air\n345 temperature are hardly visible due to clutter from the overlying time-series. The authors can\n346 also describe in the figure caption the direction of the time axis instead of using the arrows as\nPage 12 of 48\n347 it increases the clutter in the image and decreases understanding of the processes. The slope\n348 of the curves between the post-MPT, MPT, and pre-MPT periods could be described in a table\n349 or just mentioned in the text.\n351 would enable readers to better understand your plots without going back and forth between\n352 manuscript text and figure caption.\n354 the text.\n356 trend of the incoming solar radiation did not appreciably change during the past 2.6 Ma (0.0001\n359 2a–c) and simulated global mean surface temperature (Supplementary Fig. 6) have similar\n360 periodicity signals throughout the Quaternary. During the pre–MPT, the prescribed CO 2\n361 concentration (Fig. 1e) and simulated global mean surface temperature (Fig. 1h) have a gradual\n362 downward trend, intensifying across the MPT by about 2 and 3 times, respectively. On the\n363 contrary, the simulated sea ice extent accompanies an upward trend (significant at a 95%\n364 confidence level) during the pre–MPT (1.5 thousand km2/kyr), which intensifies by a factor\n365 of about 2.2 (3.4 thousand km2/kyr) during the MPT (Fig. 1g) indicating the influence of GHG\n366 and/or ice-sheet forcing on sea ice.\"\n368 references to data sources in the figure caption and main manuscript.\n370 performed by Yun et al. (2023) employing the Community Earth System Model version 1.2\n371 (CESM1.2), which is forced with transient prescribed orbital(A. Berger, 1978), GHG(Willeit\n372 et al., 2019) (including CO 2 , CH 4, and N 2 O concentration), and ice-sheet forcings(Willeit et al.,\n373 2019).”\n375 obliquity) (Fig. 1a–c) influences the intensity of summer (June–July–August; JJA) ……..”\nPage 13 of 48\n377 1e) (Supplementary Fig. 2a–c) and prescribed global ice volume in sea level equivalent(Willeit\n378 et al., 2019) (ice-sheet forcing) (Fig. 1f) (Supplementary Fig. 2j–l) show a dominant ∼41 kyr\n379 periodicity ……………………………..”\n382 series shows the prescribed orbital forcing: (a) eccentricity index(A. L. Berger, 1978), (b)\n383 precession index(A. L. Berger, 1978) (𝑒sinω, where 𝑒 is the eccentricity and ω is the longitude\n384 of the perihelion), (c) obliquity index(A. L. Berger, 1978) (units are in degrees); (d) summer\n385 (JJA) incoming solar radiation at 65°N (units are in watts per square meter), (e) prescribed CO 2\n386 concentration(Willeit et al., 2019) (units are in ppmv), (f) prescribed global ice volume in sea\n387 level equivalent(Willeit et al., 2019) (units are in meters), (g)……..\"\nPage 14 of 48\n388\n390 of the separate representations (one focuses on variability, the other on the mean-state).\n392 - unless negative values should illustrate model land sea mask or similar (see my comment\n393 below regarding land sea masks). The meaning of values below zero should then be clearly\n394 documented, and ideally using a color different than white (which currently just means \"no sea\n395 ice \").\n398 similar between all three time periods, a sideeffect of different definitions?\n400 area in the Barents Sea stems from a model land sea mask with subaerial Barents Shelf that is\n401 not shown on the plot? In particular, I assume the model considers a subaerial Barents Shelf,\n402 i.e. the white color shown there does not imply that there are ice free ocean regions, am I right?\n404 than overlaying a (high-resolution) modern land sea mask - otherwise readers will be confused\n405 regarding sea ice dynamics. In this respect, more information on the paleogeography used in\n406 the model, and potential differences between the studied time periods, should be presented\n407 somewhere in the manuscript. This is particularly relevant for continental ice sheet forcing as\n408 you argue that this forcing would control high latitude sea ice in the Northern Hemisphere.\n410 updated manuscript. As suggested, the contour value/scale is from 0 to 100% sea ice\n411 concentration. Figure 2 in the updated manuscript also shows subaerial land regions (especially\n412 subaerial Barents Shelf and closed Bering Strait) as the experiment uses a constant Last Glacial\n414 of using a constant Last Glacial Maximum (∼21ka) land-sea boundary.\nPage 15 of 48\n417 pre–MPT (1.2–2.6 Ma) periods. The shaded blue color shows the sea ice coverage over modern\n418 land-sea boundaries. The units are in percentage sea ice concentration. The grey dots show\n419 subaerial land regions (especially subaerial Barents Shelf and closed Bering Strait) as the\n420 experiment uses a constant Last Glacial Maximum land-sea boundary. The black and purple\n421 lines demarcate the sea ice extent and perennial sea ice boundary, respectively.\"\n423 boundary from the Last Glacial Maximum (∼21ka), resulting in ∼135 meters lower sea level\n424 than the present(Willeit et al., 2019). Consequently, we note that high-latitude has some\n425 unrealistic land-sea distributions, such as a closed Bering Strait and subaerial Barents Shelf\n426 (Figure 2). A closed Bering Sea would disrupt seawater transport between the Arctic and the\n428 in the Atlantic Meridional Overturning Circulation during glacial periods(Hu et al., 2012).\n430 (Supplementary Fig. 7) and the complete Pacific Ocean (Fig. 4d–f) show similar seasonal and\n431 perennial sea ice variability akin to NH sea ice variability (Fig. 3) possibly because the\n432 simulated NH meridional heat is transported equatorward mainly through the atmosphere for\n433 high obliquity and high atmospheric CO 2 levels(Yun et al., 2023). Additionally, the presence\n434 of a subaerial Barents Shelf would reduce sea ice formation during glacial-interglacial periods,\n435 influence surface albedo changes(Cao et al., 2019), and could warm the high-latitude(Cao et\n436 al., 2019). Then again, using a volume-conserving ocean model with time-evolving ocean\n437 volume and sea level would increase computational costs and may result in unidentified biases\n438 and ambiguities that studies have not fully explored(Yun et al., 2023).\"\n440 (above 80% sea ice concentration) and MPT/post-MPT (above 90% sea ice concentration). We\nPage 16 of 48\n441 reason that during the pre-MPT period, the sea ice is less expanded than in the post-MPT period\n442 because of relatively higher pre-MPT greenhouse gasses concentration and global mean surface\n443 temperature. Therefore, we used slightly different definitions and argued that, in general,\n444 perennial sea ice covered the high-latitude (north of 70°N latitude) regions throughout the\n446 \"above 90%\" in the definition, as shown in Supplementary Figure 4 in the revised manuscript.\n449 for the past 2.6 Ma. The purple dashed line shows the average sea ice concentration at 72°N\n450 during the post–MPT, MPT, and pre–MPT periods. The units are in percentage sea ice\n451 concentration. We notice that sea ice concentration in the high-latitude region (north of 70°N\n452 latitude), such as those at 72°N and 76°N, is covered with more than 80% sea ice concentration\n453 during the pre-MPT and more than 90% sea ice concentration during the MPT and post-MPT\n454 periods.\"\n455\n458 in the caption for Figure 3.\nPage 17 of 48\n461 standardized variance) during…..\"\n463 redundant.\n465 have clarified the domain analyzed in our manuscript earlier in the result section.\n467 ice concentration) in the Northern Hemisphere (NH; 35°N–90°N) ……..\"\n468\n470 circles in Fig. 4a-c to make it easier for the reader to identify the progressive shift in sea ice\n471 extent.\n473 in the revised manuscript. As suggested, we have shown more latitude circles in Figure 5a-c\n474 and presented Figure 5d-i more clearly.\n477 periodicity-variance spectrum of sea ice in those major sea ice variability regions, as shown in\n479 to identify the progressive shift in sea ice extent.\n481 authors can just refer to the indices in the title and keep only one legend relative to the region\n482 (using solid lines only instead of both solid and dashed lines).\nPage 18 of 19\n486 \"Fig. 8: Cross-correlation of summer incoming solar radiation at 65°N, surface air\n487 temperature, and sea ice with obliquity and precession indices (orbital parameters) from 0–\n488 2.6 Ma. (a) ….\"\n489\n491 coefficients are misleading.\nPage 19 of 19"}]

Decision letter
Decision Letter

Dear Professor An,

We are pleased to inform you that your manuscript, "Northern Hemisphere sea ice variability in a transient CGCM simulation of the past 2.6 Ma", has been accepted for publication in Nature Communications.

Thank you for your patience during the review process and for your efforts in addressing the reviewers' comments. We believe the manuscript has been significantly improved and is now suitable for publication.

Best regards,

Editorial Team
Nature Communications