Scorer

Compare RNA 3D structures using circular statistics on torsion angles.

This module implements a scoring pipeline that computes signed (MCD) and unsigned (MCQ) circular differences between backbone torsion angles of a target and model RNA structure, then derives a composite similarity score in [0, 1].

Terminology

MCD (Mean Circular Deviation): Signed angular difference, range [-pi, pi]. MCQ (Mean Circular Quality): Unsigned (absolute) angular difference, range [0, pi].

bootstrap_ci(data, statistic_fn, n_bootstrap, alpha=0.05)

Compute a bootstrap confidence interval for a circular statistic.

Parameters:

Name	Type	Description	Default
data		Input data (list of floats).	required
statistic_fn		Callable that takes a list and returns a scalar.	required
n_bootstrap		Number of bootstrap resamples.	required
alpha		Significance level (default 0.05 for 95% CI).	0.05

Returns:

Type	Description
	Tuple of (lower_bound, upper_bound) for the confidence interval.

Source code in src/rnapolis/scorer.py

def bootstrap_ci(data, statistic_fn, n_bootstrap, alpha=0.05):
    """Compute a bootstrap confidence interval for a circular statistic.

    Args:
        data: Input data (list of floats).
        statistic_fn: Callable that takes a list and returns a scalar.
        n_bootstrap: Number of bootstrap resamples.
        alpha: Significance level (default 0.05 for 95% CI).

    Returns:
        Tuple of (lower_bound, upper_bound) for the confidence interval.
    """
    statistics = []
    for _ in range(n_bootstrap):
        sample = np.random.choice(data, size=len(data), replace=True).tolist()
        statistics.append(statistic_fn(sample))
    lower = np.percentile(statistics, 100 * alpha / 2)
    upper = np.percentile(statistics, 100 - 100 * alpha / 2)
    return float(lower), float(upper)

circular_mad(data, median)

Compute the circular median absolute deviation.

Parameters:

Name	Type	Description	Default
data		Angular values in radians.	required
median		Circular median in radians.	required

Returns:

Type	Description
	The median of absolute wrapped deviations from median.

Source code in src/rnapolis/scorer.py

def circular_mad(data, median):
    """Compute the circular median absolute deviation.

    Args:
        data: Angular values in radians.
        median: Circular median in radians.

    Returns:
        The median of absolute wrapped deviations from ``median``.
    """
    data = np.asarray(data)
    deviations = np.abs(np.arctan2(np.sin(data - median), np.cos(data - median)))
    return float(np.median(deviations))

compute_mcd_mcq(target, model)

Compute signed (MCD) and unsigned (MCQ) angular differences between two structures.

For each residue matched by chain, number, and name, wraps the angular difference of each backbone torsion angle into [-pi, pi] (MCD) and takes its absolute value (MCQ).

If a residue at position i in the target does not match the same position in the model, a lookup by residue identity is used as a fallback.

Parameters:

Name	Type	Description	Default
target		DataFrame of torsion angles for the reference structure.	required
model		DataFrame of torsion angles for the model structure.	required

Returns:

Type	Description
	Tuple of (signed_diffs, unsigned_diffs) where each is a list of
	floats in radians.

Source code in src/rnapolis/scorer.py

def compute_mcd_mcq(target, model):
    """Compute signed (MCD) and unsigned (MCQ) angular differences between two structures.

    For each residue matched by chain, number, and name, wraps the angular
    difference of each backbone torsion angle into [-pi, pi] (MCD) and takes
    its absolute value (MCQ).

    If a residue at position *i* in the target does not match the same position
    in the model, a lookup by residue identity is used as a fallback.

    Args:
        target: DataFrame of torsion angles for the reference structure.
        model: DataFrame of torsion angles for the model structure.

    Returns:
        Tuple of (signed_diffs, unsigned_diffs) where each is a list of
        floats in radians.
    """
    signed_diffs = []
    unsigned_diffs = []
    model_lookup = _build_model_lookup(model)

    for i in range(len(target)):
        target_row = target.iloc[i]
        chain_id = target_row["chain_id"]
        residue_number = target_row["residue_number"]
        residue_name = target_row["residue_name"]

        if (
            i < len(model)
            and model.iloc[i]["chain_id"] == chain_id
            and model.iloc[i]["residue_number"] == residue_number
            and model.iloc[i]["residue_name"] == residue_name
        ):
            s, u = _compute_angle_diffs(target_row, model.iloc[i])
            signed_diffs.extend(s)
            unsigned_diffs.extend(u)
        else:
            key = (chain_id, residue_number, residue_name)
            if key in model_lookup:
                j = model_lookup[key]
                s, u = _compute_angle_diffs(target_row, model.iloc[j])
                signed_diffs.extend(s)
                unsigned_diffs.extend(u)

    return signed_diffs, unsigned_diffs

