Relaxation

See also deposited-data accessors on NMRStarEntry (relaxation, order_parameters).

RelaxationProfile

RelaxationProfile

RelaxationProfile(table, sequence, entry_id=None, entity_id=None, entry=None)

Per-residue R1/R2/hetNOE for one entity, aligned to its sequence, with RelaxDB motional labels.

Attributes:

Name	Type	Description
table	DataFrame	One row per sequence position (1-indexed). Columns: Seq_ID, residue, R1, R1_err, R2, R2_err, NOE, NOE_err, R2_R1, R2_R1_err, has_data (plus scaled_R2_R1_pred and label once those steps are run).
sequence	str	One-letter sequence the data is aligned to.
entry_id, entity_id	identifiers carried for reference.

Source code in makeshift/relaxation/relax_profile.py

def __init__(self, table, sequence, entry_id=None, entity_id=None, entry=None):
    self.table = table
    self.sequence = sequence
    self.entry_id = entry_id
    self.entity_id = entity_id
    self.entry = entry
    self.scale_factor = None

from_bmrb classmethod

from_bmrb(bmrb_id, entity_id=None, sequence=None, peaklist=None, **fetch_kw)

Fetch a BMRB entry and build a profile from its deposited relaxation.

Source code in makeshift/relaxation/relax_profile.py

@classmethod
def from_bmrb(cls, bmrb_id, entity_id=None, sequence=None, peaklist=None,
              **fetch_kw):
    """Fetch a BMRB entry and build a profile from its deposited relaxation."""
    entry = NMRStarEntry.from_bmrb(bmrb_id, **fetch_kw)
    return cls.from_entry(entry, entity_id=entity_id, sequence=sequence,
                          peaklist=peaklist)

from_entry classmethod

from_entry(entry, entity_id=None, sequence=None, peaklist=None)

Build from an already-parsed NMRStarEntry. Pulls R1 (from T1), R2 (from T2), and hetNOE, converting times to rates using each list's units tag.

peaklist (a PeakList, a DataFrame with a Seq_ID column, or an iterable of Seq_IDs) marks which residues have an assigned amide peak; if None, the entry's own backbone amide peaks are derived via PeakList. Residues with no assigned peak are flagged "missing" (.); residues with a peak but no relaxation data are flagged "no data" (x).

Source code in makeshift/relaxation/relax_profile.py

@classmethod
def from_entry(cls, entry, entity_id=None, sequence=None, peaklist=None):
    """
    Build from an already-parsed NMRStarEntry.
    Pulls R1 (from T1), R2 (from T2), and hetNOE, converting times to rates
    using each list's units tag.

    `peaklist` (a PeakList, a DataFrame with a Seq_ID column, or an iterable
    of Seq_IDs) marks which residues have an assigned amide peak; if None,
    the entry's own backbone amide peaks are derived via PeakList. Residues
    with no assigned peak are flagged "missing" (.); residues with a peak but
    no relaxation data are flagged "no data" (x).
    """
    if sequence is None:
        seqs = entry.sequences()
        if entity_id is not None:
            sequence = entry.sequences(entity_id=entity_id)
        elif len(seqs):
            sequence = seqs["Polymer_seq_one_letter_code"].iloc[0]
            entity_id = seqs["ID"].iloc[0]
        if not isinstance(sequence, str) or not sequence:
            raise ValueError(
                "could not resolve a sequence; pass sequence=... explicitly"
            )

    n = len(sequence)
    table = pd.DataFrame({
        "Seq_ID": np.arange(1, n + 1),
        "residue": list(sequence),
    })

    for kind, col in (("T1", "R1"), ("T2", "R2")):
        df = entry.relaxation(kind)
        units = cls._list_units(entry, kind)
        table[col], table[f"{col}_err"] = cls._align_rate(
            df, n, units, kind, entry_id=entry.entry_id)

    noe = entry.relaxation("NOE")
    table["NOE"], table["NOE_err"] = cls._align_plain(noe, n)

    table["R2_R1"] = table["R2"] / table["R1"]
    table["R2_R1_err"] = table["R2_R1"] * np.sqrt(
        (table["R2_err"] / table["R2"]) ** 2
        + (table["R1_err"] / table["R1"]) ** 2)
    table["has_data"] = table[["R1", "R2"]].notna().all(axis=1)

    if peaklist is None:
        peaklist = cls._entry_peaklist(entry, entity_id)
    present = cls._peaklist_seqids(peaklist)
    table["has_HN"] = table["Seq_ID"].isin(present)

    return cls(table, sequence,
               entry_id=getattr(entry, "entry_id", None),
               entity_id=entity_id, entry=entry)

add_rigid_prediction

