Mass transfer in the pathogenesis and treatment of eosinophilic esophagitis

Eosinophilic esophagitis (EoE) is an inflammatory disease of the esophagus, characterized by penetration of eosinophils into the esophageal wall. Once thought to be a rare condition, incidents of EoE are becoming more common. As EoE has risen to prominence, significant effort has been undertaken to successfully diagnose and treat the disease. Most studies of EoE have utilized a typical approach of biological assessment. Although this approach is successful, it may neglect fundamental biophysics that impact the progression and treatment of the disease. The use of engineering methods to create a physics-based assessment of EoE defines a biophysical model to answer fundamental questions concerning EoE development, which can be used in the treatment of the disease. The methods and techniques are common in traditional engineering applications and are used to obtain insight into a physical system through mathematical application. This work examines a purported relationship between eosinophil concentration and esophageal morphology. The pathogenesis of EoE using conservation principles is performed. Finally, a model for delivery of hydrophobic drugs from polymeric micelles is presented, which vessels can be utilized for delivery of drugs into the epithelium in the esophagus.

MASS TRANSFER IN THE PATHOGENESIS AND TREATMENT OF EOSINOPHILIC ESOPHAGITIS by David Dean Shaw A dissertation submitted to the faculty of The University of Utah in partial fulfillment of the requirements for the degree of Doctor of Philosophy Department of Chemical Engineering The University of Utah August 2017 Copyright © David Dean Shaw 2017 All Rights Reserved The University of Utah Graduate School STATEMENT OF DISSERTATION APPROVAL The dissertation of David Dean Shaw has been approved by the following supervisory committee members: , Chair John McLennan 28 APR 2017 Date Approved Leonard F. Pease III , Member 28 APR 2017 Date Approved Swomitra Mohanty , Member 28 APR 2017 Date Approved , Member Gerald J. Gleich 28 APR 2017 Date Approved Kathryn A. Peterson , Member 28 APR 2017 Date Approved , Member Richard L. Raun 28 APR 2017 Date Approved and by the Department/College/School of , Chair/Dean of Milind Deo Chemical Engineering and by David B. Kieda, Dean of The Graduate School. ABSTRACT Eosinophilic esophagitis (EoE) is an inflammatory disease of the esophagus, characterized by penetration of eosinophils into the esophageal wall. Once thought to be a rare condition, incidents of EoE are becoming more common. As EoE has risen to prominence, significant effort has been undertaken to successfully diagnose and treat the disease. Most studies of EoE have utilized a typical approach of biological assessment. Although this approach is successful, it may neglect fundamental biophysics that impact the progression and treatment of the disease. The use of engineering methods to create a physics-based assessment of EoE defines a biophysical model to answer fundamental questions concerning EoE development, which can be used in the treatment of the disease. The methods and techniques are common in traditional engineering applications and are used to obtain insight into a physical system through mathematical application. This work examines a purported relationship between eosinophil concentration and esophageal morphology. The pathogenesis of EoE using conservation principles is performed. Finally, a model for delivery of hydrophobic drugs from polymeric micelles is presented, which vessels can be utilized for delivery of drugs into the epithelium in the esophagus. For Emily, Eleanor, and Virginia TABLE OF CONTENTS ABSTRACT....................................................................................................................... iii  ACKNOWLEDGEMENTS .............................................................................................. vii Chapters 1. INTRODUCTION ........................................................................................................ 1  1.1 References...................................................................................................... 8  2. MASS TRANSFER OF ANTIGEN FLOWING DOWN AN IDEALIZED ESOPHAGUS ............................................................................................................. 10  2.1 Introduction.................................................................................................. 10  2.2 Methods ....................................................................................................... 11  2.3 Results and Analysis .................................................................................... 17  2.4 Conclusion ................................................................................................... 19  2.5 References.................................................................................................... 34  3. A LINEAR MODEL OF EOSINOPHIL GENERATION IN THE ESOPHAGUS ... 35 3.1 Abstract ........................................................................................................ 35  3.2 Introduction and Background ...................................................................... 35  3.3 Methods ....................................................................................................... 36  3.4 Results and Analysis .................................................................................... 43  3.5 Conclusion ................................................................................................... 46  3.6 References.................................................................................................... 56  4. DIFFERENTIAL MODEL OF EOSINOPHIL GENERATION IN THE ESOPHAGUS ............................................................................................................. 58  4.1 Abstract ........................................................................................................ 58  4.2 Introduction and Background ...................................................................... 58  4.3 Methods ....................................................................................................... 59  4.3.1 Esophagus Model Definition ........................................................... 60  4.3.2 Linear Stability Analysis.................................................................. 62  4.3.3 Nondeterministic Computational Modeling .................................... 68  4.3.4 Plackett-Burman Sensitivity Analysis ............................................. 71  4.4 Results and Analysis .................................................................................... 74  4.5 Conclusion ................................................................................................... 76  4.6 References.................................................................................................... 94  5. RELEASE OF PHARMACEUTICAL COCKTAILS FROM SMALL POLYMER MICELLES CONTROLLED BY THERMODYNAMICS........................................ 96  5.1 Abstract ........................................................................................................ 96  5.2 Introduction and Background ...................................................................... 97  5.3 Methods ....................................................................................................... 99  5.4 Results and Analysis .................................................................................. 105  5.5 Conclusions................................................................................................ 110  5.6 References.................................................................................................. 120  6. CONCLUSIONS AND FUTURE WORK ............................................................... 123  vi vii ACKNOWLEDGEMENTS This work has be supported by financial support from a University of Utah TCG grant, TCIP grants from the Utah Governor's Office of Economic Development, the American Partnership for Eosinophilic Disorders (APFED), and NSF CBET (CBET-1125490). I am grateful to all of those with whom I have had the pleasure to work during this and other related projects. I am truly grateful to the members of my committee, who have provided me extensive personal and professional guidance and taught me a great deal about both scientific research and life in general. I would specifically like to thank Dr. John McLennan, the chairman of my committee. As the individual who facilitated the completion of this journey, I will be eternally grateful. I want to also thank Dr. Richard L. Raun, my friend and mentor for so many years. He has shown me, by his example, what a good scientist (and person) should be. I would also be neglectful if I did not thank Dr. Leonard F. Pease III. Without his willingness to introduce me to the topics involved in this work and his enthusiasm I would not be in the position I am in today. In this journey nobody has been more important to me than the members of my family. I would like to thank my parents, whose love and guidance gave me the rudder to guide me through to calm waters. Most importantly, I wish to thank my loving and very supportive wife, Emily, and my two wonderful children, Ellie and Ginny, who provide unending inspiration, love, and comedy at every turn. 1 CHAPTER 1 INTRODUCTION Eosinophil granulocytes, or eosinophils, are a type of white blood cells characterized by coarse granules within their cytoplasm. Eosinophils are a key element in combating multicellular parasites and are also key participants in allergic responses.1 The esophagus, unlike other segments of the gastrointestinal tract, is normally lacking eosinophils.2 Thus, a discovery of esophageal eosinophilia signifies an uncharacteristic condition. One such cause of eosinophilia is eosinophilic esophagitis (EoE), an inflammatory disease affecting the esophagus. Once considered a rare condition, studies demonstrate that eosinophilic esophagitis, like other immunologic disorders, is becoming more prominent.3-5 The disease has been reported in every populated continent except Africa6 and incidents have risen over several decades - affecting more than 300,000 patients per year.7 This increased occurrence is likely attributed to both increased occurrence and growing recognition of the disease. Eosinophilic esophagitis commonly affects children and adults, with a predilection for young males.5 Symptoms are characterized by heartburn, failure to thrive, epigastic or chest burn, dysphagia (recurrent in adults), and food impaction (recurrent in adults). In many cases, emergency intervention is required during incidents of food impaction.6,8 In addition, the disease affects the epithelial layer and the endoscopically observed 2 morphology of the esophagus, resulting in plaques or exudates, longitudinal segmentation (i.e., periodic rings), and furrows (see Figure 1.1) along the esophagus.9-12 Based on a patient's medical history for dysphagia, an esophagogastroduodenoscopy (EGD) is collected to obtain esophageal biopsies to diagnose EoE. Biopsies are collected at several locations throughout the esophagus to overcome the "spotty manner" of the disease.13 The resultant biopsies are examined for eosinophil concentration (see Figure 1.2). Biopsy examinations resulting in 15-20 eosinophils per high power field (1 hpf ≈ 0.22 mm2) satisfies diagnostic recommendations.14,15 However, elevated eosinophil concentrations in the esophagus are not specific to a particular disease. EoE requires the elimination of other causes of esophageal eosinophilia, specifically gastroesophageal reflux disease (GERD), which may display eosinophilia albeit at much lower eosinophil densities. EoE is typically distinguished from GERD by persistent esophageal eosinophilia, despite acid neutralization therapy.6 Eosinophil granulocytes are categorized by the presence of small particles, or granules, in their cytoplasm (see Figure 2.3). Each eosinophil granule comprises a core and matrix. These granules contain several cationic proteins within the cellular matrix. The granule matrix contains eosinophil peroxidase (EPO), eosinophil cationic protein (ECP), and eosinophil-derived neurotoxin (EDN), while a crystalloid core contains major basic proteins MBP-1 and MBP-2.6,16,17 Each of these proteins are cytotoxic to most tissues, making eosinophils effective combatants against biological parasites.17 Eosinophil infiltration into the esophagus is characterized by the release of granule proteins from activated eosinophils to both the epithelial and sub-mucosal layers.18,19 3 Activation is induced by an immune stimulus and once activated, the eosinophils experience degranulation, releasing granule proteins into the surrounding tissue.20-22 The immune stimulus has been shown to be an allergic response, specifically from ingested foods. Modifying a patent's diet has shown EoE symptoms can be managed and even resolved and EoE-associated inflammation allowed to heal.23 Eliminating foods containing milk, egg, soy, wheat, nuts, and fish - foods commonly associated with allergic responses - have been effective at reducing EoE symptoms.24 Alternatively, or in addition to a modified diet, if EoE symptoms do not subside, common drug therapies such as corticosteroids and glucocorticoids may reduce EoE-induced inflammation. With the growth of diagnostic incidences, a clinical description of EoE is well established.6 Yet a biophysical description of EoE from a systems biology perspective has not been significantly examined. Biophysics is the application of the physical sciences to a biological system. EoE research has primarily been concerned with the interaction between the various systems of a cell and how these interactions are regulated or treated. The application of a biophysical study is vital for understanding the origins of EoE, which in turn provides avenues for treatments and for assessing improved efficacy of those proposed treatments. The overall objective of this research study is to use engineering techniques to answer fundamental questions concerning EoE development and treatment. The methods and techniques are common in traditional engineering applications (e.g., fluid dynamics, heat transfer, etc.) and are used to obtain insight into a biophysical system through mathematical application. Chapter 2 of this work examines if a purported relationship exists between the concentration of eosinophils and esophageal morphology through the 4 transport of antigens into the epithelium. Assessment of the pathogenesis of EoE utilizing conservation equations is performed in Chapters 3 and 4. Finally, a model for delivery of hydrophobic drugs from polymeric micelles is presented in Chapter 5. These micelles are proven vessels for delivery of hydrophobic drugs for cancer treatment and may be used in the delivery of hydrophobic drugs (e.g., a topically active corticosteroid) for the treatment of EoE in the esophagus. 5 Figure 1.1 - Endoscopy images of (a) normal esophagus, (b) esophagus with furrows, (c) esophagus with rings, and (d) esophagus with plaques or exudate (images courtesy of Fox, et al., Gastrointestinal Endoscopy, 2002; 56:260-270. and Remedios, et al., Drugs March 26 2011, 71(5), 527-40). 6 Figure 1.2 - Comparison between a diseased tissue specimen with eosinophil infiltration and a normal esophageal tissue specimen. Pathology slides courtesy of Hedieh Saffari, University of Utah, 2013. 7 Figure 1.3 - Intact eosinophil prior to degranulation showing the bilobed nucleus with granule proteins in the granules and an intact cytoplasm. Image courtesy of Hedieh Saffari, University of Utah, 2013. 8 1.1 References (1) Barbara, Y. Wheater's functional histology : a text and colour atlas; 4th ed.; Churchill Livingstone: Edinburgh ; New York, 2000. (2) Mishra, A.; Hogan, S. P.; Lee, J. J.; Foster, P. S.; Rothenberg, M. E. Journal of Clinical Investigation 1999, 103, 1719. (3) Noel, R. J.; Putnam, P. E.; Rothenberg, M. E. New England Journal of Medicine 2004, 351, 940. (4) Fox, V. L.; Nurko, S.; Furuta, G. T. Gastrointestinal Endoscopy 2002, 56, (5) Gonsalves, N.; Kahrilas, P. J. Neurogastroenterology & Motility 2009, 21, (6) Rothenberg, M. E. Gastroenterology 2009, 137, 1238. (7) Pease, L. F.; Gleich, G. J.; Hedieh, S., Personal Communication. 260. 1017. (8) Mackenzie, S. H.; Go, M.; Chadwick, B.; Thomas, K.; Fang, J.; Kuwada, S.; Lamphier, S.; Hilden, K.; Peterson, K. Alimentary Pharmacology & Therapeutics 2008, 28, 1140. (9) Remedios, M.; Jones, D.; Kerlin, P. Drugs 2011, 71, 527. (10) Gonsalves, N.; Policarpio-Nicolas, M.; Zhang, Q.; Rao, M. S.; Hirano, I. Gastrointestinal Endoscopy 2006, 64, 313. (11) Croese, J.; Fairley, S. K.; Masson, J. W.; Chong, A. K. H.; Whitaker, D. A.; Kanowski, P. A.; Walker, N. I. Gastrointestinal Endoscopy 2003, 58, 516. (12) Potter, J. W.; Saeian, K.; Staff, D.; Massey, B. T.; Komorowski, R. A.; Shaker, R.; Hogan, W. J. Gastrointestinal Endoscopy 2004, 59, 355. (13) Saffari, H.; Peterson, K. A.; Fang, J. C.; Teman, C.; Gleich, G. J.; Pease Iii, L. F. Journal of Allergy and Clinical Immunology 2012, 130, 798. (14) Furuta, G. T.; Liacouras, C. A.; Collins, M. H.; Gupta, S. K.; Justinich, C.; Putnam, P. E.; Bonis, P.; Hassall, E.; Straumann, A.; Rothenberg, M. E. Gastroenterology 2007, 133, 1342. (15) Rodrigo, S.; Abboud, G.; Oh, D.; DeMeester, S. R.; Hagen, J.; Lipham, J.; DeMeester, T. R.; Chandrasoma, P. American Journal Gastroenteroology 2008, 103, 435. 9 (16) 54, 656. Peters, M. S.; Rodriguez, M.; Gleich, G. J. Laboratory Investigation 1986, (17) Rothenberg, M. E.; Hogan, S. P. Annual Review of Immunology 2006, 24, 147. (18) Kephart, G. M.; Alexander, J. A.; Arora, A. S.; Romero, Y.; Smyrk, T. C.; Talley, N. J.; Kita, H. American Journal of Gastroenterology 2009, 105, 298. (19) Justinich, C. J.; Ricci, A., Jr.; Kalafus, D. A.; Treem, W. R.; Hyams, J. S.; Kreutzer, D. L. Journal of Pediatric Gastroenterology 1997, 25, 194. (20) Gleich, G. J.; Adolphson, C. R. Advanced Immunology 1986, 39, 177. (21) Erjefält, J. S.; Andersson, M.; Greiff, L.; Korsgren, M.; Gizycki, M.; Jeffery, P. K.; Persson, C. G. A. Journal of Allergy and Clinical Immunology 1998, 102, 286. (22) Melo, R. C. N.; Spencer, L. A.; Perez, S. A. C.; Neves, J. S.; Bafford, S. P.; Morgan, E. S.; Dvorak, A. M.; Weller, P. F. Lab Invest 2009, 89, 769. (23) Markowitz, J. E.; Spergel, J. M.; Ruchelli, E.; Liacouras, C. A. American Journal of Gastroenterology 2003, 98, 777. (24) Gonsalves, N.; Yang, G. Y.; Doerfler, B.; Ritz, S.; Ditto, A. M.; Hirano, I. Gastroenterology 2012, 142, 1451. 10 CHAPTER 2 MASS TRANSFER OF ANTIGEN FLOWING DOWN AN IDEALIZED ESOPHAGUS In a study completed by Saffari, et al.,1,2 eosinophil density in esophagectomy images appeared to be a function of the curvature of the esophagus. To assess the physical cause, analysis using computational modeling is performed. Two idealized geometries with varying degrees of esophageal undulations are examined. Solutions to the models are calculated using the commercial program ANSYS® Fluent®. We hypothesize that the observations from Saffari, et al. are correct and eosinophil density does vary as a function of curvature within the esophagus in accord with Saffari's observations. The result of this analysis will link a biophysical model of esophageal morphology to physically measured results. 2.1 Introduction Within the esophagus, the distribution of eosinophils in eosinophilic esophagitis (EoE) has been shown to be patchy. In a study completed by Saffari, et al.,1,3 eosinophil density was mapped around the circumferential location of the esophagus from esophagectomy images. Those results suggest a general trend that eosinophil density varies as curvature of the esophagus varies. A weak relationship can be seen in Figure 11 2.1, where the curvature from one of the esophageal biopsies is plotted next to its corresponding eosinophil density count as a function of the perimeter location. Although the relationship may be weak (in part because of difficulty aligning estimates of curvature with Saffari's data), it is believed that luminal flow in the esophagus is a primary driver of antigen mass transfer in the esophagus. It is expected that luminal flow is influenced by the shape of the esophagus. Thus, a study of the conditions controlling mass transfer into the epithelium is required to connect the morphological condition of the esophagus and physical data. 2.2 Methods To determine if curvature has an impact on mass transfer, a fluid model is used to calculate antigen transport into the esophagus. The concentration of granule proteins in the epithelium is expected to be directly proportional to the concentration of antigen in the epithelium. To perform a curvature analysis on the esophagus, a mode of antigen transport into the esophagus needs to be recognized. One common mode of transport is a continual flow of a luminal film along the epithelium of the esophagus as opposed to a swallowed or bolus flow through the esophagus. The thickness of the luminal film has not been directly measured, thus an estimate is made to broadly cover the range of plausible film thicknesses. It is known the esophagus lacks an adherent mucus layer, but is constantly bathed by swallowed saliva and secretions from submucosal glands.4,5 Measurements of the gastrointestinal mucus layer 12 suggest a thickness of approximately 70 µm in the corpus of the gastrointestinal system.4 Any estimate on the lumen film thickness will likely be within this order of magnitude. A thin film model is presented in Figure 2.2. This model assumes a constant flow of saliva and secretion flow down the esophagus. Additionally, this film is present in quiescent flow (i.e., between boli) and comprised primarily of saliva. It is assumed the esophageal wall acts as a no slip surface and the fluid remains Newtonian at modest flow rates and concentrations. The film thickness is assumed to be small compared to the perimeter and length of the esophagus; this assumption simplifies the model such that only r-directional velocities are pertinent. Applying the simplifications above, the general momentum conservation equation for the system in Figure 2.2 becomes, 1 ∆ ∆ where, , , , , , and , (2.1) are the velocity change in the radial-direction, location along the esophagus, viscosity of the fluid, the pressure along the esophagus, density, and the acceleration due to gravity. Thus, the only forces acting on the fluid are gravity and a differential pressure from the proximal to the distal ends of the esophagus. Using the appropriate boundary conditions on the system, eq (2.1) is solved for the velocity, 3 1 2 where, , (2.2) 13 ∆ ∆ 3 (2.3) . Eq (2.3) is the average velocity in the luminal film and is governed primarily by the pressure and gravity forces. To solve for the thickness of the film, it is assumed the variation along the curvature is sufficiently small, thus the average velocity is related to the volumetric flow rate by, (2.4) , where , , and H are the volumetric flow rate, perimeter of the esophagus, and thickness of the lumenal film. The lumenal film thickness is then determined by combining eqs (2.3) and (2.4) and solving for the film thickness, . ∆ ∆ (2.5) To simplify eq. (2.5), the ratio of the Bond and Capillary numbers is used, Bo Ca ̅ (2.6) . Substituting the average velocity in eq (2.3) and simplifying yields Bo Ca 3 ∆ ∆ . (2.7) The Bond and Capillary ratio gives the ratio of the body forces to viscous forces in a system. Applying eq (2.7) to eq (2.5) results in the equation for the thickness of the film, Bo Ca . (2.8) 14 This equation gives the thickness of the film as a function of the Bo/Ca ratio, which represents the ratio of body forces to viscous forces acting on the fluid flow in the esophagus. In this formulation, the numeric value of the Bo/Ca ratio is dependent on what forces are applied to the system. The numerical value of the ratio would be large as the viscous forces become large and small for large body forces (e.g., large pressure differential). The thickness of the fluid layer flowing within the esophagus, given in eq (2.8), is primarily a function of two possible scenarios. One potential scenario utilizes gravitational force on the liquid without a pressure gradient, for instance when no pressures are exerted through swallowing. The other is the application of gravitational force and a differential pressure along the esophagus (e.g., flow under the influence of a bolus). Both scenarios are used to constrain the thickness of the luminal film to a maximum and minimum thickness. As the film thickness increases, the model becomes a multidimensional system. Curvature and length of the esophagus must be accounted for. An analytical solution is not easily obtained to successfully model the fluid flow and mass transport down the esophagus unless the film thickness is small. Since the film thickness is not known experimentally, a parametric analysis of the definition given in eq (2.8) is used to define a plausible range of luminal fluid thickness for the model. A customized computer code could be written, but the commercial Computational Fluid Dynamics (CFD) program ANSYS® Fluent® is selected due to the program's flexibility and ease of use for the purposes and intent of this analysis. 15 The start of all computational methods is the generation of an applicable geometry that can tie the model to a physically realistic system. The cross-sectional geometry in Figures 2.3a and 2.3b was prepared as an idealization of an esophageal shape. These shapes were selected based on geometrical simplicity, focusing on curvature differences, while simultaneously approximately representing some observed esophagectomy images (see Figure 2.3c). Recognizing the symmetry present in the idealized geometry in Figures 2.3a and 2.4b, only a fraction of the geometry must be modeled. Figure 2.4 shows the simplification applied to the geometry. Figure 2.5 shows a view of the mesh for one of the models used in this assessment along the negative z-axis. Throughout this analysis, mesh skewness was kept below 0.85 to avoid divergence in the solution of the model. The meshed geometry was imported into Fluent® and boundary conditions set to simulate the expected mass transfer configuration (see Figure 2.6). The epithelium was defined as a porous media with porosity and permeability. In Fluent®, a permeable membrane is defined by the porosity and an associated frictional factor. The flow of material in the permeable membrane (i.e., epithelium) is simplified as a fluid flow through a packed bed. Assuming the fluid film has laminar flow with a velocity sufficiently low for Re<10, the Blake-Kozeny or Kozeny-Carman equation6 can be used to define the bed friction factor as, 75 1 ̅ . (2.9) The friction factor from eq (2.9) is a function of the porosity, , particle diameter, , average velocity, ̅ , and liquid's properties. The porosity of a healthy epithelium has 16 been measured to be approximately 0.43.7 The diameter of epithelial cells are assumed spherical, and measured to be around 0.4 µm.8 The liquid in the film has the same general properties as water, but also contains a casein mole fraction of 0.03 to simulate the amount of casein present in milk.9 Casein is a family name for related phosphorproteins commonly found in mammalian milk.10 Casein was chosen due to its history of inciting allergic responses.10 The epithelium transport of casein was approximated using Stokes-Einstein diffusion coefficients and experimental approximations provided by Bergman, et al.11 There was no species diffusive flux along all epithelial boundaries. The flux was set to zero due to Fluent®'s inability to simply specify the mass transfer flux along a boundary without producing potentially erroneous results of negative concentrations. Three main variables are assessed to determine the effect of the esophagus curvature on casein concentration within the epithelium. The most significant variable is the curvature of the idealized esophagus, shown in Figures 2.3a and 2.3b. Fluid film thickness is relatively unknown and is also varied based on an assessment generated from eq (2.8). Finally, the other potential major impact parameter is the porosity of the epithelium, which will impact the concentration profile of casein in the epithelium. A summary of values and boundary conditions is given in Table 2.1. The dependence of curvature of the esophagus is evaluated by using minimum and maximum film thicknesses from eq (2.8) and by varying the degree of curvature in the idealized esophagus (shown in Figures 2.3a and 2.3b). A comparison between the concentration of casein in the peaks and valleys of the epithelium in each scenario is 17 performed to characterize the effect of curvature on casein concentration in the epithelium. 2.3 Results and Analysis The thickness of the luminal film was determined through a parametric analysis of the Bo/Ca ratio defined in eq (2.8). A volumetric flow between 0.3 and 10 mL/min (from salivary rates12,13) was used with and without the effect of a pressure differential. Using this method, the film is thinnest, between 16 and 70 µm, when a pressure gradient of 80 Torr14 is applied to the esophagus and thickest, between 29 and 123 µm, when the pressure gradient is neglected, similar to the 70 µm from literature.4 Figure 2.7 shows the applicable range of fluid film thicknesses. These maximum and minimum values are used in the development of the Fluent model. Using the simplified geometries of the esophagus and the parameters (see Table 2.2) in Section 2.2, the solutions are calculated under steady state conditions. The steady state model is used to produce a conservative result for the concentration profile within the epithelium. Figures 2.8 through 2.11 show the regions of significant change within the model for each condition. Figure 2.8 shows the mass fraction of casein in an epithelium with small undulations down the length of the esophagus for a "thick" fluid layer of 80 µm. Figure 2.8a shows the result for low porosity and Figure 2.8b for high porosity. Likewise, Figure 2.9 shows the mass fraction of casein in an epithelium with small undulations down the length of the esophagus for a "thin" fluid layer of 16 µm. Figure 2.9a shows the result for low porosity and Figure 2.9b for high porosity. The differences in these results are quite 18 pronounced. The porosity appears to have little impact on the shape and form of the concentration profile, as opposed to the fluid film thickness (see Table 2.2). Figure 2.10 shows the mass fraction of casein in an epithelium with large undulations down the length of the esophagus for the thick fluid layer. Figure 2.10a shows the result for low porosity and Figure 2.10b for high porosity. Similarly, Figure 2.11 shows the mass fraction of casein in an epithelium with small undulations down the length of the esophagus for the thin fluid layer. The differences in these results are also quite pronounced and different from the small undulations seen in Figures 2.8 and 2.9. For all conditions, the porosity appears to have little impact on the mass transfer, with the fluid film thickness having a large impact on the mass fraction profile. Comparing the different curvatures for each of the conditions shows the curvature has a significant impact on the shape of the mass fraction profile in the epithelium. This can be seen in Figure 2.12, where a decrease in mass fraction at the extrema of the curvature is present. This dependence is most likely due to forces causing excessive shear flow in those regions due to the fluid-surface contact. Table 2.2 shows a comparison of integrated concentrations along the center-line of the epithelium at each peak and valley of the undulations. The difference in integrated mass fraction between the models with small undulations is no more than 7.4%, while the models with large undulations are at least double the low result. This indicates a strong dependence between mass fraction concentration profile and curvature. Interestingly, the majority of a mass fraction change is in the proximal region of the esophagus. It is fair to assume the concentration of antigen, like casein, is proportional to eosinophil concentration within the epithelium, but this phenomenon is not directly 19 measured during the analyses of esophagectomy slices.15 The result from the model is most likely due to the boundary conditions and assumptions used to treat the epithelium, specifically, the distal boundary condition being assigned an insulated boundary and the epithelium treatment as a porous media similar to a packed-bed. I would say that you likely want to make the case that the model is representative. Although the model is an idealized representation of physical reality, the overall result is the same, with mass transport showing significant differences under different curvatures. 2.4 Conclusion The computational analysis shown here is an assessment of a physics-based model calculated through the commercial computational package, Fluent®, to describe the concentration of casein within the epithelium as a function of length and radial diameter. The model provides an adequate representation of the physical system, yielding valuable information about the diseased state. The model significantly demonstrates a difference in concentration at the peaks and valleys of the esophagus. This result is observed by Saffari, et al.,15 and appears be a cause of excessive concentrations of eosinophils. This result shows that the curvature of the esophagus drives the accumulation of eosinophils in specific areas of the esophagus, accelerating a diseased state. To eliminate, or even reduce, the curvature of the esophagus should be explored to minimize the impact of curvature on a diseased esophagus. This could be achieved through the use of a stint-like object that would provide structural support to the esophagus. 20 Figure 2.1 - Plot of the eosinophil count collected by Saffari, et al., versus a measurement of curvature for a given circumferential location around the esophagus. Data courtesy of Hedieh Saffari, University of Utah, 2013 21 PO dv r  0 dz vr  0 Pf Figure 2.2 - Luminal film thickness model used to estimate the film thickness, H, for analysis. 22 Figure 2.3 - Cross-sections of the esophagus in an idealized geometry with (a) large undulations, (b) for small undulations for use in computational modeling, and (c) micrograph section of the esophagus; pathology slides courtesy of Hedieh Saffari, University of Utah, 2013. 23 Figure 2.4 - Application of symmetry from the idealized model to the CFD model. 24 Figure 2.5 - Mesh of the model. The mesh quality parsed the geometry into acceptable sections. The yellow section in this image represents the fluid layer, while the grey section represents the epithelial layer. 25 Figure 2.6 - The basic boundary conditions for each of the computational models in FLUENT®. In the model, there are two materials. The first is the fluid running down the esophagus. The second is the epithelium, shown here as a permeable membrane. The exterior boundaries of the model are set in FLUENT ® as insulated, meaning all of the conserved property's ( ) fluxes are zero (e.g., mass fraction of casein, etc.). The main driving force for mass transfer is the boundary conditions at the distal and proximal regions of the esophagus. 26 Figure 2.7 - Range of lumen film thickness as a function of the Bond to Capillary number ratio. To generate the curves shown, nominal values of 977 kg⁄m3 , 9 10 Pa∙s, and 9.81 m⁄s2 for density, viscosity, and acceleration due to gravity are used, with 5 10 m3 ⁄s and 0.025 m for minimum volumetric flow rate (Q) and esophageal perimeter and 5 10 esophageal perimeter. m3 ⁄s and 0.06 m for maximum volumetric flow rate (Q) and 27 Figure 2.8 - The mass fraction of casein in an epithelium as a function of the location within the epithelium with small undulations down the length of the esophagus for the thick fluid layer for (a) low porosity (0.215) and (b) high porosity (0.86). The geometries shown are identical, with the only variation being the porosity, thus a significant difference in mass fraction along the length of the epithelium is clearly seen. 28 Figure 2.9 - The mass fraction of casein as a function of the location within the epithelium with small undulations down the length of the esophagus for the thin fluid layer for (a) low porosity (0.215) and (b) high porosity (0.86). The geometries shown are identical, with the only variation being the porosity, thus a significant difference in mass fraction along the length of the epithelium is clearly seen. 29 Figure 2.10 - The mass fraction of casein as a function of the location within the epithelium with large undulations down the length of the esophagus for the thick fluid layer for (a) low porosity (0.215) and (b) high porosity (0.86). The geometries shown are identical, with the only variation being the porosity, thus a significant difference in mass fraction along the length of the epithelium is clearly seen. 30 Figure 2.11 - The mass fraction of casein as a function of the location within the epithelium with large undulations down the length of the esophagus for the thin fluid layer for (a) low porosity (0.215) and (b) high porosity (0.86). The geometries shown are identical, with the only variation being the porosity, thus a significant difference in mass fraction along the length of the epithelium is clearly seen. 31 Figure 2.12 - The mass fraction of casein in an epithelium with large undulations down the length of the esophagus for the thick fluid layer. The change in mass fraction profile is clearly seen, showing a dependence of curvature on mass transfer. 32 Table 2.1 - Summary of assumptions used in the model. Assumption/Parameter Fluid flow down the esophagus is laminar Uniform porosity throughout the epithelium The epithelium is treated as a permeable membrane Model Solver: Pressure Based Zero mass flux along external Diffusion Coefficient Casein Molecular Weight Casein Density Casein Inlet Mole Fraction Mesh Skewness Porosity Density Viscosity Acceleration due to gravity Volumetric Flow Rate Pressure Difference Particle Diameter Value N/A N/A N/A N/A N/A 1.56 10 m2 ⁄s 23 1033 0.03 <0.85 Low: 0.215 High: 0.86 977 kg⁄m3 10 Pa∙s 9.81 m⁄s2 Min: 5 10 m3 ⁄s Max: 5 10 m3 ⁄s 80 Torr 4 10 m 9 33 Table 2.2 - Comparison between the center-line integrated mass fractions of casein within the epithelium. r‐direction  Undulation  Thickness  Porosity  Type  High  Small  Thick  Low  High  Small  Thin  Low  High  Large  Thick  Low  High  Large  Thin  Low  Valley Integrated  Mass Fraction  Peak Integrated  Mass Fraction  Absolute %  Difference  3.51E‐03  3.35E‐03  3.86E‐04  3.76E‐04  1.91E‐03  1.90E‐03  1.27E‐04  1.26E‐04  3.68E‐03  3.51E‐03  4.15E‐04  4.06E‐04  2.43E‐03  2.41E‐03  1.11E‐04  1.10E‐04  4.62%  4.56%  6.99%  7.39%  21.40%  21.16%  14.41%  14.55%  34 2.5 References (1) Saffari, H.; Peterson, K. A.; Fang, J. C.; Teman, C.; Gleich, G. J.; Pease III, L. F. Journal of Allergy and Clinical Immunology 2012, 130, 798. (2) Saffari, H., University of Utah, 2013. (3) Gonsalves, N.; Yang, G. Y.; Doerfler, B.; Ritz, S.; Ditto, A. M.; Hirano, I. Gastroenterology 2012, 142, 1451. (4) Dixon, J.; Strugala, V.; Griffin, S. M.; Welfare, M. R.; Dettmar, P. W.; Allen, A.; Pearson, J. P. American Journal of Gastroenterology 2001, 96, 2575. (5) Orlando, R. C. GI Motility Online 2006. (6) Seader, J. D.; Henley, E. J. Separation Process Principles; New York : Wiley: New York, 2006. (7) Baber, K. Master, Universität Stuttgart, 2009. (8) Dalton, B. A.; McFarland, G. A.; Steele, J. G. Journal of Biomedical Materials Research 2001, 56, 83. (9) (10) 51, 412. Ribadeau-Dumas, B.; Grappin, R. Lait 1989, 69, 357. Docena, G. H.; Fernandez, R.; Chirdo, F. G.; Fossati, C. A. Allergy 1996, (11) Bergman, T. L.; Lavine, A. S.; Incropera, F. P.; DeWitt, D. P. Fundamentals of Heat and Mass Transfer; 7th ed.; Hoboken, NJ : John Wiley: Hoboken, NJ, 2011. (12) Moritsuka, M.; Kitasako, Y.; Burrow, M. F.; Ikeda, M.; Tagami, J.; Nomura, S. Journal of Dentistry 2006, 34, 716. (13) Navazesh, M.; Kumar, S. K. S. JADA 2008, 139, 35S. (14) Hansen, M.; Saffari, H.; Peterson, K. A.; Gleich, G. J.; Pease III, L. F. (Unsubmitted Draft) 2013. (15) Saffari, H., University of Utah, 2013. 35 CHAPTER 3 A LINEAR MODEL OF EOSINOPHIL GENERATION IN THE ESOPHAGUS 3.1 Abstract Temporal evolution of granule protein concentrations and their putative steady state from a biophysical perspective characterizes the pathogenesis of EoE. It is hypothesized that conditions leading to an EoE diseased esophagus can be predicted by analyzing a reaction-diffusion model based on an interacting protein network. This study examines the trends that govern the evolution of granule protein concentrations in the epithelial layer of the esophagus. By using a mass transfer version of the lumped capacitance method, a determination is obtained that granule protein continues to grow until a steady state is reached. The most significant factors governing granule protein evolution in the epithelium yield a series of conditions defining the diseased state. These factors/conditions may suggest treatment pathways for EoE patients. 3.2 Introduction and Background A mathematical (systems biology) model for eosinophilic esophagitis (EoE) is developed using fundamental transport principles. The pathogenesis of EoE has been analyzed extensively and is a relatively complex biological process. EoE begins with the introduction of an antigen within the esophagus, causing the activation of T-helper cells 36 (Th2), releasing the cytokine interleukin proteins IL-5 and IL-13.1 IL-13 promotes epithelial cell repair and causes the production of eotaxin-3 from the epithelial cells.2,3 Both IL-5 and eotaxin-3 promotes the activation of eosinophils, releasing major basic proteins (MBP-1 and MBP-2), eosinophil peroxidase (EPO), eosinophil cationic protein (ECP), and eosinophil-driven neurotoxin (EDN).3,4 When activated, eosinophils also produce transforming growth factor beta (TGF-β), which continues to activate additional epithelial cells and contributes to hyperplasia and fibrosis.5,6 Figure 3.1 is a graphical representation of the pathogenesis process as described. To generate a mathematical description of the pathogenesis process may be difficult and may result in a model that is too unwieldy to be solved or usable. A simplified biological process is adopted to generate a tractable model. A system of linear mass transport and reaction equations are derived using mass transfer first principles. The intent is to determine if the mathematical system described has any inherent instability, which can clarify the pertinent parameters for a diseased state. 3.3 Methods In assessing the pathogenesis process outlined in Figure 3.1 and determining the significant interactions affecting EoE, certain protein-cell interactions can be neglected or integrated into other reactive processes to make a tractable model. The emphasis of the model development is to create a systems biology description, which emphasizes the evolution of eosinophil granule protein generation. Currently, interactions only influencing eosinophil activation are pertinent. This approach is taken due to the significance of eosinophil granule protein production in defining the EoE diseased state. 37 The transport of the significant species' concentrations within the esophageal epithelium is assumed to be only a function of time and not position. This assumption implies that the concentration gradients in the epithelium are negligible. This is a common practice in heat transfer when the dimensionless quantities are of the appropriate magnitude. This approach is often called a lumped capacitance method and is considered a simplified transient conduction model.7 This simplified transient approach performs an integral balance of species around a region of interest (i.e., the epithelial layer). The flux of species crossing each boundary is examined, in addition to accumulation within the region to give the species within the region of interest throughout time. Using a simplified one-dimensional esophagus, shown in Figure 3.2, the lumped capacitance method is applied by performing an integral balance around the epithelial layer. The balance for the antigen concentration yields the transport equation , where , , , , , and (3.1) are the concentration of the antigen in the epithelium, the molar flux into or out of the control volume, the surface area of the epithelium, the reaction coefficient of the antigen in the epithelium, and the volume of the epithelium, respectively. First order reactions are utilized throughout the derivation based upon observations from Nobuhisa, et al.8 and Humbles, et al.9 The molar flux into the system is defined using the mass transfer coefficient formulation,7,10 which becomes ≡ where , , , and , (3.2) are the mass transfer coefficient in the luminal film, concentration of the antigen at the surface of the film, and the concentration of the antigen at the luminal film-epithelium interface, respectively. 38 Through the same approach, the flux out of the control volume of the epithelium layer is defined as ≡ Similarly, , , , and . , (3.3) are the mass transfer coefficient in the submucosa, concentration of the antigen at the submucosa-muscular interface, and the concentration of the antigen at the submucosa-epithelium interface, respectively. The concentration at both sides of the epithelium can be related to the concentration of the antigen in the epithelium through the partition coefficient,11 also known as the "Kvalue" or "K-factor", and relates the concentration of a component in one phase to another phase across an interface. In mathematical terms, this is written as ≡ where, and (3.4) are the mole fraction in phase one and the mole fraction in phase two. Applying the definition in eq (3.4) to the interfaces on both sides of the epithelium gives , (3.5) and , (3.6) Thus, the interfacial concentration in the luminal film and submucosa can be related to the concentration of antigen in the epithelium. Combining eqs (3.1)-(3.3), (3.5), (3.6) and simplifying gives, 39 , , , , , . (3.7) Following this pattern described above and assuming that the generation of eotaxin-3 is solely attributed to antigen concentration and the generation of granule proteins is solely attributed to interaction with eotaxin-3, the transport of eotaxin-3 and granule protein is defined as , , , , , (3.8) and , , , , , (3.9) . Eqs (3.7), (3.8), and (3.9) are a series of coupled ordinary differential equations (ODEs) and can be solved analytically in a sequential manner. To facilitate the solution process, the above equations are simplified and rewritten to become, (3.10) , (3.11) , and The values of . (3.12) , (3.13) are defined as , , 40 , , , , , , , (3.14) , (3.15) , (3.16) , (3.17) , , , , , , , (3.18) , (3.19) and (3.20) , Solving the coupled ODEs in eqs (3.10)-(3.12) yields the following solution, Γ, 1 , Γ, 1 Γ, 1 Γ Γ , (3.21) 1 , Γ , , (3.22) 1 (3.23) Γ , 1 . Defining α as 1 , (3.24) The coefficients in eqs (3.21)-(3.23) become Γ, , (3.25) 41 1 Γ, , Γ, (3.26) , (3.27) Γ, (3.28) , Γ Γ , , (3.29) , , (3.30) and Γ . , (3.31) The solutions given above in eqs (3.21), (3.22), and (3.23) are the time-dependent concentration of the antigen, eotaxin-3, and granule protein in the epithelium. The solutions to the transient equations above are variants of the increasing form of the exponential decay equation. In the increasing form of the exponential decay equation, there exists an asymptotic solution as the independent variable goes to infinity. This implies that there is an asymptote as time increases for each of the derived equations. Taking the limit of eqs (3.21), (3.22), and (3.23) as → ∞ yields the steady state solution for each of the components as lim Γ, , → lim → and Γ, Γ, , (3.32) (3.33) 42 lim → Γ, Γ , Γ . , (3.34) The long-term, steady state, solutions given in eqs (3.32), (3.33), and (3.34) allow for a parametric assessment of the input parameters to those equations. Understanding that granule protein concentration from the eosinophils is the primary factor driving the diseased state, only eq (3.34) is specifically examined across a range of plausible inputs. These factors are examined to assess the long-term trends of granule protein concentration leading to a diseased state. The primary independent variable for the evaluation of eq (3.34) is the antigen concentration in the liquid along the epithelial wall, , with the dependent variable being the concentration of granule protein in the epithelium, lim . To better gauge granule → protein response to the physical inputs, the inputs are estimated using engineering first principles and fundamental assumptions. The physical inputs deemed to be significant, e.g., reaction rate of all the parameters, mass transfer coefficients, and concentration of each species, are varied to determine how the granule protein concentration changes as the physical inputs are perturbed. The blood flow removing biological constituents was assumed to be sufficient to rapidly remove the species within the submucosa. Thus, , , and , are assumed to be zero. The mass transfer coefficients are estimated by calculating the diffusion coefficients using the Stokes-Einstein equation and the heat and mass transfer analogy for a tube bank. This protocol was summarized by Incropera, et al.7 for diffusivity. The limiting case uses the transient mass transfer conduction model for the Biot number (Bim≈0.1). A conservative approach to mass transfer along the interfaces of the epithelium is implemented by assuming an uninterrupted flow of species across the 43 epithelial interfaces, thus the thermodynamic factors (i.e., the values) in eqs (3.14),(3.16), and (3.19) are assumed to be unity. The reaction rates are estimated using data provided by Nobuhisa, et al.8 and Humbles, et al.9 using the binary reaction rate model as recommended by Fogler.12 The initial antigen concentration at the luminal film surface was estimated using the properties of casein in milk.13,14 Casein is a family name for related phosphorproteins commonly found in mammalian milk. Casein was chosen due to its history of triggering allergic responses. The volume and area of the epithelium are calculated using a typical adult esophagus, approximately 25 cm in length,15 and thickness of the epithelial layer of approximately 3 mm.16 The physical values of the epithelial layer are kept constant in this analysis, while attention is focused on the physical inputs believed to influence granule protein mass transfer in the epithelium. A summary of parameters and assumptions is given in Tables 3.1 and 3.2. 3.4 Results and Analysis The concentration of granule protein is calculated using eqs (3.32), (3.33), and (3.34). To calculate the influence of each parameter on the generation of granule protein, the solutions are calculated as a function of the antigen concentration in the luminal film ( ). A parametric analysis is performed using the physical values first estimated from the information given in Section 3.3 and then perturbing the parameter values by ±50%. The parameters examined in this work focus on reaction rates, concentrations of material at the luminal film interface, and mass transport coefficients in the luminal film and submucosa. The mass transport coefficients in eqs (3.12) through (3.20) are calculated 44 using diffusion coefficients and the lumped capacitance mass transfer analogies indicated in Section 3.3. The parametric analysis varies each physical input independently to reveal significant trends. Figures 3.3 through 3.5 show the concentration of granule protein in the epithelium as a function of the antigen concentration in the luminal film when these physical inputs are varied. In general, each perturbation of the inputs had an impact on the concentration of granule protein in the epithelium. Only one parameter, the mass transfer coefficient of antigen in the submucosa, did not significantly change the concentration of the granule protein in the submucosa and is not shown in this calculation. Figure 3.3 shows an increase in granule protein concentration within the epithelium as the concentrations of eotaxin-3 and granule protein are increased within the luminal film, where the intercept along the granule protein concentration in the figures is more impacted by the concentration of eotaxin-3 in the luminal film (Figure 3.3a). Figure 3.4 shows an increasing trend in granule protein concentration within the epithelium as the mass transfer coefficients for eotaxin-3 and antigen are increased in the luminal film. Conversely, a decrease in granule protein concentration shows an increase in mass transport of granule protein in the fluid. Figure 3.5 shows an increase for granule protein concentration within the epithelium as the mass transfer coefficient for eotaxin-3 in the submucosa is increased. A decrease in granule protein concentration is observed when an increase in mass transport of granule protein in the submucosa occurs. Figure 3.6 shows an increasing trend in granule protein concentration within the epithelium as the reaction coefficients for eotaxin-3 and granule protein are increased, as opposed to a 45 decreasing trend in granule protein concentration with an increase in reaction coefficient for antigen. In all of the figures cited, a positive, increasing, linear trend is observed. Depending on which input is perturbed, either a change in the granule protein concentration intercept (i.e., y-intercept) or the slope is experienced. The significance of difference in granule protein intercept locations is highly dependent on the influence the input has on the equations. An input that influences the slope and intercept would have the most impact on the mass transport of granule protein in the epithelium. It is clear from the plots that the mass transport physical inputs significantly impact the granule protein concentration, specifically, the mass transport parameters of the antigen and granule protein in the luminal film. Additionally, all of the reaction rates significantly influence the concentration of granule protein in the epithelium. Although the perturbations in the physical input's effect on the intercept and slope differ in magnitude and severity, all of the trends are linear and increasing as the concentration of antigen is increased. This result is correct mathematically, but imposes very specific physical assumptions on the system. It would be expected that the concentration of granule protein should, at some point, have a decreasing trend when a specific key parameter is sufficiently small. This is shown to not be correct when examining the equations as the system approaches a steady state. Additionally, the primary assumption used in the lumped capacitance method is valid in the ranges of small Biot numbers (Bim<0.1). This assumption was based on limited information concerning the mass transfer and reaction of granule protein in the 46 epithelium. This model is an acceptable first treatment of a biophysical mass transport model driving eosinophilic esophagitis, but a more complete model should be pursued. 3.5 Conclusion The analysis described here yields a first assessment of a biophysics-based mathematical model to describe the concentration of granule protein within the epithelium of the esophagus. The equations used in this model yield an asymptotic solution as the limit of time is taken to infinity, giving a series of steady state equations. From these steady state equations, a parametric analysis was performed to explain granule protein concentration trends by perturbing the key physical inputs influencing the diseased condition. The concentration of granule protein in the epithelium linearly increases using this steady state model. The mass transfer of the antigen and granule protein along with all of the reaction rates has a significant impact on the slope and intercept of the granule protein-antigen linear relationship. It is believed this is primarily due to the linear interaction of these parameters in the transport equations. The model also presents a continually increasing steady state model of the concentration of granule protein. Even when some parameters decrease the transport (e.g., mass transfer of granule protein in the luminal film), the slope is never negative. Valuable information is obtained through this first effort in defining a mathematical representation of the pathogenesis of EoE. The model does indicate that the mass transfer coefficients and reaction rates have a significant impact on long-term growth of 47 eosinophils in the epithelium. From this evaluation, the main driving force behind EoE are the magnitudes of the mass transfer coefficients and reaction rates. Key components driving the pathogenesis of a diseased state will eventually diffuse out of the epithelium when an adverse stimulus is no longer present. The equations given make the prediction that a continual increase will occur and never reach a state of decreasing concentrations. Also, the linearized transport assumption used in the lumped capacitance method may not be adequate to fully describe the diseased state. A more detailed model must be considered to explore the stability of the diseased state and determine which parameters have the highest impact. 48 Figure 3.1 - The basic pathogenesis of EoE showing the relation between antigen/allergic response and eosinophil granule protein release (courtesy of Rothenberg, Gastroenterology, 2009, 137, 1238-1249). This model is the basis of this study and simplifications to this pathogenesis model are made to generate a mathematical model using engineering first principles. 49 Fluid Cil Epithelium Submucosa Muscular Cif Cie Cis CV Ci∞ Figure 3.2 - A notional one-dimensional model for transient conduction in the esophagus with the concentration profiles of a generic biological component, i, shown in the image. A control volume (CV) is made around the epithelial layer of the esophagus, where a species balance is performed. The lumped capacitance method assumes the concentration of the biological component, i, is constant in the epithelial layer. 50 Figure 3.3 - Concentration of granule protein in the epithelium as a function of the antigen concentration at the luminal film as → ∞ for (a) concentration of eotaxin-3 in the luminal film and (b) concentration of granule protein in the luminal film. The plots show increasing trends in granule protein concentration within the epithelium as the concentrations of eotaxin-3 and granule protein are increased. 51 Figure 3.4 - Concentration of granule protein in the epithelium as a function of the antigen concentration at the luminal film as → ∞ for (a) mass transfer coefficient of eotaxin-3 in the luminal film, (b) mass transfer coefficient antigen in the luminal film, and (c) mass transfer coefficient of the antigen in the luminal film. The plots show positive trends in granule protein concentration within the epithelium as the mass transfer coefficients for eotaxin-3 and antigen are increased in the luminal film, as opposed to a decrease in granule protein concentration with an increase in mass transport of granule protein in the fluid. 52 Figure 3.5 - Concentration of granule protein in the epithelium as a function of the antigen concentration at the luminal film surface as → ∞ for (a) mass transfer coefficient of eotaxin-3 in the submucosa and (b) mass transfer coefficient of granule protein in the submucosa. A positive trend exists for granule protein concentration within the epithelium as the mass transfer coefficient for eotaxin-3 in the submucosa is increased, while there is a decrease in granule protein concentration with an increase in mass transport of granule protein in the submucosa. There is no significant impact on granule protein concentration for the mass transfer of the antigen in the submucosa and is thus not shown here. 53 Figure 3.6 - Concentration of granule protein in the epithelium as a function of the antigen concentration at the liquid surface as → ∞ for (a) concentration of eotaxin-3 in the luminal film edge, (b) concentration of granule protein in the luminal film edge antigen reaction rate, and (c) diffusion coefficient of the antigen in the epithelium. The plots show positive trends in granule protein concentration within the epithelium as the reaction coefficients for eotaxin-3 and granule protein are increased, as opposed to a decrease in granule protein concentration with an increase in reaction coefficient for antigen. 54 Table 3.1 - The calculated values of the physical parameters used in the model. The diffusion coefficients were calculated using Stokes-Einstein's equation and the reaction rates were estimated using data from Nobuhisa, et al.8 and Humbles, et al.9 Parameter Diffusion Coef (cm2/s) Reaction Rate (1/s) Species Antigen Calculated Value Eotaxin-3 4.993x10 Eosinophils 6.773x10 Antigen/Eotaxin-3 5.794x10 Eotaxin-3/Eosinophils 1.786x10 -6 1.564x10 -6 -6 -7 -5 55 Table 3.2 - A summary of the assumptions made in the development of the model. Assumptions One-Dimensional in the r-direction Constant properties throughout geometry Linear concentration changes within the volume Biot number is <<0.1 Reactions are all first order Epithelium is a constant matrix of material Constant mass fraction of 0.383 for casein entering the system at the boundary Concentrations of components from the boundaries go to zero Time independent 56 3.6 References (1) Spergel, J. M. Genome Medicine 2010, 2. (2) Allahverdian, S.; Harada, N.; Singhera, G. K.; Knight, D. A.; Dorscheid, D. R. American Journal of Respiratory Cell and Molecular Biology 2008, 38, 153. (3) Rothenberg, M. E. Gastroenterology 2009, 137, 1238. (4) Blanchard, C.; Mingler, M. K.; Vicario, M.; Abonia, J. P.; Wu, Y. Y.; Lu, T. X.; Collins, M. H.; Putnam, P. E.; Wells, S. I.; Rothenberg, M. E. The Journal of Allergy and Clinical Immunology 2007, 120, 1292. (5) Phipps, S.; Ying, S.; Wangoo, A.; Ong, Y.-E.; Levi-Schaffer, F.; Kay, A. B. Journal of Immunology 2002, 169, 4604. (6) Gharaee-Kermani, M.; Phan, S. H. International Journal of Molecular Medicine 1998, 1, 43. (7) Bergman, T. L.; Lavine, A. S.; Incropera, F. P.; DeWitt, D. P. Fundamentals of Heat and Mass Transfer; 7th ed.; Hoboken, NJ : John Wiley: Hoboken, NJ, 2011. (8) Terada, N.; Hamano, N.; Kim, W.; Hirai, K.; Nakajima, T.; Yamada, H.; Kawasaki, H.; Yamashita, T.; Kishi, H.; Nomura, T.; Numata, T.; Yoshie, O.; Konno, A. American Journal of Respiratory and Critical Care Medicine 2001, 164, 575. (9) Humbles, A. A.; Conroy, D. M.; Marleau, S.; Rankin, S. M.; Palframan, R. T.; Proudfoot, A. E.; Wells, T. N.; Li, D.; Jeffery, P. K.; Griffiths-Johnson, D. A.; Williams, T. J.; Jose, P. J. Journal of Experimental Medicine 1997, 186, 601. (10) Taylor, R.; Krishna, R. Multicomponent Mass Transfer; New York : Wiley: New York, 1993. (11) Seader, J. D.; Henley, E. J. Separation Process Principles; New York : Wiley: New York, 2006. (12) Fogler, H. S. Elements of Chemical Reaction Engineering; 4th ed., 2006. (13) Ribadeau-Dumas, B.; Grappin, R. Lait 1989, 69, 357. (14) Goff, H. D. In Dairy Science and Technology Education Series; University of Guelph https://www.uoguelph.ca/foodscience/book-page/milk-proteins 2015; Vol. 2015. (15) Kuo, B.; Urma, D. GI Motility Online 2006. 57 (16) Giuli, R. Barrett's Esophagus: Columnar Lined Esophagus : 250 Questions, 250 Answers; J. Libbey Eurotext, 2003. 58 CHAPTER 4 DIFFERENTIAL MODEL OF EOSINOPHIL GENERATION IN THE ESOPHAGUS 4.1 Abstract The evolution of eosinophilic esophagitis in the epithelium is heavily dependent on the physical transport parameters dictating the pathogenesis of the disease. It is hypothesized that the conditions leading to a disease (versus normal) state of an esophagus can be predicted by analyzing a differential reaction-diffusion model for EoE based on an interacting protein network. This study examines the influence of the physical parameters governing both the mass transfer and reaction network on granule protein concentrations within the epithelium. Through the use of computational experiments, continual growth trends occur as physical parameter magnitudes increase. A Plackett-Burman screening design of experiments is completed to identify significant parameters governing mass transfer evolution. The significant parameters governing granule protein evolution in the epithelium suggest treatment pathways for EoE patients. 4.2 Introduction and Background In the pursuit of a biophysical description for eosinophilic esophagitis (EoE), a mathematical (systems biology) model is developed using fundamental transport principles, similar to what was performed in Chapter 3. Following the same approach to 59 pathogenesis, as in Chapter 3, EoE begins with the introduction of an antigen within the esophagus, causing the activation of T-helper cells (Th2), releasing the cytokine interleukin proteins IL-5 and IL-13.1 IL-13 promotes epithelial cell repair and causes the production of eotaxin-3 from the epithelial cell.2,3 Both IL-5 and eotaxin-3 promote the activation of eosinophils, releasing major basic proteins (MBP-1 and MBP-2), eosinophil peroxidase (EPO), eosinophil cationic protein (ECP), and eosinophil-driven neurotoxin (EDN).3,4 When activated, eosinophils also produce transforming growth factor beta (TGF-β), which continues to activate additional epithelial cells and contributes to hyperplasia and fibrosis.5,6 Figure 4.1 is a graphical representation of the pathogenesis process as just described. Here, a slightly more complex mass transport process is approached. Using the principles of differential mass balances, a system of linear partial differential equations, describing the mass transport and reaction equations, is derived. These model equations are analyzed to determine inherent instability. The equations are also assessed using computational experimentation using a pseudo Monte Carlo method and design of experiments to identify significant trends in the estimated parameters. The parameteric assessment identifies the factors most significantly impacting the concentration of granule proteins within the epithelium. 4.3 Methods A model to describe the migration of granule proteins from eosinophils is achieved by using a differential approach to the mathematical model for species conservation of granule protein generation. The previous analysis, in Chapter 3, made significant 60 assumptions in the system, most significantly, neglecting nonlinear variability in concentration throughout the esophagus and the relation between antigen, eotaxin, and granule protein concentrations. Relaxing these assumptions provides a convenient solution, but no significant insight is gained about a diseased condition. Through the species conservation approach, the conditions where the system would be stable or unstable are reached by performing a classical linear stability analysis.7 Additionally, a nondeterministic modeling approach, where variable input parameters are used to generate pseudo-random output results, is used to examine the parameter trends. The Latin Hypercube Sampling (LHS)8 method is used to generate a pseudo-random range of input parameters, generating an array of pseudo-random outputs. In utilizing a nondeterministic approach, trends can be determined that lead to the effect parameters will have on the output of the chemical diffusion reaction system. A significant weakness to the nondeterministic approach is predictions are only as reliable as the prior information used. If no data are present, a model calibration cannot be reached and the model is a prediction of the expected outcome defined by the model equations. To determine the sensitivity of the model inputs to the given outputs, a Plackett-Burman9 sensitivity analysis is executed. The Plackett-Burman analysis is a general method that can be utilized to provide insight as to which parameters significantly contribute to a diseased esophagus. 4.3.1 Esophagus Model Definition To generate a more sophisticated model, the linear assumptions outlined in Chapter 3 are removed from the derivation. Following the same pathogenesis process outlined in Chapter 3 with the antigen-eotaxin-3-granule protein interaction simplification, the 61 equations are derived for the region of the epithelial layer in the esophagus. EoE primarily is manifest in the epithelium and only this region is examined in this assessment. Applying a general species mass balance to a control volume around the epithelium (see Figure 4.2), the system is described by, 1 , (4.1) for the antigen, 1 , (4.2) for the eotaxin-3, 1 and for the granule protein, where , , , (4.3) , , and , are mass fraction, diffusion coefficient, reaction rate, radial position, and time, respectively. Assuming a constant concentration of antigen for a given length of time, in addition to a linear flux of species out of the epithelium for all species into the blood stream, the boundary conditions become, Δ (4.4) | 0 and ,0 , , (4.5) where is the maximum radial thickness of the esophagus. Eqs (4.1) through (4.3) define the transport in the epithelium, while the conditions defined in eqs (4.4) and (4.5) give 62 the initial and boundary conditions resulting in antigen exposure and the flux of species due to blood flow at the epithelial wall. 4.3.2 Linear Stability Analysis The analytical solutions to eqs (4.1) through (4.3) are not easily obtainable. Alternative methods to assess the system's behavior can be obtained by using methods from hydrodynamic stability.7,10,11 The key concept with hydrodynamic stability is examining a stationary system under an infinitesimal perturbation. Under this perturbation, the system will exhibit specific characteristics. When the system is disrupted, three characteristic behaviors may occur. If a perturbed system reacts with a diminishing response for all inputs the system is characterized as stable. If the amplitude of a system's response increases for all inputs, such that a progressive departure from the stationary state occurs, the system is characterized as unstable. When the system's response diminishes for some inputs and departs from the stationary state for other inputs, the system is considered marginally stable.11 Application of hydrodynamic stability methods can be applied to a wide variety of applications. Although the physics involved in defining a system may be different from hydrodynamics, the mathematics are similar and the same methods can be applied.7 The most fundamental approach to hydrodynamic stability is to use a linear or weakly nonlinear stability analysis.7 In linear stability analysis, the equations defining the physics of a system are linearized for small perturbations.7 To begin this process for assessing EoE disease, the solutions to eqs (4.1), (4.2), and (4.3) are assumed to be a linear combination of a steady state solution and a perturbed solution, 63 . (4.6) Typically, the stability of a system is assessed against the steady state solution under small perturbations.7 Through this assessment, the growth or decay in the perturbation can be ascertained. For a complete analysis to be performed, the analytical solution for the steady state condition is recommended. Due to the coupled nature of the equations used to describe EoE in Section 4.3.1, an analytical form for the steady state solution does exist. Although the analytical steady state solution cannot be utilized for a complete analysis, the equations can be simplified to examine only the unstable components in eq (4.6). This approach is ideal when the profile under small perturbations (e.g., concentration profile) is of no concern and the analysis is focused under conditions the equations would be stable or unstable. Thus, to determine the perturbed state, only the perturbed operator is required.10 The perturbation analysis begins by nondimensionalizing eqs (4.1), (4.2), and (4.3) to give, for antigen, 1 ∗ ∗ ∗ ∗ ∗ , (4.7) for eotaxin-3, 1 ∗ ∗ ∗ ∗ , ∗ (4.8) and for granule protein, 1 ∗ ∗ ∗ ∗ ∗ . (4.9) 64 where, and are the characteristic time and length, respectively. Similarly, the steady state equations are defined for antigen by, 1 0 ∗ ∗ ∗ , ∗ (4.10) for eotaxin-3 by, 1 0 ∗ ∗ ∗ , ∗ (4.11) and for granule protein by, 1 0 ∗ ∗ ∗ . ∗ (4.12) The linear combination in eq (4.6) is inserted into eqs (4.7), (4.8), and (4.9), and steady state terms are eliminated by subtracting the steady state equations eqs (4.10), (4.11), and (4.12) to yield for antigen, 1 ∗ ∗ ∗ ∗ ∗ , (4.13) for eotaxin-3, 1 ∗ ∗ ∗ ∗ , ∗ (4.14) and for granule protein, 1 ∗ ∗ ∗ ∗ ∗ . (4.15) Eqs (4.13), (4.14), and (4.15) are the perturbed state of the EoE system defined in Section 4.3.1. 65 To examine the stability, a complex growth operator is selected to identify normal modes in the system. Traditionally, the complex operator is assumed to be a linear superposition of imaginary exponential equations. This is not a practical approach with these equations due to the cylindrical coordinate system from the esophagus. Through the Bessel function of the first kind,12 the perturbation for the i-th species is given by, , where ,… (4.16) , , is the wave vector, pointing in the direction of the wave propagation, and is the growth rate of the wave. Using the relation given in Deen,12 (4.17) 0, to evaluate the cylindrical coordinate derivative of the Bessel function of the first kind in eq (4.16) yields the solution, 1 , where , , (4.18) is the magnitude of the wave vector, or generally, the wave number.13 Only the magnitude of the wave vector is considered since the direction of propagation does not matter in this analysis. Inserting eq (4.18) into eqs (4.13), (4.14), (4.15), and simplifying gives the general perturbed equations, , , , , , , , , , , . (4.17) 66 Recognizing eq (4.17) forms a system of algebraic equations, it is rewritten in matrix notation, , , , 0 (4.18) 0 , 0 , . , 0 For eq (4.18) to have nontrivial solutions, the condition 0 0 0 (4.19) 0 0 must be satisfied. From eq (4.19), the growth rate is defined through the characteristic equation, (4.20) 0. Thus, solving for the growth rate gives, . (4.21) 67 The stability of the system is defined through the magnitude of the growth rate. The growth rate can be defined in terms of real and imaginary components. For this assessment, only the magnitude of the real component is of importance. If instability is present in the system, the growth rate must be positive for at least one wave number. Stability is present when the real component of the growth rate is negative for all wave numbers. In eq (4.21), all of the possible solutions to the growth rate are negative, meaning the system is unconditionally stable. This result holds true for multiple dimensions. Only a negative diffusion coefficient for any one of the species examined through this model would cause instability, resulting in an accumulation of eosinophils in the esophagus, leading to a diseased state. The physical phenomena of negative diffusion is possible in nature and is a result of interparticle interactions or various thermodynamic forces.14,15 These typically result in either spinodal decomposition16 or first order phase transition by nucleation.14 In the case of spinodal decomposition, like species are highly attracted to one another. This attraction is significant enough to cause a flux in the direction up the gradient, giving a diffusion coefficient which is effectively negative. Spinodial decomposition typically is manifest macroscopically as packets of "concentration waves".16 In the case of nucleation, phase concentration is reduced in the immediate location of a growing particle.14 The components feeding the nucleation process flow through a normal diffusion process, down the concentration gradient, but a negative concentration gradient occurs in the global frame of reference.14 Nucleation results in particle growth and coalescence, which may occur in all phases of material. 68 In either process of negative diffusion, no significant research has proven eosinophils, or any of the other biophysical components in eosinophilic esophagitis, resemble either spinodal decomposition or nucleation characteristics. Without evidence of negative diffusion, a parametric study of eqs (4.1), (4.2), and (4.3) is sought to better understand the influence of the physical parameters on mass transfer. 4.3.3 Nondeterministic Computational Modeling In assessing the pathogenesis process outlined in Figure 4.1 and determining the significant interactions affecting EoE, certain protein-cell interactions can be neglected or integrated into other processes to make a tractable and discernable model. The emphasis of the mathematical model development is to create a systems biology description, which emphasizes the evolution of eosinophil granule protein generation, thus interactions only influencing eosinophil activation are pertinent. With an analytical solution being impractical, if not impossible to obtain, a computational solution is needed. Generally speaking, any description of a physical system may be summarized in the form of, , where is the vector of input variables, methods, the function ,…, ∈ , (4.22) is the solution via analytical or numerical is the model describing the system, and is the input variable space. In scientific investigation, the ability to understand the relationship between the system's response, , and a given set of inputs, , is needed to make informed decisions about the system.17,18 69 If the relationship between inputs and outputs are explored experimentally, the resulting tests may be costly and potentially prohibitive to execute. Significant research has been performed to explore more efficient approaches to minimize the burden of physical experimental studies while maximizing the information gain. Thus, alternative approaches using models to explore a system's relationships are more desirable due to the limitations of physical experiments. With the development of more efficient computers, numerical methods have become a critical component in scientific exploration. A model given by eq (4.22) is often difficult to solve analytically and one must typically resort to computational methods to determine an approximate solution. Through in silico experimentation, these relationships can be examined and valuable information can be obtained easily, but can quickly become burdensome if the system is too complex. Consequently, a primary objective in computer experimentation is to find a system definition that is simple enough to compute a solution while simultaneously providing an adequate description of the physics. In silico experiments are very powerful because large numbers of input variables can be explored simultaneously. A significant challenge with the numerical solution to eq (4.22) is that the solution is deterministic, in which the solution for a given set of inputs will always produce a specific set of outputs. To overcome the challenge of determinism within a numerically determined solution, additional processing of the input variables can be performed to create a nondeterministic algorithm. The most common approach to processing of the input variables is to introduce randomness into the input variables.19 The process of introducing 70 randomness can be performed in several ways, but can give significant information as to the basic characteristics of the system. The process of in silico experimentation provides unique abilities and requires special attention to the experimental design approach to produce beneficial information. In this work, a Latin Hypercube Sampling (LHS) algorithm is used to randomized the inputs. LHS was proposed by McKay, Beckman, and Conover in 19798 and is considered as the seminal paper in the design of computer experiments.17 If the inputs into the computational model are not precisely known, but a range of values can be postulated, either through statistical interpretation or estimation, LHS is an ideal method. LHS divides the domain of all plausible variables into n strata of equally probable and uniform intervals. The variable intervals are arranged to create a multidimensional cube (a.k.a., a Latin hypercube) of equally probable interval intersections. The interval intersections are selected pseudo-randomly to ensure only one intersection is selected for all rows and columns defining the variables' plausible range. With the intersection now chosen, discrete values are then selected at random within each intersection. As an example, assume the LHS method is applied to two variables. The first step in LHS divides the range of each variable into 5 intervals. The strata of variables are then arranged with one variable ( ) along the x-axis and the other ( ) along the y-axis of the domain. Intersections are selected pseudo-randomly to ensure no other intersections can be chosen in the same row and column of the current selection. With each selected interval, a discrete value is used for computation. Figure 4.3 shows an example of this process. 71 LHS is utilized to examine parameter correlations from eqs (4.1), (4.2), and (4.3). A time-dependent solution is calculated computationally using the finite difference methods outlined by Patankar20 and Ferziger.21 In the calculation, the esophagus is exposed to a constant concentration of antigen for 30 minutes. The diffusivities and reaction rates in eqs (4.1), (4.2), and (4.3) are varied giving four different parameters to compare against each concentration. Similar assumptions are used to determine initial parameter values from Chapter 3. The mass transfer coefficients are estimated by calculating the diffusion coefficients using the Stokes-Einstein equation and the heat and mass transfer analogy for a tube bank as summarized by Incropera, et al.22 for diffusivity. The reaction rates are estimated using data provided by Nobuhisa, et al.23 and Humbles, et al.24 to a linear reaction rate model recommended by Fogler.25 The initial antigen concentration at the luminal film surface was estimated using the properties of casein in milk.26,27 The calculated values (shown in Table 4.1) were varied significantly by several orders of magnitude to yield a wide range of potential responses. Table 4.2 summarizes the assumptions made in the development of the model. 4.3.4 Plackett-Burman Sensitivity Analysis A weakness of the LHS algorithm is it does not intuitively determine the impact of each input parameter on the overall output of the concentration profiles. To determine the significance of each parameter on the outcome of the model, a different computer experiment method is needed, in which the primary purpose is to "screen-out" the important parameters, or factors, from the less important ones. This type of experiment 72 is typically used for physical experiments, but can also be utilized for computational experiments, where the screened results can give better direction in future physical experimentation. A common screening design is the full factorial experiment design. This design gives a complete picture of the importance of each factor in addition to the influence on interacting factors. The primary disadvantage of a sophisticated experimental design like a full factorial approach is that as one increases the total number of factors analyzed, the total number of experimental runs (i.e., observations) will increase as number of effects and , where is the is the number of levels for each factor. Thus, for an experiment with two levels (i.e., a "high" and a "low" value) and six factors, the total number of observations would need to be 2 64, which is not ideal in a physical or computational experiment.28 Additional designs may be performed that will decrease the overall number of observations and are more ideal in an experiment where there are several factors being considered. One such design is the Plackett-Burman screening design. Originally developed in 1946 by statisticians Robin L. Plackett and J.P. Burman,9 it is a screening method that efficiently identifies the active factors of a system using as few experimental runs as possible. For two factor levels, the minimum number of observations required must be a multiple of four instead of a power of two for full and partial factorial designs. A Plackett-Burman experimental design is used to identify main factors within a system. Due to the properties of a Plackett-Burman design, the main factors are estimated independently of each other, thus interacting effects between factors are not considered. Therefore, the Plackett-Burman design is used to study main effects. Since this study is 73 only concerned with main responses and not interactions, a two-level Plackett-Burman experiment is used for all of the variable parameters in the model. Table 4.3 shows the Plackett-Burman experimental design for this work. To assess the contribution of each parameter in the Plackett-Burman design, the contributions of each factor are evaluated in terms of a parameter called the contrast. The contrast is a variable that can be used to compare two or more factors by examining their differences with zero. Letting denote a contrast, the algebraic expression can be written as, , (4.23) which is subject to the condition that the coefficients, , sum to zero. Values of the contrast are obtained for each of the factors at a given response (e.g., integrated concentration, mean value, maximum value, etc.). To compare the contrast of each factor, the sum of squares for the contrast is calculated, i.e., . where is the sum of squares for the contrast, and (4.24) is the number of observations. The sum of squares for the contrast for each factor and response is compared to determine the significance of each factor in the system.28 Since the magnitudes of each concentration differ significantly, the responses for the sum of squares for the contrast are normalized against the largest magnitude of the sum of squares. Through comparing the normalized values for the sum of squares of the contrast, the Plackett-Burman design "screens-out" parameters that are important against parameters that do not impact the physical system under consideration. 74 Using the same mean parameters and solution method given in Section 4.3.3, the time-dependent solution is determined. To provide the positive and negative conditions needed for the Plackett-Burman approach shown in Table 4.3, the parameters were varied ±50%. The model was calculated for each observation and concentration profiles were collected. The integrated and maximum concentration for each species is used as the system response quantity. 4.4 Results and Analysis From the equations describing EoE transport and reaction in the epithelium, the mass fraction of each species is calculated over the length of the epithelium. The LHS analysis is performed by varying from the nominal parameters calculated (see Table 4.1) and repeating the calculation 500 times with an array of pseudo-random parameters. Figure 4.4 shows one example of a profile from that array. To compare the relative outcomes as a function of the mass fraction of the species, each profile is integrated and compared to the parameters used to create the output. Figures 4.5 through 4.10 show the various integrated profiles as a function of the parameter value used to generate the profile. A strong correlation exists between the antigen diffusion and the integrated mass fractions of antigen (Figure 4.5a), eotaxin-3 (Figure 4.7a), and granule proteins (i.e., eosinophils) (Figure 4.9a). This is manifest by a discernible logarithmic growth trend in each of the figures. Antigen has the most significant correlation, although this is primarily due to the nature of the boundary condition, changes in mass fraction profiles, and the antigen diffusion coefficient. 75 The diffusion coefficients for eotaxin-3 and eosinophils do not have a significant impact on the model's outputs. The pathogenesis is linear, requiring an initial prime mover, the antigen, to start the process. Without the transport of antigen into the epithelium, no other species materialize. Thus, a strong correlation should be expected between the antigen diffusion and the mass fraction profiles of the other species in the epithelium. There is a strong correlation between the reaction rate for antigen-eotaxin-3 reaction rate, in eqs (4.1) and (4.2), and the mass fraction of eotaxin-3 (see Figure 4.8a). Additionally, a three-way correlation also occurs between the rate for antigen-eotaxin-3 reaction and the eotaxin-3-eosinophil reaction rate, in eqs (4.2) and (4.3), and the eosinophil mass fraction within the epithelium (see Figure 4.10). This is manifest in a more linear growth, showing constant growth in the mass fractions as reaction rates increase. From the LHS assessment, the increase in the magnitude of specific parameters resulted in the emergence of a consistent and continual growth pattern. Thus, it is important to understand what physical parameters have a significant impact on the mass fraction profiles within the epithelium. Performing the Plackett-Burman screening test on the system, the sum of squares of the contrast for each parameter can be identified and compared. Figures 4.11 and 4.12 show the normalized sum of squares of the contrast for each species and physical parameter. The normalized maximum mass fraction within the epithelium is shown in Figure 4.11. The maximum mass fraction is not shown for antigen because the maximum value occurs at the boundary condition and is unchanged 76 throughout the analysis. The Plackett-Burman analysis shows that the only significant responses in the pathogenesis of EoE are reaction rate coefficients ( and ) for the mass fractions for eotaxin-3 and granule proteins in the epithelium. The normalized integrated mass fraction is shown in Figure 4.12. This result differs slightly from the maximum concentration, which allows for the inclusion of the antigen mass fraction responses. Significantly, the main diffusion response occurs only for antigen diffusivity for each of the species, with the most significant response being between the integrated antigen mass fractions. The reaction rates also have a significant impact on the integrated mass fractions for eotaxin-3 and granule proteins, similar to the maximum concentration outcome. These results correlate with the LHS result shown in Figures 4.5 through 4.10, but with the benefit of a weighting scale on which the significance of each parameter may be placed. 4.5 Conclusion The analysis performed here provides a differential assessment for a conservationbased biophysical model to describe the concentration of eosinophils within the epithelium of the esophagus. A linear stability analysis was performed on the equations describing the transport of the main components, antigen, eotaxin-3, and granule proteins from the eosinophils. Linear stability analysis showed that the equations describing the system were unconditionally stable. This stability implied a constant growth of concentrations within the epithelium without the interaction from alternative phenomena, like thermodynamic forces. 77 Without evidence indicating thermodynamic phenomena, like spinodal decomposition or interparticle interactions, it is safe to assume increasing parameter values would cause increasing concentrations of granule proteins within the epithelium. A parametric analysis to yield parameter growth interactions on the mass fraction outputs was completed using LHS, a pseudo Monte Carlo method of parameter exploration. The LHS analysis showed key parameters, which highly influence the mass fraction of species within the epithelium. Antigen diffusion appeared to be the only transportspecific parameter affecting the mass fraction of all the species within the epithelium, due to the fact that antigen is a marker to promote the response and activation of eotaxin-3 and granule proteins in the pathogenesis process. This phenomenon resulted in a logarithmic increase in concentration as the diffusion coefficient is increased. Reaction rates have a discernable impact on the mass fraction profiles of eotaxin-3 and the granule protein, which resulted in a linear increase as the reaction rates are increased. To quantify the significance of parameter influence on the transport equations, a Plackett-Burman screening test was completed. Plackett-Burman screening tests provide a method to determine the significant parameters within an experiment, thus screening trivial information from the significant factors. Performing a computational factorial experiment using the Plackett-Burman method determined the significant parameters affecting the diseased state of EoE. The reaction rates of antigen-eotaxin-3 and eotaxin3-granule protein were the most significant factors driving the presence of granule proteins within the epithelium. The results from the Plackett-Burman analysis reflect what is observed in the LHS study. If the equations used to describe the pathogenesis are correct and are stable, a 78 continual increase in the concentration of granule protein is expected. Thus, to prevent a diseased condition, one must determine a way to retard or halt the activation of eosinophils from markers like eotaxin-3 or determine a way to prevent the transport of antigen into the epithelium. 79 Figure 4.1 - The basic pathogenesis of EoE showing the relation between antigen/allergic response and eosinophil granule protein release (courtesy of Rothenberg, Gastroenterology, 2009, 137, 1238-1249). This model is the basis of this study and simplifications to this pathogenesis model are made to generate a mathematical model using engineering first principles. 80 Ji  hi i  Ji  0 Figure 4.2 - A notional one-dimensional model for transport of species in the epithelium of the esophagus with the expected concentration profile of a generic biological component, i, shown in the image. A control volume (CV) is made around the epithelial layer of the esophagus, where a species balance is performed. 81 Randomly Selected Interval Intersections 1 Randomly Selected Value in Intersections 1 v v Parameter Interval 0 0 Parameter Interval u 1 0 0 Figure 4.3 - Example of LHS with u for two variables 1 82 Figure 4.4 - The mass fraction of (a) antigen, (b) eotaxin-3, and (c) eosinophils as a function of the distance into the epithelium. 83 Figure 4.5 - Integrated antigen mass fraction as a function of the (a) antigen diffusion coefficient, (b) eotaxin diffusion coefficient, and (c) eosinophil diffusion coefficient. A logarithmic growth pattern can be observed in the mass fraction as a function of the antigen diffusion coefficient. The impact of this parameter is due to the serial nature of the mass transport equations derived in the pathogenesis of EoE. The high correlation is also due to the fixed boundary condition of the antigen coming into the epithelium. The lack of correlation in (b) and (c) indicates a poor correlation between the parameter and the output of the model. The similarity of these images to a "scatter" plot implies no correlation. 84 Figure 4.6 - Integrated antigen concentration as a function of the (a) the antigeneotaxin-3 reaction rate (i.e., ) and (b) the eotaxin-3-eosinophil reaction rate (i.e., ). No significant correlation is observed, although a slight decrease can be seen with the antigen-eotaxin-3 reaction rate, which is not significant enough to be considered relevant. The lack of correlation in (a) and (b) indicates a poor correlation between the parameter and the output of the model. The similarity of these images to a "scatter" plot implies no correlation. 85 Figure 4.7 - Integrated eotaxin-3 mass fraction as a function of the (a) antigen diffusion coefficient, (b) eotaxin-3 diffusion coefficient, and (c) eosinophil diffusion coefficient. A logarithmic growth pattern can be observed in the mass fraction as a function of the antigen diffusion coefficient. The impact of this parameter is due to the serial nature of the mass transport equations derived in the pathogenesis of EoE. The lack of correlation in (b) and (c) indicates a poor correlation between the parameter and the output of the model. The similarity of these images to a "scatter" plot implies no correlation. 86 Figure 4.8 - Integrated eotaxin-3 concentration as a function of the (a) the antigeneotaxin-3 reaction rate (i.e., ) and (b) the eotaxin-3-eosinophil reaction rate (i.e., ). A significant correlation is observed for the antigen-eotaxin-3 reaction rate, which manifests in a linearly increasing integrated mass fraction of eotaxin-3. The lack of correlation in (b) indicates a poor correlation between the parameter and the output of the model. The similarity of this image to a "scatter" plot implies no correlation. 87 Figure 4.9 - Integrated eosinophil mass fraction as a function of the (a) antigen diffusion coefficient, (b) eotaxin-3 diffusion coefficient, and (c) eosinophil diffusion coefficient. A logarithmic growth pattern can be observed in the mass fraction as a function of the antigen diffusion coefficient. The impact of this parameter is due to the serial nature of the mass transport equations derived in the pathogenesis of EoE. The lack of correlation in (b) and (c) indicates a poor correlation between the parameter and the output of the model. The similarity of these images to a "scatter" plot implies no correlation. 88 Figure 4.10 - Integrated eosinophil mass fraction vs. the reaction rate of (a) the antigen-eotaxin-3 reaction rate (i.e., (i.e., ) and (b) the eotaxin-3-eosinophil reaction rate ). A significant correlation is observed for the antigen-eotaxin-3 and the eotaxin- 3-eosinophil reaction rate, which manifests in a linearly increasing integrated mass fraction of eosinophils (i.e., granule proteins) within the epithelium. 89 Normalized Sum of Squares of the Contrast Impact of Model Parameters on Maximum Concentration  Using  Plackett‐Burman Method 1 0.8 ωAntigen ωEO3 ωEosinophils 0.6 0.4 0.2 0 Diff of Antigen Diff of  Eotaxin Diff of  Eosinophils Rxn Coef Antigen Rxn Coef Eosinophils Parameter in P‐B Analysis Figure 4.11 - Summary of the sum of squares of the contrast for the maximum mass fraction for each parameter within the epithelium. The most significant parameters using this response quantity are the reaction rates for the antigen-eotaxin-3 (i.e., eotaxin-3-eosinophil (i.e., ). ) and the 90 Normalized Sum of Squares of the Contrast Impact of Model Parameters on Integrated Concentration  Using  Plackett‐Burman Method 1 0.8 ωAntigen ωEO3 ωEosinophils 0.6 0.4 0.2 0 Diff of Antigen Diff of  Eotaxin Diff of  Eosinophils Rxn Coef Antigen Rxn Coef Eosinophils Parameter in P‐B Analysis Figure 4.12 - Summary of the sum of squares of the contrast for the integrated mass fraction for each parameter within the epithelium. The most significant parameters using this response quantity are the diffusion of antigen and the reaction rates for the antigeneotaxin-3 (i.e., ) and the eotaxin-3-eosinophil (i.e., ). 91 Table 4.1 - The calculated values of the physical parameters used in the model. The diffusion coefficients were calculated using Stokes-Einstein's equation and the reaction rates were estimated using data from Nobuhisa, et al.23 and Humbles, et al.24 Parameter Diffusion Coef (cm2/s) Reaction Rate (1/s) Species Antigen Calculated Value 1.564x10 Eotaxin-3 4.993x10 Eosinophils 6.773x10 Antigen/Eotaxin-3 5.794x10 Eotaxin-3/Eosinophils 1.786x10 -6 -6 -6 -7 -5 92 Table 4.2 - A summary of the assumptions made in the development of the computational model. Assumptions One-Dimensional in the r-direction Constant properties throughout geometry Flux between elements in discrete geometry is first order an changes minimally over time Reactions are all first order Epithelium is a constant matrix of material Constant concentration of antigen at the boundary Time dependent, calculated profiles for 30 minutes of constant exposure Flux out of the epithelium assumed to be analogous to mass transfer in a tube bank with flux estimated as 5.26 10 m/s for antigen, 6.73 10 m/s for eotaxin-3, and 8.51 10 m/s for granule protein. 93 Table 4.3 - The Plackett-Burman design used to analyze the sensitivity of parameters in the Eosinophilic Esophagitis computational model. A high value is denoted as a " " symbol and a low value is denoted as a " " symbol in the table. For this study, only eight observations are required. Input Variables Observation A (DAntigen) B (DEO3) C (DEosinophils ) D (k Antigen) E (k Eosinophils ) Number 1 2 3 4 5 6 7 8 + + + + + + + + + + + + + + + + + + + + 94 4.6 References (1) Spergel, J. M. Genome Medicine 2010, 2. (2) Allahverdian, S.; Harada, N.; Singhera, G. K.; Knight, D. A.; Dorscheid, D. R. American Journal of Respiratory Cell and Molecular Biology 2008, 38, 153. (3) Rothenberg, M. E. Gastroenterology 2009, 137, 1238. (4) Blanchard, C.; Mingler, M. K.; Vicario, M.; Abonia, J. P.; Wu, Y. Y.; Lu, T. X.; Collins, M. H.; Putnam, P. E.; Wells, S. I.; Rothenberg, M. E. The Journal of Allergy and Clinical Immunology 2007, 120, 1292. (5) Phipps, S.; Ying, S.; Wangoo, A.; Ong, Y.-E.; Levi-Schaffer, F.; Kay, A. B. Journal of Immunology 2002, 169, 4604. (6) Gharaee-Kermani, M.; Phan, S. H. International Journal of Molecular Medicine 1998, 1, 43. (7) Drazin, P. G. Introduction to Hydrodynamic Stability; Cambridge University Press, 2002. (8) McKay, M. D.; Beckman, R. J.; Conover, W. J. Technometrics 1979, 21, (9) Plackett, R. L.; Burman, J. P. Biometrika 1946, 33, 305. 239. (10) Press, 2004. Drazin, P. G.; Reid, W. H. Hydrodynamic Stability; Cambridge University (11) Chandrasekhar, Clarendon Press, 1961. S. Hydrodynamic and Hydromagnetic Stability; (12) Deen, W. M. Analysis of Transport Phenomena; 2nd ed.; New York : Oxford University Press: New York, 2012. (13) Hata, S.; Nakao, H.; Mikhailov, A. S. Progress of Theoretical and Experimental Physics 2014, 2014. (14) Karpov, V. G. Physical Review Letters 1995, 75, 2702. (15) Argyrakis, P.; Chumak, A. A.; Maragakis, M.; Tsakiris, N. Physical Review B 2009, 80, 104203. (16) Binder, K.; Fratzl, P. In Phase Transformations in Materials; Wiley-VCH Verlag GmbH & Co. KGaA: 2005, p 409. 95 (17) Fang, K. T.; Li, R.; Sudjianto, A. Design and Modeling for Computer Experiments; CRC Press, 2005. (18) Santner, T. J.; Williams, B. J.; Notz, W. I. The Design and Analysis of Computer Experiments; Springer New York, 2013. (19) Oberkampf, W. L.; Roy, C. J. Verification and Validation in Scientific Computing; Cambridge University Press, 2010. (20) Patankar, S. Numerical Heat Transfer and Fluid Flow; Taylor & Francis, 1980. (21) Ferziger, J. H.; Perić, M. Computational Methods for Fluid Dynamics; Berlin, 1999. (22) Bergman, T. L.; Lavine, A. S.; Incropera, F. P.; DeWitt, D. P. Fundamentals of Heat and Mass Transfer; 7th ed.; Hoboken, NJ : John Wiley: Hoboken, NJ, 2011. (23) Terada, N.; Hamano, N.; Kim, W.; Hirai, K.; Nakajima, T.; Yamada, H.; Kawasaki, H.; Yamashita, T.; Kishi, H.; Nomura, T.; Numata, T.; Yoshie, O.; Konno, A. American Journal of Respiratory and Critical Care Medicine 2001, 164, 575. (24) Humbles, A. A.; Conroy, D. M.; Marleau, S.; Rankin, S. M.; Palframan, R. T.; Proudfoot, A. E.; Wells, T. N.; Li, D.; Jeffery, P. K.; Griffiths-Johnson, D. A.; Williams, T. J.; Jose, P. J. Journal of Experimental Medicine 1997, 186, 601. (25) Fogler, H. S. Elements of Chemical Reaction Engineering; 4th ed., 2006. (26) Ribadeau-Dumas, B.; Grappin, R. Lait 1989, 69, 357. (27) Goff, H. D. In Dairy Science and Technology Education Series; University of Guelph https://www.uoguelph.ca/foodscience/book-page/milk-proteins 2015; Vol. 2015. (28) Wheeler, D. J. Understanding Industrial Experimentation; Statistical Process Controls, Incorporated, 1988. 96 CHAPTER 5 RELEASE OF PHARMACEUTICAL COCKTAILS FROM SMALL POLYMER MICELLES CONTROLLED BY THERMODYNAMICS 5.1 Abstract Most emerging anti-cancer agents remain sparingly soluble in water, so they require carriers such as poly(ethylene glycol)-block-poly(D,L-lactic acid) (PEG-b-PLA) micelles to enhance their delivery. Shin, et al. used these micelles as vectors for combinations of three hydrophobic cancer drugs (paclitaxel, 17-allylamino-17-demethoxygeldanamycin, and rapamycin) to empirically determine that the half-life of agent release is a function of the oil-water partition coefficient. This is in contrast to their assertion that the cocktail release rate is diffusion controlled. We determine this unexpected, yet important, result using a lumped capacitance method. The proposed model provides excellent agreement with experimental release profiles provided by Shin, et al. indicating that the release rate is not controlled by simple diffusion. The model finds that the half-life is a function of the micelle radius, external mass transfer coefficient, and lumped partition coefficient. This explains the apparent dependence of the release rate profile on the oil-water partition coefficient. These findings simplify design and enable tuning the release of multidrug mixtures from micelles, which will be increasingly important since many new drug candidates are hydrophobic macromolecules. 97 5.2 Introduction and Background Advances in combinatorial chemistry and high throughput screening have accelerated the discovery of therapeutic agents. Although the use of these methods allows rapid evaluation of a wide array of drug candidates, a growing number of these molecules remain poorly water-soluble with limited bioavailability.1 Effective therapies must be soluble, release appropriately, and have excellent biodistribution.2 Therefore, significant research has focused on the development of carriers to enhance the solubility of poorly water-soluble pharmaceutical agents. Further complicating matters is an emerging trend in cancer therapy to combine pharmaceuticals to target multiple tumorigenic pathways to overcome drug resistance and tumor heterogeneity.3-5 Delivery of these potent yet poorly water-soluble drugs may be enhanced by low-molecular-weight surfactants (e.g., synthetic Kolliphor EL may be incorporated to stabilize nonpolar emulsions) to improve solubilization.2,6 However, the use of these excipients may induce cytotoxicity and accentuate drug instabilities.7 Improved encapsulation methods remain essential for effective delivery of combinatorial drug treatments. Polymeric micelles remain promising because they are versatile and have been used as drug carriers for decades.2,8-16 These spheroidal structures form by self-assembly of amphiphilic block copolymers in aqueous environments into a hydrophobic core and a hydrophilic shell (see Figure 5.1).2,17 The micelle's core serves as a micro-reservoir to encapsulate hydrophobic molecules, while the hydrophilic shell interfaces with biological media. Variations in micelle-forming chemistry facilitate customization of solubilization, biocompatibility, and biodegradability of micelles.18,19 98 In an intriguing experiment, Shin, et al. examined the solubilization of three hydrophobic cancer drugs (paclitaxel (PTX), 17-allylamino-17-demethoxygeldanamycin (17-AAG), and rapamycin (RAP)) in poly(ethylene glycol)-block-poly(D,L-lactic acid) (PEG-b-PLA) micelles.3,19 Paclitaxel stabilizes microtubule polymers,20 a major constituent of cell mitosis, thus inhibiting the mitotic process of cancer cells. 17-AAG, also known as Tanespimycin, is an Hsp90 inhibitor, preventing proper folding of proteins and destabilizing the proteins required for tumor growth.21 Rapamycin is an mTOR inhibitor, preventing cell growth and cell cycle progression through its ability to minimize growth factor stimuli.22 By combining these agents, Shin, et al. demonstrated the effectiveness of the drug cocktail, showing the agents together reduce tumor growth more rapidly than an individual drug treatment.19 In their experiment, they specifically evaluated the release rate of these agents from the micelle both individually and in combination. Each combination of drug loadings decreased agent solubility compared to their individual solubilities. For example, simultaneous release of PTX, 17-AAG, and RAP was slower in combination than independently, as measured by the half-life, t1/2, of drug release. Remarkably, Shin, et al. also asserted that t1/2 values correlated well with the logarithm of the oil-water partition coefficient. However, for nondegrading polymeric matrices, solute transport is typically assumed to be a diffusion-driven process, with transport typically defined via Fick's laws of diffusion.17,23,24 The timescale for diffusiondriven systems is given by the square of the radius divided by the diffusivity coefficient.25 Here the diffusivity is 1.3×10-12 m2/s, calculated from the Stokes-Einstein equation given a micelle diameter of 40 nm and viscosity of 8.5×10-4 Pa.s at 300 K.26 For systems 99 explored by Shin, et al., this characteristic time is on the order of 1 microsecond-orders of magnitude smaller than experimentally measured t1/2 values that are on the order of tens of hours. In this article, we resolve this inconsistency by developing a mechanistic model to explain the experimental results presented by Shin, et al. From first principles, we derive a two-parameter model that clarifies the mechanisms governing the half-life and allows straightforward fitting to experimental data. 5.3 Methods Figure 5.2 depicts a concentration profile for an agent (e.g., any of the three defined by Shin, et al.) within the micelle in four regions: the core, the hydrophobic shell (e.g., PLA), the hydrophilic shell (e.g., PEG), and an external aqueous media. Within each region, the concentration is constant (justified below) albeit time dependent. Partition coefficients (e.g., local thermodynamic equilibrium coefficient), also known as Kvalues,27 relate concentrations between regions. An integral balance around the exterior of the micelle gives the rate of agents leaving the micelle, , (5.1) where ntot represents the moles of agent in all three internal layers, kc is the external mass transfer coefficient, As is the surface area of the micelle, C3(t) is the concentration at the surface of the micelle as a function of time t, and C∞ is the concentration far from the micelle. The total moles of any individual agent within the spherical micelle may be determined by integration, 100 , 4 , . (5.2) Each region of constant composition must be integrated separately 4 (5.3) , yielding 4 3 (5.4) . (5.5) . where Vm is the micelle volume and K is a lumped partition coefficient, defined as 1 ≡ 1 . (5.6) The overbar denotes scaling with respect to distance from the center of the micelle (e.g., ⁄ , / ). Taking the derivative of eq (5.5) and substituting into eq (5.1) gives the agent concentration at the surface of the micelle as a function of time, . (5.7) Eq (5.7) is an ordinary differential equation that may be solved analytically by separation of variables. With the initial condition of C3(to)=C30, the solution to eq 5.7 becomes . (5.8) 101 This equation shows that the characteristic time for agent release is . 3 (5.9) Substituting eq (5.8) into eq (5.9) gives the total moles of drug in the micelle as a function of time, (5.10) . The interfacial concentration is related to the concentration of the core via . (5.11) The total percent released, P, from the micelle as a function of time is defined as ≡1 . (5.12) Substituting eqs (5.10) and 5.11 into eq (5.12) yields 1 1 , (5.13) which gives the percent drug release from the micelle as a function of the mass transfer coefficient, the ratio of central and infinite drug concentrations, and the lumped thermodynamic factor (i.e., kc, C∞/Coo, and K, respectively). This equation is used to analyze the experimental data reported by Shin, et al. Additional physical insight may be obtained from eq (5.13). The steady state concentration is →∞ 1 . (5.14) Likewise, the half-life is calculated by determining when half this concentration is reached via eq (5.13) as 102 ⁄ ln 2 ln 2 3 . (5.15) This expression shows that the half-life depends on the partition coefficients through K as was anticipated by Shin, et al., even though simple diffusion does not govern. With time as the only independent variable and the percent release as the dependent variable, only two parameters remain available to match each data set: a concentration ratio C∞/(kIII kII kI C30), and the characteristic time, defined in eq (5.9), herein termed the lumped parameter. To estimate each of these parameters, a nonlinear regression algorithm from the Statistics Toolbox® in the MATLAB® software tool was used. The regression algorithm uses an iterative least squares method to estimate the parameters as reported in Table 5.1 for each data set presented by Shin, et al. This model will be compared to well-known solutions based on Fick's second law. In its most complete formulation, the driving force, di, for mass transfer of each species, i, is its chemical potential, i. Taylor and Krishna26 show that the gradient in chemical potential is proportional to the gradient in species mole fraction, xi, , , where is the gas constant, number of species present, and is the absolute temperature, (5.16) is the pressure, is the is the thermodynamic factor.26 For ideal solutions, the thermodynamic factor becomes unity. For a pseudo-binary system, the diffusive flux vector becomes , (5.17) 103 where Di and Γi represent the binary diffusion coefficient and thermodynamic factor.26 Substitution into species conservation in the absence of convection and reaction (neither of which is likely significant for small micelles) gives 1 (5.18) in spherical coordinates. From eq (5.18), the time scale for diffusion may be expressed in terms of the thermodynamic factor as (5.19) . This expression resolves to the usual length squared over diffusivity ratio mentioned previously when thermodynamic nonideality remains negligible (i.e., 1). Smith, et al.27 show that the partition coefficient is related to the activity coefficient by taking the ratio of activity coefficients in each phase, ≡ , (5.20) where γi(1) and γi(2) are the activity coefficients of chemical species i in phase one and phase two, respectively. Likewise, the thermodynamic factor is related to the partition coefficient through, Γ 1 . (5.21) Thus, the partition coefficients at each interface shown in Figure 5.2 and the thermodynamic factor in eqs (5.20) and (5.19), respectively, may be estimated from an activity coefficient calculation. Typically, experimental phase equilibria data are not available to adequately define the activity coefficient. Semi-empirical methods are usually used to estimate the activity 104 coefficient. One of the best methods to estimate for liquid-liquid systems is the UNIFAC method.26,28-30 The UNIFAC method represents each pharmaceutical agent in terms of "structural units" or "subgroups."26-31 The UNIFAC method is a frequently used method when no thermodynamic data are available, typical of drug delivery formulations. For the purposes of this study, each agent was divided into subgroups as given in the 2011 UNIFAC group contribution by Wittig, et al.31 This generated a list of functional groups (see Table 5.1). To determine the value of the lumped parameter, K from eq (5.6), the relative thicknesses of each layer in Figure 5.2 must be estimated. Although the micelle diameter was reported as 36.9-43.8 nm in Table 1 of Shin, et al., h/R and H/R remain unknown but may be approximated using Kuhn's freely jointed chain model. The stretched out chain length is l=Nmlm=Nklk, where Nm is the number of mers determined in Table 5.3, lm is the length of a mer (0.148 nm for PEG32 and 0.295 nm for PLA33), Nk is the number of Kuhn lengths, and lk is the Kuhn length. The Kuhn length is 0.60 nm for PEG32 and 6.07lm12.0lm=1.79-3.54 nm for PLA.33,34 The root-mean-square end-to-end distance is then given as N1/2lk, which with Table 5.3 is 2.90 nm for PEG and 3.41-4.79 nm for PLA. The relative thicknesses then becomes (5.22) and (5.23) . Substitution gives ⁄ 0.583 0.712 and ⁄ 0.843 0.867. 105 Using eq (5.9) to describe mass transfer from an approximately 40 nm diameter stationary sphere in 300 K water gives a mass transfer coefficient of 6.36×10-5 m/s, the value used in this analysis. 5.4 Results and Analysis This article presents a simple model to describe the release of pharmaceutical agents from micelles - individually and in combination. This model, derived above, uses the lumped capacitance method instead of Fick's second law to describe drug concentration profiles (see Figure 5.2). The remainder of this article compares the model to experimental data, describes model parameters, evaluates the assumptions on which the model rests, and discusses when the model applies. Figures 5.3 through 5.5 quantitatively compare the model to experimental data from Shin, et al. The model (see eq (5.13)) fits the data well for each pharmaceutical agent individually and in combination. Each curve in Figures 5.3 through 5.5, calculated using eq (5.13), generally falls within the error bounds of each measurement. Occasionally (e.g., Figure 5.4c and Figure 5.5d), the curves falls slightly outside the experimental error bound. The quality of the fit remains independent of the number and combination of pharmaceutical agent within the center of the micelle. Combined, these data suggest that the proposed model accurately describes the delivery of a drug cocktail from a micelle. The model has only two parameters - as given in eq (5.13). The first parameter is a concentration ratio that vanishes when the solution far from the micelle does not contain the agent. Table 5.3 shows that this parameter is often nonzero, indicating some of the agent has "leaked" into the surrounding solution prior to the start of experimental 106 measurement. Although the value is small, the model accounts for this discrepancy and is not biased when fitting the data. The second parameter is the lumped time scale given in eq (5.9). This parameter depends on the micelle radius, the external mass transfer coefficient, and the lumped partition coefficient. The micelle radius remains readily tunable through the polymer composition and the relative size of hydrophilic and hydrophobic blocks. This has been reviewed by Poon and Andelman.35 Larger micelles increase the time scale and decrease the rate of agent release. The external mass transfer coefficient may be calculated from first principles as summarized by Deen25 or estimated using a heat and mass transfer analogy as summarized by Incropera, et al.36 For a micelle in quiescent solution, the Sherwood number is two. Recognizing this, the mass transfer coefficient may be determined directly (kc=ShΓiDi/(2R)). Where there is significant external flow relative to the micelle, Incropera, et al. 36 provide correlations for the Sherwood number as a function of the flow velocity. As the external velocity increases, the timescale for release shrinks and the rate of release goes up. This may be important in the circulatory system, for example, where the velocity of a micelle in the blood stream would vary considerably. Considering the circulatory system, these correlations combined with eq (5.13) provide an opportunity (outside the scope of this article) to estimate the rate of agent release in each portion of the blood stream from first principles. The final factor is the lumped partition coefficient. This factor may be determined directly from experimental data or estimated from thermodynamic correlations. Table 5.3 shows the characteristic time from eq (5.13) fitted from the experimental data from Shin, et al. 107 Alternatively, the lumped parameter may be estimated directly from the UNIFAC model. As explained above, the activity coefficients may be estimated by selecting UNIFAC subgroup parameters in the 2011 UNIFAC parameter update30 as shown in Table 5.2.31 Table 5.3 shows that the characteristic time may be predicted within a factor of 2 to 4 for highly hydrophobic agents, like PTX and RAP using the UNIFAC method to predict the lumped partition coefficient. There is significantly poorer agreement for the less hydrophobic drug, 17-AAG. This agreement remains far better than estimates from simple diffusion, which overpredict mass transfer by more than six orders of magnitude. This indicates that the UNIFAC method most likely overpredicts the partition coefficients. The half-life for each data set can be determined by using the empirical model fits and inserting estimated parameters into eq (5.13). The calculated half-lives are also shown in Table 5.3 and compared to the half-lives given by Shin, et al. Comparison between the model and the experiment demonstrates less than a 16% difference in all but one case. This is also evidenced qualitatively from the plots in Figure 5.3 through Figure 5.5. The most significant aspect of the model is that it does not use any transport properties that describe diffusion within the micelle, thus signifying that mass transfer within the micelle is not diffusion controlled. Indeed, none of the fitted parameters involve factors associated directly with diffusion within the micelle, although mass transfer outside of the micelle remains an essential factor in agent release. This result contrasts with the assertion given by Shin, et al., that delivery is purely diffusion driven. Additionally, the characteristic times show similar characteristics to the LogP calculation presented by Shin, et al. The LogP factor is a measure of the relative affinity of a 108 compound in nonionized water versus an immiscible model lipid solvent, which in this application is the pharmaceutical agent, at equilibrium. The value of LogP can be utilized as a gauge of the relative toxicity of a drug within the body.37 Thus, LogP values can be calculated to give a relative measure of the hydrophobicity of a molecule and by extension the relative toxicity of a pharmaceutical agent. The model presented here makes effective use of the "lumped capacitance method," commonly used in heat transfer analyses36 to provide an original and highly accurate solution to the problem of estimating drug release rate from hydrophobic/hydrophilic micelles. Formally, this transient mass transfer method applies only when convective mass transfer away from the micelle is significantly slower than diffusive mass transfer within the micelle. Such gradients would only persist outside the micelle.36 Traditionally, this requirement has been expressed via the mass transfer Biot number, Bim, Bi , (5.24) given here with the binary thermodynamic factor. A small Bim indicates that mass transfer is much greater within the micelle than outside of the micelle, and the lumped capacitance method is only valid when Bim<0.1.36 In contrast, the Bim for this system is approximated by using the mass transfer coefficient for a spherical micelle in quiescent media, using the Stokes-Einstein diffusivity approximation, and the thermodynamic coefficient from experimental data, yielding a Bim much larger than the valid range of less than 0.1. Thus, using common approximations for the system under assessment the lumped capacitance method should not be valid since Bim>>0.1. 109 However, the thermodynamics of phase separation resolve this apparent discrepancy. Using the UNIFAC method, the activity coefficients for each of the pharmaceutical agents in water may be calculated as a function of the mole fraction, shown in Figure 5.6. From Figure 5.6, ln(γi) of each component is initially positive and rapidly decreases as mole fraction increases. A negative value for ln(γi) denotes a region where there are negative deviations from ideal behavior (where adhesive forces exceed cohesive forces).27 The rate of mass transfer would thus be attenuated in this region of non-ideality and any agent with negative deviation would resist transport into an aqueous phase. In addition to this region of a negative deviation, spinodal decomposition occurs at low concentrations before the local minima of ln(γi) is reached (see Figure 5.6). Unstable phase separation occurs within this region of spinodal decomposition and suggests the possibility of hydrophobic molecule clustering on the outer edge of the micelle to form in essence a separate phase prior to separating from the micelle. As the derivative of ln(γi) goes to zero (i.e., the minima) the spinodal limit is approached, as shown by the solid black line in Figure 5.6. Therefore, mass transfer within the micelle is rapid, until the hydrophobic drug molecules reach the micelle surface where molecular forces effectively delay release into the water phase, satisfying the intent of the lumped capacitance method, even where the Bim restriction does not formally govern. Given the success of using internal diffusion to describe release rates of other systems, the question arises of when this lumped capacitance model should be employed as opposed to an internal diffusion model. The ratio of the characteristic time scale for 110 internal diffusion (eq (5.19)) to the characteristic time scale for the proposed model provides some insight: 3 Γ , (5.25) when diffusional forces within the micelle are much faster than other mass transfer processes, or when micelles are particularly small, the lumped capacitance method presented herein is valid for estimation of drug delivery rates. If diffusion is slow or the micelle is large, significant concentration gradients would exist violating the constant composition approximation required for our model. 5.5 Conclusions A lumped capacitance model similar to that often used in heat transfer calculations has been applied to the problem of estimating the rate of drug delivery from hydrophobic/hydrophilic micelles. The proposed lumped capacitance model for drug delivery shows excellent agreement with the experimental data of Shin, et al.19 The results demonstrate clearly that external mass transfer and not diffusion governs release from small micelles. 111 Figure 5.1 - Conceptual model of a PEG-b-PLA micelle containing a three-drug micelle prior to release. 112 Figure 5.2 - A conceptual one-dimensional representation of the expected concentration profile of a pharmaceutical agent within a micelle. The highest concentration of the pharmaceutical agent is located in the core of the micelle, shown as Co. The concentration decreases step-wise as distance from the center of the micelle increases. Concentrations within the hydrophobic PLA shell from h to H and the hydrophilic PEG shell from H to the micelle radius, R, are C1 and C2, respectively. The micelle surface concentration, located immediately outside the micelle, is represented as C3, with a concentration of C∞ far from the micelle surface. The partition coefficient at each interface of the micelle arises from the concentration ratio of adjacent phases (shown as kI, kII, and kIII). 113 Figure 5.3 - The fitted curves using eq (14) for 17-AAG loaded (a) individually, (b) with PTX, (c) with RAP, and (d) with PTX and RAP against experimental data provided by Shin, et al. 114 Figure 5.4 - The fitted curves using eq (14) for PTX loaded (a) with RAP, (b) with 17AAG, and (c) with RAP and 17-AAG against experimental data provided by Shin, et al. 115 Figure 5.5 - The fitted curves using eq (14) RAP loaded (a) individually, (b) with PTX, (c) with 17-AAG, and (d) with PTX and 17-AAG against experimental data provided by Shin, et al. 116 3 2 ln(γi) 1 0 ‐1 ‐2 ‐3 0 0.25 0.5 0.75 1 Mole Fraction, xi Figure 5.6 - The calculated natural logarithm of the activity coefficient (ln(γi)) of each pharmaceutical agent (i.e, 17-AAG (green solid line), PTX (blue dashed line), and RAP(red dotted line)) in water through the UNIFAC method as a function of the mole fraction. The spinodal limit is shown (solid black line) demarking where spinodal decomposition occurs. 117 Table 5.1 - Model-fitted parameters using the experimental data reported by Shin, et al. for in vitro release of pharmaceutical agents from using eq (5.13). Agents PTX 17-AAG RAP PTX 17-AAG RAP 17-AAG RAP PTX RAP PTX 17-AAG a Concentration Ratio Characteristic time (hr)b C∞/(kIII kII kI Coo) a NA NAa 0.056 2.172 0.109 12.28 0.208 7.230 0.060 2.512 0.071 12.61 0.060 2.601 0.036 16.26 0.072 10.08 2.10.10-9 18.23 5.78.10-4 12.17 0.059 3.508 Half-life (hr)c NAa 1.505 8.511 5.011 1.741 8.743 1.803 11.27 6.984 12.63 16.07 2.432 Not applicable due to excessive drug precipitation and no data available for fitting. Defined in eq (5.9). c Defined in eq (5.15). b 118 Table 5.2 - Selected UNIFAC subgroups used in the activity coefficient calculation. Main Group 1 Subgroup Number 1 CH3 6 1 6 0 0 0 1 2 CH2 2 2 12 95 0 0 1 3 CH 3 10 15 0 1 0 1 4 C 3 1 2 0 0 0 15 3 ACH 16 3 0 0 0 0 19 4 ACCH2 0 1 0 0 0 0 34 5 OH 3 1 3 1 1 0 36 7 H2O 0 0 0 0 0 1 42 9 CH3CO 2 2 1 0 23 0 33 9 CH2CO 4 2 4 0 0 0 59 13 CH3O 0 2 3 0 0 0 60 13 CH2O 0 0 0 95 0 0 61 13 CHO 5 1 2 0 22 0 72 15 CH2NH 1 1 1 0 0 0 Subgroup PTX 17-AAG RAP PEG PLA Water 119 Table 5.3 - Experimental values reported by Shin, et al. for in vitro release of agents compared to the fitted and calculated parameters. Agent(s) PTX Fitted UNIFAC Model Concentration Experimental Ratio Half-life LogPb Characteristic Characteristic Half-life (hr) Time (hr) Time (hr) (hr) C∞/(kIII kII kI Coo) a a a a NA NA NA 3.0 NA NAa 17-AAG 0.056 1.32 1.505 1.3 1.55 2.16.10-2 RAP 0.109 8.52 8.511 5.8 8.77 3.48 PTX 17-AAG 0.208 0.060 5.01 1.74 5.011 1.741 3.0 1.3 5.16 1.79 2.65 2.32.10-2 RAP 0.071 8.73 8.743 5.8 9.01 3.42 17-AAG 0.060 1.80 1.803 1.3 1.86 2.17.10-2 RAP PTX 0.036 0.072 10.05 6.00 11.27 6.984 5.8 3.0 11.6 7.20 3.19 2.73 2.10.10-9 5.78.10-4 0.059 13.93 9.20 2.52 12.63 16.07 2.432 5.8 3.0 1.3 13.0 8.70 2.51 3.17 2.72 2.28.10-2 RAP PTX 17-AAG 120 5.6 References (1) Lipinski, C. A. Journal of Pharmacological and Toxicological Methods 2000, 44, 235. (2) Aliabadi, H. M.; Lavasanifar, A. Expert Opinion on Drug Delivery 2006, 3, 139. (3) Shin, H.-C.; Cho, H.; Lai, T. C.; Kozak, K. R.; Kolesar, J. M.; Kwon, G. S. Journal of Controlled Release 2012, 163, 93. (4) Grant, S. Journal of Clinical Investigation 2008, 118, 3003. (5) Roforth, M. M.; Tan, C. Anticancer Drugs 2008, 19, 681. (6) 2012, 10. Savjani, K. T.; Gajjar, A. K.; Savjani, J. K. ISRN Pharmaceutics 2012, (7) Shin, H.-C.; Alani, A. W. G.; Rao, D. A.; Rockich, N. C.; Kwon, G. S. Journal of Controlled Release 2009, 140, 294. (8) Lavasanifar, A.; Samuel, J.; Kwon, G. S. Advanced Drug Delivery Reviews 2002, 54, 169. (9) Adams, M. L.; Lavasanifar, A.; Kwon, G. S. Journal of Pharmaceutical Sciences 2003, 92, 1343. (10) Kwon, G. S.; Okano, T. Pharmaceutical Research 1999, 16, 597. (11) Kataoka, K.; Harada, A.; Nagasaki, Y. Advanced Drug Delivery Reviews 2001, 47, 113. (12) Kazunori, K.; Glenn S, K.; Masayuki, Y.; Teruo, O.; Yasuhisa, S. Journal of Controlled Release 1993, 24, 119. (13) 1998, 15. Kwon, G. S. Critical Reviews™ in Therapeutic Drug Carrier Systems (14) 20, 357. Kwon, G. S. Critical Reviews in Therapeutic Drug Carrier Systems 2003, (15) Torchilin, V. P. Journal of Controlled Release 2001, 73, 137. (16) Torchilin, V. P. Advanced Drug Delivery Reviews 2002, 54, 235. 121 (17) Arifin, D. Y.; Lee, L. Y.; Wang, C.-H. Advanced Drug Delivery Reviews 2006, 58, 1274. (18) Kwon, G. S.; Okano, T. Advanced Drug Delivery Reviews 1996, 21, 107. (19) Shin, H.-C.; Alani, A. W. G.; Cho, H.; Bae, Y.; Kolesar, J. M.; Kwon, G. S. Molecular Pharmaceutics 2011, 8, 1257. (20) Jordan, M. A.; Wilson, L. Nature Reviews. Cancer 2004, 4, 253. (21) Solit, D. B.; Chiosis, G. Drug Discovery Today 2008, 13, 38. (22) Dancey, J. E. Cancer Biology & Therapy 2006, 5, 1065. (23) 364, 328. (24) Siepmann, J.; Siepmann, F. International Journal of Pharmaceutics 2008, Fu, Y.; Kao, W. J. Expert Opinion on Drug Delivery 2010, 7, 429. (25) Deen, W. M. Analysis of Transport Phenomena; 2nd ed.; New York : Oxford University Press: New York, 2012. (26) Taylor, R.; Krishna, R. Multicomponent Mass Transfer; New York : Wiley: New York, 1993. (27) Smith, J. M.; Van Ness, H. C.; Abbott, M. M. Introduction to Chemical Engineering Thermodynamics; 5th ed.; New York : McGraw-Hill: New York, 1996. (28) Fredenslund, A.; Jones, R. L.; Prausnitz, J. M. AIChE Journal 1975, 21, 1086. (29) Seader, J. D.; Henley, E. J. Separation Process Principles; New York : Wiley: New York, 2006. (30) Kang, J. W.; Diky, V.; Chirico, R. D.; Magee, J. W.; Muzny, C. D.; Abdulagatov, I.; Kazakov, A. F.; Frenkel, M. Fluid Phase Equilibria 2011, 309, 68. (31) Wittig, R.; Lohmann, J.; Gmehling, J. Industrial and Engineering Chemistry Research 2003, 42, 183. (32) Russel, W. B.; Saville, D. A.; Schowalter, W. R. Colloidal Dispersions; Cambridge University Press, 1992. (33) Grijpma, D. W.; Penning, J. P.; Pennings, A. J. Colloid Polymer Science 1994, 272, 1068. 122 (34) Henton, D. E.; Gruber, P.; Lunt, J.; Randall, J. Natural Fibers, Biopolymers, and Biocomposites, Taylor & Francis, Boca Raton, FL 2005, 527. (35) Poon, W. C. K.; Andelman, D. Soft Condensed Matter Physics in Molecular and Cell Biology; New York : Taylor & Francis: New York, 2006. (36) Bergman, T. L.; Lavine, A. S.; Incropera, F. P.; DeWitt, D. P. Fundamentals of Heat and Mass Transfer; 7th ed.; Hoboken, NJ : John Wiley: Hoboken, NJ, 2011. (37) Cronin, M. Current Computer-aided Drug Design 2006, 2, 405. 123 CHAPTER 6 CONCLUSIONS AND FUTURE WORK The general purpose of the study was to use engineering techniques to answer fundamental questions concerning EoE development and treatment. These methods and techniques are common in traditional engineering applications (e.g., fluid dynamics, heat transfer, etc.) and are used to obtain insight into a physical system through mathematical application. By applying these engineering methods, insight is developed into the physical parameters and conditions contributing to a diseased state. First, an examination of the relationship between concentration and the curvature of the esophagus was performed. The transport of casein was modeled using an idealized esophageal geometry; casein concentration profile in the epithelium was observed to be significantly impacted by curvature of the esophagus. This result agrees with experimental observations and if esophageal curvature can be reduced, plaque and regions of elevated eosinophil concentrations may be minimized. A species conservation model was developed to represent a simplified pathogenesis model of EoE. These equations were examined using linear stability analysis, which concluded transport described by these equations is unconditionally stable. Thus, a parametric assessment of the conservation equations was performed. The parametric assessment determined the reactive interactions between biomarker components (e.g., 124 eotaxin-3) and eosinophils are highly significant. To better understand and control EoE, continued work should investigate the rates in which these components react with each other. If the reactions can be slowed or even disabled, EoE could be conceivably eliminated. In addition, a model was developed to describe the delivery of polymer micelles loaded with hydrophobic drugs. This model can be used to predict drug delivery time scales for a wide range of pharmaceutics, which can be used to treat EoE. Using the model, drug release can be tailored such that the micelles can be applied directly to the esophagus for a timed release to the affected areas, leading to more efficient and effective treatment of EoE. The use of numerical and analytical analyses is an important foundation to study the development of EoE. This foundation was used here, and to build upon this effort, a more detailed model of the esophagus may be developed. Ideally, this model would account for fluid-epithelial interaction and absorption of antigen into the esophagus. By accounting for the highly variable geometry of the esophagus, and the flow of mucin down the esophagus, a more detailed map of eosinophil density may be possible. That type of information could be directly related to patient measurements, leading to an absorption-diffusion-reaction model of EoE. The implications of this approach can give an individual a model that can be used to best treat EoE for a wide range of patient scenarios and conditions, which can be examined to determine the precise cause of a patient's individual response to an EoE response. To identify long-term trends of the disease and under what conditions one could be diagnosed with EoE, additional effort should be undertaken to carefully examine the 125 equations that define the pathogenesis. The model presented in this work was simplified to yield general trends. This resulted in unconditionally stable conditions, but it reveals important factors governing EoE development, in particular, the interactions and rates of reaction between activation components and eosinophils that govern the mass transfer within the epithelium. Naturally, different approaches and models may result in a system of equations that will clarify the applicable physical parameters that could lead to a diseased state. Collecting EoE data will be important for validating computational any models and appreciating transport phenomenon within the esophagus. Currently, the data available are not well suited from modeling to identify a coherent pathogenic model of EoE. It is believed that through the development of biophysics-based models and the application to biological systems, better predictions of diseases, like EoE, can be determined and methods of treatment can be more readily defined. This effort is opening a new avenue in the assessment and treatment of the EoE disease state and points to methods for better understanding of disease pathogenesis.