BindCraft2 Protein Binder Design: How to Generate and Score Target-Specific Binders?

Designing a protein binder that recognizes one target with high affinity and high specificity is a core need in biology and medicine, and doing this from a target structure alone, without long rounds of experimental optimization, would remove a major bottleneck. Identifying high-affinity candidates remains a bottleneck because binder design typically requires generating and virtually screening thousands of designs to recover a few promising hits [5]. This means that any pipeline that can raise the fraction of computationally generated designs that turn out to bind in the wet lab has direct practical value: fewer designs need to be synthesized and tested, which saves time and reagent cost.

Two families of computational approaches address this problem. Structure-based approaches such as RFDiffusion and AlphaProteo can generate binder backbones conditioned on a target structure [5]. A separate family, hallucination-based methods, optimizes a structure directly inside a structure-prediction network rather than sampling from a generative backbone model. BindCraft is a hallucination-based binder design method that leverages the AlphaFold2 network [5]. This distinction matters for a builder because it determines what software and weights are needed: a hallucination method needs a differentiable structure predictor and a gradient-based optimization loop, while a backbone-diffusion method needs a trained diffusion model plus a separate sequence design step.

The question of how good these binders actually are once synthesized is separate from how the pipeline works, and the evidence on this point is uneven across papers. BindCraft has reported experimental success rates ranging from 10% to 100% [7], which is a wide range that depends heavily on the target class. Confidence scores from structure-prediction models, particularly interface pTM (ipTM), correlate with binding strength [5], but the correlation of confidence metrics such as pLDDT, ipTM, and ipAE with experimental binding outcomes has been validated on only limited datasets and remains poorly understood across targets and design regimes [6]. A builder should therefore not treat any single confidence number as proof of binding, only as a filter that has shown some correlation in some settings.

Because of this gap between in silico confidence and wet-lab outcome, benchmarking efforts and scoring-function improvements have become active areas of work alongside the generative pipelines themselves. BindEnergyCraft was built specifically to give a denser, more informative optimization signal than ipTM [5], and ProtDBench was built to give a standardized, multi-verifier evaluation protocol for comparing generated binders [6]. Together these two threads, better optimization objectives and better benchmarks, are what a developer needs to understand before deciding which pipeline and which scoring function to adopt for a given target.

BindCraft is an open-source and automated pipeline for de novo protein binder design [7]. It leverages the weights of AlphaFold2 to generate protein binders from target structures [7]. This means the user does not train a new model: instead the pipeline runs an optimization procedure inside AlphaFold2's existing computation graph, treating parts of the input (the binder sequence or structure) as free variables and everything else (the target structure) as fixed.

The pipeline has been applied broadly. BindCraft was applied to cell-surface receptors, common allergens, de novo designed proteins, and multi-domain nucleases such as CRISPR-Cas9 [7]. In several of these applications, the designed binders were shown to have functional effects beyond simple binding. Designed binders were reported to reduce IgE binding to birch allergen, modulate Cas9 gene-editing activity, reduce the cytotoxicity of a foodborne bacterial enterotoxin, and redirect adeno-associated virus capsids for targeted gene delivery [7]. In one demonstration, BindCraft generated binders with nanomolar affinity without high-throughput screening or experimental optimization, including for targets without known binding sites [7].

Independent groups have since stress-tested BindCraft on specific target classes to see how well these claims hold up outside the original paper. For MDM2, BindCraft generated 70 unique peptides, 15 were synthesized, and 7 showed specific binding with nanomolar affinities [12]. Competition assays confirmed site-specific binding of the MDM2 peptides at the intended target site [12], which is an important check because binding that is not site-specific would not be useful for most applications. For WDR5, six of nine candidates bound the MYC binding WBM site with submicromolar affinity [12], and rational chemical modification improved the potency of one WDR5 binder six-fold to a KD of 39 nM [12]. This shows that BindCraft outputs can serve as a starting point for further medicinal-chemistry optimization, not just a final answer. However, none of the tested BindCraft-generated peptides for PD-1 and PD-L1 showed detectable binding [12], which is a direct counter-example showing that success is target-dependent and not guaranteed. BindCraft can generate high-affinity peptides solely from a target structure [12], but the PD-1/PD-L1 result shows this capability does not transfer uniformly to every target.