add_rigid_prediction(pdb=None, source='auto', config=None, chain=None, noe_cut=0.65, **fetch_kw)

Run HYDRONMR on a structure and scale its rigid R2/R1 (T1_over_T2) to the observed data, so elevated R2/R1 stands out as exchange.

pdb may be a local file, a 4-character PDB id (fetched from RCSB), or a UniProt accession (fetched from AlphaFold DB) — pass source= to force one.

If pdb is None, the entry's own deposited PDB code is used when it cites one; otherwise this raises (makeshift does not predict structure).

The scale factor is fit by least squares on ordered residues (hetNOE above noe_cut where available), mirroring classify.fit_R2_rigid. Adds scaled_R2_R1_pred and NOE_pred; residues outside the modeled region keep NaN.

Source code in makeshift/relaxation/relax_profile.py

def add_rigid_prediction(self, pdb=None, source="auto", config=None,
                         chain=None, noe_cut=0.65, **fetch_kw):
    """
    Run HYDRONMR on a structure and scale its rigid R2/R1 (T1_over_T2) to the
    observed data, so elevated R2/R1 stands out as exchange.

    `pdb` may be a local file, a 4-character PDB id (fetched from RCSB), or a
    UniProt accession (fetched from AlphaFold DB) — pass `source=` to force
    one. 

    If `pdb` is None, the entry's own deposited PDB code is used when it
    cites one; otherwise this raises (makeshift does not predict structure).

    The scale factor is fit by least squares on ordered residues (hetNOE
    above `noe_cut` where available), mirroring classify.fit_R2_rigid. Adds
    `scaled_R2_R1_pred` and `NOE_pred`; residues outside the modeled region
    keep NaN.
    """
    from ..hydronmr import run as run_hydronmr
    from ..utils.structures import fetch_structure

    if pdb is None:
        pdb_ids = self.entry.get_pdb_ids() if self.entry is not None else []
        af_ids = self.entry.get_alphafold_ids() if self.entry is not None else []
        if source == "rcsb":
            if not pdb_ids:
                raise ValueError("entry cites no PDB; pass pdb=<PDB id | path>")
            pdb = pdb_ids[0]
        elif source == "afdb":
            if not af_ids:
                raise ValueError("entry cites no AlphaFold/UniProt accession; "
                                 "pass pdb=<UniProt accession | path>")
            pdb = af_ids[0]
        elif pdb_ids:                    # source == "auto": prefer deposited PDB
            pdb, source = pdb_ids[0], "rcsb"
        elif af_ids:                     # else fall back to AlphaFold
            pdb, source = af_ids[0], "afdb"
        else:
            raise ValueError(
                "no structure given and the entry cites no PDB or "
                "AlphaFold/UniProt accession; pass pdb=<path | PDB id | "
                "UniProt accession> (experimental or predicted) to enable "
                "exchange labeling"
            )
        print(f"  no pdb given; using {source} structure {pdb}")

    pdb_path = fetch_structure(pdb, source=source, **fetch_kw)
    result = run_hydronmr(pdb_path, config_path=config) if config \
        else run_hydronmr(pdb_path)
    hydro = result.to_dataframe()
    if chain is not None:
        hydro = hydro[hydro["chain"] == chain]
    hydro = hydro.rename(columns={"seqpos": "Seq_ID"})

    t = self.table.merge(
        hydro[["Seq_ID", "T1_over_T2", "NOE"]].rename(
            columns={"T1_over_T2": "_pred_ratio", "NOE": "NOE_pred"}),
        on="Seq_ID", how="left")

    ordered = (t["_pred_ratio"].notna() & t["R2_R1"].notna()
               & ((t["NOE"] > noe_cut) | t["NOE"].isna()))
    n_match = int((t["_pred_ratio"].notna() & t["R2_R1"].notna()).sum())
    pred = t.loc[ordered, "_pred_ratio"]
    obs = t.loc[ordered, "R2_R1"]
    self.scale_factor = float((obs * pred).sum() / (pred ** 2).sum())

    t["scaled_R2_R1_pred"] = self.scale_factor * t["_pred_ratio"]
    t = t.drop(columns="_pred_ratio")
    self.table = t
    print(f"  HYDRONMR: {n_match} residues matched structure to data, "
          f"scale factor {self.scale_factor:.3f} "
          f"({int(ordered.sum())} ordered residues used)")
    return self

label

label(rex_n_std=1.0, noe_cut=0.65)

Assign a label token to every residue and return the label string.

Requires add_rigid_prediction first for the exchange (^) call, which flags residues whose R2/R1 exceeds the rigid prediction by more than rex_n_std standard deviations of that excess across modeled residues. ps-ns motion (v) is hetNOE <= noe_cut.

