Physicochemical QSAR Analysis of Passive Permeability Across Caco-2 Monolayers

Journal of Pharmaceutical Sciences xxx (2018) 1-9

Contents lists available at ScienceDirect

Journal of Pharmaceutical Sciences journal homepage: www.jpharmsci.org

Drug DiscoveryeDevelopment Interface

Physicochemical QSAR Analysis of Passive Permeability Across Caco-2 Monolayers Kiril Lanevskij 1, 2, *, Remigijus Didziapetris 1, 2 1 2

V
sI˛„Auk
stieji algoritmai“, A.Mickevi
ciaus 29, LT-08117 Vilnius, Lithuania ACD/Labs, Inc., 8 King Street East, Toronto, Ontario M5C 1B5, Canada

a r t i c l e i n f o

a b s t r a c t

Article history: Received 17 September 2018 Revised 5 October 2018 Accepted 5 October 2018

Caco-2 cell line is frequently used as a simplified in vitro model of intestinal absorption. In this study, a database of 1366 Caco-2 permeability coefficients (Pe) for 768 diverse drugs and drug-like compounds was compiled from public sources. The collected data represent permeation rates measured at varying experimental conditions (pH from 4.0 to 8.0, and stirring rates from 0 to >1000 rpm) that presumably account for passive diffusion across mucosal epithelium. These data were subjected to multistep nonlinear regression analysis using a minimal set of physicochemical descriptors (octanol-water log D, pKa, hydrogen bonding potential, and molecular size). The model was constructed in a mechanistic manner incorporating the following components: (i) a hydrodynamic equation of size- and chargespecific along with nonspecific diffusion across the paracellular pathway; (ii) transcellular diffusion represented by thermodynamic membrane/water partitioning ratio; (iii) stirring-dependent limit of maximum achievable permeability due to the presence of unstirred water layer. The obtained model demonstrates good accuracy of log Pe predictions with a residual mean square error <0.5 log units for all training and validation sets. Given its robust performance and straightforward interpretation in terms of simple physicochemical properties, the proposed model may serve as a valuable tool to guide drug discovery efforts toward readily absorbable compounds. © 2018 Published by Elsevier Inc. on behalf of the American Pharmacists Association.

Keywords: absorption Caco-2 cells computational ADME drug transport in silico modeling nonlinear regression passive diffusion permeability coefficient physicochemical properties QSAR

Introduction Intestinal absorption and oral bioavailability are key properties of any new drug candidates because oral dosage is the preferred route of drug administration.1 In vivo experimental determination of these parameters is not plausible until the latest stages of preclinical development; therefore, various in vitro epithelial cell lines have long been used to substitute these data.2-4 One of the most

Abbreviations used: Daq, aqueous diffusion coefficient; Do/w, pH-dependent octanol/water distribution coefficient; N, number of compounds; NHA, number of Hbond acceptors; NHD, number of H-bond donors; NLS, nonlinear least squares; Po/w, octanol/water partitioning coefficient of neutral species; pKa, acidic ionization constant; Ppara, paracellular permeability; Pe, effective permeability coefficient; Ptrans, transcellular permeability; PUWL, permeability of the unstirred water layer; R, average pore size; RMSE, residual mean square error; Vx, McGowan characteristic volume; Dj, transepithelial potential drop; ε/d, ratios of epithelium porosity to the tortuosity factor of the pores. Conflicts of interest: None. This article contains supplementary material available from the authors by request or via the Internet at https://doi.org/10.1016/j.xphs.2018.10.006. * Correspondence to: Kiril Lanevskij (Telephone: þ370 5 262 3408). E-mail address: [email protected] (K. Lanevskij).

popular surrogate end points is human colonic adenocarcinoma (Caco-2) cell line.2 Caco-2 cells grown on porous filter supports spontaneously differentiate into polarized monolayers with microvillous “brush border” on the apical side, tight junctions between adjacent cells, and enzyme expression levels that are sufficiently well representative of in vivo enterocytes.5,6 Drug molecules may cross the epithelial barrier by a variety of different mechanisms.2 The compounds may passively diffuse across the monolayer, either directly through the cells (transcellular route) or through the water-filled pores in tight junctions interconnecting the cells (paracellular route). Certain chemicals that are substrates of various transporter proteins expressed in the mucosal epithelium may also undergo facilitated diffusion or be actively carried against the concentration gradient either in the apical to basolateral (absorptive) or in the basolateral to apical (secretive) direction.7 Finally, some chemicals may permeate by receptormediated endocytosis or transcytosis, although this pathway is more common to peptides and is of less importance to small molecules.2,7-9 The overall permeation rate observed for a compound is the cumulative outcome of all involved mechanisms, but passive permeability has a distinctive role because it is an intrinsic property common to all molecules, governed by their physicochemical

https://doi.org/10.1016/j.xphs.2018.10.006 0022-3549/© 2018 Published by Elsevier Inc. on behalf of the American Pharmacists Association.

2

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

Table 1 Summary of Recent QSAR Models of In Vitro Permeabilities Across Cell Monolayers Study

Cell Line

Methods

Descriptors

N

R2

RMSE

Sherer et al., 201212 Wang et al., 201611

LLC-PK1, Caco-2 Caco-2

RF MLR, PLS, SVM, boosting

MOE 2D descriptors, atom pair parameters MOE 2D and 3D descriptors

Over et al., 201613

MDCK Caco-2 (macrocycles)

0.52, 0.64 0.79-0.83 0.75-0.81 0.75 0.69 0.81-0.84

0.20, 0.62 0.31-0.34 0.31-0.36 0.36 0.38 0.41-0.50

Fredlund et al., 201714 This work

