{
  "skill": "alterlab-molfeat",
  "evals": [
    {
      "id": "ecfp-qsar-pipeline",
      "prompt": "Featurize my training SMILES with ECFP (radius 3, 2048 bits) and build a scikit-learn pipeline that trains a random forest regressor for a QSAR model. Save the featurizer config so I can redeploy it.",
      "expected_output": "Invokes alterlab-molfeat: wraps FPCalculator('ecfp', radius=3, fpSize=2048) in a MoleculeTransformer(n_jobs=-1), composes a sklearn Pipeline with the transformer and a RandomForestRegressor, fits on SMILES, and persists the config with transformer.to_state_yaml_file for reproducibility.",
      "assertions": [
        { "type": "should_trigger", "value": true },
        { "type": "behavior", "value": "Uses molfeat FPCalculator + MoleculeTransformer in a scikit-learn pipeline and saves the featurizer state YAML." }
      ]
    },
    {
      "id": "compare-representations",
      "prompt": "I want to compare ECFP, MACCS, RDKit 2D descriptors, and ChemBERTa embeddings as input representations for the same property-prediction task. Set up all four featurizers.",
      "expected_output": "Invokes alterlab-molfeat: instantiates FPCalculator('ecfp'), FPCalculator('maccs'), and RDKitDescriptors2D() wrapped in MoleculeTransformer (with a dtype so each returns an array), plus PretrainedHFTransformer(kind='ChemBERTa-77M-MLM', notation='smiles') used directly (not wrapped), and produces feature matrices to benchmark the representations head-to-head.",
      "assertions": [
        { "type": "should_trigger", "value": true },
        { "type": "output_contains", "value": "ChemBERTa" }
      ]
    },
    {
      "id": "chemberta-embeddings",
      "prompt": "Generate ChemBERTa embeddings for my 50k SMILES for transfer learning, and cache them so I don't recompute.",
      "expected_output": "Invokes alterlab-molfeat: uses PretrainedHFTransformer(kind='ChemBERTa-77M-MLM', notation='smiles') (molfeat[transformer]; PretrainedMolTransformer is only the abstract base class) to produce 384-dim embeddings, processes large input in chunks, and caches the resulting array (e.g. pickle) to avoid recomputation on reruns.",
      "assertions": [
        { "type": "should_trigger", "value": true },
        { "type": "behavior", "value": "Uses molfeat's pretrained ChemBERTa transformer for embeddings and caches the result for reuse." }
      ]
    },
    {
      "id": "concat-featurizers-screen",
      "prompt": "For virtual screening I want a combined feature vector that concatenates MACCS and ECFP fingerprints, then trains a classifier to rank a 1M compound library.",
      "expected_output": "Invokes alterlab-molfeat: uses FeatConcat(['maccs', 'ecfp'], params={'ecfp': {'length': 2048}}) (FeatConcat takes fingerprint names / FPVecTransformer objects, not FPCalculator) to build the combined 167+2048 representation, featurizes train and library with it, fits a classifier, and ranks the library by predicted probability.",
      "assertions": [
        { "type": "should_trigger", "value": true },
        { "type": "behavior", "value": "Uses FeatConcat to combine MACCS and ECFP into one feature vector for the screening classifier." }
      ]
    },
    {
      "id": "near-miss-pytdc",
      "prompt": "I don't have a dataset yet. Get me a standardized lipophilicity benchmark with scaffold splits so I have labeled SMILES to model.",
      "expected_output": "Does NOT invoke this skill; defers to alterlab-pytdc. The user needs a curated AI-ready dataset with proper splits (Therapeutics Data Commons), not featurization of molecules they already have; molfeat starts once labeled SMILES exist.",
      "assertions": [
        { "type": "should_not_trigger", "value": true },
        { "type": "output_contains", "value": "alterlab-pytdc" }
      ]
    },
    {
      "id": "near-miss-deepchem",
      "prompt": "Train an end-to-end graph convolutional network for molecular property prediction with built-in data loaders, splitters, and model classes in one framework.",
      "expected_output": "Does NOT invoke this skill; defers to alterlab-deepchem. The user wants a full deep-learning framework with bundled GCN model classes, data loaders and splitters end-to-end, whereas molfeat only produces feature vectors and stops at the featurization boundary.",
      "assertions": [
        { "type": "should_not_trigger", "value": true },
        { "type": "output_contains", "value": "alterlab-deepchem" }
      ]
    }
  ]
}

# Molfeat API Reference

## Core Modules

Molfeat is organized into several key modules that provide different aspects of molecular featurization:

- **`molfeat.store`** - Manages model loading, listing, and registration
- **`molfeat.calc`** - Provides calculators for single-molecule featurization
- **`molfeat.trans`** - Offers scikit-learn compatible transformers for batch processing
- **`molfeat.utils`** - Utility functions for data handling
- **`molfeat.viz`** - Visualization tools for molecular features

---

## molfeat.calc - Calculators

Calculators are callable objects that convert individual molecules into feature vectors. They accept either RDKit `Chem.Mol` objects or SMILES strings as input.

### SerializableCalculator (Base Class)

Base abstract class for all calculators. When subclassing, must implement:
- `__call__()` - Required method for featurization
- `__len__()` - Optional, returns output length
- `columns` - Optional property, returns feature names
- `batch_compute()` - Optional, for efficient batch processing

**State Management Methods:**
- `to_state_json()` - Save calculator state as JSON
- `to_state_yaml()` - Save calculator state as YAML
- `from_state_dict()` - Load calculator from state dictionary
- `to_state_dict()` - Export calculator state as dictionary

### FPCalculator

Computes molecular fingerprints. Supports 15+ fingerprint methods.

**Supported Fingerprint Types:**

**Structural Fingerprints:**
- `ecfp` - Extended-connectivity fingerprints (circular)
- `fcfp` - Functional-class fingerprints
- `rdkit` - RDKit topological fingerprints
- `maccs` - MACCS keys (166-bit structural keys)
- `avalon` - Avalon fingerprints
- `pattern` - Pattern fingerprints
- `layered` - Layered fingerprints

**Atom-based Fingerprints:**
- `atompair` - Atom pair fingerprints
- `atompair-count` - Counted atom pairs
- `topological` - Topological torsion fingerprints
- `topological-count` - Counted topological torsions

**Specialized Fingerprints:**
- `map4` - MinHashed atom-pair fingerprint up to 4 bonds
- `secfp` - SMILES extended connectivity fingerprint
- `erg` - Extended reduced graphs
- `estate` - Electrotopological state indices

**Parameters:**
- `method` (str) - Fingerprint type name
- `radius` (int) - Radius for circular fingerprints (default: 3)
- `fpSize` (int) - Fingerprint size (default: 2048)
- `includeChirality` (bool) - Include chirality information
- `counting` (bool) - Use count vectors instead of binary

