Protein Descriptors for Machine Learning in Drug Discovery

When you’re modeling protein-ligand binding or reactivity, one of the first decisions is what to actually compute as a protein feature. The options range from a single vector derived from the amino acid sequence to a full ensemble of structural and interaction features extracted from molecular dynamics. Picking the wrong level of detail either throws away information the model needs, or buries it in noise the model can’t use.

The descriptors below split into three fundamentally different categories depending on what inputs they require and how they respond to target movement:

• Invariant (sequence-only): Invariant per protein and computed once from sequence inputs.
• State-dependent (structure/MD): Conformational state-dependent; computed per structure or trajectory frame.
• Complex-dependent (interaction): Dependent on the protein-ligand structural interface; calculated from bound complex snapshots.

Features are structured hierarchically: sequence descriptors form the invariant base layer, while structure and interaction descriptors are computed per snapshot and aggregated over the MD ensemble (using mean, standard deviation, percentiles, and state occupancies).

1. Sequence-Based Descriptors

Sequence features are computationally inexpensive and remain invariant across conformations.

• AAindex / Z-scales: Residue-level physicochemical properties, including hydrophobicity, polarity, and volume, represented as compressed scales.
• BioVec / ProtVec: Distributed representation embeddings capturing local sequence context and motifs.
• PSSM: Position-Specific Scoring Matrices defining evolutionary conservation from multiple alignments.
• AA Composition (AAComp), CTD (CTDC/CTDT/CTDD), PAAC: Global composition, property transitions, and sequence-order correlations.
• Deep Protein Language Models (PLMs): Numerical representations from pretrained transformer models (such as ESM, ProtBERT, or UniRep) that capture structural and biophysical information directly from raw sequences.

Pretrained PLM embeddings are utilized for:

• Cross-species and cross-target transfer learning.
• Classification of substrate and inhibitor specificity.
• Modeling targets with sparse experimental datasets.

Mutation indicators are binary vectors encoding specific residue substitutions in active or allosteric sites, serving as another invariant sequence-level feature.

2. Structure-Based Descriptors

Structure-based features characterize target conformations and vary across molecular dynamics trajectories.

• Secondary structure occupancy (such as loop regions surrounding the pocket).
• Backbone dihedral angle ($\phi$/$\psi$) distributions.
• Distances characterizing access channels and gate loop configurations.
• Active site pocket geometry, including volume and bottleneck radius.

3. Complex-Based Descriptors

Complex-based descriptors characterize the protein-ligand interface and require coordinate data of the bound complex.

• Interaction Fingerprints (IFPs): Binary representation of hydrogen bonds, $\pi$-stacking, cation-$\pi$, hydrophobic contacts, and salt bridges (using tools such as PLIP, BINANA, or ODDT).
• PLEC Fingerprints (protein-ligand extended connectivity): Encode local atom environments for protein-ligand contacts within a cutoff. Useful for virtual screening, affinity prediction, and QSAR on docked or MD complexes, but not applicable to sequence-only modeling.

[!NOTE] Sequence-to-ligand generation frameworks (such as Prot2Drug) demonstrate that sequence-only drug-target interaction (DTI) pipelines can generalize to uncharacterized proteins.

Sequence Descriptor Characterization using Propy3

The propy3 package calculates classical sequence descriptors:

• Amino Acid Composition (AAComp) - 20 descriptors: Fraction of each amino acid in the sequence. Global composition, no order information.
• CTD Descriptors (Composition, Transition, Distribution) - 147 descriptors:
  • CTDC (Composition) - 21: Physicochemical class composition including hydrophobicity, volume, and polarity.
  • CTDT (Transition) - 21: Transition frequencies between property classes.
  • CTDD (Distribution) - 105: Positions of class occurrences along the sequence (1%, 25%, 50%, 75%, and 100% position thresholds).
• Pseudo Amino Acid Composition (PAAC) - 30 descriptors: Adds correlation factors on top of composition to capture local sequence-order effects.

Feature Selection for Binding Models

Sequence-based feature relevance varies by modeling objective:

The combination of PAAC, CTDC, and CTDT forms a robust sequence baseline.

Ensemble and MD-Derived Descriptors

