Loading deepchem/feat/__init__.py +1 −0 Original line number Diff line number Diff line Loading @@ -26,3 +26,4 @@ from deepchem.feat.atomic_coordinates import AtomicConvFeaturizer from deepchem.feat.adjacency_fingerprints import AdjacencyFingerprint from deepchem.feat.smiles_featurizers import SmilesToSeq, SmilesToImage from deepchem.feat.materials_featurizers import ElementPropertyFingerprint, SineCoulombMatrix, StructureGraphFeaturizer from deepchem.feat.contact_fingerprints import ContactCircularFingerprint deepchem/feat/atomic_coordinates.py +28 −56 Original line number Diff line number Diff line Loading @@ -227,49 +227,54 @@ class AtomicConvFeaturizer(ComplexFeaturizer): # TODO(rbharath): extend to more fragments if len(frag_num_atoms) != 2: raise ValueError("Currently only supports two fragments") self.frag1_num_atoms = frag_num_atoms[0] self.frag2_num_atoms = frag_num_atoms[0] self.frag_num_atoms = frag_num_atoms self.complex_num_atoms = sum(frag_num_atoms) self.max_num_neighbors = max_num_neighbors self.neighbor_cutoff = neighbor_cutoff self.strip_hydrogens = strip_hydrogens def _featurize(self, mol_pdb_file, protein_pdb_file): def _featurize(self, molecular_complex): """Featurize a single molecular complex. Parameters ---------- molecular_complex: Object Some representation of a molecular complex. """ frag_coords = [] frag_mols = [] try: frag1_coords, frag1_mol = rdkit_util.load_molecule( mol_pdb_file, is_protein=False, sanitize=True, add_hydrogens=False) frag2_coords, frag2_mol = rdkit_util.load_molecule( protein_pdb_file, is_protein=True, sanitize=True, add_hydrogens=False) fragments = rdkit_util.load_complex( molecular_complex, add_hydrogens=False) except MoleculeLoadException: # Currently handles loading failures by returning None # TODO: Is there a better handling procedure? logging.warning("Some molecules cannot be loaded by Rdkit. Skipping") logging.warning("This molecule cannot be loaded by Rdkit. Returning None") return None system_mol = rdkit_util.merge_molecules([frag1_mol, frag2_mol]) mols = [frag[1] for frag in fragments] system_mol = rdkit_util.merge_molecules(mols) system_coords = rdkit_util.get_xyz_from_mol(system_mol) frag1_coords, frag1_mol = self._strip_hydrogens(frag1_coords, frag1_mol) frag2_coords, frag2_mol = self._strip_hydrogens(frag2_coords, frag2_mol) system_coords, system_mol = self._strip_hydrogens(system_coords, system_mol) if self.strip_hydrogens: fragments = [ rdkit_util.strip_hydrogens(frag[0], frag[1]) for frag in fragments ] system_coords, system_mol = rdkit_util.strip_hydrogens( system_coords, system_mol) try: frag1_coords, frag1_neighbor_list, frag1_z = self.featurize_mol( frag1_coords, frag1_mol, self.frag1_num_atoms) frag2_coords, frag2_neighbor_list, frag2_z = self.featurize_mol( frag2_coords, frag2_mol, self.frag2_num_atoms) frag_inputs = [ self.featurize_mol(frag[0], frag[1], frag_num_atoms) for (frag, frag_num_atoms) in zip(fragments, self.frag_num_atoms) ] system_coords, system_neighbor_list, system_z = self.featurize_mol( system_coords, system_mol, self.complex_num_atoms) system_outputs = self.featurize_mol(system_coords, system_mol, self.complex_num_atoms) except ValueError as e: logging.warning( "max_atoms was set too low. Some complexes too large and skipped") return None return frag1_coords, frag1_neighbor_list, frag1_z, frag2_coords, frag2_neighbor_list, frag2_z, \ system_coords, system_neighbor_list, system_z return frag_outputs, system_outputs def get_Z_matrix(self, mol, max_atoms): if len(mol.GetAtoms()) > max_atoms: Loading @@ -287,39 +292,6 @@ class AtomicConvFeaturizer(ComplexFeaturizer): coords = pad_array(coords, (max_num_atoms, 3)) return coords, neighbor_list, z def _strip_hydrogens(self, coords, mol): class MoleculeShim(object): """ Shim of a Molecule which supports #GetAtoms() """ def __init__(self, atoms): self.atoms = [AtomShim(x) for x in atoms] def GetAtoms(self): return self.atoms class AtomShim(object): def __init__(self, atomic_num): self.atomic_num = atomic_num def GetAtomicNum(self): return self.atomic_num if not self.strip_hydrogens: return coords, mol indexes_to_keep = [] atomic_numbers = [] for index, atom in enumerate(mol.GetAtoms()): if atom.GetAtomicNum() != 1: indexes_to_keep.append(index) atomic_numbers.append(atom.GetAtomicNum()) mol = MoleculeShim(atomic_numbers) coords = coords[indexes_to_keep] return coords, mol ############################# Deprecation warning for old name of AtomicConvFeaturizer ############################### Loading deepchem/feat/contact_fingerprints.py +59 −38 Original line number Diff line number Diff line Loading @@ -3,35 +3,40 @@ Topological fingerprints for macromolecular structures. """ import numpy as np import logging import itertools from deepchem.utils.hash_utils import hash_ecfp from deepchem.feat import ComplexFeaturizer from deepchem.utils import rdkit_util from deepchem.utils.hash_utils import vectorize from deepchem.utils.voxel_utils import voxelize from deepchem.utils.voxel_utils import convert_atom_to_voxel from deepchem.utils.rdkit_util import compute_all_ecfp from deepchem.utils.rdkit_util import compute_pairwise_distances from deepchem.utils.rdkit_util import MoleculeLoadException logger = logging.getLogger(__name__) def featurize_binding_pocket_ecfp(protein_xyz, protein, ligand_xyz, ligand, def featurize_contacts_ecfp(frag1, frag2, pairwise_distances=None, cutoff=4.5, ecfp_degree=2): """Computes ECFP dicts for ligand and binding pocket of the protein. """Computes ECFP dicts for pairwise interaction between two molecular fragments. Parameters ---------- protein_xyz: np.ndarray Of shape (N_protein_atoms, 3) protein: rdkit.rdchem.Mol Contains more metadata. ligand_xyz: np.ndarray Of shape (N_ligand_atoms, 3) ligand: rdkit.rdchem.Mol Contains more metadata frag1: Tuple A tuple of (coords, mol) returned by `rdkit_util.load_molecule`. frag2: Tuple A tuple of (coords, mol) returned by `rdkit_util.load_molecule`. #protein_xyz: np.ndarray # Of shape (N_protein_atoms, 3) #protein: rdkit.rdchem.Mol # Contains more metadata. #ligand_xyz: np.ndarray # Of shape (N_ligand_atoms, 3) #ligand: rdkit.rdchem.Mol # Contains more metadata pairwise_distances: np.ndarray Array of pairwise protein-ligand distances (Angstroms) cutoff: float Loading @@ -40,15 +45,26 @@ def featurize_binding_pocket_ecfp(protein_xyz, ECFP radius """ if pairwise_distances is None: pairwise_distances = compute_pairwise_distances(protein_xyz, ligand_xyz) #pairwise_distances = compute_pairwise_distances(protein_xyz, ligand_xyz) pairwise_distances = compute_pairwise_distances(frag1[0], frag2[0]) # contacts is of form (x_coords, y_coords), a tuple of 2 lists contacts = np.nonzero((pairwise_distances < cutoff)) protein_atoms = set([int(c) for c in contacts[0].tolist()]) #protein_atoms = set([int(c) for c in contacts[0].tolist()]) # contacts[0] is the x_coords, that is the frag1 atoms that have # nonzero contact. frag1_atoms = set([int(c) for c in contacts[0].tolist()]) # contacts[1] is the y_coords, the frag2 atoms with nonzero contacts frag2_atoms = set([int(c) for c in contacts[1].tolist()]) protein_ecfp_dict = compute_all_ecfp( protein, indices=protein_atoms, degree=ecfp_degree) ligand_ecfp_dict = compute_all_ecfp(ligand, degree=ecfp_degree) #protein_ecfp_dict = compute_all_ecfp( # protein, indices=protein_atoms, degree=ecfp_degree) frag1_ecfp_dict = compute_all_ecfp( frag1[1], indices=frag1_atoms, degree=ecfp_degree) #ligand_ecfp_dict = compute_all_ecfp(ligand, degree=ecfp_degree) frag2_ecfp_dict = compute_all_ecfp(frag2[1], indices=frag2_atoms, degree=ecfp_degree) return (protein_ecfp_dict, ligand_ecfp_dict) #return (protein_ecfp_dict, ligand_ecfp_dict) return (frag1_ecfp_dict, frag2_ecfp_dict) class ContactCircularFingerprint(ComplexFeaturizer): Loading Loading @@ -84,7 +100,7 @@ class ContactCircularFingerprint(ComplexFeaturizer): self.size = size def _featurize_complex(self, mol, protein): def _featurize_complex(self, molecular_complex): """ Compute featurization for a single mol/protein complex Loading @@ -98,25 +114,30 @@ class ContactCircularFingerprint(ComplexFeaturizer): Parameters ---------- mol: object Representation of the molecule protein: object Representation of the protein molecular_complex: Object Some representation of a molecular complex. """ (lig_xyz, lig_rdk), (prot_xyz, prot_rdk) = mol, protein distances = compute_pairwise_distances(prot_xyz, lig_xyz) return [ vectorize( try: fragments = rdkit_util.load_complex(molecular_complex, add_hydrogens=False) except MoleculeLoadException: logging.warning("This molecule cannot be loaded by Rdkit. Returning None") return None pairwise_features = [] for (frag1, frag2) in itertools.combinations(fragments, 2): # Get coordinates distances = compute_pairwise_distances(frag1[0], frag2[0]) vector = [vectorize( hash_ecfp, feature_dict=ecfp_dict, size=self.size) for ecfp_dict in featurize_binding_pocket_ecfp( prot_xyz, prot_rdk, lig_xyz, lig_rdk, for ecfp_dict in featurize_contacts_ecfp( frag1, frag2, distances, cutoff=self.cutoff, ecfp_degree=self.radius) ] ecfp_degree=self.radius)] pairwise_features.append(vector) return np.array(pairwise_features) class ContactCircularVoxelizer(ComplexFeaturizer): """Computes ECFP fingerprints on a voxel grid. Loading deepchem/trans/tests/test_transformers.py +21 −0 Original line number Diff line number Diff line Loading @@ -238,6 +238,27 @@ class TestTransformers(unittest.TestCase): # Check that untransform does the right thing. np.testing.assert_allclose(normalization_transformer.untransform(y_t), y) def test_normalization_batch_one(self): """This test is based on a user bug report. When `batch_size==1` and `n_tasks > 1`, `NormalizationTransformer.untransform` had a weird edge case where it would return a larger output that input. """ n_samples = 11 n_features = 5 n_tasks = 2 train_dataset = dc.data.NumpyDataset(np.random.rand(n_samples, n_features), np.random.rand(n_samples, n_tasks)) transformer = dc.trans.NormalizationTransformer( transform_y=True, dataset=train_dataset, move_mean=True) # batch of size 1 shape batch = np.random.rand(1, 2, 1) output = transformer.untransform(batch) assert len(batch) == len(output) def test_X_normalization_transformer(self): """Tests normalization transformer.""" solubility_dataset = load_solubility_data() Loading deepchem/trans/transformers.py +2 −2 Original line number Diff line number Diff line Loading @@ -499,8 +499,8 @@ class NormalizationTransformer(Transformer): # Get the reversed shape of z: (..., n_tasks, batch_size) z_shape.reverse() # Find the task dimension of z for dim in z_shape: if dim != n_tasks and dim == 1: for ind, dim in enumerate(z_shape): if ind < (len(z_shape) - 1) and dim == 1: # Prevent broadcasting on wrong dimension y_stds = np.expand_dims(y_stds, -1) y_means = np.expand_dims(y_means, -1) Loading Loading
deepchem/feat/__init__.py +1 −0 Original line number Diff line number Diff line Loading @@ -26,3 +26,4 @@ from deepchem.feat.atomic_coordinates import AtomicConvFeaturizer from deepchem.feat.adjacency_fingerprints import AdjacencyFingerprint from deepchem.feat.smiles_featurizers import SmilesToSeq, SmilesToImage from deepchem.feat.materials_featurizers import ElementPropertyFingerprint, SineCoulombMatrix, StructureGraphFeaturizer from deepchem.feat.contact_fingerprints import ContactCircularFingerprint
deepchem/feat/atomic_coordinates.py +28 −56 Original line number Diff line number Diff line Loading @@ -227,49 +227,54 @@ class AtomicConvFeaturizer(ComplexFeaturizer): # TODO(rbharath): extend to more fragments if len(frag_num_atoms) != 2: raise ValueError("Currently only supports two fragments") self.frag1_num_atoms = frag_num_atoms[0] self.frag2_num_atoms = frag_num_atoms[0] self.frag_num_atoms = frag_num_atoms self.complex_num_atoms = sum(frag_num_atoms) self.max_num_neighbors = max_num_neighbors self.neighbor_cutoff = neighbor_cutoff self.strip_hydrogens = strip_hydrogens def _featurize(self, mol_pdb_file, protein_pdb_file): def _featurize(self, molecular_complex): """Featurize a single molecular complex. Parameters ---------- molecular_complex: Object Some representation of a molecular complex. """ frag_coords = [] frag_mols = [] try: frag1_coords, frag1_mol = rdkit_util.load_molecule( mol_pdb_file, is_protein=False, sanitize=True, add_hydrogens=False) frag2_coords, frag2_mol = rdkit_util.load_molecule( protein_pdb_file, is_protein=True, sanitize=True, add_hydrogens=False) fragments = rdkit_util.load_complex( molecular_complex, add_hydrogens=False) except MoleculeLoadException: # Currently handles loading failures by returning None # TODO: Is there a better handling procedure? logging.warning("Some molecules cannot be loaded by Rdkit. Skipping") logging.warning("This molecule cannot be loaded by Rdkit. Returning None") return None system_mol = rdkit_util.merge_molecules([frag1_mol, frag2_mol]) mols = [frag[1] for frag in fragments] system_mol = rdkit_util.merge_molecules(mols) system_coords = rdkit_util.get_xyz_from_mol(system_mol) frag1_coords, frag1_mol = self._strip_hydrogens(frag1_coords, frag1_mol) frag2_coords, frag2_mol = self._strip_hydrogens(frag2_coords, frag2_mol) system_coords, system_mol = self._strip_hydrogens(system_coords, system_mol) if self.strip_hydrogens: fragments = [ rdkit_util.strip_hydrogens(frag[0], frag[1]) for frag in fragments ] system_coords, system_mol = rdkit_util.strip_hydrogens( system_coords, system_mol) try: frag1_coords, frag1_neighbor_list, frag1_z = self.featurize_mol( frag1_coords, frag1_mol, self.frag1_num_atoms) frag2_coords, frag2_neighbor_list, frag2_z = self.featurize_mol( frag2_coords, frag2_mol, self.frag2_num_atoms) frag_inputs = [ self.featurize_mol(frag[0], frag[1], frag_num_atoms) for (frag, frag_num_atoms) in zip(fragments, self.frag_num_atoms) ] system_coords, system_neighbor_list, system_z = self.featurize_mol( system_coords, system_mol, self.complex_num_atoms) system_outputs = self.featurize_mol(system_coords, system_mol, self.complex_num_atoms) except ValueError as e: logging.warning( "max_atoms was set too low. Some complexes too large and skipped") return None return frag1_coords, frag1_neighbor_list, frag1_z, frag2_coords, frag2_neighbor_list, frag2_z, \ system_coords, system_neighbor_list, system_z return frag_outputs, system_outputs def get_Z_matrix(self, mol, max_atoms): if len(mol.GetAtoms()) > max_atoms: Loading @@ -287,39 +292,6 @@ class AtomicConvFeaturizer(ComplexFeaturizer): coords = pad_array(coords, (max_num_atoms, 3)) return coords, neighbor_list, z def _strip_hydrogens(self, coords, mol): class MoleculeShim(object): """ Shim of a Molecule which supports #GetAtoms() """ def __init__(self, atoms): self.atoms = [AtomShim(x) for x in atoms] def GetAtoms(self): return self.atoms class AtomShim(object): def __init__(self, atomic_num): self.atomic_num = atomic_num def GetAtomicNum(self): return self.atomic_num if not self.strip_hydrogens: return coords, mol indexes_to_keep = [] atomic_numbers = [] for index, atom in enumerate(mol.GetAtoms()): if atom.GetAtomicNum() != 1: indexes_to_keep.append(index) atomic_numbers.append(atom.GetAtomicNum()) mol = MoleculeShim(atomic_numbers) coords = coords[indexes_to_keep] return coords, mol ############################# Deprecation warning for old name of AtomicConvFeaturizer ############################### Loading
deepchem/feat/contact_fingerprints.py +59 −38 Original line number Diff line number Diff line Loading @@ -3,35 +3,40 @@ Topological fingerprints for macromolecular structures. """ import numpy as np import logging import itertools from deepchem.utils.hash_utils import hash_ecfp from deepchem.feat import ComplexFeaturizer from deepchem.utils import rdkit_util from deepchem.utils.hash_utils import vectorize from deepchem.utils.voxel_utils import voxelize from deepchem.utils.voxel_utils import convert_atom_to_voxel from deepchem.utils.rdkit_util import compute_all_ecfp from deepchem.utils.rdkit_util import compute_pairwise_distances from deepchem.utils.rdkit_util import MoleculeLoadException logger = logging.getLogger(__name__) def featurize_binding_pocket_ecfp(protein_xyz, protein, ligand_xyz, ligand, def featurize_contacts_ecfp(frag1, frag2, pairwise_distances=None, cutoff=4.5, ecfp_degree=2): """Computes ECFP dicts for ligand and binding pocket of the protein. """Computes ECFP dicts for pairwise interaction between two molecular fragments. Parameters ---------- protein_xyz: np.ndarray Of shape (N_protein_atoms, 3) protein: rdkit.rdchem.Mol Contains more metadata. ligand_xyz: np.ndarray Of shape (N_ligand_atoms, 3) ligand: rdkit.rdchem.Mol Contains more metadata frag1: Tuple A tuple of (coords, mol) returned by `rdkit_util.load_molecule`. frag2: Tuple A tuple of (coords, mol) returned by `rdkit_util.load_molecule`. #protein_xyz: np.ndarray # Of shape (N_protein_atoms, 3) #protein: rdkit.rdchem.Mol # Contains more metadata. #ligand_xyz: np.ndarray # Of shape (N_ligand_atoms, 3) #ligand: rdkit.rdchem.Mol # Contains more metadata pairwise_distances: np.ndarray Array of pairwise protein-ligand distances (Angstroms) cutoff: float Loading @@ -40,15 +45,26 @@ def featurize_binding_pocket_ecfp(protein_xyz, ECFP radius """ if pairwise_distances is None: pairwise_distances = compute_pairwise_distances(protein_xyz, ligand_xyz) #pairwise_distances = compute_pairwise_distances(protein_xyz, ligand_xyz) pairwise_distances = compute_pairwise_distances(frag1[0], frag2[0]) # contacts is of form (x_coords, y_coords), a tuple of 2 lists contacts = np.nonzero((pairwise_distances < cutoff)) protein_atoms = set([int(c) for c in contacts[0].tolist()]) #protein_atoms = set([int(c) for c in contacts[0].tolist()]) # contacts[0] is the x_coords, that is the frag1 atoms that have # nonzero contact. frag1_atoms = set([int(c) for c in contacts[0].tolist()]) # contacts[1] is the y_coords, the frag2 atoms with nonzero contacts frag2_atoms = set([int(c) for c in contacts[1].tolist()]) protein_ecfp_dict = compute_all_ecfp( protein, indices=protein_atoms, degree=ecfp_degree) ligand_ecfp_dict = compute_all_ecfp(ligand, degree=ecfp_degree) #protein_ecfp_dict = compute_all_ecfp( # protein, indices=protein_atoms, degree=ecfp_degree) frag1_ecfp_dict = compute_all_ecfp( frag1[1], indices=frag1_atoms, degree=ecfp_degree) #ligand_ecfp_dict = compute_all_ecfp(ligand, degree=ecfp_degree) frag2_ecfp_dict = compute_all_ecfp(frag2[1], indices=frag2_atoms, degree=ecfp_degree) return (protein_ecfp_dict, ligand_ecfp_dict) #return (protein_ecfp_dict, ligand_ecfp_dict) return (frag1_ecfp_dict, frag2_ecfp_dict) class ContactCircularFingerprint(ComplexFeaturizer): Loading Loading @@ -84,7 +100,7 @@ class ContactCircularFingerprint(ComplexFeaturizer): self.size = size def _featurize_complex(self, mol, protein): def _featurize_complex(self, molecular_complex): """ Compute featurization for a single mol/protein complex Loading @@ -98,25 +114,30 @@ class ContactCircularFingerprint(ComplexFeaturizer): Parameters ---------- mol: object Representation of the molecule protein: object Representation of the protein molecular_complex: Object Some representation of a molecular complex. """ (lig_xyz, lig_rdk), (prot_xyz, prot_rdk) = mol, protein distances = compute_pairwise_distances(prot_xyz, lig_xyz) return [ vectorize( try: fragments = rdkit_util.load_complex(molecular_complex, add_hydrogens=False) except MoleculeLoadException: logging.warning("This molecule cannot be loaded by Rdkit. Returning None") return None pairwise_features = [] for (frag1, frag2) in itertools.combinations(fragments, 2): # Get coordinates distances = compute_pairwise_distances(frag1[0], frag2[0]) vector = [vectorize( hash_ecfp, feature_dict=ecfp_dict, size=self.size) for ecfp_dict in featurize_binding_pocket_ecfp( prot_xyz, prot_rdk, lig_xyz, lig_rdk, for ecfp_dict in featurize_contacts_ecfp( frag1, frag2, distances, cutoff=self.cutoff, ecfp_degree=self.radius) ] ecfp_degree=self.radius)] pairwise_features.append(vector) return np.array(pairwise_features) class ContactCircularVoxelizer(ComplexFeaturizer): """Computes ECFP fingerprints on a voxel grid. Loading
deepchem/trans/tests/test_transformers.py +21 −0 Original line number Diff line number Diff line Loading @@ -238,6 +238,27 @@ class TestTransformers(unittest.TestCase): # Check that untransform does the right thing. np.testing.assert_allclose(normalization_transformer.untransform(y_t), y) def test_normalization_batch_one(self): """This test is based on a user bug report. When `batch_size==1` and `n_tasks > 1`, `NormalizationTransformer.untransform` had a weird edge case where it would return a larger output that input. """ n_samples = 11 n_features = 5 n_tasks = 2 train_dataset = dc.data.NumpyDataset(np.random.rand(n_samples, n_features), np.random.rand(n_samples, n_tasks)) transformer = dc.trans.NormalizationTransformer( transform_y=True, dataset=train_dataset, move_mean=True) # batch of size 1 shape batch = np.random.rand(1, 2, 1) output = transformer.untransform(batch) assert len(batch) == len(output) def test_X_normalization_transformer(self): """Tests normalization transformer.""" solubility_dataset = load_solubility_data() Loading
deepchem/trans/transformers.py +2 −2 Original line number Diff line number Diff line Loading @@ -499,8 +499,8 @@ class NormalizationTransformer(Transformer): # Get the reversed shape of z: (..., n_tasks, batch_size) z_shape.reverse() # Find the task dimension of z for dim in z_shape: if dim != n_tasks and dim == 1: for ind, dim in enumerate(z_shape): if ind < (len(z_shape) - 1) and dim == 1: # Prevent broadcasting on wrong dimension y_stds = np.expand_dims(y_stds, -1) y_means = np.expand_dims(y_means, -1) Loading