At larger scale, the Human Bindome embeds the experimentally benchmarked BindCraft method in an accelerated, parallelized framework with automated domain-level target selection [13]. Using this framework, 306,146 binder candidates were generated for 8,296 human proteins, covering 40.9% of the full proteome [13]. Each Human Bindome candidate includes a defined sequence, a predicted binder-target structure model, and in silico confidence metrics [13], which gives a builder a template for what a candidate record should contain when running a pipeline at scale.

BindCraft hallucinates binder structures with AlphaFold2 and performs gradient-based optimization of the interface predicted TM-score (ipTM), among other losses [5]. In practice this means the binder sequence (and implicitly its structure, as predicted by AlphaFold2) is iteratively updated by gradient descent so as to increase the ipTM score AlphaFold2 assigns to the predicted binder-target complex, alongside other loss terms. This is a fundamentally different mechanism from a diffusion model that samples a backbone and then designs a sequence for it: here the structure-prediction network itself acts as the objective function being optimized.

The choice of ipTM as the main optimization target has a known weakness. ipTM does not directly reflect the statistical likelihood of the full binder-target interaction [5]. This matters because the goal of binder design is to find sequences that are likely to actually fold into a stable, correctly-docked complex, and a score that does not track this likelihood can push the optimizer toward structures that look confident to the score but are not necessarily physically favorable complexes.

The specific mechanical reason for this weakness lies in how ipTM is computed. Because ipTM computes a maximum over target residue indices, it produces sparse gradients concentrated on only a small subset of interface residue pairs [5]. In gradient-based optimization, a sparse gradient means only a few positions in the structure receive a meaningful update signal at each step, while the majority of the interface is not directly shaped by the objective. This can lead to designs that satisfy the score at a few residue pairs without the rest of the interface being well-formed, which is one plausible reason for structural clashes or poor overall complex quality in generated designs.

This mechanistic diagnosis, made explicit in the BindEnergyCraft paper, is what motivated a search for an alternative scoring function that gives a denser gradient signal across the whole interface rather than a sparse one. The rest of the optimization framework, the AlphaFront2-based hallucination loop and the general concept of gradient-based structure and sequence updates, is retained; only the interface score being optimized needs to change.

pTMEnergy reinterprets AlphaFold predicted alignment error (pAE) logits as an energy-based model over binder-target complexes [5]. Rather than taking a single scalar summary (as ipTM does, via a max operation over residue pairs), pTMEnergy treats the full distribution of pAE logits as defining an energy landscape, giving a score that in principle draws on information from across the whole predicted interface rather than a small subset of it.

pTMEnergy provides a dense, differentiable signal intended to reflect the likelihood of a folded complex under AlphaFold's learned distribution [5]. Because it is differentiable, it can be substituted directly into the same gradient-based optimization loop that BindCraft already uses, without needing a different kind of optimizer. This directly targets the sparse-gradient weakness described above: a dense signal means more of the interface residue pairs contribute a meaningful gradient at each optimization step, in contrast to the max-based ipTM computation.

A key implementation detail for a builder is that pTMEnergy is derived from the same pAE outputs used to compute ipTM and requires no additional model calls or architecture modifications [5]. This means adopting pTMEnergy does not require retraining AlphaFold2 or adding new network components: the pAE logits AlphaFold2 already outputs are simply processed differently to produce the new score. This is an important practical point because it keeps the compute cost of the modified pipeline essentially the same as the original BindCraft pipeline, since no extra forward or backward passes through a different model are needed.

The resulting method, BindEnergyCraft, which retains BindCraft's optimization framework while replacing ipTM with pTMEnergy, achieved higher in silico binder success rates and reduced structural clashes than BindCraft, RFDiffusion, and ESM3 across multiple challenging targets [5]. This is an in silico comparison, evaluated across the challenging targets used in that paper's benchmark; the source claims do not report a wet-lab head-to-head between BindEnergyCraft and BindCraft, RFDiffusion, or ESM3 on the same synthesized candidates, so a builder should read the improvement as a benchmark-level, computational finding rather than a confirmed experimental advantage over these specific alternatives.

