StollLab / chiLife

Python package for site-directed spin labeling of proteins
https://stolllab.github.io/chiLife/
GNU General Public License v3.0
9 stars 0 forks source link

IntrinsicLabel() with a single atom has no spin_centers set #159

Closed Mishakolok closed 1 month ago

Mishakolok commented 1 month ago

In the IntrinsicLabel.py file, there is a method spin_centers that calculates the spin center of rotamers in an ensemble. The method uses a check if np.any(self.spin_idx): to determine if there are any spin indices present. However, this check fails when spin_idx contains only a single element [0], which can occur when an IntrinsicLabel has only one atom, such as a metal ion.


def __init__(self, res, atom_selection, spin_atoms=None, name='IntrinsicLabel'):

        self._selection = atom_selection
        self._coords = atom_selection.positions.copy()[None, :, :]
        self.coords = self._coords
        self.atom_names = atom_selection.names.copy()
        self.atom_types = atom_selection.types
        self.name = name
        self.label = res
        self.res = res

        self.chain = "".join(set(atom_selection.segids))
        self.weights = np.array([1])
        if isinstance(spin_atoms, (list, tuple, np.ndarray)):
            self.spin_atoms = np.asarray(spin_atoms)
            self.spin_weights = np.ones(len(spin_atoms)) / len(spin_atoms)
        elif isinstance(spin_atoms, dict):
            self.spin_atoms = np.asarray([key for key in spin_atoms])
            self.spin_weights = np.asarray([val for val in spin_atoms.values()])
        elif isinstance(spin_atoms, str):
            self.spin_atoms = np.asarray([spin_atoms])
            self.spin_weights = np.asarray([1])
        else:
            raise RuntimeError("spin_atoms must contain a string, a list/tuple/array of strings, a dict")

        self.site = atom_selection.select_atoms(f'name {self.spin_atoms[np.argmax(self.spin_weights)]}').resnums[0]
        sa_mask = np.isin(self.atom_names, self.spin_atoms,)

        self.spin_idx = np.argwhere(sa_mask)
        self.spin_idx.shape = -1

        self.atoms = [FreeAtom(atom.name, atom.type, idx, atom.resname, atom.resnum, atom.position)
                      for idx, atom in enumerate(self._selection.atoms)]

    @property
    def spin_centers(self):
        """get the spin center of the rotamers in the ensemble"""
        if np.any(self.spin_idx):
            spin_centers = np.average(self._coords[:, self.spin_idx, :], weights=self.spin_weights, axis=1)
        else:
            spin_centers = np.array([])
        return np.atleast_2d(np.squeeze(spin_centers))

If an IntrinsicLabel has only one atom (like a metal ion) the spin_idx becomes just [0], which is fine by itself. However, the np.any in this case returns False (since 0 is False), which leads to the spin_centers set to []. This leads to problems when calculating distance distributions with use_spin_centers = True. I think it could simply fixed with changing the check to something like if len(self.spin_idx) > 0

mtessmer commented 1 month ago

Hi @Mishakolok, Good catch! I will get on this asap. Thank you!