Source code in makeshift/relaxation/relax_profile.py

def label(self, rex_n_std=1.0, noe_cut=0.65):
    """
    Assign a label token to every residue and return the label string.

    Requires `add_rigid_prediction` first for the exchange (`^`) call, which
    flags residues whose R2/R1 exceeds the rigid prediction by more than
    `rex_n_std` standard deviations of that excess across modeled residues.
    ps-ns motion (`v`) is hetNOE <= `noe_cut`.
    """
    t = self.table
    have_pred = "scaled_R2_R1_pred" in t.columns

    excess = pd.Series(np.nan, index=t.index)
    rex_mask = pd.Series(False, index=t.index)
    if have_pred:
        modeled = t["scaled_R2_R1_pred"].notna() & t["R2_R1"].notna()
        excess[modeled] = t.loc[modeled, "R2_R1"] - t.loc[modeled, "scaled_R2_R1_pred"]
        thresh = excess[modeled].mean() + rex_n_std * excess[modeled].std()
        rex_mask = excess > thresh
    else:
        warnings.warn(
            "no rigid prediction set; exchange (^) cannot be called. "
            "trying to run add_rigid_prediction(pdb) first with default "
            "parameters .", UserWarning)

    psns_mask = t["NOE"] <= noe_cut

    labels = []
    for i, row in t.iterrows():
        if row["residue"] == "P":
            labels.append(PROLINE)
        elif not row["has_data"]:
            if row.get("has_HN", False):
                labels.append(NODATA)
            else:
                labels.append(MISSING)
        else:
            rex = bool(rex_mask.get(i, False))
            psns = bool(psns_mask.get(i, False))
            labels.append(BOTH if (rex and psns) else
                          REX if rex else
                          PSNS if psns else ORDERED)

    t["label"] = labels
    self.table = t
    return "".join(labels)

plot

plot(data_type='R2_R1', ax=None, figsize=(6, 1.5))

Plot a relaxation observable along the sequence with motion labels: orange = exchange (^/b) blue = ps-ns (v) black = ordered (A) purple P = proline red star = missing gray = no data The scaled rigid prediction is overlaid for R2_R1. Requires label() first.

Source code in makeshift/relaxation/relax_profile.py

def plot(self, data_type="R2_R1", ax=None, figsize=(6, 1.5)):
    """
    Plot a relaxation observable along the sequence with 
    motion labels:
        orange = exchange (^/b)
        blue = ps-ns (v)
        black = ordered (A)
        purple P = proline
        red star = missing 
        gray = no data
    The scaled rigid prediction is overlaid for R2_R1. 
    Requires `label()` first.
    """
    import matplotlib.pyplot as plt

    t = self.table
    if "label" not in t.columns:
        raise ValueError("call label() before plot()")
    if data_type not in t.columns:
        raise ValueError(f"no column {data_type!r} in table")

    if ax is None:
        _, ax = plt.subplots(figsize=figsize)
    x = t["Seq_ID"].values
    y = t[data_type].values
    ax.plot(x, y, color="black", lw=0.5, zorder=5)

    if data_type == "R2_R1" and "scaled_R2_R1_pred" in t.columns:
        ax.plot(x, t["scaled_R2_R1_pred"].values, color="grey", zorder=5)

    ymin, ymax = ax.get_ylim()
    p_pos = ymin + 0.05 * (ymax - ymin)
    star_pos = ymin + 0.9 * (ymax - ymin)
    err = t.get(f"{data_type}_err")

    colors = {REX: "tab:orange", BOTH: "tab:orange",
              PSNS: "tab:blue", ORDERED: "black", TERMINUS: "grey"}
    for _, row in t.iterrows():
        j, lab = row["Seq_ID"], row["label"]
        if lab == PROLINE:
            ax.axvline(j, color="tab:purple", alpha=0.5, lw=0.5)
            ax.text(j - 0.5, p_pos, "P", color="tab:purple",
                    fontsize=5, weight="bold")
        elif lab == MISSING:
            ax.axvline(j, color="tab:red", alpha=0.5, lw=0.5)
            ax.scatter([j], [star_pos], marker="*", color="tab:red")
        elif lab == NODATA:
            ax.axvline(j, color="tab:grey", alpha=0.5, lw=0.5)
        elif lab in colors and lab != TERMINUS and pd.notna(row[data_type]):
            e = err.loc[row.name] if err is not None else None
            ax.errorbar(j, row[data_type], yerr=e, fmt=".",
                        color=colors[lab], zorder=10)
    ax.set_xlim(0, len(self.sequence) + 1)
    ax.set_xlabel("Residue")
    return ax