AI for Climate: Weather Prediction and Carbon Tracking

Every June, the India Meteorological Department's forecast of monsoon onset over Kerala moves markets, sowing decisions, and reservoir releases for a country where roughly 70% of annual rainfall arrives in a single four-month window. IMD's operational forecast comes from NCUM, the NCMRWF Unified Model run at the National Centre for Medium Range Weather Forecasting in Noida, a descendant of the UK Met Office's Unified Model, executed on a supercomputer at a global grid spacing of roughly 12 km, refined further for the Indian domain. In November 2023, Google DeepMind published a competing approach in Science: a neural network called GraphCast that was trained once on four decades of historical weather data and then, on a single Google Cloud TPU v4 chip, produced a 10-day global forecast in under a minute, beating the European Centre's flagship physics-based model (IFS HRES) on about 90% of 1,380 tested variable-and-lead-time combinations (Lam et al., 2023). Two systems, same target, radically different machinery: one integrates the equations of fluid dynamics forward in time on a supercomputer; the other learned an approximation of that integration from data and applies it almost for free. This chapter builds both machines from first principles, and then does the same for the sibling problem of tracking carbon and methane from orbit, where a very similar tension between physical measurement and learned inference shows up again.

The physics problem underneath the forecast

Weather is, formally, an initial value problem for a system of partial differential equations, the primitive equations, which describe how temperature, pressure, humidity, and the wind vector evolve over a rotating, stratified fluid (the atmosphere) under gravity, the Coriolis force, and pressure-gradient forces. There is no closed-form solution. Operational forecasting discretizes the atmosphere into a 3D grid, typically tens of vertical pressure levels stacked over a global latitude-longitude mesh, and steps the equations forward in small time increments using finite-difference or spectral methods. NCUM at ~12 km horizontal spacing with dozens of vertical levels has on the order of hundreds of millions of grid points (roughly 5-6 million horizontal points × ~70 vertical levels ≈ 390 million), each carrying several state variables, updated every few minutes of simulated time to produce a multi-day forecast.

That time-stepping is not free to choose. Explicit numerical schemes for wave-like PDEs are constrained by the Courant-Friedrichs-Lewy (CFL) condition: information cannot be allowed to physically cross more than roughly one grid cell per time step, or the scheme becomes numerically unstable and blows up. Formally, for a characteristic wave speed v, grid spacing Δx, and time step Δt, stability requires vΔt/Δx ≤ C for some Courant number C close to 1. Halve the grid spacing and, to keep the left-hand side from exceeding the same bound, you must halve the time step too. This single constraint is why every resolution upgrade a national weather service buys is also a multi-year supercomputer procurement.

Worked example: what doubling resolution actually costs

Suppose IMD halves its horizontal grid spacing, from 12 km to 6 km, keeping the number of vertical levels fixed. Trace the compute cost step by step.

Step 1, spatial grid points. Halving the spacing in both the east-west and north-south directions quadruples the number of horizontal grid points for the same physical domain: (Δx/2) → 4× the columns.

Step 2, time step. By the CFL condition above, halving Δx with the wave speed v unchanged forces Δt to halve as well, so covering the same forecast horizon (say, a 10-day run) now takes twice as many time steps.

Step 3, total work. Total computational work is proportional to (grid points touched per step) × (number of steps): 4 × 2 = 8. Doubling horizontal resolution costs roughly eight times the compute, not two, and this factor doesn't depend on whether you start from 12 km or 3 km — it's a property of the doubling itself, driven jointly by the geometry of a 2D grid and the CFL constraint on the time step. This is the standard back-of-envelope reasoning meteorological agencies use when they justify supercomputer upgrades, and it is exactly why an incremental resolution gain (say 12 km to 10 km) is cheap while a full doubling is a generational hardware decision.

Learning the time-evolution operator instead of solving it

GraphCast and its contemporaries, Huawei Cloud's Pangu-Weather (Bi et al., 2023, Nature) and NVIDIA's FourCastNet (Pathak et al., 2022), sidestep the CFL bottleneck entirely by not integrating any differential equation at inference time. Instead they learn a single function f that maps the atmospheric state at time t directly to the state at t + Δt, trained by supervised learning on ERA5, the European Centre's 40-year reanalysis dataset that reconstructs the historical atmosphere on a fine grid by blending observations with a physical model. GraphCast specifically trains on 0.25° resolution ERA5 data (about 28 km at the equator, roughly 721 × 1440 grid points) spanning 1979 to 2017, learning to predict 6-hour steps.

