8.5 Paleoceanography

Paleoceanography reconstructs past ocean conditions from geochemical proxies preserved in marine sediments, corals, and ice cores. It reveals how the ocean responded to past climate perturbations driven by orbital forcing, tectonic changes, and carbon cycle events, providing essential context for understanding future climate trajectories.

Geochemical Proxy Records

Marine sediment cores contain microfossils and geochemical signals that record past ocean conditions. The primary proxies are:

δ¹&sup8;O from Benthic Foraminifera

The oxygen isotope ratio in calcite shells of bottom-dwelling foraminifera records both deep-water temperature and global ice volume:

$$\delta^{18}O = \left(\frac{(^{18}O/^{16}O)_{\text{sample}}}{(^{18}O/^{16}O)_{\text{VPDB}}} - 1\right) \times 1000 \text{ (‰)}$$

Higher δ¹&sup8;O indicates colder temperatures and/or larger ice sheets (preferential evaporation of lighter ¹&sup6;O). The LR04 benthic stack provides the standard reference curve for the past 5.3 Myr.

Mg/Ca Paleothermometry

The magnesium-to-calcium ratio in foraminiferal calcite is temperature-dependent, following an exponential relationship:

$$\text{Mg/Ca} = B \cdot \exp(A \cdot T)$$

Typically A ≈ 0.09°C²¹, B ≈ 0.38 mmol/mol for planktonic species. Combined with δ¹&sup8;O, Mg/Ca allows separation of temperature and ice volume signals.

U'K'37 Alkenone SST

Alkenones are long-chain unsaturated ketones produced by coccolithophores. Their degree of unsaturation correlates with growth temperature:

$$U^{K'}_{37} = \frac{[C_{37:2}]}{[C_{37:2}] + [C_{37:3}]}, \quad T = \frac{U^{K'}_{37} - 0.044}{0.033}$$

Provides SST estimates independent of carbonate chemistry. Valid range approximately 1–28°C. Calibrated against modern core-top samples.

Milankovitch Orbital Cycles

Variations in Earth's orbital parameters modulate the distribution of solar radiation, pacing the glacial-interglacial cycles of the Quaternary:

Eccentricity (~100 kyr & ~400 kyr)

Orbit shape varies from nearly circular (e ≈ 0.005) to slightly elliptical (e ≈ 0.058). Modulates the amplitude of precession effects. Dominates the δ¹&sup8;O record for the past 1 Myr.

Obliquity (~41 kyr)

Axial tilt varies between 22.1° and 24.5° (currently 23.44°). Controls seasonality and high-latitude insolation. Dominated ice age cycles before 1 Ma (41 kyr world).

Precession (~23 kyr)

Wobble of Earth's axis determines which hemisphere has summer at perihelion. Controls the seasonal distribution of insolation. Climatic precession parameter: $e\sin(\varpi)$.

$$Q_{\text{summer 65°N}} \approx Q_0\left[1 + e\sin\varpi + \epsilon\cos\varpi + \frac{\Delta\varepsilon}{\varepsilon_0}\right]$$

Summer insolation at 65°N is the canonical pacemaker of ice ages (Hays, Imbrie & Shackleton, 1976)

Marine Isotope Stages and Cenozoic Climate

The benthic δ¹&sup8;O record defines numbered Marine Isotope Stages (MIS). Odd numbers denote interglacials (warm), even numbers denote glacials (cold):

MIS 1

Current interglacial (Holocene, 0–11.7 ka). δ¹&sup8;O ≈ 3.2‰.

MIS 2

Last Glacial Maximum (~20 ka). δ¹&sup8;O ≈ 5.0‰. Sea level −120 m.

MIS 5e

Last interglacial (~125 ka). ~1–2°C warmer than pre-industrial. Sea level +6–9 m higher.

MIS 11

Interglacial ~400 ka. Long warm period. Orbital analog to modern Holocene.

Over the longer Cenozoic timescale (66 Ma), the deep ocean has cooled from ~12°C in the early Eocene to ~2°C today, with major transitions including the Eocene-Oligocene glaciation (34 Ma) and the Mid-Pleistocene Transition (~1 Ma, switch from 41 kyr to 100 kyr cyclicity).