Caco-2 (“intrinsic”) Caco-2

Cross-prediction using Caco-2 model PLS, RF, SVM logP/logD, PSA, NHD, NHA, charges, structural and topological parameters Same with exp. logD PLS, RF, SVM AZ descriptor set Nonlinear least squares log Do/w, pKa, NHD, Vx

15,791, 313 1017a 255b 298c 220 214c 214c 284b 497a 442b 427c

0.82-0.93 No data 0.80 0.77 0.53

0.36-0.43 0.45 0.44 0.49 0.47

a b c

Training set. Test set. External validation set.

characteristics. When its rate is sufficiently high, passive diffusion may account for the major part of total permeability even with simultaneous occurrence of carrier-mediated processes.7,8 Given the vast numbers of new chemical entities explored in modern drug discovery projects, performing even simpler in vitro transport experiments for all proposed molecules is not feasible, and there is significant demand for reliable computational approaches allowing to estimate absorption potential of new compounds before their synthesis.10 A detailed review of previous computational studies dealing with prediction of permeabilities of chemicals in Caco-2 or other cell lines can be found in a publication by Wang et al.,11 whereas several recent works in this field have been summarized in Table 1 of this article. Wang and coauthors themselves performed a comprehensive evaluation of several machine learning techniques using a training set of more than 1000 entries and a selection of more than 300 two- and three-dimensional molecular descriptors. Overall, they obtained consistent results for both training and validation sets with a residual mean square error (RMSE) <0.4 log units, providing a good indication of the accuracy of fit achievable by the current state-of-art statistical methods. Over et al.13 focused on a specific class of macrocyclic molecules, which may exhibit regio-, stereo-, and conformation-specific trends in permeability. Their models relating Caco-2 permeability to common physicochemical properties and other characteristics describing the composition and the topology of the macrocycles were able to reach predictivity levels similar to those reported in the study by Wang et al. 11 when experimental log Do/w was used as a descriptor instead of calculated values. Fredlund et al.14 specifically targeted their efforts at passive permeability and proposed a so-called “intrinsic” transport measure, expressed as Caco-2 permeability determined in the presence of an inhibitor cocktail of 3 major efflux transporters. In silico analysis based on >2500 in-house compounds yielded a model with hydrogen bond donor properties, lipophilicity, and positive charge among the most important descriptors. Finally, a slightly older work by Sherer et al.12 is notable for the size of the data set. Although primarily based on a different cell line (LLC-PK1), with more than 15,000 proprietary compounds, this is one of the largest data sets in cell permeability studies described in the literature up to date. Here, the authors also aimed to study passive permeability and excluded compounds with bidirectional permeability ratios outside the 0.5-2.0 interval. Such prefiltering and the overall breadth of the covered chemical space allowed observing a clear parabolic trend in the relationship between permeation rates and chromatographically measured log Do/w values.

Many QSAR studies use purely statistical approaches aiming to identify the best combination of descriptors and modeling methods that maximizes the predictive performance characteristics without going into deeper detail about the involved transport mechanisms. However, because a significant part of variation in compound permeabilities across epithelial cells is determined by nonspecific diffusion processes, it is beneficial to separate the effects of passive and active transport routes and to study passive diffusion pathways from a mechanistic point of view.7 Even if slightly less accurate, such models have an advantage in their clear interpretability and ability to guide drug discovery efforts toward the molecules with favorable physicochemical profiles. In the field of cellular absorptive systems, this principle has been adopted by Avdeef and coworkers who developed in combo approach incorporating a biophysical hydrodynamic model of paracellular pathway and experimentally determined permeabilities in PAMPA assay as a model of transcellular pathway.15,16 In our previous works, we used a similar mechanistic approach and performed ionization-specific physicochemical analysis of in vivo absorption in human jejunum,17 and in situ permeability across the blood-brain barrier.18 The objective of our present study was to extend this approach to in vitro Caco-2 system, and to obtain a reliable and interpretable in silico model of transepithelial permeation rates that would be suitable for evaluating the passive absorption potential of new chemical entities, validating experimental measurements, and as a reference permeability estimate for studying active transport processes.

Theory A system of equations for calculating Caco-2 permeability coefficients was constructed on the basis of common theoretical knowledge of absorptive processes. The overall resistance to diffusion, that is, the inverse of permeability, can be calculated as a sum of resistances of individual components of the system19:

1 1 1 ¼ þ Pe PUWL Ppara þ Ptrans

(1)

Here, Pe is the effective permeability coefficient, PUWL is the permeability of the unstirred water layer (UWL), and the final term corresponds to permeability of the cell layer, which in turn can be split into the contributions of paracellular (Ppara) and transcellular (Ptrans) pathways.

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

UWL Permeability The UWL resistance depends on the actual thickness of the aqueous boundary layer, which can be altered by changing the efficiency of stirring. The relationship between stirring rate and maximum achievable UWL-controlled permeability is usually described by the following equation19,20:

PUWL ¼ Kna

(2)

Ppara ¼

1=3 Cpara1 ,V x ,F

þ


r  R

3

" , fo þ

X i

kzi fi , 1  ekzi

#

1=3 Cpara2 ,V x

(7)

Transcellular Permeability

Here K is an empirically fitted constant that reflects diffusivity of the solute, kinematic viscosity, and geometric factors; n is the stirring rate; the exponent a reflects the type of hydrodynamics and is expected to have values between 0.5 and 1.0.20 We used a slightly modified equation, assuming linear dependence (a ¼ 1) for simplicity, and incorporating an additional term Ko corresponding to the minimal UWL permeability that would be observed in the absence of stirring:

According to the modified Fick's law, transmembrane permeability of the molecule can be expressed by the product of its kinetic diffusion coefficient Dk and membrane/water partitioning ratio Kd, corrected by membrane thickness h19:

PUWL ¼ Ko þ Kn

The variance in Dk can be modeled by McGowan volume (Vx) in the same way as shown above for Daq, while h, as a property of the membrane rather than the drug itself, may be incorporated into the constant term Ctrans. Finally, the transmembrane partitioning coefficient Kd can be approximated by the partitioning ratio in a reference end point, such as the frequently used octanol-water system. The intrinsic permeability of neutral species could then be described by Equation 9:

(3)

Paracellular Permeability Ppara term was expressed by a theoretical model of molecular size-restricted diffusion through aqueous pores in tight junctions originally devised by Adson et al.,21,22 and modified by Avdeef and Tam23:

Ppara ¼


ε d

1

,Daq ,F


r  R

" , fo þ

X i

#
ε kzi fi , ,Daq þ d 2 1  ekzi

(4)

Here (ε/d)x terms are the ratios of epithelium porosity to the tortuosity factor of the pores, Daq is the the aqueous diffusion coefficient of the molecule, and F is the Renkin molecular sieving function21,22 that decreases with increasing molecular size, approaching zero as the molecule's hydrodynamic radius (r) approaches the average radius of the pore (R):

F



r 3
r 5  r 2 r ¼ 1 1  2:104 þ 2:09  0:95 R R R R R


r 

(5)

kzi The 1e kzi term in Equation 4 corresponds to the electrochemical energy function that has distinct values for ionic forms (i) with different net charges (fi denote their fractions in the solution) and is equal to 1 for uncharged species (fo). The net charges are denoted by zi, and k is a constant reflecting the electrochemical characteristics of the monolayer, in particulardthe transepithelial potential drop.21,22 The original model of Adson accounts for only one type of pores characterized by (ε/d)1, whereas the modification from the study by Avdeef and Tam23 assumes the presence of an additional low-capacity size- and charge-unspecific pathway (ε/d)2. The latter term proves to be useful to account for nonzero permeability of large hydrophilic monolayer integrity markers, such as oligosaccharides, or Lucifer Yellow. Finally, according to the Stokes-Einstein equation, the aqueous diffusion coefficient (Daq) is inversely proportional to the molecule's hydrodynamic radius (r),24 which can be approximated from the cubic root of its McGowan characteristic volume that can be easily calculated from structure25:

r¼

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 3Vx 166 4p

Aggregating all constant terms together yields Equation 7:

(6)

Ptrans ¼

  1 1=3 D K ¼ Ctrans V x Kd h k d

(8)

log K od ¼ co þ c1 logPo=w þ c2 NHD

(9)

Here cx values correspond to model coefficients, log Po/w is the octanol/water partitioning ratio, and NHD stands for the number of hydrogen bond donors in the moleculedan additional descriptor to correct for excessive hydrogen bonding potential of octanol. The choice of the hydrogen bonding parameter will be further discussed in the Results section. It is also important to note that permeability-lipophilicity dependence may not necessarily be linear throughout the entire scale of parameter values. This can be related to hydrophobic entrapment of highly lipophilic compounds in the membrane and may be expressed by a bilinear model first introduced by Kubinyi.26 We also incorporated this effect into our model by inclusion of an additional term with a and b coefficients standing for the second (falling) log Po/w slope and the inflection point, respectively:


log K od ¼ co þ c1 logPo=w þ a log 1 þ 10logPo=w b þ c2 NHD (10) To account for the effect of ionization, there are 2 straightforward approaches. One obvious choice is to substitute log Po/w of neutral species with pH-dependent distribution ratio log Do/w that represents the overall partitioning coefficient as a sum of partitioning coefficients (Pi) of all ionic species in the solution in accordance with thermodynamic partitioning theory of electrolytes:

logDo=w

X ¼ log fi ,Pi

!

X fi ,10logPo þDi ¼ log

i

! (11)

i

Another plausible approach would be to apply the equation of the same form directly to the analyzed Kd coefficient and to determine Di values specific to the system under consideration:

X logKd ¼ log fi ,K id i

!

X o fi ,10log K d þDi ¼ log i

! (12)

4

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

The second approach is preferable in general case because Di values in any particular end point may differ from those observed in octanol-water system. This is especially true for tissue distribution end points, where factors other than membrane permeation may play significant roles, and the overall pattern of ionization dependence may be completely different.27 In absorptive systems, a qualitatively similar profile of Di could be expected, but precise numerical values may still be somewhat different. For example, in our previous studies focusing on in vivo systems (rodent blood-brain barrier penetration,18 and human intestinal absorption17), the estimated differences between neutral and charged species were only about 2-3 log units compared to 3-4 log units in octanol-water system.19 However, these are normalized values assuming one-to-one correspondence between log P/log D and diffusion rates, whereas on real-world data sets, fitted c1 coefficients from Equation 10 are usually significantly less than unity (e.g., 0.5-0.7 in Refs.17,18), and the differences between systems become much less pronounced, falling to the level of experimental uncertainty. Preliminary analysis of our current Caco-2 permeability data set yielded very similar predictivities of models based on Equations 11 and 12. In these circumstances, opposite to our previous works, we used the first approach because inferring the effect of ionization from log Do/w has several practical advantages:  Reducing the complexity of calculations and the number of degrees of freedom in the model  Straightforward ability to evaluate permeation rates of compounds with multiple ionizable groups, whereas obtaining accurate end pointespecific Di estimates for multiprotic acids or bases is problematic due to the lack of reliable data  Most new publications dealing with characterization of novel compounds measure lipophilicity in terms of log Do/w. Using the same quantity as a descriptor in the model opens broader opportunities of using in combo approach3 to improve the accuracy of calculations by replacing predicted parameter values by experimentally determined ones. Therefore, the final equation used for calculating transcellular permeability looks as follows: 1=3

