{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3D Molecular Descriptors for ML/QSAR\n",
    "\n",
    "This notebook demonstrates calculating 3D molecular descriptors from Auto3D-generated conformers for:\n",
    "\n",
    "- **QSAR modeling** - Quantitative Structure-Activity Relationships\n",
    "- **Machine learning** - Features for property prediction\n",
    "- **Molecular similarity** - 3D pharmacophore fingerprints\n",
    "- **Drug design** - Shape and property optimization\n",
    "\n",
    "## Descriptor Categories\n",
    "\n",
    "| Type | Description | Examples |\n",
    "|------|-------------|----------|\n",
    "| **Constitutional** | Atom/bond counts | MW, nAtoms, nRotBonds |\n",
    "| **Topological** | 2D connectivity | TPSA, Chi indices |\n",
    "| **Geometrical** | 3D shape/size | RadiusOfGyration, Asphericity |\n",
    "| **Electronic** | Charge distribution | Dipole, partial charges |\n",
    "| **Quantum** | From NNP energies | HOMO-LUMO gap, atomization energy |\n",
    "\n",
    "The key advantage of 3D descriptors is capturing **conformationally-dependent properties**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import tempfile\n",
    "from pathlib import Path\n",
    "\n",
    "import numpy as np\n",
    "import Auto3D\n",
    "from Auto3D import Auto3DOptions, main\n",
    "from rdkit import Chem\n",
    "from rdkit.Chem import AllChem, Descriptors, Descriptors3D\n",
    "from rdkit.Chem import rdMolDescriptors, Lipinski\n",
    "from rdkit.Chem import rdFreeSASA\n",
    "from rdkit.ML.Descriptors import MoleculeDescriptors\n",
    "import pandas as pd\n",
    "\n",
    "print(f\"Auto3D version: {Auto3D.__version__}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Generate 3D Structures\n",
    "\n",
    "First, generate optimized 3D structures using Auto3D."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Drug molecules for descriptor calculation\n",
    "drug_dataset = {\n",
    "    \"aspirin\": \"CC(=O)OC1=CC=CC=C1C(=O)O\",\n",
    "    \"ibuprofen\": \"CC(C)CC1=CC=C(C=C1)C(C)C(=O)O\",\n",
    "    \"acetaminophen\": \"CC(=O)NC1=CC=C(C=C1)O\",\n",
    "    \"caffeine\": \"CN1C=NC2=C1C(=O)N(C(=O)N2C)C\",\n",
    "    \"naproxen\": \"COC1=CC2=CC(C(C)C(=O)O)=CC=C2C=C1\",\n",
    "    \"diazepam\": \"CN1C(=O)CN=C(C2=CC=CC=C2)C3=CC=CC(Cl)=C31\",\n",
    "    \"fluoxetine\": \"CNCCC(C1=CC=CC=C1)OC2=CC=C(C=C2)C(F)(F)F\",\n",
    "    \"atorvastatin\": \"CC(C)C1=C(C(=C(N1CCC(CC(CC(=O)O)O)O)C2=CC=C(C=C2)F)C3=CC=CC=C3)C(=O)NC4=CC=CC=C4\",\n",
    "    \"metformin\": \"CN(C)C(=N)NC(=N)N\",\n",
    "    \"omeprazole\": \"COC1=CC2=NC(CS(=O)C3=NC4=CC=CC=C4N3)=NC2=CC1OC\",\n",
    "}\n",
    "\n",
    "with tempfile.NamedTemporaryFile(mode='w', suffix='.smi', delete=False) as f:\n",
    "    for name, smi in drug_dataset.items():\n",
    "        f.write(f\"{smi} {name}\\n\")\n",
    "    input_file = f.name\n",
    "\n",
    "print(f\"Dataset: {len(drug_dataset)} molecules\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Generate 3D structures\n",
    "if __name__ == \"__main__\":\n",
    "    config = Auto3DOptions(\n",
    "        path=input_file,\n",
    "        k=1,                        # Lowest energy conformer\n",
    "        optimizing_engine=\"AIMNET\",\n",
    "        use_gpu=True,\n",
    "    )\n",
    "    \n",
    "    output_sdf = main(config)\n",
    "    print(f\"Output: {output_sdf}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. 2D Descriptors (Constitutional & Topological)\n",
    "\n",
    "These are calculated from the molecular graph (don't require 3D)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_2d_descriptors(mol):\n",
    "    \"\"\"\n",
    "    Calculate key 2D molecular descriptors.\n",
    "    \"\"\"\n",
    "    return {\n",
    "        # Constitutional\n",
    "        \"MW\": Descriptors.MolWt(mol),\n",
    "        \"nHeavyAtoms\": mol.GetNumHeavyAtoms(),\n",
    "        \"nRotBonds\": Lipinski.NumRotatableBonds(mol),\n",
    "        \"nRings\": rdMolDescriptors.CalcNumRings(mol),\n",
    "        \"nAromaticRings\": rdMolDescriptors.CalcNumAromaticRings(mol),\n",
    "        \"nHeteroatoms\": rdMolDescriptors.CalcNumHeteroatoms(mol),\n",
    "        \n",
    "        # Lipinski\n",
    "        \"HBD\": Lipinski.NumHDonors(mol),\n",
    "        \"HBA\": Lipinski.NumHAcceptors(mol),\n",
    "        \"LogP\": Descriptors.MolLogP(mol),\n",
    "        \"TPSA\": Descriptors.TPSA(mol),\n",
    "        \n",
    "        # Electronic\n",
    "        \"FormalCharge\": Chem.GetFormalCharge(mol),\n",
    "        \n",
    "        # Complexity\n",
    "        \"FractionCSP3\": rdMolDescriptors.CalcFractionCSP3(mol),\n",
    "        \"NumAliphaticRings\": rdMolDescriptors.CalcNumAliphaticRings(mol),\n",
    "    }\n",
    "\n",
    "\n",
    "# Calculate for all molecules\n",
    "if 'output_sdf' in dir():\n",
    "    mols = list(Chem.SDMolSupplier(output_sdf, removeHs=False))\n",
    "    \n",
    "    data_2d = []\n",
    "    for mol in mols:\n",
    "        if mol is None:\n",
    "            continue\n",
    "        name = mol.GetProp(\"_Name\")\n",
    "        desc = calculate_2d_descriptors(mol)\n",
    "        desc[\"Name\"] = name\n",
    "        data_2d.append(desc)\n",
    "    \n",
    "    df_2d = pd.DataFrame(data_2d)\n",
    "    df_2d = df_2d.set_index(\"Name\")\n",
    "    \n",
    "    print(\"2D Descriptors:\")\n",
    "    print(df_2d.round(2).to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. 3D Descriptors (Geometrical)\n",
    "\n",
    "These require 3D coordinates and describe molecular shape."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_3d_descriptors(mol):\n",
    "    \"\"\"\n",
    "    Calculate 3D geometrical descriptors.\n",
    "    \n",
    "    These depend on the conformer and capture shape information.\n",
    "    \"\"\"\n",
    "    descriptors = {}\n",
    "    \n",
    "    try:\n",
    "        # Shape descriptors\n",
    "        descriptors[\"RadiusOfGyration\"] = rdMolDescriptors.CalcRadiusOfGyration(mol)\n",
    "        descriptors[\"Asphericity\"] = rdMolDescriptors.CalcAsphericity(mol)\n",
    "        descriptors[\"Eccentricity\"] = rdMolDescriptors.CalcEccentricity(mol)\n",
    "        descriptors[\"InertialShapeFactor\"] = rdMolDescriptors.CalcInertialShapeFactor(mol)\n",
    "        descriptors[\"SpherocityIndex\"] = rdMolDescriptors.CalcSpherocityIndex(mol)\n",
    "        \n",
    "        # Principal moments of inertia\n",
    "        pmi = rdMolDescriptors.CalcPMI1(mol), rdMolDescriptors.CalcPMI2(mol), rdMolDescriptors.CalcPMI3(mol)\n",
    "        descriptors[\"PMI1\"] = pmi[0]\n",
    "        descriptors[\"PMI2\"] = pmi[1]\n",
    "        descriptors[\"PMI3\"] = pmi[2]\n",
    "        \n",
    "        # Normalized PMI ratios (useful for shape classification)\n",
    "        if pmi[2] > 0:\n",
    "            descriptors[\"NPR1\"] = pmi[0] / pmi[2]  # Rod-likeness\n",
    "            descriptors[\"NPR2\"] = pmi[1] / pmi[2]  # Disk-likeness\n",
    "        \n",
    "        # Plane of best fit\n",
    "        descriptors[\"PBF\"] = rdMolDescriptors.CalcPBF(mol)\n",
    "        \n",
    "    except Exception as e:\n",
    "        print(f\"Error calculating 3D descriptors: {e}\")\n",
    "    \n",
    "    return descriptors\n",
    "\n",
    "\n",
    "if 'output_sdf' in dir():\n",
    "    data_3d = []\n",
    "    for mol in mols:\n",
    "        if mol is None:\n",
    "            continue\n",
    "        name = mol.GetProp(\"_Name\")\n",
    "        desc = calculate_3d_descriptors(mol)\n",
    "        desc[\"Name\"] = name\n",
    "        data_3d.append(desc)\n",
    "    \n",
    "    df_3d = pd.DataFrame(data_3d)\n",
    "    df_3d = df_3d.set_index(\"Name\")\n",
    "    \n",
    "    print(\"\\n3D Shape Descriptors:\")\n",
    "    print(df_3d.round(3).to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Surface Area Descriptors\n",
    "\n",
    "Solvent-accessible surface area (SASA) is important for solubility and binding."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_surface_descriptors(mol):\n",
    "    \"\"\"\n",
    "    Calculate surface area descriptors.\n",
    "    \"\"\"\n",
    "    descriptors = {}\n",
    "    \n",
    "    try:\n",
    "        # Need to compute radii first\n",
    "        radii = rdFreeSASA.classifyAtoms(mol)\n",
    "        \n",
    "        # Total SASA\n",
    "        sasa = rdFreeSASA.CalcSASA(mol, radii)\n",
    "        descriptors[\"SASA_total\"] = sasa\n",
    "        \n",
    "        # LabuteASA (approximate)\n",
    "        descriptors[\"LabuteASA\"] = rdMolDescriptors.CalcLabuteASA(mol)\n",
    "        \n",
    "        # PEOE VSA descriptors (partial charge weighted)\n",
    "        peoe_vsa = rdMolDescriptors.PEOE_VSA_(mol)\n",
    "        descriptors[\"PEOE_VSA_sum\"] = sum(peoe_vsa)\n",
    "        \n",
    "        # SMR VSA descriptors (MR weighted)\n",
    "        smr_vsa = rdMolDescriptors.SMR_VSA_(mol)\n",
    "        descriptors[\"SMR_VSA_sum\"] = sum(smr_vsa)\n",
    "        \n",
    "    except Exception as e:\n",
    "        print(f\"Error calculating surface descriptors: {e}\")\n",
    "    \n",
    "    return descriptors\n",
    "\n",
    "\n",
    "if 'output_sdf' in dir():\n",
    "    data_surf = []\n",
    "    for mol in mols:\n",
    "        if mol is None:\n",
    "            continue\n",
    "        name = mol.GetProp(\"_Name\")\n",
    "        desc = calculate_surface_descriptors(mol)\n",
    "        desc[\"Name\"] = name\n",
    "        data_surf.append(desc)\n",
    "    \n",
    "    df_surf = pd.DataFrame(data_surf)\n",
    "    df_surf = df_surf.set_index(\"Name\")\n",
    "    \n",
    "    print(\"\\nSurface Area Descriptors:\")\n",
    "    print(df_surf.round(1).to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Energy-Based Descriptors from Auto3D\n",
    "\n",
    "The neural network potentials provide quantum-mechanical-like properties."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_energy_descriptors(mol):\n",
    "    \"\"\"\n",
    "    Extract energy-based descriptors from Auto3D output.\n",
    "    \"\"\"\n",
    "    descriptors = {}\n",
    "    \n",
    "    # Electronic energy from NNP\n",
    "    if mol.HasProp(\"E_hartree\"):\n",
    "        e_hartree = float(mol.GetProp(\"E_hartree\"))\n",
    "        descriptors[\"E_hartree\"] = e_hartree\n",
    "        \n",
    "        # Atomization energy (E / n_atoms) - rough measure of stability\n",
    "        n_atoms = mol.GetNumAtoms()\n",
    "        descriptors[\"E_per_atom\"] = e_hartree * 627.509 / n_atoms  # kcal/mol per atom\n",
    "    \n",
    "    # Relative energy (if multiple conformers)\n",
    "    if mol.HasProp(\"E_rel\"):\n",
    "        descriptors[\"E_rel\"] = float(mol.GetProp(\"E_rel\"))\n",
    "    \n",
    "    # Thermodynamic properties (if calculated)\n",
    "    if mol.HasProp(\"G_hartree\"):\n",
    "        descriptors[\"G_hartree\"] = float(mol.GetProp(\"G_hartree\"))\n",
    "    if mol.HasProp(\"H_hartree\"):\n",
    "        descriptors[\"H_hartree\"] = float(mol.GetProp(\"H_hartree\"))\n",
    "    \n",
    "    return descriptors\n",
    "\n",
    "\n",
    "if 'output_sdf' in dir():\n",
    "    data_energy = []\n",
    "    for mol in mols:\n",
    "        if mol is None:\n",
    "            continue\n",
    "        name = mol.GetProp(\"_Name\")\n",
    "        desc = calculate_energy_descriptors(mol)\n",
    "        desc[\"Name\"] = name\n",
    "        data_energy.append(desc)\n",
    "    \n",
    "    df_energy = pd.DataFrame(data_energy)\n",
    "    df_energy = df_energy.set_index(\"Name\")\n",
    "    \n",
    "    print(\"\\nEnergy-Based Descriptors:\")\n",
    "    print(df_energy.round(4).to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Distance-Based Descriptors\n",
    "\n",
    "Pharmacophore-relevant distances from 3D structures."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_distance_descriptors(mol):\n",
    "    \"\"\"\n",
    "    Calculate distance-based descriptors.\n",
    "    \"\"\"\n",
    "    conf = mol.GetConformer()\n",
    "    n_atoms = mol.GetNumAtoms()\n",
    "    \n",
    "    # Get all pairwise distances\n",
    "    distances = []\n",
    "    for i in range(n_atoms):\n",
    "        for j in range(i+1, n_atoms):\n",
    "            pos_i = conf.GetAtomPosition(i)\n",
    "            pos_j = conf.GetAtomPosition(j)\n",
    "            dist = pos_i.Distance(pos_j)\n",
    "            distances.append(dist)\n",
    "    \n",
    "    distances = np.array(distances)\n",
    "    \n",
    "    return {\n",
    "        \"MaxDist\": distances.max(),           # Molecular diameter\n",
    "        \"MinDist\": distances.min(),           # Closest atoms\n",
    "        \"MeanDist\": distances.mean(),         # Average distance\n",
    "        \"StdDist\": distances.std(),           # Distance variability\n",
    "    }\n",
    "\n",
    "\n",
    "if 'output_sdf' in dir():\n",
    "    data_dist = []\n",
    "    for mol in mols:\n",
    "        if mol is None:\n",
    "            continue\n",
    "        name = mol.GetProp(\"_Name\")\n",
    "        desc = calculate_distance_descriptors(mol)\n",
    "        desc[\"Name\"] = name\n",
    "        data_dist.append(desc)\n",
    "    \n",
    "    df_dist = pd.DataFrame(data_dist)\n",
    "    df_dist = df_dist.set_index(\"Name\")\n",
    "    \n",
    "    print(\"\\nDistance-Based Descriptors (Å):\")\n",
    "    print(df_dist.round(2).to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7. Combined Descriptor Table\n",
    "\n",
    "Merge all descriptors for ML modeling."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if 'df_2d' in dir():\n",
    "    # Combine all descriptors\n",
    "    df_all = df_2d.copy()\n",
    "    \n",
    "    if 'df_3d' in dir():\n",
    "        df_all = df_all.join(df_3d, how='left')\n",
    "    if 'df_surf' in dir():\n",
    "        df_all = df_all.join(df_surf, how='left')\n",
    "    if 'df_dist' in dir():\n",
    "        df_all = df_all.join(df_dist, how='left')\n",
    "    if 'df_energy' in dir():\n",
    "        df_all = df_all.join(df_energy, how='left')\n",
    "    \n",
    "    print(f\"\\nComplete Descriptor Matrix: {df_all.shape[0]} molecules × {df_all.shape[1]} descriptors\")\n",
    "    print(\"\\nDescriptor columns:\")\n",
    "    for i, col in enumerate(df_all.columns):\n",
    "        print(f\"  {i+1:2d}. {col}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Export to CSV for ML\n",
    "if 'df_all' in dir():\n",
    "    output_csv = \"/tmp/molecular_descriptors.csv\"\n",
    "    df_all.to_csv(output_csv)\n",
    "    print(f\"\\nExported descriptors to: {output_csv}\")\n",
    "    \n",
    "    # Show summary statistics\n",
    "    print(\"\\nDescriptor Statistics:\")\n",
    "    print(df_all.describe().round(2).T[[\"mean\", \"std\", \"min\", \"max\"]].to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 8. Molecular Fingerprints\n",
    "\n",
    "3D pharmacophore fingerprints capture spatial arrangement of features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from rdkit.Chem import AllChem\n",
    "from rdkit.Chem.Pharm2D import Gobbi_Pharm2D, Generate\n",
    "\n",
    "def calculate_fingerprints(mol):\n",
    "    \"\"\"\n",
    "    Calculate various molecular fingerprints.\n",
    "    \"\"\"\n",
    "    fps = {}\n",
    "    \n",
    "    # Morgan (ECFP-like) - 2D\n",
    "    fp_morgan = AllChem.GetMorganFingerprintAsBitVect(mol, radius=2, nBits=1024)\n",
    "    fps[\"Morgan_OnBits\"] = fp_morgan.GetNumOnBits()\n",
    "    \n",
    "    # RDKit fingerprint - 2D\n",
    "    fp_rdkit = Chem.RDKFingerprint(mol)\n",
    "    fps[\"RDKit_OnBits\"] = fp_rdkit.GetNumOnBits()\n",
    "    \n",
    "    # 2D Pharmacophore\n",
    "    try:\n",
    "        factory = Gobbi_Pharm2D.factory\n",
    "        fp_pharm = Generate.Gen2DFingerprint(mol, factory)\n",
    "        fps[\"Pharm2D_OnBits\"] = fp_pharm.GetNumOnBits()\n",
    "    except:\n",
    "        fps[\"Pharm2D_OnBits\"] = np.nan\n",
    "    \n",
    "    return fps\n",
    "\n",
    "\n",
    "if 'output_sdf' in dir():\n",
    "    data_fp = []\n",
    "    for mol in mols:\n",
    "        if mol is None:\n",
    "            continue\n",
    "        name = mol.GetProp(\"_Name\")\n",
    "        fps = calculate_fingerprints(mol)\n",
    "        fps[\"Name\"] = name\n",
    "        data_fp.append(fps)\n",
    "    \n",
    "    df_fp = pd.DataFrame(data_fp).set_index(\"Name\")\n",
    "    \n",
    "    print(\"\\nFingerprint Summary:\")\n",
    "    print(df_fp.to_string())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 9. Descriptor Selection for QSAR\n",
    "\n",
    "Best practices for selecting descriptors for ML models."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def select_descriptors(df, variance_threshold=0.01, correlation_threshold=0.95):\n",
    "    \"\"\"\n",
    "    Select descriptors for ML modeling.\n",
    "    \n",
    "    Removes:\n",
    "    - Constant or near-constant descriptors\n",
    "    - Highly correlated descriptors\n",
    "    \"\"\"\n",
    "    # Remove columns with NaN\n",
    "    df_clean = df.dropna(axis=1)\n",
    "    \n",
    "    # Remove low variance\n",
    "    variances = df_clean.var()\n",
    "    high_var = variances[variances > variance_threshold].index\n",
    "    df_clean = df_clean[high_var]\n",
    "    \n",
    "    # Remove highly correlated (keep first of each pair)\n",
    "    corr_matrix = df_clean.corr().abs()\n",
    "    upper = corr_matrix.where(np.triu(np.ones(corr_matrix.shape), k=1).astype(bool))\n",
    "    to_drop = [col for col in upper.columns if any(upper[col] > correlation_threshold)]\n",
    "    df_selected = df_clean.drop(columns=to_drop)\n",
    "    \n",
    "    return df_selected, {\n",
    "        \"original\": len(df.columns),\n",
    "        \"after_nan\": len(df_clean.columns) + len(to_drop),\n",
    "        \"after_variance\": len(high_var),\n",
    "        \"after_correlation\": len(df_selected.columns),\n",
    "        \"removed_correlated\": to_drop\n",
    "    }\n",
    "\n",
    "\n",
    "if 'df_all' in dir():\n",
    "    df_selected, info = select_descriptors(df_all)\n",
    "    \n",
    "    print(\"\\nDescriptor Selection:\")\n",
    "    print(f\"  Original: {info['original']} descriptors\")\n",
    "    print(f\"  After removing NaN columns: {info['after_nan']}\")\n",
    "    print(f\"  After variance filter: {info['after_variance']}\")\n",
    "    print(f\"  After correlation filter: {info['after_correlation']}\")\n",
    "    print(f\"\\nRemoved (highly correlated): {info['removed_correlated']}\")\n",
    "    print(f\"\\nSelected descriptors: {list(df_selected.columns)}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 10. Best Practices\n",
    "\n",
    "### For QSAR/ML\n",
    "1. **Use 3D descriptors** when conformation matters (binding, permeability)\n",
    "2. **Boltzmann-average** for flexible molecules\n",
    "3. **Standardize** descriptors before modeling (z-score)\n",
    "4. **Remove collinear** descriptors to avoid overfitting\n",
    "\n",
    "### Descriptor Selection\n",
    "- Start with comprehensive set (~200 RDKit descriptors)\n",
    "- Remove constant/near-constant\n",
    "- Remove highly correlated (keep more interpretable)\n",
    "- Use feature importance from random forest for final selection\n",
    "\n",
    "### 3D-Specific Considerations\n",
    "- Shape descriptors (PMI, asphericity) are conformer-dependent\n",
    "- Use lowest-energy conformer or Boltzmann average\n",
    "- Energy-based descriptors from NNPs can capture electronic effects"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Cleanup\n",
    "if 'input_file' in dir() and os.path.exists(input_file):\n",
    "    os.unlink(input_file)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary\n",
    "\n",
    "This tutorial covered:\n",
    "\n",
    "1. **2D descriptors** - MW, LogP, TPSA, counts\n",
    "2. **3D shape descriptors** - RadiusOfGyration, Asphericity, PMI\n",
    "3. **Surface descriptors** - SASA, VSA-based descriptors\n",
    "4. **Energy descriptors** - from Auto3D NNP calculations\n",
    "5. **Distance descriptors** - molecular dimensions\n",
    "6. **Fingerprints** - Morgan, RDKit, pharmacophore\n",
    "7. **Descriptor selection** - removing correlated features\n",
    "\n",
    "Key workflow:\n",
    "1. Generate 3D with Auto3D\n",
    "2. Calculate descriptors with RDKit\n",
    "3. Filter and select features\n",
    "4. Use for QSAR/ML modeling"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.10.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}