www.gusucode.com > simbio 案例源码程序 matlab代码 > simbio/EstimatingTheBioavailabilityOfADrugExample.m
%% Estimating the Bioavailability of a Drug % In this example, you will use the parameter estimation capabilities of % SimBiology(TM), to calculate |F|, the bioavailability, of the drug % ondansetron. You will calculate |F| by fitting a model of absorption and % excretion of the drug to experimental data tracking drug concentration % over time. % % This example requires Optimization Toolbox(TM). % Copyright 2015-2016 The MathWorks, Inc. %% Background % Most drugs must be absorbed into the bloodstream in order to become % active. An intravenous (IV) administration of a drug is one way to % achieve this. However, it is impractical or impossible in many cases. %% % When a drug is not given by IV, it follows some other route into the % bloodstream, such as absorption through the mucous membranes of the GI % tract or mouth. Drugs administered through a route other than IV % administration are generally not completely absorbed. Some portion of the % drug is directly eliminated and never reaches the bloodstream. %% % The percentage of drug absorbed is the bioavailability of the drug. % Bioavailability is one of the most important pharmacokinetic properties % of a drug. It is useful when calculating safe dosages for non-IV % routes of administration. Bioavailability is calculated relative to % an IV administration. When administered intravenously, a drug has 100% % bioavailability. Other routes of administration tend to reduce the % amoutn of drug that reaches the blood stream. %% Modeling Bioavailability % Bioavailability can be modeled using one of several approaches. In this % example, you use a model with a GI compartment and a blood plasma % compartment. Oral administration is modeled by a dose event in the GI % compartment. IV adminsitration is modeled by a dose event in the blood % plasma compartment. %% % The example models the drug leaving the GI compartment in two ways. The available % fraction of the drug is absorbed into the bloodstream. The remainder is % directly eliminated. The total rate of elimination, |ka|, is % divided into absorption, |ka_Central|, and direct elimination, |Cl_Oral|. The % bioavailability, |F|, connects total elimination with |ka_Central| and % |Cl_Oral| via two initial assignment rules. %% % % ka_Central = F*ka % Cl_Oral = (1-F)*ka % %% % The drug is eliminated from the |Blood_Plasma| compartment through first-order % kinetics, at a rate determined by the parameter |Cl_Central|. %% % % <<../bioavailability_model.png>> % %% % Load the project that contains the model |m1|. sbioloadproject(fullfile(matlabroot,'examples','simbio','Bioavailability.sbproj'),'m1'); %% Format of the Data for Estimating Bioavailability % You can estimate bioavailability by comparing intrapatient measurements % of drug concentration under different dosing conditions. For instance, a % patient receives an IV dose on day 1, then receives an oral dose on day 2. On both days, we can measure the blood plasma % concentration of the drug over some period of time. %% % Such data allow us to estimate the bioavailability, as well as % other parameters of the model. Intrapatient time courses were generated % for the drug ondansetron, reported in [2] and reproduced in [1]. %% % Load the data, which is a table. load(fullfile(matlabroot,'examples','simbio','ondansetron_data.mat')); %% % Convert the data to a |groupedData| object because the fitting function % |sbiofit| requires it to be a |groupedData| object. %% gd = groupedData(ondansetron_data); %% % Display the data. %% gd %% % The data have variables for time, drug concentration, grouping % information, IV, and oral dose amounts. Group 1 % contains the data for the IV time course. Group 2 contains the data for % the oral time course. |NaN| in the Drug column means no measurement was % made at that time. |NaN| in one of the dosing columns means no dose was % given through that route at that time. %% % Plot the pharmacokinetic profiles of the oral dose and IV administration. %% plot(gd.Time(gd.Group==1),gd.Drug(gd.Group==1),'Marker','+') hold on plot(gd.Time(gd.Group==2),gd.Drug(gd.Group==2),'Marker','x') legend({'8 mg IV','8 mg Oral'}) xlabel('Time (hour)') ylabel('Concentration (milligram/liter)') %% % Notice there is a lag phase in the oral dose of about an hour while the % drug is absorbed from the GI tract into the bloodstream. %% Fitting the Data % Estimate the following four parameters of the model: %% % % * Total forward rate out of the dose compartment, |ka| % * Clearance from the |Blood_Plasma| compartment, |clearance| % * Volume of the |Blood_Plasma| compartment % * Bioavailability of the orally administered drug, |F| %% % Set the initial values of these parameters and specify the log transform % for all parameters using an |estimatedInfo| object. init = [1 1 2 .8]; estimated_parameters = estimatedInfo({'log(ka)','log(clearance)',... 'log(Blood_Plasma)','logit(F)'},'InitialValue',init); %% % Because |ka|, |clearance|, and |Blood_Plasma| are positive physical % quantities, log transforming reflects the underlying physical constraint % and generally improves fitting. This example uses a logit transform on % |F| % because it is a quantity constrained between 0 and 1. The logit transform % takes the interval of 0 to 1 and transforms it by taking the log-odds of % |F| (treating |F| as a probability). For a few drugs, like % theophyline, constraining |F| between 0 and 1 is inappropriate because oral % bioavailability can be greater than 1 for drugs with unusual absorption % or metabolism mechanisms. %% % Next, map the response data to the corresponding model component. In the % model, the plasma drug concentration is represented by % |Blood_Plasma.Drug_Central|. The corresponding concentration data is the % |Drug| variable of the groupedData object |gd|. responseMap = {'Blood_Plasma.Drug_Central = Drug'}; %% % Create the dose objects required by |sbiofit| to handle the dosing % information. First, create the IV dose targeting |Drug_Central| and the % oral dose targeting |Dose_Central|. iv_dose = sbiodose('IV','TargetName','Drug_Central'); oral_dose = sbiodose('Oral','TargetName','Drug_Oral'); %% % Use these dose objects as template doses to generate an array of dose % objects from the dosing data variables |IV| and |Oral|. doses_for_fit = createDoses(gd,{'IV','Oral'},'',[iv_dose, oral_dose]); %% % Estimate parameters using |sbiofit|. opts = optimoptions('lsqnonlin','Display','final'); results = sbiofit(m1, gd,responseMap,estimated_parameters,doses_for_fit,... 'lsqnonlin',opts,[],'pooled',true); %% Interpreting Results % First, check if the fit is successful. plot(results) %% % Overall, the results seem to be a good fit. However, they do not capture % a distribution phase over the first hour. It might be possible to improve % the fit by adding another compartment, but more data would be required to % justify such an increase in model complexity. %% % When satisfied with the model fit, you can draw conclusions about % the estimated parameters. Display the parameters stored in the results % object. results.ParameterEstimates %% % The parameter |F| is the bioavailability. The result indicates % that ondansetron has approximately a 66% bioavailability. This estimate in line % with the literature reports that oral administration of ondansetron in % the 2-24 milligram range has a 60% bioavailability [1,2]. %% % |Blood_Plasma| is the volume of distribution. This result is reasonably close to % the 160 liter Vd reported for ondansetron [1]. The estimated clearance % is 45.4 L/hr. %% % |ka| does not map directly onto a widely reported % pharmacokinetic parameter. Consider it from two perspectives. We % can say that 66% of the drug is available, and that the available drug % has an absorption parameter of 0.4905/hr. Or, we can say that drug % clearance from the GI compartment is 0.7402/hr, and 66% of the drug % cleared from the GI tract is absorbed into the bloodstream. %% Generalizing This Approach % |lsqnonlin|, as well as several other optimization algorithms supported % by |sbiofit|, are local algorithms. Local algorithms are subject to the % possibility of finding a result that is not the best result over all % possible parameter choices. Because local algorithms do not guarantee % convergence to the globally best fit, when fitting PK models, restarting % the fit with different initial conditions multiple times is a good % practice. Alternatively, |sbiofit| supports several global methods, such % as particle swarm, or genetic algorithm optimization. Verifying that a % fit is of sufficient quality is an important step before drawing % inferences from the values of the parameters. %% % This example uses data that was the mean time course of several patients. % When fitting a model with data from more patients, some parameters might % be the same between patients, some not. Such requirements introduce the % need for hierarchical modeling. You can perform hierarchical modeling can % by configuring the |CategoryVariableName| flag of |EstimatedInfo| object. %% References % % # Roila, Fausto, and Albano Del Favero. "Ondansetron clinical pharmacokinetics." Clinical Pharmacokinetics 29.2 (1995): 95-109. % # Colthup, P. V., and J. L. Palmer. "The determination in plasma and pharmacokinetics of ondansetron." European Journal of Cancer & Clinical Oncology 25 (1988): S71-4. %