Ptrans ¼ V x

,10c0 þc1



b

logDo=w þc2 NHD þa log 1þ10logDo =w

 (13)

Data and Methods Data Compilation The observed values of permeability coefficients of drugs, druglike compounds, and other small molecules were collected from the original publications starting from early works of Artursson’s group28,29 up to the most recent studies, where Caco-2 permeability assay is typically performed as a part of ADME profile evaluation for congeneric compound series. Data from earlier QSAR studies or other compilations were not included in our database unless these records could be verified in the original sources. All permeability coefficients were recorded together with experimental conditions at which they were obtained, namely pH of the medium on the apical side, and stirring rate (in rpm) if stirring was used. Therefore, the database contains multiple entries for some compounds, representing either different combinations of pH and stirring rate, or the same conditions, but measured in different laboratories. Data Verification

indeed represent the rate of passive diffusion by paracellular or transcellular route and are not affected by carrier-mediated influx or efflux processes. In those cases when the original article provided permeabilities in both directionsdapical to basolateral (Pab) and basolateral to apical (Pba)dthe potential involvement of active processes could be determined from Pba/Pab ratios. Accordingly, entries with Pba/Pab < 0.5 (influx) or Pba/Pab > 2.0 (efflux) were excluded from the database. When such details were not available, decisions were made after preliminary analysisdcompounds demonstrating large deviations from basic physicochemical trends and known as substrates of P-gp, BCRP, MRPs, PepT1, MCT1, OCTs, or other transporters expressed in Caco-2 cells were also excluded from modeling. Furthermore, several suspicious entries were removed on account of other issues, such as big discrepancies with other results for the same compound at similar conditions, large inconsistencies with in vivo absorption efficiency, evidence of extremely leaky epithelium (unusually high Pe values for large hydrophilic molecules), and so on. The final curated data set consisted of 939 entries for 529 diverse chemicals covering the pH range from 4.0 to 8.0 and stirring rates from 0 to 1090 rpm. Preparation for Modeling To avoid artifacts that could potentially arise from residual active transport effects or poor physicochemical property predictions, for model development purposes, we extracted about a half of the database subjectively classified as higher quality data. These mostly represented marketed drugs and other small molecules with known experimental log Po/w and/or pKa values and considered to be reasonably free of unwanted effects. The remaining data collected from more recent studies with less detailed assay descriptions were reserved for model validation purposes. Although this test set also contained some entries for known drugs, the overall chemical space covered by these data was shifted toward more novel scaffolds. Finally, another portion of experimental data was collected during model development stage. This subset was based exclusively on articles dealing with optimization of novel compound series published in the last few years, which had undergone slightly less stringent criteria for inclusion in the database. In particular, some publications provided permeability coefficients without full account for assay conditions used in measurements. If stirring rate was not indicated, it was assumed that the experiments had been conducted without stirring. In cases

Table 2 Optimized Coefficients of Physicochemical Descriptors in the Nonlinear Least Squares Model of Caco-2 Permeability Expressed as Pe, 106 cm/s Coefficient UWL permeability Ko K Paracellular permeability Cpara1 Cpara2 R

k Transcellular permeability c0 c1 c2

a

Because this study focuses on passive transport across mucosal epithelium, it is important to ensure that the collected values

b

Meaning

Value

UWL-limited permeability at no stirring Effect of stirring rate increase by 1 rpm

43 0.4

Scaling factor for charge- and sizeselective pores Scaling factor for low-capacity nonselective pores Average pore radius, Å Electrochemical energy factor

160

Scaling factor for Ptrans Effect of lipophilicity before optimum Effect of hydrogen bond donor count Hydrophobic entrapment (slope ¼ c1 þ a) Bilinear inflection point

0.2 5.0 1.8 2.0 0.5 0.4 1.0 3.0

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

when even pH information was missing, such entries were included only if full assay description could be retrieved from earlier publications by the same group. It was then presumed that the same protocol was used in later studies (typically, pH ¼ 7.4, no stirring). The resulting data set consisting of 427 additional entries was used as an external validation set with an intention to evaluate the potential utility of the obtained predictive model in drug discovery projects. The entire Caco-2 permeability data set used in the present study is available as Supplementary Information. Model Development The predictive model was obtained by nonlinear least squares fitting of training set data to the system of Equations 1, 3, 7, and 13 in several steps. First, the upper limits of permeability imposed by UWL were estimated using data from the study by Karlsson and Artursson.28 The contribution of size- and charge-dependent paracellular pathway was evaluated on a set of 45 data points for small

hydrophilic molecules known to be transported by paracellular route30,31 and several other entries where preliminary analysis indicated prevalence of paracellular over transcellular pathway. The optimal inflection point value in the log Do/w versus log Pe relationship was estimated by visual inspection of property dependence plots (e.g., see those presented in the Results and Discussion section) as it could not be fitted reliably due to the low proportion of highly lipophilic molecules in the data set. The remaining transcellular permeability coefficients were fitted using the entire training set, accounting for PUWL and Ppara values calculated as outlined previously. Software Database manipulations and descriptor calculations were performed using ACD/Algorithm Builder32 and ACD/Percepta33 programs. Subsequent model development and evaluation was carried out in R statistical environment.34

Distribution by log Do/w

a

Observed vs. predicted log Pe

d

Training Set

Training Set 4

140

3

log Pe (observed)

120

N

100 80 60

2 1 0

-2

40 20

-1

-1

0

< - -4 - -3 - -2 - -1 - 0

1

2

3

4

