Ligands: SMILES from a structure

/sites treats a bound ligand like a motif: pick it, and the page cuts it out with its pocket and draws its orbitals where it sits in the protein. Before that, the page works out which molecule it is. A structure file gives a three-letter code and heavy atoms, usually without hydrogens or bond orders. The other pages of the app start from a SMILES, so the ligand needs one.

#Four routes to the SMILES, each checked against the atoms

How the geometry is read.

• A bond is any pair closer than the sum of the covalent radii (Cordero 2008) plus 0.45 Å.
• An atom needs a double bond when it is flat with three neighbours, has a wide angle with two, or sits in a flat five- or six-membered ring.
• Each such atom gets exactly one double bond, from the neighbour at the shortest double-bond length. This is a matching, the same step that turns an aromatic SMILES into a Kekulé structure (Baber & Hodgkin 1992).
• Textbook tautomer preferences break ties that the lengths cannot: amino over imino on an aromatic ring, amide over imidic acid, lactam over lactim, ketone over enol.
• A porphyrin core gets two C=N and two N–H, or N⁻ on a metal.

How the routes are compared. Two routes give the same molecule when their heavy-atom graphs match and every atom carries the same number of hydrogens once both are neutralised. The comparison also checks each @/@@ in a template against the handedness in 3D. The report shows every route, its SMILES and whether it agrees.

Which route is used. The first that has the same atoms as the file supplies the bond orders and hydrogens, in the order: file dictionary, CCD, name. Where none does, the geometry is used:

• a wrong code: one test file calls p-xylene "TMP", which is thymidine monophosphate;
• a "LIG" from a modelling program;
• a model with atoms missing.

The geometry route was tested on 104 ligand instances against the CCD (78/78 complete instances right) and on 156 held-out instances against the files' own dictionaries (151/153). See the experiment: results and forecast.

#Bonds to the protein

Some ligands are bonded to the protein. Atoms of other residues within bonding distance, or declared by the depositors (LINK, _struct_conn), are kept as links:

• an N-linked sugar: NAG C1 to Asn ND2;
• an acyl-enzyme: a carbonyl C to Ser OG;
• a chromophore thioester: p-coumaric acid C1 to Cys69 SG of photoactive yellow protein;
• haem c: two thioethers, CAB and CAC to Cys SG.

A link counts as a neighbour when the geometry is read. In the SMILES the protein atom is replaced by H. In the orbital fragment that H points along the bond. The dictionary's leaving atom (O1 of the sugar, the OH of the acid) is then missing for a reason, and the report says so.

If a dictionary differs from the atoms only next to a link, the link changed the chemistry. Haem c is not haem b: its two former vinyl carbons are sp³. There the geometry is used.

#Metals

Metals are not part of the SMILES. A metal-bound N or S that would carry an H is drawn deprotonated: the porphyrin N⁻ on iron, a thiolate, a sulfonamide N⁻ on zinc. That is how the orbital layer treats a metal complex, as ligand anion plus metal ion.

With the pocket switched on, a zinc that binds the ligand comes into the fragment with its protein ligands, and the ligand's donor atoms join them as dative bonds. Mg, Ca, Na and K come in as ions. Other metals are listed, not computed.

#The bound fragment

The fragment holds three things:

• the ligand: bond orders and hydrogens from the chosen route, on the deposited coordinates, protonated at the chosen pH with the app's site pKa model;
• with "pocket": the side chains of every residue within 4.5 Å (link-H on Cβ, as the ligand scope of a motif);
• a metal it binds, as above.

Glycines have no side chain to cut out. Their contacts, like every contact, are listed in the pocket table with the polar pairs within 3.5 Å.

The orbital layer is the unchanged /api/model layer. Typical cost: 300–360 atoms and 0.3–0.5 s on the server for a drug-sized ligand with its pocket.

What is drawn first: the quick view.

• Orbitals: the ligand's lone pairs, where it gives or takes hydrogen bonds; with a zinc in the fragment, the donor pairs and the zinc's open site instead.
• π: the HOMO of each of the ligand's π systems, as on /molviewer: Hückel coefficients on the Hartree–Fock radial p functions, amber/teal for the two signs.

