Attribute Correlation Modeling for Synthetic Spatial Data

A synthetic feature is more than a geometry. Every generated parcel, sensor, or trip carries a vector of attributes — income, elevation, temperature, dwell time, land-use class — and the value a synthetic dataset delivers depends almost entirely on whether the relationships between those attributes survive generation. This page is part of Spatial Distribution & Pattern Generation: where point process simulation models decide where features land, attribute correlation modeling decides what they carry and how those values move together. Sample each attribute from its own marginal in isolation and you get a table where every column is individually plausible and every row is nonsense — high income beside zero elevation beside tropical temperature at a polar latitude. The joint distribution, not the marginals, is what a downstream regressor, classifier, or graph model actually learned from real data, and it is the thing most naive generators silently destroy.

Problem Framing: The Independent-Column Fallacy

The failure this page solves is the independent-column dataset — a synthetic table where every attribute matches its real marginal distribution and yet the rows are jointly impossible. It is the single most common defect in attribute synthesis because the naive approach is so tempting: fit a distribution to each column, draw from each independently, done. The marginals validate perfectly. The dataset is still useless.

Consider a synthetic dataset of urban buildings with attributes for floor count, footprint area, assessed value, and year built. Real skyscrapers have many floors, large footprints, high value, and recent construction, all together. Sample each column independently and you produce single-story buildings with skyscraper valuations, or hundred-floor towers built in 1890 on a ten-square-metre footprint. Any model trained to predict value from structure — the whole reason someone wanted synthetic training data — learns from correlations that no longer exist.

The defects show up in specific, testable ways:

• Destroyed pairwise correlation. The Pearson or Spearman correlation between two attributes collapses toward zero. A regression coefficient that was strongly positive in real data becomes insignificant in synthetic data, and models trained on the synthetic set fail to generalize back to real inputs.
• Broken conditional structure. Even when a global correlation survives, the conditional relationship — how attribute B behaves given a specific band of attribute A — is lost. Land-use class should sharply constrain building height; independent sampling lets industrial parcels sprout residential towers.
• Severed spatial coupling. Attributes that should track location drift free of it. Temperature stops falling with latitude and elevation; income stops clustering by neighbourhood. The attribute table is internally plausible but spatially incoherent.
• Silent tail decoupling. Marginals can match while joint tail behaviour — the co-occurrence of extreme values — is wrong. Two attributes that spike together in real floods or heatwaves spike independently in the synthetic set, so any tail-risk model is miscalibrated.

Each defect passes a per-column check and fails the moment a consumer looks at two columns at once. The remedy is to model the marginals and the dependence structure separately and explicitly, then recombine them — the discipline this page formalizes.

Prerequisites & Toolchain

Attribute correlation modeling is numerically light but statistically exacting. Pin majors so the linear-algebra and rank-transform behaviour does not shift between environments:

Copy
python >= 3.10
numpy == 1.26.*
scipy == 1.11.*        # multivariate_normal, rankdata, norm, spearmanr
pandas == 2.*
geopandas == 0.14.*    # attribute join back onto geometry
pyproj == 3.*          # metric CRS for spatial-covariate fields

Two environment decisions dominate correctness. First, decide up front whether your target dependence is linear or rank-based. Pearson correlation is not invariant under the marginal transforms a copula applies; Spearman rank correlation is. If you fit a correlation matrix on raw values but reconstruct through arbitrary marginals, the realized correlation will drift from the target. Work in rank space (Spearman) end to end, or apply the closed-form Pearson-to-Spearman adjustment for the Gaussian copula, and document which you chose. Second, keep the correlation matrix positive semi-definite. Empirical rank-correlation matrices estimated column-pairwise are frequently not PSD, and the Cholesky factorization at the heart of the Gaussian copula will fail on them. Plan to project any estimated matrix to the nearest PSD matrix before it enters the sampler.

Core Concept: Marginals, Copulas, and Conditional Coupling

The governing idea comes from Sklar’s theorem: any multivariate distribution can be decomposed into its one-dimensional marginals and a copula that encodes the dependence between them, and the two pieces can be modeled independently. This separation is exactly what fixes the independent-column fallacy — you fit each marginal to its real column shape, fit the copula to the real dependence structure, and recombine them so both survive.