The architecture is an encoder-processor-decoder graph neural network. The encoder takes the regular latitude-longitude grid and embeds each grid point's variables (temperature, wind components, humidity, geopotential, at multiple pressure levels) onto the nodes of a separate mesh graph, typically a refined icosahedron, using a small per-node multilayer perceptron. The processor then runs several rounds of message passing on that mesh: each node updates its embedding by aggregating messages from its graph neighbors through a learned function. This is the part doing the actual work of approximating atmospheric dynamics — advection, pressure-gradient forcing, diffusion — as a purely local, learned interaction repeated many times, since real weather physics is itself local (a parcel of air is influenced directly only by its immediate surroundings; distant influence has to propagate step by step). The decoder maps the processed mesh embeddings back onto the original latitude-longitude grid, producing the predicted state 6 hours later. To get a 10-day forecast, the model runs autoregressively: the 6-hour output is fed back in as the next input, 40 times over (10 days × 24 hours ÷ 6 hours per step = 40 steps). The diagram below shows this loop.

Worked example: one message-passing round, traced by hand

To see what a "learned local physics" update actually computes, strip the processor down to four grid cells arranged in a square (northwest, northeast, southeast, southwest), each holding a single scalar (a toy temperature field), connected in a cycle: NW-NE, NE-SE, SE-SW, SW-NW. A minimal message-passing update sends a message along each edge proportional to the difference between neighboring values, then sums the incoming messages into each node, exactly the structure real GNN layers use, just with a single scalar weight instead of a learned neural network:

x = {"NW": 20.0, "NE": 22.0, "SW": 19.0, "SE": 21.0}
neighbors = {
    "NW": ["NE", "SW"],
    "NE": ["NW", "SE"],
    "SW": ["SE", "NW"],
    "SE": ["NE", "SW"],
}
w = 0.1  # shared edge weight

new_x = {}
for node in x:
    incoming = sum(w * (x[nbr] - x[node]) for nbr in neighbors[node])
    new_x[node] = x[node] + incoming

print(new_x)
print("total before:", sum(x.values()), " total after:", sum(new_x.values()))

Tracing NW by hand: its neighbors are NE (22.0) and SW (19.0), so its incoming message is 0.1×(22-20) + 0.1×(19-20) = 0.2 - 0.1 = 0.1, giving a new value of 20.1. Tracing NE similarly: neighbors NW (20.0) and SE (21.0) give 0.1×(20-22) + 0.1×(21-22) = -0.2 - 0.1 = -0.3, so NE drops to 21.7. The same arithmetic for SW and SE gives 19.3 and 20.9. Running the code prints {'NW': 20.1, 'NE': 21.7, 'SW': 19.3, 'SE': 20.9}, with total before: 82.0 total after: 82.0. The total is exactly conserved, because every edge contributes a message and its exact negative to the two nodes it connects, so the sum over the whole graph cancels regardless of the edge weight's value. This is not a coincidence and not something the toy network was told to do: it falls straight out of the update rule's structure, the same way real diffusion conserves heat. Production models like GraphCast use a learned vector-valued message function instead of a single scalar weight, and the messages needn't be symmetric or purely diffusive (real atmospheric transport also advects, which is directional), but the underlying idea, that each round of message passing is a local, learnable approximation of one increment of physical evolution, is exactly this mechanism scaled up across roughly 41,000 mesh nodes (GraphCast's processor runs on an icosahedron refined 6 times, giving 40,962 nodes) and 16 rounds of message passing per step.

Common misconception: satellites measure concentration, not emissions

A student who has seen satellite images of smog or algorithms that classify photographs will naturally assume that a "carbon-tracking satellite" watches a coal plant and reads off its emission rate the way a speed camera reads a car's velocity. It doesn't, and the distinction matters. What an instrument like NASA's OCO-2 or ESA's Sentinel-5P (carrying the TROPOMI spectrometer) measures is a concentration: the total column of CO₂ or CH₄ molecules that sunlight passed through on its way from the sun, down to the ground, and back up to the sensor, inferred from how much light was absorbed at wavelengths those gas molecules are known to absorb (near-infrared bands around 1.6 and 2.06 microns for CO₂; comparable bands for methane). That is a snapshot of how much gas is sitting in a column of air at one moment, not a rate of release. Emission is a flux, mass per unit time, and going from a concentration field to a flux requires a second step: an atmospheric transport model that accounts for wind speed and direction to work out how much gas would have had to be released, and how fast, to produce the observed pile-up downwind of a source. This inversion step is exactly analogous to the weather forecasting problem run in reverse: instead of pushing a known atmospheric state forward through the equations of motion, you push an unknown source term forward through the same kind of transport model and adjust it until the simulated concentration field matches what the satellite actually saw. Confusing "the satellite detected elevated methane over this facility" with "the satellite measured this facility's emission rate in tonnes per year" is the single most common error in reporting on this technology, and it is the gap that companies like GHGSat and the Climate TRACE coalition are explicitly built to close.

Worked example: scoring a pixel with a matched filter

