Construction Cost Prediction

This repository is a working environment for estimating construction cost per square metre from a mix of structured economic indicators and geospatial observations. The core idea is simple, but the modeling problem is not: construction cost is shaped by national macro conditions, local infrastructure, geography, climate risk, and visible patterns of development on the ground. A useful system has to combine those signals without overreacting to noise, country-specific quirks, or one-off shocks.

The project is built around a competition-style dataset spanning quarterly observations in Japan and the Philippines. Each row describes a location and period, then ties that tabular record to satellite composites. The tabular side captures context such as GDP, inflation, accessibility to transport links, and hazard-related variables. The raster side contributes Sentinel-2 multispectral imagery and VIIRS nighttime lights, which together provide a rough picture of land cover, built density, and economic activity. The goal is not just to fit a leaderboard target, but to understand which parts of that signal are genuinely informative and which parts are artifacts of the benchmark.

In practice, this repo is where we assemble that pipeline end to end: ingest raw competition files, normalize keys and timestamps, derive satellite features, store intermediate artifacts in DuckDB and Parquet, and experiment with multimodal models. Some experiments are deliberately lightweight and tabular-first; others bring in learned remote-sensing representations such as CLAY embeddings. The emphasis is on building a workflow that makes it easy to move between exploratory analysis, feature generation, and modeling without rewriting the whole stack each time.

What We Are Trying To Accomplish

The immediate task is to predict cost in USD per square metre as accurately as possible for unseen location-quarter pairs. The broader objective is more ambitious: build a workflow for cost estimation that is spatially aware, economically grounded, and resilient to messy real-world data. A model that performs well here should be able to learn differences between countries, pick up changes over time, and make use of remotely sensed context where tabular data alone is too coarse.

That makes this project a useful test bed for several questions:

• How much predictive power comes from macroeconomic and infrastructure variables before imagery is even introduced?
• When do satellite-derived features add signal, and when do they mainly act as a noisy proxy for variables we already know?
• How should temporal structure be handled so that validation reflects deployment conditions instead of leaking future information backward?
• Can learned image embeddings improve on hand-crafted spectral features for a small, heterogeneous regression dataset?

The repository is intended to support that investigation rather than lock us into one modeling recipe.

Competition Context And Limitations

One important lesson from this benchmark is that it should not be treated as a clean measure of machine learning quality. The target is closely related to publicly available government statistics, which means strong leaderboard performance may come from reverse-engineering the data-generating process rather than from better modeling. That does not make the task useless, but it does change what a good result means.

For this reason, the project treats the competition as both a prediction exercise and a data-forensics exercise. We want to understand where the benchmark is informative, where it is fragile, and where it may reward reconstruction of external administrative data more than generalizable learning. The remote-sensing component is still valuable as a research direction, especially for cases where public cost statistics are incomplete, delayed, or inconsistent. But any conclusions drawn from the leaderboard need to be filtered through the possibility of target leakage, country-specific data revisions, and shortcuts that would not survive in a production forecasting setting.

Repository Workflow

The working flow in this repository is:

1. Set up the environment.
2. Ingest the tabular and raster data into a normalized workspace.
3. Generate engineered satellite features and, if needed, learned raster embeddings.
4. Use the resulting tables for exploration, validation design, and model experiments.

Intermediate artifacts are written to inferred/, including:

• base_table.parquet for the normalized tabular backbone
• satellite_features.parquet for derived raster statistics
• workspace.duckdb for joined analysis tables
• sentinel2_embeddings.parquet when CLAY embeddings are generated
• implisat_features.parquet when ImpliSat scene-complexity features are generated (see below)

The primary working surface is a set of marimo notebooks under notebooks/, covering ingestion, exploratory analysis, and embedding generation. For the regression benchmark notebook, feature views and model preprocessing are separate concepts: the embedding views expose raw CLAY or SatMAE embedding columns, while individual pipelines decide whether to apply PCA, clustering, or a native mixed-type model fit.

Data Shape

At a high level, each observation combines:

• Geographic identity and quarter-level timing
• Economic variables such as deflated GDP and CPI-like signals
• Infrastructure and geographic accessibility indicators
• Hazard and climate descriptors
• A Sentinel-2 composite linked by filename
• A VIIRS nighttime-light composite linked by filename
• A scalar regression target for the training split

