Commit a97636b7 authored by Bharath Ramsundar's avatar Bharath Ramsundar
Browse files

Features

parent 4407e4ce
Loading
Loading
Loading
Loading
+27 −10
Original line number Diff line number Diff line
@@ -14,7 +14,13 @@ logger = logging.getLogger(__name__)

class AtomicCoordinates(MolecularFeaturizer):
  """
  Nx3 matrix of Cartesian coordinates [Angstrom]
  Nx3 matrix of Cartesian coordinates [Bohr]

  RDKit stores atomic coordinates in Angstrom. Atomic unit of length
  is the bohr (1 bohr = 0.529177 Angstrom). Converting units makes
  gradient calculation consistent with most QM software packages.

  TODO(rbharath): Add option for angstrom computation.
  """
  name = ['atomic_coordinates']

@@ -26,6 +32,10 @@ class AtomicCoordinates(MolecularFeaturizer):
    ----------
    mol : RDKit Mol
          Molecule.

    Returns
    -------
    A Numpy ndarray of shape `(N,3)` in Bohr.
    """
    N = mol.GetNumAtoms()
    coords = np.zeros((N, 3))
@@ -43,7 +53,6 @@ class AtomicCoordinates(MolecularFeaturizer):
      coords[atom, 1] = coords_in_bohr[atom].y
      coords[atom, 2] = coords_in_bohr[atom].z

    coords = [coords]
    return coords


@@ -145,13 +154,23 @@ class NeighborListAtomicCoordinates(MolecularFeaturizer):


class NeighborListComplexAtomicCoordinates(ComplexFeaturizer):
  """
  Adjacency list of neighbors for protein-ligand complexes in 3-space.
  """Featurizes a molecular complex as coordinates and adjacency list of neighbors.
  
  Neighbors dtermined by user-dfined distance cutoff.
  Computes the 3D coordinates of the system and an adjacency list of
  neighbors for protein-ligand complexes in 3-space. Neighbors are
  determined by user-defined distance cutoff.
  """

  def __init__(self, max_num_neighbors=None, neighbor_cutoff=4):
    """Initialize this featurizer.

    Parameters
    ----------
    max_num_neighbors: int, optional
      set the maximum number of neighbors
    neighbor_cutoff: int, optional
      The number of neighbors to store in the neighbor list.
    """
    if neighbor_cutoff <= 0:
      raise ValueError("neighbor_cutoff must be positive value.")
    if max_num_neighbors is not None:
@@ -169,10 +188,8 @@ class NeighborListComplexAtomicCoordinates(ComplexFeaturizer):

    Parameters
    ----------
    mol_pdb_file: Str 
      Filename for ligand pdb file. 
    protein_pdb_file: Str 
      Filename for protein pdb file. 
    molecular_complex: Object
      Some representation of a molecular complex.
    """
    mol_coords, ob_mol = rdkit_util.load_molecule(mol_pdb_file)
    protein_coords, protein_mol = rdkit_util.load_molecule(protein_pdb_file)
+3 −0
Original line number Diff line number Diff line
"""
Feature calculations.
"""
import os
import logging
import types
import numpy as np
import multiprocessing
from typing import Iterable, Union, Dict, Any
import logging
from deepchem.utils import rdkit_util

logger = logging.getLogger(__name__)

+5 −16
Original line number Diff line number Diff line
@@ -29,14 +29,6 @@ def featurize_contacts_ecfp(frag1,
    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
@@ -45,7 +37,6 @@ def featurize_contacts_ecfp(frag1,
    ECFP radius
  """
  if pairwise_distances is None:
    #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))
@@ -56,14 +47,10 @@ def featurize_contacts_ecfp(frag1,
  # 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)
  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 (frag1_ecfp_dict, frag2_ecfp_dict)


@@ -121,9 +108,10 @@ class ContactCircularFingerprint(ComplexFeaturizer):
      fragments = rdkit_util.load_complex(molecular_complex, add_hydrogens=False)

    except MoleculeLoadException:
      logging.warning("This molecule cannot be loaded by Rdkit. Returning None")
      logger.warning("This molecule cannot be loaded by Rdkit. Returning None")
      return None
    pairwise_features = []
    # We compute pairwise contact fingerprints
    for (frag1, frag2) in itertools.combinations(fragments, 2):
      # Get coordinates
      distances = compute_pairwise_distances(frag1[0], frag2[0])
@@ -135,9 +123,10 @@ class ContactCircularFingerprint(ComplexFeaturizer):
            distances,
            cutoff=self.cutoff,
            ecfp_degree=self.radius)]
      pairwise_features.append(vector)
      pairwise_features += vector
    
    return np.array(pairwise_features)
    pairwise_features = np.concatenate(pairwise_features)
    return pairwise_features

class ContactCircularVoxelizer(ComplexFeaturizer):
  """Computes ECFP fingerprints on a voxel grid.
+18 −0
Original line number Diff line number Diff line
# This example featurizes a small subset of the PDBBind v2015 protein-ligand core data with the AtomicConvFeaturizer
import deepchem as dc
import logging
from deepchem.molnet.load_function.pdbbind_datasets import get_pdbbind_molecular_complex_files

complex_files = get_pdbbind_molecular_complex_files(subset="core", version="v2015", interactions="protein-ligand", load_binding_pocket=False)
ligand_files = [complex[1] for complex in complex_files]
core_subset = ligand_files[:2]

featurizer = dc.feat.AtomicCoordinates()
features = featurizer.featurize(core_subset)
print("features.shape")
print(features.shape)
print("features[0].shape")
print(features[0].shape)
print("features[1].shape")
print(features[1].shape)
assert features.shape == (2,)
Loading