Clashfinder

AtomType

Bases: Enum

Simple atom type enum used to map atoms to approximate vdW radii.

Source code in src/rnapolis/clashfinder.py

class AtomType(Enum):
    """Simple atom type enum used to map atoms to approximate vdW radii."""

    C = "C"
    N = "N"
    O = "O"
    P = "P"

    @cached_property
    def radius(self) -> float:
        """Return an approximate van der Waals radius for this atom type."""
        if self.value == "C":
            return CARBON_RADIUS
        elif self.value == "N":
            return NITROGEN_RADIUS
        elif self.value == "O":
            return OXYGEN_RADIUS
        elif self.value == "P":
            return PHOSPHORUS_RADIUS
        raise RuntimeError(f"Unknown atom type: {self}")

    def matches(self, atom: Atom) -> bool:
        """
        Check whether a given atom belongs to this atom type.

        Args:
            atom: Atom to test.

        Returns:
            ``True`` if the atom name starts with this atom type symbol,
            ``False`` otherwise.
        """
        return atom.name.strip().startswith(self.value)

radius cached property

Return an approximate van der Waals radius for this atom type.

matches(atom)

Check whether a given atom belongs to this atom type.

Parameters:

Name	Type	Description	Default
atom	Atom	Atom to test.	required

Returns:

Type	Description
bool	True if the atom name starts with this atom type symbol,
bool	False otherwise.

Source code in src/rnapolis/clashfinder.py

def matches(self, atom: Atom) -> bool:
    """
    Check whether a given atom belongs to this atom type.

    Args:
        atom: Atom to test.

    Returns:
        ``True`` if the atom name starts with this atom type symbol,
        ``False`` otherwise.
    """
    return atom.name.strip().startswith(self.value)

classify_clash(atom_i, atom_j, occupancy)

Classify a clash into a simple text label.

Currently this function detects O3' clashes with phosphate atoms and labels them as "O3'"; all other clashes are left unclassified.

Parameters:

Name	Type	Description	Default
atom_i	Atom	First atom in the clashing pair.	required
atom_j	Atom	Second atom in the clashing pair.	required
occupancy	float	Sum of occupancies for the two atoms.	required

Returns:

Type	Description
Optional[str]	Short label for the clash type, or None if it does not match any known category.

Source code in src/rnapolis/clashfinder.py

def classify_clash(atom_i: Atom, atom_j: Atom, occupancy: float) -> Optional[str]:
    """
    Classify a clash into a simple text label.

    Currently this function detects O3' clashes with phosphate atoms
    and labels them as ``"O3'"``; all other clashes are left unclassified.

    Args:
        atom_i: First atom in the clashing pair.
        atom_j: Second atom in the clashing pair.
        occupancy: Sum of occupancies for the two atoms.

    Returns:
        Short label for the clash type, or ``None`` if it does not match any known category.
    """
    if atom_i.name == "O3'" and atom_j.name in (
        "OP1",
        "OP2",
        "OP3",
        "O1P",
        "O2P",
        "O3P",
    ):
        return "O3'"
    return None

find_clashes(residues, ignore_occupancy, ignore_autoclashes, nucleic_acid_only, require_same_atom_name, enable_molprobity_mode)

Find steric clashes between atoms in a list of residues.

The function collects carbon, nitrogen, oxygen and phosphorus atoms, builds a KDTree and looks for pairs that are closer than the sum of their vdW radii (optionally using a MolProbity-like cutoff).

Parameters:

Name	Type	Description	Default
residues	List[Residue3D]	List of 3D residues to analyse.	required
ignore_occupancy	bool	If True, ignore atom occupancies when reporting clashes.	required
ignore_autoclashes	bool	If True, skip clashes within the same residue.	required
nucleic_acid_only	bool	If True, consider only residues marked as nucleotides.	required
require_same_atom_name	bool	If True, report only clashes of atoms with the same name.	required
enable_molprobity_mode	bool	If True, add a 0.5 Å margin as in MolProbity.	required

Returns:

Name	Type	Description
list	list[tuple[tuple[Residue3D, Atom], tuple[Residue3D, Atom], float]]	A list of clashes. Each clash is represented as a tuple: ((res_i, atom_i), (res_j, atom_j), occupancy_sum), where: res_i, res_j (Residue3D): clashing residues, atom_i, atom_j (Atom): clashing atoms, occupancy_sum (float): sum of occupancies for the atom pair.