A Gaussian copula is the workhorse. You take the target correlation matrix $R$, form its Cholesky factor $L$ so that $LL^{\top} = R$, and draw correlated standard-normal vectors $z = L\xi$ from independent standard normals $\xi$. Applying the normal CDF $\Phi$ column-wise maps those to correlated uniforms $u = \Phi(z)$. The final step, the probability integral transform (PIT), pushes each uniform column back through the inverse CDF of that attribute’s own marginal: $x_j = F_j^{-1}(u_j)$. The result carries each attribute’s real marginal shape and the prescribed dependence — the two halves Sklar’s theorem let us model apart. The drill-down on preserving copula correlations in synthetic attributes covers the rank-space corrections and PSD repair that keep the realized matrix on target when marginals are heavy-tailed or discrete.

The Gaussian copula has one well-known limitation: it has no tail dependence. Extreme values in different attributes become asymptotically independent, which understates joint-extreme co-occurrence. When the consumer cares about co-extremes — joint flood depth and wind speed, simultaneous demand spikes — reach for a Student-t copula (symmetric tail dependence) or a Clayton/Gumbel copula (asymmetric), fit only their dependence parameter, and keep the same PIT machinery.

Conditional sampling handles the case where one attribute deterministically or near-deterministically constrains another. Rather than forcing an implausible relationship through a single global correlation, factor the joint into a chain of conditionals — sample land-use class first, then draw height from the height distribution conditioned on that class. This composes cleanly with the copula: model the continuous attributes jointly with a copula, condition that block on the discrete class drawn first.

Spatial-covariate coupling is the step that makes this a spatial problem rather than a generic tabular one. Real attributes co-vary with location: temperature with latitude and elevation, income with neighbourhood, soil moisture with a hydrological gradient. Encode that by conditioning the copula’s latent draws on covariate fields evaluated at each feature’s coordinates — add a location-driven mean shift to the latent normal before the PIT, or make the marginal parameters themselves functions of position. This is where attribute synthesis meets geometry, and why it lives beside the generators that place the points. For richer, learned dependence that a fixed copula cannot express, GAN-based spatial generation can synthesize geometry and attributes jointly, at the cost of the reproducibility guarantees a copula gives for free.

Step-by-Step Implementation

The pipeline runs as four auditable stages. Each block assumes the pinned toolchain and an explicit integer RNG seed.

Step 1 — Fit marginals and the rank-correlation target from real data. Estimate each column’s empirical CDF and the Spearman rank-correlation matrix. Keep them as the two separable artifacts Sklar’s theorem promises.

pythonCopy
import numpy as np
from scipy import stats

def fit_targets(real: np.ndarray):
    """real: (n_samples, n_attrs) array of real attribute values."""
    # Spearman rank correlation is invariant to the marginal transforms below.
    rank_corr, _ = stats.spearmanr(real)
    rank_corr = np.atleast_2d(rank_corr)
    # Store per-column sorted values as an empirical inverse-CDF lookup.
    marginals = [np.sort(real[:, j]) for j in range(real.shape[1])]
    return rank_corr, marginals

Step 2 — Repair the correlation matrix to positive semi-definite. A pairwise-estimated matrix is often not PSD. Project it to the nearest PSD matrix by clipping negative eigenvalues, then renormalize the diagonal to 1.

pythonCopy
def nearest_psd_corr(R: np.ndarray) -> np.ndarray:
    R = (R + R.T) / 2.0                       # enforce symmetry
    vals, vecs = np.linalg.eigh(R)
    vals = np.clip(vals, 1e-8, None)          # kill negative eigenvalues
    R_psd = vecs @ np.diag(vals) @ vecs.T
    d = np.sqrt(np.diag(R_psd))               # renormalize to unit diagonal
    return R_psd / np.outer(d, d)

Step 3 — Sample the Gaussian copula and apply the probability integral transform. Correlate independent normals through the Cholesky factor, map to uniforms, then push each column through its empirical inverse CDF.

pythonCopy
def sample_copula(rank_corr, marginals, n, rng_seed: int) -> np.ndarray:
    rng = np.random.default_rng(rng_seed)     # PCG64, explicit + reproducible
    # Adjust Spearman target to the Pearson correlation the Gaussian copula needs.
    pearson_target = 2.0 * np.sin(np.pi * rank_corr / 6.0)
    R = nearest_psd_corr(pearson_target)
    L = np.linalg.cholesky(R)                 # R = L @ L.T
    xi = rng.standard_normal((n, R.shape[0]))
    z = xi @ L.T                              # correlated standard normals
    u = stats.norm.cdf(z)                     # correlated uniforms in (0, 1)
    out = np.empty_like(u)
    for j, col in enumerate(marginals):       # PIT: inverse empirical CDF per column
        out[:, j] = np.quantile(col, u[:, j], method="linear")
    return out