**Usage:**
```python
from molfeat.calc import FPCalculator

# Create fingerprint calculator
calc = FPCalculator("ecfp", radius=3, fpSize=2048)

# Compute fingerprint for single molecule
fp = calc("CCO")  # Returns numpy array

# Get fingerprint length
length = len(calc)  # 2048

# Get feature names
names = calc.columns
```

**Common Fingerprint Dimensions:**
- MACCS: 167 dimensions
- ECFP (default): 2048 dimensions
- MAP4 (default): 1024 dimensions

### Descriptor Calculators

**RDKitDescriptors2D**
Computes 2D molecular descriptors using RDKit.

```python
from molfeat.calc import RDKitDescriptors2D

calc = RDKitDescriptors2D()
descriptors = calc("CCO")  # Returns 200+ descriptors
```

**RDKitDescriptors3D**
Computes 3D molecular descriptors (requires conformer generation).

**MordredDescriptors**
Calculates over 1800 molecular descriptors using Mordred.

```python
from molfeat.calc import MordredDescriptors

calc = MordredDescriptors()
descriptors = calc("CCO")
```

### Pharmacophore Calculators

**Pharmacophore2D**
RDKit's 2D pharmacophore fingerprint generation.

**Pharmacophore3D**
Consensus pharmacophore fingerprints from multiple conformers.

**CATS**
Computes Chemically Advanced Template Search (CATS) descriptors - pharmacophore point pair distributions.
(The class is `CATS`; there is no `CATSCalculator`.)

**Parameters** (`CATS(max_dist=None, bins=None, scale="raw", use_3d_distances=False)`):
- `use_3d_distances` - use 3D conformer distances instead of topological (2D) distances
- `max_dist`, `bins` - distance range and bins for the pair distributions
- `scale` - Scaling mode: "raw", "num", or "count"

```python
from molfeat.calc import CATS

calc = CATS(scale="raw")
cats = calc("CCO")  # 189 values with the default 2D settings (molfeat 1.0)
```

### Shape Descriptors

**USRDescriptors**
Ultrafast shape recognition descriptors (multiple variants).

**ElectroShapeDescriptors**
Electrostatic shape descriptors combining shape, chirality, and electrostatics.

### Graph-Based Calculators

**ScaffoldKeyCalculator**
Computes 40+ scaffold-based molecular properties.

**AtomCalculator**
Atom-level featurization for graph neural networks.

**BondCalculator**
Bond-level featurization for graph neural networks.

### Utility Function

**get_calculator()**
Factory function to instantiate calculators by name.

```python
from molfeat.calc import get_calculator

# Instantiate any calculator by name
calc = get_calculator("ecfp", radius=3)
calc = get_calculator("maccs")
calc = get_calculator("desc2D")
```

Raises `ValueError` for unsupported featurizers.

---

## molfeat.trans - Transformers

Transformers wrap calculators into complete featurization pipelines for batch processing.

### MoleculeTransformer

Scikit-learn compatible transformer for batch molecular featurization.

**Key Parameters:**
- `featurizer` - Calculator or featurizer to use
- `n_jobs` (int) - Number of parallel jobs (-1 for all cores)
- `dtype` - Output data type (numpy float32/64, torch tensors)
- `verbose` (bool) - Enable verbose logging
- `ignore_errors` is **not** a constructor argument: pass it to the call —
  `feats, ids = transformer(mols, ignore_errors=True)` returns only the valid rows plus their
  input indices, while `transformer.transform(mols, ignore_errors=True)` keeps `None` placeholders.
  Without `dtype`, calls return a list of per-molecule arrays.

**Essential Methods:**
- `transform(mols)` - Processes batches and returns representations
- `_transform(mol)` - Handles individual molecule featurization
- `__call__(mols)` - Convenience wrapper around transform()
- `preprocess(mol)` - Prepares input molecules (not automatically applied)
- `to_state_yaml_file(path)` - Save transformer configuration
- `from_state_yaml_file(path)` - Load transformer configuration

**Usage:**
```python
from molfeat.calc import FPCalculator
from molfeat.trans import MoleculeTransformer
import datamol as dm

# Load molecules
smiles = dm.data.freesolv().sample(100).smiles.values

# Create transformer
calc = FPCalculator("ecfp")
transformer = MoleculeTransformer(calc, n_jobs=-1)

# Featurize batch
features = transformer(smiles)  # Returns numpy array (100, 2048)

# Save configuration
transformer.to_state_yaml_file("ecfp_config.yml")

# Reload
transformer = MoleculeTransformer.from_state_yaml_file("ecfp_config.yml")
```

**Performance:** Set `n_jobs=-1` for near-linear speedup on multi-core CPUs; the win grows with batch size and per-molecule cost (largest for descriptor sets like Mordred, negligible for tiny inputs where process spawn overhead dominates).

### FeatConcat

Concatenates multiple featurizers into unified representations.

```python
import numpy as np
from molfeat.trans import FeatConcat

# Combine multiple fingerprints. FeatConcat only accepts FPVecTransformer objects or
# fingerprint names (it rejects FPCalculator instances) and is itself callable.
concat = FeatConcat(
    ["maccs", "ecfp"],
    params={"ecfp": {"length": 2048}},   # FPVecTransformer's ecfp default length is 2000
    dtype=np.float32,
)
features = concat(smiles)   # (n_mols, 2215) = 167 + 2048; concat.length == 2215
```

### PretrainedMolTransformer

Subclass of `MoleculeTransformer` for pre-trained deep learning models.

**Unique Features:**
- `_embed()` - Batched inference for neural networks
- `_convert()` - Transforms SMILES/molecules into model-compatible formats (e.g. SELFIES strings for language models)
- Integrated caching system for efficient storage

It is an abstract base class; instantiate a concrete subclass (molfeat 1.x ships
`PretrainedHFTransformer`, `FCDTransformer`, `CheMeleonTransformer`, `MolJEPATransformer`;
the DGL GIN and Graphormer transformers were removed in 1.0).

**Usage:**
```python
import numpy as np
from molfeat.trans.pretrained import PretrainedHFTransformer

# Load pretrained model (Hugging Face backbone; needs molfeat[transformer])
transformer = PretrainedHFTransformer(kind="ChemBERTa-77M-MLM", notation="smiles", dtype=np.float32)

# Generate embeddings -> (n_mols, 384)
embeddings = transformer(smiles)
```

### PrecomputedMolTransformer

Transformer for cached/precomputed features.

---

## molfeat.store - Model Store

Manages featurizer discovery, loading, and registration.

### ModelStore

Central hub for accessing available featurizers.

**Key Methods:**
- `available_models` - Property listing all available featurizers
- `search(name=None, **kwargs)` - Search for specific featurizers
- `load(name, **kwargs)` - Load a featurizer by name
- `register(name, card)` - Register custom featurizer