compute_score(metrics)

Compute a composite similarity score from circular statistics.

Combines three weighted sub-scores:

• Fit: How close the mean and median MCQ are to zero.
• Concentration: How tightly the MCQ distribution is concentrated.
• Uniformity: How concentrated the MCD distribution is around zero, indicating lack of systematic bias.

Parameters:

Name	Type	Description	Default
metrics		Dictionary containing at minimum the keys mcq, medcq, rmcq, circular_mad_mcq, p_watson_dmcq, p_sim_test_dmcq, mcd, medcd, rmcd, circular_mad_mcd, p_rayleigh_dmcd, p_wilcoxon_dmcd.	required

Returns:

Type	Description
	Similarity score in [0, 1], where 1.0 means identical structures.

Source code in src/rnapolis/scorer.py

def compute_score(metrics):
    """Compute a composite similarity score from circular statistics.

    Combines three weighted sub-scores:

    - **Fit**: How close the mean and median MCQ are to zero.
    - **Concentration**: How tightly the MCQ distribution is concentrated.
    - **Uniformity**: How concentrated the MCD distribution is around zero,
      indicating lack of systematic bias.

    Args:
        metrics: Dictionary containing at minimum the keys ``mcq``, ``medcq``,
            ``rmcq``, ``circular_mad_mcq``, ``p_watson_dmcq``, ``p_sim_test_dmcq``,
            ``mcd``, ``medcd``, ``rmcd``, ``circular_mad_mcd``,
            ``p_rayleigh_dmcd``, ``p_wilcoxon_dmcd``.

    Returns:
        Similarity score in [0, 1], where 1.0 means identical structures.
    """
    wf = WEIGHTS_FIT
    wc = WEIGHTS_CONCENTRATION
    wu = WEIGHTS_UNIFORMITY
    ws = WEIGHTS_COMPOSITE

    fit = wf[0] * max(0, 1 - metrics["mcq"] / MCQ_NORMALIZATION) + wf[1] * max(
        0, 1 - metrics["medcq"] / MCQ_NORMALIZATION
    )

    concentration = (
        wc[0] * metrics["rmcq"]
        + wc[1] * max(0, 1 - metrics["circular_mad_mcq"])
        + wc[2] * (1 - metrics["p_watson_dmcq"])
        + wc[3] * (1 - metrics["p_sim_test_dmcq"])
    )

    uniformity = (
        wu[0] * max(0, 1 - (abs(metrics["mcd"]) / MCD_NORMALIZATION))
        + wu[1] * max(0, 1 - (abs(metrics["medcd"]) / MCD_NORMALIZATION))
        + wu[2] * metrics["rmcd"]
        + wu[3] * max(0, 1 - (metrics["circular_mad_mcd"] / MCD_NORMALIZATION))
        + wu[4] * (1 - metrics["p_rayleigh_dmcd"])
        + wu[5] * (1 - metrics["p_wilcoxon_dmcd"])
    )

    composite = ws[0] * fit + ws[1] * concentration + ws[2] * uniformity

    return composite

evaluate_similarity(target_path, model_path, n_bootstrap=10000, visualize_on=False, output_dir='.')

Evaluate the similarity between two RNA structures.

Parses both structures, computes torsion angle differences, derives circular statistics with bootstrap confidence intervals, and returns a results dictionary including a composite score.

Parameters:

