Documentation
Theory & Documentation
Admetshiny is an R package that provides an interactive Shiny application and a toolbox of R functions for the management, calculation, filtering, visualization and exploratory analysis of molecular descriptors and ADMET (Absorption, Distribution, Metabolism, Excretion and Toxicity) properties of small molecules.
The package structure is organized into 24 R source files. The NAMESPACE declare 31 exports, the description 10 hard imports and 16 suggest, but we’ve recommended install all the packages with the guide script in the website.
Two working modules
CDK & Webchem: a four step pipeline where: 1. obtain SMILES via PubChem identifiers, manual entry or CSV upload; 2. Calculate CDK descriptors; 3. apply drug-likeness filters; 4. visualize the results.
ADMET Master Manager: a four step wizard: 1. Upload & Preview; 2. Map Columns; 3. Filter; 4. Plots
Both modules feed their final filtered dataset into the shared Report tab, which uses rmarkdown to render a comprehensive analysis document in HTML, PDF or Word format.
Data Standardization Schema
The application uses a single canonical column schema regardless of the original data source.
| Field code | Standard column name | Type |
|---|---|---|
| SMILES | CanonicalSMILES | string |
| Name | Name | string |
| LogP | LogP | numeric |
| WLOGP | WLOGP | numeric |
| TPSA | TPSA | numeric |
| HBD | #H-bond donors | numeric |
| HBA | #H bond acceptors | numeric |
| Rotatable Bonds | #Rotatable bonds | numeric |
| Molar Refractivity | MR | numeric |
| Heavy Atoms | #Heavy atoms | numeric |
| Aromatic Heavy Atoms | #Aromatic heavy atoms | numeric |
| GI Absorption | GI absorption | categorical |
| Gi Absorption_num | GI_absorption_num | numeric |
| BBB Permeant | BBB permeant | categorical |
| BBB Permeant_num | BBB_num | numeric |
| Pgp Substrate | Pgp substrate | categorical |
| Pgp Substrate_num | Pgp_num | numeric |
| LogS | LogS | numeric |
| LogD | LogD | numeric |
The mapADMETColumns algorithm
The function executes seven sequential steps:
Rename columns: For each user column pair in the mapping, the field code is looked up in code to standard. If a target name already exists and points to a different column, the old column is dropped first to prevent duplicated name confusion.
Numeric coercion: The standard numeric columns (
MW, LogP, WLOGP, TPSA, HBD, HBA, RB, MR, heavy atoms, aromatics atoms, LogS, LogD, GI_absorption_num, BBB_num_Pgp_num) are coerced with suppress warnings.ADMET numeric to categorical conversion: If the user supplied a numeric column mapped to GI Absorption and another, the corresponding categorical column is derived using a threshold of 0.5
if ("GI_absorption_num" %in% names(data)) {
v <- data$GI_absorption_num
data[["GI absorption"]] <- ifelse(
is.na(v), NA,
ifelse(v >= 0.5, "High", "Low"))
data- Optional CDK calculation for missing descriptors: When the function
calculate_cdk = TRUE, and a CanonicalSMILES column exists, the function computes the set of missing standard descriptors viasetdiff(), MAPS each to its CDK short code. For more information about to this function, you could read the technical documentation. - WLOGP/LogP precedence logic: WLOGP is the BOILED-Egg’s original calibration axis (Daina and Zoete 2016a). When the user has explicitly mapped a column to WLOGP but has not mapped another column to LogP, the generic LogP column is silently set to the WLOGP values:
## ----- 5. If WLOGP present and LogP not mapped, set LogP = WLOGP -----
## WLOGP is preferred for BOILED-Egg (official calibration), but only
## if the user didn't explicitly map a different column to LogP.
if ("WLOGP" %in% names(data) && !"LogP" %in% names(data)) {
data$LogP <- data$WLOGP
}
History note. An earlier version of this step (Issue 7 in worklog.md) used the simpler condition "WLOGP"
%in% names(data) and silently overwrote a user-mapped LogP column. The current guarded condition
respects an explicit LogP mapping.- Compute violation columns: The standardized data is passed through
computeViolationColumns()to add the five “#Violations” columns used by the filters. - Compute ADMET properties: If any of GI absorption, BBB permeant or Pgp substrate is missing,
computeADMETProperties()is called. To honour user-provided categorical columns, the function temporarily stashes any pre-existing categorical ADMET column runscomputeADMETProperties(), then restores the stashed values.
Physicochemical Descriptors
Molecular Weight (MW)
The molecular weight is the sum of the atomic weights of every atom in the molecule, weighted by the count of each atom type:
\[\mathrm{MW} \;=\; \sum_{i=1}^{N} A_i \times \mathrm{MW}_i,\]
That’s equation define the next functions, where N is the number of distinct atom types Ai is the count of atoms of type and MWi is the monoisotopic or average atomic weight of atom type taken from the IUPAC standard atomic weights (Roskoski 2020).
LogP (Lipophilicity)
The partition coefficient LogP between n-octanol and water is estimated by several independent methods. ADMETShiny exposes the two atomistic fragmental predictors ALogP and WLOGP, together with Consensus LogP that aggregates five distinct methods (Filimonov et al. 2014).
At the case to ALogP , is calculated with Development Kit (CDK), that’s implementation of uses an atomistic fragmental method: every atom is assigned to a typed fragment whose contribution has been regression fitted on a curated training set. The LogP is the sum of fragment contributions:
\[ \mathrm{LogP} \;=\; \sum_{i=1}^{N} f_i \times n_i,\]
Topological Polar Surface Area (TPSA)
The TPSA is computed by Ertl’s fragment based method as the sum of polar atom type surface contributions over 43 predefined atom types:
\[\mathrm{WLOGP} \;=\; \sum_{i=1}^{N} f_i \times n_i,\]
where the count of polar atom type is defined by elements, hybridization, and number of attached hydrogens.
Molar Refractivity (MR)
Molar refractivity is related to the refractive index n, the molecular weight MW, and the density d of the pure compound by the Lorenz-Lorentz equation:
\[\mathrm{MR} \;=\; \frac{n^2 - 1}{n^2 + 2} \times \frac{\mathrm{MW}}{d}\]
Because the refractive index and density are not generally available for arbitrary molecules, ADMETShiny uses the Wildman-Crippen fragmental estimator:
\[\mathrm{MR} \;=\; \sum_{i=1}^{N} r_i \times n_i,\]
Hydrogen Bond Donors (HBD)
The HBD count is the number of N-H and O-H bonds in the molecule, that’s the total number hydrogens bound to a nitrogen or oxygen atom.
\[\mathrm{HBD} \;=\; \sum_{i} (\text{atom}_i \in \{\mathrm{N}, \mathrm{O}\}) \times h_i,\]
Hydrogen Bond Acceptors (HBA)
The HBA count is the number of nitrogen and oxygen atoms in the molecule with Lipinski’s definition.
\[\mathrm{HBA} \;=\; \sum_{i} (\text{atom}_i \in \{\mathrm{N}, \mathrm{O}\}).\]
Rotatable Bonds (RB)
A rotatable bond is any single non-ring non-terminal bond between two non-hydrogen atoms:
\[\mathrm{RB} = \sum_{(i,j)\in E} \mathbb{I}(b_{ij} = 1) \cdot \mathbb{I}(\text{ring}_{ij} = 0) \cdot \mathbb{I}(\deg(i) > 1) \cdot \mathbb{I}(\deg(j) > 1)\]
Heavy Atoms and Aromatic Heavy Atoms
The heavy-atom count is the total number of non-hydrogen atoms:
\[\mathrm{HA} = \sum_{i} \mathbb{I}(\text{atom}_i \neq \mathrm{H})\]
The heavy-atom count is the total number of non-hydrogen atoms. The aromatic heavy-atom count is the subset of heavy atoms that belong to an aromatic ring:
Drug-Likeness Filters
Each drug-likeness filter is formalized as an indicator function combination. For each rule we define a violation indicator. V=0 denotes as compound that satisfies the rule, V=1 a violation.
Lipinski Rule of Five
The Lipinski violation count is the number of threshold superations over the four canonical physicochemical properties:
\[V_{\mathrm{Lipinski}} = \mathbb{I}(\mathrm{MW} > 500) + \mathbb{I}(\mathrm{LogP} > 5) + \mathbb{I}(\mathrm{HBA} > 10) + \mathbb{I}(\mathrm{HBD} > 5)\]
A compound is considered drug-like when Vlipinski = 0. A single violation (V=1) is still considered acceptable, following Lipinski’s original “one violation allowed” convention.
Veber Filter
The Veber rule predicts good oral bioavailability from molecular flexibility and polarity. The original Veber criterion uses only RB > 10 and TPSA > 140 branches; ADMETShiny additionally requires the H-bond pair count to exceed 12 in the TPSA based branch, matching the form used by SwissADME (Daina et al. 2017).
\[V_{\mathrm{Veber}} = \Ind\Bigl(\mathrm{RB} > 10 \lor \bigl(\mathrm{TPSA} > 140 \land (\mathrm{HBA} + \mathrm{HBD}) > 12\bigr)\Bigr)\]
Ghose Filter
The Ghose filter defines a drug-like chemical space by four simultaneous range constraints:
\[\begin{aligned} V_{\mathrm{Ghose}} = {} & \mathbf{1}(\mathrm{MW} < 160 \lor \mathrm{MW} > 480) \\ & + \mathbf{1}(\mathrm{MR} < 40 \lor \mathrm{MR} > 130) \\ & + \mathbf{1}(\mathrm{LogP} < -0.4 \lor \mathrm{LogP} > 5.6) \\ & + \mathbf{1}(\mathrm{HA} < 20 \lor \mathrm{HA} > 70) \end{aligned}\]
Each of the four properties must lie inside its permissible range for a compound to satisfy the rule.
Egan Filter
The Egan filter predicts oral absorption from lipophilicity and polarity alone:
\[V_{\mathrm{Egan}} = \mathbb{I}(\mathrm{TPSA} > 131.6) + \mathbb{I}(\mathrm{LogP} > 5.88)\]
A compound satisfies the Egan rule when both indicators are zero.
Muegge Filter
The Muegge “Lead-like” filter adds a pharmacophore points constraint PP = HBA + HBD.
\[\def\Ind{\operatorname{Ind}} \begin{aligned} V_{\mathrm{Muegge}} = {} & \Ind(\mathrm{MW} < 200 \lor \mathrm{MW} > 600) \\ & + \Ind(\mathrm{LogP} < -2 \lor \mathrm{LogP} > 5) \\ & + \Ind(\mathrm{HBA} > 10) \\ & + \Ind(\mathrm{HBD} > 5) \\ & + \Ind(\mathrm{RB} > 15) \\ & + \Ind(\mathrm{TPSA} > 150) \\ & + \Ind(\mathrm{PP} < 4) \end{aligned}\]
The last term enforces the presence of at least four pharmacophore points H-bond donor or acceptor atoms, excluding molecules that lack the minimal polar functionality required for specific target recognition.
BOILED-Egg Model
The BOILED-Egg model classifies small molecules for two ADME endpoints human intestinal absorption (HIA) and blood-brain barrier (BBB) permeation from physicochemical descriptors only: \[TPSA\] on the x-axis and \[LogP\] on the y-axis.
A molecule whose point \[(x,y) = (TPSA, LogP) \] lies inside the HIA ellipse is predicted to be passively absorbed by the gastrointestinal tract; a molecule inside the BBB ellipse “yolk” is predicted to cross the blood-brain barrier.
Official WLOGP Polygons
Each ellipse is fully defined by scalar parameters: the two focus points \[ F_1 = (x_1,y_1)\] and \[ F_2 = (x_2, y_2) \] and the major axis length. A point \[ P= (x,y)\] lies inside the ellipse if and only if the sum of its distances to the two points does not exceed d:
\[d(P, F_1) + d(P, F_2) \;\leq\; d,\]
where the Euclidean point to focus distances are:
\[\begin{aligned} d(P, F_1) &= \sqrt{(x - x_1)^2 + (y - y_1)^2} \\ d(P, F_2) &= \sqrt{(x - x_2)^2 + (y - y_2)^2} \end{aligned} \tag{1}\]
The constants defining the official WLOGP ellipses for HIA+ and BBB+ are calibrated on the SwissADME training set (Daina and Zoete 2016b; Daina et al. 2017). We had use this when WLOGP are available for ADMETShiny.
Point in Polygon Test approximation
For runtime classification, the analytical ellipse is discrete into a closed polygon of n boundary points. The query point P is then decided by the ray-casting algorithm as a horizontal ray is cast from P in the +x direction axis, and the number of polygon edges crossed by the ray is counted. The parity of the count gives the classification:
\[\mathrm{inside}(P) = \left( \sum_{e \in E_{\mathrm{poly}}} \mathbb{I}(\text{ray crosses } e) \right) \bmod 2,\]
That definition in the practice in ADMETShiny delegate this test to sp::point.in.polygon() when available, with the ray-casting formulation above as a pure R fallback.
LogP Source selection and trained ellipses
Because the official ellipses were calibrated against WLOGP, their application to data produced with a different LogP estimator can produce systematic misclassification. ADMETShiny therefore selects the polygon set as a function of the user-supplied LogP source:
\[\mathrm{polygons} = \begin{cases} \text{WLOGP official} & \text{if } \texttt{logp\_source} = \text{"WLOGP"}, \\ \text{ALogP-trained} & \text{otherwise}. \end{cases} \tag{2}\]
When the input data contain ALogP but no WLOGP, ADMETShiny uses ALogP trained polygons. These polygons are re-derived by ADMETShiny following the Monte-Carlo optimization protocol using a CDK training set of 439 molecules (Daina and Zoete 2016b).
Monte-Carlo optimization
We have considered the Daina & Zoete methodology for this optimization. This is performed in the normalized (TPSA, ALogP) space. Normalization uses the empirical mean and standard deviation of the training set.
\[\begin{aligned} x' &= \frac{x - \mu_x}{\sigma_x} & y' &= \frac{y - \mu_y}{\sigma_y} \end{aligned}\]
Each ellipse is parameterized by a five-dimensional vectors. The optimization maximizes a regularized with Matthews Correlation Coefficient:
\[\mathrm{Score} = \mathrm{MCC} - 0.05 \times S\] The MCC itself is the standard balanced classification metric,
\[\mathrm{MCC} = \frac{\mathrm{TP} \times \mathrm{TN} - \mathrm{FP} \times \mathrm{FN}}{\sqrt{(\mathrm{TP}+\mathrm{FP})(\mathrm{TP}+\mathrm{FN})(\mathrm{TN}+\mathrm{FP})(\mathrm{TN}+\mathrm{FN})}} \tag{3}\]
Validation data and accuracy
The re-trained ALogP ellipses were validated against the CDK training set of 439 molecules. The classification performance is summarized in the next table:
| Model | MCC | Accuracy | TP | FN | FP | TN |
|---|---|---|---|---|---|---|
| HIA | 0.986 | 99.5% | 348 | 2 | 0 | 89 |
| BBB | 0.977 | 98.9% | 230 | 5 | 0 | 204 |
The HIA ellipse correctly classifies 348 of 350 absorbed molecules and all 89 non absorbed molecules. The BBB ellipse correctly classifies 230 of 235 permeant molecules and all 204 non permeant molecules. Both ellipses achieve zero false positives, which is desirable for a screening tool where over prediction of absorption & permeation would generate costly false leads.
P-glycoprotein Substrate Prediction & Random Forest Architecture
P-glycoprotein (P-gp, ABCB1) substrate status is predicted by a Random Forest of K= 100 binary decision trees. Each tree is trained independently on a bootstrap sample of the training set; at each internal node, the split is chosen over a random subset of features. Each leaf node stores the fraction of training set samples in that leaf that are P-gp substrates:
\[P(\mathbf{x}) = \dfrac{\#\{\text{substrates in leaf reached by } \mathbf{x}\}}{\#\{\text{samples in leaf reached by } \mathbf{x}\}} \tag{4}\]
The final binary classification applies a threshold of 0.5 to the averaged probability:
\[\text{P-gp substrate} = \begin{cases} \text{"Yes"} & \text{if } P(\text{substrate} \mid \mathbf{x}) \ge 0.5, \\ \text{"No"} & \text{otherwise}. \end{cases} \tag{5}\]
Training Data, methodology & model performance
The model was trained on 882 experimental P-gp (ABCB1) compounds extracted from Metrabase. The class distribution is 474 substrates (positive class) and 408 non-substrates (negative class). The nine physicochemical descriptors of equation are used as features. Their means and standard deviations in the training set are listed in the next table.
| Descriptor | Mean | SD |
|---|---|---|
| MW | 426.05 | 168.47 |
| LogP | 3.28 | 2.22 |
| TPSA | 89.63 | 59.68 |
| HBD | 1.91 | 1.83 |
| HBA | 5.71 | 3.65 |
| RB | 5.72 | 4.54 |
| Heavy Atoms | 30.28 | 12.06 |
| Aromatic Atoms | 11.17 | 7.52 |
| MR | 115.31 | 43.87 |
The model performance was evaluated by five fold stratified cross validation.
| Metric | Random Forest | SD |
|---|---|---|
| Accuracy | 0.694 | 0.025 |
| MCC | 0.383 | 0.050 |
| Sensitivity | 0.700 | - |
| Specificity | 0.591 | - |
Additional Metrics
Composite Drug-Likeness Score
ADMETShiny combines the available drug-likeness filters into a single composite score. For a compound evaluated against n rules, the score is the percentage of rules with zero violations:
\[\mathrm{Score} = \frac{1}{n_{\mathrm{rules}}} \sum_{i=1}^{n_{\mathrm{rules}}} \mathbf{1}(V_i = 0) \times 100 \tag{6}\]
The score is discreted into three bands:
- Excellent — \(\mathrm{Score} \ge 80\).
- Acceptable — \(60 \le \mathrm{Score} < 80\).
- Poor — \(\mathrm{Score} < 60\).
Column Mapping System (ADMET Master Manager)
The ADMET Master module accepts arbitrary CSV or Excel inputs whose columns names may differ from the canonical ADMETShiny schema. The function mapADMETColumns takes a user supplied mapping vector m that maps each canonical name to the corresponding source column name, and applies a deterministic column rename:
\[\mathrm{data}_{\mathrm{standard}} = \operatorname{rename}(\mathrm{data}_{\mathrm{raw}}, \mathbf{m}) \tag{7}\]
After column mapping, ADMET Master converts numerical endpoint predictions into the categorical labels used by the visualization modules. For the HIA endpoint, the conversion is:
\[\text{GI absorption} = \begin{cases} \text{"High"} & \text{if } p_{\mathrm{HIA}} \ge 0.5, \\ \text{"Low"} & \text{if } p_{\mathrm{HIA}} < 0.5. \end{cases} \tag{8}\]