For XK2 in 1HVR that is 6 lone pairs and 5 π systems: about 60 000 triangles, drawn in about 0.2 s. "Every orbital of the bound molecule" (184 orbitals, about 630 000 triangles) is one choice away.

The other settings.

• π field: π density, the LUMO (lilac/violet: empty), the interference part or the potential.
• "π of the bound molecule only": switch it off to include the pocket's side chains.

A motif keeps its own settings, so moving between a ligand and a motif never mixes them up.

#The free molecule

The SMILES of the free molecule opens on the other pages, which build their own conformer: /molviewer, /molecule, /pka-density, /catalogue, and /name-to-smiles with the name. It is:

• the dictionary's molecule, with leaving atoms back and without metals (protoporphyrin IX for HEM), when a dictionary matched;
• otherwise the molecule as bound, with each link replaced by H (mesoporphyrin for haem c).

#Finding structures by ligand: GL1X

The panel "Find structures by ligand" searches GL1X's ligand index from your browser:

• search: GET https://gl1x.com/api/bio/ligands?q=… lists components;
• binding sites: GET /api/bio/ligands/<code> lists the binding sites of one;
• open: a binding site loads that MolNFT record here, on that ligand: /sites?gl1x=<contract>/<token>&ligand=HEM:A:154.

"Open on GL1X" goes back to the record.

#Surveys: every ligand at once (17.0)

On a lab server (see the lab server) the page offers two surveys:

• "Survey every ligand in this file", under the ligand list;
• "Survey these N binding sites", under a component's binding sites in the GL1X index.

Each ligand is read as above, without drawing it: which route gave its SMILES and whether the others agree, its links, metals and pocket. With "with π orbitals" its bound fragment's π systems come too, each with the Hückel x of its HOMO and LUMO. A table per component comes first:

• how many instances;
• how many distinct molecules;
• the routes, links and metals;
• the commonest pocket residues.

A row opens its ligand, and the survey downloads as CSV or JSON.

Molecules. Two instances are the same molecule when they have the same heavy-atom graph and the same number of hydrogens, wherever those sit. The 18 haems of 19HC are read in six tautomeric forms: which two pyrrole nitrogens carry the charge is a choice of the matching, not chemistry. They count as one molecule, as in InChI's mobile-H layer.

POST /api/batch takes jobs (each a structure as for /api/site: structure_ref, gl1x, example or the file, optionally with ligands) and the options orbitals, ph, pocket, pocket_radius, classes. It answers with rows, consensus and errors. The same survey runs on a folder of files without the server:

python3 tools/lab/batch_ligands.py /data/structures --ccd components.cif.gz --orbitals --out survey

A local dictionary. --ccd-file components.cif.gz gives the server the PDB's whole Chemical Component Dictionary. The dictionary route then needs neither the browser nor data.rcsb.org. Each entry's CACTVS SMILES is used, or its bond table when it has none.

#The API

POST /api/site takes the same structure fields as for a motif (see Protein sites), plus:

Every answer, list_only included, now lists ligands (key, comp, class, heavy_atoms, name, dictionary, center). With ligand, the answer holds:

• ligand: the report, with routes, chosen, links, metal_contacts, pocket, polar_contacts, notes and warnings;
• site, with ligand_protonation;
• model, in the protein's frame;
• free_molecule: smiles, smiles_bound and pages.

curl -s https://app.lcserve.com/api/site -H 'Content-Type: application/json' \
  -d '{"example":"1hvr","ligand":"XK2:A:263"}' | jq '.ligand.chosen, .free_molecule.pages.molviewer'

Errors (HTTP 400, code): ligand_not_found, ligand_not_molecular (a water or a lone metal ion), site_too_large.

#Limits and what is not done

• Stereo. Tetrahedral stereo is read from 3D; E/Z of double bonds is not written.
• Charges. The SMILES of the free molecule is neutral; charges at a pH come from the pKa model, on the page that opens it.
• Corrins and other unusual macrocycles may come out wrong from the geometry. Use the dictionary route; the report shows when it disagrees.
• Incomplete models. Atoms missing from the model are named, not rebuilt.
• Multi-residue ligands. An oligosaccharide is read residue by residue, each link capped.