Point-source detection, used by systems like NASA's EMIT imaging spectrometer and airborne surveys of methane plumes, works by scanning each pixel of a hyperspectral image (dozens to hundreds of narrow wavelength bands per pixel) for a spectral shape that matches methane's known absorption signature, against a noisy, spatially correlated background. The standard tool, a cluster-tuned matched filter, closely following the method used by Thorpe et al. (2013) for airborne methane mapping, scores each pixel by

mf(x) = t^T Σ⁻¹ (x - μ) / sqrt(t^T Σ⁻¹ t)

where x is the pixel's observed radiance vector, μ and Σ are the background's mean radiance and covariance across bands (estimated from thousands of "clean" pixels elsewhere in the scene), and t is the target's differential absorption signature, how much a unit increase in methane concentration would depress the radiance in each band. Reduce this to two bands, one strongly absorbed by methane and one clean reference band, so the whole computation can be checked by hand.

import numpy as np

Sigma = np.array([[4, 2],
                   [2, 3]])         # background covariance, 2 bands
mu    = np.array([50, 48])         # background mean radiance
t     = np.array([-0.8, -0.1])     # target absorption signature
x     = np.array([46, 47.5])       # observed pixel radiance

Sigma_inv   = np.linalg.inv(Sigma)
diff        = x - mu
numerator   = t @ Sigma_inv @ diff
denominator = np.sqrt(t @ Sigma_inv @ t)
mf_score    = numerator / denominator

print(f"numerator   = {numerator:.3f}")
print(f"denominator = {denominator:.4f}")
print(f"mf score    = {mf_score:.3f}")

Working it by hand: the determinant of Σ is 4×3 - 2×2 = 8, so Σ⁻¹ = (1/8)×[[3,-2],[-2,4]] = [[0.375,-0.25],[-0.25,0.5]]. The deviation from background is x - μ = [-4, -0.5]. Multiplying, Σ⁻¹(x-μ) = [0.375×(-4) + (-0.25)×(-0.5), -0.25×(-4) + 0.5×(-0.5)] = [-1.375, 0.75]. Dotting with t gives the numerator: (-0.8)×(-1.375) + (-0.1)×0.75 = 1.1 - 0.075 = 1.025. For the denominator, Σ⁻¹t = [0.375×(-0.8)+(-0.25)×(-0.1), -0.25×(-0.8)+0.5×(-0.1)] = [-0.275, 0.15], and t·Σ⁻¹t = (-0.8)×(-0.275)+(-0.1)×0.15 = 0.22-0.015 = 0.205, whose square root is ≈0.4528. The code prints numerator = 1.025, denominator = 0.4528, mf score = 2.264. A score of about 2.26 background-noise units above what an uncorrelated pixel would show is well past the threshold operational systems use to flag a pixel for follow-up, exactly the kind of signal that triggers a human analyst or a downstream inversion pipeline to look closer at that facility.

From a flagged pixel to a tonnage estimate

Once a plume is delineated in an image, converting it to an emission rate uses the integrated mass enhancement (IME) method (Varon et al., 2018): sum the excess gas mass across every plume pixel to get IME (kilograms), then estimate Q ≈ IME × U_eff / L, where U_eff is an effective wind speed at plume height and L is the plume's length along the wind direction. This is the transport-model inversion from the previous section collapsed into a single closed-form estimate for a compact point source, and it is the step that turns "this facility's pixels look anomalous" into "this facility is estimated to be releasing approximately this many tonnes of methane per hour," the number that actually matters for regulation, carbon markets, and coalitions like Climate TRACE that aggregate such estimates across hundreds of millions of individual emitting assets worldwide.

Active recall

Attempt each question before reading its answer.

1. Why is it wrong to describe GraphCast as "an image classifier that looks at satellite photos and predicts the weather"?
2. If IMD upgrades again from 6 km to 3 km resolution, keeping vertical levels fixed, how many times more compute does the CFL-based reasoning predict, relative to the 6 km run?
3. In the matched filter worked example, suppose recalibration raises the assumed target signature's first component from -0.8 to -1.2 (leaving the second band at -0.1, and Σ, μ, x unchanged). Recompute the matched filter score. Does it rise in proportion to the stronger assumed signature?
4. In the message-passing worked example, change only the SW-NW edge weight from 0.1 to 0.3, keeping the other three edges at 0.1. Recompute all four updated node values. Is the total still conserved, and why?
5. Why are satellite-based "emissions maps" actually the output of two separate models chained together, and what would go wrong if you skipped the second one?
6. Forecast skill degrades with lead time even for a perfect model with a tiny error in the initial condition. Name the mathematical property responsible and the operational technique meteorological centers use to quantify the resulting uncertainty.

Answers.

