Loading deepchem/feat/__init__.py +9 −0 Original line number Diff line number Diff line """ Making it easy to import in classes. """ # flake8: noqa # base classes for featurizers from deepchem.feat.base_classes import Featurizer from deepchem.feat.base_classes import MolecularFeaturizer from deepchem.feat.base_classes import MaterialStructureFeaturizer from deepchem.feat.base_classes import MaterialCompositionFeaturizer from deepchem.feat.base_classes import ComplexFeaturizer from deepchem.feat.base_classes import UserDefinedFeaturizer from deepchem.feat.graph_features import ConvMolFeaturizer from deepchem.feat.graph_features import WeaveFeaturizer from deepchem.feat.fingerprints import CircularFingerprint Loading @@ -22,6 +26,11 @@ from deepchem.feat.atomic_coordinates import AtomicCoordinates from deepchem.feat.atomic_coordinates import NeighborListComplexAtomicCoordinates from deepchem.feat.adjacency_fingerprints import AdjacencyFingerprint from deepchem.feat.smiles_featurizers import SmilesToSeq, SmilesToImage # molecule featurizers from deepchem.feat.molecule_featurizers import MolGraphConvFeaturizer # material featurizers from deepchem.feat.material_featurizers import ElementPropertyFingerprint from deepchem.feat.material_featurizers import SineCoulombMatrix from deepchem.feat.material_featurizers import CGCNNFeaturizer deepchem/feat/molecule_featurizers/__init__.py 0 → 100644 +2 −0 Original line number Diff line number Diff line # flake8: noqa from deepchem.feat.molecule_featurizers.mol_graph_conv_featurizer import MolGraphConvFeaturizer deepchem/feat/molecule_featurizers/mol_graph_conv_featurizer.py 0 → 100644 +205 −0 Original line number Diff line number Diff line from typing import List, Optional, Sequence, Tuple, Union import numpy as np from deepchem.utils.typing import RDKitAtom, RDKitBond, RDKitMol from deepchem.utils.graph_conv_utils import get_atom_type_one_hot, get_atomic_number, \ construct_hydrogen_bonding_info, get_atom_hydrogen_bonding_one_hot, \ get_atom_is_in_aromatic_one_hot, get_atom_hybridization_one_hot, \ get_atom_total_num_Hs, get_atom_chirality_one_hot, get_atom_formal_charge, \ get_atom_partial_charge, get_atom_ring_size_one_hot, get_bond_type_one_hot, \ get_bond_is_in_same_ring_one_hot, get_bond_graph_distance_one_hot, \ get_bond_euclidean_distance from deepchem.feat.base_classes import MolecularFeaturizer from deepchem.feat.graph_data import GraphData def constrcut_atom_feature( atom: RDKitAtom, use_mpnn_style: bool, hydrogen_bonding: List[Tuple[int, str]], chiral_center: Optional[List[Tuple[int, str]]] = None, sssr: Optional[Sequence] = None) -> List[Union[int, float]]: """TODO: add docstring""" # common feature atom_type = get_atom_type_one_hot(atom) aromatic = get_atom_is_in_aromatic_one_hot(atom) hybridization = get_atom_hybridization_one_hot(atom) acceptor_donor_one_hot = get_atom_hydrogen_bonding_one_hot( atom, hydrogen_bonding) if use_mpnn_style: # MPNN style atom vecotor atomic_number = get_atomic_number(atom) num_Hs = get_atom_total_num_Hs(atom) return atom_type + atomic_number + acceptor_donor_one_hot + aromatic + \ hybridization + num_Hs # Weave style atom vector if sssr is None or chiral_center is None: raise ValueError("Must set the values to `sssr` and `chiral_center`.") chirality = get_atom_chirality_one_hot(atom, chiral_center) formal_charge = get_atom_formal_charge(atom) partial_charge = get_atom_partial_charge(atom) ring_size = get_atom_ring_size_one_hot(atom, sssr) return atom_type + chirality + formal_charge + partial_charge + \ ring_size + hybridization + acceptor_donor_one_hot + aromatic def construct_bond_feature( bond: RDKitBond, use_mpnn_style: bool, graph_dist_matrix: Optional[np.ndarray] = None, euclidean_dist_matrix: Optional[np.ndarray] = None, ) -> List[Union[int, float]]: """TODO: add docstring""" # common feature bond_type = get_bond_type_one_hot(bond) if use_mpnn_style: # MPNN style bond vecotor if euclidean_dist_matrix is None: raise ValueError("Must set the value to `euclidean_dist_matrix`.") euclidean_distance = get_bond_euclidean_distance(bond, euclidean_dist_matrix) return bond_type + euclidean_distance # Weave style atom vector if graph_dist_matrix is None: raise ValueError("Must set the value to `graph_dist_matrix`.") graph_distance = get_bond_graph_distance_one_hot(bond, graph_dist_matrix) same_ring = get_bond_is_in_same_ring_one_hot(bond) return bond_type + graph_distance + same_ring class MolGraphConvFeaturizer(MolecularFeaturizer): """This class is a featurizer of gerneral graph convolution networks for molecules. The default featurization is based on WeaveNet style edge and node annotation. TODO: add more docstrings. Examples ------- >>> smiles = ["C1CCC1", "C1=CC=CN=C1"] >>> featurizer = MolGraphConvFeaturizer() >>> out = featurizer.featurize(smiles) >>> type(out[0]) <class 'deepchem.feat.graph_data.GraphData'> """ def __init__(self, add_self_loop: bool = False, use_mpnn_style: bool = False): """ Paramters --------- add_self_loop: bool, default False TODO: Docstring use_mpnn_style: bool, default False TODO: Docstring """ self.add_self_loop = add_self_loop self.use_mpnn_style = use_mpnn_style def _featurize(self, mol: RDKitMol) -> GraphData: """Calculate molecule graph features from RDKit mol object. Parametrs --------- mol: rdkit.Chem.rdchem.Mol RDKit mol object. Returns ------- graph: GraphData A molecule graph with some features. """ try: from rdkit import Chem from rdkit.Chem import rdmolops, AllChem except ModuleNotFoundError: raise ValueError("This method requires RDKit to be installed.") # construct atom and bond features hydrogen_bonding = construct_hydrogen_bonding_info(mol) if self.use_mpnn_style: # MPNN style # compute 3D coordinate. Sometimes, this operation raise Error mol_for_coord = AllChem.AddHs(mol) conf_id = AllChem.EmbedMolecule(mol_for_coord) mol_for_coord = AllChem.RemoveHs(mol_for_coord) dist_matrix = rdmolops.Get3DDistanceMatrix(mol_for_coord, confId=conf_id) # construct atom (node) feature atom_features = np.array( [ constrcut_atom_feature(atom, self.use_mpnn_style, hydrogen_bonding) for atom in mol.GetAtoms() ], dtype=np.float, ) # construct edge (bond) information src, dist, bond_features = [], [], [] for bond in mol.GetBonds(): # add edge list considering a directed graph start, end = bond.GetBeginAtomIdx(), bond.GetEndAtomIdx() src += [start, end] dist += [end, start] bond_features += 2 * [ construct_bond_feature( bond, self.use_mpnn_style, euclidean_dist_matrix=dist_matrix) ] if self.add_self_loop: src += [i for i in range(mol.GetNumAtoms())] dist += [i for i in range(mol.GetNumAtoms())] bond_fea_length = len(bond_features[0]) bond_features += 2 * [[0 for _ in range(bond_fea_length)]] return GraphData( node_features=atom_features, edge_index=np.array([src, dist], dtype=np.int), edge_features=np.array(bond_features, dtype=np.float)) # Weave style # compute partial charges AllChem.ComputeGasteigerCharges(mol) dist_matrix = Chem.GetDistanceMatrix(mol) chiral_center = Chem.FindMolChiralCenters(mol) sssr = Chem.GetSymmSSSR(mol) # construct atom (node) feature atom_features = np.array( [ constrcut_atom_feature(atom, self.use_mpnn_style, hydrogen_bonding, chiral_center, sssr) for atom in mol.GetAtoms() ], dtype=np.float, ) # construct edge (bond) information src, dist, bond_features = [], [], [] for bond in mol.GetBonds(): # add edge list considering a directed graph start, end = bond.GetBeginAtomIdx(), bond.GetEndAtomIdx() src += [start, end] dist += [end, start] bond_features += 2 * [ construct_bond_feature( bond, self.use_mpnn_style, graph_dist_matrix=dist_matrix) ] if self.add_self_loop: src += [i for i in range(mol.GetNumAtoms())] dist += [i for i in range(mol.GetNumAtoms())] bond_fea_length = len(bond_features[0]) bond_features += 2 * [[0 for _ in range(bond_fea_length)]] return GraphData( node_features=atom_features, edge_index=np.array([src, dist], dtype=np.int), edge_features=np.array(bond_features, dtype=np.float)) deepchem/feat/tests/test_mol_graph_conv_featurizer.py 0 → 100644 +42 −0 Original line number Diff line number Diff line import unittest from deepchem.feat import MolGraphConvFeaturizer # TODO: Add more test cases class TestMolGraphConvFeaturizer(unittest.TestCase): def test_default_featurizer(self): smiles = ["C1=CC=CN=C1", "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C"] featurizer = MolGraphConvFeaturizer() graph_feat = featurizer.featurize(smiles) assert len(graph_feat) == 2 # assert "C1=CC=CN=C1" assert graph_feat[0].num_nodes == 6 assert graph_feat[0].num_node_features == 25 assert graph_feat[0].num_edges == 12 assert graph_feat[0].num_edge_features == 13 # assert "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C" assert graph_feat[1].num_nodes == 22 assert graph_feat[1].num_node_features == 25 assert graph_feat[1].num_edges == 44 assert graph_feat[1].num_edge_features == 13 def test_mpnn_style_featurizer(self): smiles = ["C1=CC=CN=C1", "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C"] featurizer = MolGraphConvFeaturizer(use_mpnn_style=True) graph_feat = featurizer.featurize(smiles) assert len(graph_feat) == 2 # assert "C1=CC=CN=C1" assert graph_feat[0].num_nodes == 6 assert graph_feat[0].num_node_features == 17 assert graph_feat[0].num_edges == 12 assert graph_feat[0].num_edge_features == 5 # assert "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C" assert graph_feat[1].num_nodes == 22 assert graph_feat[1].num_node_features == 17 assert graph_feat[1].num_edges == 44 assert graph_feat[1].num_edge_features == 5 deepchem/utils/graph_conv_utils.py 0 → 100644 +468 −0 Original line number Diff line number Diff line """ Utilities for constructing node features or bond features. Some functions are based on chainer-chemistry or dgl-lifesci. Repositories: - https://github.com/chainer/chainer-chemistry - https://github.com/awslabs/dgl-lifesci """ import os import logging from typing import List, Union, Sequence, Tuple import numpy as np from deepchem.utils.typing import RDKitAtom, RDKitBond, RDKitMol logger = logging.getLogger(__name__) DEFAULT_ATOM_TYPE_SET = [ "C", "N", "O", "F", "P", "S", "Br", "I", ] DEFAULT_HYBRIDIZATION_SET = ["SP1", "SP2", "SP3"] DEFAULT_RING_SIZE_SET = [3, 4, 5, 6, 7, 8] DEFAULT_BOND_TYPE_SET = ["SINGLE", "DOUBLE", "TRIPLE", "AROMATIC"] DEFAULT_GRAPH_DISTANCE_SET = [1, 2, 3, 4, 5, 6, 7] class _ChemicalFeaturesFactory: """This is a singleton class for RDKit base features.""" _instance = None @classmethod def get_instance(cls): try: from rdkit import RDConfig from rdkit.Chem import ChemicalFeatures except ModuleNotFoundError: raise ValueError("This class requires RDKit to be installed.") if not cls._instance: fdefName = os.path.join(RDConfig.RDDataDir, 'BaseFeatures.fdef') cls._instance = ChemicalFeatures.BuildFeatureFactory(fdefName) return cls._instance def one_hot_encode(val: Union[int, str], allowable_set: Union[List[str], List[int]], include_unknown_set: bool = False) -> List[int]: """One hot encoder for elements of a provided set. Examples -------- >>> one_hot_encode("a", ["a", "b", "c"]) [1, 0, 0] >>> one_hot_encode(2, [0, 1, 2]) [0, 0, 1] >>> one_hot_encode(3, [0, 1, 2]) [0, 0, 0] >>> one_hot_encode(3, [0, 1, 2], True) [0, 0, 0, 1] Parameters ---------- val: int or str The value must be present in `allowable_set`. allowable_set: List[int] or List[str] List of allowable quantities. include_unknown_set: bool, default False If true, the index of all values not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] An one hot vector of val. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. Raises ------ `ValueError` if include_unknown_set is False and `val` is not in `allowable_set`. """ if include_unknown_set is False: if val not in allowable_set: logger.warning("input {0} not in allowable set {1}:".format( val, allowable_set)) if include_unknown_set is False: one_hot_legnth = len(allowable_set) else: one_hot_legnth = len(allowable_set) + 1 one_hot = [0 for _ in range(one_hot_legnth)] try: one_hot[allowable_set.index(val)] = 1 except: if include_unknown_set: # If include_unknown_set is True, set the last index is 1. one_hot[-1] = 1 else: pass return one_hot ################################################################# # atom (node) featurization ################################################################# def get_atom_type_one_hot(atom: RDKitAtom, allowable_set: List[str] = DEFAULT_ATOM_TYPE_SET, include_unknown_set: bool = True) -> List[int]: """Get an one hot feature of an atom type. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object allowable_set: List[str] The atom types to consider. The default set is `["C", "N", "O", "F", "P", "S", "Br", "I"]`. include_unknown_set: bool, default True If true, the index of all atom not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] An one hot vector of atom types. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ return one_hot_encode(atom.GetSymbol(), allowable_set, include_unknown_set) def get_atomic_number(atom: RDKitAtom) -> List[int]: """Get an atomic number of an atom. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of the atomic number. """ return [atom.GetAtomicNum()] def construct_hydrogen_bonding_info(mol: RDKitMol) -> List[Tuple[int, str]]: """Construct hydrogen bonding infos about a molecule. Paramters --------- mol: rdkit.Chem.rdchem.Mol RDKit mol object Returns ------- List[Tuple[int, str]] A list of tuple `(atom_index, hydrogen_bonding_type)`. The `hydrogen_bonding_type` value is "Acceptor" or "Donor". """ factory = _ChemicalFeaturesFactory.get_instance() feats = factory.GetFeaturesForMol(mol) hydrogen_bonding = [] for f in feats: hydrogen_bonding.append((f.GetAtomIds()[0], f.GetFamily())) return hydrogen_bonding def get_atom_hydrogen_bonding_one_hot( atom: RDKitAtom, hydrogen_bonding: List[Tuple[int, str]]) -> List[int]: """Get an one hot feat about whether an atom accepts electrons or donates electrons. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object hydrogen_bonding: List[Tuple[int, str]] The return value of `construct_hydrogen_bonding_info`. The value is a list of tuple `(atom_index, hydrogen_bonding)` like (1, "Acceptor"). Returns ------- List[int] A one hot vector of the ring size type. The first element indicates "Donor", and the second element indicates "Acceptor". """ one_hot = [0, 0] atom_idx = atom.GetIdx for hydrogen_bonding_tuple in hydrogen_bonding: if hydrogen_bonding_tuple[0] == atom_idx: if hydrogen_bonding_tuple[1] == "Donor": one_hot[0] = 1 elif hydrogen_bonding_tuple[1] == "Acceptor": one_hot[1] = 1 return one_hot def get_atom_is_in_aromatic_one_hot(atom: RDKitAtom) -> List[int]: """Get ans one hot feature about whether an atom is in aromatic system or not. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of whether an atom is in aromatic system or not. """ return [int(atom.GetIsAromatic())] def get_atom_hybridization_one_hot( atom: RDKitAtom, allowable_set: List[str] = DEFAULT_HYBRIDIZATION_SET, include_unknown_set: bool = False) -> List[int]: """Get an one hot feature of hybridization type. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object allowable_set: List[str] The hybridization types to consider. The default set is `["SP1", "SP2", "SP3"]` include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] An one hot vector of the hybridization type. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ return one_hot_encode( str(atom.GetHybridization()), allowable_set, include_unknown_set) def get_atom_total_num_Hs(atom: RDKitAtom) -> List[int]: """Get the number of hydrogen which an atom has. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of the number of hydrogen which an atom has. """ return [atom.GetTotalNumHs()] def get_atom_chirality_one_hot( atom: RDKitAtom, chiral_center: List[Tuple[int, str]]) -> List[int]: """Get an one hot feature about an atom chirality type. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object chiral_center: List[Tuple[int, str]] The return value of `Chem.FindMolChiralCenters(mol)`. The value is a list of tuple `(atom_index, chirality)` like (1, 'S'). Returns ------- List[int] A one hot vector of the chirality type. The first element indicates "R", and the second element indicates "S". """ one_hot = [0, 0] atom_idx = atom.GetIdx() for chiral_tuple in chiral_center: if chiral_tuple[0] == atom_idx: if chiral_tuple[1] == "R": one_hot[0] = 1 elif chiral_tuple[1] == "S": one_hot[1] = 1 return one_hot def get_atom_formal_charge(atom: RDKitAtom) -> List[int]: """Get a formal charge of an atom. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of the formal charge. """ return [atom.GetFormalCharge()] def get_atom_partial_charge(atom: RDKitAtom) -> List[float]: """Get a partial charge of an atom. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[float] A vector of the parital charge. Notes ----- Before using this function, you must calculate `GasteigerCharge` like `AllChem.ComputeGasteigerCharges(mol)`. """ gasteiger_charge = atom.GetProp('_GasteigerCharge') if gasteiger_charge in ['-nan', 'nan', '-inf', 'inf']: gasteiger_charge = 0 return [float(gasteiger_charge)] def get_atom_ring_size_one_hot(atom: RDKitAtom, sssr: Sequence, allowable_set: List[int] = DEFAULT_RING_SIZE_SET, include_unknown_set: bool = False) -> List[int]: """Get an one hot feature about the ring size if an atom is in a ring. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object sssr: Sequence The return value of `Chem.GetSymmSSSR(mol)`. The value is a sequence of rings. allowable_set: List[int] The ring size types to consider. The default set is `["SINGLE", "DOUBLE", "TRIPLE", "AROMATIC"]`. include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] A one hot vector of the ring size type. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ one_hot = [0 for _ in range(len(allowable_set))] atom_index = atom.GetIdx() if atom.IsInRing(): for ring in sssr: ring = list(ring) if atom_index in ring: ring_size = len(ring) try: one_hot[DEFAULT_RING_SIZE_SET.index(ring_size)] = 1 except: pass return one_hot ################################################################# # bond (edge) featurization ################################################################# def get_bond_type_one_hot(bond: RDKitBond, allowable_set: List[str] = DEFAULT_BOND_TYPE_SET, include_unknown_set: bool = False) -> List[int]: """Get an one hot feature of bond type. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object allowable_set: List[str] The bond types to consider. The default set is `["SINGLE", "DOUBLE", "TRIPLE", "AROMATIC"]`. include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] A one hot vector of the bond type. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ return one_hot_encode( str(bond.GetBondType()), allowable_set, include_unknown_set) def get_bond_is_in_same_ring_one_hot(bond: RDKitBond) -> List[int]: """Get an one hot feature about whether atoms of a bond is in the same ring or not. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object Returns ------- List[int] A one hot vector of whether a bond is in the same ring or not. """ return [int(bond.IsInRing())] def get_bond_graph_distance_one_hot( bond: RDKitBond, graph_dist_matrix: np.ndarray, allowable_set: List[int] = DEFAULT_GRAPH_DISTANCE_SET, include_unknown_set: bool = True) -> List[int]: """Get an one hot feature of graph distance. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object graph_dist_matrix: np.ndarray The return value of `Chem.GetDistanceMatrix(mol)`. The shape is `(num_atoms, num_atoms)`. allowable_set: List[str] The graph distance types to consider. The default set is `[1, 2, ..., 7]`. include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] A one hot vector of the graph distance. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ graph_dist = graph_dist_matrix[bond.GetBeginAtomIdx(), bond.GetEndAtomIdx()] return one_hot_encode(graph_dist, allowable_set, include_unknown_set) def get_bond_euclidean_distance( bond: RDKitBond, euclidean_dist_matrix: np.ndarray) -> List[float]: """Get an one hot feature of euclidean distance. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object euclidean_dist_matrix: np.ndarray The return value of `Chem.GetDistanceMatrix(mol)`. The shape is `(num_atoms, num_atoms)`. Returns ------- List[float] A vector of the euclidean distance. """ euclidean_dist = euclidean_dist_matrix[bond.GetBeginAtomIdx(), bond.GetEndAtomIdx()] return [euclidean_dist] Loading
deepchem/feat/__init__.py +9 −0 Original line number Diff line number Diff line """ Making it easy to import in classes. """ # flake8: noqa # base classes for featurizers from deepchem.feat.base_classes import Featurizer from deepchem.feat.base_classes import MolecularFeaturizer from deepchem.feat.base_classes import MaterialStructureFeaturizer from deepchem.feat.base_classes import MaterialCompositionFeaturizer from deepchem.feat.base_classes import ComplexFeaturizer from deepchem.feat.base_classes import UserDefinedFeaturizer from deepchem.feat.graph_features import ConvMolFeaturizer from deepchem.feat.graph_features import WeaveFeaturizer from deepchem.feat.fingerprints import CircularFingerprint Loading @@ -22,6 +26,11 @@ from deepchem.feat.atomic_coordinates import AtomicCoordinates from deepchem.feat.atomic_coordinates import NeighborListComplexAtomicCoordinates from deepchem.feat.adjacency_fingerprints import AdjacencyFingerprint from deepchem.feat.smiles_featurizers import SmilesToSeq, SmilesToImage # molecule featurizers from deepchem.feat.molecule_featurizers import MolGraphConvFeaturizer # material featurizers from deepchem.feat.material_featurizers import ElementPropertyFingerprint from deepchem.feat.material_featurizers import SineCoulombMatrix from deepchem.feat.material_featurizers import CGCNNFeaturizer
deepchem/feat/molecule_featurizers/__init__.py 0 → 100644 +2 −0 Original line number Diff line number Diff line # flake8: noqa from deepchem.feat.molecule_featurizers.mol_graph_conv_featurizer import MolGraphConvFeaturizer
deepchem/feat/molecule_featurizers/mol_graph_conv_featurizer.py 0 → 100644 +205 −0 Original line number Diff line number Diff line from typing import List, Optional, Sequence, Tuple, Union import numpy as np from deepchem.utils.typing import RDKitAtom, RDKitBond, RDKitMol from deepchem.utils.graph_conv_utils import get_atom_type_one_hot, get_atomic_number, \ construct_hydrogen_bonding_info, get_atom_hydrogen_bonding_one_hot, \ get_atom_is_in_aromatic_one_hot, get_atom_hybridization_one_hot, \ get_atom_total_num_Hs, get_atom_chirality_one_hot, get_atom_formal_charge, \ get_atom_partial_charge, get_atom_ring_size_one_hot, get_bond_type_one_hot, \ get_bond_is_in_same_ring_one_hot, get_bond_graph_distance_one_hot, \ get_bond_euclidean_distance from deepchem.feat.base_classes import MolecularFeaturizer from deepchem.feat.graph_data import GraphData def constrcut_atom_feature( atom: RDKitAtom, use_mpnn_style: bool, hydrogen_bonding: List[Tuple[int, str]], chiral_center: Optional[List[Tuple[int, str]]] = None, sssr: Optional[Sequence] = None) -> List[Union[int, float]]: """TODO: add docstring""" # common feature atom_type = get_atom_type_one_hot(atom) aromatic = get_atom_is_in_aromatic_one_hot(atom) hybridization = get_atom_hybridization_one_hot(atom) acceptor_donor_one_hot = get_atom_hydrogen_bonding_one_hot( atom, hydrogen_bonding) if use_mpnn_style: # MPNN style atom vecotor atomic_number = get_atomic_number(atom) num_Hs = get_atom_total_num_Hs(atom) return atom_type + atomic_number + acceptor_donor_one_hot + aromatic + \ hybridization + num_Hs # Weave style atom vector if sssr is None or chiral_center is None: raise ValueError("Must set the values to `sssr` and `chiral_center`.") chirality = get_atom_chirality_one_hot(atom, chiral_center) formal_charge = get_atom_formal_charge(atom) partial_charge = get_atom_partial_charge(atom) ring_size = get_atom_ring_size_one_hot(atom, sssr) return atom_type + chirality + formal_charge + partial_charge + \ ring_size + hybridization + acceptor_donor_one_hot + aromatic def construct_bond_feature( bond: RDKitBond, use_mpnn_style: bool, graph_dist_matrix: Optional[np.ndarray] = None, euclidean_dist_matrix: Optional[np.ndarray] = None, ) -> List[Union[int, float]]: """TODO: add docstring""" # common feature bond_type = get_bond_type_one_hot(bond) if use_mpnn_style: # MPNN style bond vecotor if euclidean_dist_matrix is None: raise ValueError("Must set the value to `euclidean_dist_matrix`.") euclidean_distance = get_bond_euclidean_distance(bond, euclidean_dist_matrix) return bond_type + euclidean_distance # Weave style atom vector if graph_dist_matrix is None: raise ValueError("Must set the value to `graph_dist_matrix`.") graph_distance = get_bond_graph_distance_one_hot(bond, graph_dist_matrix) same_ring = get_bond_is_in_same_ring_one_hot(bond) return bond_type + graph_distance + same_ring class MolGraphConvFeaturizer(MolecularFeaturizer): """This class is a featurizer of gerneral graph convolution networks for molecules. The default featurization is based on WeaveNet style edge and node annotation. TODO: add more docstrings. Examples ------- >>> smiles = ["C1CCC1", "C1=CC=CN=C1"] >>> featurizer = MolGraphConvFeaturizer() >>> out = featurizer.featurize(smiles) >>> type(out[0]) <class 'deepchem.feat.graph_data.GraphData'> """ def __init__(self, add_self_loop: bool = False, use_mpnn_style: bool = False): """ Paramters --------- add_self_loop: bool, default False TODO: Docstring use_mpnn_style: bool, default False TODO: Docstring """ self.add_self_loop = add_self_loop self.use_mpnn_style = use_mpnn_style def _featurize(self, mol: RDKitMol) -> GraphData: """Calculate molecule graph features from RDKit mol object. Parametrs --------- mol: rdkit.Chem.rdchem.Mol RDKit mol object. Returns ------- graph: GraphData A molecule graph with some features. """ try: from rdkit import Chem from rdkit.Chem import rdmolops, AllChem except ModuleNotFoundError: raise ValueError("This method requires RDKit to be installed.") # construct atom and bond features hydrogen_bonding = construct_hydrogen_bonding_info(mol) if self.use_mpnn_style: # MPNN style # compute 3D coordinate. Sometimes, this operation raise Error mol_for_coord = AllChem.AddHs(mol) conf_id = AllChem.EmbedMolecule(mol_for_coord) mol_for_coord = AllChem.RemoveHs(mol_for_coord) dist_matrix = rdmolops.Get3DDistanceMatrix(mol_for_coord, confId=conf_id) # construct atom (node) feature atom_features = np.array( [ constrcut_atom_feature(atom, self.use_mpnn_style, hydrogen_bonding) for atom in mol.GetAtoms() ], dtype=np.float, ) # construct edge (bond) information src, dist, bond_features = [], [], [] for bond in mol.GetBonds(): # add edge list considering a directed graph start, end = bond.GetBeginAtomIdx(), bond.GetEndAtomIdx() src += [start, end] dist += [end, start] bond_features += 2 * [ construct_bond_feature( bond, self.use_mpnn_style, euclidean_dist_matrix=dist_matrix) ] if self.add_self_loop: src += [i for i in range(mol.GetNumAtoms())] dist += [i for i in range(mol.GetNumAtoms())] bond_fea_length = len(bond_features[0]) bond_features += 2 * [[0 for _ in range(bond_fea_length)]] return GraphData( node_features=atom_features, edge_index=np.array([src, dist], dtype=np.int), edge_features=np.array(bond_features, dtype=np.float)) # Weave style # compute partial charges AllChem.ComputeGasteigerCharges(mol) dist_matrix = Chem.GetDistanceMatrix(mol) chiral_center = Chem.FindMolChiralCenters(mol) sssr = Chem.GetSymmSSSR(mol) # construct atom (node) feature atom_features = np.array( [ constrcut_atom_feature(atom, self.use_mpnn_style, hydrogen_bonding, chiral_center, sssr) for atom in mol.GetAtoms() ], dtype=np.float, ) # construct edge (bond) information src, dist, bond_features = [], [], [] for bond in mol.GetBonds(): # add edge list considering a directed graph start, end = bond.GetBeginAtomIdx(), bond.GetEndAtomIdx() src += [start, end] dist += [end, start] bond_features += 2 * [ construct_bond_feature( bond, self.use_mpnn_style, graph_dist_matrix=dist_matrix) ] if self.add_self_loop: src += [i for i in range(mol.GetNumAtoms())] dist += [i for i in range(mol.GetNumAtoms())] bond_fea_length = len(bond_features[0]) bond_features += 2 * [[0 for _ in range(bond_fea_length)]] return GraphData( node_features=atom_features, edge_index=np.array([src, dist], dtype=np.int), edge_features=np.array(bond_features, dtype=np.float))
deepchem/feat/tests/test_mol_graph_conv_featurizer.py 0 → 100644 +42 −0 Original line number Diff line number Diff line import unittest from deepchem.feat import MolGraphConvFeaturizer # TODO: Add more test cases class TestMolGraphConvFeaturizer(unittest.TestCase): def test_default_featurizer(self): smiles = ["C1=CC=CN=C1", "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C"] featurizer = MolGraphConvFeaturizer() graph_feat = featurizer.featurize(smiles) assert len(graph_feat) == 2 # assert "C1=CC=CN=C1" assert graph_feat[0].num_nodes == 6 assert graph_feat[0].num_node_features == 25 assert graph_feat[0].num_edges == 12 assert graph_feat[0].num_edge_features == 13 # assert "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C" assert graph_feat[1].num_nodes == 22 assert graph_feat[1].num_node_features == 25 assert graph_feat[1].num_edges == 44 assert graph_feat[1].num_edge_features == 13 def test_mpnn_style_featurizer(self): smiles = ["C1=CC=CN=C1", "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C"] featurizer = MolGraphConvFeaturizer(use_mpnn_style=True) graph_feat = featurizer.featurize(smiles) assert len(graph_feat) == 2 # assert "C1=CC=CN=C1" assert graph_feat[0].num_nodes == 6 assert graph_feat[0].num_node_features == 17 assert graph_feat[0].num_edges == 12 assert graph_feat[0].num_edge_features == 5 # assert "O=C(NCc1cc(OC)c(O)cc1)CCCC/C=C/C(C)C" assert graph_feat[1].num_nodes == 22 assert graph_feat[1].num_node_features == 17 assert graph_feat[1].num_edges == 44 assert graph_feat[1].num_edge_features == 5
deepchem/utils/graph_conv_utils.py 0 → 100644 +468 −0 Original line number Diff line number Diff line """ Utilities for constructing node features or bond features. Some functions are based on chainer-chemistry or dgl-lifesci. Repositories: - https://github.com/chainer/chainer-chemistry - https://github.com/awslabs/dgl-lifesci """ import os import logging from typing import List, Union, Sequence, Tuple import numpy as np from deepchem.utils.typing import RDKitAtom, RDKitBond, RDKitMol logger = logging.getLogger(__name__) DEFAULT_ATOM_TYPE_SET = [ "C", "N", "O", "F", "P", "S", "Br", "I", ] DEFAULT_HYBRIDIZATION_SET = ["SP1", "SP2", "SP3"] DEFAULT_RING_SIZE_SET = [3, 4, 5, 6, 7, 8] DEFAULT_BOND_TYPE_SET = ["SINGLE", "DOUBLE", "TRIPLE", "AROMATIC"] DEFAULT_GRAPH_DISTANCE_SET = [1, 2, 3, 4, 5, 6, 7] class _ChemicalFeaturesFactory: """This is a singleton class for RDKit base features.""" _instance = None @classmethod def get_instance(cls): try: from rdkit import RDConfig from rdkit.Chem import ChemicalFeatures except ModuleNotFoundError: raise ValueError("This class requires RDKit to be installed.") if not cls._instance: fdefName = os.path.join(RDConfig.RDDataDir, 'BaseFeatures.fdef') cls._instance = ChemicalFeatures.BuildFeatureFactory(fdefName) return cls._instance def one_hot_encode(val: Union[int, str], allowable_set: Union[List[str], List[int]], include_unknown_set: bool = False) -> List[int]: """One hot encoder for elements of a provided set. Examples -------- >>> one_hot_encode("a", ["a", "b", "c"]) [1, 0, 0] >>> one_hot_encode(2, [0, 1, 2]) [0, 0, 1] >>> one_hot_encode(3, [0, 1, 2]) [0, 0, 0] >>> one_hot_encode(3, [0, 1, 2], True) [0, 0, 0, 1] Parameters ---------- val: int or str The value must be present in `allowable_set`. allowable_set: List[int] or List[str] List of allowable quantities. include_unknown_set: bool, default False If true, the index of all values not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] An one hot vector of val. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. Raises ------ `ValueError` if include_unknown_set is False and `val` is not in `allowable_set`. """ if include_unknown_set is False: if val not in allowable_set: logger.warning("input {0} not in allowable set {1}:".format( val, allowable_set)) if include_unknown_set is False: one_hot_legnth = len(allowable_set) else: one_hot_legnth = len(allowable_set) + 1 one_hot = [0 for _ in range(one_hot_legnth)] try: one_hot[allowable_set.index(val)] = 1 except: if include_unknown_set: # If include_unknown_set is True, set the last index is 1. one_hot[-1] = 1 else: pass return one_hot ################################################################# # atom (node) featurization ################################################################# def get_atom_type_one_hot(atom: RDKitAtom, allowable_set: List[str] = DEFAULT_ATOM_TYPE_SET, include_unknown_set: bool = True) -> List[int]: """Get an one hot feature of an atom type. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object allowable_set: List[str] The atom types to consider. The default set is `["C", "N", "O", "F", "P", "S", "Br", "I"]`. include_unknown_set: bool, default True If true, the index of all atom not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] An one hot vector of atom types. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ return one_hot_encode(atom.GetSymbol(), allowable_set, include_unknown_set) def get_atomic_number(atom: RDKitAtom) -> List[int]: """Get an atomic number of an atom. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of the atomic number. """ return [atom.GetAtomicNum()] def construct_hydrogen_bonding_info(mol: RDKitMol) -> List[Tuple[int, str]]: """Construct hydrogen bonding infos about a molecule. Paramters --------- mol: rdkit.Chem.rdchem.Mol RDKit mol object Returns ------- List[Tuple[int, str]] A list of tuple `(atom_index, hydrogen_bonding_type)`. The `hydrogen_bonding_type` value is "Acceptor" or "Donor". """ factory = _ChemicalFeaturesFactory.get_instance() feats = factory.GetFeaturesForMol(mol) hydrogen_bonding = [] for f in feats: hydrogen_bonding.append((f.GetAtomIds()[0], f.GetFamily())) return hydrogen_bonding def get_atom_hydrogen_bonding_one_hot( atom: RDKitAtom, hydrogen_bonding: List[Tuple[int, str]]) -> List[int]: """Get an one hot feat about whether an atom accepts electrons or donates electrons. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object hydrogen_bonding: List[Tuple[int, str]] The return value of `construct_hydrogen_bonding_info`. The value is a list of tuple `(atom_index, hydrogen_bonding)` like (1, "Acceptor"). Returns ------- List[int] A one hot vector of the ring size type. The first element indicates "Donor", and the second element indicates "Acceptor". """ one_hot = [0, 0] atom_idx = atom.GetIdx for hydrogen_bonding_tuple in hydrogen_bonding: if hydrogen_bonding_tuple[0] == atom_idx: if hydrogen_bonding_tuple[1] == "Donor": one_hot[0] = 1 elif hydrogen_bonding_tuple[1] == "Acceptor": one_hot[1] = 1 return one_hot def get_atom_is_in_aromatic_one_hot(atom: RDKitAtom) -> List[int]: """Get ans one hot feature about whether an atom is in aromatic system or not. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of whether an atom is in aromatic system or not. """ return [int(atom.GetIsAromatic())] def get_atom_hybridization_one_hot( atom: RDKitAtom, allowable_set: List[str] = DEFAULT_HYBRIDIZATION_SET, include_unknown_set: bool = False) -> List[int]: """Get an one hot feature of hybridization type. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object allowable_set: List[str] The hybridization types to consider. The default set is `["SP1", "SP2", "SP3"]` include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] An one hot vector of the hybridization type. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ return one_hot_encode( str(atom.GetHybridization()), allowable_set, include_unknown_set) def get_atom_total_num_Hs(atom: RDKitAtom) -> List[int]: """Get the number of hydrogen which an atom has. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of the number of hydrogen which an atom has. """ return [atom.GetTotalNumHs()] def get_atom_chirality_one_hot( atom: RDKitAtom, chiral_center: List[Tuple[int, str]]) -> List[int]: """Get an one hot feature about an atom chirality type. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object chiral_center: List[Tuple[int, str]] The return value of `Chem.FindMolChiralCenters(mol)`. The value is a list of tuple `(atom_index, chirality)` like (1, 'S'). Returns ------- List[int] A one hot vector of the chirality type. The first element indicates "R", and the second element indicates "S". """ one_hot = [0, 0] atom_idx = atom.GetIdx() for chiral_tuple in chiral_center: if chiral_tuple[0] == atom_idx: if chiral_tuple[1] == "R": one_hot[0] = 1 elif chiral_tuple[1] == "S": one_hot[1] = 1 return one_hot def get_atom_formal_charge(atom: RDKitAtom) -> List[int]: """Get a formal charge of an atom. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[int] A vector of the formal charge. """ return [atom.GetFormalCharge()] def get_atom_partial_charge(atom: RDKitAtom) -> List[float]: """Get a partial charge of an atom. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object Returns ------- List[float] A vector of the parital charge. Notes ----- Before using this function, you must calculate `GasteigerCharge` like `AllChem.ComputeGasteigerCharges(mol)`. """ gasteiger_charge = atom.GetProp('_GasteigerCharge') if gasteiger_charge in ['-nan', 'nan', '-inf', 'inf']: gasteiger_charge = 0 return [float(gasteiger_charge)] def get_atom_ring_size_one_hot(atom: RDKitAtom, sssr: Sequence, allowable_set: List[int] = DEFAULT_RING_SIZE_SET, include_unknown_set: bool = False) -> List[int]: """Get an one hot feature about the ring size if an atom is in a ring. Paramters --------- atom: rdkit.Chem.rdchem.Atom RDKit atom object sssr: Sequence The return value of `Chem.GetSymmSSSR(mol)`. The value is a sequence of rings. allowable_set: List[int] The ring size types to consider. The default set is `["SINGLE", "DOUBLE", "TRIPLE", "AROMATIC"]`. include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] A one hot vector of the ring size type. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ one_hot = [0 for _ in range(len(allowable_set))] atom_index = atom.GetIdx() if atom.IsInRing(): for ring in sssr: ring = list(ring) if atom_index in ring: ring_size = len(ring) try: one_hot[DEFAULT_RING_SIZE_SET.index(ring_size)] = 1 except: pass return one_hot ################################################################# # bond (edge) featurization ################################################################# def get_bond_type_one_hot(bond: RDKitBond, allowable_set: List[str] = DEFAULT_BOND_TYPE_SET, include_unknown_set: bool = False) -> List[int]: """Get an one hot feature of bond type. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object allowable_set: List[str] The bond types to consider. The default set is `["SINGLE", "DOUBLE", "TRIPLE", "AROMATIC"]`. include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] A one hot vector of the bond type. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ return one_hot_encode( str(bond.GetBondType()), allowable_set, include_unknown_set) def get_bond_is_in_same_ring_one_hot(bond: RDKitBond) -> List[int]: """Get an one hot feature about whether atoms of a bond is in the same ring or not. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object Returns ------- List[int] A one hot vector of whether a bond is in the same ring or not. """ return [int(bond.IsInRing())] def get_bond_graph_distance_one_hot( bond: RDKitBond, graph_dist_matrix: np.ndarray, allowable_set: List[int] = DEFAULT_GRAPH_DISTANCE_SET, include_unknown_set: bool = True) -> List[int]: """Get an one hot feature of graph distance. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object graph_dist_matrix: np.ndarray The return value of `Chem.GetDistanceMatrix(mol)`. The shape is `(num_atoms, num_atoms)`. allowable_set: List[str] The graph distance types to consider. The default set is `[1, 2, ..., 7]`. include_unknown_set: bool, default False If true, the index of all types not in `allowable_set` is `len(allowable_set)`. Returns ------- List[int] A one hot vector of the graph distance. If `include_unknown_set` is False, the length is `len(allowable_set)`. If `include_unknown_set` is True, the length is `len(allowable_set) + 1`. """ graph_dist = graph_dist_matrix[bond.GetBeginAtomIdx(), bond.GetEndAtomIdx()] return one_hot_encode(graph_dist, allowable_set, include_unknown_set) def get_bond_euclidean_distance( bond: RDKitBond, euclidean_dist_matrix: np.ndarray) -> List[float]: """Get an one hot feature of euclidean distance. Paramters --------- bond: rdkit.Chem.rdchem.Bond RDKit bond object euclidean_dist_matrix: np.ndarray The return value of `Chem.GetDistanceMatrix(mol)`. The shape is `(num_atoms, num_atoms)`. Returns ------- List[float] A vector of the euclidean distance. """ euclidean_dist = euclidean_dist_matrix[bond.GetBeginAtomIdx(), bond.GetEndAtomIdx()] return [euclidean_dist]