Conformational and interaction descriptors must be calculated per trajectory frame and aggregated statistically over the ensemble:

Active-Site Geometry

• Distances from ligand reactive groups to active site centers, such as a metal atom or active site residue.
• Pocket residue planarity and ligand centroid alignment.
• Pocket volume and bottleneck radius (calculated using CAVER).

Access and Dynamics

• Loop opening distances around target pockets.
• Residue root-mean-square fluctuations (RMSF).
• Secondary structure transitions (calculated via DSSP).

Contact Interfaces

• Hydrogen bonds, salt bridges, $\pi$-interactions, and hydrophobic contacts (using PLIP/BINANA/ODDT).
• Frame-level PLEC profiles.
• Residue-ligand contact occupancies and water-mediated network persistence.

Surface and Electrostatics

• Solvent-accessible surface area (SASA) of the pocket and bound ligand.
• Pocket electrostatic potential fields (calculated via APBS).

Dynamic Trajectory Metrics

• Fraction of trajectory frames meeting a defined catalytic cutoff, such as reaction site distance less than 4.5 Angstroms.
• Dwell times in open or closed conformation states.
• Markov State Model (MSM) state probabilities.

Feature Blueprint

1. Data Inputs:
   • Coordinates: Topology file (top.pdb) and trajectory file (traj.xtc).
   • Atom selections: Pocket reference center (such as a catalytic metal or pocket centroid C$\alpha$), ligand atom group, and pocket residues (within 6-8 Å of the pocket).
2. Feature Blocks:
   • 3D Shape: 3D Zernike moments computed on the pocket surface, aggregated over the trajectory (mean $\pm$ standard deviation).
   • Geometry: Distance from catalytic center to reactive atom, ligand centroid to pocket plane, pocket volume and bottleneck radius, and solvent-accessible surface area (SASA), aggregated over frames.
   • Interactions: PLEC profiles and PLIP contact classes, aggregated by state occupancy.
3. Output: A unified feature matrix combining continuous metrics with fingerprint occupancies.
4. Data Preprocessing:
   • Thin trajectories to 10-50 ps intervals.
   • Standardize continuous features ($z$-score scaling).
   • Reduce high-dimensional shape vectors using Principal Component Analysis (PCA).

Code: computing basic geometry features with MDAnalysis

Script to compute reference distances and plane projections:

import MDAnalysis as mda
import numpy as np
import pandas as pd

u = mda.Universe("top.pdb", "traj.xtc")

# Select reference atom in the active site (such as pocket center CA)
ref_atom = u.select_atoms("protein and resid 102 and name CA")[0]
lig = u.select_atoms("resname LIG")
# Select specific reacting ligand atom
react_atom = lig.select_atoms("name C1")[0]

# Select pocket residues defining a local plane
plane_atoms = u.select_atoms("protein and name CA and (resid 90 or resid 95 or resid 100)")

rows = []
for ts in u.trajectory[::10]:  # thin trajectory
    # 1. Reference distance
    d_ref = np.linalg.norm(ref_atom.position - react_atom.position)
    
    # 2. Ligand-to-pocket-plane distance
    X = plane_atoms.positions - plane_atoms.positions.mean(0)
    _, _, vh = np.linalg.svd(X, full_matrices=False)
    n = vh[-1]  # normal vector of the plane
    p0 = plane_atoms.positions.mean(0)
    
    lc = lig.center_of_mass()
    d_plane = abs(np.dot(lc - p0, n))
    
    rows.append({
        "frame": ts.frame, 
        "d_ref": d_ref, 
        "d_plane": d_plane
    })

pd.DataFrame(rows).to_csv("md_geom.csv", index=False)

Shape descriptors: Zernike moments vs ReOmokage feature points

For machine learning workflows, Zernike moments provide fixed-length feature vectors. For structural alignment, ReOmokage-style landmarks preserve absolute spatial relationships.

Takeaways

• Sequence descriptors form the invariant baseline of the feature space.
• Incorporate structure and interaction descriptors only when coordinate data are available, calculating them frame-by-frame.
• Aggregate descriptors statistically over trajectories; single-frame coordinates do not represent dynamic systems.
• Protein-ligand interaction fingerprints (PLEC, IFPs) capture transient contact states that static models omit.