4

log Pe (predicted)

e

Test Set

Test Set 3 2

log Pe (observed)

120 100

N

3

y = 0.9589x + 0.0393 R² = 0.8039

140

80

60

1 0

-2

40 20

-1

< - -4 - -3 - -2 - -1 - 0

1

2

3

4

0

1

2

3

-1

-2

0

y = 0.8665x + 0.1998 R² = 0.7651

-3

>

log Do/w

log Pe (predicted)

f

External Set

External Set 2.5

140

2

log Pe (observed)

120

100

N

2

-3

>

log Do/w

c

1

-2

0

b

5

80 60 40

-1.5

20

1.5 1 0.5 -1

0 -0.5 0 -0.5

-1

0 < - -4 - -3 - -2 - -1 - 0

log Do/w

1

2

3

4

>

-1.5

0.5

1

1.5

2

2.5

y = 0.9019x - 0.1411 R² = 0.5301

log Pe (predicted)

Figure 1. Distribution of octanol/water log D values (a-c) and comparison of experimental versus predicted log Pe (d-f) for training, test, and external validation set compounds.

6

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

Results and Discussion

Transcellular Permeability: Lipophilicity and H-Bonding

The fitted model coefficients from Equations 3, 7 and 13 are presented in Table 2, whereas the populations of data subsets by log Do/w values, as well as the relationships between predicted and observed log Pe for all training and validation sets are displayed in Figure 1. The individual components of the model equation are discussed in the following.

As expected, the transcellular diffusion rate is dominated by lipophilicity. The hydrogen bond donating potential also has a significant contribution with inclusion of each additional mobile hydrogen atom (OH or NH) leading to an average 2.5 times decrease in Pe values. When hydrogen acceptor count (NHA) descriptor was added to the model, it did not significantly increase the predictive power of the model, and the absolute value of its fitted coefficient was less than 0.1. This observation is consistent with recent findings of Over et al.13 who showed that even for large macrocyclic molecules, the relationship of log Pe with log Do/w becomes evident when the compounds are split into groups by NHD values, while the trend is much less clear for NHA. We tried using other hydrogen bonding descriptors, such as topological polar surface area,37 or Abraham's solvation parameters,38 but these also did not enhance predictivity, so that the simple H-donor count was preferred. The estimated log Do/w dependence profiles for groups of compounds with the same NHD values are depicted in Figure 3. It is quite clear that permeability increases with lipophilicity only up to an optimum point (transformed into a plateau at small NHD due to UWL limit) and then decreases at a similar rate. The overall effect of H-bonding (in terms of NHD) observed here is very similar compared to our previous work dealing with in vivo jejunal permeability.17 The lipophilicity slope (c1) is somewhat lower in the present study, but the respective values are not directly comparable due to the differences in the employed ionization model. Importantly, the bell-shaped curves analogous to those presented in Figure 3 could not be obtained in the former work because the molecules within a very wide range of lipophilicities tend to have 100% absorption, and thus, fraction-absorbed data did not reveal the full pattern of lipophilicity-permeability relationship.

UWL Permeability The rate of diffusion through unstirred water layer was analyzed using the data from the study by Karlsson and Artursson,28 where permeabilities of several lipophilic compounds were determined at a variety of stirring rates (y) from 0 to 1090 rpm. The correlation between y and Pe is shown in Figure 2. Except for several outliers (excluded from the plot), a linear relationship is observed indicating that in the absence of stirring, Pe barely exceeds 40 106 cm/s, but this PUWL limit is almost doubled when 100 rpm stirring is applied. Paracellular Permeability The parameter values reported in Table 2 correspond to an average pore radius R ¼ 5.0 Å, and transepithelial potential drop Dj ¼ 48 mV. It should be noted that in a global QSAR model parameterized on data from multiple sources like the one presented here, these parameters should not be interpreted with any special meaning. In fact, each particular cell line used in a certain laboratory can be characterized by its own combination of R, Dj, and (ε/d)x values (the latter ratios are incorporated into Cpara constants), so that the values reported here merely reflect some “average” characteristics that work well for an aggregated set of interlaboratory data. Nevertheless, it is worth pointing out that the obtained values are consistent with the results of Sugano et al.,24 who reported R ¼ 5.6 Å for in vivo human intestine model, given that leakiness of colonic epithelium as in Caco-2 cells is likely not significantly different from jejunum.35 Moreover, these values are in good agreement with many individual cell lines in the published literature (e.g.,21,31,36), and the average potential drop of 43 mV as summarized by Avdeef35 is very close to Dj ¼ 48 mV obtained in the present study.

Ionization In theory, in vitro cell lines such as Caco-2 could provide great potential for studying the influence of the compound's ionization state on its membrane permeation rate since contrary to in vivo experiments, here one can arbitrarily adjust pH of the medium and measure permeabilities of the same molecule throughout the entire range of physiologically relevant pH values. However, in practice, quite few such studies can be found in the literature. In

PUWL vs. stirring rate 600

Testosterone Warfarin

500

PUWL, 10-6 cm/s

Antipyrine Alprenolol

400

Propranolol 300

Linear (Overall)

200 100

y = 0.3986x + 43.252 R² = 0.8589

0 0

200

400

600

800

1000

1200

ν, rpm Figure 2. Relationship between stirring rate (n) and UWL-limited Caco-2 permeability for a selection of compounds from the study by Karlsson and Artursson.28

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