Because confidence scores from structure-prediction models are known to be imperfect proxies for real binding, and because their reliability is poorly characterized across targets, a standardized benchmark for comparing binder-design methods and scoring functions is needed. ProtDBench is one such benchmark, built specifically to evaluate generated binder backbones in a reproducible, multi-verifier way.

ProtDBench evaluates generated binder backbones by sampling multiple sequences for each backbone, applying a verifier and filter to each sequence, and computing the per-sequence success rate [6]. This procedure separates two design stages that are often conflated: the backbone generation step (which fixes the fold) and the sequence design step (which fixes the actual amino acid identities threaded onto that fold). By sampling several sequences per backbone, the benchmark can distinguish backbones that are easy to design good sequences for from backbones that are only occasionally successful.

ProtDBench defines the passing backbone set as the backbones for which at least one sampled sequence passes the filter [6]. This is a lenient, best-of-many criterion at the backbone level: a backbone counts as a success if it produced at least one good sequence, even if most sampled sequences for it failed. ProtDBench clusters only passing backbones by structural similarity using FoldSeek or TM-score [6], which lets the benchmark report not just a raw success count but also how structurally diverse the successful designs are, avoiding the case where all reported successes are near-duplicates of the same fold.

A central finding of this benchmark work is a caution against over-trusting any single verifier. ProtDBench reports substantial verifier-dependent bias and limited agreement among structure-prediction models used as evaluation verifiers [6]. This directly reinforces the earlier point that ipTM and related confidence metrics have validated correlation with real binding only on limited datasets [6]: if different structure-prediction models used as verifiers disagree substantially with each other, then a success rate computed with one verifier (for example AlphaFold2-based ipTM) cannot be assumed to hold when a different verifier (for example AlphaFold3) is used instead. For a builder, this means reporting results from a single verifier without disclosing which one was used, or without cross-checking against a second verifier, risks overstating confidence in a design.

Not all binder-design pipelines are hallucination-based. A separate, structure-based workflow illustrates both the promise and the practical difficulty of getting from a computational design to a validated binder. A hotspot-driven binder-design workflow integrated RFdiffusion, ProteinMPNN, AlphaFold, molecular dynamics simulations, and manual structural assessment [3]. This workflow combines a backbone-diffusion step (RFdiffusion), a sequence-design step (ProteinMPNN), a structure-confidence check (AlphaFold), a physical-simulation step (molecular dynamics), and human review, which is a heavier and more manually-supervised pipeline than the automated, single-network BindCraft loop.

The hotspot-driven workflow generated four protein designs that were subjected to experimental validation [3]. This is a small number of designs relative to the

Binder design pipelines increasingly draw on protein language models (pLMs) as a complementary or alternative sequence-generation tool. Protein language models (pLMs) are pre-trained on large-scale evolutionary sequence data and adapt large language models for biological sequences [17]. This is a distinct mechanism from AlphaFold2-based hallucination or RFdiffusion-based backbone generation: a pLM learns patterns directly from large sequence databases rather than from structural training data, and it is typically used to propose or score sequences rather than 3D coordinates.

Several specific architectures exist within this family, and a builder choosing among them should know their structural differences. ESM2 is a BERT-style encoder transformer model, while ProtGPT2 and ProGen2 are GPT-style decoder transformer models and ProtT5 is an encoder-decoder transformer model [17]. An encoder-only model like ESM2 is typically used to score or embed an existing sequence, a decoder-only model like ProtGPT2 or ProGen2 is typically used to generate new sequences autoregressively, and an encoder-decoder model like ProtT5 can be used for tasks that map one sequence representation to another. None of the cited claims describe pLMs being substituted directly into the BindCraft or BindEnergyCraft optimization loop, so a builder should treat pLM-based sequence generation as a separate, complementary tool for proposing or filtering candidate sequences, not as a documented drop-in replacement for AlphaFold2-based hallucination in these specific pipelines.

The central limitation for anyone building on this literature is that confidence scores are not proof of binding. The correlation of confidence metrics such as pLDDT, ipTM, and ipAE with experimental binding outcomes has been validated on only limited datasets and remains poorly understood across targets and design regimes [6]. This means that a high ipTM or pTMEnergy score should be read as a filter that increases the chance of finding a hit, not as a guarantee, and any reported success rate is tied to the specific target class and dataset it was measured on.