Mid-Pleistocene Transition (MPT)

Around 1.2–0.7 Ma, ice age cycles switched from 41 kyr (obliquity-paced) to ~100 kyr (eccentricity-paced) without any change in orbital forcing. Proposed mechanisms include progressive CO' decline, Antarctic ice sheet expansion, and regolith removal exposing crystalline bedrock in North America and Scandinavia.

Terminations

Rapid deglaciations (Terminations I–V) are among the most dramatic climate transitions. Termination I (19–11 ka) raised sea level by ~120 m and increased CO' from 190 to 270 ppm. Meltwater pulses (MWP-1A: ~14.5 ka, ~20 m SLR in <500 years) caused rapid AMOC disruptions.

Paleocene-Eocene Thermal Maximum (PETM)

The PETM (~56 Ma) was a rapid global warming event caused by massive carbon release over ~10,000–20,000 years. It serves as the best geological analog for modern anthropogenic CO' emissions:

Characteristics

Global warming of +5–8°C. Carbon isotope excursion: δ¹³C dropped by −3 to −4‰ (light carbon input). Duration: ~170,000 years total, with peak warmth lasting ~100,000 years. Estimated 3,000–10,000 GtC released.

Ocean Impacts

Benthic foraminifera mass extinction (30–50% of species). Ocean acidification (CCD shoaled by 2 km). Deep ocean temperatures reached ~20°C. Reduced ocean overturning. Recovery took >100,000 years via silicate weathering feedback.

Glacial-Interglacial CO' Changes

Ice core records from Antarctica show that atmospheric CO' varied between ~180 ppm (glacials) and ~280 ppm (interglacials) over the past 800,000 years, tightly coupled with temperature and sea level. The ocean is the primary CO' reservoir mediating these changes:

$$\Delta pCO_2 \approx -10 \frac{\text{ppm}}{°C} \cdot \Delta T_{\text{ocean}} + \Delta pCO_{2,\text{biology}} + \Delta pCO_{2,\text{circulation}}$$

Mechanisms for glacial CO' drawdown include: increased ocean solubility (colder water absorbs more CO'), increased biological pump efficiency (iron fertilization from dust), reduced ocean ventilation (less deep water formation), and expanded sea ice (reduced gas exchange). No single mechanism accounts for the full 80–100 ppm change.

Python: Milankovitch Cycles, δ¹&sup8;O Analysis & Power Spectrum

Python

!/usr/bin/env python3

script.py104 lines

#!/usr/bin/env python3
"""paleoceanography.py - Orbital parameters, d18O, spectral analysis"""
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import periodogram

# --- 1. Milankovitch orbital parameters (simplified Berger 1978) ---
fig, axes = plt.subplots(2, 2, figsize=(14, 10))

t_kyr = np.linspace(-800, 0, 8000)  # 800 kyr BP to present

# Simplified orbital elements
eccentricity = 0.028 + 0.012 * np.cos(2*np.pi*t_kyr/413) + \
               0.008 * np.cos(2*np.pi*t_kyr/95) + \
               0.006 * np.cos(2*np.pi*t_kyr/124)

obliquity = 23.3 + 0.7 * np.cos(2*np.pi*t_kyr/41) + \
            0.3 * np.cos(2*np.pi*t_kyr/39.7)

precession = eccentricity * np.sin(2*np.pi*t_kyr/23 + \
             0.3*np.cos(2*np.pi*t_kyr/19))

# Summer insolation at 65N (simplified)
Q_65N = 480 + 20*np.cos(2*np.pi*t_kyr/41) + \
        30*precession/np.max(np.abs(precession))

ax1 = axes[0, 0]
ax1.plot(t_kyr, eccentricity, 'r-', lw=1, label='Eccentricity')
ax1.set_ylabel('Eccentricity'); ax1.set_title('Orbital Parameters')
ax1t = ax1.twinx()
ax1t.plot(t_kyr, obliquity, 'b-', lw=1, label='Obliquity')
ax1t.set_ylabel('Obliquity (deg)', color='blue')
ax1.set_xlabel('Time (kyr BP)'); ax1.grid(True, alpha=0.3)
ax1.legend(loc='upper left', fontsize=8)
ax1t.legend(loc='upper right', fontsize=8)