Source code in src/rnapolis/clashfinder.py

def find_clashes(
    residues: List[Residue3D],
    ignore_occupancy: bool,
    ignore_autoclashes: bool,
    nucleic_acid_only: bool,
    require_same_atom_name: bool,
    enable_molprobity_mode: bool,
) -> list[tuple[tuple[Residue3D, Atom], tuple[Residue3D, Atom], float]]:
    """
    Find steric clashes between atoms in a list of residues.

    The function collects carbon, nitrogen, oxygen and phosphorus atoms,
    builds a KDTree and looks for pairs that are closer than the sum of
    their vdW radii (optionally using a MolProbity-like cutoff).

    Args:
        residues: List of 3D residues to analyse.
        ignore_occupancy: If True, ignore atom occupancies when reporting clashes.
        ignore_autoclashes: If True, skip clashes within the same residue.
        nucleic_acid_only: If True, consider only residues marked as nucleotides.
        require_same_atom_name: If True, report only clashes of atoms with the same name.
        enable_molprobity_mode: If True, add a 0.5 Å margin as in MolProbity.

    Returns:
        list:
            A list of clashes. Each clash is represented as a tuple:

            - **((res_i, atom_i), (res_j, atom_j), occupancy_sum)**, where:
              - **res_i**, **res_j** (Residue3D): clashing residues,
              - **atom_i**, **atom_j** (Atom): clashing atoms,
              - **occupancy_sum** (float): sum of occupancies for the atom pair.
    """

    reference_residues = []
    reference_atoms = []
    coordinates = []

    for residue in residues:
        if (
            nucleic_acid_only is True and residue.is_nucleotide
        ) or nucleic_acid_only is False:
            for atom in residue.atoms:
                if any([atom_type.matches(atom) for atom_type in AtomType]):
                    reference_residues.append(residue)
                    reference_atoms.append(atom)
                    coordinates.append(atom.coordinates)

    if len(coordinates) < 2:
        return []

    kdtree = KDTree(coordinates)
    result = []
    max_radius = max([atom_type.radius for atom_type in AtomType])
    molprobity_factor = 0.5 if enable_molprobity_mode is True else 0.0

    for i, j in kdtree.query_pairs(2.0 * max_radius + molprobity_factor):
        ai: Atom = reference_atoms[i]
        aj: Atom = reference_atoms[j]
        ri, rj = reference_residues[i], reference_residues[j]

        if ignore_autoclashes is True and ri == rj:
            continue
        if require_same_atom_name is True and ai.name != aj.name:
            continue

        distance = np.linalg.norm(ai.coordinates - aj.coordinates)
        sum_vdw_radii = AtomType[ai.name[0]].radius + AtomType[aj.name[0]].radius
        if distance > sum_vdw_radii + molprobity_factor:
            continue

        sum_occupancies = (ai.occupancy or 1.0) + (aj.occupancy or 1.0)
        if ignore_occupancy is True or math.isclose(sum_occupancies, 1.0):
            result.append(((ri, ai), (rj, aj), sum_occupancies))

    return result

main()

Command-line entry point for the clash finder tool.

It reads a PDB or mmCIF file, detects atomic clashes using the chosen options, prints a human-readable summary and optionally writes a CSV report with detailed clash information.

Source code in src/rnapolis/clashfinder.py