Name	Type	Description	Default
target_path		Path to the reference structure (.pdb or .cif).	required
model_path		Path to the model structure (.pdb or .cif).	required
n_bootstrap		Number of bootstrap resamples for confidence intervals.	10000
visualize_on		If True, save polar plots to output_dir.	False
output_dir		Directory for output files (plots, CSV).	'.'

Returns:

Type	Description
	Dictionary of circular statistics and a score key with the
	composite similarity score in [0, 1].

Source code in src/rnapolis/scorer.py

def evaluate_similarity(
    target_path, model_path, n_bootstrap=10000, visualize_on=False, output_dir="."
):
    """Evaluate the similarity between two RNA structures.

    Parses both structures, computes torsion angle differences, derives
    circular statistics with bootstrap confidence intervals, and returns
    a results dictionary including a composite score.

    Args:
        target_path: Path to the reference structure (``.pdb`` or ``.cif``).
        model_path: Path to the model structure (``.pdb`` or ``.cif``).
        n_bootstrap: Number of bootstrap resamples for confidence intervals.
        visualize_on: If ``True``, save polar plots to ``output_dir``.
        output_dir: Directory for output files (plots, CSV).

    Returns:
        Dictionary of circular statistics and a ``score`` key with the
        composite similarity score in [0, 1].
    """
    target_torsion, model_torsion = _parse_structures(target_path, model_path)

    signed_diffs, unsigned_diffs = compute_mcd_mcq(target_torsion, model_torsion)

    if visualize_on:
        unsigned_diffs_arr = np.array(unsigned_diffs)
        signed_diffs_arr = np.array(signed_diffs)
        visualize(
            unsigned_diffs_arr,
            os.path.join(output_dir, "dmcq.svg"),
        )
        visualize(
            signed_diffs_arr,
            os.path.join(output_dir, "dmcd.svg"),
        )

    results = _compute_statistics(signed_diffs, unsigned_diffs, n_bootstrap)
    results["score"] = compute_score(results)

    return results

main()

CLI entry point for RNA structure scoring.

Source code in src/rnapolis/scorer.py

def main():
    """CLI entry point for RNA structure scoring."""
    parser = argparse.ArgumentParser(
        description="Compare two RNA 3D structures using circular statistics"
    )
    parser.add_argument(
        "--target", type=str, help="path to reference structure", required=True
    )
    parser.add_argument(
        "--model", type=str, help="path to model structure", required=True
    )
    parser.add_argument(
        "--bootstrap-reps",
        type=int,
        help="number of bootstrap resamples",
        required=False,
        default=10000,
    )
    parser.add_argument(
        "--visualize",
        action="store_true",
        help="save polar plots to output directory",
        required=False,
        default=False,
    )
    parser.add_argument(
        "--output-dir",
        type=str,
        help="output directory for CSV and plots",
        required=False,
        default=".",
    )

    args = parser.parse_args()

    results = evaluate_similarity(
        args.target,
        args.model,
        args.bootstrap_reps,
        args.visualize,
        args.output_dir,
    )

    results_list = [[str(key) for key in results.keys()]] + [list(results.values())]
    save_csv(os.path.join(args.output_dir, "result.csv"), results_list)
    print(results["score"])

monte_carlo_mad_test(data, observed_mad, use_full_circle, n_simulations)

Test concentration by comparing observed MAD to random circular samples.

Generates n_simulations random samples from a uniform circular distribution and counts how often their MAD is at most as small as the observed_mad.

Parameters:

Name	Type	Description	Default
data		Original angular data (used only for its length).	required
observed_mad		The observed circular MAD to compare against.	required
use_full_circle		If True, sample from [-pi, pi]; otherwise [0, pi].	required
n_simulations		Number of Monte Carlo iterations.	required

Returns:

Type	Description
	Proportion of random samples with MAD <= observed_mad (p-value).

Source code in src/rnapolis/scorer.py