ipTM itself carries a known mechanistic weakness. ipTM does not directly reflect the statistical likelihood of the full binder-target interaction [5], and because it computes a maximum over target residue indices, it produces sparse gradients concentrated on only a small subset of interface residue pairs [5]. This is a specific, mechanistic reason to prefer pTMEnergy for the optimization objective when reproducing BindCraft-style pipelines, though the comparison between BindEnergyCraft and BindCraft on real wet-lab targets has not been reported in the cited claims, only the in silico comparison across benchmark targets [5].

Different benchmark and verifier choices also complicate any comparison across papers. ProtDBench reports substantial verifier-dependent bias and limited agreement among structure-prediction models used as evaluation verifiers [6], so success rates computed with different verifiers, or reported by different papers using different pipelines, cannot be directly compared without knowing which verifier was used. Similarly, results on specific target classes are not uniform: the same general method (BindCraft) succeeded strongly for MDM2 and WDR5 [12] but produced no detectable binders for PD-1 and PD-L1 [12], and a separate hotspot-driven workflow using entirely different tools (RFdiffusion, ProteinMPNN, AlphaFold, molecular dynamics [3]) succeeded for only one of four designs against the LEDGF PWWP domain [3], with that one success complicated by unwanted homodimerization [3]. None of these results should be extrapolated to targets or conditions not directly tested.

Finally, identifying high-affinity candidates remains a bottleneck because binder design typically requires generating and virtually screening thousands of designs to recover a few promising hits [5]. Even with an improved scoring function like pTMEnergy, or a large-scale generation effort like the Human Bindome covering 40.9% of the human proteome [13], the fundamental need to filter a large candidate pool down to a small number of wet-lab-testable designs has not been eliminated by any of the cited work, only partially mitigated. A developer should budget for this filtering step and should not expect any single in silico score to remove the need for experimental validation.

1. Set up the base BindCraft environment. Install the open-source and automated BindCraft pipeline for de novo protein binder design [7], which leverages the weights of AlphaFold2 to generate protein binders from target structures [7]. You will need AlphaFold2 model weights available locally and a working JAX/AlphaFold2 inference stack, since the pipeline runs gradient-based optimization inside this network rather than calling an external API.
2. Prepare the target structure input. Obtain a PDB structure file for your target protein. If you plan to reproduce a hotspot-driven design instead, note that a separate published workflow integrated RFdiffusion, ProteinMPNN, AlphaFold, molecular dynamics simulations, and manual structural assessment [3] as an alternative to single-network hallucination; decide up front which of the two families you are building, since they need different software stacks.
3. Choose the optimization objective. The default BindCraft loss is a hallucination-based, gradient-based optimization of the interface predicted TM-score (ipTM), among other losses [5]. If you want the improved variant, replace ipTM with pTMEnergy, which reinterprets AlphaFold predicted alignment error (pAE) logits as an energy-based model over binder-target complexes [5] and requires no additional model calls or architecture modifications [5], since it is derived from pAE outputs you already compute.
4. Run the optimization loop. Iterate gradient updates on the binder sequence/structure to increase the chosen score (ipTM or pTMEnergy) against the fixed target. Keep in mind that ipTM produces sparse gradients concentrated on only a small subset of interface residue pairs [5], so if you observe clashes or poorly-formed interfaces outside a few contact residues, this is a known failure mode of ipTM-based optimization, and pTMEnergy is the documented mitigation [5].
5. Generate a large candidate pool. Because identifying high-affinity candidates remains a bottleneck and typically requires generating and virtually screening thousands of designs to recover a few promising hits [5], plan your compute budget for a batch of many candidates per target, not a handful.
6. Score and filter candidates with structure-prediction models. Structure-prediction models such as AlphaFold2 and AlphaFold3 enable binding evaluation without crystal structures [5], and confidence scores, particularly ipTM, correlate with binding strength [5]. Use these scores as an initial filter, but remember this correlation has been validated on only limited datasets [6].
7. Benchmark with a multi-verifier, per-sequence protocol. Follow the ProtDBench approach: sample multiple sequences per generated backbone, apply a verifier and filter to each sequence, and compute the per-sequence success rate [6]. Define your passing backbone set as backbones with at least one sequence that passes the filter [6], and cluster only the passing backbones by structural similarity using FoldSeek or TM-score [6] to check design diversity.
8. Cross-check with more than one verifier. Because ProtDBench reports substantial verifier-dependent bias and limited agreement among structure-prediction models used as verifiers [6], run your filter with at least two different structure-prediction models and report both results rather than relying on one.
9. Select a shortlist for synthesis. From your filtered, clustered candidates, pick a diverse shortlist for wet-lab testing. Expect target-dependent variance: independent tests found 7 of 15 synthesized MDM2 peptides bound specifically [12], 6 of 9 WDR5 candidates bound submicromolar [12], but none of the tested PD-1/PD-L1 peptides bound at all [12], so treat your in silico shortlist as a hypothesis to be tested, not a guarantee.
10. Validate experimentally. Use a biophysical binding assay (as used to confirm low-micromolar affinity for a LEDGF PWWP domain binder [3]) and, where relevant, a competition assay to confirm site-specificity, as was done for the MDM2 peptides [3][12]. Watch for unexpected behavior such as homodimerization, which in one case interfered with binding in solution despite the protein still being co-crystallizable with its target at 2.1 Å resolution [3].
11. Report results with full conditions. When publishing or logging success rates, always state which verifier, which benchmark, and which target class produced the number, since success rates ranging from 10% to 100% [7] and the verifier-disagreement finding [6] mean an unqualified success-rate number is not informative on its own.