**Usage:**
```python
from molfeat.store.modelstore import ModelStore

# Initialize store
store = ModelStore()

# List all available models
all_models = store.available_models
print(f"Found {len(all_models)} featurizers")

# Search for specific model
results = store.search(name="ChemBERTa-77M-MLM")
if results:
    model_card = results[0]

    # View usage information
    model_card.usage()

    # Load the model
    transformer = model_card.load()

# Direct loading
transformer = store.load("ChemBERTa-77M-MLM")
```

**ModelCard Attributes:**
- `name` - Model identifier
- `description` - Model description
- `version` - Model version
- `authors` - Model authors
- `tags` - Categorization tags
- `usage()` - Display usage examples
- `load(**kwargs)` - Load the model

---

## Common Patterns

### Error Handling

```python
# Enable error tolerance at call time (not in the constructor)
featurizer = MoleculeTransformer(calc, n_jobs=-1, verbose=True, dtype=np.float32)

# Only valid molecules are returned, with their input positions
features, valid_ids = featurizer(smiles_with_errors, ignore_errors=True)
```

### Data Type Control

```python
# NumPy float32 (default)
features = transformer(smiles, enforce_dtype=True)

# PyTorch tensors
import torch
transformer = MoleculeTransformer(calc, dtype=torch.float32)
features = transformer(smiles)
```

### Persistence and Reproducibility

```python
# Save transformer state
transformer.to_state_yaml_file("config.yml")
transformer.to_state_json_file("config.json")

# Load from saved state
transformer = MoleculeTransformer.from_state_yaml_file("config.yml")
transformer = MoleculeTransformer.from_state_json_file("config.json")
```

### Preprocessing

```python
# Manual preprocessing
mol = transformer.preprocess("CCO")

# Transform with preprocessing
features = transformer.transform(smiles_list)
```

---

## Integration Examples

### Scikit-learn Pipeline

```python
from sklearn.pipeline import Pipeline
from sklearn.ensemble import RandomForestClassifier
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator

# Create pipeline
pipeline = Pipeline([
    ('featurizer', MoleculeTransformer(FPCalculator("ecfp"))),
    ('classifier', RandomForestClassifier())
])

# Fit and predict
pipeline.fit(smiles_train, y_train)
predictions = pipeline.predict(smiles_test)
```

### PyTorch Integration

```python
import torch
from torch.utils.data import Dataset, DataLoader
from molfeat.trans import MoleculeTransformer

class MoleculeDataset(Dataset):
    def __init__(self, smiles, labels, transformer):
        self.smiles = smiles
        self.labels = labels
        self.transformer = transformer

    def __len__(self):
        return len(self.smiles)

    def __getitem__(self, idx):
        features = self.transformer(self.smiles[idx])
        return torch.tensor(features), torch.tensor(self.labels[idx])

# Create dataset and dataloader
transformer = MoleculeTransformer(FPCalculator("ecfp"))
dataset = MoleculeDataset(smiles, labels, transformer)
loader = DataLoader(dataset, batch_size=32)
```

---

## Performance Tips

1. **Parallelization**: Use `n_jobs=-1` to utilize all CPU cores
2. **Batch Processing**: Process multiple molecules at once instead of loops
3. **Caching**: Leverage built-in caching for pretrained models
4. **Data Types**: Use float32 instead of float64 when precision allows
5. **Error Handling**: Call with `ignore_errors=True` (it returns `(features, valid_ids)`) for large datasets with potential invalid molecules

# Available Featurizers in Molfeat

This document provides a comprehensive catalog of all featurizers available in molfeat, organized by category.

> **molfeat 1.0 (2026-09) removed the DGL-based pretrained GNNs (GIN variants, JTVAE) and Graphormer.**
> Their model-store cards still list them, but loading fails in 1.x; they are kept below only as
> legacy entries (usable with `molfeat<1`, Python ≤ 3.10). 1.0 added the foundation-model
> featurizers `CheMeleonTransformer` and `MolJEPATransformer` (see "Foundation Models").

## Transformer-Based Language Models

Pre-trained transformer models for molecular embeddings using SMILES/SELFIES representations.

### RoBERTa-style Models
- **Roberta-Zinc480M-102M** - RoBERTa masked language model trained on ~480M SMILES strings from ZINC database
- **ChemBERTa-77M-MLM** - Masked language model based on RoBERTa trained on 77M PubChem compounds
- **ChemBERTa-77M-MTR** - Multitask regression version trained on PubChem compounds

### GPT-style Autoregressive Models
- **GPT2-Zinc480M-87M** - GPT-2 autoregressive language model trained on ~480M SMILES from ZINC
- **ChemGPT-1.2B** - Large transformer (1.2B parameters) pretrained on PubChem10M
- **ChemGPT-19M** - Medium transformer (19M parameters) pretrained on PubChem10M
- **ChemGPT-4.7M** - Small transformer (4.7M parameters) pretrained on PubChem10M

### Specialized Transformer Models
- **MolT5** - Self-supervised framework for molecule captioning and text-based generation

## Foundation Models (molfeat 1.x)

- **CheMeleonTransformer** — CheMeleon descriptor-foundation-model fingerprints (2,048-d); weights downloaded from the authors' Zenodo record and MD5-verified (`from molfeat.trans.pretrained import CheMeleonTransformer`)
- **MolJEPATransformer** — Mol-JEPA embeddings from the authors' Hugging Face checkpoint (CC BY-NC 4.0, custom model code); requires `trust_remote_code=True` and `accept_noncommercial_license=True`; uses the `transformer` + `pyg` extras

## Graph Neural Networks (GNNs) — legacy, removed in molfeat 1.0

Pre-trained graph neural network models operating on molecular graph structures.

### GIN (Graph Isomorphism Network) Variants
All pre-trained on ChEMBL molecules with different objectives:
- **gin-supervised-masking** - Supervised with node masking objective
- **gin-supervised-infomax** - Supervised with graph-level mutual information maximization
- **gin-supervised-edgepred** - Supervised with edge prediction objective
- **gin-supervised-contextpred** - Supervised with context prediction objective

### Other Graph-Based Models
- **JTVAE_zinc_no_kl** - Junction-tree VAE for molecule generation (trained on ZINC)
- **Graphormer-pcqm4mv2** - Graph transformer pretrained on PCQM4Mv2 quantum chemistry dataset for HOMO-LUMO gap prediction

## Molecular Descriptors

Calculators for physico-chemical properties and molecular characteristics.

### 2D Descriptors
- **desc2D** / **rdkit2D** - 200+ RDKit 2D molecular descriptors including:
  - Molecular weight, logP, TPSA
  - H-bond donors/acceptors
  - Rotatable bonds
  - Ring counts and aromaticity
  - Molecular complexity metrics

### 3D Descriptors
- **desc3D** / **rdkit3D** - RDKit 3D molecular descriptors (requires conformer generation)
  - Inertial moments
  - PMI (Principal Moments of Inertia) ratios
  - Asphericity, eccentricity
  - Radius of gyration