def main():
    """
    Command-line entry point for the clash finder tool.

    It reads a PDB or mmCIF file, detects atomic clashes using the
    chosen options, prints a human-readable summary and optionally
    writes a CSV report with detailed clash information.
    """
    parser = argparse.ArgumentParser()
    parser.add_argument("input", help="Path to PDB or mmCIF file")
    parser.add_argument(
        "--ignore-occupancy",
        help="By default clashes are reported if atoms' occupancies are not equal to 1.0, but you can ignore this check. If you ignore this check, any pair of atoms too close to each other will be reported regardless of their occupancy",
        action="store_true",
    )
    parser.add_argument(
        "--nucleic-acid-only",
        help="By default all kind of clashes will be found, but you can focus only on nucleic acids chains",
        action="store_true",
    )
    parser.add_argument(
        "--ignore-autoclashes",
        help="By default clashes will be reported even in scope of the same residue, but you can disable this behaviour",
        action="store_true",
    )
    parser.add_argument(
        "--require-same-atom-name",
        help="By default any two clashing atoms are reported (e.g. P vs OP1), but when this is set, the program will report only clashes of the same atom name (e.g. OP1 vs OP1)",
        action="store_true",
    )
    parser.add_argument(
        "--enable-molprobity-mode",
        help="By default this tool will report any *strong* clash, i.e., when two atoms are closer than their sum of vdW radii; when this option is set, additional 0.5A is added as in MolProbity, so that more clashes are detected, but some of them might be *weak*",
        action="store_true",
    )
    parser.add_argument("--csv", help="Store result in CSV format")
    args = parser.parse_args()

    with open(args.input) as f:
        structure3d = read_3d_structure(f, 1)

    clashing_chains = {}
    max_occupancy_residues = {}
    max_occupancy_chains = {}

    clashes = find_clashes(
        structure3d.residues,
        args.ignore_occupancy,
        args.ignore_autoclashes,
        args.nucleic_acid_only,
        args.require_same_atom_name,
        args.enable_molprobity_mode,
    )

    if clashes:
        for pi, pj, occupancy in clashes:
            ri, ai = pi
            rj, aj = pj

            chain_key = (ri.chain, rj.chain)
            residue_key = (ri, rj)
            if chain_key not in clashing_chains:
                clashing_chains[chain_key] = {}
            if residue_key not in clashing_chains[chain_key]:
                clashing_chains[chain_key][residue_key] = set()
            clashing_chains[chain_key][residue_key].add((ai, aj, occupancy))

            max_occupancy_residues[(ri, rj)] = max(
                [max_occupancy_residues.get((ri, rj), 0.0), occupancy]
            )
            max_occupancy_chains[(ri.chain, rj.chain)] = max(
                [max_occupancy_residues.get((ri.chain, rj.chain), 0.0), occupancy]
            )

    if clashing_chains:
        for ci, cj in sorted(clashing_chains):
            if ci == cj:
                print(
                    f"Clashes found in chain {ci} with maximum occupancy sum equal to {max_occupancy_chains[(ci, cj)]}"
                )
            else:
                print(
                    f"Clashes found between chains {ci} and {cj} with maximum occupancy sum equal to {max_occupancy_chains[(ci, cj)]}"
                )
            for ri, rj in clashing_chains[(ci, cj)]:
                if ri == rj:
                    print(
                        f"    Clashes found in residue {ri} with maximum occupancy sum equal to {max_occupancy_residues[(ri, rj)]}"
                    )
                else:
                    print(
                        f"    Clashes found between residues {ri} and {rj} with maximum occupancy sum equal to {max_occupancy_residues[(ri, rj)]}"
                    )
                for ai, aj, occupancy in sorted(clashing_chains[(ci, cj)][(ri, rj)]):
                    print(
                        f"        Clashes found between atoms {ai.name} and {aj.name} with occupancy sum of {occupancy}"
                    )

        if args.csv:
            metadata = read_metadata(args.input, ["exptl", "refine"])

            with open(args.csv, "w") as f:
                writer = csv.writer(f)
                writer.writerow(
                    [
                        "Filename",
                        "Experimental method",
                        "Resolution",
                        "Atom 1",
                        "Atom 2",
                        "Occupancy sum",
                        "Classification",
                    ]
                )

                for ci, cj in sorted(clashing_chains):
                    for ri, rj in clashing_chains[(ci, cj)]:
                        for ai, aj, occupancy in sorted(
                            clashing_chains[(ci, cj)][(ri, rj)]
                        ):
                            writer.writerow(
                                [
                                    f"{os.path.splitext(os.path.basename(args.input))[0]}",
                                    metadata["exptl"][0]["method"],
                                    metadata["refine"][0]["ls_d_res_high"],
                                    f"{ri} {ai.name}",
                                    f"{rj} {aj.name}",
                                    occupancy,
                                    classify_clash(ai, aj, occupancy),
                                ]
                            )