Computers and Chemical Engineering 28 (2004) 2407–2419
Neural network-based prediction and optimization of estradiol release from ethylene–vinyl acetate membranes Laurent Simon∗ , Maria Fernandes a
Otto H. York Department of Chemical Engineering, New Jersey Institute of Technology, Newark, NJ 07102, USA Received 11 December 2003; received in revised form 24 May 2004; accepted 9 June 2004 Available online 7 July 2004
Abstract Drug-delivery systems, with predictable delivery rates, were designed using an artificial neural network-based optimization algorithm. A two-chamber diffusion cell was used to study the permeation of estradiol through ethylene–vinyl acetate copolymer membrane. The explanatory variables were the vinyl acetate (VA) content of the membrane, poly(ethylene glycol) (PEG)—solvent-composition, and membrane thickness. After deriving a neural network model to predict estradiol delivery rates as a function of these input variables, a constrained optimization procedure was applied to estimate the membrane/vehicle properties necessary to achieve a prescribed dosage. The results compared adequately well with experimental data with 71% of the data agreeing within one standard deviation. Input sensitivity analysis showed that at specific VA levels, drug delivery was more sensitive to changes in PEG compositions. The non-uniqueness of the inversion method and the accuracy of the procedure were investigated using neural network-based two-dimensional contour plots. The methodology proposed could be used to design customized polymer-based drug-delivery systems that meet specific end-user requirements. © 2004 Elsevier Ltd. All rights reserved. Keywords: Neural networks; Constrained optimizations; Drug delivery rate
1. Introduction A major challenge in the pharmaceutical industry is that current drug-delivery methods, such as tablets, injections and sprays, are inefficient in administering long-term consistent drug release to a target site. Certain drugs lose effectiveness over a period of time, requiring multiple injections to maintain the necessary therapeutic level. The plasma drug concentration may fall below the minimum effective concentration (MEC) or rise above the minimum toxic concentration (MTC) in a dosing cycle. This work investigates a new mathematical formulation that ensures accurate, cost-effective polymer-based delivery of drugs to target sites. Although optimal experimental designs (e.g. sequential simplex delivery method or response surface experiments) are helpful tools in the development of drug-delivery devices, they are costly and provide only limited predictive capabilities. Therefore, it is imperative to
∗
Corresponding author. Tel.: +1 973 596 5263; fax: +1 973 596 8436. E-mail address:
[email protected] (L. Simon).
0098-1354/$ – see front matter © 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.compchemeng.2004.06.002
investigate the use and benefits of new mathematical frameworks to enhance understanding of controlled drug release. One of the outcomes of the Genome project is the influx of new drugs being developed and evaluated to treat human diseases. Polymer-based delivery devices are currently being studied for delivery of these new drugs. Such systems can be engineered to release the medicinal agent to the target site in a controlled fashion. Recent advances in polymer science also contribute to the continued use of polymer-based drug delivery systems. These developments include nanoporous materials with photosensitive pores (Flanagan, 2003) and temperature-sensitive polymer membranes (Grassi, Yuk, & Cho, 1999). Innovative material, currently under development, has the potential to pave the way for sophisticated drug/polymer systems that have the ability to respond to changes in the environment (Bettini, Colombo, & Peppas, 1995; Siegel & Pitt, 1995). Biodegradable polymers, based on polyactic acid and polyanhydrides, are employed for the long time-release of drugs from implanted patches and other drug-delivery systems (Kocova El Arini & Leuenberger, 1995). Hydrophilic water-soluble polymers are preferred for oral administration of the drug. Koizumi, Ritthidej,
2408
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
and Phaechamud (2001) proposed a mechanistic model for drug release from chitosan coated tablets. In their applications, spherical core tablets, containing a model drug propranolol hydrochloride, were investigated. The rate of drug release was controlled by diffusion through the inner space and the polymer coating of the pills. In view of this rapid development and growing complexity of drug delivery systems, a new methodology is proposed to integrate experimental data, prior system knowledge, and artificial neural networks (ANNs), to assist in the development of novel drug-release device technologies. Controlled drug-release-based technology has been studied extensively (Langer, 1998) and exploited by the pharmaceutical industry with major breakthroughs. However, even when the release of drugs from a polymer matrix occurs mainly by diffusion, it is still difficult to achieve the desired release profile. One reason is that the release rate is a function of several factors, such as the characteristics of the excipients added to the polymeric drug device to control the drug release pattern (Shah & Zatz, 1992) and the conditions of the external environment (e.g. pH and temperature). In some cases, a mathematical model can be derived using Fick’s second law; however, the appropriate drug formulation is not explicit in such models. Current estimation procedures are iterative and time-consuming since the effect of pertinent contributing factors on the drug delivery rate is often evaluated in a qualitative manner, which precludes the use of any mathematical theory of optimization. As a result, numerous experiments are necessary to approximate the drug formulation required to achieve the desired release profile. Further refinement of these procedures is often based on trial and error approaches or experimental optimization schemes (e.g. sequential simplex methods). These techniques are cost-prohibitive and often result in sub-optimum drug-delivery devices. Berkland, King, Kim, and Pack (2002) agree that the ability to tailor drug release profiles is difficult. In their investigation, using poly(d,l-lactide-co-glycolide) (PLG) microspheres, they identified (by testing different systems) appropriate mixtures of uniform microspheres that provided constant release of rhodamine and piroxicam for 8 and 14 days, respectively (Berkland et al., 2002). This work is based on the use of ANN models to predict the drug-release rates as a function of key membrane and solvent properties. The reason to use an artificial neural network-based model is that the transport of pharmaceutics through membranes is complex and influenced by a host of factors such as additives, polymer functionality, porosity, film casting conditions, temperature, and even the type of membrane. Since the relationship between these factors and drug delivery rate is not well defined, it is customary to write the diffusion coefficient as a nonlinear function of the systems property (e.g. polymer compositions, excipients, etc.). For instance, Siepmann, Lecomte, and Bodmeier (1999) investigated the effect of the composition of diffusion-controlled release devices (type and amount
of plasticizer, type of polymer) on the drug diffusivity, and ultimately on the release kinetics. Their goal was to achieve desired release profiles by manipulating the drug diffusion coefficient. Since only a few researchers have studied the effect of plasticizer on the drug diffusion in a quantitative manner, Siepmann, Lecomte, and Bodmeier (1999) used exponential functions to describe the relationships between the diffusion coefficient of theophylline and the plasticizer level of various polymer–plasticizer systems. Since the synergistic effect of the properties studied in the current work on drug-release rate has not been described from first-principle modeling, an ANN function was used. The ANN function was then inverted to estimate the necessary solvent/membrane properties required to achieve a desired delivery rate. The use of a new Mathematica® (Wolfram Research, Inc.) add-on package is also described. Since the package provides symbolic representation of ANN models, contour plots can easily be drawn based on the generated model instead of polynomial approximations, routinely used for this type of analysis. A new method for ANN-based input sensitivity analysis is also provided. This paper is divided into three sections: (1) A Materials and methods section underlying the experimental procedure and concepts involved; (2) A Results and discussions section addressing: (a) Release kinetics due to the effects of vinyl acetate (VA) contents of ethylene–vinyl acetate (EVA) copolymer membrane, membrane thickness, and poly(ethylene glycol) (PEG) (solvent) composition; (b) Development of an ANN model to predict release rates as a function of the variables mentioned above; (c) Design of drug-delivery devices to meet specific end-user requirements based on an inversion of the ANN model; and (d) Local input sensitivity analysis; (3) Conclusion.
2. Materials and methods 2.1. Materials Estradiol and PEG 400 were purchased from Sigma– Aldrich (St. Louis, MO). The 3M Company generously donated the ethylene–vinyl acetate membranes. 2.2. Determination of estradiol saturation concentration A drug suspension in different PEG 400 solution concentrations was placed in the donor cell to guarantee that the drug concentration, following a zero-order kinetics, would remain constant in the solvent. The drug equilibrium solubilities, in different concentrations of PEG 400 solutions and at a temperature of 37 ◦ C, were determined by adding an excess amount of estradiol to 3 ml of PEG 400-saline solutions: 20, 30 and 40%. Samples were produced in duplicates and placed in a water bath at 37 ◦ C for 3 days with occasional mixing. The saturated solutions were centrifuged and the supernatant liquid collected. The concentration was
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
determined by measuring the absorbance at 290 nm using a UV spectrophotometer. 2.3. Drug release through EVA membranes A two-chamber diffusion cell was used to study the permeation of estradiol through EVA membranes. Each half-cell contained 3 ml of solution with an effective diffusional area of 0.785 cm2 . The EVA membrane was clamped between the donor and receiver cells and the system was kept at a temperature of 37 ◦ C by placing the apparatus in a temperature-controlled water bath. A pump was used to circulate water through the device. Saturated drug-in-solvent solutions, to be used in the experiments, were prepared, in advance, to allow some time for equilibration. Three milliliters of the equilibrated drug/PEG solution were placed in the donor cell compartment. The receptor cell compartment was filled with 3 ml of the same concentration of PEG solution. At specific time intervals, the receiver cell content was removed for analysis and replaced with the same volume of fresh PEG solution. The concentrations of the collected samples were determined by measuring the sample’s absorbance in a spectrophotometer at 290 nm. 2.4. Computational methods 2.4.1. Drug permeation through the EVA membrane The following equation was used to calculate the cumulative amount of drug released per unit area (Shin & Lee, 2002): Q = P(Cd − Cr )t,
(1)
where P is the permeability coefficient, t the time, Cd and Cr are the drug concentrations in the donor and receiver cells, respectively. Since the drug in the donor compartment is above the drug saturation concentration (Ceq ) and the drug concentration in the receptor solution is maintained under sink condition (i.e. Cr Ceq ), Eq. (1) becomes: Q = PCeq t
(2)
2409
ANN models offer a better filtering capacity and generally perform better than empirical models. Artificial neural network models have emerged as an accepted tool for modeling nonlinear processes. Different architectures of these models have been implemented in the literature (Eikens, Karim, & Simon, 2001; Gao & Loney, 2001; Simon & Karim, 2001; Simon, Karim, & Schreiweis, 1998). Artificial neural networks consist of a set of weights, w, a set of basis functions, f(·), and a set of parameters, p (biases, centers, etc.). The general equation can be written as: yˆ = fNN (x, w, p)
(4)
where yˆ denotes the NN estimation (i.e., the permeation rate Q/t), x is the input vector (i.e., VA, Th, and PEG) and fNN is the network function. Developing a supervised ANN model involves three phases: A training phase in which parameters w and p are adjusted in order to minimize a performance criterion, such as the squared error between the experimental or measured output and the output predicted by the network. This phase requires a training set consisting of known input and output patterns. Examples, representative of the problem, are propagated through a series of layers. The network learns the underlying relationship between input and output variables. A recall phase, in which the network is exposed to the input patterns seen in the training phase, is then used. Adjustments, such as the number of the layers, can be made to enhance the robustness and reliability of the system. The weights and biases are computed. This step is followed by a testing phase, in which new input patterns (not part of the training set) are fed through the network with fixed weights and biases. The performance of the system is obtained by comparing experimental and predicted network outputs. Fig. 1 shows the structure of a neural network with one hidden layer and two nodes (h1 and h2 ). The choice of the number of nodes is explained in the Results and Discussion section. A “bias” unit is used to help the training phase. The bias term is equivalent to a y-intercept in a linear regression model. The components of w are w2h1 , w2h2 , w3h1 , w3h2 , w4h1 , and w4h2 . The components of p are p1h1 , p1h2 , and p2o . Optimum values of the weights (wopt )
The drug permeation rate per unit area is then: Q = PCeq t
(3)
Its value was estimated by plotting Cd cum V/A as a function of time. The diffusion area is denoted by A; V is the receiver cell volume (3 ml) and Cd cum the cumulative receptor cell concentration. 2.4.2. Artificial neural network prediction of the permeation rate A feed-forward ANN model was developed to predict the drug permeation rate as a function of vinyl acetate contents (VA) of the ethylene–vinyl acetate (EVA) copolymer membrane, membrane thickness (Th), and PEG 400 compositions (PEG). When dealing with noisy and incomplete data,
1 p1h2
h1
p1h1 w 2h1
Σσ
w h1o
VA
w 2h2 w 3h1
h2
PEGw 3h2
Σσ
w 4h1 Th
w h2o
1 p2o
Σ
[Q/t ]hat
w 4h2
Fig. 1. Neural network architecture using one hidden layer and two nodes. Weight wihj connects the ith input to the jth hidden neuron; weight whio connects output from the ith hidden layer to the output; p1hj connects the unit bias to the jth hidden node; p2o connects the unit bias to the output node. The weighted inputs are summed up and passed through a nonlinear activation function sigma.
2410
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
and biases (popt ) are obtained and used when solving the inverse problem (Section 2.4.4). 2.4.3. Input sensitivity analysis Neural network-based sensitivity analysis is usually implemented to assess the causal importance of individual input variables. Different applied methods have been reported in the literature (Caldwell, 1995; Caldwell, 1996; Howes & Crook, 1999; Karim, Hodge, & Simon, 2003; Simon & Karim, 2002). The technique used by Simon and Karim (2002) can be summarized as follows: (1) find an optimum ANN model to represent the input–output mapping; (2) perturb the variable of interest of the input vector while keeping the remaining components at their mean values; and (3) compute the output of the network. This method was based on previous work by Principe, Euliano, and Lefebvre (2000), who used the following expression to calculate the sensitivity of input k on the network output: p O ¯ ip )2 p=1 i=1 (yip − y Sk = (5) σk2 where yip is the ith output of the optimum network for the pth pattern before input perturbation (the other input variables are kept at their main values), y¯ ip the ith output of the optimum network for the pth pattern after input pertubation, O stands for the number of network output (O = 1), P is the number of patterns, and σk2 represents the variance of the input perturbation. Although this method yields successful results, Eq. (5) can be viewed as a discretized form of the equation for input sensitivity (i.e. squared sensitivity) and is useful for cases in which the input–output mapping function is too complex. One of the advantages of using Eq. (5) is that irrelevant input variables (variables that do not affect the outcome of a process) can be eliminated (Simon & Karim, 2002). As a result, data collection is reduced and a better insight into the underlying relationships between input and output variables can be gained. However, in processes, such as drug delivery, researchers engaged in data analysis and interpretation may be interested in knowing whether a small change in an input variable results in a relatively large change in the outcome (e.g. drug-release rate) for values of the remaining input variables other than the average quantities. The approach used in the present work makes use of the definition of input sensitivity and can help researchers to take into account the other input values. The effect of an input χ (component of a input vector χ) ¯ on the output ζ such that ζ = f(χ) ¯ is the partial derivative of the function f with respect to the independent variable of interest: ∂f S= (6) ∂χ χ0 where χ0 is the remaining components. The methodology proposed in this work uses Eq. (6) to perform the neural network-based input sensitivity analysis. As a result, an
ANN model with a minimum number of nodes and hidden layers, still capable of explaining changes in estradiol permeation based on vinyl acetate (VA) content of the membrane, poly(ethylene glycol) (PEG)—solvent-composition, and membrane thickness will be sought. Contrary to previous work, an explicit expression will be given for the network model to show such dependence. 2.4.4. Solution of the inverse problem The inverse problem consists of finding the VA content of the membrane, the membrane thickness, and the percentage of PEG (i.e. the input vector xi ) necessary to achieve a prescribed drug permeation rate (yset i ). A constrained optimization scheme, implemented in Mathematica® (Wolfram Research, Inc.), was used to solve the problem. “NMinimize” is a function that provides four different methods for global optimization. The “Random search” option was selected. The optimization problem consists of finding the optimal input vector parameters xi min such that the objective function F defined as: F(xi ) = [yset i − fNN (xi , wopt , popt )]2
(7)
subject to xi ∈ Ψ
(8)
is minimized: xi min = Min F(xi )
(9)
where xi min is a member of an admissible solution set Ψ . The search space Ψ is made up of the data used in the training and testing sets. The constraint defined by Eq. (8) forces the algorithm to explore only solutions that can be tested immediately under the experimental conditions studied. The weights wopt and biases popt correspond to the optimum network parameters and remain fixed when solving the inverse problem. Because of the intrinsic non-uniqueness of inversion methods, two-dimensional contour plots were employed to analyze the characteristics of the optimized ANN model and to estimate the accuracy of the procedure outlined above. These contours are often based on polynomial equations from a multivariate regression analysis. However, since a parsimonious network structure is derived, it is possible to use the network model directly to draw the contours. The two approaches to design customized polymer-based drug-delivery systems that meet specific end-user requirements were compared.
3. Results and discussion 3.1. Release kinetics Twenty-four samples were generated by varying PEG 400 volume fractions (20, 30, 40, and 50%), VA contents of the
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
2411
14.00
12.00 Drug release rate (ug/cm^2.h)
Membrane thickness: 0.00508 cm
PEG = 20% PEG = 30% PEG = 40% PEG = 50% PEG = 50% (NN) PEG = 40% (NN) PEG = 30% (NN) PEG = 20% (NN)
10.00
8.00
6.00
4.00
2.00
0.00 0
5
10
15
20
Vinyl acetate content (%)
Fig. 2. Experimental (—) and predicted (- - -) drug-release rate per area at a membrane thickness of 0.00508 cm. PEG compositions are 20, 30, 40, 50%; VA contents are 9 and 19%.
membrane (9 and 19%), and membrane thickness (0.00508, 0.00762, 0.01016 cm). The release rates, corresponding to the various combinations, are shown in Figs. 2–6. These results suggest that the permeation rate rose with an increase in VA content, and PEG composition and a decrease in membrane thickness. These findings agree with research led by Shin and Byun (1996). Through their investigations of the permeability of ethinylestradiol through EVA copolymer membranes, they reported that an increase in VA comonomer content increased the drug-release rate and permeability coefficient (Shin & Byun, 1996). This result may be due to an increase in the diffusivity of the drug. Their conclusion was reached after considering previous works which had proved that the introduction of VA comonomers
to a highly crystalline polyethylene decreased the crystallinity of the system, which, in turn, increased the diffusivity of polymers (see Shin & Byun, 1996 for a complete reference). It should be noted that if the diffusion coefficients are of the same order of magnitude for the EVA copolymer membranes studied, the increase in the permeation rate may be the result of an increase in the partition of the drug into the EVA membrane (Shin & Byun, 1996). In this work, estradiol equilibrium solubility was found to increase with the addition of PEG 400: 51, 180, 600, 1000 g/cm3 for 20, 30, 40, and 50% PEG, respectively, which resulted in an increase of the drug-release rate. This trend can be predicted by the analysis of Eq. (3).
7.00 PEG = PEG = PEG = PEG = PEG = PEG = PEG = PEG =
Drug release rate (ug/cm^2.h)
6.00
5.00
4.00
20% 30% 40% 50% 50% 40% 30% 20%
Membrane thickness: 0.00762 cm (NN) (NN) (NN) (NN)
3.00
2.00
1.00
0.00 0
5
10
15
20
Vinyl acetate content (%)
Fig. 3. Experimental (—) and predicted (- - -) drug-release rate per area at a membrane thickness of 0.00762 cm. PEG compositions are 20, 30, 40, 50%; VA contents are 9 and 19%.
2412
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419 4.50
Drug release rate (ug/cm^2.h)
Membrane thickness: 0.01016 cm
PEG = 20% PEG = 30% PEG = 40% PEG = 50% PEG = 50% (NN) PEG = 40% (NN) PEG = 30% (NN) PEG = 20% (NN)
4.00 3.50 3.00 2.50 2.00 1.50 1.00 0.50 0.00 0
5
10
15
20
Vinyl acetate content (%)
Fig. 4. Experimental (—) and predicted (- - -) drug-release rate per area at a membrane thickness of 0.01016 cm. PEG compositions are 20, 30, 40, 50%; VA contents are 9 and 19%.
Shin and Byun noted that the permeation rate was inversely proportional to the membrane thickness (Shin & Byun, 1996). From Figs. 5 and 6, it is evident that such a trend exists (although not quite linear for all PEG and VA compositions) and can be used to control drug permeation. It is known that an increase in membrane thickness implies an increase in the distance the drug has to travel through the membrane, thereby reducing the diffusion of the drug.
3.2. Artificial neural network model The data generated for the kinetics study were used for the training and testing of the ANN. These sets were pre-processed by initially randomizing the rows of data corresponding to different input and output vectors. This step is important because similar input data, grouped together, may considerably affect the training phase given
7.00
6.00
Drug release rate (ug/cm^2.h)
Vinyl acetate: 9%
PEG = 20% PEG = 30% PEG = 40% PEG = 50% PEG = 50% (NN) PEG = 40% (NN) PEG = 30% (NN) PEG = 20% (NN)
5.00
4.00
3.00
2.00
1.00
0.00 0
0.002
0.004
0.006
0.008
0.01
0.012
Menbrane thickness (cm)
Fig. 5. Experimental (—) and predicted (- - -) drug-release rate per area at a VA content of 9%. PEG compositions are 20, 30, 40, 50%; membrane thicknesses are 0.00508, 0, 0.00762, 0.01016 cm.
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
2413
14.00
12.00
Drug release rate (ug/cm^2.h)
Vinyl acetate: 19%
PEG = 20% PEG = 30% PEG = 40% PEG = 50% PEG = 50% (NN) PEG = 40% (NN) PEG = 30% (NN) PEG = 20% (NN)
10.00
8.00
6.00
4.00
2.00
0.00 0
0.002
0.004
0.006
0.008
0.01
0.012
Membrane thickness (cm)
Fig. 6. Experimental (—) and predicted (- - -) drug-release rate per area at a VA content of 19%. PEG compositions are 20, 30, 40, 50%; membrane thicknesses are 0.00508, 0, 0.00762, 0.01016 cm.
that they are unlikely to be fully representative of the range of possible values. The input and output patterns were scaled to lie in the interval [0.1, 0.9] for efficient processing of the network because of the sigmoid activation function, later used in the training phase. Seventy-five percent of the data were used for tutorial purposes (training set) and the remaining (testing set) were used to test the performance of the model learned during the training phase. When dealing with a large data set, several algorithms, such as “optimal brain damage” (Cun, Denker, & Solla, 1990) and “optimal brain surgeon” (Hassibi & Stork, 1993) can be employed to decide on the optimum number of nodes to use in the hidden layer. However, with a relatively small data set, the choices become limited. Since it is considered good statistical practice that the number of adjustable parameters be less than the number of data points, two hidden nodes were chosen. A feed-forward network, consisting of one hidden layer, two neurons and a sigmoid-type activation function, was built with three inputs and one output (often denoted as a 3:2:1 network). The maximum number of iterations to train the network was fixed at thirty. The Neural Networks—“Train and analyze neural networks to fit your data” (Jonas Sjöberg, 2002)—add-on Mathematica® package from Wolfram Research, Inc., was used in this study. This work is the first attempt to take advantage of the capability of this software in engineering controlled drug-release devices. Artificial neural network models are usually considered implicit black-box predictors. Subramanian, Yajnik, and
Murthy (2003) investigated the use of these models as alternatives to multiple regression analysis in optimizing formulation parameters. They observed that ANN models were more accurate than polynomial equations (commonly used to generate contour plots) in predicting percentage drug entrapment (PDE). These network models can also handle more input variables and process a larger number of experimental data (Subramanian et al., 2003). However, in their work, data analysis such as relationship between independent variables and the PDE, and second-order main effects of both drug/lipid ratio and volume of hydration, was performed using a reduced polynomial equation. This work, by taking advantage of Mathematica® ’s mixed numeric–symbolic environment and the Neural Network add-on package, extends the usefulness of ANN models beyond that of black-box predictors. An explicit mathematical representation of the input–output mapping is first derived and later used for developing a new method for implementing neural network-based input sensitivity analysis, and solving drug formulation problems using contour plots and function minimization. The steps involved in building the ANN models are: • Data pre-processing. The data (randomized) was initially read from a text file using the Mathematica® command “ReadList”. They were then normalized from 0.1 to 0.9. • Seventy-five percent of the data was used for training, and the rest for testing. • A feed-forward network with two neurons in the hidden layer was initialized using the Mathematica® command: “InitializeFeedForwardNet”.
2414
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
Model
Model
0.8
0.5
0.4 0.6 0.3 0.4 0.2 0.2
0.1
0.1
0.2
0.3
0.4
0.5
Data 0.2
0.4
0.6
0.8
Data
Fig. 7. Neural network predicted output (ˆq/t) vs. the normalized output (q/t). The clustering of points around line y = x indicates a good fit. The training set was used.
Fig. 8. Neural network predicted output (ˆq/t) vs. the normalized output (q/t). The clustering of points around line y = x indicates a good fit. The testing set was used.
• The network was trained using the Mathematica® command “NeuralFit”. A symbolic function was generated at the end of the training procedure. • The network was tested by using the symbolic function with the testing set as argument. • If the performance of the testing set was not satisfactory, the training phase was repeated.
function is:
Learning took place after only 14 iterations, at which point the minimum squared error no longer diminished—therefore, preventing overtraining. The optimized neural network
qˆ 1386.56 = 0.10 + 23.04−2.06peg+1.58th−15.34va t 1+e 0.58 + 1 + e1.55−4.38peg+1.98th−16.28va
(10)
where peg, th, va, and qˆ /t are the normalized PEG volume faction, membrane thickness, vinyl acetate content, and drug permeation rate per unit area, respectively. The results of the training phase are summarized in Fig. 7. The clustering of points around the regression line y = x shows that the
Table 1 Estimation of the VA content, PEG composition and membrane thickness using a specified permeation rate per unit area Sample
VA (%)
TH (cm)
PEG (%)
Q/t (g/cm2 h)
S.D. (g/cm2 h)
VAopt (%)
THopt (cm)
PEGopt (%)
(Q/t)opt (g/cm2 h)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
9 9 19 19 9 19 19 19 9 19 19 9 9 19 9 9 9 9
0.00508 0.00508 0.00508 0.00508 0.00762 0.01016 0.00762 0.00762 0.01016 0.00762 0.00762 0.00762 0.00762 0.01016 0.00762 0.00508 0.00508 0.01016
30 40 20 40 40 20 30 20 50 40 50 20 30 40 50 20 50 40
1.11 3.18 2.04 6.47 3.19 1.03 2.76 1.77 2.01 3.81 5.79 0.72 1.06 2.00 4.29 1.00 5.54 1.44
0.104 0.236 0.205 0.299 0.207 0.017 0.247 0.188 0.120 0.195 0.412 0.095 0.097 0.178 0.543 0.150 0.592 0.158
9 9 9 19 9 19 9 9 9 19 19 9 19 9 9 9 9 19
0.00762 0.00508 0.00762 0.00508 0.00508 0.01016 0.01016 0.00508 0.00762 0.00508 0.00762 0.00762 0.01016 0.00762 0.00762 0.00508 0.00508 0.01016
30 40 40 40 40 20 50 30 40 30 50 20 20 40 50 20 50 30
1.06 3.18 3.19 6.47 3.18 1.03 2.01 1.11 3.19 4.18 5.79 0.72 1.03 3.19 4.29 1.00 5.54 1.08
The training set was used.
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
2415
Table 2 Estimation of the VA content, PEG composition and membrane thickness using a specified permeation rate per unit area Sample
VA (%)
TH (cm)
PEG (%)
Q/t (g/cm2 h)
S.D. (g/cm2 h)
VAopt (%)
THopt (cm)
PEGopt (%)
19 20 21 22 23 24
19 9 19 9 19 19
0.00508 0.01016 0.01016 0.01016 0.01016 0.00508
30 30 50 20 30 50
4.18 0.99 3.91 0.56 1.08 11.48
0.155 0.102 0.257 0.026 0.127 0.479
9 9 19 9 19 19
0.00762 0.00508 0.00508 0.01016 0.01016 0.00508
50 20 30 20 20 50
(Q/t)opt (g/cm2 h) 4.29 1.00 4.18 0.56 1.03 11.48
The testing set was used.
model is accurate in predicting the system output. The network performed equally well when unlabeled testing input data (not part of the training set) were fed into the network (Fig. 8). These results are also shown in Figs. 2–6, in which the network predictions and experimental data are plotted. Eq. (10) will be used in the rest of the work. 3.3. Inversion of the artificial neural network model Although artificial neural networks have been implemented to predict drug-release profiles (Lim, Quek, & Peh, 2003), additional research is needed to investigate the
feasibility of ANNs in designing effective pharmaceutical formulations based on user-defined drug-release rate targets. Recent contributions include the work of Takahara, Takayama, and Nagai (1997) and Takayama, Fujikawa, Obata, and Morishita (2003). Using the methodology outlined in Section 2.4.4, the denormalized results summarized in Table 1 were obtained for the training data. The first column “Sample” numbers the data points; the second, third, and fourth columns denote the input set, and the fifth column represents the experimental output (Q/t): permeation rate per diffusion area. The sixth column is the standard deviation of the output.
Fig. 9. PEG sensitivity analysis at different vinyl acetate contents and membrane thicknesses.
2416
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
Fig. 10. Membrane thickness sensitivity analysis at different vinyl acetate contents and PEG compositions.
Using the scaled target output (yset i = qˆ /t), Eq. (10), and the optimization procedure outlined by Eqs. (7)-(9), an optimum scaled input vector was calculated xi min = [vaopt , thopt , pegopt ]. Columns seven through nine represent the denormalized optimum input vector [VAopt , THopt , PEGopt ]. For comparison, the last column, which corresponds to the actual output of [VAopt , THopt , PEGopt ], is added. The bold values correspond to cases in which the methodology predicted the experimental values (taking into account standard deviations). Eleven out of 18 such cases (61%) were observed. The underlined values correspond to cases where the output estimated by the inverse problem corresponds exactly to the experimental value. These results can be explained by the fact that several combinations of the inputs gave a similar output. A more complex network architecture (additional nodes) would probably predict 100% of the input vector. This would cause the network to memorize the input and output patterns in the training sets, but testing would prove unsuccessful. This point becomes clearer when we examine Table 2, corresponding to the testing set. One-hundred percent of the input data were predicted correctly (taking into account standard deviations).
This method gives satisfactory results, although the inverse problem performed better on the testing than on the training sets. This example illustrates a very important point. Although overparametrization and overtraining would eventually lead to a more accurate prediction of the training data, the network may lose the ability to predict unseen patterns. Tutorial (training) data are important insofar as they help to achieve a model that best describes the underlying process. Such an empirical model can then be used, in various applications, to perform system analysis, design and performance evaluations. With a large number of hidden units, it is possible to model spurious relationships, which do not reflect true dependencies among the variables. This will lead to meaningless results and erroneous process decisions. The methodology proposed could be used to design accurate polymer-based drug-delivery systems to meet specific end-user requirements. The predictive capability of the methodology allows for a drug-delivery device performance assessment without having recourse to a large number of experiments. Such a black-box modeling approach, although it does not substitute for a first-principle physical
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
knowledge of the system, can provide a starting point for future investigation and discussions. Contour plots were drawn to solve the inverse problem and evaluate the results of the procedure outlined above. The command “contour (Z)” from Matlab (The MathWorks, Natick, MA) was used to draw a contour plot of matrix Z with contour lines and their values chosen automatically. After using the functions “linspace” and “meshgrid” to create a set of points at values that spanned the data range, Eq. (10) (denormalized form) was used at 9 and 19% vinyl acetate content. The result is shown in Fig. 9. The contour plots agree, very well, with the best solutions of the minimization procedure. Consider, for instance, the first two points of Table 1. The neural network prediction for the input vectors:[9, 0.00762, 30], and [9, 0.00508, 40], representing the percentage of vinyl acetate, membrane thickness, and PEG compositions, are 1.06 and 3.18 g/cm2 , respectively. From the “%VA = 9” plot (Fig. 9), the non-uniqueness of these two problems is evident. A release rate of 1.1 g/cm2 can be obtained by using several combinations of PEG concentrations and membrane thicknesses, including the ones given by the minimization method. The same is true for
2417
a release rate of 3.2 mg/cm2 . Both the contour plots and the minimization procedure can assist in the design of customized polymer-based drug-delivery systems that meet specific end-user requirements. However, for pharmaceutical applications, the contour plots provide the clinician with several viable alternatives that can be considered. Therefore, least-squares optimization routines are not necessary, unless one is dealing with a high dimensional input feature space. 3.4. Local input sensitivity analysis Previously published methods calculate average sensitivity analysis based on the ANN output. The proposed methodology used the ANN model directly and allowed calculation of input sensitivities when the other inputs were kept at levels other than their mean values. To compare, for instance, how sensitive the drug-release rate was to changes in PEG and TH, derivatives of Eq. (10) with respect to the variables: PEG, and TH, respectively, were calculated. The sensitivity function, for a particular input (Sin ), was evaluated for four different combinations of the other input variables. Four functions for PEG sensitivity (SPEG ) were
Fig. 11. Contour plots of the influence of PEG compositions, and membrane thickness on the drug-release rate per area at different VA contents.
2418
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419
obtained: SPEG at low vinyl acetate composition and low membrane thickness (LOW VA–LOW TH), SPEG at high vinyl acetate composition and high membrane thickness (HIGH VA–HIGH TH), SPEG at low vinyl acetate composition and high membrane thickness (LOW VA–HIGH TH), and SPEG at high vinyl acetate composition and low membrane thickness (HIGH VA–LOW TH). The results are shown in Fig. 10. As expected, SPEG is higher for HIGH VA–LOW TH. Also, an increase in PEG composition is accompanied by an increase in the drug-release rate (positive correlation). When compared to STH at different combinations of VA and PEG (Fig. 11), the drug-release rate is overall more sensitive to a change in PEG compositions for the levels investigated. Also, a negative correlation exists between drug-release rate and membrane thickness. From Fig. 11, the system responds better to a change in membrane thickness for at HIGH VA–HIGH PEG. This finding is in good agreement with Figs. 5 and 6 which show a more drastic change in the drug-release rate at 19% VA and 50% PEG.
4. Conclusion Franz diffusion cells were used to study the permeation of estradiol through ethylene–vinyl acetate copolymer membrane. Twenty-four samples were obtained by varying PEG 400 volume fractions (20, 30, 40, and 50%), VA contents of the membrane (9 and 19%), and membrane thickness (TH) (0.00508, 0, 0.00762, 0.01016 cm). The results show that the permeation rate rose with an increase in VA content and PEG composition and a decrease in membrane thickness. An artificial neural network, with two hidden nodes, was developed and validated to predict estradiol release rates (per unit area) as a function of the explanatory variables. A high correlation between experimental and predicted drug delivery rates was achieved for the training (75% of the data) and testing sets (25% of the data). A constrained optimization procedure, implemented in Mathematica® (Wolfram Research, Inc.), was applied to calculate the optimum VA content, PEG volume fraction, and membrane thickness for a given delivery rate. Because of the non-unique character of this inverse problem, it was necessary to define an admissible solution space Ψ , made up of data used in the training and testing phases. Sixty-one percent of the training data agreed with experimental results while 100% of the testing data were correctly predicted. Input sensitivity analysis showed that at specific VA levels, drug delivery was more sensitive to changes in PEG compositions than to changes in membrane thickness. The non-uniqueness of the inversion method and the accuracy of the procedure were studied using neural network-based two-dimensional contour plots. Using the methodology proposed, accurate drug-delivery devices can be tailored to meet specific drug-release requirements with the least amount of experiments.
References Berkland, C., King, M., Kim, K., & Pack, D. W. (2002). Precise control of PLG microsphere size provides enhanced control of drug-release rate. Journal of Controlled Release, 82(1), 137–147. Bettini, R., Colombo, P., & Peppas, N. A. (1995). Solubility effects on drug transport through pH-sensitive, swelling-controlled release systems: Transport of theophylline and metoclopramide monohydrochloride. Journal of Controlled Release, 37(1–2), 105–111. Caldwell, R. B. (1995). Three methods of neural network sensitivity analysis for input variable reduction: A case study in forecasting the S&P 500 index (Part 1). NeuroVe$t Journal, 3(6), 22–25. Caldwell, R. B. (1996). Three methods of neural network sensitivity analysis for input variable reduction: A case study in forecasting the S&P 500 index (part 2). NeuroVe$t Journal, 4(1), 16–22. Cun, Y. L., Denker, J. S., Solla, S. A. (1990). Optimal brain damage. In D. S. Touretzky (Ed.), Advances in neural information processing systems (Part 2). San Mateo, CA: Morgan Kaufmann Publishers, Inc. Eikens, B., Karim, M. N., Simon, L. (2001). Combining neural networks with first principle models for bioprocess modeling. In I. M. Mujtaba, M. A. Hussain (Eds.), Applications of neural networks and other learning technologies in process engineering. London: Imperial College Press. Flanagan, N. (2003). Fast-dissolving technologies, pulsing tablets, and transdermal means of drug transport. Genetic Engineering News, 23(4), 30–31. Gao, L., & Loney, N. W. (2001). Evolutionary polymorphic neural network in chemical process modeling: Computers. Computers and Chemical Engineering, 25(11–12), 1403–1410. Grassi, M., Yuk, S. H., & Cho, S. H. (1999). Modeling of solute transport across a temperature-sensitive polymer membrane. Journal of Membrane Science, 152(2), 241–249. Hassibi, B., Stork, D. G. (1993). Second order derivatives for network pruning: Optimal brain surgeon. In D. S. Lippman, J. E. Moody, D. S. Touretzky (Eds.), Advances in neural information processing systems (Part 5). San Mateo, CA: Morgan Kaufmann Publishers, Inc. Howes, P., & Crook, N. T. (1999). Using input parameter influences to support the decisions of feed-forward neural networks. Neurocomputing, 24(1–3), 191–206. Karim, M. N., Hodge, D., & Simon, L. (2003). Data-based modeling and analysis of bioprocesses: Some real experiences. Biotechnology Progress, 19(5), 1591–1605. Kocova El Arini, S., & Leuenberger, H. (1995). Modeling of drug release from polymer matrices: Effect of drug loading. International Journal of Pharmaceutics, 121(2), 141–148. Koizumi, T., Ritthidej, G. C., & Phaechamud, T. (2001). Mechanistic modeling of drug release from chitosan coated tablets. Journal of Controlled Release, 70(3), 277–284. Langer, R. (1998). Drug delivery and targeting. Nature, 392(6679 Suppl.), 5–10. Lim, C. P., Quek, S. S., & Peh, K. K. (2003). Prediction of drug release profiles using an intelligent learning system: an experimental study in transdermal iontophoresis. Journal of Pharmaceutical and Biomedical Analysis, 31(1), 159–168. Principe, J. C., Euliano, N. R., Lefebvre, W. C. (2000). Neural and adaptive systems: fundamentals through simulations. New York: John Wiley and Sons. Shah, P. S., & Zatz, J. L. (1992). Plasticization of cellulose esters used in the coating of sustained release solid dosage forms. Drug Development and Industrial Pharmacy, 18(18), 1759–1772. Shin, S. C., & Byun, S. Y. (1996). Controlled release of ethinylestradiol from ethylene–vinyl acetate membrane. International Journal of Pharmaceutics, 137(1), 95–102. Shin, S.-C., & Lee, H.-J. (2002). Controlled release of triprolidine using ethylene–vinyl acetate membrane and matrix systems. European Journal of Pharmaceutics and Biopharmaceutics, 54(2), 201–206.
L. Simon, M. Fernandes / Computers and Chemical Engineering 28 (2004) 2407–2419 Siegel, R. A., & Pitt, C. G. (1995). A strategy for oscillatory drug release: General scheme and simplified theory. Journal of Controlled Release, 33(1), 173–188. Siepmann, J., Lecomte, F., & Bodmeier, R. (1999). Diffusion-controlled drug delivery systems: calculation of the required composition to achieve desired release profiles. Journal of Controlled Release, 60(2), 379–389. Simon, L., & Karim, M. N. (2001). Probabilistic neural networks using Bayasian decision strategies and a modified Gompertz model for growth phase classification in the batch culture of Bacillus subtilis. Biochemical Engineering Journal, 7(1), 41–48. Simon, L., Karim, M. N., & Schreiweis, A. (1998). Prediction and classification of different phases in a fermentation using neural networks. Biotechnology Techniques, 12(4), 301–304.
2419
Simon, L., & Karim, M. N. (2002). Control of starvation-induced apoptosis in CHO cell cultures. Biotechnology and Bioengineering, 78(6), 645– 657. Subramanian, N., Yajnik, A., & Murthy, R. S. R. (2003). Artificial Neural Network as an alternative to multiple regression analysis in optimizing formulation parameters of cytarabine liposomes. AAPS PharmSciTech, 5(1), 1–9 (Article 4). Takahara, J., Takayama, K., & Nagai, T. (1997). Multi-objective simultaneous optimization technique based on an artificial neural network in sustained release formulations. Journal of Controlled Release, 49(1), 11–20. Takayama, K., Fujikawa, M., Obata, Y., & Morishita, M. (2003). Neural network-based optimization of drug formulations. Advanced Drug Delivery Reviews, 55(9), 1217–1231.