1. An image classifier maps pixels to a label drawn from a fixed, small set (cat, dog, stop sign) using spatial texture cues, with no obligation to respect any physical law. GraphCast's inputs and outputs are the same physical state variables (temperature, wind components, humidity, geopotential at dozens of pressure levels) on the same grid; it is trained to predict the next physical state given the current one, i.e. to approximate the time-evolution operator of the atmosphere, and it is applied autoregressively exactly the way a numerical integrator is. It never sees a photograph and its output is a full 3D physical field, not a label.

2. Still 8×. The derivation (quadrupled grid points from halving spacing in two horizontal dimensions, doubled time steps from the CFL condition) never referenced the starting resolution, only the fact of a doubling, so the same factor applies at any starting point.

3. Recomputing with t = [-1.2, -0.1]: Σ⁻¹t changes because t changed, so both numerator and denominator change. The new numerator is t·Σ⁻¹(x-μ) = (-1.2)×(-1.375)+(-0.1)×0.75 = 1.65-0.075 = 1.575, and the new denominator, from t·Σ⁻¹t, works out to ≈0.6964, giving mf ≈ 1.575/0.6964 ≈ 2.262. The score barely moved (2.264 to 2.262) even though the assumed signature strength in the dominant band grew 50%. That's because the matched filter score is exactly invariant to a uniform rescaling of the entire target vector t (scaling numerator and denominator by the same factor cancels), and this change, though only applied to one component, was dominated by that already-dominant band, so it behaved almost like a uniform rescaling. The lesson: matched filter detection confidence depends on the target's spectral shape relative to background noise correlations, not on getting its absolute magnitude exactly right.

4. With only the SW-NW weight raised to 0.3: NW's incoming message is now 0.1×(22-20) + 0.3×(19-20) = 0.2-0.3 = -0.1, giving 19.9. NE is unaffected by the change (neighbors NW and SE, both still weight 0.1), so it still becomes 21.7. SE is also unaffected (neighbors NE, SW at weight 0.1), staying 20.9. SW's incoming message becomes 0.1×(21-19) + 0.3×(20-19) = 0.2+0.3 = 0.5, giving 19.5. New total: 19.9+21.7+19.5+20.9 = 82.0, still conserved. This holds because conservation only requires that the weight used for a message from i to j equal the weight used for the return message from j to i for that same edge; it never required all edges to share one weight. Any per-edge (but symmetric) weighting still conserves the total, which is why real GNN-based weather models can freely learn different message strengths for different neighbor pairs without breaking mass or energy consistency, as long as the update rule keeps that pairwise symmetry.

5. A concentration map is the output of a forward radiative-transfer / spectroscopy model (or a matched filter, for point-source detection) applied to raw satellite radiance. An emission map additionally requires a transport-inversion model (or the closed-form IME formula for a compact plume) that uses wind data to convert an observed pile-up of gas into a release rate. Skipping the second step and reporting concentration numbers as if they were emission rates would conflate, for example, a facility that emits modestly into locally still air (and so shows a large accumulated concentration) with one that emits heavily into a strong crosswind (which disperses the plume and can show a smaller concentration signal despite releasing more gas overall).

6. Sensitive dependence on initial conditions, the property identified by Edward Lorenz (1963) in his study of a simplified convection model, popularly the "butterfly effect": in a chaotic system, any nonzero error in the initial state grows over time, roughly exponentially, at a rate set by the system's leading Lyapunov exponent, which is why forecast skill decays with lead time even for a flawless model. Operational centers handle this with ensemble forecasting: running many forecasts (physical or, increasingly, neural, as GraphCast also supports) from slightly perturbed initial conditions and treating the spread across the ensemble as a quantitative estimate of forecast uncertainty, rather than trusting any single run's output as certain beyond a few days.

Think About It

Think about this: How would you explain ai for climate: weather prediction and carbon tracking to a friend who has never seen a computer? What real-world analogy would you use? Imagine you had to build a system using these concepts — what would be your first step? Try this: before moving on, write down three things you learned and one question you still have.

Practice Exercises

Now it is time to practice! Complete these challenges to solidify your understanding:

• Exercise 1: Write a short program that demonstrates the core concept from this chapter. Test it with at least 3 different inputs.
• Exercise 2: Find a real-world example where ai for climate: weather prediction and carbon tracking is used in an Indian company (like TCS, Infosys, Flipkart, or ISRO). Write a paragraph explaining the connection.
• Exercise 3: Create a mind-map connecting ai for climate: weather prediction and carbon tracking to at least 3 other topics you have studied.

Key Takeaways — Summary and Recap

Let us recap what we covered: the core ideas behind ai for climate: weather prediction and carbon tracking, how they connect to real-world applications, and why they matter for your journey in computer science. Remember these key points as you move forward. For competitive exam preparation (CBSE, JEE, BITSAT), focus on understanding the WHY behind each concept, not just the WHAT.