This dissertation presents the results of statistical analyses of large climate datasets from two time intervals – the 20th century instrumental record and the proxy record of the last deglaciation – in order to understand the forcings and mechanisms of past climate variability. A longstanding question in climate dynamics concerns the origin of Pacific decadal variability (PDV). This issue is addressed by calculating the Southern Hemisphere equivalent of the Pacific Decadal Oscillation (PDO) index from Pacific sea surface temperature (SST) anomalies over the 20th century, which is found to be similar to its Northern Hemisphere counterpart. The Northern and Southern PDO indices both exhibit pronounced seasonality in autocorrelation with interannual persistence of winter SST anomalies despite their absence during the intervening summer, suggesting a role for reemergence. These two indices can be reasonably well reproduced using a first-order autoregressive model forced by the El Niño-Southern Oscillation, suggesting that Pacific decadal variability on both sides of the equator may largely be a reddened response to tropical forcing. Over 100 high-resolution proxy climate time series spanning parts or all of the Last Glacial Maximum (LGM) to Holocene interval are analyzed to characterize spatiotemporal patterns of glacial-interglacial climate change. Peak glacial and interglacial conditions occurred statistically synchronously between the Northern (22.1 ± 4.3 ka and 8.0 ± 3.2 ka) and Southern (22.3 ± 3.6 ka and 7.4 ± 3.7 ka) Hemispheres, suggesting that the hemispheres were synchronized by greenhouse gases, local insolation, and/or ocean circulation. Global cooling during the LGM was likely ≥4.9°C, but only ~0.6°C during the Younger Dryas. Younger Dryas climate anomalies exhibit a general hemispheric seesaw pattern with the largest negative anomalies in the high northern latitudes, mixed sign anomalies in the low latitudes, and modest positive anomalies in the high southern latitudes, consistent with an ocean circulation driver of this event. Empirical Orthogonal Function analysis of 71 records from 19-11 ka indicates that 72% of deglacial climate variability can be explained by two modes. The first mode (61% of the variance) shows a globally near-uniform signal and its associated principal component is strongly correlated with ice-core records of atmospheric CO₂. The second mode (11% of the variance) displays a more variable spatial pattern and its principal component parallels variations in Atlantic Meridional Overturning Circulation (AMOC) strength. Averaging 77 calibrated proxy temperature records over the last deglaciation indicates that global mean temperature was highly correlated and varied nearly in phase with the rise in CO₂, which differs from ice-core studies suggesting Antarctic temperature led CO₂. This result thus suggests a primary role for CO₂ in driving deglacial warming and the global deglaciation. Northern and Southern Hemisphere mean temperature time series both bear the imprint of CO₂ forcing, but differ largely in response to the millennial-scale seesawing of heat between the hemispheres related to AMOC variability. An analysis of deglacial temperature variability by latitude indicates that deglacial warming began near-synchronously throughout the Southern Hemisphere and tropics at ~19 ka and was coeval with abrupt cooling in the northern extratropics. These worldwide temperature changes may have been driven by a collapse of the AMOC related to the 19-kyr meltwater pulse. As this meltwater event has been tied to initial melting of Northern Hemisphere ice sheets in response to boreal summer insolation forcing, a classic Milankovitch trigger is implicated for the last deglaciation.