### Comprehensive Descriptor Sets
- **mordred** - Over 1800 molecular descriptors covering:
  - Constitutional descriptors
  - Topological indices
  - Connectivity indices
  - Information content
  - 2D/3D autocorrelations
  - WHIM descriptors
  - GETAWAY descriptors
  - And many more

### Electrotopological Descriptors
- **estate** - Electrotopological state (E-State) indices encoding:
  - Atomic environment information
  - Electronic and topological properties
  - Heteroatom contributions

## Molecular Fingerprints

Binary or count-based fixed-length vectors representing molecular substructures.

### Circular Fingerprints (ECFP-style)
- **ecfp** / **ecfp:2** / **ecfp:4** / **ecfp:6** - Extended-connectivity fingerprints
  - Radius variants (2, 4, 6 correspond to diameter)
  - Default: radius=3, 2048 bits
  - Most popular for similarity searching
- **ecfp-count** - Count version of ECFP (non-binary)
- **fcfp** / **fcfp-count** - Functional-class circular fingerprints
  - Similar to ECFP but uses functional groups
  - Better for pharmacophore-based similarity

### Path-Based Fingerprints
- **rdkit** - RDKit topological fingerprints based on linear paths
- **pattern** - Pattern fingerprints (similar to MACCS but automated)
- **layered** - Layered fingerprints with multiple substructure layers

### Key-Based Fingerprints
- **maccs** - MACCS keys (166-bit structural keys)
  - Fixed set of predefined substructures
  - Good for scaffold hopping
  - Fast computation
- **avalon** - Avalon fingerprints
  - Similar to MACCS but more features
  - Optimized for similarity searching

### Atom-Pair Fingerprints
- **atompair** - Atom pair fingerprints
  - Encodes pairs of atoms and distance between them
  - Good for 3D similarity
- **atompair-count** - Count version of atom pairs

### Topological Torsion Fingerprints
- **topological** - Topological torsion fingerprints
  - Encodes sequences of 4 connected atoms
  - Captures local topology
- **topological-count** - Count version of topological torsions

### MinHashed Fingerprints
- **map4** - MinHashed Atom-Pair fingerprint up to 4 bonds
  - Combines atom-pair and ECFP concepts
  - Default: 2048 dimensions in molfeat 1.0 (`FPCalculator("map4")`)
  - Needs the `map4` package from https://github.com/reymond-group/map4 (not on PyPI)
  - Fast and efficient for large datasets
- **secfp** - SMILES Extended Connectivity Fingerprint
  - Operates directly on SMILES strings
  - Captures both substructure and atom-pair information

### Extended Reduced Graph
- **erg** - Extended Reduced Graph
  - Uses pharmacophoric points instead of atoms
  - Reduces graph complexity while preserving key features

## Pharmacophore Descriptors

Features based on pharmacologically relevant functional groups and their spatial relationships.

### CATS (Chemically Advanced Template Search)
- **cats2D** - 2D CATS descriptors
  - Pharmacophore point pair distributions
  - Distance based on shortest path
  - Calculator class `molfeat.calc.CATS` (189 values with default 2D settings)
- **cats3D** - 3D CATS descriptors
  - Euclidean distance based
  - Requires conformer generation
- **cats2D_pharm** / **cats3D_pharm** - Pharmacophore variants

### Gobbi Pharmacophores
- **pharm2D-gobbi** - 2D pharmacophore fingerprints (`Pharmacophore2D(factory="gobbi")`; there is no `FPCalculator("gobbi2D")`)
  - 8 pharmacophore feature types:
    - Hydrophobic
    - Aromatic
    - H-bond acceptor
    - H-bond donor
    - Positive ionizable
    - Negative ionizable
    - Lumped hydrophobe
  - Good for virtual screening

### Pmapper Pharmacophores
- **pmapper2D** - 2D pharmacophore signatures
- **pmapper3D** - 3D pharmacophore signatures
  - High-dimensional pharmacophore descriptors
  - Useful for QSAR and similarity searching

## Shape Descriptors

Descriptors capturing 3D molecular shape and electrostatic properties.

### USR (Ultrafast Shape Recognition)
- **usr** - Basic USR descriptors
  - 12 dimensions encoding shape distribution
  - Extremely fast computation
- **usrcat** - USR with pharmacophoric constraints
  - 60 dimensions (12 per feature type)
  - Combines shape and pharmacophore information

### Electrostatic Shape
- **electroshape** - ElectroShape descriptors
  - Combines molecular shape, chirality, and electrostatics
  - Useful for protein-ligand docking predictions

## Scaffold-Based Descriptors

Descriptors based on molecular scaffolds and core structures.

### Scaffold Keys
- **scaffoldkeys** - Scaffold key calculator
  - 40+ scaffold-based properties
  - Bioisosteric scaffold representation
  - Captures core structural features

## Graph Featurizers for GNN Input

Atom and bond-level features for constructing graph representations for Graph Neural Networks.

### Atom-Level Features
- **atom-onehot** - One-hot encoded atom features
- **atom-default** - Default atom featurization including:
  - Atomic number
  - Degree, formal charge
  - Hybridization
  - Aromaticity
  - Number of hydrogen atoms

### Bond-Level Features
- **bond-onehot** - One-hot encoded bond features
- **bond-default** - Default bond featurization including:
  - Bond type (single, double, triple, aromatic)
  - Conjugation
  - Ring membership
  - Stereochemistry

## Integrated Pretrained Model Collections

Molfeat integrates models from various sources:

### HuggingFace Models
Access to transformer models through HuggingFace hub:
- ChemBERTa variants
- ChemGPT variants
- MolT5
- Custom uploaded models

### DGL-LifeSci Models (legacy — removed in molfeat 1.0)
Pre-trained GNN models from DGL-Life:
- GIN variants with different pre-training tasks
- AttentiveFP models
- MPNN models

### FCD (Fréchet ChemNet Distance)
- **fcd** - Pre-trained CNN for molecular generation evaluation

### Graphormer Models (legacy — removed in molfeat 1.0)
- Graph transformers from Microsoft Research
- Pre-trained on quantum chemistry datasets

## Usage Notes

### Choosing a Featurizer

**For traditional ML (Random Forest, SVM, etc.):**
- Start with **ecfp** or **maccs** fingerprints
- Try **desc2D** for interpretable models
- Use **FeatConcat** to combine multiple fingerprints

**For deep learning:**
- Use **ChemBERTa** or **ChemGPT** for transformer embeddings
- Use **CheMeleonTransformer** for descriptor-pretrained foundation-model fingerprints
- (GIN / Graphormer embeddings are legacy: only available with `molfeat<1`)

