1. What is PyCDFt?
PyCDFt is a small, free, open-source Python library that does one job: it corrects and refines climate-model data so it matches real-world observations more closely. The technique is CDF-t (Cumulative Distribution Function transform), and PyCDFt is a fast, clean implementation of it.
Two things make it practical:
- Works on grids. Correct one time series, or a whole field of many locations, in a single call.
- Built for scale. Vectorized and chunked, so thousands of grid points × decades of data fit on a laptop or a supercomputer.
Think of PyCDFt as a translator between two languages: the one spoken by a climate model (coarse, and a bit off) and the one spoken by real observations (fine, and trustworthy). PyCDFt learns the "dictionary" from a past period, then uses it to rewrite the model's data for another period — typically the future — in a way that is both realistic and consistent with the model's projected change.
The library has two layers:
| Layer | What it does | Where it lives |
|---|---|---|
| Core | The CDF-t math in numpy. Always installed; the stable, tested heart. | pycdft |
| GPU backend | Optional. The same math in torch, which can be faster on large grids (see §8). | pycdft.gpu |
This page is the plain-language companion to README.md (quick
start, API, comparison tables). Mathematical details and design choices
are in docs/design.md; the GPU backend is in docs/gpu.md.
2. The problem PyCDFt solves
Global climate models (GCMs) simulate the whole planet, but they are:
- Coarse — a grid cell may be ~100 km wide, too big for local decisions (a farm, a city, a solar plant).
- Biased — they often run systematically too warm, too dry or too cloudy compared with measurements.
Two techniques address this:
- Bias correction: remove the model's systematic error by comparing it with observations over a shared calibration period.
- Downscaling: turn coarse model data into fine, location-specific data.
CDF-t does both at once, and it also preserves how the model's climate changes between the calibration period and the target period.
A simple method says "the model is 2 °C too warm, so subtract 2 °C." But the bias itself can differ across the distribution, and the climate shifts over time. CDF-t carries the model's evolution (how its distribution shifts and widens) onto the observation scale. A mean shift (§5.3) is used only as a first step; the quantile transfer handles the rest.
2.1 The inputs, in standard vocabulary
The API takes exactly three inputs and returns one output:
| PyCDFt name | Standard term | What it is |
|---|---|---|
obs | observations, calibration period (HR) | High-resolution observations — the reference the model is corrected against. |
cm_hist | model, calibration period (LR) | Coarse model data for the same calibration period as obs. |
cm_future | model, target period (LR) | Coarse model data you want to refine. The name says "future", but it can be any target period, including a held-out period for validation. |
output of transform(cm_future) | downscaled output (HR) | The corrected, high-resolution estimate for the target period. |
The model (LR) is the input; the observations (HR) are the reference. You calibrate on a past model/observation pair, then apply the mapping to the target-period model data. When you work with a field of locations, each location pairs its observations with the model series chosen for it. To validate, hold out a later period, downscale it, and compare with that period's observations — the library only does the mapping, the evaluation is up to you.
3. Installation
You need Python 3.9 or newer. Install from the repository:
pip install . # core only (numpy)
pip install -e ".[test]" # development install with the test tools(The quotes keep shells such as zsh from interpreting the brackets.) The core has no dependencies beyond numpy. Optional extras:
| Extra | What it adds |
|---|---|
[gpu] | torch, for the GPU backend (pycdft.gpu). For a specific build (ROCm, CPU-only) install that torch first — see §8. |
[examples] | Matplotlib, only to run the example scripts and save figures. |
[cross] | ibicus, only for the optional cross-validation tests. |
[test] | pytest, for running the test suite. |
[dev] | pytest + pyflakes, for development (tests and lint). |
4. Quick start
The whole library in a few lines (time is the first axis):
import numpy as np
from pycdft import CDFt
rng = np.random.default_rng(7)
obs = rng.normal(0.0, 1.0, 1000) # observations (calibration period)
cm_hist = rng.normal(3.0, 1.0, 1000) # model (calibration period)
cm_future = rng.normal(4.0, 1.0, 1000) # model (target period)
model = CDFt(n_quantiles=1001, tail="linear", delta_shift="additive")
ds = model.fit_transform(obs, cm_hist, cm_future) # downscaled output
The 2-D field form (each column an independent location), the one-shot
cdft_transform, the fit/transform
split and the histogram path are in the README quick start.
5. How CDF-t works
CDF-t was published by Michelangeli, Vrac & Loukos in 2009. It corrects the model's target-period distribution so that it looks like the observations would have looked had they experienced the same climate change as the model.
5.1 The building blocks: CDF and quantile function
A CDF, written \(F(x)\), answers: "what fraction of the values in a dataset are smaller than or equal to \(x\)?" It is a number between 0 and 1 that tells you how far up the ranking \(x\) sits.
- \(F(x) = 0\): \(x\) lies below every value in the dataset.
- \(F(x) = 0.5\): \(x\) is the median.
- \(F(x) = 1\): \(x\) is at or above the largest value.
The inverse of a CDF is the quantile function \(q(p)\): given a fraction \(p\) (say 0.5), it returns the value at that position in the ranking. Notation, as in the R CDFt package and ibicus:
| Symbol | Meaning |
|---|---|
| \(F_m^{cal}\), \(q_m^{cal}\) | model (LR) CDF and quantile function, calibration period — this is cm_hist. |
| \(F_o^{cal}\), \(q_o^{cal}\) | observation (HR) CDF and quantile function, calibration period — this is obs. |
| \(F_m^{tar}\), \(q_m^{tar}\) | model (LR) CDF and quantile function, target period — this is cm_future. |
| \(F_o^{tar}\) | target-period observation CDF — the HR distribution we are after (the output). |
| \(m_{cal}\), \(o_{cal}\) | model and observation series in the calibration period (used in the delta shift) |
Superscripts: \(cal\) = calibration period, \(tar\) = target period. Subscripts: \(o\) = observations (HR), \(m\) = model (LR). Below, \((F_m^{cal})^{-1}\) and \(q_m^{cal}\) mean the same thing. Model series are taken after the delta shift of §5.3.
5.2 The CDF-t formulas
CDF-t estimates the target observation CDF from three known ones:
\[ F_o^{tar} = F_o^{cal} \circ (F_m^{cal})^{-1} \circ F_m^{tar} \]Read right to left: where does a value rank in the target-period model (\(F_m^{tar}\)), which calibration-period model value has that rank (\((F_m^{cal})^{-1}\)), and how often does a calibration-period observation fall below that value (\(F_o^{cal}\))?
The downscaled value for a target-period model value \(x\) is found by inverting that CDF at the rank of \(x\):
\[ y(x) = q_m^{tar}\Big( F_m^{cal}\big( q_o^{cal}\big( F_m^{tar}(x) \big) \big) \Big) \]- \(p = F_m^{tar}(x)\) — where does \(x\) rank in the target-period model? A fraction from 0 to 1.
- \(q_o^{cal}(p)\) — which calibration-period observation has that rank? An observation value.
- \(F_m^{cal}(\cdot)\) — where does that value rank among the calibration-period model values? A fraction.
- \(q_m^{tar}(\cdot)\) — which target-period model value has that rank? That is the downscaled output \(y\).
Steps 2 and 3 are the calibration part: they compare observations with the model over the calibration period, and this is where the bias is corrected. Steps 1 and 4 use the target-period model distribution, which is how the model's change is carried through. If the model had no bias (steps 2–3 do nothing), the output equals the input.
A useful check: if the target-period model distribution equals the
calibration one (\(F_m^{tar} = F_m^{cal}\)), the formula collapses to
classical quantile mapping, \(y = q_o^{cal}(F_m^{cal}(x))\). So when the
model does not change, CDF-t reduces to ordinary quantile mapping
(verified in tests/test_core.py).
5.3 The delta shift
Before the transfer, PyCDFt (like the reference implementations) can re-centre the model so its average matches the observations' average. The same shift, computed from the calibration period, is applied to both the model's calibration and target series. Three modes:
\[ \begin{aligned} \text{additive:} & \quad x \to x + \big( \mathrm{mean}(o_{cal}) - \mathrm{mean}(m_{cal}) \big) \\ \text{multiplicative:} & \quad x \to x \times \frac{\mathrm{mean}(o_{cal})}{\mathrm{mean}(m_{cal})} \\ \text{none:} & \quad \text{(no change)} \end{aligned} \]The model usually has a persistent average error relative to the calibration-period observations. "Additive" adds or subtracts a constant to remove it. "Multiplicative" scales by a ratio — better for quantities that must stay positive, such as rainfall or sunlight. "None" shifts nothing. Because the same shift is applied to both periods, the average bias is removed but the model's change between them is kept.
5.4 How PyCDFt evaluates it (the grid)
Rather than running the formula for every single value (slow), PyCDFt evaluates it on a fixed set of probabilities \(p_k = \mathrm{linspace}(0, 1, Q)\) (default \(Q = 1001\)) and then interpolates. Three quantities are computed:
\[ \begin{aligned} R_k &= F_m^{cal}\big( q_o^{cal}(p_k) \big) && \text{(calibration only — independent of the target period)} \\ x_{\text{grid},k} &= q_m^{tar}(p_k) && \text{(target period)} \\ T_k &= q_m^{tar}(R_k) && \text{(target period)} \end{aligned} \]The output for any target-period model value \(x\) is then a piecewise-linear interpolation of \(T\) against \(x_{\text{grid}}\):
\[ y(x) = \mathrm{interp}(x,\, x_{\text{grid}},\, T) \]Building a lookup table once is much faster than recomputing everything per value. And since \(R\) depends only on the calibration period, you can fit it once and reuse it for many target periods, scenarios or ensemble members.
6. Core concepts you should know
6.1 Tails — what happens beyond the data
Real data is finite, so the transfer curve has a start and an end. The question is what to do for values outside that range. PyCDFt offers two choices:
| Option | Behaviour | When to use |
|---|---|---|
tail="clip" | Hold the constant end value. | When the output must never leave the observed range. |
tail="linear" (default) |
Extend with a straight line whose slope is measured near the edge and capped. | When extremes matter and you want a smooth, sensible extension. |
In linear mode, PyCDFt fits a straight line to the last 30 points of the transfer, then limits its steepness to at most twice the median slope in the middle. The extrapolation stays realistic and bounded.
6.2 Missing data (NaN)
Missing values are handled in four cases:
- NaNs scattered in the calibration data are ignored (excluded from counts and averages).
- NaNs in the target-period input pass through as NaN, at those positions only.
- If one location's calibration data is entirely NaN, that location maps to NaN instead of crashing the run.
- An entirely NaN target field gives an entirely NaN output.
For comparison, ibicus does not support missing data and the R
CDFt package drops NaNs. In PyCDFt the output stays aligned
with the input and one bad location does not affect the others.
6.3 Memory-friendly design
- Calibration data keeps its original type (float32 stays float32) to save memory; internal math uses float64 for accuracy.
- Calculations run in small chunks of grid cells so intermediate results stay under a memory limit (default 2 GiB).
- The caller's arrays are never modified.
- For very long records, histogram mode (
use_histograms=True) builds a compact summary instead of keeping every value, trading a little accuracy for much less memory.
7. API reference
7.1 Entry points
| Name | Use |
|---|---|
CDFt(...) | High-level class with .fit(obs, cm_hist), .transform(cm_future) and .fit_transform(obs, cm_hist, cm_future). |
cdft_transform(obs, cm_hist, cm_future, ...) | One-shot functional form. |
calibrate, apply_cdft | Low-level building blocks (below). |
7.2 CDFt options
| Parameter | Default | Meaning |
|---|---|---|
n_quantiles | 1001 | Number of probability points \(Q\) on the transfer grid (§5.4). |
tail | "linear" | Handling of values beyond the data range: "linear" or "clip" (§6.1). |
delta_shift | "additive" | How to remove the model's average bias: "additive", "multiplicative" or "none" (§5.3). |
use_histograms | False | Compact histogram path for very long records (§6.3). |
n_bins | 2000 | Number of bins in histogram mode. |
workers | 1 | Threads for parallel per-location work. |
7.3 Low-level building blocks
For custom pipelines, pycdft.core exposes
calibrate/apply_cdft and pycdft.math
the individual kernels; see the README's "Low-level API" for a runnable
example.
calibrate takes the model first and the observations second,
the reverse of CDFt.fit(obs, cm_hist). Confirm against the
signature in your installed version before relying on it.
The calibration object is target-independent: fit it once, then reuse it for any number of target periods, scenarios or ensemble members.
8. The GPU backend
The optional pycdft.gpu backend does the same math
in torch, so it can run on a graphics card (NVIDIA CUDA, AMD ROCm or Apple
MPS) and can be faster on large grids. It mirrors the core interface and
adds a device argument:
from pycdft.gpu import CDFt
ds = CDFt(device="cuda").fit_transform(obs, cm_hist, cm_future)
Device names follow torch: "cuda" (also used by ROCm builds),
"mps" and "cpu". Install the torch build that
matches your hardware — see docs/gpu.md.
The numpy core never imports the GPU backend. Without torch,
pycdft.gpu is still importable but calls raise a clear error.
It was validated against the numpy backend on an AMD MI210 (ROCm),
matching to within floating-point precision. Validation on CUDA and MPS is
not reported here.
9. How PyCDFt compares to other tools
PyCDFt follows the conventions of the R CDFt package and the
Python package ibicus, so it should agree closely with
both in the middle of the data. The differences are mostly in the
extremes (tails) and the exact grid. The comparison table, the measured
agreement and the analysis of where the implementations differ are in
the README and in docs/design.md §4.