Orbital Forcing Formula
Earthβs climate-relevant orbital forcing arises from the gravitational interplay of all eight planets. Their orbital and rotational cycles synchronise over a common Solar System Resonance Cycle of 8H = 2,682,536 years, and every climate-relevant cycle on Earth therefore lands at an integer divisor of 8H. Spectral analysis of the LR04 benthic δ¹βΈO stack (Lisiecki & Raymo 2005) confirms this structure empirically and yields an explicit predictive formula.
Orbital forcing is not climate. The formula on this page captures the orbital-forcing component only. Joint OLS fit on LR04 explains ~24 % of the observed variance (RΒ² = 0.238); the remaining ~76 % comes from non-orbital climate-system response β ice-sheet hysteresis, COβ and carbon-cycle feedbacks, internal variability. The model takes no position on those components. Orbital cycles are the clock that sets the timing of glacial/interglacial transitions; the magnitude of the observed climate signal is dominated by Earth-system feedbacks, not orbital forcing directly. Every prediction here β including the next-glaciation forward projection β describes when the orbital clock makes a phase transition possible, not when surface climate necessarily follows.
This page summarises the empirical climate findings; the full technical record (methodology, pre-registration, reproducibility scripts) lives in Milankovitch Evidence on GitHub (doc 17)Β .
Five headline findings
-
Every significant LR04 climate peak sits at an integer divisor of 8H. The formula has 26 such integers β 20 detected above the 3Γ median significance threshold in full LR04, plus 6 more present in pre-MPT data. 25 of the 26 have clean physical interpretations as standard celestial-mechanics beats (k+g_j, k+s_j, g_jβg_k, s_jβs_k) or direct planet apsidal/nodal cycles from the Solar System Resonance Cycle period tableΒ .
-
Mars dominates the per-planet climate fingerprint. Two exclusive direct matches in LR04 full (n=35 Mars apsidal, n=53 Mars eccentricity cycle) and three more exclusive matches in pre-MPT (n=16 Mars Axial, n=21 Mars Obliquity, n=53 confirmed). Marsβs near-resonance with Earthβs apsidal eigenmode (gβ β 17.37 β³/yr, gβ β 17.92 β³/yr) produces the cleanest direct climate signal of any planet.
-
The 100-kyr glacial cycle is an inclination-side eigenmode beat, not direct eccentricity forcing. The empirical centroid sits at the Mercury-Mars sββsβ nodal beat at 107.3 kyr β a planet-pair orbital-plane coupling. Direct eccentricity attribution faces specific failure modes (405-kyr term essentially absent in post-MPT LR04, no bispectral 95k+125k coupling). Vindicates Muller & MacDonald 1997βs βinclination, not eccentricityβ framing with a specific eigenmode label they didnβt have.
-
Pre-MPT and post-MPT differ in climate sensitivity, not orbital forcing. Orbital forcing is essentially stationary over 5.3 Myr. LR04βs volatility growth from past to present reflects climate-system response changing β Northern Hemisphere ice-sheet establishment and the Mid-Pleistocene Transition crossing a hysteresis threshold.
-
Forward projection: the next natural glaciation peak is predicted in ~38,000 years (~40,000 AD). The strongest predicted glaciation in the next quarter-million years sits at ~194,500 years from now. Consistent with the classical Berger & Loutre 2002 prediction within model uncertainty. Anthropogenic COβ may delay the next natural glaciation by 50+ kyr (Ganopolski et al. 2016) β not modelled here.
The 8H Orbital Forcing Formula
Definition
with:
- t = age in kyr BP (positive = past, negative = future, t = 0 β 2000 AD)
- 8H = 2,682,536 years (Solar System Resonance Cycle)
- N = the 26 active integer divisors listed below
- a_n, b_n = OLS-fitted coefficients (amplitude = β(a_nΒ² + b_nΒ²))
- C(t) = normalised δ¹βΈO proxy (positive = colder/glacial, negative = warmer/interglacial)
Joint OLS fit on full LR04 (T = 5,320 kyr) gives RΒ² = 0.238, condition number 1.6 (all 26 candidates Rayleigh-resolvable β no collinearity), and past-200-kyr local RΒ² = 0.320.
Three components dominate by amplitude: obliquity (n=65) at 0.275, MarsβJupiter eccentricity beat (n=28) at 0.238, MercuryβMars nodal beat (n=25) at 0.213.
The 26 active integer divisors
Each n corresponds to a specific celestial-mechanics quantity. Letters refer to Laskar 2004 secular eigenfrequencies (g_j apsidal, s_j nodal); k = Earthβs general precession in longitude (50.29β³/yr). Planet labels on g_j / s_j follow Bergerβs convention; see Eigenfrequencies Β§βThe eigenmodes are realβ¦β for the modelβs frame.
| n | period (kyr) | amp | physical interpretation |
|---|---|---|---|
| 7 | 383.2 | 0.110 | gββgβ long eccentricity (~405k) |
| 9 | 298.1 | 0.111 | gββgβ eccentricity beat / Mercury Axial = AscNode = 8H/9 |
| 12 | 223.5 | 0.104 | sβ βsβ nodal beat / Uranus AscNode = 8H/12 |
| 14 | 191.6 | 0.061 | gββgβ eccentricity beat |
| 16 | 167.7 | 0.082 | Mars Axial precession = 8H/16 (model direct) |
| 18 | 149.0 | 0.122 | sββsβ nodal beat |
| 20 | 134.1 | 0.103 | gββgβ eccentricity beat |
| 21 | 127.7 | 0.084 | Mars Obliquity / Jupiter Axial = 8H/21 (model direct) |
| 22 | 121.9 | 0.095 | sββsβ nodal / gββgβ |
| 25 | 107.3 | 0.213 | sββsβ Mercury-Mars nodal (100-kyr centroid) |
| 28 | 95.8 | 0.238 | gββgβ eccentricity beat (Berger 95k) |
| 30 | 89.4 | 0.124 | gββgβ eccentricity beat |
| 31 | 86.5 | 0.155 | gββgβ eccentricity beat |
| 35 | 76.6 | 0.126 | Mars apsidal = 8H/35 (model direct) |
| 38 | 70.6 | 0.104 | sββsβ nodal beat |
| 39 | 68.8 | 0.164 | sβ βsβ Earth nodal precession |
| 48 | 55.9 | 0.101 | sββsβ nodal beat |
| 50 | 53.7 | 0.156 | gββgβ eccentricity beat |
| 53 | 50.6 | 0.138 | Mars Eccentricity cycle = 8H/53 (model direct) |
| 65 | 41.27 | 0.275 | k+sβ Earth obliquity (Berger 41k) β dominant peak |
| 66 | 40.6 | 0.043 | obliquity-band arithmetic-mean cycle length (cycle-counting artifact, near-zero at full LR04 resolution) |
| 68 | 39.5 | 0.162 | k+sβ obliquity sub-peak |
| 73 | 36.8 | 0.092 | 2|sβ| / gββsβ |
| 76 | 35.3 | 0.091 | gββsβ apsidalβnodal beat |
| 113 | 23.74 | 0.086 | k+gβ climatic precession sub-peak (Berger 23.7k) |
| 120 | 22.35 | 0.104 | k+gβ climatic precession sub-peak = H/15 (Berger 22.4k) |
Six integers correspond directly to specific planet cycles from the Solar System Resonance Cycle period tableΒ : n=9 (Mercury Axial), n=12 (Uranus AscNode), n=16 (Mars Axial / Jupiter Obliquity), n=21 (Mars Obliquity / Jupiter Axial β βMarsβJupiter Axial-Obliquity Swapβ), n=35 (Mars Perihelion ecliptic), n=53 (Mars Eccentricity cycle).
The other 20 integers are eigenmode beats between planet pairs. All 26 land on integer divisors of 8H because 8H is the natural synchronisation period of the solar system.
In-app Explorer
Try it interactively. The 3d simulationΒ ships with an Orbital Forcing Formula Explorer modal (Tools menu) that plots the formulaβs 26-component reconstruction directly on top of the LR04 stack. Five tabs span the full 5.3 Myr record, post-MPT, the last 200 kyr, and a 250-kyr forward projection. The x-axis runs past β future (today on the left edge of the projection tab); years are labelled BC / AD with today = 2000 AD.
Closure: no orphan peaks off the 8H lattice
The strongest possible test of the 8H framework is its closure: do all significant spectral peaks in LR04 land on integer divisors of 8H, or are there orphan peaks at non-integer positions that would imply forcing from outside the planetary-eigenmode framework? A peak at, say, n = 43.5 with non-negligible amplitude would falsify the framework.
We performed this closure test by fitting all 200 integer divisors of 8H jointly to LR04 (joint OLS, RΒ² = 0.443) and scanning the residual for non-integer power. The result:
- Residual amplitude at integer positions: 0.000 (machine zero, by construction)
- Noise floor at random non-integer positions: 0.029 median, 0.127 (95th percentile)
- 14 residual peaks above noise threshold sitting >0.3 from any integer β every single one between two adjacent integer divisors in or near the formula
The biggest βorphanβ (n = 65.45, period 40.99 kyr, amp 0.544) sits between integers 65 (k+sβ Earth obliquity, 41.27 kyr) and 66 β the cycle-counting arithmetic mean position of the non-stationary obliquity band, exactly the Jensenβs-inequality phenomenon documented in doc 17 Β§6.6Β . The other 13 orphans similarly sit between adjacent integer beats β they are the expected fingerprint of cycle-length dispersion (LR04 obliquity cycles vary 31β59 kyr) and spectral leakage between close integer signals.
Crucially: zero orphan peaks land in βemptyβ regions of the 8H lattice β no peaks at positions like n = 12.7, n = 42.3, n = 80, or n = 140 that would suggest genuinely off-lattice forcing. The 8H integer-divisor framework captures the full frequency structure of LR04βs significant spectral content.
This is the third independent empirical confirmation of the framework, alongside the 405-kyr-absence test and the bispectral coupling absence (see Eccentricity attribution has specific empirical headwinds above).
Per-Planet Contributions
For each climate peak, cross-referencing against the full doc 55 period table (8 planets Γ 6 cycle types) reveals which planets contribute directly versus via eigenmode beats.
Mars dominates
LR04 full window (T = 5,320 kyr):
| Planet | Exact direct matches | Total (with Β±1 near) | Exclusive |
|---|---|---|---|
| Mars | 2 | 5 | n=35, n=53 β both unique to Mars |
| Jupiter | 0 | 5 | none (shared with Mars, Earth, Saturn, Uranus) |
| Saturn | 0 | 4 | none |
| Mercury | 2 | 3 | n=9 (Axial = AscNode by Cassini lock) |
| Venus | 0 | 3 | none (via eigenmode beats) |
| Earth | 0 | 3 | none (Fibonacci H/n positions empty, see below) |
| Uranus | 1 | 2 | n=12 (AscNode) |
| Neptune | 0 | 0 | none direct in LR04 full |
Marsβs two exclusive matches (n=35 Mars ecliptic perihelion = 8H/35 at 76.6 kyr; n=53 Mars eccentricity cycle = 8H/53 at 50.6 kyr) are doc-55 entries unique to Mars. Probability of two such exclusive matches by chance is roughly (6/46)Β² β 1.7%.
Neptune and Uranus appear only pre-MPT
The pre-MPT window (1,200β3,000 kyr BP) reveals three additional peaks not in LR04 full that involve outer planets:
- n = 14 (191.6 kyr) β gββgβ Venus-Neptune eccentricity beat
- n = 30 (89.4 kyr) β gββgβ Earth-Uranus eccentricity beat
- n = 38 (70.6 kyr) β sββsβ Neptune-Earth nodal beat
Neptuneβs eigenmodes are the slowest (gβ β 0.67 β³/yr, |sβ| β 0.69 β³/yr); beats with inner planets produce relatively long periods (~70β200 kyr) with small amplitudes. Pre-MPT, ice sheets were smaller and less hysteretic β slow outer-planet beats registered above the noise floor. Post-MPT large ice sheets impose a ~100-kyr hysteretic response that dominates and filters out the slower outer-planet signals.
Cross-planet obliquity validation
The modelβs obliquity-period predictions for the inner solar system match three independent published references to within 2.5 % with zero free parameters:
| Planet | Published period | Reference | Model H/n | Deviation |
|---|---|---|---|---|
| Mercury | 895,000 yr | Bills & Comstock 2005 | 8H/3 = 894,179 yr | +0.09 % |
| Earth | 41,000 yr | Laskar 2004; Berger 1978 | H/8 = 41,915 yr | +2.2 % |
| Mars | 124,800 yr | Ward 1973; Laskar 2004 | 8H/21 = 127,740 yr | +2.4 % |
Mercuryβs 0.09 % match against the independent Cassini-state forced-obliquity calculation (Bills & Comstock used Cassini-state dynamical theory, corroborated by Yseboodt & Margot 2006, Peale 2005, Bois & Rambaux 2007) is the modelβs tightest cross-validation against non-Holistic published references.
The 100-kyr Glacial Cycle
The 100-kyr glacial cycle is paleoclimatologyβs most-debated feature. The data point clearly to an inclination-side / orbital-plane origin via planet-pair eigenmode beats rather than direct eccentricity forcing.
What the spectrum shows
LR04βs 100-kyr band is a broad single peak spanning ~80β125 kyr β not the split structure (95k + 125k) that direct eccentricity forcing would produce. The energy-weighted centroid sits at n = 25 (= 107.3 kyr in 8H units), with strong adjacent contributions at n = 28 (95.8 kyr) and weaker contributions at n = 22 (~122 kyr).
Empirical centroid: Mercury-Mars nodal beat
The dominant 100-kyr-band integer n = 25 corresponds to the sβ β sβ Mercury-Mars nodal beat in eigenmode terms (predicted 25.11, error 0.11). This is a nodal (orbital-plane) beat between two inner rocky planets, not an eccentricity beat. Consistent in spirit with Muller & MacDonald 1997βs βinclination, not eccentricityβ framing, with a specific eigenmode identification not previously stated.
Eccentricity attribution has specific empirical headwinds
Direct-eccentricity attribution faces three discriminating empirical failure modes that the inclination-side framework does not:
| Test | Eccentricity attribution | Inclination-side framework |
|---|---|---|
| 405-kyr term presence (Bergerβs strongest predicted eccentricity beat) | Predicts dominant; observed essentially absent in post-MPT LR04 (amplitude ratio 0.12) | No contradiction |
| Eccentricity-beat phase coupling (95k+125k bispectrum) | Predicts strong; not detected (bicoherence 0.507 < null-95 0.555 β replicates M-M 1997) | No contradiction |
| 100-kyr centroid identification | Predicts 95β99 kyr eccentricity beats; observed centroid at 107 kyr is a nodal beat | Mercury-Mars nodal (sββsβ) at n=25 = 107 kyr β matches |
| Berger climatic-precession 6-peak spectrum | Not predicted (these are k+g_j sub-peaks, not eccentricity-beat periods) | All 6 peaks match H/n divisors to < 0.4 % |
Theoretical proposal: H/3 second-obliquity component
Within the inclination-side family, the model proposes Earthβs intrinsic inclination precession (H/3 = ~111,772 years = 111.77 kyr) produces a second obliquity component, contributing alongside the Mercury-Mars planet-pair nodal beat. The mechanism (no dust required):
Inclination precession (H/3) β second obliquity component at H/3 β standard Milankovitch insolation forcing β ice sheets
Every step after βsecond obliquity componentβ is standard Milankovitch physics; the mechanism needs no new forcing.
At T = 1.2 Myr the Rayleigh resolution at P = 110 kyr is ΞP β 10 kyr, so 95k / 100k / 112k are spectrally collinear. H/3 lies within one Rayleigh element of the empirical 107-kyr centroid and remains a candidate within the inclination-side family β but cannot be singled out from the Mercury-Mars nodal beat at this data length. The empirical signal is consistent with this family; the specific H/3 second-obliquity contribution is theoretical.
Earthβs intrinsic Fibonacci H/n positions are empty in the climate spectrum
The modelβs intrinsic Earth precession periods (Fibonacci divisors of H β H/3, H/5, H/8, H/13, H/16) all sit at near-zero amplitude in LR04 (0.014 to 0.078, none above 0.085). This is structurally consistent with eigenfrequencies Β§βThe eigenmodes are realβ¦β: the Fibonacci H/n periods are Earthβs intrinsic geometric precession in its own orbital frame; climate observes Earth in the heliocentric frame where planetary-induced orbital-plane motion adds in, shifting the observed signal to integers adjacent to the Fibonacci anchors (Bergerβs k+g_j and k+s_j sub-peaks structure).
Pre-MPT vs Post-MPT
Same forcing, different response
Three structural transitions are visible in the full LR04 record:
| Approximate epoch | Transition | Effect on climate spectrum |
|---|---|---|
| ~2.7 Ma BC | Late Pliocene cooling onset (iNHG) | Northern Hemisphere ice sheets establish |
| 2.7 β 1.2 Ma BC | Early Pleistocene β41-kyr worldβ | Obliquity dominates |
| ~1.2 β 0.7 Ma BC | Mid-Pleistocene Transition (MPT) | Shift to β100-kyr worldβ: ice-sheet hysteresis crosses threshold |
The orbital forcing β planetary eigenmodes, Earthβs precession, obliquity oscillation β is essentially constant over the past 5 Myr. What changed is the climate systemβs sensitivity: long-term COβ decline, growing Northern Hemisphere ice sheets, and ice-sheet hysteresis crossing the Willeit threshold (~1 Ma).
The 8H formula has stationary amplitudes, so the reconstruction shows the same oscillation amplitude in pre-MPT as post-MPT. LR04βs actual amplitude grows from past to present. The mismatch is the expected outcome: the formula captures orbital forcing; LR04 amplitude variation reflects changing climate sensitivity.
MPT amplitude growth pattern
Windowed amplitude analysis across the MPT (pre-MPT 1,500β2,500 kyr vs post-MPT 0β1,000 kyr) gives concrete post/pre growth ratios:
| Band | Preβpost ratio | Direction |
|---|---|---|
| 41-kyr obliquity | 0.72Γ | shrank |
| 100-kyr band (n = 25, 28) | 1.64Γ | grew |
| 23.7-kyr climatic-precession (n = 113) | 2.19Γ | grew most |
| 405-kyr eccentricity (n = 7) | 0.34Γ | shrank |
The 41-kyr peak decreased β consistent with the standard βice-sheet saturation silences the obliquity pacemakerβ framing (Willeit 2019). The 100-kyr and 23-kyr bands both grew, with the 23-kyr growth concentrated in one specific sub-peak (k+gβ ) rather than spread across the precession band.
Forward Projection
The next 250,000 years
Evaluating C(t) at negative t (future):
Holocene check at t = 0: C(0) β β0.44 (normalised) β interglacial β (the orbital signal is well below zero, but already past its deepest interglacial value)
Phase-transition timeline and predicted glacial maxima:
| Years from now | AD date | C(t) | Phase |
|---|---|---|---|
| 0 (today) | 2000 AD | β0.44 | Holocene interglacial |
| ~5,700 | ~7,700 AD | β0.56 | Peak orbital interglacial warmth β sustained cooling begins after this |
| ~32,300 | ~34,300 AD | 0.00 | Orbital signal crosses zero (interglacial β glacial-favoring) |
| 38,000 | ~40,000 AD | +0.25 | next natural glaciation onset |
| 81,500 | ~83,500 AD | +0.45 | moderate |
| 110,500 | ~112,500 AD | +0.32 | mild |
| 154,500 | ~156,500 AD | +0.67 | strong |
| 194,500 | ~196,500 AD | +0.92 | strongest in window |
The orbital signal places peak interglacial warmth around ~7,700 AD β after that, the formula predicts sustained orbital-driven cooling, with the signal becoming glacial-favoring around ~34,300 AD and reaching the first glacial maximum at ~40,000 AD. The 9-kyr LGM offset in past validation (predicted 29 kyr BP vs. observed ~20 kyr BP) suggests the surface-temperature response can lag the orbital clock by several kyr; the actual peak-warmth surface date may therefore differ from the ~7,700 AD orbital peak. The formula does not include anthropogenic forcing.
Validation against past glacials
| Formula prediction | Actual paleoclimate event | Discrepancy |
|---|---|---|
| Glacial maximum at 29 kyr BP | Last Glacial Maximum β 20 kyr BP | +9 kyr early |
| Glacial maximum at 138 kyr BP | MIS 6 β 140 kyr BP | β2 kyr |
The LGM offset of ~9 kyr is expected: pure orbital forcing predicts the insolation drive, not the ice-sheet response. Ice sheets carry memory; peak ice volume lags peak orbital cooling by several kyr (standard climate physics). The MIS 6 match within 2 kyr is excellent.
Pacing shift: the next 250 kyr looks more like a 41-kyr world
A striking feature of the forward projection is that the intervals between predicted glacial peaks are not ~100 kyr β they cluster near the obliquity band (~40 kyr):
| Interval | Years between peaks |
|---|---|
| Glacial onset β 2nd | 43.1 kyr |
| 2nd β 3rd | 29.2 kyr |
| 3rd β 4th | 43.8 kyr |
| 4th β 5th | 40.4 kyr |
| 5th β 6th | 48.5 kyr |
| Mean | ~41 kyr |
By contrast the past ~700 kyr in LR04 was dominated by ~100-kyr pacing (LGM β MIS 6 β MIS 8 β MIS 10 β MIS 11, mean interval ~100 kyr). The forward projectionβs ~40-kyr orbital pacing therefore resembles the pre-MPT β41-kyr worldβ (Early Pleistocene, ~2.7β1.2 Ma) more than the late-Pleistocene 100-kyr regime.
Why: Earthβs eccentricity is heading into a long-term minimum (the Venus-Jupiter gββgβ ~400-kyr envelope), suppressing the climatic-precession amplitude (eΒ·sin Ο) and letting the obliquity band (k+sβ at 41 kyr) show through cleanly; the 100-kyr-band inclination-side beats happen to be partially out of phase for the next ~250 kyr. This is the same orbital prediction, from a different angle, as the famous Berger & Loutre 2002Β βexceptionally long interglacial aheadβ paper.
Climate response is a separate question. Whether Earthβs surface climate actually tracks this obliquity-band orbital pacing depends on ice-sheet hysteresis. Three scenarios: (A) ice sheets shrink below the Willeit 2019 threshold and the climate response re-couples to the orbital clock (a βreturn to 41-kyr worldβ); (B) post-MPT hysteresis persists and the climate keeps producing ~100-kyr cycles even with weaker orbital triggers (βskippedβ or subdued glacials); (C) anthropogenic COβ overrides orbital pacing entirely for the near-term (Ganopolski 2016). The 8H formula speaks to the orbital half only.
Empirical analogue: the late Pliocene (8H ago)
The 8H formula is exactly periodic at 8H β every component is an integer divisor of 8H, so C(t) = C(t + 8H) by construction (verified to floating-point precision). This means the orbital signal in the next 250 kyr is byte-identical to the orbital signal in 2.43β2.68 Ma BC, making the late-Pliocene LR04 record a direct empirical analogue for how Earthβs climate responded to the same forcing weβre about to enter.
| Metric | Late Pliocene (2.43β2.68 Ma BC) | Post-MPT (0β250 kyr BP) |
|---|---|---|
| LR04 correlation with formula C(t) | r = 0.760 | r = 0.675 |
| LR04 normalised range | [β1.62, +2.34] | [β4.18, +2.89] |
| Amplification factor (LR04 std / orbital std) | 2.4Γ | 4.0Γ |
| Glacial-peak intervals (kyr) | 32, 46, 35, 44 (mean 39.2) | 44, 25, 22, 31, 45, 20, 18 (mean 29.3) |
| Dominant spectral period | 41.8 kyr (amp 1.04) | 41.8 kyr (amp 0.55) |
Three findings:
- The orbital signal is empirically 41-kyr-paced in the late-Pliocene window β LR04 intervals cluster at the obliquity period, spectral peak at 41.8 kyr is dominant. This confirms the next 250 kyr will look 41-kyr-paced at the orbital level.
- The climate response was cleaner and more orbital-coupled in the late Pliocene. Higher correlation (0.76 vs 0.68) and lower amplification (2.4Γ vs 4.0Γ) β without large hysteretic ice sheets, climate tracked the orbital signal more directly.
- The climate envelope reachable by orbital forcing alone is modest. Late-Pliocene normalised range is narrower than post-MPT β if ice sheets shrink below the Willeit threshold and climate re-couples to the orbital clock, natural climate amplitude over the next 250 kyr would resemble the late Pliocene more than the late Pleistocene.
This empirical anchor strengthens Scenario A above: the late Pliocene happened, the LR04 record exists, and it shows what an orbital-coupled climate response to this exact orbital signal looks like β moderate-amplitude, 41-kyr-paced, less extreme than post-MPT swings. Late-Pliocene global mean temperature was ~2β3 Β°C warmer than pre-industrial, comparable to projected anthropogenic warming β so Scenarios A and C may converge if anthropogenic COβ pushes Earth toward a Pliocene-like baseline.
Limitations
The formula captures orbital forcing only β not ice-sheet hysteresis, carbon-cycle feedbacks, or anthropogenic COβ. It tells you when the orbital clock makes glaciation possible; the actual ice-volume response depends on the full climate system. The ~76 % of LR04 variance beyond the 26-component fit lives in ice-volume dynamics, carbon cycle, and other internal feedbacks.
Anthropogenic COβ may delay the next natural glaciation by 50+ kyr in moderate-emission scenarios (Ganopolski et al. 2016Β ). This is not modelled here.
Methodology highlights
| Topic | Detail |
|---|---|
| Spectral methods | MTM (Thomson multitaper, NW=3, K=5); Lomb-Scargle for Cheng2016βs irregular sampling; Hinich bispectrum for the eccentricity-beat coupling test; multi-component OLS amplitude fit for the formula; single-component OLS scan for the integer-divisor spectrum |
| Rayleigh resolution | ΞP β PΒ²/T β the fundamental spectral-resolution constraint, independent of method. At T = 5,320 kyr (full LR04) and P = 41 kyr, ΞP β 0.32 kyr (Berger sub-peaks resolvable). At T = 1.2 Myr and P = 110 kyr, ΞP β 10 kyr (95k/100k/112k spectrally collinear). |
| Pre-registration | Every empirical test had its data source, method, parameters, and verdict rules locked in writing before any data analysis β so that methodology issues couldnβt be rationalised away in the verdict. |
| Bootstrap | 1,000 resampling iterations for amplitude confidence intervals. |
For the full methodology, reproducibility scripts, and per-test verdict tables, see Milankovitch Evidence (GitHub doc 17)Β .
Data sources
- LR04 stack β Lisiecki, L. E. & Raymo, M. E. (2005), Paleoceanography 20, PA1003. 57 globally distributed benthic δ¹βΈO records, 0β5,320 kyr BP. Used for primary spectral analysis.
- Cheng 2016 Asian speleothem composite β Cheng, H. et al. (2016), Nature 534, 640. U-Th-dated Asian Monsoon δ¹βΈO, 0β640 kyr BP. Used as the non-tuned chronology-bias control (places the dominant peak at the same FFT bin as orbitally-tuned LR04, refuting a systematic ~10 % dating offset).
See also
- Eigenfrequencies β the secular eigenmodes (g_j, s_j) used throughout this page, with the model-frame vs Berger-attribution comparison
- Obliquity & Inclination β Earthβs two-component obliquity construction (axial tilt + inclination tilt)
- Eccentricity β Earthβs eccentricity oscillation and the 100-kyr problem statement
- Why Earth Is Special β Earthβs reference frame duality
- Supporting Evidence β broader validation across multiple model domains
- Scientific Background β historical context for the 100-kyr problem
- Milankovitch Evidence on GitHub (doc 17)Β β full technical record including methodology, scripts, per-test verdicts
- Solar System Resonance Cycle period table on GitHub (doc 55)Β β the 8 planets Γ 6 cycle types reference table used for direct-match cross-referencing