Pseudocode sketch of the core hallucination loop with the pTMEnergy option:

load target_structure from PDB file
initialize binder_sequence randomly
for step in range(num_steps):
    complex_pred = alphafold2_forward(binder_sequence, target_structure)
    if objective == "ipTM":
        loss = -complex_pred.ipTM
    elif objective == "pTMEnergy":
        loss = -pTMEnergy(complex_pred.pae_logits)
    grads = backprop(loss, wrt=binder_sequence)
    binder_sequence = update(binder_sequence, grads)
return binder_sequence, complex_pred.ipTM, complex_pred.pTMEnergy

The script below implements, on our own SQLite data, a simplified version of the two evaluation ideas central to this note: (1) a ProtDBench-style per-sequence success-rate and passing-backbone-set computation over generated candidates (claims R and S), using structural similarity clustering of passing designs via a TM-score-like proxy (claim T), and (2) a comparison of two scoring functions analogous to ipTM versus pTMEnergy (claims B, I, J, M, N, O, P), computed from our own `forecasts` table columns (q10, q50, q90, p_up) standing in for confidence-style scores, since our data does not contain real AlphaFold pAE logits or ipTM values. Because our schema has no binder-design tables, we treat each row of `forecasts` as one candidate design for one symbol/target, `p_up` as the ipTM-style sparse confidence score, and a dense combination of q10/q50/q90 spread as the pTMEnergy-style dense score (this substitution is clearly commented in the code; it is a demonstration of the scoring-comparison method, not a claim that our market data is protein data). The baseline to beat is the plain ipTM-style filter's success rate; the pTMEnergy-style dense score is expected to give a more stable, less sparse ranking, which we check by comparing how many distinct residue-pair-proxy features (q10, q50, q90) actually vary versus how much of the score change is driven by a single feature (sparse vs dense gradient proxy check, claim J). The script also builds a small ProtDBench-style report: per-candidate success rate, passing-backbone set size, and a toy clustering by rounding scores (proxy for FoldSeek/TM-score clustering, claim T). It prints its numbers and compares against a random-baseline success rate. Everything runs offline against `data.sqlite` (or QOURAT_DB), finishes in well under three minutes on CPU, and raises no exceptions if the tables exist (it creates a small synthetic dataset if the file has no data, so it is runnable standalone too).

import os
import sqlite3
import sys
import time
import numpy as np
import pandas as pd

# ---------------------------------------------------------------------------
# Config
# ---------------------------------------------------------------------------
DB_PATH = os.environ.get("QOURAT_DB", "data.sqlite")
RNG = np.random.default_rng(0)


def get_connection(path):
    # Opens (or creates) the sqlite file. No network access needed.
    conn = sqlite3.connect(path)
    return conn