Sentinel-2 contributes multispectral bands that support feature engineering from reflectance patterns and spectral indices. VIIRS adds a compact proxy for human activity intensity and settlement structure. In this repo, those raster sources are converted into model-friendly tabular features and can also be passed through pretrained remote-sensing encoders.

Validation Modes

The choice of validation method is something we found to be quite defining for this project and data given. The competition does not look like a clean independent and identically distributed (IID) regression task in which each row can be treated as exchangeable with every other row. Locations and time periods recur. Japan and the Philippines do not cover the same history in the same way. And some years, most obviously 2020 (COVID), carry shocks that are large enough to change the interpretation of an error by themselves. A convincing result therefore has to say not only which model was used, but include failure modes of the generalizations was being tested.

For that reason the benchmark stack exposes three validation modes. Country-year holdout is the main one and is implemented as leave-one-group-out over the derived country-year column. Each fold withholds an entire country-year slice at once and asks the model to reconstruct it from the remaining slices. This is the most informative stress test for the benchmark because it preserves the structured nature of the task: the model still sees the same overall location universe and the same broad time range, but it cannot lean on direct exposure to the held-out macro regime within that country-year. In practice this is where differences between Japan and the Philippines matter most. The two countries do not contribute interchangeable observations, and the coverage of unusual years is not symmetric, so a method that looks stable under easier splits can degrade sharply once a whole country-year block is removed.

The second mode, group-location, is a 5-fold GroupKFold over location. This is a deliberately different question. Instead of asking whether the model can infer missing cells inside a familiar spatiotemporal table, it asks whether the model can transfer to locations it has never seen at all. That is a harder spatial extrapolation problem than the competition itself appears to require, but it could still be useful for comparison. It tells us how much of the score is coming from reusable geographic and macroeconomic structure, and how much is coming from memorized location-specific effects that would disappear the moment the location "cheat" is removed.

The third mode, kfold, is ordinary shuffled 5-fold cross-validation. It is included because it is convenient for rapid iteration and rough ranking of pipelines when we just want some signal. But it is the weakest of the three as a substantive estimate. Random row-wise splits allow location identity and time-specific structure to leak across train and validation folds, which makes the reported scores easier than the benchmark geometry suggests. In this repository it is treated as a permissive diagnostic, not as the primary number to anchor a conclusion on.

The reason we mainly focus on country-year is that it is most aligned with structure in the data that we are actually trying to understand. The difficulty in this competition turned out to be about whether a model can separate stable location effects, country-level macro shifts, and imagery-conditioned residual information without borrowing the answer from a nearly identical slice of time. A country-year holdout does not solve every evaluation problem, but it is the cleanest available compromise between realism and credibility in a dataset whose main difficulty comes from structured missing cells rather than from classical out-of-sample independence.

Environment Setup

Express Exploration

You can use the Benchmark Log Explorer and EDA notebooks to investigate results stored in this repo, without cloning it or installing anything on your machine. Just open a terminal and type:

nix run github:janinge/ccp-ml

Assuming that you have Nix installed, this will fetch the repo and its dependencies automatically and open a notebook browser in your browser.

Competition Data

Download evaluation_dataset.zip and train_dataset.zip from Solafune.

After downloading, either:

1. Extract the contents of both zip files into a data/ directory at the project root.
2. Create a symlink here named data that points to a directory containing the extracted contents of both zip files.

Software Dependencies

Recommended Setup (using Nix w/included uv)

1. Clone the repository using git clone https://github.com/janinge/ccp-ml.git && cd ccp-ml.
2. Enter the dev shell with nix develop (or enable direnv allow).
3. Launch the notebook launcher: uv run marimo run notebooks

This lists all notebooks in notebooks/ and lets you start them from the web interface.

Without Nix

You can still run this repo with uv or pip if you do not use Nix.

1. Clone the repository.
2. Create/activate a Python 3.13 environment.
3. Install either of these, then run:
   • uv: uv sync (plus groups below)
   • pip: pip install -e . (plus extras below)

pip can work, if you're lucky:

• no lockfile parity with uv.lock (version drift is likely)
• geospatial dependencies must be managed manually on your OS
• higher chance of wheel/build failures (especially EO/GPU stack)