Step 4 — Couple attributes to location. Shift the latent normal by a covariate field evaluated at each feature’s coordinates before the PIT, so a target attribute tracks latitude, elevation, or a neighbourhood field while keeping its marginal shape.

pythonCopy
def spatial_shift(z: np.ndarray, coords: np.ndarray, attr_idx: int,
                  beta: float) -> np.ndarray:
    """Add a location-driven mean shift to one latent column before PIT.
    coords: (n, 2) metric x/y; beta scales the covariate's latent influence."""
    covariate = (coords[:, 1] - coords[:, 1].mean()) / coords[:, 1].std()
    z = z.copy()
    z[:, attr_idx] += beta * covariate        # e.g. temperature falls with northing
    return z

Validation & Testing

An attribute synthesizer is only trustworthy if a CI gate proves the joint structure survived — not just the marginals. Gate on both, because a per-column check passes on exactly the datasets this page exists to reject.

pythonCopy
import numpy as np
from scipy import stats

def gate_attributes(real: np.ndarray, synth: np.ndarray,
                    corr_tol: float = 0.05, ks_tol: float = 0.05) -> None:
    # 1. Marginals: every column matches its real shape (KS test per attribute).
    for j in range(real.shape[1]):
        d, _ = stats.ks_2samp(real[:, j], synth[:, j])
        assert d < 0.15, f"marginal {j} drifted: KS D={d:.3f}"

    # 2. Dependence: the rank-correlation matrices agree in Frobenius norm.
    r_real, _ = stats.spearmanr(real)
    r_synth, _ = stats.spearmanr(synth)
    frob = np.linalg.norm(np.atleast_2d(r_real) - np.atleast_2d(r_synth))
    assert frob < corr_tol * real.shape[1], f"correlation drift: ||dR||_F={frob:.3f}"

    # 3. Conditional structure survives on a representative slice.
    hi = real[:, 0] > np.median(real[:, 0])
    hi_s = synth[:, 0] > np.median(synth[:, 0])
    m_real = real[hi, 1].mean()
    m_synth = synth[hi_s, 1].mean()
    assert abs(m_real - m_synth) / (abs(m_real) + 1e-9) < 0.1, "conditional mean drift"

The Frobenius-norm check on the full correlation matrix is the load-bearing gate: it catches the independent-column failure directly, because a decorrelated synthetic set produces a near-zero synthetic matrix and a large distance from the real one. Wire these as build-failing checks through the CI/CD integration pipeline. For a more principled distributional distance than pairwise correlation — one sensitive to higher moments the correlation matrix ignores — layer in the summary statistics the realism metrics evaluation stage formalizes, and gate a spatial-autocorrelation statistic to confirm the location coupling actually took.

Performance & Scale Considerations

Copula sampling is cheap per row but the surrounding matrix and join work needs care at continental feature counts.

• Cholesky once, sample many. Factor the correlation matrix a single time and reuse L for every batch; the O(d³) factorization is trivial for tens of attributes but must never sit inside the per-row loop. Draw the latent normals in vectorized batches, not row by row.
• Vectorize the PIT. np.quantile over a whole column is orders of magnitude faster than a Python-level inverse-CDF call per value. For very large marginals, precompute a fixed-resolution CDF lookup table and interpolate.
• Chunk the spatial join, not the sampler. The attribute draw is embarrassingly parallel; the expensive step at scale is joining attributes back onto geometry and evaluating covariate fields. Tile the extent and process tiles independently, deriving each tile’s RNG seed from its integer coordinates plus a pipeline version hash so results stay reproducible — the same deterministic-tiling discipline the async execution for large grids layer relies on.
• Stream, don’t accumulate. Write each tile’s attribute table to partitioned Parquet and release it rather than holding the full multi-attribute frame in memory. Peak heap then stays flat as feature count climbs.

Failure Modes & Troubleshooting

Frequently Asked Questions

Related

• Point Process Simulation Models — places the features whose attributes this stage populates and couples to location.
• GAN-Based Spatial Data Generation — learns geometry and attribute dependence jointly when a fixed copula is too rigid.
• Preserving Copula Correlations in Synthetic Attributes — rank-space corrections and PSD repair that keep the realized matrix on target.
• Realism Metrics & Evaluation — distributional distances that gate joint fidelity beyond pairwise correlation.
• Spatial Distribution & Pattern Generation — the parent area tying placement, attributes, and density together.