most publications, at most 2 pH points are used to determine permeation rates. Yet, a few comprehensive sets of pH-dependent Caco-2 permeability measurements for several marketed drugs can be found in the works of Palm et al.39 and Neuhoff et al.40,41 In Figure 4, experimental pH-dependent log Pe values for these drugs are plotted against their log Do/w, and a clear sigmoid relationship emerges for a group of compounds. A relatively hydrophilic drug cimetidine forms an inflection in the bottom part of the curve because at lower pH, a significant contribution to its permeation rate comes from the paracellular pathway. More hydrophobic molecules of alfentanil (base) and indomethacin (acid) diffuse by transcellular route, and in the latter case, permeability apparently reaches the UWL limit. Other compounds significantly deviate from this sigmoid, having permeabilities much higher than could be anticipated from their physicochemical characteristics. This is not surprising for salicylic acid, which is a known substrate of several proton-dependent transporters.41 However, the outliers also include 2 b-blockers, atenolol and metoprolol, which are sometimes referred as markers of passive paracellular and transcellular permeability, respectively.42 Yet, recent studies indicate that bblocker class drugs tend to interact with organic cation transporters. Some of these compounds inhibit cellular uptake of standard substrates of several OCTs and OCTNs,43 and atenolol itself exhibits saturable uptake by OCT1.44,45 There is less such evidence for metoprolol, but it has comparable permeability to a similar drug propranolol despite being 1.5 units less lipophilic and having almost identical pKa and hydrogen bonding pattern.8 Therefore, available data indicate that metoprolol is very likely to cross the cell layer by the means of facilitated diffusion, at least at conditions close to physiological. The authors of the aforementioned publications39-41 also attempted to quantify the contributions of ionized species to the overall diffusion rates. Because available data are scarce and subject to experimental uncertainty, precise quantification is difficult, and values obtained for different compounds vary

significantly. Still, the reported results suggest that neutral species can be transported about 102 to 104 times faster than cations or anions, so that in a global model, inferring the ionization effect from octanol/water log D ratio is a reasonable approximation. Model Performance The statistical characteristics listed in Table 2 indicate that the obtained model can predict passive permeability coefficients of drug-like molecules across Caco-2 monolayers with RMSE about 0.5 log units. The accuracy of predictions depends on the quality of underlying physicochemical parameters, particularly log Do/w, and the slight difference between training and test sets is in part related to the fact that less experimental physicochemical characteristics were known for test set compounds. A relatively low R2 value is obtained for external validation set due to a much narrower variation range of log Do/w and log Pe values, as the external set contains quite few hydrophilic poorly permeable chemicals, which can be clearly seen by comparing Figure 1c to Figures 1a and b. This is understandable, as drug discovery projects usually aim for readily absorbable molecules. It is also evident that many compounds are predicted to have UWL-limited permeability, whereas their experimental values exhibit substantial variability and form a distinct vertical line in Figure 1f. This can be explained by a variety of reasons, such as inaccurate log Do/w predictions for novel scaffolds, unaccounted active transport effects, or even issues with experimental data (e.g., when stirring conditions are not reported, but Pe exceeds PUWL at n ¼ 0). In any case, the overall RMSE of predictions is on par with training and internal test set compounds indicating that the model is not overfitted. Interestingly, in a recent evaluation of interlaboratory variability in experimental Caco-2 data, it was estimated that the standard deviation of published measurements for the same compound can reach almost 0.6 log units.46 These data account for all possible transport mechanisms,

NHD = 1

NHD = 2

3

3

2.5

2.5 2

1.5

log Pe

log Pe

2

1 0.5

-0.5

1.5 1 0.5 0

-2

0 -2

0

2

4

6

-0.5

0

-1

log Do/w

3

2.5

2.5

2

2

1.5

1.5

1 0.5 0 -0.5 0

2

6

4

6

log Do/w

4

6

-2

1 0 -0.5 0 -1

-1.5

-1.5

log Do/w

4

0.5

-1 -2

2

NHD = 4

3

log Pe

log Pe

NHD = 3

-2

7

-2

2

log Do/w

Figure 3. Bilinear relationship between Caco-2 permeability coefficients and octanol/water log D for subsets of our database grouped by number of H-Donor sites (NHD) in the molecules. The dotted lines indicate the UWL limit at n ¼ 200 rpm, and solid curves represent model predictions at these conditions.

8

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9

log Pe vs. log Do/w Cimetidine

3

Alfentanil Indomethacin

2

Atenolol

log Pe

Metoprolol

1

Salicylic acid

0 -4

-3

-2

-1

0

1

2

3

4

-1

-2

log Do/w Figure 4. Relationship between Caco-2 permeability determined at multiple pH points, and octanol/water log D for passively transported compounds (filled points) and compounds that presumably exhibit active influx across mucosal epithelium (blank points). Experimental data taken from the studies by Palm et al. and Neuhoff et al.39-41

and such large variability can be caused by, for example, differences in expression levels of transporter proteins, so that better data consistency can be expected for passively diffusing chemicals. Still, in light of these results, the predictive power of the passive permeation model presented in the present study can be considered close to the level of interlaboratory variation. In Combo Approach As already mentioned in Theory section, usage of interpretable physicochemical descriptors makes it possible to enhance the accuracy of predictions by substituting calculated input parameters (such as log Do/w) with experimentally measured values. Moreover, this in combo approach provides opportunities to guide research efforts toward more permeable regions of chemical space, as one can easily simulate what effect structural modifications would incur on the molecule's absorptive potential by performing the calculations using custom values of lipophilicity, H-bonding, molecular size, and indicating the relevant experimental conditions (pH and stirring rate).

