Mathematical and Numerical Methods for Modeling Drug Delivery to the Posterior Segment of the Eye
The high incidence of diseases that affect the posterior segment of the eye (PSE)-here intended as composed by sclera, choroid and retina-prompts for establishing effective and well tolerated therapies. Topical application (instillation of drops) and systemic assumption remain the most widespread drug administration routes. However, the drug achieved levels are not therapeutically sufficient therapeutic, since in the first case the drug is mainly washed away by different pathways (aqueous humor, systemic adsorption, tears) and in the second reaches the PSE target in a minimal fraction. Intravitreal injections and biodegradable episcleral implants have emerged in the last decades as alternative, more effective, administration routes. Whilst these techniques already offer significant improvements, much space is still open to research. The difficulty in delivering drugs to the PSE via the intravitreal (IV) or episcleral (EP) route resides in the several physical and dynamical barriers-including blood retinal barrier, clearance from choriocapillaries and lymphatics-which hinder the passage of the drug molecules. In combination, and in support of clinical experiments, mathematical and computational methods can be used to simulate drug concentration levels in the tissues upon IV or EP administration. In this respect, this short review gives an update of the most recent (from year 2000 to present time) developments on mathematical models for drug delivery to the PSE. We specifically focus our attention on physiological modeling works that include spatial dependency. Our review work is organized in short sections accompanied by detailed tables. We discuss descriptions of the PSE morphology, considering "anatomically accurate" as well as reduced models and we analyze the biophysical phenomena included in the examined models. We present the numerical techniques adopted to solve the resulting systems of partial differential equations and we deal with the delicate issue of parameter choice and validation of the results. Eventually, we examine the main results and scientific achievements of the considered models. Our conclusions point out the absence of a "standard mathematical model" and highlight a significant scattering of the results obtained from the different authors.
Drug delivery, Posterior segment of the eye, Episcleral implant, Intravitreal injection, Mathematical model
ASE: Anterior Segment of the Eye; BCs: Boundary Conditions; EP: Episcleral; FEM: Finite Element Method; IOP: Intraocular Pressure; IV: Intravitreal; NS: Navier Stokes Equations; PSE: Posterior Segment of the Eye; RPE: Retinal Pigmented Epithelium
The unique anatomy and physiology of the eye makes drug delivery to the PSE - here intended as composed by sclera, choroid and retina - a highly challenging task . In the past, drugs targeting the PSE were administered for the most part via topical or systemic routes, but severe limitations and inefficiency sources were associated with these therapeutic approaches. Intravitreal injections and (juxta-) episcleral implants have more recently emerged as alternative routes to achieve sustained therapeutic drug levels . The penetration rate of a certain drug through these latter routes is influenced by its solubility, hydro/lipophilicity, charge, degree of ionization and molecular size . Furthermore, drug concentration in the PSE tissues is not homogeneous but distributes in a complex way according to several concurrent biophysical phenomena as convection, diffusion, losses through anterior and posterior pathways, clearance mechanisms through blood and lymphatic vessels and degradation processes (see also Table 1 of . It has also been observed  that only a weak correlation exists between the actual local drug concentration in the PSE, which determines the pharmacological response, and the mean concentration values that experiments can monitor in the more accessible aqueous chamber and vitreous body under in vivo conditions. In this perspective, mathematical and computational models can offer important tools of analysis. Previous reviews specifically focusing on mathematical modeling of drug delivery to the PSE can be found in [4,5]. We refer to [2,6-9] for up-to-date, comprehensive, reviews on physiological and clinical aspects of this topic.
The present short review gives an update of the most recent (from year 2000 on) developments on mathematical modeling of drug delivery to the PSE. We specifically focus our attention on physiological modeling works that include spatial dependency in the interested domains. This choice stems from the above observations about the inhomogeneous distribution of drug in the complex PSE structures, which definitely prompts for including spatial dependency in the models. We do not include mathematical models tailored to study drug distribution in the sole vitreous, for which a specific and relatively vast literature exists, see [10,11]. In our analysis, we have identified 11 papers which fulfill the above described criteria, and whose content is synthetically reviewed in table 1, table 2 and table 3. Each table is focused on a particular aspect of mathematical modeling of drug delivery to the PSE: representation of the PSE geometry, modeling of the relevant biophysical phenomena and discretization techniques adopted to obtain the numerical solution of the model on a computer.
We organize the review as follows. In Sect. 2, we analyze and compare the geometrical modelling choices; in Sect. 3, we present the biophysical phenomena included in the models; in Sect. 4, we deal with the numerical techniques adopted to solve the models, and in Sect. 5 we discuss the main results and scientific achievements of the models. Eventually, in Sect. 6 we draw the conclusions.
Morphometric information on the PSE geometry has been traditionally obtained from dissection of specimina and, more recently, from digital imaging, both in animals and humans. In establishing a geometrical model for the PSE, important issues must be faced: i) which are the relevant anatomical structures to consider as mathematical domains of the model? ii) What are anatomically correct measures of the included volumes? (for example, the volume of the vitreous chamber considered in Jooybar, et al.  is equal to 4.4 ml , where as in the original model of , the volume of the vitreous is about 2.2 ml); iii) what are the significant clearance pathways (for example, blood vessels, periocular structures or lymphatic vessels) for a certain drug administration mode? iv) how to represent un-resolved structures via boundary conditions on the external surfaces? How to deal with internal interface/coupling conditions? In addition, one should consider the fact that the PSE has unique geometrical features in each subject. In this respect, the common procedure consists in extensive literature searches and data comparisons so to obtain models based on a representative, average, geometry for a certain species (human, rabbit, mouse, monkey.
In table 1, we provide a comparative analysis between the geometrical models of the PSE adopted in the considered mathematical works (being for the totality 3D or 1D models, except for the 2D model of . The most part of the papers considers geometrical parameters corresponding to the human eye, possibly partially derived from other species (rabbit eye, see also the comparative analysis carried out in . Moreover, for the most part, these recent models consider, at least, vitreous and retina as separated domains. Choroid and sclera are sometimes lumped together in a unique domain, as in . Anterior structures are described with different levels of precision, ranging from the "anatomically accurate" representation of  to the representation of mimicking boundary conditions of . Interface conditions between adjacent layers l and m are, in their most general form, expressed by the permeation relation
Dl being the diffusion coefficient in layer l with outward unit normal vector nl, Cl and Cm the drug concentrations in layers l and m, respectively, Llm the porous permeability of the membrane between layers l and m and Plm the partition coefficient keeping into account the affinity of the drug in the two layers. As a general remark, it is worth noticing that, while in  it is emphasized that fine geometric details are not needed to capture the salient features of drug distribution and only the proportions of different domains of the eye are considered to be important, more accurate analyses may need the inclusion of structures responsible of specific pathways (as for example the periocular space in  or the gap of Petit in .
Biophysical Phenomena Modeling
In table 2, we provide a comparative analysis between the biophysical phenomena included in the considered mathematical models. All the models include diffusion mechanisms for the drug and the majority also includes porous/fluid filtration, modeled for the most part by Darcy equations, but also by Navier-Stokes equations  and by Darcy-Brinkman equations . Simulation results show that the filtration field has almost a negligible effect for low-weight drugs [17,18] , whilst it becomes significant for high-weight drugs. On the contrary, the contribution of RPE pumping always appear important, especially when the retina is the therapeutic target. The RPE pumping effect, included in a very limited selection of papers, is modeled as a fictitious advective field or through different inward/outward permeabilities of the RPE membrane [13,19,20]. The inclusion of the significant clearance pathways for a certain drug administration mode is rather varied. Choriocapillaries are often treated as a perfect sink [16,21] or they are modeled via porous jump equations . Retinal blood vessels are only considered in  and, via a mixture equation approach, in . Generally, the fluid filtration is assumed to be steady state (except for)  and not influenced by the drug transport.
Numerical Approaches and Validation Techniques
In table 3, we provide a comparative analysis between the molecules considered in each model, the mathematical description of the drug source, the numerical approaches and parameter estimates used in the considered models for drug delivery to the PSE. For the most part, low-weight molecules are considered, with a predominance of the marker fluorescein on which the most part of experimental data are based. Only a few papers address different molecules with therapeutic significance, as the anti-VEGF macromolecular drug ranibizumab (IgG1 Fab) in . The drug source is represented for IV injection as a bolus of prescribed concentration of spherical shape [15,21] or time evolving shape  and for EP plug as a prescribed concentration  or prescribed time dependent concentration  at the scleral boundary. In Ferreira, et al. , the drug implant is inserted in the vitreous.
The inclusion in the model of the spatial dependency leads to large systems of partial differential equations, whose solution with analytical techniques is unfeasible, if not impossible. Thus, the use of numerical methods is in order. Spatial discretization can be based on finite elements (FEMs) [13,19] or finite volume techniques , generally implemented in commercial packages or on finite differences . Different discretization schemes are adopted for time marching, including implicit/semi-implicit finite difference schemes [13,14] and adaptive Runge-Kutta schemes . Three-dimensional problems give rise to large size systems, with a number of degrees of freedom ranging from about 200,000 nodes  to 321,000 nodes  and may require the use of parallel computing (as in the simulations carried out on a quad processor in . Special care must be taken when simulating diffusion/convection/reaction equations for drug delivery. Convection dominates when considering high weight/scarcely diffusing molecules or modeling RPE pumping actions as a fictitious velocity field. Reaction, representing metabolic terms/degradation effects, may also become dominant and induce instabilities and/or unphysical negative concentration fields . Best fit between the simulated drug levels and experimental/in vivo data have been performed to assess the model validity and fit critical parameters (using for the most part fluorescein data as those from Palestine & Brubaker)  in Balachandran & Barocas  or those from  in Park, et al. ; Gd-DTPA data were used in . Parameters estimated from data fitting of fluorescein have been often used also for different drugs/delivery routes. Diffusivity, permeability and RPE pumping parameters have been computed from data fitting [13,16,18], from analysis of TEM images  or from analytical relations/correlations between different layers . It should be observed that, whilst the actual value of some of these parameters will not probably affect too much the results when dealing with low weight molecules, when considering larger/heavier molecules the parameter values can have a more significant impact .
In table 4, we present the main results obtained in the considered mathematical models. Simulations have been performed to study concentration profiles in different conditions. Sensitivity studies have been carried out to address the significance of: i) certain biophysical mechanisms (active transport in the retina) [13,20], clearance rates ; ii) technical details in the drug administration (needle angle) , implant position ; different boundary conditions (vitreous outflow) ; implications of interspecies variability ; pathological conditions (inflammation of the retina), posterior vitreous detachment .
With the purpose of comparison of the different models when computing a similar quantity, in figure1 we show drug peak (and mean) concentration value and time-to-peak in the retina and choroid. We represent the data where available. The results are all obtained for low-weight molecules comparable to fluorescein. The drug is delivered via EP plug or EP injection. It is evident from the graphs a significant scattering of the results, of which several factors may be responsible. What appears more evident from our review work is that:
1. There is no general agreement on the detail level required for the geometrical description of the anatomical structures. Different choices are carried out, depending on the considered drug delivery route (IV injection, episcleral implant and even with the same approach with the precise location of the drug source), species (man, rabbit, monkey, rat) and on the intraspecies variability
2. The characteristics of the considered drug (hydrophilicity/lipophilicity, molecular weight and size, charge) can impact the choice of including or not in the model certain bio-physical phenomena (for example, the clearance routes and the permeability across membranes as the RPE) for a specific drug. Moreover they can also influence, by themselves, the geometrical detail needed in the description of a certain structure
3. The value of the parameters strongly influence the model results. This is clearly demonstrated by the results obtained by the sensitivity analyses performed in several papers. The parameter values obtained from experiments or from data fitting on the model are often very different and thus do not represent valid benchmarks
4. In general, experimental data are very sparse and significantly differ from each other. This lack of "golden standard experimental data" makes more difficult to dispose of well assessed data for model validation purposes
5. The mathematical representation itself of the phenomena (for example, representing RPE pumping effects by a fictitious convective field or by an altered permeability) impacts on the results. However, again, there is no "gold standard" for carrying out a comparison.
Conclusions and Perspectives
According to the World Health Organization, in 2010, more than 13% of the cases of blindness were the result of diseases of the PSE . It is thus of high importance to dispose of efficient and well tolerated ways to deliver drugs to this region. Mathematical models can provide tools useful to enlarge the knowledge on this topic, supporting and completing clinical experiments. In general, much space is still open when it comes to address therapeutically relevant drugs, especially macromolecules (anti-angiogenic antibodies, oligonucleotides, growth factors, and trans-gene expression)  and prodrugs improving drug bio-availability) . For these drugs, parameter estimation may need reviewing along with reconsideration of the significant biophysical mechanisms. To cite one example, the convection field can become a determinant of the drug fate, and contrast/enhance RPE pumping actions. A better understanding of the impact of aging/pathological conditions is necessary to develop more effective therapies so that drug concentration is maintained within the therapeutic window at the target site for the desired period of time (so far the effect of the liquefaction of the vitreous, breakdown of blood retinal barrier and elevated intraocular pressure have been considered in Ferreira, et al.) . This aspect should also include the effect of surgical procedures, for example the increased drug clearance due to vitrectomy  or the modified fluid-dynamics due to the use of vitreous substitutes .
- Macwan JS, Hirani A, Pathak Y (2016) Challenges in ocular pharmacokinetics and drug delivery. In: Nano-Biomaterials For Ophthalmic Drug Delivery, Springer, 593-607.
- Gaudana R, Megha Barot, Ashaben Patel, et al. (2013) Barriers for posterior segment ocular drug delivery. Treatise on Ocular Drug Delivery 1: 68-95.
- El Sanharawi M, Kowalczuk L, Touchard E, et al. (2010) Protein delivery for retinal diseases: from basic considerations to clinical applications. Prog Retin Eye Res 29: 443-465.
- Ranta V, Urtti A (2006) Transscleral drug delivery to the posterior eye: prospects of pharmacokinetic modeling. Adv Drug Deliv Rev 58: 1164-1181.
- Missel P (2012a) Computer modeling for ocular drug delivery. In: Ocular Drug Delivery Systems: Barriers and Application of Nanoparticulate Systems. CRC Press, 59-92.
- Janoria KG, Gunda S, Boddu SH, et al. (2007) Novel approaches to retinal drug delivery. Expert Opin Drug Deliv 4: 371-388.
- Gaudana R, Ananthula HK, Parenky A (2010) Ocular drug delivery. AAPS J 12: 348-360.
- Jaffe G, Ashton P (2006) Intraocular Drug Delivery. In: Pearson A, Taylor & Francis Group, New York.
- Hughes PM, Olejnik O, Chang-Lin JE, et al. (2005) Topical and systemic drug delivery to the posterior segments. Adv Drug Deliver Rev 57: 2010-2032.
- Stay M, J X, Randolph TW, et al. (2003) Computer simulation of convective and diffusive transport of controlled-release drugs in the vitreous humor. Pharm Res 20: 96-102.
- Tojo K (2004) A Pharmacokinetic Model for Ocular Drug Delivery. Chem Pharm Bull (Tokyo) 52: 1290-1294.
- Jooybar E, Abdekhodaie MJ, Farhadi F, et al. (2014) Computational modeling of drug distribution in the posterior segment of the eye: effects of device variables and positions. Math Biosci 255: 11-20.
- Balachandran R, BarocasV (2008) Computer modeling of drug delivery to the posterior eye: effect of active transport and loss to choroidal blood flow. Pharmaceut Res 25: 2685-2696.
- Narasimhan A, Vishnampet R (2012) Effect of choroidal blood flow on transscleral retinal drug delivery using a porous medium model. International Journal of Heat and Mass Transfer 55: 5665-5672.
- Missel P (2012b) Simulating intravitreal injections in anatomically accurate models for rabbit, monkey, and human eyes. Pharmaceut Res 29: 3251-3272.
- Kim H, Lizak MJ, Tansey G, et al. (2005) Study of ocular transport of drugs released from an intravitreal implant using magnetic resonance imaging. Ann Biomed Eng 33: 150-164.
- Kavousanakis M, Kalogeropoulos N, Dimitrios T Hatziavramidis (2014) Computational modeling of drug delivery to the posterior eye. Chemical Engineering Science 108: 203-212.
- Mac Gabhann F, Demetriades AM, Deering T, et al. (2007) Protein transport to choroid and retina following periocular injection: theoretical and experimental study. Ann Biomed Eng 35: 615-630.
- Kotha S, Murtomaki L (2014) Virtual pharmacokinetic model of human eye. Math Biosci 253: 11-18.
- Causin P, Malgaroli F (2016) Mathematical assessment of drug build-up in the posterior eye following transscleral delivery. Journal of Mathematics in Industry 6: 9.
- Park J, Bungay PM, Lutz RJ, et al. (2005) Evaluation of coupled convective-diffusive transport of drugs administered by intravitreal injection and controlled release implant. J Control Release 105: 279-295.
- Ferreira JA, de Oliveira P, da Silva PM, et al. (2017) Mathematics of aging: Diseases of the posterior segment of the eye. Computers & Mathematics with Applications 73: 11-26.
- Palestine AG, Brubaker RF (1981) Pharmacokinetics of fluorescein in the vitreous. Invest Ophthalmol Vis Sci 21: 542-549.
- Araie M, Maurice D (1991) The loss of fluorescein, fluorescein glucuronide and fluorescein isothiocyanate dextran from the vitreous by the anterior and retinal pathways. Exp Eye Res 52: 27-39.
- Pascolini D, Mariotti S (2012) Global estimates of visual impairment: 2010. Br J Ophthalmol 96: 614-618.
- Isakova K, Pralits O, Repetto R, et al. (2014) A model for the linear stability of the interface between aqueous humor and vitreous substitutes after vitreoretinal surgery. Physics of Fluids 26: 124101.
- Ninawe P, Hatziavramidis D, Parulekar SJ (2010) Delivery of drug macromolecules from thermally responsive gel implants to the posterior eye. Chemical Engineering Science 65: 5170-5177.
Paola Causin, Department of Mathematics, Università degli Studi di Milano via Saldini 50, 20133 Milano, Italy.
© 2017 Causin P, et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.