def ensure_tables(conn):
    # Creates the required tables if they do not exist yet, so the script
    # is runnable even against a fresh empty file. Schema matches the spec.
    cur = conn.cursor()
    cur.execute("""CREATE TABLE IF NOT EXISTS bars (
        symbol TEXT, tf TEXT, ts TEXT, open REAL, high REAL, low REAL,
        close REAL, volume REAL)""")
    cur.execute("""CREATE TABLE IF NOT EXISTS forecasts (
        symbol TEXT, horizon TEXT, made_at TEXT, q10 REAL, q50 REAL,
        q90 REAL, p_up REAL)""")
    cur.execute("""CREATE TABLE IF NOT EXISTS trades (
        symbol TEXT, ts TEXT, price REAL, size REAL, side TEXT)""")
    cur.execute("""CREATE TABLE IF NOT EXISTS book (
        symbol TEXT, ts TEXT, level INTEGER, bid_price REAL,
        bid_size REAL, ask_price REAL, ask_size REAL)""")
    conn.commit()


def maybe_seed_forecasts(conn):
    # If the forecasts table is empty, populate it with a small synthetic
    # set of "candidate designs" so the pipeline below has something to
    # score. This stands in for generated binder candidates: each row is
    # one candidate for one "target" (symbol), analogous to one sampled
    # sequence for one backbone in ProtDBench (claim R).
    df = pd.read_sql_query("SELECT COUNT(*) AS n FROM forecasts", conn)
    if df["n"].iloc[0] > 0:
        return
    symbols = ["AAPL", "MSFT", "BTC-USD"]
    rows = []
    for sym in symbols:
        for i in range(60):
            q50 = RNG.normal(0, 1)
            spread = abs(RNG.normal(0.5, 0.2)) + 0.01
            q10 = q50 - spread
            q90 = q50 + spread
            p_up = float(np.clip(RNG.normal(0.5, 0.2), 0, 1))
            rows.append((sym, "1d", f"2024-01-{(i % 28) + 1:02d}T00:00:00",
                         q10, q50, q90, p_up))
    conn.executemany(
        "INSERT INTO forecasts (symbol, horizon, made_at, q10, q50, q90, p_up) "
        "VALUES (?,?,?,?,?,?,?)", rows)
    conn.commit()


# ---------------------------------------------------------------------------
# Scoring functions
# ---------------------------------------------------------------------------

def sparse_score(row):
    # Proxy for ipTM-style scoring: BindCraft optimizes ipTM, a score that
    # computes a max over target residue indices and gives sparse
    # gradients concentrated on a small subset of interface pairs
    # (claims B and J). Here we proxy this by taking the single most
    # extreme of the three quantile-derived features, i.e. a max-like
    # operation over a small set of "residue pair" proxies (q10, q50, q90).
    feats = np.array([abs(row["q10"]), abs(row["q50"]), abs(row["p_up"] - 0.5)])
    return float(feats.max())


def dense_score(row):
    # Proxy for pTMEnergy-style scoring: pTMEnergy is derived from the same
    # underlying outputs (here: q10, q50, q90, p_up) but combines them into
    # a dense signal rather than taking a single max (claims M, N, O).
    # We use a simple mean of the same features, so every feature
    # contributes to the score, not just the largest one.
    feats = np.array([abs(row["q10"]), abs(row["q50"]), abs(row["q90"]),
                       abs(row["p_up"] - 0.5)])
    return float(feats.mean())


def gradient_sparsity(row):
    # Checks how much of the sparse score is driven by a single feature
    # versus spread across features, as a numeric illustration of claim J
    # (sparse gradients concentrated on a small subset of interface pairs).
    # Returns the fraction of total feature magnitude held by the largest
    # single feature: close to 1.0 means very sparse, close to 1/n means
    # evenly dense.
    feats = np.array([abs(row["q10"]), abs(row["q50"]), abs(row["q90"]),
                       abs(row["p_up"] - 0.5)])
    total = feats.sum()
    if total = threshold
    per_seq_success_rate = float(df["passes"].mean())

    backbone_groups = df.groupby("symbol")
    passing_backbones = []
    for symbol, group in backbone_groups:
        if group["passes"].any():
            passing_backbones.append(symbol)
    passing_backbone_set = set(passing_backbones)

    return per_seq_success_rate, passing_backbone_set, df