Dependency Groups

Base install covers common notebook + benchmark usage. Add groups as needed:

• Earth observation I/O/features: --group eo (or .[eo])
• UI widgets used in benchmark/EDA notebooks: --group widgets (or .[widgets])
• Model family extras (CatBoost/TabPFN/EVoC): --group models (or .[models])
• Embedding stack (CLAY/SatMAE/UMAP tooling): --group embeddings (or .[embeddings])
• PyTorch runtime on Linux GPU hosts: --group gpu (or .[gpu])

Typical installs:

• Full local research stack (CPU): uv sync --group eo --group widgets --group models --group embeddings
• Linux GPU host: NIXPKGS_ALLOW_UNFREE=1 nix develop .#gpu --impure && uv sync --group eo --group widgets --group models --group embeddings --group gpu

Notebook/script group dependencies:

• notebooks/ingest.py: base + eo
• notebooks/eda.py: base + eo + widgets
• notebooks/regression_benchmarks.py: base + widgets + models (+gpu when running NAM)
• notebooks/catboost_shap_explorer.py: base + models
• notebooks/clay_embeddings.py: base + eo + embeddings + gpu
• notebooks/embedding_geometry_explorer.py: base + embeddings
• notebooks/benchmark_log_explorer.py: base
• transforms/implisat_embeddings.py: base + eo + gpu
• transforms/satmae_embeddings.py: base + eo + embeddings + gpu
• experiments/regression/run_benchmarks.py: base (+models for CatBoost/TabPFN pipelines, +gpu for NAM)

Authentication Tokens For Model Downloads

If you use pretrained model downloads, set these environment variables first:

# Hugging Face: improves reliability and download rate limits for model pulls
export HF_TOKEN="hf_..."

# TabPFN: easiest path is to accept the license at https://ux.priorlabs.ai/
# then copy your token and export it:
export TABPFN_TOKEN="..."

With TABPFN_TOKEN set, TabPFN can fetch gated pretrained checkpoints without interactive prompts.

CLI Benchmark Runs

To run benchmark suites without marimo:

• python -m experiments.regression.run_benchmarks list
• python -m experiments.regression.run_benchmarks run-all --cv-mode all --tabpfn-mode full --nam-mode full --tune-catboost --tune-random-forest --n-trials 10

Included Tooling

• Core analysis: marimo, duckdb, pandas, numpy, pyarrow, scikit-learn
• Notebook and visualization support: Mainly Vega-Altair and VegaFusion
• Image processing for "earth observation": rasterio, rioxarray, xarray, tifffile
• Modeling extras: EVōC, TabPFN, and CLAY

ImpliSat Scene Complexity Features

One direction this project explores is whether the compressibility of a satellite image can serve as a meaningful proxy for scene complexity and, by extension, construction activity. The intuition is that images containing fine-grained urban structure and dense built fabric are harder to compress accurately than images of homogeneous terrain.

This is operationalized through a framework called ImpliSat (Cho et al., 2025), which compresses Sentinel-2 multispectral imagery using Implicit Neural Representations (INRs). Rather than storing pixel values directly, an INR trains a small neural network to map spatial coordinates to pixel intensities, encoding the image as a set of network weights. ImpliSat extends this with a hypernetwork that conditions each band's representation on its ground sampling distance (GSD) — the physical size of one pixel on the ground, which ranges from 10 m to 60 m across Sentinel-2 bands — and a unique band identifier. This conditioning is applied through Fourier-frequency modulation, which lets the model adapt the spatial frequency content of the representation to match each band's native resolution rather than treating all bands as equivalent.

The paper's primary contribution is a compression framework. The feature engineering here repurposes it differently: instead of using the trained network weights as scene descriptors, we measure how well the network manages to reconstruct each image patch, quantified by PSNR (Peak Signal-to-Noise Ratio, a standard image-quality metric measuring reconstruction fidelity in decibels). A low PSNR indicates the network struggled to represent the patch — typically because it contains high-frequency detail such as dense rooftops, road networks, or industrial structures. A high PSNR indicates a smooth, easily representable scene. Extracting PSNR over a 4 × 4 spatial grid and summarizing it as mean, standard deviation, minimum, and maximum per band produces twelve scalars that together capture both average scene complexity and its spatial distribution across the image.