# --- 2. Synthetic benthic d18O record ---
ax2 = axes[0, 1]
# Ice volume response (lagged, nonlinear)
ice_response = np.zeros_like(t_kyr)
for i in range(1, len(t_kyr)):
    dt_kyr = t_kyr[1] - t_kyr[0]
    forcing = -Q_65N[i] + 480
    # Slow growth, fast deglaciation (sawtooth)
    if forcing > 0:
        ice_response[i] = ice_response[i-1] + 0.003 * forcing * dt_kyr
    else:
        ice_response[i] = ice_response[i-1] + 0.008 * forcing * dt_kyr
    ice_response[i] = np.clip(ice_response[i], -1.5, 2.0)

d18O = 3.3 + 0.8 * ice_response + 0.2 * np.random.randn(len(t_kyr))

ax2.plot(t_kyr, d18O, 'b-', lw=0.5)
ax2.fill_between(t_kyr, d18O, 3.5,
                 where=d18O > 3.5, alpha=0.3, color='lightblue',
                 label='Glacial')
ax2.fill_between(t_kyr, d18O, 3.5,
                 where=d18O < 3.5, alpha=0.3, color='orange',
                 label='Interglacial')
ax2.invert_yaxis()
ax2.set_xlabel('Time (kyr BP)'); ax2.set_ylabel('d18O (permil)')
ax2.set_title('Synthetic Benthic d18O Record')
ax2.legend(fontsize=8); ax2.grid(True, alpha=0.3)

# --- 3. Power spectrum of d18O ---
ax3 = axes[1, 0]
# Compute periodogram
fs = 1.0 / abs(t_kyr[1] - t_kyr[0])  # samples per kyr
freqs, power = periodogram(d18O, fs=fs)
periods = 1.0 / freqs[1:]
power = power[1:]

ax3.semilogy(periods, power, 'b-', lw=1)
# Mark Milankovitch periods
for p, label, color in [(100, '100 kyr', 'red'),
                         (41, '41 kyr', 'green'),
                         (23, '23 kyr', 'orange')]:
    ax3.axvline(p, color=color, ls='--', alpha=0.7, label=label)
ax3.set_xlabel('Period (kyr)'); ax3.set_ylabel('Power')
ax3.set_title('Spectral Analysis of d18O')
ax3.set_xlim(10, 500); ax3.legend(); ax3.grid(True, alpha=0.3)

# --- 4. CO2-temperature relationship ---
ax4 = axes[1, 1]
# Simplified ice-core CO2 vs temperature
co2_glacial = 180 + 15*np.random.randn(20)
co2_interglacial = 280 + 10*np.random.randn(20)
dT_glacial = -5 + 0.8*np.random.randn(20)
dT_interglacial = 1 + 0.5*np.random.randn(20)

ax4.scatter(co2_glacial, dT_glacial, c='blue', s=50, alpha=0.7,
            label='Glacial periods')
ax4.scatter(co2_interglacial, dT_interglacial, c='red', s=50,
            alpha=0.7, label='Interglacial periods')
# Modern
ax4.scatter([420], [1.2], c='black', s=200, marker='*',
            zorder=5, label='Modern (2024)')
ax4.set_xlabel('CO2 (ppm)'); ax4.set_ylabel('Temperature anomaly (C)')
ax4.set_title('CO2-Temperature Relationship (800 kyr)')
ax4.legend(); ax4.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('paleoceanography.png', dpi=150)
plt.show()

Click Run to execute the Python code

Code will be executed with Python 3 on the server

Fortran: Milankovitch Orbital Elements (Berger 1978 Solution)

Fortran

Compute Milankovitch orbital elements using simplified

program.f9078 lines

program milankovitch_orbits
  ! Compute Milankovitch orbital elements using simplified
  ! Berger (1978) trigonometric series
  implicit none
  integer, parameter :: dp = selected_real_kind(15)
  real(dp), parameter :: pi = 3.141592653589793_dp
  real(dp), parameter :: deg2rad = pi / 180.0_dp