**For similarity searching:**
- **ecfp** - General purpose, most popular
- **maccs** - Fast, good for scaffold hopping
- **map4** - Efficient for large-scale searches
- **usr** / **usrcat** - 3D shape similarity

**For pharmacophore-based approaches:**
- **fcfp** - Functional group based
- **cats2D/3D** - Pharmacophore pair distributions
- **pharm2D-gobbi** - Explicit pharmacophore features

**For interpretability:**
- **desc2D** / **mordred** - Named descriptors
- **maccs** - Interpretable substructure keys
- **scaffoldkeys** - Scaffold-based features

### Model Dependencies

Some featurizers require optional dependencies:

- **Transformers** (ChemBERTa, ChemGPT, MolT5): `uv pip install "molfeat[transformer]"`
- **Mordred**: `uv pip install "molfeat[mordred]"`
- **FCD**: `uv pip install "molfeat[fcd]"`
- **Mol-JEPA**: `uv pip install "molfeat[transformer,pyg]"`
- **MAP4**: install the `map4` package from https://github.com/reymond-group/map4 (no molfeat extra)
- **All dependencies**: `uv pip install "molfeat[all]"`
- DGL (gin-*, jtvae) and Graphormer models: removed in molfeat 1.0 (no `dgl`/`graphormer` extras)

### Accessing All Available Models

```python
from molfeat.store.modelstore import ModelStore

store = ModelStore()
all_models = store.available_models

# Print all available featurizers
for model in all_models:
    print(f"{model.name}: {model.description}")

# Search for specific types
transformers = [m for m in all_models if "transformer" in m.tags]
gnn_models = [m for m in all_models if "gnn" in m.tags]
fingerprints = [m for m in all_models if "fingerprint" in m.tags]
```

## Performance Characteristics

### Computational Speed (relative)
**Fastest:**
- maccs
- ecfp
- rdkit fingerprints
- usr

**Medium:**
- desc2D
- cats2D
- Most fingerprints

**Slower:**
- mordred (1800+ descriptors)
- desc3D (requires conformer generation)
- 3D descriptors in general

**Slowest (first run):**
- Pretrained models (ChemBERTa, ChemGPT, GIN)
- Note: Subsequent runs benefit from caching

### Dimensionality

**Low (< 200 dims):**
- maccs (167)
- usr (12)
- usrcat (60)

**Medium (200-2000 dims):**
- desc2D (223 in molfeat 1.0)
- ecfp (2048 default, configurable)
- map4 (2048 default in molfeat 1.0)

**High (> 2000 dims):**
- mordred (1800+)
- Concatenated fingerprints
- Some transformer embeddings

**Variable:**
- Transformer models (e.g. 384 for ChemBERTa-77M-MLM; larger models 768-1024)
- CheMeleon (2048); Mol-JEPA (512 for the default `cls` output)

# Molfeat Usage Examples

This document provides practical examples for common molfeat use cases.

## Installation

```bash
uv pip install molfeat                 # 1.0.0 (2026-09), Python >= 3.11

# With all optional dependencies
uv pip install "molfeat[all]"

# With specific dependencies
uv pip install "molfeat[transformer]"  # For ChemBERTa, ChemGPT, MolT5
uv pip install "molfeat[mordred]"      # For Mordred descriptors
# molfeat 1.0 removed the DGL GIN and Graphormer models (no dgl/graphormer extras)
```

---

## Quick Start

### Basic Featurization Workflow

```python
import datamol as dm
from molfeat.calc import FPCalculator
from molfeat.trans import MoleculeTransformer

# Load sample data
data = dm.data.freesolv().sample(100).smiles.values

# Single molecule featurization
calc = FPCalculator("ecfp")
features_single = calc(data[0])
print(f"Single molecule features shape: {features_single.shape}")
# Output: (2048,)

# Batch featurization with parallelization
transformer = MoleculeTransformer(calc, n_jobs=-1)
features_batch = transformer(data)
print(f"Batch features shape: {features_batch.shape}")
# Output: (100, 2048)
```

---

## Calculator Examples

### Fingerprint Calculators

```python
from molfeat.calc import FPCalculator

# ECFP (Extended-Connectivity Fingerprints)
ecfp = FPCalculator("ecfp", radius=3, fpSize=2048)
fp = ecfp("CCO")  # Ethanol
print(f"ECFP shape: {fp.shape}")  # (2048,)

# MACCS keys
maccs = FPCalculator("maccs")
fp = maccs("c1ccccc1")  # Benzene
print(f"MACCS shape: {fp.shape}")  # (167,)

# Count-based fingerprints
ecfp_count = FPCalculator("ecfp-count", radius=3)
fp_count = ecfp_count("CC(C)CC(C)C")  # Non-binary counts

# MAP4 fingerprints
map4 = FPCalculator("map4")
fp = map4("CC(=O)Oc1ccccc1C(=O)O")  # Aspirin
```

### Descriptor Calculators

```python
from molfeat.calc import RDKitDescriptors2D, MordredDescriptors

# RDKit 2D descriptors (200+ properties)
desc2d = RDKitDescriptors2D()
descriptors = desc2d("CCO")
print(f"Number of 2D descriptors: {len(descriptors)}")

# Get descriptor names
names = desc2d.columns
print(f"First 5 descriptors: {names[:5]}")

# Mordred descriptors (1800+ properties)
mordred = MordredDescriptors()
descriptors = mordred("c1ccccc1O")  # Phenol
print(f"Mordred descriptors: {len(descriptors)}")
```

### Pharmacophore Calculators

```python
from molfeat.calc import CATS

# 2D CATS descriptors (topological distances)
cats = CATS(scale="raw")
descriptors = cats("CC(C)Cc1ccc(C)cc1C")  # Cymene
print(f"CATS descriptors: {descriptors.shape}")  # (189,) with default settings

# 3D CATS descriptors (Euclidean distances; requires a conformer)
cats3d = CATS(use_3d_distances=True, scale="num")
```

---

## Transformer Examples

### Basic Transformer Usage

```python
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator
import datamol as dm

# Prepare data
smiles_list = [
    "CCO",
    "CC(=O)O",
    "c1ccccc1",
    "CC(C)O",
    "CCCC"
]

# Create transformer
calc = FPCalculator("ecfp")
transformer = MoleculeTransformer(calc, n_jobs=-1)

# Transform molecules
features = transformer(smiles_list)
print(f"Features shape: {features.shape}")  # (5, 2048)
```

### Error Handling

```python
# Handle invalid SMILES gracefully
smiles_with_errors = [
    "CCO",           # Valid
    "invalid",       # Invalid
    "CC(=O)O",       # Valid
    "xyz123",        # Invalid
]

transformer = MoleculeTransformer(
    FPCalculator("ecfp"),
    n_jobs=-1,
    verbose=True,           # Log errors
)

# ignore_errors belongs to the call/transform, not the constructor
features, valid_ids = transformer(smiles_with_errors, ignore_errors=True)
print(valid_ids)  # [0, 2] — only the valid molecules are returned

# transform() keeps placeholders instead
features = transformer.transform(smiles_with_errors, ignore_errors=True)
print(features)  # [array(...), None, array(...), None]
```