def cluster_passing_backbones(df, passing_backbone_set):
    # Proxy for claim T: ProtDBench clusters only passing backbones by
    # structural similarity using FoldSeek or TM-score. We do not have
    # real structures, so we use a simple proxy: cluster passing backbones
    # by rounding their mean score to one decimal place, as a stand-in
    # for a structural-similarity bucket.
    passing_df = df[df["symbol"].isin(passing_backbone_set)]
    if len(passing_df) == 0:
        return {}
    means = passing_df.groupby("symbol")["score"].mean()
    clusters = {}
    for symbol, m in means.items():
        key = round(m, 1)
        clusters.setdefault(key, []).append(symbol)
    return clusters


# ---------------------------------------------------------------------------
# Baseline: random scoring, to check the two methods beat pure chance
# ---------------------------------------------------------------------------

def random_baseline_success_rate(n_rows, threshold, n_trials=200):
    # A simple baseline: if scores were pure random noise in [0, 1.5],
    # what fraction would pass the threshold? This is the number to beat.
    rates = []
    for _ in range(n_trials):
        fake_scores = RNG.uniform(0, 1.5, size=n_rows)
        rates.append(float((fake_scores >= threshold).mean()))
    return float(np.mean(rates))


def main():
    start = time.time()
    conn = get_connection(DB_PATH)
    ensure_tables(conn)
    maybe_seed_forecasts(conn)

    df = pd.read_sql_query(
        "SELECT symbol, horizon, made_at, q10, q50, q90, p_up FROM forecasts",
        conn)
    conn.close()

    if len(df) == 0:
        print("No forecast rows found even after seeding, aborting.")
        sys.exit(1)

    threshold = 0.5  # fixed filter threshold for both scoring functions

    sparse_rate, sparse_passing, sparse_df = protdbench_style_eval(
        df, sparse_score, threshold)
    dense_rate, dense_passing, dense_df = protdbench_style_eval(
        df, dense_score, threshold)

    sparsity_values = df.apply(gradient_sparsity, axis=1)
    mean_sparsity = float(sparsity_values.mean())

    sparse_clusters = cluster_passing_backbones(sparse_df, sparse_passing)
    dense_clusters = cluster_passing_backbones(dense_df, dense_passing)

    baseline_rate = random_baseline_success_rate(len(df), threshold)

    print("=== BindCraft-style scoring comparison (proxy data) ===")
    print(f"Rows (candidate designs) evaluated: {len(df)}")
    print(f"Backbones (symbols) evaluated: {df['symbol'].nunique()}")
    print()
    print("-- ipTM-style (sparse, max-based) scoring --")
    print(f"Per-sequence success rate: {sparse_rate:.3f}")
    print(f"Passing backbone set size: {len(sparse_passing)} "
          f"of {df['symbol'].nunique()}")
    print(f"Structural-similarity clusters (proxy): {sparse_clusters}")
    print()
    print("-- pTMEnergy-style (dense, mean-based) scoring --")
    print(f"Per-sequence success rate: {dense_rate:.3f}")
    print(f"Passing backbone set size: {len(dense_passing)} "
          f"of {df['symbol'].nunique()}")
    print(f"Structural-similarity clusters (proxy): {dense_clusters}")
    print()
    print(f"Mean gradient sparsity of ipTM-style score (1.0 = fully "
          f"sparse, 0.25 = fully dense over 4 features): {mean_sparsity:.3f}")
    print()
    print(f"Random baseline success rate to beat: {baseline_rate:.3f}")
    print()
    if dense_rate >= sparse_rate and dense_rate > baseline_rate:
        print("RESULT: dense (pTMEnergy-style) scoring matched or beat "
              "sparse (ipTM-style) scoring, and both beat the random baseline.")
    else:
        print("RESULT: dense scoring did not clearly beat sparse scoring "
              "on this data; inspect thresholds and feature scaling.")

    elapsed = time.time() - start
    print(f"\nElapsed time: {elapsed:.2f} seconds")


if __name__ == "__main__":
    main()