reach UWL permeability limits if no carrier-mediated efflux occurs. NHD  4 is highly detrimental for permeability, and in this case, compounds with log Do/w outside of the optimum range are likely to exhibit low permeation rates (Pe < 1 106 cm/s), unless the molecules are small enough for a significant contribution of paracellular route to be observed. Very high lipophilicities may result in lower permeabilities due to the hydrophobic entrapment in the membranes of the epithelial cells or other factors, such as adsorption onto the experimental device.4 The results of internal validation and external validation indicate that the derived model can estimate passive permeability in Caco-2 monolayers close to the level of interlaboratory variability, especially when measured values of physicochemical input parameters (primarily log Do/w) are available. Strong deviations from model predictions may signify involvement of carrier-mediated processes or other confounding factors. Therefore, the presented model may not only provide a reliable tool for early estimation of passive absorption potential of new chemicals but also aid in validation of experimental permeability measurements, as well as identification and quantification of the contributions of influx or efflux processes.

Conclusion A mechanistic model of drugs' absorption in Caco-2 monolayers has been developed in the present study relating passive permeability by paracellular and transcellular routes to key physicochemical parameters of molecules, such as lipophilicity, ionization, hydrogen bond donating potential, and molecular size. Paracellular diffusion was described by a theoretical hydrodynamic model accounting for 2 types of pores in the tight junctions between the cells. According to the fitted coefficients, mono-cations diffuse through charge-selective pores about 2 times faster, and monoanions about 3 times slower than neutral molecules. It has also been demonstrated that transcellular diffusion, which accounts for >80% of total variation in permeation rates, can be successfully approximated by a bilinear relationship with octanol/water distribution ratio (log Do/w) corrected for excess hydrogen bonding. According to the proposed model, the optimal log Do/w for absorption lies in the interval between 2 and 4. The molecules with lipophilicities in the optimal range and NHD  1 can be expected to

Acknowledgments This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors. Supporting information available: Caco-2 permeability data set used in the present study along with predictions performed by the described model. References 1. Lin L, Wong H. Predicting oral drug absorption: mini review on physiologicallybased pharmacokinetic models. Pharmaceutics. 2017;9(4):41. https://doi.org/ 10.3390/pharmaceutics9040041. 2. Artursson P, Palm K, Luthman K. Caco-2 monolayers in experimental and theoretical predictions of drug transport. Adv Drug Deliv Rev. 2001;46(1-3):2743. 3. van de Waterbeemd H. Which in vitro screens guide the prediction of oral absorption and volume of distribution? Basic Clin Pharmacol Toxicol. 2005;96(3):162-166.

K. Lanevskij, R. Didziapetris / Journal of Pharmaceutical Sciences xxx (2018) 1-9 4. Volpe DA. Variability in caco-2 and MDCK cell-based intestinal permeability assays. J Pharm Sci. 2008;97(2):712-725. 5. Hubatsch I, Ragnarsson EGE, Artursson P. Determination of drug permeability and prediction of drug absorption in caco-2 monolayers. Nat Protoc. 2007;2(9): 2111-2119.  pez-Expo sito I, et al., eds. 6. Lea T. Caco-2 cell line. In: Verhoeckx K, Cotter P, Lo The Impact of Food Bioactives on Health: In Vitro and Ex Vivo Models. Cham, CH: Springer; 2015:103-111. 7. Sugano K, Kansy M, Artursson P, et al. Coexistence of passive and carriermediated processes in drug transport. Nat Rev Drug Discov. 2010;9(8):597-614. 8. Smith D, Artursson P, Avdeef A, et al. Passive lipoidal diffusion and carriermediated cell uptake are both important mechanisms of membrane permeation in drug disposition. Mol Pharm. 2014;11(6):1727-1738. €mer SD. Quantitative aspects of drug permeation across in vitro and in vivo 9. Kra barriers. Eur J Pharm Sci. 2016;87:30-46. 10. van de Waterbeemd H, Gifford E. ADMET in silico modelling: towards prediction paradise? Nat Rev Drug Discov. 2003;2(3):192-204. 11. Wang N-N, Dong J, Deng Y-H, et al. ADME properties evaluation in drug discovery: prediction of caco-2 cell permeability using a combination of NSGA-II and boosting. J Chem Inf Model. 2016;56(4):763-773. 12. Sherer EC, Verras A, Madeira M, et al. QSAR prediction of passive permeability in the LLC-PK1 cell line: trends in molecular properties and cross-prediction of caco-2 permeabilities. Mol Inform. 2012;31(3-4):231-245. 13. Over B, Matsson P, Tyrchan C, et al. Structural and conformational determinants of macrocycle cell permeability. Nat Chem Biol. 2016;12(12):1065-1074. 14. Fredlund L, Winiwarter S, Hilgendorf C. In vitro intrinsic permeability: a transporter-independent measure of caco-2 cell permeability in drug design and development. Mol Pharm. 2017;14(5):1601-1609. € J, Tavelin S. Caco-2 15. Avdeef A, Artursson P, Neuhoff S, Lazorova L, Gråsjo permeability of weakly basic drugs predicted with the double-sink PAMPA pKa(flux) method. Eur J Pharm Sci. 2005;24(4):333-349. 16. Avdeef A, Bendels S, Di L, et al. PAMPA–critical factors for better predictions of absorption. J Pharm Sci. 2007;96(11):2893-2909. 17. Reynolds DP, Lanevskij K, Japertas P, Didziapetris R, Petrauskas A. Ionization-specific analysis of human intestinal absorption. J Pharm Sci. 2009;98(11):4039-4054. 18. Lanevskij K, Japertas P, Didziapetris R, Petrauskas A. Ionization-specific prediction of blood-brain permeability. J Pharm Sci. 2009;98(1):122-134. 19. Avdeef A. Absoption and Drug Development. 1st ed. New York: Wiley-Interscience; 2003. € nkko €nen J. Analysis of unstirred water layer in 20. Korjamo T, Heikkinen AT, Mo in vitro permeability experiments. J Pharm Sci. 2009;98(12):4469-4479. 21. Adson A, Raub TJ, Burton PS, et al. Quantitative approaches to delineate paracellular diffusion in cultured epithelial cell monolayers. J Pharm Sci. 1994;83(11):1529-1536. 22. Adson A, Burton PS, Raub TJ, Barsuhn CL, Audus KL, Ho NF. Passive diffusion of weak organic electrolytes across caco-2 cell monolayers: uncoupling the contributions of hydrodynamic, transcellular, and paracellular barriers. J Pharm Sci. 1995;84(10):1197-1204. 23. Avdeef A, Tam KY. How well can the caco-2/Madin-Darby canine kidney models predict effective human jejunal permeability? J Med Chem. 2010;53(9):3566-3584. 24. Sugano K, Takata N, Machida M, Saitoh K, Terada K. Prediction of passive intestinal absorption using bio-mimetic artificial membrane permeation assay and the paracellular pathway model. Int J Pharm. 2002;241(2):241-251. 25. Abraham MH, McGowan JC. The use of characteristic volumes to measure cavity terms in reversed phase liquid-chromatography. Chromatographia. 1987;23:243-246. 26. Kubinyi H. Quantitative structure–activity relationships. 7. The bilinear model, a new model for nonlinear dependence of biological activity on hydrophobic character. J Med Chem. 1977;20:625-629.