### Concatenating Multiple Featurizers

```python
import numpy as np
from molfeat.trans import FeatConcat
from molfeat.trans.fp import FPVecTransformer

# FeatConcat accepts fingerprint names or FPVecTransformer objects (not FPCalculator)
# Combine MACCS (167) + ECFP (2048) = 2215 dimensions
concat = FeatConcat(["maccs", "ecfp"], params={"ecfp": {"length": 2048}}, dtype=np.float32)
features = concat(smiles_list)
print(f"Combined features shape: {features.shape}")  # (n, 2215)

# Triple combination with explicit transformers
triple_concat = FeatConcat([
    FPVecTransformer("maccs"),
    FPVecTransformer("ecfp", length=2048),
    FPVecTransformer("rdkit", length=2048),
], dtype=np.float32)
```

### Saving and Loading Configurations

```python
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator

# Create and save transformer
transformer = MoleculeTransformer(
    FPCalculator("ecfp", radius=3, fpSize=2048),
    n_jobs=-1
)

# Save to YAML
transformer.to_state_yaml_file("my_featurizer.yml")

# Save to JSON
transformer.to_state_json_file("my_featurizer.json")

# Load from saved state
loaded_transformer = MoleculeTransformer.from_state_yaml_file("my_featurizer.yml")

# Use loaded transformer
features = loaded_transformer(smiles_list)
```

---

## Pretrained Model Examples

### Using the ModelStore

```python
from molfeat.store.modelstore import ModelStore

# Initialize model store
store = ModelStore()

# List all available models
print(f"Total available models: {len(store.available_models)}")

# Search for specific models
chemberta_models = store.search(name="ChemBERTa")
for model in chemberta_models:
    print(f"- {model.name}: {model.description}")

# Get model information
model_card = store.search(name="ChemBERTa-77M-MLM")[0]
print(f"Model: {model_card.name}")
print(f"Version: {model_card.version}")
print(f"Authors: {model_card.authors}")

# View usage instructions
model_card.usage()

# Load model directly
transformer = store.load("ChemBERTa-77M-MLM")
```

### ChemBERTa Embeddings

```python
import numpy as np
from molfeat.trans.pretrained import PretrainedHFTransformer

# Load ChemBERTa model (needs molfeat[transformer])
chemberta = PretrainedHFTransformer(kind="ChemBERTa-77M-MLM", notation="smiles", dtype=np.float32)

# Generate embeddings
smiles = ["CCO", "CC(=O)O", "c1ccccc1"]
embeddings = chemberta(smiles)
print(f"ChemBERTa embeddings shape: {embeddings.shape}")
# Output: (3, 384) - ChemBERTa-77M-MLM has hidden size 384

# Use in ML pipeline
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(
    embeddings, labels, test_size=0.2
)

clf = RandomForestClassifier()
clf.fit(X_train, y_train)
predictions = clf.predict(X_test)
```

### ChemGPT Models

```python
# ChemGPT models were trained on SELFIES
# Small model (4.7M parameters)
chemgpt_small = PretrainedHFTransformer(kind="ChemGPT-4.7M", notation="selfies", dtype=np.float32)

# Medium model (19M parameters)
chemgpt_medium = PretrainedHFTransformer(kind="ChemGPT-19M", notation="selfies", dtype=np.float32)

# Large model (1.2B parameters)
chemgpt_large = PretrainedHFTransformer(kind="ChemGPT-1.2B", notation="selfies", dtype=np.float32)

# Generate embeddings
embeddings = chemgpt_small(smiles)
```

### Foundation-Model Embeddings (molfeat 1.x)

```python
from molfeat.trans.pretrained import CheMeleonTransformer, MolJEPATransformer

# CheMeleon: 2,048-d fingerprints; checkpoint downloaded from Zenodo and MD5-checked
chemeleon = CheMeleonTransformer()
embeddings = chemeleon(smiles)          # (3, 2048)

# Mol-JEPA: CC BY-NC 4.0 weights with custom model code — both flags are required
moljepa = MolJEPATransformer(trust_remote_code=True, accept_noncommercial_license=True)
```

The DGL GIN (`gin_supervised_*`) and Graphormer models used in older examples were removed in
molfeat 1.0; reproduce them only in a separate `molfeat<1` environment (Python ≤ 3.10).

---

## Machine Learning Integration

### Scikit-learn Pipeline

```python
from sklearn.pipeline import Pipeline
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import cross_val_score
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator

# Create ML pipeline
pipeline = Pipeline([
    ('featurizer', MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)),
    ('classifier', RandomForestClassifier(n_estimators=100))
])

# Train and evaluate
pipeline.fit(smiles_train, y_train)
predictions = pipeline.predict(smiles_test)

# Cross-validation
scores = cross_val_score(pipeline, smiles_all, y_all, cv=5)
print(f"CV scores: {scores.mean():.3f} (+/- {scores.std():.3f})")
```

### Grid Search for Hyperparameter Tuning

```python
from sklearn.model_selection import GridSearchCV
from sklearn.svm import SVC

# Define pipeline
pipeline = Pipeline([
    ('featurizer', MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)),
    ('classifier', SVC())
])

# Define parameter grid
param_grid = {
    'classifier__C': [0.1, 1, 10],
    'classifier__kernel': ['rbf', 'linear'],
    'classifier__gamma': ['scale', 'auto']
}

# Grid search
grid_search = GridSearchCV(pipeline, param_grid, cv=5, n_jobs=-1)
grid_search.fit(smiles_train, y_train)

print(f"Best parameters: {grid_search.best_params_}")
print(f"Best score: {grid_search.best_score_:.3f}")
```

### Multiple Featurizer Comparison

```python
from sklearn.metrics import roc_auc_score

# Test different featurizers
featurizers = {
    'ECFP': FPCalculator("ecfp"),
    'MACCS': FPCalculator("maccs"),
    'RDKit': FPCalculator("rdkit"),
    'Descriptors': RDKitDescriptors2D(),
}

results = {}
for name, calc in featurizers.items():
    transformer = MoleculeTransformer(calc, n_jobs=-1, dtype=np.float32)
    X_train = transformer(smiles_train)
    X_test = transformer(smiles_test)

    clf = RandomForestClassifier(n_estimators=100)
    clf.fit(X_train, y_train)

    y_pred = clf.predict_proba(X_test)[:, 1]
    auc = roc_auc_score(y_test, y_pred)
    results[name] = auc

    print(f"{name}: AUC = {auc:.3f}")
```

### PyTorch Deep Learning