def monte_carlo_mad_test(data, observed_mad, use_full_circle, n_simulations):
    """Test concentration by comparing observed MAD to random circular samples.

    Generates ``n_simulations`` random samples from a uniform circular
    distribution and counts how often their MAD is at most as small as
    the ``observed_mad``.

    Args:
        data: Original angular data (used only for its length).
        observed_mad: The observed circular MAD to compare against.
        use_full_circle: If ``True``, sample from [-pi, pi]; otherwise [0, pi].
        n_simulations: Number of Monte Carlo iterations.

    Returns:
        Proportion of random samples with MAD <= ``observed_mad`` (p-value).
    """
    lower = -math.pi if use_full_circle else 0
    count = 0
    for _ in range(n_simulations):
        random_sample = [random.uniform(lower, math.pi) for _ in range(len(data))]
        random_mad = circular_mad(
            random_sample, descriptive.circ_median(np.array(random_sample))
        )
        if random_mad <= observed_mad:
            count += 1
    return count / n_simulations

parse_file(filepath)

Parse an RNA structure file in PDB or mmCIF format.

Parameters:

Name	Type	Description	Default
filepath		Path to the input file. Must have a .pdb or .cif extension.	required

Returns:

Type	Description
	Parsed atom data as a DataFrame.

Raises:

Type	Description
ValueError	If the file extension is not .pdb or .cif.

Source code in src/rnapolis/scorer.py

def parse_file(filepath):
    """Parse an RNA structure file in PDB or mmCIF format.

    Args:
        filepath: Path to the input file. Must have a ``.pdb`` or ``.cif`` extension.

    Returns:
        Parsed atom data as a DataFrame.

    Raises:
        ValueError: If the file extension is not ``.pdb`` or ``.cif``.
    """
    if filepath.endswith(".pdb"):
        with open(filepath) as file:
            data = rna_parser.parse_pdb_atoms(file)
    elif filepath.endswith(".cif"):
        with open(filepath) as file:
            data = rna_parser.parse_cif_atoms(file)
    else:
        raise ValueError(
            f"Invalid file format for '{filepath}'. Expected .pdb or .cif extension."
        )
    return data

save_csv(result_file_path, data)

Write rows of data to a CSV file.

Parameters:

Name	Type	Description	Default
result_file_path		Path to the output CSV file.	required
data		List of rows, where each row is a list of values.	required

Source code in src/rnapolis/scorer.py

def save_csv(result_file_path, data):
    """Write rows of data to a CSV file.

    Args:
        result_file_path: Path to the output CSV file.
        data: List of rows, where each row is a list of values.
    """
    with open(result_file_path, "w", newline="") as file:
        writer = csv.writer(file)
        writer.writerows(data)

visualize(d, outfile)

Save a polar plot of angular data to an SVG file.

Plots the circular distribution of angle differences using an adaptive configuration: for highly concentrated data (resultant length > 0.9) the density layer is disabled to avoid an uninformative spike, and the rose diagram uses finer bins.

Input angles are wrapped to [0, 2pi) before plotting so that signed differences (which may be negative) do not confuse the density estimator or scatter positioning.

Parameters:

Name	Type	Description	Default
d		Array of angular values in radians.	required
outfile		Output path for the polar plot (typically .svg).	required

Source code in src/rnapolis/scorer.py

def visualize(d, outfile):
    """Save a polar plot of angular data to an SVG file.

    Plots the circular distribution of angle differences using an adaptive
    configuration: for highly concentrated data (resultant length > 0.9) the
    density layer is disabled to avoid an uninformative spike, and the rose
    diagram uses finer bins.

    Input angles are wrapped to [0, 2pi) before plotting so that signed
    differences (which may be negative) do not confuse the density estimator
    or scatter positioning.

    Args:
        d: Array of angular values in radians.
        outfile: Output path for the polar plot (typically ``.svg``).
    """
    d = d % (2 * np.pi)

    if np.allclose(d, d[0]):
        fig, ax = plt.subplots(subplot_kw={"projection": "polar"}, figsize=(6, 6))
        ax.set_title("All values equal", pad=30)
        ax.plot([0], [1], "o")
        plt.savefig(outfile, format="svg")
        plt.close()
        return

    c = Circular(d, unit="radian")
    config = _build_plot_config(c)

    fig, ax = plt.subplots(
        figsize=(10, 10),
        subplot_kw={"projection": "polar"},
        layout="constrained",
    )
    c.plot(ax, config=config)
    ax.set_title("Complete data", pad=30)
    plt.savefig(outfile, format="svg")
    plt.close()