9

27. Lanevskij K, Dapkunas J, Juska L, Japertas P, Didziapetris R. QSAR analysis of blood-brain distribution: the influence of plasma and brain tissue binding. J Pharm Sci. 2011;100(6):2147-2160. 28. Karlsson J, Artursson P. A method for the determination of cellular permeability coefficients and aqueous boundary layer thickness in monolayers of intestinal epithelial ( Caco-2) cells grown in permeable filter chambers. Int J Pharm. 1991;71(1):55-64. 29. Artursson P, Karlsson J. Correlation between oral drug absorption in humans and apparent drug permeability coefficients in human intestinal epithelial (Caco-2) cells. Biochem Biophys Res Commun. 1991;175(3):880-885. 30. Collett A, Sims E, Walker D, et al. Comparison of HT29-18-C1 and caco-2 cell lines as models for studying intestinal paracellular drug absorption. Pharm Res. 1996;13(2):216-221. 31. Knipp GT, Ho NF, Barsuhn CL, Borchardt RT. Paracellular diffusion in caco-2 cell monolayers: effect of perturbation on the transport of hydrophilic compounds that vary in charge and size. J Pharm Sci. 1997;86(10):1105-1110. 32. Japertas P, Didziapetris R, Petrauskas A. Fragmental methods in the design of new compounds. Applications of the advanced algorithm builder. Quant Struct Act Rel. 2002;21:23-37. 33. ACD/Percepta v. 2017.2.1. Toronto, Ontario, Canada: ACD/Labs, Inc.; 2017. Available at: http://www.acdlabs.com/products/percepta/. Accessed October 31, 2018. 34. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing; 2017. Available at: http://www.r-project. org/. Accessed October 31, 2018. 35. Avdeef A. Leakiness and size exclusion of paracellular channels in cultured epithelial cell monolayers-interlaboratory comparison. Pharm Res. 2010;27(3): 480-489. 36. Liang E, Chessic K, Yazdanian M. Evaluation of an accelerated caco-2 cell permeability model. J Pharm Sci. 2000;89(3):336-345. 37. Ertl P, Rohde B, Selzer P. Fast calculation of molecular polar surface area as a sum of fragment-based contribution and its application to the prediction of drug transport properties. J Med Chem. 2000;43:3714-3717. 38. Abraham MH. Application of solvation equations to chemical and biochemical processes. Pure Appl Chem. 1993;65:2503-2512. 39. Palm K, Luthman K, Ros J, Grasjo J, Artursson P. Effect of molecular charge on intestinal epithelial drug transport: pH-dependent transport of cationic drugs. J Pharmacol Exp Ther. 1999;291(2):435-443. 40. Neuhoff S, Ungell AL, Zamora I, Artursson P. pH-dependent bidirectional transport of weakly basic drugs across caco-2 monolayers: implications for drug-drug interactions. Pharm Res. 2003;20(8):1141-1148. 41. Neuhoff S, Ungell AL, Zamora I, Artursson P. pH-Dependent passive and active transport of acidic drugs across caco-2 cell monolayers. Eur J Pharm Sci. 2005;25(2-3):211-220. 42. Dixit P, Jain DK, Dumbwani J. Standardization of an ex vivo method for determination of intestinal permeability of drugs using everted rat intestine apparatus. J Pharmacol Toxicol Methods. 2012;65(1):13-17. 43. Grube M, Ameling S, Noutsias M, et al. Selective regulation of cardiac organic cation transporter novel type 2 (OCTN2) in dilated cardiomyopathy. Am J Pathol. 2011;178(6):2547-2559. 44. Mimura Y, Yasujima T, Ohta K, Inoue K, Yuasa H. Functional identification of organic cation transporter 1 as an atenolol transporter sensitive to flavonoids. Biochem Biophys Rep. 2015;2:166-171. 45. Mimura Y, Yasujima T, Ohta K, Inoue K, Yuasa H. Functional identification of plasma membrane monoamine transporter (PMAT/SLC29A4) as an atenolol transporter sensitive to flavonoids contained in apple juice. J Pharm Sci. 2017;106(9):2592-2598. 46. Lee JB, Zgair A, Taha DA, et al. Quantitative analysis of lab-to-lab variability in caco-2 permeability assays. Eur J Pharm Biopharm. 2017;114:38-42.