```python
import torch
import torch.nn as nn
from torch.utils.data import Dataset, DataLoader
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator

# Custom dataset
class MoleculeDataset(Dataset):
    def __init__(self, smiles, labels, transformer):
        self.features = transformer(smiles)
        self.labels = torch.tensor(labels, dtype=torch.float32)

    def __len__(self):
        return len(self.labels)

    def __getitem__(self, idx):
        return (
            torch.tensor(self.features[idx], dtype=torch.float32),
            self.labels[idx]
        )

# Prepare data
transformer = MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)
train_dataset = MoleculeDataset(smiles_train, y_train, transformer)
train_loader = DataLoader(train_dataset, batch_size=32, shuffle=True)

# Simple neural network
class MoleculeClassifier(nn.Module):
    def __init__(self, input_dim):
        super().__init__()
        self.network = nn.Sequential(
            nn.Linear(input_dim, 512),
            nn.ReLU(),
            nn.Dropout(0.3),
            nn.Linear(512, 256),
            nn.ReLU(),
            nn.Dropout(0.3),
            nn.Linear(256, 1),
            nn.Sigmoid()
        )

    def forward(self, x):
        return self.network(x)

# Train model
model = MoleculeClassifier(input_dim=2048)
optimizer = torch.optim.Adam(model.parameters(), lr=0.001)
criterion = nn.BCELoss()

for epoch in range(10):
    for batch_features, batch_labels in train_loader:
        optimizer.zero_grad()
        outputs = model(batch_features).squeeze()
        loss = criterion(outputs, batch_labels)
        loss.backward()
        optimizer.step()
```

---

## Advanced Usage Patterns

### Custom Preprocessing

```python
from molfeat.trans import MoleculeTransformer
import datamol as dm

class CustomTransformer(MoleculeTransformer):
    def preprocess(self, mol):
        """Custom preprocessing: standardize molecule"""
        if isinstance(mol, str):
            mol = dm.to_mol(mol)

        # Standardize
        mol = dm.standardize_mol(mol)

        # Remove salts / solvents (datamol has no remove_salts())
        mol = dm.remove_salts_solvents(mol)

        return mol

# Use custom transformer
transformer = CustomTransformer(FPCalculator("ecfp"), n_jobs=-1)
features = transformer(smiles_list)
```

### Featurization with Conformers

```python
import datamol as dm
from molfeat.calc import RDKitDescriptors3D

# Generate conformers
def prepare_3d_mol(smiles):
    mol = dm.to_mol(smiles)
    # dm.conformers.generate adds hydrogens for embedding and returns a Mol with conformers
    mol = dm.conformers.generate(mol, n_confs=1)
    return mol

# 3D descriptors
calc_3d = RDKitDescriptors3D()

smiles = "CC(C)Cc1ccc(C)cc1C"
mol_3d = prepare_3d_mol(smiles)
descriptors_3d = calc_3d(mol_3d)
```

### Parallel Batch Processing

```python
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator
import time

# Large dataset
smiles_large = load_large_dataset()  # e.g., 100,000 molecules

# Test different parallelization levels
for n_jobs in [1, 2, 4, -1]:
    transformer = MoleculeTransformer(
        FPCalculator("ecfp"),
        n_jobs=n_jobs
    )

    start = time.time()
    features = transformer(smiles_large)
    elapsed = time.time() - start

    print(f"n_jobs={n_jobs}: {elapsed:.2f}s")
```

### Caching for Expensive Operations

```python
from molfeat.trans.pretrained import PretrainedHFTransformer
import numpy as np
import pickle

# Load expensive pretrained model
transformer = PretrainedHFTransformer(kind="ChemBERTa-77M-MLM", notation="smiles", dtype=np.float32)

# Cache embeddings for reuse
cache_file = "embeddings_cache.pkl"

try:
    # Try loading cached embeddings
    with open(cache_file, "rb") as f:
        embeddings = pickle.load(f)
    print("Loaded cached embeddings")
except FileNotFoundError:
    # Compute and cache
    embeddings = transformer(smiles_list)
    with open(cache_file, "wb") as f:
        pickle.dump(embeddings, f)
    print("Computed and cached embeddings")
```

---

## Common Workflows

### Virtual Screening Workflow

```python
from molfeat.calc import FPCalculator
from sklearn.ensemble import RandomForestClassifier
import datamol as dm

# 1. Prepare training data (known actives/inactives)
train_smiles = load_training_data()
train_labels = load_training_labels()  # 1=active, 0=inactive

# 2. Featurize training set
transformer = MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)
X_train = transformer(train_smiles)

# 3. Train classifier
clf = RandomForestClassifier(n_estimators=500, n_jobs=-1)
clf.fit(X_train, train_labels)

# 4. Featurize screening library
screening_smiles = load_screening_library()  # e.g., 1M compounds
X_screen = transformer(screening_smiles)

# 5. Predict and rank
predictions = clf.predict_proba(X_screen)[:, 1]
ranked_indices = predictions.argsort()[::-1]

# 6. Get top hits
top_n = 1000
top_hits = [screening_smiles[i] for i in ranked_indices[:top_n]]
```

### QSAR Model Building

```python
from molfeat.calc import RDKitDescriptors2D
from sklearn.linear_model import Ridge
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import cross_val_score
import numpy as np

# Load QSAR dataset
smiles = load_molecules()
y = load_activity_values()  # e.g., IC50, logP

# Featurize with interpretable descriptors once (deterministic, no leakage)
transformer = MoleculeTransformer(RDKitDescriptors2D(), n_jobs=-1)
X = transformer(smiles)

# Put scaling INSIDE the pipeline so each CV fold fits its own scaler
# (fitting StandardScaler on all X before cross_val_score leaks test-fold stats)
model = Pipeline([
    ("scaler", StandardScaler()),
    ("ridge", Ridge(alpha=1.0)),
])
scores = cross_val_score(model, X, y, cv=5, scoring="r2")
print(f"R² = {scores.mean():.3f} (+/- {scores.std():.3f})")

# Fit final model on all data
model.fit(X, y)

# Interpret feature importance (coefficients live on the ridge step)
feature_names = transformer.featurizer.columns
coef = model.named_steps["ridge"].coef_
top_features_idx = np.abs(coef).argsort()[-10:][::-1]

print("Top 10 important features:")
for idx in top_features_idx:
    print(f"  {feature_names[idx]}: {coef[idx]:.3f}")
```

### Similarity Search

```python
from molfeat.calc import FPCalculator
from sklearn.metrics.pairwise import cosine_similarity
import numpy as np

# Query molecule
query_smiles = "CC(=O)Oc1ccccc1C(=O)O"  # Aspirin

# Database of molecules
database_smiles = load_molecule_database()  # Large collection

# Compute fingerprints
calc = FPCalculator("ecfp")
query_fp = calc(query_smiles).reshape(1, -1)

transformer = MoleculeTransformer(calc, n_jobs=-1)
database_fps = transformer(database_smiles)

# Compute similarity
similarities = cosine_similarity(query_fp, database_fps)[0]

# Find most similar
top_k = 10
top_indices = similarities.argsort()[-top_k:][::-1]

print(f"Top {top_k} similar molecules:")
for i, idx in enumerate(top_indices, 1):
    print(f"{i}. {database_smiles[idx]} (similarity: {similarities[idx]:.3f})")
```