To keep computation tractable, only three bands are used — one from each of Sentinel-2's three resolution tiers: B2 at 10 m, B5 at 20 m, and B1 at 60 m. Bands within the same resolution tier are expected to carry very similar spatial structure and would produce redundant features. The network architecture is also scaled down from the paper's configuration (which used larger capacity and tens of thousands of iterations on a high-end GPU), using smaller hidden dimensions and fewer iterations for each image fit.

A further efficiency choice is to compute features per location rather than per observation. The dataset contains 1 024 training rows but only 125 unique locations. Since structural change within each location is gradual over the 2019–2024 time span covered by the training data, a single representative image per location — the most recent one in the training split — is a reasonable and far cheaper approximation. This reduces computation by roughly 8× with limited expected loss in downstream predictive power.

One deliberate departure from the original paper is that the stochastic modulation vector used during image fitting is fixed per band rather than resampled on each forward pass. Random resampling improves single-image reconstruction quality, but it makes PSNR values across different images incomparable, since each location would be evaluated under a different modulation basis. Fixing the vector ensures that differences in PSNR reflect scene complexity rather than variance introduced by the fitting procedure.

The twelve features are stored in inferred/implisat_features.parquet (one row per location) and merged into the workspace by location identifier. They appear in the tabular_core_plus_implisat feature view, which combines the standard tabular baseline with these scene-complexity descriptors.

To generate these features, run:

nix develop .#gpu
uv sync --group eo --group gpu
uv run python -m transforms.implisat_embeddings

Once the parquet exists, the tabular_core_plus_implisat view is automatically available in the regression benchmark notebook and in all tabular-capable pipeline runs.

Neural Additive Models

A recurring tension in this project is that prediction quality and interpretability tend to pull in opposite directions. Tree-based models and deep networks capture the non-linear structure that construction cost data exhibits — the relationship between GDP and material prices is not linear, and neither is the effect of climate zone on build standards — but they make it difficult to reason about how individual variables are contributing. Linear models stay transparent, but they leave that non-linear structure on the table.

Neural Additive Models (NAMs) offer a middle path. Rather than learning a single monolithic function over all inputs, a NAM decomposes the prediction into a sum of independent, feature-wise terms: each input variable is routed through its own small neural network, and the final output is the sum of all per-feature contributions plus a global bias. Because each sub-network operates on exactly one feature, the learned transformation for each variable can be inspected in isolation — what shape the model has found for GDP, for CPI, for distance to the capital — without that interpretation being obscured by high-dimensional interactions. This is the neural-network analogue of generalized additive models (GAMs), a well-established class of statistical regression models that decompose predictions into a sum of smooth per-feature curves.

In the implementation included here, each per-feature network is a three-layer multilayer perceptron (MLP) — a standard feedforward neural network — with two hidden layers of 64 units and ReLU activations. ReLU (Rectified Linear Unit) is a simple non-linearity that passes positive values through unchanged and clips negative values to zero; in combination with multiple layers, it produces piecewise non-linear mappings that can approximate a wide range of smooth functions. Each sub-network outputs a single scalar contribution, and the final prediction is the sum of all scalar outputs plus the bias. Training uses the Adam optimizer — an adaptive gradient method that adjusts the effective learning rate per parameter — with mean-squared-error loss, mini-batches of 32 samples, and data shuffled at each epoch so that the gradient estimates do not depend on presentation order.

During development, an alternative activation function called ExU (Exponential Units) was also evaluated. ExU aims to produce smoother and more expressive feature transformations than ReLU by applying an exponential non-linearity, but that expressiveness comes at a cost: small shifts in input can produce large shifts in output, making the network highly sensitive to feature scaling and to the choice of learning rate. On a dataset of this size — a few hundred training rows per country — that sensitivity was enough to destabilize training and produce worse generalization than the simpler ReLU architecture, even after careful tuning. The ReLU variant was therefore retained as the final model. The episode is a useful reminder that activation functions designed to improve representational capacity in large-data regimes do not necessarily transfer to small, heterogeneous tabular datasets.

