Climate Forcing Attribution

import matplotlib.pyplot as plt
import numpy as np
from scipy import stats

np.random.seed(456)

# Climate data
co2 = np.linspace(280, 420, 70) + np.random.normal(0, 3, 70)
temp = -0.6 + 0.013 * (co2 - 280) + np.random.normal(0, 0.12, 70)

slope, intercept, r_value, p_value, std_err = stats.linregress(co2, temp)
x_fit = np.linspace(275, 425, 100)
y_fit = slope * x_fit + intercept

# Prediction interval
n = len(co2)
se = np.sqrt(np.sum((temp - (slope * co2 + intercept))**2) / (n - 2))
ci_95 = 1.96 * se * np.sqrt(1 + 1/n + (x_fit - np.mean(co2))**2 / np.sum((co2 - np.mean(co2))**2))

fig, ax = plt.subplots(figsize=(10, 7), facecolor='#0a0a0f')
ax.set_facecolor('#0a0a0f')

# Prediction band
ax.fill_between(x_fit, y_fit - ci_95, y_fit + ci_95, color='#F54927', alpha=0.1, linewidth=0)
ax.fill_between(x_fit, y_fit - ci_95*0.5, y_fit + ci_95*0.5, color='#F54927', alpha=0.15, linewidth=0)

ax.plot(x_fit, y_fit, color='#F54927', linewidth=2.5, zorder=3)
ax.scatter(co2, temp, c='#27D3F5', s=60, alpha=0.8, edgecolors='white', linewidths=0.4, zorder=4)

# Zero line
ax.axhline(y=0, color='#444444', linestyle='-', linewidth=1, alpha=0.5)

# Equation
eq_text = f'ΔT = {slope:.4f} × CO₂ {intercept:+.2f}\nR² = {r_value**2:.3f}'
ax.text(0.03, 0.97, eq_text, transform=ax.transAxes, fontsize=10,
        color='#aaaaaa', ha='left', va='top', family='monospace')

for spine in ['top', 'right']:
    ax.spines[spine].set_visible(False)
for spine in ['bottom', 'left']:
    ax.spines[spine].set_color('#333333')

ax.set_xlabel('CO₂ Concentration (ppm)', fontsize=12, color='white', fontweight='500', labelpad=10)
ax.set_ylabel('Temperature Anomaly (°C)', fontsize=12, color='white', fontweight='500', labelpad=10)
ax.set_title('Greenhouse Gas Forcing', fontsize=15, color='white', fontweight='bold', pad=20, loc='left')
ax.tick_params(colors='#666666', labelsize=10, length=0)

plt.tight_layout()
plt.show()