real(dp) :: t_kyr, ecc, obl, prec, omega_bar
  real(dp) :: Q_65N, S0, lat_rad
  real(dp) :: ecc_terms(5,3), obl_terms(5,3), prec_terms(5,3)
  integer :: i, j, nsteps

! Simplified Berger (1978) series coefficients
  ! Format: amplitude, period(kyr), phase(deg)
  ! Eccentricity
  ecc_terms(1,:) = (/ 0.01102_dp, 413.0_dp, 148.0_dp /)
  ecc_terms(2,:) = (/ 0.00873_dp,  95.0_dp, 292.0_dp /)
  ecc_terms(3,:) = (/ 0.00749_dp, 124.0_dp,  53.0_dp /)
  ecc_terms(4,:) = (/ 0.00536_dp,  99.6_dp, 151.0_dp /)
  ecc_terms(5,:) = (/ 0.00416_dp, 131.0_dp, 280.0_dp /)

! Obliquity (degrees)
  obl_terms(1,:) = (/ 0.5983_dp, 41.00_dp, 251.9_dp /)
  obl_terms(2,:) = (/ 0.0637_dp, 39.73_dp, 280.4_dp /)
  obl_terms(3,:) = (/ 0.0505_dp, 40.40_dp, 128.3_dp /)
  obl_terms(4,:) = (/ 0.0413_dp, 53.61_dp, 292.7_dp /)
  obl_terms(5,:) = (/ 0.0235_dp, 54.68_dp, 116.2_dp /)

! Precession (e*sin(omega_bar))
  prec_terms(1,:) = (/ 0.01860_dp, 23.72_dp, 165.2_dp /)
  prec_terms(2,:) = (/ 0.01558_dp, 22.43_dp, 190.1_dp /)
  prec_terms(3,:) = (/ 0.01402_dp, 24.02_dp,  94.5_dp /)
  prec_terms(4,:) = (/ 0.00760_dp, 18.99_dp, 323.0_dp /)
  prec_terms(5,:) = (/ 0.00714_dp, 19.16_dp, 143.0_dp /)

S0 = 1361.0_dp          ! Solar constant (W/m^2)
  lat_rad = 65.0_dp * deg2rad  ! 65 degrees N

nsteps = 800  ! 800 kyr BP to present
  print *, '=== Milankovitch Orbital Elements ==='
  print *, 'Time(kyr)  Eccentricity  Obliquity(deg)  Precession  Q_65N(W/m2)'

do i = 0, nsteps, 10
    t_kyr = -real(i, dp)

! Compute eccentricity
    ecc = 0.0280_dp  ! mean value
    do j = 1, 5
      ecc = ecc + ecc_terms(j,1) * &
            cos(2.0_dp*pi*t_kyr/ecc_terms(j,2) + ecc_terms(j,3)*deg2rad)
    end do
    ecc = max(ecc, 0.0_dp)

! Compute obliquity
    obl = 23.320_dp  ! mean value (degrees)
    do j = 1, 5
      obl = obl + obl_terms(j,1) * &
            cos(2.0_dp*pi*t_kyr/obl_terms(j,2) + obl_terms(j,3)*deg2rad)
    end do

! Compute precession parameter (e*sin(omega_bar))
    prec = 0.0_dp
    do j = 1, 5
      prec = prec + prec_terms(j,1) * &
             sin(2.0_dp*pi*t_kyr/prec_terms(j,2) + prec_terms(j,3)*deg2rad)
    end do

! Summer solstice insolation at 65N (simplified)
    Q_65N = S0 / pi * (sin(lat_rad) * sin(obl*deg2rad) * pi + &
            cos(lat_rad) * cos(obl*deg2rad)) * &
            (1.0_dp + 2.0_dp * prec)

print '(F8.0, F14.6, F14.3, F14.6, F12.1)', &
          t_kyr, ecc, obl, prec, Q_65N
  end do

end program milankovitch_orbits

Click Run to execute the Fortran code