---

## Troubleshooting

### Handling Invalid Molecules

```python
# Use ignore_errors at call time to skip invalid molecules
transformer = MoleculeTransformer(FPCalculator("ecfp"), verbose=True, dtype=np.float32)

# Only valid rows are returned, plus their positions in the input
valid_features, valid_ids = transformer(smiles_list, ignore_errors=True)
valid_smiles = [smiles_list[i] for i in valid_ids]
```

### Memory Management for Large Datasets

```python
# Process in chunks for very large datasets
def featurize_in_chunks(smiles_list, transformer, chunk_size=10000):
    all_features = []

    for i in range(0, len(smiles_list), chunk_size):
        chunk = smiles_list[i:i+chunk_size]
        features = transformer(chunk)
        all_features.append(features)
        print(f"Processed {i+len(chunk)}/{len(smiles_list)}")

    return np.vstack(all_features)

# Use with large dataset
features = featurize_in_chunks(large_smiles_list, transformer)
```

### Reproducibility

```python
import random
import numpy as np
import torch

# Set all random seeds
def set_seed(seed=42):
    random.seed(seed)
    np.random.seed(seed)
    torch.manual_seed(seed)
    torch.cuda.manual_seed_all(seed)

set_seed(42)

# Save exact configuration
transformer.to_state_yaml_file("config.yml")

# Document version
import molfeat
print(f"molfeat version: {molfeat.__version__}")
```

# Molfeat Workflows and Advanced Patterns

Copy-ready recipes for end-to-end molfeat workflows, ModelStore discovery, and advanced usage patterns extracted from the skill body.

## Common Workflows

### Building a QSAR Model

```python
from molfeat.trans import MoleculeTransformer
from molfeat.calc import FPCalculator
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import cross_val_score

# Featurize molecules
transformer = MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)
X = transformer(smiles_train)

# Train model
model = RandomForestRegressor(n_estimators=100)
scores = cross_val_score(model, X, y_train, cv=5)
print(f"R² = {scores.mean():.3f}")

# Save configuration for deployment
transformer.to_state_yaml_file("production_featurizer.yml")
```

### Virtual Screening Pipeline

```python
from sklearn.ensemble import RandomForestClassifier

# Train on known actives/inactives
transformer = MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)
X_train = transformer(train_smiles)
clf = RandomForestClassifier(n_estimators=500)
clf.fit(X_train, train_labels)

# Screen large library
X_screen = transformer(screening_library)  # e.g., 1M compounds
predictions = clf.predict_proba(X_screen)[:, 1]

# Rank and select top hits
top_indices = predictions.argsort()[::-1][:1000]
top_hits = [screening_library[i] for i in top_indices]
```

### Similarity Search

```python
from sklearn.metrics.pairwise import cosine_similarity

# Query molecule
calc = FPCalculator("ecfp")
query_fp = calc(query_smiles).reshape(1, -1)

# Database fingerprints
transformer = MoleculeTransformer(calc, n_jobs=-1)
database_fps = transformer(database_smiles)

# Compute similarity
similarities = cosine_similarity(query_fp, database_fps)[0]
top_similar = similarities.argsort()[-10:][::-1]
```

### Scikit-learn Pipeline Integration

```python
from sklearn.pipeline import Pipeline
from sklearn.ensemble import RandomForestClassifier

# Create end-to-end pipeline
pipeline = Pipeline([
    ('featurizer', MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1)),
    ('classifier', RandomForestClassifier(n_estimators=100))
])

# Train and predict directly on SMILES
pipeline.fit(smiles_train, y_train)
predictions = pipeline.predict(smiles_test)
```

### Comparing Multiple Featurizers

```python
import numpy as np
from molfeat.trans.pretrained import PretrainedHFTransformer

featurizers = {
    'ECFP': MoleculeTransformer(FPCalculator("ecfp"), n_jobs=-1, dtype=np.float32),
    'MACCS': MoleculeTransformer(FPCalculator("maccs"), n_jobs=-1, dtype=np.float32),
    'Descriptors': MoleculeTransformer(RDKitDescriptors2D(), n_jobs=-1, dtype=np.float32),
    # pretrained models are already transformers — don't wrap them in MoleculeTransformer
    'ChemBERTa': PretrainedHFTransformer(kind="ChemBERTa-77M-MLM", notation="smiles", dtype=np.float32),
}

results = {}
for name, transformer in featurizers.items():
    X = transformer(smiles)
    # Evaluate with your ML model
    score = evaluate_model(X, y)
    results[name] = score
```

## Discovering Available Featurizers

Use the ModelStore to explore all available featurizers:

```python
from molfeat.store.modelstore import ModelStore

store = ModelStore()

# List all available models
all_models = store.available_models
print(f"Total featurizers: {len(all_models)}")

# Search for specific models
chemberta_models = store.search(name="ChemBERTa")
for model in chemberta_models:
    print(f"- {model.name}: {model.description}")

# Get usage information
model_card = store.search(name="ChemBERTa-77M-MLM")[0]
model_card.usage()  # Display usage examples

# Load model
transformer = store.load("ChemBERTa-77M-MLM")
```

## Advanced Features

### Custom Preprocessing

```python
class CustomTransformer(MoleculeTransformer):
    def preprocess(self, mol):
        """Custom preprocessing pipeline"""
        if isinstance(mol, str):
            mol = dm.to_mol(mol)
        mol = dm.standardize_mol(mol)
        mol = dm.remove_salts_solvents(mol)   # datamol has no remove_salts()
        return mol

transformer = CustomTransformer(FPCalculator("ecfp"), n_jobs=-1)
```

### Batch Processing Large Datasets

```python
def featurize_in_chunks(smiles_list, transformer, chunk_size=10000):
    """Process large datasets in chunks to manage memory"""
    all_features = []
    for i in range(0, len(smiles_list), chunk_size):
        chunk = smiles_list[i:i+chunk_size]
        features = transformer(chunk)
        all_features.append(features)
    return np.vstack(all_features)
```

### Caching Expensive Embeddings

```python
import pickle

cache_file = "embeddings_cache.pkl"
transformer = PretrainedHFTransformer(kind="ChemBERTa-77M-MLM", notation="smiles", dtype=np.float32)

try:
    with open(cache_file, "rb") as f:
        embeddings = pickle.load(f)
except FileNotFoundError:
    embeddings = transformer(smiles_list)
    with open(cache_file, "wb") as f:
        pickle.dump(embeddings, f)
```