Within the regression benchmark, NAM runs alongside the other tabular-capable pipelines and supports all eight feature views, including satellite-derived and embedding-augmented variants. Two modes are exposed in the benchmark notebook: Fast, which runs 20 training epochs for rapid iteration, and Full, which runs 50 epochs for a more converged result. The NAM mode is off by default because training a separate network per feature makes it considerably slower than tree-based models, in particular on the higher-dimensional embedding views where one sub-network is instantiated per embedding dimension.

TabPFN

Almost every other model in this stack is fit the usual way: parameters are estimated from the training rows, optionally tuned on held-out folds, and the resulting weights or trees are then applied to test rows. TabPFN (Hollmann et al., 2025) is structurally different. It is a transformer trained once, in advance, on a very large collection of synthetic supervised tasks drawn from a prior over plausible tabular data-generating processes. At prediction time the model is not retrained; instead the entire training set is concatenated with the test rows and passed through the network in a single forward pass, with attention serving as the mechanism that conditions test predictions on the training examples in context. Architecturally this is a prior-data fitted network: training-time meta-learning is used to approximate Bayesian inference over the synthetic prior, and inference-time prediction reduces to one evaluation of that approximation on a new task. The second-version release extends the original tabular classifier to regression and to mixed numeric–categorical inputs, which is what makes it applicable here.

This setup raises a question that is not obvious from architecture diagrams alone: what does the input have to look like for TabPFN to be competitive, and is the answer different from what works for the rest of the stack? Two configurations are evaluated. The first restricts the input to numeric columns only, with missing values median-imputed and remaining columns standardized, producing a dense matrix in the format the original TabPFN paper used. The second passes the selected feature view through largely unchanged and supplies the positions of categorical columns to the model, relying on the v2 preprocessor to handle mixed types natively. Both configurations are run across the same set of feature views as the other pipelines so that the choice between encoded-numeric and native-categorical inputs can be separated from the choice between tabular-only and embedding-augmented inputs.

The aggregate picture is that TabPFN reaches the same performance band as the well-tuned tree ensembles without any tuning of its own, and is occasionally a few hundredths of an RMSLE point behind the best tree baseline rather than ahead of it. That is a defensible result for a model with no per-task fitting, but it is not the dominant performance reported on standard tabular benchmarks. The variance across folds is what makes it less attractive in practice: TabPFN is competitive on the folds drawn from the better-sampled country and degrades sharply on the folds drawn from the country with fewer locations, sparser observations, and a different cost regime. The model is therefore not uniformly weaker than the ensembles; it is weaker specifically where the in-context prior has to extrapolate across a structural shift inside the dataset, and the absence of a built-in mechanism for additive location or time effects shows up as a long tail of large errors.

The shape of the input matters at least as much as the choice of model. Raw embedding views with several hundred dimensions are a poor fit for TabPFN regardless of configuration: the attention context is finite and the extra dimensions consume capacity that tree-based methods would simply ignore. PCA-reduced and cluster-summarized variants of the same embeddings recover most of the lost ground, which is why those are the embedding views routed into the TabPFN pipelines at all. On tabular views, the mixed-type configuration consistently outperforms the numeric-only one, suggesting that the v2 categorical handling carries real information and that flattening categoricals into one-hot indicators is the wrong move for this model. What ultimately throws TabPFN off in this dataset is less feature engineering than the geometry of the split itself: training and evaluation share the same locations and quarters but in mutually exclusive cells, and a model whose only inductive bias comes from an architecture-level synthetic prior cannot replace the explicit row- and column-effect structure that this matrix-completion-shaped task seems to reward. TabPFN earns its place as a strong zero-tuning reference point, but on this benchmark it is most useful as a calibration signal for the other pipelines rather than as the model that closes the remaining gap.

CLAY Checkpoint

If you want to generate CLAY embeddings, keep the weights outside the repository and point the environment variable at a local checkpoint path:

export CLAY_CHECKPOINT_PATH="$HOME/.cache/clay/clay-v1.5.ckpt"

The current code expects the CLAY v1.5 checkpoint layout.

Project Stance

This repository is not organized around a claim that one architecture has already won. It is organized around a harder question: what can actually be learned about construction cost from the combination of macro data and earth observation data once shortcuts, leakage, and benchmark oddities are taken seriously? That is the standard the experiments here are meant to meet.