Code will be compiled with gfortran and executed on the server

Paleo Ocean Circulation and Gateway Events

Tectonic changes in ocean gateways have fundamentally altered global ocean circulation and climate. The opening and closing of seaways controlled the routing of heat and salt by ocean currents:

Drake Passage Opening (~34 Ma)

Opening of the Southern Ocean gateway between South America and Antarctica enabled the Antarctic Circumpolar Current (ACC). This thermally isolated Antarctica, allowing major ice sheet growth at the Eocene-Oligocene boundary. The $\delta^{18}O$ record shows a +1.5‰ shift marking this transition.

Closure of Panama (~3 Ma)

The Central American Seaway closure separated Atlantic and Pacific surface waters, strengthening the Gulf Stream and enhancing moisture supply to northern latitudes. This may have promoted Northern Hemisphere glaciation by increasing snowfall on Greenland and Scandinavia.

Closure of Tethys/Indonesian Throughflow

The progressive closure of the Tethys Ocean and restriction of the Indonesian Throughflow altered Indo-Pacific warm pool structure and monsoon circulation. The Messinian Salinity Crisis (5.96–5.33 Ma) desiccated the Mediterranean when its connection to the Atlantic was severed.

AMOC Variability

Paleo-proxy evidence ($\epsilon_{Nd}$, $^{231}$Pa/$^{230}$Th, $\delta^{13}C$) shows that Atlantic Meridional Overturning Circulation (AMOC) weakened or collapsed during Heinrich events and the Younger Dryas (12.9–11.7 ka), causing abrupt cooling in the North Atlantic.

Dansgaard-Oeschger Events

During the last glacial period (100–15 ka), Greenland ice cores record 25 abrupt warming events (D-O events) of 8–16°C within decades, linked to rapid AMOC resumptions. The associated bipolar seesaw transported heat between hemispheres, with Southern Ocean warming during Northern Hemisphere stadials.

Paleoceanographic Methods and Age Control

Reconstructing past ocean conditions requires precise age models for sediment cores, obtained through multiple dating methods:

$$\text{Age} = \frac{-1}{\lambda}\ln\left(\frac{^{14}C_{\text{sample}}}{^{14}C_{\text{modern}}}\right) \quad \text{(radiocarbon, valid to ~50 ka)}$$

Magnetostratigraphy

Geomagnetic polarity reversals recorded in sediment provide absolute age tie-points. The Brunhes-Matuyama reversal at 780 ka is a widely used marker. Particularly useful for Cenozoic timescales.

Biostratigraphy

First/last appearance datums of microfossil species (foraminifera, nannofossils, diatoms, radiolaria) provide chronological markers. Orbital tuning of cyclical climate proxies to astronomical solutions refines age models to ~5 kyr precision.

Tephra Layers

Volcanic ash layers provide isochronous (time-equivalent) markers across wide regions. Geochemical fingerprinting links tephra in different cores. Important for high-resolution correlations in the late Quaternary.

Core Logging

X-ray fluorescence (XRF) core scanning provides continuous elemental records at mm-scale resolution. Fe/Ca, Ti/Ca, and Ba/Al ratios trace terrigenous input, productivity, and redox conditions without destroying the core.

Key Concepts Summary

Proxy Records

δ¹&sup8;O records ice volume + temperature; Mg/Ca gives independent SST; U'K'37 alkenones provide organic-based SST estimates. Combined proxies separate multiple signals.

Milankovitch Cycles

Eccentricity (100/400 kyr), obliquity (41 kyr), and precession (23 kyr) pace ice ages. 65°N summer insolation is the primary trigger. 100 kyr cycle dominant for past 1 Myr.

PETM

Best analog for modern CO' release: +5–8°C warming, ocean acidification, benthic extinction. Recovery took >100 kyr. Current CO' release rate is ~10× faster than PETM.

Glacial CO' Cycles

CO' varied 180–280 ppm over ice ages. Ocean solubility, biological pump, and circulation changes controlled the partitioning. Now at >420 ppm (unprecedented in 3+ Myr).

← 8.4 Sea Level Rise Part 9: Marine Resources →