Probabilistic Physiologically Based Pharmacokinetic Model for Penicillin G in Milk From Dairy Cows Following Intramammary or Intramuscular Administrations

Probabilistic Physiologically Based Pharmacokinetic Model for Penicillin G in Milk From Dairy... Abstract Penicillin remains one of the most frequently identified violative drug residues in food-producing animals. The predominant violations of penicillin were found in cull dairy cows. In the United States, procaine penicillin G is approved to be used in dairy cows through intramuscular (IM) and intramammary (IMM) administrations. Physiologically based pharmacokinetic (PBPK) models are useful tools to predict withdrawal intervals and tissue residues of drugs in food animals to ensure food safety, especially for extralabel drug use due to the scarcity of experimental data after extralabel administrations. Currently, no PBPK model is available to predict penicillin concentrations in milk. A population PBPK model with a physiologically based compartment for the mammary gland was established for penicillin G in dairy cows. The model predicted the tissue and milk residues well based on comparison with data from previous pharmacokinetic studies. The predicted milk discard interval of procaine penicillin G administered at 10 times the label dose for 3 repeated IM administrations was 182 h, and 122 h at 4 times the label dose after 3 repeated IMM infusions. Predicted results showed that even 4 times label dose did not lead to violative tissue residues in healthy dairy cows with IMM infusions. The predominant violations found in cull dairy cows may be caused by altered pharmacokinetics due to mastitis, other diseases, and/or interactions with other drugs, which have impacts on penicillin distribution and elimination. The current PBPK model can help predict milk discard interval for penicillin following extralabel use through IM and IMM administrations. PBPK modeling, dairy cows, food safety, drug residue, Food Animal Residue Avoidance Databank (FARAD), penicillin INTRODUCTION Penicillin is one of the top 3 most commonly detected violative residues in animal-derived foods reported by the National Residue Program of United States Department of Agriculture (USDA) from 2014 to 2016 (USDA, 2015, 2017a,b). This is mainly because penicillin G is one of the most widely used antimicrobials and tested regularly (FDA, 2013; Portis et al., 2012; Vogel et al., 2001). Extralabel use of penicillin G in food-producing animals is very common (Chiesa et al., 2006) due to reduced sensitivity to penicillin in bacteria. The extralabel (also called off-label) use of veterinary drugs means to use a drug in a manner that is different from the approved labeling under the supervision of a veterinarian, and extralabel use may only occur if prescribed by a veterinarian (FDA, 2017). Consumption of beef or pork products containing violative penicillin residues can lead to anaphylactic reactions for sensitive individuals (Dayan, 1993; Gomes and Demoly, 2005; Raison-Peyron et al., 2001). In addition, violative residues of penicillin G in milk could interfere with starter cultures for fermented dairy products (Payne et al., 2006; Suhren, 1996) and may result in suspension of the producer’s permit or certification (NRC, 1999). Drug residues above the regulatory safe level in animal-derived products challenge the global food safety (Baynes and Riviere, 2014; Baynes et al., 2016; Paige et al., 1997). The U.S. Food and Drug Administration (FDA) has established a zero tolerance for penicillin in milk products (FDA, 2017). In actual practice, FDA uses the “safe level” of 5 ng/ml for penicillin residues in milk. The USDA National Residue Program reported that 135 out of 213 (63%) and 153 out of 216 (71%) violative samples in cull dairy cows were found with penicillin G violations in 2015 and 2016, respectively (USDA, 2017a,b). The predominant violations of penicillin found in cull dairy cows lead to the concerns of potential drug residues in the cow-derived food products, including milk. In the United States, procaine penicillin G is approved to be used in dairy cows through intramuscular (IM) and intramammary (IMM) administrations. The approved doses of procaine penicillin G are at a daily dose of 6600 IU/kg of body weight (6.5 mg/kg) for no more than 7 consecutive days through the IM injection (Papich et al., 1993), and 100 000 IU (99 mg) into each affected quarter of udders every 12 h for no more than 3 doses (United States Pharmacopeia, 2003). The procaine penicillin G is slowly absorbed after IM administration (Papich, 1987). The major metabolite of penicillin G, penicilloic acid, does not have the antimicrobial activity and accounts for about 16%–30% of an IM dose of penicillin G in humans (DrugBank, 2018; Wishart et al., 2006). Unlike in humans, the metabolism of penicillin G in food-producing animals is not well characterized, and the elimination of penicillin G is mainly through renal secretion (Papich and Riviere, 2009). The determination of withdrawal intervals (WDIs) for extralabel use of veterinary drugs is important to avoid the violations of tissue residues. The WDI or milk discard time is the time for drug residues in the edible tissues or milk to deplete below concentrations that are considered safe for human consumptions (FDA, 2012). The milk produced before the milk discard time must be discarded to avoid drug residues in milk (EMA, 2000). According to the Pasteurized Milk Ordinance (FDA, 2015a), milk violation tests carried out by FDA are based on the samples from bulk milk tank. In the United States, the 99th percentile tolerance limit with 95% confidence is used to determine the withdrawal period and milk discard time (FDA, 2006, 2016). In Europe, it is recommended to determine the withdrawal period with the upper one-sided 95% tolerance limit for the drug residue below the maximum residue limit with 95% confidence (EMA, 2016). Due to the lack of pharmacokinetic data after extralabel use of veterinary drugs, it is difficult to use statistical tolerance methods or classic pharmacokinetic methods to determine the WDIs after extralabel administrations. The physiologically based pharmacokinetic (PBPK) model is a useful tool to predict tissue residues and WDIs of veterinary drugs for food safety assessment (Henri et al., 2016a; Lin et al., 2016b). PBPK modeling is a mechanism-based approach to simulate the absorption, distribution, metabolism, and elimination of chemicals (Andersen, 2003; Lin et al., 2016a; WHO, 2010) and to predict the concentrations in the edible tissues and products beyond times when data are not available. By changing the related parameter values, the PBPK model can be extrapolated to different species and therapeutic scenarios. There is no whole-body PBPK model for predicting the milk discard intervals of drugs after extralabel use. Recently, we published a multiroute population PBPK model for penicillin G in beef cattle (Li et al., 2017). On the basis of cattle model, in the present study we developed a PBPK model with physiologically based compartment of mammary glands to predict milk concentrations and discard intervals after IM or IMM administration of penicillin G at extralabel doses in dairy cows. The systemic absorption from the mammary gland and the distribution of penicillin G were considered. The volume change by milk secretion and episodic milking process were also included. Probabilistic analysis using Monte Carlo simulation was applied to address the uncertainties of parameter values. This PBPK model contributes to the field of toxicological science and biological modeling by establishing a quantitative tool that allows to conduct drug residue safety assessment of cow-derived food products and creating a new physiologically based approach to simulate milk production and episodic milking processes, which can be applied to other species, including rodents, humans, and goats. MATERIALS AND METHODS Data source for model calibration The pharmacokinetic data used in the calibration and evaluation of the PBPK model were acquired from the Comparative Pharmacokinetic Database from the Food Animal Residue Avoidance Databank (FARAD). FARAD is supported by USDA National Institute of Food and Agriculture with an overarching goal to provide veterinary practitioners the most current and accurate information to facilitate the production of safe foods of animal origin through the prevention and mitigation of violative chemical residues in food animal products (Craigmill et al., 2006; Riviere et al., 1986, 2017). Pharmacokinetic data in dairy cows after IM or IMM administration of procaine penicillin G were selected. The formulation of aqueous suspension was used for IM injection of procaine penicillin G. For IMM injection, the formulation of procaine penicillin G in peanut oil, which is more relevant to the current clinical formulations, was used. Brief description and key information of selected pharmacokinetic studies are summarized in Table 1. The graphic pharmacokinetic data were extracted from selected studies using WebPlotDigitizer (version 3.10, https://automeris.io/WebPlotDigitizer/; last accessed March 29, 2018.). Table 1. Summary of Pharmacokinetic Studies of Penicillin G in Dairy Cows Used for Calibration and Evaluation of the PBPK Model PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg All data used for model calibration and evaluation were from healthy cows. The abbreviations for routes: IM, intramuscular injection; IMM, intramammary infusion. The abbreviations for matrices: P, plasma; M, milk. The abbreviations for assays: HPLC, high performance liquid chromatography; MM, microbiological methods. NA, not available or not applicable. Table 1. Summary of Pharmacokinetic Studies of Penicillin G in Dairy Cows Used for Calibration and Evaluation of the PBPK Model PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg All data used for model calibration and evaluation were from healthy cows. The abbreviations for routes: IM, intramuscular injection; IMM, intramammary infusion. The abbreviations for matrices: P, plasma; M, milk. The abbreviations for assays: HPLC, high performance liquid chromatography; MM, microbiological methods. NA, not available or not applicable. Model structure The PBPK model for penicillin G with the mammary gland compartment was established based on our pervious generic PBPK model for penicillin G in beef cattle (Li et al., 2017). The model consisted of 8 compartments corresponding to different tissues in the body, including liver, kidney, muscle, fat, lung, the rest of body and udder connected by the circulating blood system (Figure 1A). For food safety assessment purposes, all the major edible tissues, including liver, kidney, muscle, and fat, were included. The udder was considered as a compartment, because the model considered the absorption of penicillin G from milk to the systemic circulation and the excretion of penicillin G into the milk. The lung was also included as a compartment as it is the therapeutic target for penicillin G. Each compartment was defined by a tissue weight and tissue blood flow rate. The compartments for urine and feces were established as virtual excretory compartments without volume changes. The flow-limited model, which performed well for the previous PBPK model of penicillin G, was applied in the current model. Figure 1. View largeDownload slide A schematic diagram of the physiologically based pharmacokinetic (PBPK) model for penicillin G in dairy cows and the 1- and 2-compartment udder models. A, The PBPK model for penicillin G in dairy cows. Two different administration routes including intramuscular (IM) and intramammary (IMM) administrations are presented in the model. The dash line for the compartment of udder indicates that the volume changes since the preceding milking and this could be simulated using the Hill-Langmuir equation. B, The 1-compartment model for the udder. The milk space represents the combined milk space including alveoli, ducts, and the cistern. C, The 2-compartment model for the udder. Figure 1. View largeDownload slide A schematic diagram of the physiologically based pharmacokinetic (PBPK) model for penicillin G in dairy cows and the 1- and 2-compartment udder models. A, The PBPK model for penicillin G in dairy cows. Two different administration routes including intramuscular (IM) and intramammary (IMM) administrations are presented in the model. The dash line for the compartment of udder indicates that the volume changes since the preceding milking and this could be simulated using the Hill-Langmuir equation. B, The 1-compartment model for the udder. The milk space represents the combined milk space including alveoli, ducts, and the cistern. C, The 2-compartment model for the udder. Model calibration and parameterization The current model was based on the published PBPK model for penicillin G in beef cattle, and the physiologically based mammary gland compartment and IMM infusion were added to the model. The physiological parameters related to dairy cows were acquired from previous research. The average values of blood flow to the udder (QUC), the tissue volume fractions of the udder (VUC), venous blood volume (VvenC), and arterial blood volume (VartC) were obtained or calculated based on previous experimental studies (Campbell and Marshall, 2016; Gionbelli et al., 2015). For chemical-specific parameters (eg, partition coefficients) of penicillin G, the original values used in the cattle model (Li et al., 2017) were applied in the current dairy cow model. The parameters related to the udder compartment and IMM infusions, including VmaxC, Km, and KabC were fitted with pharmacokinetic data chosen for calibration (Table 1). VmaxC and Km were optimized by fitting with calibration datasets following IM administrations (Hogh and Rasmussen, 1966; Van OS et al. 1974), and KabC was estimated by fitting to the dataset after IMM administration (Mercer et al., 1974). The values of all physiological parameters and chemical-specific parameters used in the PBPK model are shown in Table 2. Table 2. Parameter Values Used in the PBPK Model for Penicillin G in Milk From Dairy Cows Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 “Model fitted” indicates that the parameter values were determined by fitting the model with calibration pharmacokinetic data. Table 2. Parameter Values Used in the PBPK Model for Penicillin G in Milk From Dairy Cows Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 “Model fitted” indicates that the parameter values were determined by fitting the model with calibration pharmacokinetic data. Establishment and validation of the udder compartment The compartments representing the milk and udder were established according to the physiologically based mathematical model described by Whittem et al. (2012). In brief, the Hill-Langmuir equation was used to describe the volume-time relationship for milk production. The concentration change of penicillin G due to milk production and episodic milk emptying was considered. Both the 1-compartment and 2-compartment models were created and validated using previously reported data (Davis et al., 1998). The 1-compartment model was established considering the milk space in the udder as 1 compartment with the function of milk secretion and storage (Figure 1B). The 2-compartment model divided the milk space into the alveoli, which serve the function of milk secretion and storage, and the cistern, which only has the function for milk storage (Figure 1C). For the 2-compartment model, the distribution of penicillin G between alveoli and cistern compartments was incorporated. Berkeley Madonna (Version 8.3.23.0; University of California at Berkeley, California) was used to develop the PBPK model with the mammary gland compartment and run all simulations. Codes of the model are provided in the Supplementary Data and will be available on our website (http://iccm.k-state.edu/; last accessed March 29, 2018). As the model was constructed based on the previous PBPK model, only key and new mathematical equations are described in detail in following paragraphs. For all the other equations, details have been provided in the Supplementary Data and the previously published PBPK model for penicillin G in cattle and swine (Li et al., 2017). For the 1-compartment udder model, the volume change of the milk space in udder was modeled the same as the volume change of the milk production. The milk volume in the udder compartment depends on the episodic milking and the milk production between milkings. The milk production was described using the Hill-Langmuir equation (equation 1) reported in the previous milk model (Whittem et al., 2012). According to the previous experimental data (Davis et al., 1998), the equation describes the milk production rises rapidly in a linear trend at the beginning, and the rate of milk secretion slows down with the increased hydrostatic pressure in the udder compartment. The volume of residual milk in udder can be determined by the fraction of maximum volume of the milk space (equation 2). Previous experimental research reported the residual milk fractions are in the range of 4%–14% (Carruthers et al., 1993; Isaksson and Arnarp, 1988; Knight et al., 1994), here 4% was used as the value for the residual fraction as reported in the previous model (Whittem et al., 2012). The total volume of milk space should account for the residual volume after milking and the volume of milk secreted after the preceding milking (equation 3). Vmilk=tm×VMmilksptm+s, (1) Vresidue=VMmilksp×F, (2) Vmilksp=Vresidue+Vmilk, (3) where Vmilk is the volume of milk produced since the preceding milking (l); tm is the time after the preceding milking (h); VMmilksp is the maximum volume of the milk space including spaces in alveoli, ducts, and the cistern (l); s is the time for Vmilk reaching half of VMmilksp (h); Vmilksp is the volume of the milk space (l); Vresidue is the residual volume of the milk after milking (l); F is the ratio of the residual volume to the maximum volume of the milk space (unitless). Compared with the 1-compartment model, the 2-compartment model is more physiologically based, and can be used to simulate the drug distribution and movement in different parts of udder compartment. However, the 2-compartment model is relatively complex, includes more parameters and requires robust datasets to validate. Briefly, the volume changes in the total milk space and the alveolar compartment were both simulated with the Hill-Langmuir equation (equation 4). The volume in cistern compartment was determined by subtracting the volume of alveoli from the total milk space (equation 5). The residual volume of milk in the milk space of udder was mostly retained in the alveoli (Vetharaniam et al., 2003), and the residual volume in cistern was assumed as zero immediately after milking (Knight et al., 1994). Valv= Valvresidue+tm×VMalvtm+salv, (4) Vcis=Vmilksp-Valv, (5) where Valv is the volume of milk in the alveoli (l); tm is the time after the preceding milking (h); VMalv is the maximum volume of the alveoli (l); salv is the time for the volume of milk in alveoli reaching half of VMalv (h); Valvresidue is the residual volume of the alveoli after milking (l); Vcis is the volume of milk in the cistern (l). Both the 1-compartment and 2-compartment udder models were validated (Figs. 2A and 2B) based on published experimental data (Davis et al., 1998). As available pharmacokinetic data in milk are mostly sparse data and to avoid over parameterization, the 1-compartment udder model was used in the final PBPK model. The code describing the 2-compartment udder model is provided in the code in the Supplementary Data and can be used in the future when needed and when more robust datasets are available. Figure 2. View largeDownload slide Validation and application of the milk models. Model predictions for milk volumes (solid line) and observed data (black squares, circles, and diamonds) are compared in panel A for the 1-compartment model and in panel B for the 2-compartment model. The observed data are from previous study by Davis et al. (1998). In panel C, the milk production volumes by different milking intervals were simulated using 1-compartment milking model based on the parameter values reported in Whittem et al. (2012) (VMmilksp = 26.8 l, s = 23.2 h, F = 0.04). Figure 2. View largeDownload slide Validation and application of the milk models. Model predictions for milk volumes (solid line) and observed data (black squares, circles, and diamonds) are compared in panel A for the 1-compartment model and in panel B for the 2-compartment model. The observed data are from previous study by Davis et al. (1998). In panel C, the milk production volumes by different milking intervals were simulated using 1-compartment milking model based on the parameter values reported in Whittem et al. (2012) (VMmilksp = 26.8 l, s = 23.2 h, F = 0.04). The previous milk model assumed no absorption and mass transfer of drugs between the udder compartment and systemic circulation (Whittem et al., 2012). In order to incorporate the udder compartment into a whole-body PBPK model, the absorption of penicillin G from the udder to the systemic circulation and the excretion of penicillin G from plasma to milk were considered so that the current model could be applied to both IMM and IM administrations. The absorption and excretion processes of penicillin G involve both the passive diffusion and active carrier-mediated transport (Rasmussen, 1959; Schadewinkel-Scherkl et al., 1993; Ziv and Sulman, 1974). On the basis of previously reported experimental data (Schadewinkel-Scherkl et al., 1993), the carrier-mediated transport accounts for around 25% of the absorption of penicillin G at mammary gland (Supplementary Figure 1A) and for over 75% of the penicillin G excretion (Supplementary Figure 1B). Due to very limited experimental data for both processes, the parameters related to the kinetics of the absorption and excretion need to be determined by fitting with pharmacokinetic data. In order to reduce uncertainties, the excretion process was considered as the active carrier-mediated transport only, and the absorption process was described as only through the passive diffusion (equation 6). The Michaelis-Menten equation was used for the carrier-mediated transport of penicillin excreted into milk (equation 7), and the first-order linear equation was applied for the passive diffusion. The rate of change of penicillin G in the udder compartment was described using mass balance differential equations as previously described (DeWoskin et al., 2013; Leavens et al., 2012; Lin et al., 2016b). Only the penicillin G not bound to plasma proteins was considered as freely available for excretion. The mass balance equations of penicillin G in the udder compartment are described below. RU=QU (CAfree−CVU) – Rmilkex+Kab×Amilk (6) Rmilkex =Vmax×CVUCVU+Km, (7) CVU=CU/PU, (8) where RU is the rate of change in the concentration of penicillin G in the udder (mg/h); QU represents the volume of blood flow to the udder per hour (l/h); CAfree is the arterial blood concentration of penicillin G not bound with plasma proteins (mg/l); CVU is the concentration of penicillin G in the udder venous blood (mg/l); Rmilkex is the rate for penicillin excreted into milk (mg/h); Kab is the rate constant of penicillin G absorbed from milk (/h); Amilk is the amount of penicillin G in milk (mg); Vmax represents the maximum rate for penicillin G excretion from udder tissue to milk (mg/h); Km is the concentration of penicillin G to achieve the half Vmax (mg/l); CU is the penicillin G concentration in the udder tissue (mg/l); PU is the udder tissue: plasma partition coefficient (unitless). The penicillin input through IMM administration, excretion from plasma to milk, penicillin absorption from milk, and elimination through the episodic milking process were considered for the change of penicillin concentrations in milk (equation 9). The concentration of penicillin in milk was simulated by dividing amount of penicillin in milk by the milk volume (equations 10 and 11). Rmilk=Rinputimm+Rmilkex−Kab*Amilk−Rmilking, (9) d/dt(Amilk)=Rmilk, (10) Cmilk=Amilk/Vmilk, (11) where Rmilk is the rate of change in the concentration of penicillin G in milk (mg/h); Rinputimm stands for the input rate of penicillin through IMM administration (mg/h); Rmilkex is the rate for penicillin excreted into milk (mg/h); Kab is the rate constant of penicillin G absorbed from milk (/h); Amilk is the amount of penicillin G in milk (mg); Rmilking is the rate of change in the concentration of penicillin G by milking process (mg/h); Vmilk is the volume of milk produced since the preceding milking (l); Cmilk is the penicillin G concentration in milk (mg/l). The milking intervals and frequencies have impacts on the milk production (Gehring and Smith, 2006) and lead to the change of drug concentrations in milk. Briefly, higher milk production could lead to further dilution of drug concentrations in milk. In pharmacokinetic studies for penicillin G in milk, many employed unequal intervals such as 10 + 14 h or 15 + 9 h intervals for milking instead of equal 12-h milking intervals. In order to better fit the experimental data, the array equations were used (Supplemenatry equations 1 and 2). The IM and IMM administrations were incorporated into this PBPK model because they are the label routes for procaine penicillin G in dairy cows. The IM administration was simulated using a 2-compartment dissolution model based on the approach used to simulate the absorption of long-acting oxytetracycline (Lin et al., 2015). Detailed description of the dissolution model and the IM injection was included in our previous papers (Li et al., 2017; Lin et al., 2015). The multiple-dose therapeutic scenario was achieved by using the conditional operator of “IF…THEN…ELSE…” to create a control factor (Supplementary equations 3 and 7). The volume of milk was simulated using array functions, which helped to simulate the milking processes as the discrete steps. However, the amount of penicillin G in the milk and the udder was simulated using continuous functions. The DELAY functions were also applied here to help overcome the false peaks due to the function differences. For more details, please refer to the Supplementary Data. Evaluation and sensitivity analysis The performance of the PBPK model was evaluated by comparing model simulations with experimental pharmacokinetic data of milk and plasma not used in the model calibration (Table 1). On the basis of World Health Organization (WHO) guidelines (WHO, 2010), if the results from simulations matched the measured kinetic profiles and were generally within a factor of 2 of the measured values, the model was considered reasonable and validated. These criteria are based on the considerations that calibration datasets and evaluation datasets are obtained under different conditions (eg, different experimental animals/human subjects, laboratories, and detection methods), so some level of discordance is to be expected (WHO, 2010). The goodness of fit between log-transformed values of observed and predicted plasma and milk concentrations were further analyzed with linear regression and the determination coefficient (R2) was calculated for both calibration and evaluation results. The goodness of fit was evaluated using the determination coefficient (R2) for linear regression analysis and the mean absolute percentage error (MAPE) value. The R2 is an indicator for the overall model simulation performance, and a model simulation with a R2 value equal to or higher than 0.75 is considered acceptable. The analysis for MAPE was also carried out based on the previously reported methods (Cheng et al., 2016; Lin et al., 2017). The MAPE value lower than 50% was considered as the acceptable prediction. A local sensitivity analysis was performed for a discrete time point of 24 h to determine which parameters were most influential on the 24-h area under the curve (AUC) of the plasma, liver, kidney, and milk concentrations of penicillin G. Briefly, each parameter was increased by 1% and the corresponding 24-h AUC of penicillin concentrations was computed. Normalized sensitivity coefficient (NSC) was calculated using equations reported previously (Lin et al., 2011; Yoon et al., 2009). We conduced local sensitivity analysis rather than global sensitivity analysis because local sensitivity analysis is sufficient to identify highly sensitive parameters on our selected dose metrics of interest based on the multiple relevant studies (Craigmill, 2003; Leavens et al., 2014; Zeng et al., 2017; Zhu et al., 2017). Global sensitivity analysis could be done in the future if the goal is to study interaction effects between parameters (Lumen et al., 2015; McNally et al., 2011). Monte Carlo analysis Monte Carlo simulation obtains numerical results based on the repeated random sampling of parameter values from specified distribution of each parameter. This method has been used in the applications of PBPK modeling to estimate drug tissue residues and WDIs in the area of food safety (Buur et al., 2006; Yang et al., 2012; Zeng et al., 2017). In the current PBPK model for penicillin G in milk, Monte Carlo simulation was also applied to estimate the effects of parameter uncertainties and intraspecies variability of dairy cows on milk penicillin G concentrations. For these simulations, the relatively sensitive parameters, which were identified from the local sensitivity analysis (Supplementary Table 1), distributed randomly around the mean values were specified in Table 3. Log-normal distributions of model parameters were assumed for all chemical-specific parameters, such as partition coefficients, absorption rate constants, and elimination rate constants. The lognormal transformation for parameters in Monte Carlo simulations was achieved by using Supplementary equations 8–10. Physiological parameters, including body weights, cardiac outputs, and fractions of blood flows and tissue volumes were assumed to be normally distributed (Clewell III et al., 2000; Li et al., 2017; Shankaran et al., 2013; Sterner et al., 2013; Tan et al., 2006; Yang et al., 2015). Table 3. Values and Distributions of Parameters Used in the Monte Carlo Analysis for the PBPK Model of Dairy Cows Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 A pound sign (#) indicates the CV was calculated based on previous experimental data. Default values reported by previous models were used for CVs of the other parameter values. The default values of CVs for physiological parameters were 0.3, and for chemical-specific parameters were 0.2. The 5th and 95th percentiles of each parameter were calculated as the lower and upper bounds to represent 95% confidence interval. Please refer to Table 2 for the description of each parameter. Table 3. Values and Distributions of Parameters Used in the Monte Carlo Analysis for the PBPK Model of Dairy Cows Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 A pound sign (#) indicates the CV was calculated based on previous experimental data. Default values reported by previous models were used for CVs of the other parameter values. The default values of CVs for physiological parameters were 0.3, and for chemical-specific parameters were 0.2. The 5th and 95th percentiles of each parameter were calculated as the lower and upper bounds to represent 95% confidence interval. Please refer to Table 2 for the description of each parameter. The detail of Monte Carlo simulation in Berkeley Madonna was introduced in our previous generic PBPK model for penicillin G in beef cattle (Li et al., 2017), and the code of the current model is available in the Supplementary Data. Briefly, model parameters were varied randomly around the values (central tendencies) used or estimated in model calibration by their probabilistic distributions (variability). The values of sensitive physiological and chemical-specific parameters were randomly selected based on their distributions, it is necessary to use adjustment factors to ensure the plausibility and mass balance for the PBPK model (Covington et al., 2007; Gelman et al., 1996; Sterner et al., 2013). The Monte Carlo simulation was set up to batch run for 1000 times with model parameters randomly selected from the defined distributions. Application of the model for milk production and tissue concentration prediction The model of milk production can be used as a stand-alone program to predict the milk production based on the different milking intervals. The 12, 12-h interval, 8, 16-h interval, and 8, 8, 8-h interval were used to simulate and compare the milk productions with different milking intervals. The PBPK model for penicillin G with the mammary gland compartment can be used to predict the systemic absorption and distribution of penicillin G after the IMM injection and to predict the milk concentration of penicillin G following IM injections. For IMM infusions, the label dose 100 000 IU (99 mg per quarter of udders) and 2 commonly used extralabel doses (2× label dose, 198 mg per quarter; 4× label dose, 396 mg per quarter) to 2 of the 4 quarters were simulated for 3 repeated doses with 12-h intervals. For IM administrations, the label dose 6000 IU/kg (6.5 mg/kg) and 2 commonly used extralabel doses (5× label dose, 32.5 mg/kg; 10× label dose, 65 mg/kg) were simulated for 3 repeated doses with 24-h intervals. A schematic illustrating the dosing regimens is provided in the Supplementary Figure 2. The commonly used extralabel doses for IMM infusions were based on the phone call records from the FARAD call center, and for IM administrations were obtained from the previous report (Payne et al., 2006). Infusions of penicillin G for all 4 quarters are only for the worst scenario of mastitis, and infusions for 2 quarters are common (Vilim et al., 1979). Here, IMM infusions for 2 quarters were simulated. These doses were applied as input for the current PBPK model to predict the concentrations of penicillin G in the milk, plasma, liver, and kidney. The Monte Carlo analysis was also used to simulate the distribution of penicillin concentrations based on the distribution of parameter values. Different therapeutic scenarios were analyzed using Monte Carlo simulations. Each simulation was run 1000 times. The mean value, 1st and 99th percentile values of simulated results were calculated and plotted. The extended meat WDI and the milk discard interval were determined when the values of 99th percentiles of the target tissue and milk concentrations, respectively, fell below the tolerance or the safe level. The tolerance of penicillin G for edible tissues in cattle is 0.05 µg/g (FDA, 2013). The safe level of penicillin G for milk is 5 ng/ml (FDA, 2005, 2015b) and the tolerance is 0 by U.S. FDA (Brynes, 2005; FDA, 2017). The limits of detection (LODs) for penicillin G in the milk from current methods are lower than or equal to 1 ng/ml, such as 1 ng/ml (Holstege et al., 2002), 0.25 ng/ml (Liu et al., 2011), and 0.23 ng/ml (Huang et al., 2013). In addition, we compared the model-predicted milk discard intervals from Monte Carlo simulations between 1000 and 10 000 iterations based on a 3 daily IM administration paradigm at the labeled dose level. RESULTS Validation and Application of the Milk Model The 1- and 2-compartment milk models were used to simulate the volume of milk produced. Due to the lack of published data for the volumes of milk produced by different milking intervals, the model was only calibrated with experimental data from a study that determined the pattern of milk accumulation in the milk space of udder over time (Davis et al., 1998), and the parameter values were based on the parameter values from the previous milk model (VMmilksp = 26.8 l, s = 23.2 h, F = 0.04) (Whittem et al., 2012). There were no other data available for the model evaluation. The 1-compartment model predicted the total milk volume in the milk space well (Figure 2A). The simulated results from 2-compartment model were also in good agreement with milk volumes measured in alveoli, cistern, and total milk space in udder (Figure 2B). Both the 1- and 2-compartment milk models could be used as the udder compartment for the PBPK model. The milk model could also be applied to simulate the impact of different milking intervals on the daily milk production (Figure 2C represents the results from the 1-compartment model). The results showed that the 12-h milking intervals produced slightly higher milk volumes compared with 6–18 h alternative intervals. Milking 4 times per day achieved the highest daily milk yields compared with milking 2 or 3 times a day, which is in accordance with previous reported research (Allen et al., 1986; Erdman and Varner, 1995). However, milking more than twice a day does not fit the circadian rhythm and the real-life working shifts well and requires more labor (Lyons et al., 2014). Model Calibration The PBPK model with the 1-compartment milk model was used to simulate penicillin G concentrations in plasma and milk after different therapeutic regimens. Model predictions were compared with observed concentrations in dairy cows exposed to penicillin G through IMM infusion or IM injection (results are shown in Figs. 3 and 4, respectively). Overall, the model well simulated the kinetic profiles of penicillin G in the milk and plasma in dairy cows through both IM and IMM administrations. In particular, the model can accurately predict the penicillin G residues in the milk and plasma even at the later time points, which are important for residue monitoring and the determination of times when the concentrations of residues in the milk fall below or near the tolerance of penicillin G. Figure 3. View largeDownload slide Calibration of the physiologically based pharmacokinetic model using pharmacokinetic data through the intramammary (IMM) infusion. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from dairy cows exposed to procaine penicillin G via single IMM infusion (396 mg per quarter for all quarters). Experimental data are from the previous study (Mercer et al. 1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 3. View largeDownload slide Calibration of the physiologically based pharmacokinetic model using pharmacokinetic data through the intramammary (IMM) infusion. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from dairy cows exposed to procaine penicillin G via single IMM infusion (396 mg per quarter for all quarters). Experimental data are from the previous study (Mercer et al. 1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 4. View largeDownload slide Calibration of the model using the pharmacokinetic data through the intramuscular (IM) injection. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from daily cows exposed to procaine penicillin G via single IM injection (10 mg/kg, A, B; 40 mg/kg, C, D). Experimental data of panels A and B are from the study carried out by Hogh and Rasmussen (1966), and panels C and D are from the pharmacokinetic data reported by Van OS et al. (1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 4. View largeDownload slide Calibration of the model using the pharmacokinetic data through the intramuscular (IM) injection. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from daily cows exposed to procaine penicillin G via single IM injection (10 mg/kg, A, B; 40 mg/kg, C, D). Experimental data of panels A and B are from the study carried out by Hogh and Rasmussen (1966), and panels C and D are from the pharmacokinetic data reported by Van OS et al. (1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Model Evaluation The pharmacokinetic datasets not used for the model calibration were applied to evaluate the model performance. Measured concentrations of penicillin G in the milk and plasma of dairy cows after IMM and IM administrations for single or multiple doses were compared with the model predictions (Figs. 5 and 6). Overall, the model predicted the concentrations of penicillin G adequately at different time points. For the IMM administration, the evaluation results shown in Figure 5 indicate the model well simulated the penicillin G concentrations in milk following different milking intervals. The R2 values of the regression analyses between log-transformed values of measured and simulated concentrations of penicillin G in the milk and plasma from both calibration and evaluation datasets were all higher than 0.75 (Supplementary Figure 3). For the evaluation datasets, the value of overall R2 for the IMM administration was 0.99 (Supplementary Figure 3C), and for the IM administration was 0.92 (Supplementary Figure 3D). The R2 for data in milk was 0.81, and for plasma was 0.97 for the IM administration. All the MAPE values for calibration and evaluation results were lower than 40% (Supplementary Figure 4), which indicates that the model predictions are acceptable. In particular, the MAPE values for plasma for different simulations were generally lower than or near 20% indicating good predictions of plasma concentrations following either IM or IMM administrations. Overall, the PBPK model adequately simulates the measured results in the milk and plasma from the independent studies of dairy cows via both IMM and IM administrations. Figure 5. View largeDownload slide Evaluation of the model using intramammary (IMM) data with different milking intervals. Comparison of model predictions (solid line) and observed data (blue circles) for concentrations of penicillin G in the milk of dairy cows exposed to procaine penicillin G via IMM infusions for repeated 3 doses (1898 mg per quarter for all quarters, A, B; 100 mg per quarter for 2 quarters, C; 100 mg per quarter for all quarters, D) is shown. Experimental data for A, B (mean) are from the study reported by Knappstein et al. (2003), for panel C are from the previous pharmacokinetic study (Vilim et al. 1979), and for panel D are from the study performed by Vilim et al. (1980). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 5. View largeDownload slide Evaluation of the model using intramammary (IMM) data with different milking intervals. Comparison of model predictions (solid line) and observed data (blue circles) for concentrations of penicillin G in the milk of dairy cows exposed to procaine penicillin G via IMM infusions for repeated 3 doses (1898 mg per quarter for all quarters, A, B; 100 mg per quarter for 2 quarters, C; 100 mg per quarter for all quarters, D) is shown. Experimental data for A, B (mean) are from the study reported by Knappstein et al. (2003), for panel C are from the previous pharmacokinetic study (Vilim et al. 1979), and for panel D are from the study performed by Vilim et al. (1980). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 6. View largeDownload slide Evaluation of the model using intramuscular (IM) data with different milking intervals. Comparison of model predictions (solid line and/or dotted line) and observed data (blue squares or circles) for concentrations of penicillin G in the plasma and milk of dairy cows exposed to procaine penicillin G via 1 single IM injection (11 mg/kg, A, B; 2.2 mg/kg, D; 4.4 mg/kg, E; 11 mg/kg, F), and IM injection for repeated 4 doses (6.5 mg/kg, C) is shown. Experimental data for A, B (mean) are from the study reported by Randall et al. (1954), for panel C are from reported pharmacokinetic data (FDA, 2010), and for panels D–F are from the previous study by Edwards et al. (1966). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 6. View largeDownload slide Evaluation of the model using intramuscular (IM) data with different milking intervals. Comparison of model predictions (solid line and/or dotted line) and observed data (blue squares or circles) for concentrations of penicillin G in the plasma and milk of dairy cows exposed to procaine penicillin G via 1 single IM injection (11 mg/kg, A, B; 2.2 mg/kg, D; 4.4 mg/kg, E; 11 mg/kg, F), and IM injection for repeated 4 doses (6.5 mg/kg, C) is shown. Experimental data for A, B (mean) are from the study reported by Randall et al. (1954), for panel C are from reported pharmacokinetic data (FDA, 2010), and for panels D–F are from the previous study by Edwards et al. (1966). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Sensitivity Analysis The sensitivity analysis was applied to screen the sensitive parameters in the PBPK model. Thirty-four model parameters in the PBPK model for dairy cows were analyzed using the method of local sensitivity analysis. Results of the local sensitivity analysis based on 1% variation of the parameter values are shown in Supplementary Table 1. Only parameters with at least 1 absolute value of NSC > 0.15 are shown in the table. For physiological parameters, only the body weight (BW), the cardiac output (QCC), and the blood flow to kidney (QKC) were the relatively sensitive parameters for the model. The results indicate that the PBPK model was not that sensitive to physiological parameters. The rate constant of penicillin G absorbed from milk back into the udder tissue (KabC) was the sensitive parameter for AUCs of the plasma, liver and kidney, and the fraction of undissolved procaine penicillin G (Frac) was sensitive for all AUCs. The AUC of liver was also highly sensitive to the liver partition coefficient (PL) with NSC value of 1.00. The AUC of kidney was highly sensitive to the urine elimination rate constant and the kidney partition coefficient with NSC values of −0.90 and 1.00, respectively. For AUC of milk, it was sensitive to the body weight (BW), the maximum volume of the milk space in udder (VMmilksp), the time for Vmilk reaching half of VMmilksp (s), the ratio of the residual volume to the maximum volume of the milk space (F), and the maximum rate of penicillin G excretion in the mammary gland (VmaxC; body weight-scaled) with NSC values of 0.99, −0.99, 0.61, −0.24, and 1.00, respectively. Probabilistic Analysis The relatively sensitive parameters listed in Supplementary Table 1 were subsequently used for probabilistic analysis based on the current PBPK model for penicillin G in dairy cows. The Monte Carlo sampling technique was applied for the probabilistic analysis. From the simulation results of the PBPK model after label and extralabel use of procaine penicillin G through IM and IMM administrations, the concentrations of penicillin G in the liver are higher than other edible tissues and plasma (Supplementary Figs. 5B and 5D). Therefore, the liver was chosen as the tissue to determine the WDIs for extralabel use of penicillin G in dairy cows. The Monte Carlo simulations showed that the WDIs after 3 repeated IMM infusions with label dose and 4× label dose in dairy cows were 98 and 122 h for milk, and were both 0 h for edible tissues, respectively (Figs. 7A–D). The predicted WDIs after 3 doses via IM injections with label dose and 10× label dose were 125 and 182 h for milk, and were 4 and 5 days for edible tissues (Figs. 7E–H). The exact model-predicted milk discard intervals are reported in this article, but in practice, the recommended milk discard intervals should be adjusted accordingly to fit within a milking schedule for a particular dairy farm. If the estimated WDI of edible tissues was a fraction of a day, it was rounded up to the next whole day. The results of the probabilistic analysis of the PBPK model for 2× label dose through IMM infusions and 5× label dose following IM administrations are available in the Supplementary Figure 6, and the results for probabilistic analysis in the plasma and kidney are provided in the Supplementary Figure 7. No considerable differences were seen between the Monte Carlo simulation results between 1000 and 10 000 iterations (Supplementary Figure 8). The label withdrawal periods were obtained from the Veterinarian’s Guide to Residue Avoidance Management (VetGRAM) of FARAD (Riviere et al., 2017). The labeled milk discard time via IMM infusions is 84 h, and for edible tissues the label withdrawal period is 4 days. The label withdrawal period of milk via IM injections is 48 h, and of edible tissues is 14 days in dairy cows. Figure 7. View largeDownload slide Probabilistic analysis for penicillin G concentrations in milk and liver via intramammary (IMM) and intramuscular (IM) administrations. IMM infusions with label dose (99 mg per quarter for 2 quarters) and 4× label dose (396 mg per quarter for 2 quarters) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (A–D). IM injections with label dose (6.5 mg/kg) and 10× label dose (65 mg/kg) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (E–H). The milking intervals were 12 h. Though zero tolerance was established for penicillin G in milk, the safe level (SL) of 5 µg/kg was used by FDA in actual practice. Tolerances (TOL) for penicillin G in edible tissues are 50 µg/kg. Figure 7. View largeDownload slide Probabilistic analysis for penicillin G concentrations in milk and liver via intramammary (IMM) and intramuscular (IM) administrations. IMM infusions with label dose (99 mg per quarter for 2 quarters) and 4× label dose (396 mg per quarter for 2 quarters) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (A–D). IM injections with label dose (6.5 mg/kg) and 10× label dose (65 mg/kg) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (E–H). The milking intervals were 12 h. Though zero tolerance was established for penicillin G in milk, the safe level (SL) of 5 µg/kg was used by FDA in actual practice. Tolerances (TOL) for penicillin G in edible tissues are 50 µg/kg. DISCUSSION In this study, a PBPK model for penicillin G in dairy cows was established for labeled routes of IM and IMM administrations. On the basis of model, probabilistic analysis was performed to predict the milk discard intervals and the extended WDIs for edible tissues after IM or IMM administration of procaine penicillin G at extralabel doses. The PBPK model for dairy cows was established and calibrated with pharmacokinetic data of penicillin G through IM and IMM administrations. As the systemic absorption and the drug excretion were both taken into account, the model can be used to predict the concentrations of penicillin G in milk and edible tissues through either the systemic or the local administration. Application of this PBPK model can help predict the extended milk discard intervals after extralabel use of penicillin G. In the field of toxicology, several PBPK models with a milk compartment are available, such as the lactational PBPK models for atrazine in rodents (Lin et al., 2013), perfluorooctanoic acid in humans (Loccisano et al., 2013), manganese in humans (Yoon et al., 2011), and flunixin in dairy cows (Leavens et al., 2014). The milk compartments in all these models were described as virtual compartments and were not physiologically based, which limits the models’ capability to extrapolate to other species. In the present model, a physiologically based approach simulating the mammary gland compartment was established based on the milk production and episodic milking processes to facilitate realistic prediction of drug concentrations in the milk and systemic distribution through the IMM infusion. The Hill-Langmuir equation was adapted from previous milk model (Whittem et al., 2012) to simulate the decrease of milk production due to the increased hydrostatic pressure (Mercer et al., 1970; Smith et al., 2004; Whittem et al., 2012). The 1- and 2-compartment milk models were both adapted to the current PBPK model for penicillin G. The 1-compartment is good enough for the sparse data available from existing pharmacokinetics studies. The 2-compartment model could help predict the drug movements within the 2 compartments of mammary glands, which is useful for the estimation of the drug concentration at the diseased state. The 2-compartment model would help design therapeutic regimens of penicillin G and other veterinary drugs for cows. Both 1- and 2-compartment models provide a basis to simulate pharmacokinetics of penicillin G in mastitis cows. The systemic and local routes of drug administrations, the milk production process including the milk producing rate and milking intervals, and the movement of penicillin G such as the episodic milking process and transfer of penicillin G between the plasma, udder tissue and milk were all considered in the current PBPK model. This model greatly improves our understanding of the pharmacokinetic property of penicillin G in dairy cows, especially in the mammary glands after both IM and IMM administrations. This physiologically based approach simulating the milk compartment can be extrapolated to other species, including rats, humans, and goats. Currently, different statistical methods, such as the “Time To Safe Concentration,” the “Safe Concentration based on Linear Regression” and the “Safe Concentration Per Milking” methods, were used to calculate the withdrawal period for milk (Henri et al., 2017b). However, these methods have various limitations, including the assumption of linearity of the log-transformed tissue/milk depletion data, which may not be true for all drugs (Chevance et al., 2017; Henri et al., 2017b). The nonlinear mixed-effects pharmacokinetic modeling has many advantages (Martin-Jimenez and Riviere, 1998) as a pharmaco-statistical approach to estimate the withdrawal period for milk (Chevance et al., 2017). If the pharmacokinetic data are not available for extralabel use of veterinary drugs, the method of half-life multipliers could be applied to calculate the extended WDIs (Riviere et al., 1998; Smith et al., 2004). However, the nonlinear mixed-effects pharmacokinetic model and the half-life multipliers method are not based on the mechanisms of drug action and the physiology, and cannot be extrapolated to different therapeutic scenarios and animal species. The current PBPK model for milk is a physiologically- and mechanistic-based model and could be extrapolated to the different species, drugs and therapeutic scenarios to predict tissue residues and milk concentrations. The current model is based on the assumption that the 4 quarters of the udder tissues are identical in dairy cows, and no milk production changes were considered for different lactating stages. As the model simulations are about 2 weeks, the potential impact of lactational stages on the milk production was not considered in this model. On the basis of available research (Novak et al., 2009; Yoon et al., 2004), the lactational stages could be incorporated into the model if needed. The parameters of the absorption and excretion processes were determined through model fitting with pharmacokinetics data. These parameters are related to the organic cation transporters and organic anion transporters expressed in the mammary glands of dairy cows (Al-Bataineh et al., 2009). However, the information on the expression of transporters related to the drug transport in mammary glands of dairy cows is limited. In the future, with more research on transporters in food-producing animals, the strategy of in vitro to in vivo extrapolation can be used to calculate the relevant parameter values for drug distribution through mammary glands. For some sensitive parameters related to the milk secretion, such as VMmilksp, s, and F, although their values were validated in an earlier study (Whittem et al., 2012), the experimentally measured values are not available and not commonly reported by pharmacokinetic studies for drugs in milk. These values, whenever available experimentally, would help reduce uncertainties for the current milk model. The binding of penicillin G to milk proteins is not considered in the current model, due to the limited information on protein binding of penicillin G in milk. The current model could be improved when more experimental data on protein binding of penicillin G in milk become available in the future. Also note that the milk is a suspension of fat droplets in an aqueous phase containing proteins, carbohydrates, and electrolytes, in which proteins can extensively bind some antibiotics, resulting in uneven distribution of the drug in different components of the milk (Ziv and Rasmussen, 1975). It has been shown that drug distribution in different fractions of milk can be very different for a systemic route of administration versus a local infusion depending on the drug lipophilicity (Ziv and Rasmussen, 1975). Thus, drug concentrations can be different between different milk samples (whole milk vs skim milk; foremilk vs hindmilk) depending on the drug lipophilicity and sampling methods. At this stage, there are no sufficient data that are granular enough to predict drug concentrations in different components of milk. Therefore, the present model may need some refinements depending on the investigated drugs (lipophilic vs hydrophilic), the dairy products of interest (whole milk vs skim milk), and the sampling methods (cisternal milk vs total milk). Additional studies are needed in order to establish a milk model considering different milk components, which will help model extrapolation to different animal species and different drugs. This is important because the milk compositions are different among different species (eg, 3.3%–5.4% fat in cow milk vs 5.3%–9.0% in buffalo vs 10.2%–21.5% in deer milk) (Claeys et al., 2014) and different drugs have different interactions with milk components. The milk discard interval for penicillin G following labeled IMM infusion predicted by the current model was in accordance with the labeled milk discard time from FDA (98 vs 84 h). The model-predicted milk discard interval following labeled IM injection was longer than the labeled milk discard time by FDA (125 vs 48 h), suggesting that the model-predicted result after labeled IM administration is relatively conservative. The inconsistence between the model-predicted results and the labeled milk discard time after IM administration may be due to different pharmacokinetic studies used in our model compared with the FDA method. The current model was based on a number of datasets from different studies using different drug formulations in different animal populations, thus the results consider the variability of a large diverse population of animals. On the other hand, the labeled milk discard time is typically obtained based on the data from a small set of animals. Overall, the present model results support the recommendations for the appropriate time to harvest milk following the use of penicillin G. On the basis of simulation results, the model-predicted WDIs in the target tissue liver after labeled IMM or IM exposure were consistently much lower than the labeled withdrawal periods. Also, the IMM infusion did not lead to the violations in edible tissues even at 4× label dose. However, major violations of penicillin G were identified in cull dairy cows by USDA. Also, a prolonged elimination of penicillin G in emergency-slaughter dairy cows after IM administration of penicillin G has been reported (Nouws and Ziv, 1978). These inconsistencies may be due to the pharmacokinetic data used for the model calibration of the current model were all from healthy dairy cows. Mastitis and other diseased conditions, as well as the interactions with other drugs in diseased animals may lead to the different distribution and elimination patterns of penicillin G. Due to the fact that pathophysiological data of mastitis in dairy cows are limited, the diseased model for dairy cows with mastitis is not available now. The diseased model is needed to more accurately simulate tissue distribution of penicillin G and to avoid tissue residue violations in cull dairy cows. Thus, the present model cannot be used to predict extended WDIs in target tissue after extralabel use of penicillin G in cull dairy cows. The PBPK model in dairy cows with mastitis would be a future direction for our research when additional data become available. In conclusion, the PBPK model for penicillin G in dairy cows successfully predicts the drug concentrations in the milk and plasma following IM or IMM administration of procaine penicillin G. This is the first reported PBPK model with a physiologically based compartment of mammary glands. By utilizing probabilistic Monte Carlo analysis, the model could help estimate a conservative milk discard time after extralabel use of penicillin G, thereby facilitating food safety assessment of cow-derived milk products. This PBPK framework can be applied to different antibiotics used for dairy cows to help avoid residue violations in the milk and reduce the amount of waste milk. This modeling strategy could also be extrapolated to other species, such as rodents, humans, and goats. SUPPLEMENTARY DATA Supplementary data are available at Toxicological Sciences online. ACKNOWLEDGMENTS The authors would like to acknowledge Dr Pritam Sidhu, Dr Yi-Hsien Cheng, Dr Dongping Zeng, and Shiqiang Jin in the Institute of Computational Comparative Medicine at Kansas State University for helpful discussions. The authors declare no conflict of interest. FUNDING United States Department of Agriculture (USDA) National Institute of Food and Agriculture (NIFA) for the Food Animal Residue Avoidance Databank (FARAD) Program (Award No. 2016-41480-25729 and 2017-41480-27310); The Kansas Bioscience Authority funds to the Institute of Computational Comparative Medicine (ICCM) at Kansas State University. REFERENCES Al-Bataineh M. M. , Van Der Merwe D. , Schultz B. D. , Gehring R. ( 2009 ). Cultured mammary epithelial monolayers (BME-UV) express functional organic anion and cation transporters . J. Vet. Pharmacol. Ther. 32 , 422 – 428 . Google Scholar CrossRef Search ADS PubMed Allen D. B. , DePeters E. J. , Laben R. C. ( 1986 ). Three times a day milking: Effects on milk production, reproductive efficiency, and udder health . J. Dairy Sci. 69 , 1441 – 1446 . Google Scholar CrossRef Search ADS PubMed Andersen M. E. ( 2003 ). Toxicokinetic modeling and its applications in chemical risk assessment . Toxicol. Lett. 138 , 9 – 27 . Google Scholar CrossRef Search ADS PubMed Baynes R. , Riviere J. E. ( 2014 ). Importance of veterinary drug residues. In Strategies for Reducing Drug and Chemical Residues in Food Animals: International Approaches to Residue Avoidance, Management, and Testing ( Baynes R. , Riviere J.E. , Eds.), pp. 1 – 8 . Wiley , New York, NY . Baynes R. E. , Dedonder K. , Kissell L. , Mzyk D. , Marmulak T. , Smith G. , Tell L. , Gehring R. , Davis J. , Riviere J. E. ( 2016 ). Health concerns and management of select veterinary drug residues . Food Chem. Toxicol. 88 , 112 – 122 . Google Scholar CrossRef Search ADS PubMed Brynes S. D. ( 2005 ). Demystifying 21 CFR Part 556–tolerances for residues of new animal drugs in food . Regul. Toxicol. Pharmacol. 42 , 324 – 327 . Google Scholar CrossRef Search ADS PubMed Buur J. , Baynes R. , Smith G. , Riviere J. ( 2006 ). Use of probabilistic modeling within a physiologically based pharmacokinetic model to predict sulfamethazine residue withdrawal times in edible tissues in swine . Antimicrob. Agents Chemother. 50 , 2344 – 2351 . Google Scholar CrossRef Search ADS PubMed Campbell J. R. , Marshall R. T. ( 2016 ). Physiology of Lactation, Dairy Production and Processing: The Science of Milk and Milk Products . Waveland Press , Long Grove, IL . Carruthers V. R. , Davis S. R. , Bryant A. M. , Henderson H. V. , Morris C. A. , Copeman P. J. A. ( 1993 ). Response of Jersey and Friesian cows to once a day milking and prediction of response based on udder characteristics and milk composition . J. Dairy Res. 60 , 1 – 11 . Google Scholar CrossRef Search ADS Cheng Y.-H. , Lin Y.-J. , You S.-H. , Yang Y.-F. , How C. M. , Tseng Y.-T. , Chen W.-Y. , Liao C.-M. ( 2016 ). Assessing exposure risks for freshwater tilapia species posed by mercury and methylmercury . Ecotoxicology 25 , 1181 – 1193 . Google Scholar CrossRef Search ADS PubMed Chevance A. , Jacques A. M. , Laurentie M. , Sanders P. , Henri J. ( 2017 ). The present and future of withdrawal period calculations for milk in the European Union: Focus on heterogeneous, nonmonotonic data . J. Vet. Pharmacol. Ther. 40 , 218 – 230 . Google Scholar CrossRef Search ADS PubMed Chiesa O. A. , Von Bredow J. , Smith M. , Heller D. , Condon R. , Thomas M. H. ( 2006 ). Bovine kidney tissue/biological fluid correlation for penicillin . J. Vet. Pharmacol. Ther. 29 , 299 – 306 . Google Scholar CrossRef Search ADS PubMed Claeys W. L. , Verraes C. , Cardoen S. , De Block J. , Huyghebaert A. , Raes K. , Dewettinck K. , Herman L. ( 2014 ). Consumption of raw or heated milk from different species: An evaluation of the nutritional and potential health benefits . Food Control 42 , 188 – 201 . Google Scholar CrossRef Search ADS Clewell H. J. III , Gentry P. R. , Covington T. R. , Gearhart J. M. ( 2000 ). Development of a physiologically based pharmacokinetic model of trichloroethylene and its metabolites for use in risk assessment . Environ. Health Perspect . 108 , 283 – 305 . Google Scholar CrossRef Search ADS PubMed Covington T. R. , Robinan Gentry P. , Van Landingham C. B. , Andersen M. E. , Kester J. E. , Clewell H. J. III ( 2007 ). The use of Markov chain Monte Carlo uncertainty analysis to support a Public Health Goal for perchloroethylene . Regul. Toxicol. Pharmacol. 47 , 1 – 18 . Google Scholar CrossRef Search ADS PubMed Craigmill A. L. ( 2003 ). A physiologically based pharmacokinetic model for oxytetracycline residues in sheep . J. Vet. Pharmacol. Ther . 26 , 55 – 63 . Google Scholar CrossRef Search ADS PubMed Craigmill A. L. , Riviere J. E. , Webb A. I. ( 2006 ). Tabulation of FARAD Comparative and Veterinary Pharmacokinetic Data . Wiley-Blackwell Publishing , Ames, IA . Davis S. R. , Farr V. C. , Copeman P. J. , Carruthers V. R. , Knight C. H. , Stelwagen K. ( 1998 ). Partitioning of milk accumulation between cisternal and alveolar compartments of the bovine udder: Relationship to production loss during once daily milking . J. Dairy Res. 65 , 1 – 8 . Google Scholar CrossRef Search ADS PubMed Dayan A. D. ( 1993 ). Allergy to antimicrobial residues in food: Assessment of the risk to man . Vet. Microbiol. 35 , 213 – 226 . Google Scholar CrossRef Search ADS PubMed DeWoskin R. S. , Sweeney L. M. , Teeguarden J. G. , Sams R. , Vandenberg J. ( 2013 ). Comparison of PBTK model and biomarker based estimates of the internal dosimetry of acrylamide . Food Chem. Toxicol. 58 , 506 – 521 . Google Scholar CrossRef Search ADS PubMed DrugBank ( 2018 ). Benzylpenicillin. Available at: https://www.drugbank.ca/drugs/DB01053. Last Accessed March 29, 2018. Edwards S. ( 1966 ). Penicillin levels in the milk following intramuscular injection . Vet. Record 78 , 583 – 585 . Google Scholar CrossRef Search ADS Erdman R. A. , Varner M. ( 1995 ). Fixed yield responses to increased milking frequency . J. Dairy Sci. 78 , 1199 – 1203 . Google Scholar CrossRef Search ADS PubMed European Medicines Agency (EMA) ( 2000 ). Guideline on Withdrawal Periods for Milk. European Medicines Agency , London, UK . Available at: http://www.ema.europa.eu/docs/en_GB/document_library/Scientific_guideline/2009/10/WC500004496.pdf. Last Accessed March 29, 2018. European Medicines Agency (EMA) ( 2016 ). Guideline on Approach Towards Harmonisation of Withdrawal Periods. European Medicines Agency , London, UK . Available at: http://www.ema.europa.eu/docs/en_GB/document_library/Scientific_guideline/2016/07/WC500210929.pdf. Last Accessed March 29, 2018. Gehring R. , Smith G. W. ( 2006 ). An overview of factors affecting the disposition of intramammary preparations used to treat bovine mastitis . J. Vet. Pharmacol. Ther. 29 , 237 – 241 . Google Scholar CrossRef Search ADS PubMed Gelman A. , Bois F. , Jiang J. ( 1996 ). Physiological pharmacokinetic analysis using population modeling and informative prior distributions . J. Am. Stat. Assoc. 91 , 1400 – 1412 . Google Scholar CrossRef Search ADS Gionbelli M. P. , Duarte M. S. , Valadares Filho S. C. , Detmann E. , Chizzotti M. L. , Rodrigues F. C. , Zanetti D. , Gionbelli T. R. S. , Machado M. G. ( 2015 ). Achieving body weight adjustments for feeding status and pregnant or non-pregnant condition in beef cows . PLoS One 10 , e0112111. Google Scholar CrossRef Search ADS PubMed Gomes E. R. , Demoly P. ( 2005 ). Epidemiology of hypersensitivity drug reactions . Curr. Opin. Allergy Clin. Immunol. 5 , 309 – 316 . Google Scholar CrossRef Search ADS PubMed Henri J. , Carrez R. , Méda B. , Laurentie M. , Sanders P. ( 2017a ). A physiologically based pharmacokinetic model for chickens exposed to feed supplemented with monensin during their lifetime . J. Vet. Pharmacol. Ther. 40 , 370 – 382 . Google Scholar CrossRef Search ADS Henri J. , Jacques A. M. , Sanders P. , Chevance A. , Laurentie M. ( 2017b ). The present and future of withdrawal period calculations for milk in the European Union: Dealing with data below the limit of quantification . J. Vet. Pharmacol. Ther. 40 , 116 – 122 . Google Scholar CrossRef Search ADS Hogh P. , Rasmussen F. ( 1966 ). Mammary excretion of penicillin fellowing intramuscular injection of Leocillin_1tnR (penethamate hydriodide BAN) and penicillin procaine in cows . Nord. Vet. Med . 18 , 545 – 554 . Holstege D. M. , Puschner B. , Whitehead G. , Galey F. D. ( 2002 ). Screening and mass spectral confirmation of β-lactam antibiotic residues in milk using LC-MS/MS . J. Agric. Food Chem. 50 , 406 – 411 . Google Scholar CrossRef Search ADS PubMed Huang X. , Chen L. , Chen M. , Yuan D. , Nong S. ( 2013 ). Sensitive monitoring of penicillin antibiotics in milk and honey treated by stir bar sorptive extraction based on monolith and LC-electrospray MS detection . J. Sep. Sci. 36 , 907 – 915 . Google Scholar CrossRef Search ADS PubMed Isaksson A. , Arnarp L. ( 1988 ). Quantitative estimation of residual milk in bovine udders—A methodological study . Acta Vet. Scand. 29 , 259 – 262 . Google Scholar PubMed Knappstein K. , Suhren G. , Walte H. ( 2003 ). Influence of milking frequency on withdrawal period after application of beta-lactam antibiotics-based drugs . Anal. Chim. Acta 483 , 241 – 249 . Google Scholar CrossRef Search ADS Knight C. H. , Hirst D. , Dewhurst R. J. ( 1994 ). Milk accumulation and distribution in the bovine udder during the interval between milkings . J. Dairy Res. 61 , 167 – 177 . Google Scholar CrossRef Search ADS PubMed Leavens T. L. , Tell L. A. , Clothier K. A. , Griffith R. W. , Baynes R. E. , Riviere J. E. ( 2012 ). Development of a physiologically based pharmacokinetic model to predict tulathromycin distribution in goats . J. Vet. Pharmacol. Ther. 35 , 121 – 131 . Google Scholar CrossRef Search ADS PubMed Leavens T. L. , Tell L. A. , Kissell L. W. , Smith G. W. , Smith D. J. , Wagner S. A. , Shelver W. L. , Wu H. , Baynes R. E. , Riviere J. E. ( 2014 ). Development of a physiologically based pharmacokinetic model for flunixin in cattle (Bos taurus) . Food Addit. Contam. Part A Chem. Anal. Control Expo. Risk Assess. 31 , 1506 – 1521 . Google Scholar CrossRef Search ADS PubMed Li M. , Gehring R. , Riviere J. E. , Lin Z. ( 2017 ). Development and application of a population physiologically based pharmacokinetic model for penicillin G in swine and cattle for food safety assessment . Food Chem. Toxicol. 107 , 74 – 87 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Fisher J. W. , Ross M. K. , Filipov N. M. ( 2011 ). A physiologically based pharmacokinetic model for atrazine and its main metabolites in the adult male C57BL/6 mouse . Toxicol. Appl. Pharmacol. 251 , 16 – 31 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Fisher J. W. , Wang R. , Ross M. K. , Filipov N. M. ( 2013 ). Estimation of placental and lactational transfer and tissue distribution of atrazine and its main metabolites in rodent dams, fetuses, and neonates with physiologically based pharmacokinetic modeling . Toxicol. Appl. Pharmacol. 273 , 140 – 158 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Gehring R. , Mochel J. P. , Lave T. , Riviere J. E. ( 2016a ). Mathematical modeling and simulation in animal health—Part II: Principles, methods, applications, and value of physiologically based pharmacokinetic modeling in veterinary medicine and food safety assessment . J. Vet. Pharmacol. Ther. 39 , 421 – 438 . Google Scholar CrossRef Search ADS Lin Z. , Jaberi-Douraki M. , He C. , Jin S. , Yang R. S. , Fisher J. W. , Riviere J. E. ( 2017 ). Performance assessment and translation of physiologically based pharmacokinetic models from acslX to Berkeley Madonna, MATLAB, and R language: Oxytetracycline and gold nanoparticles as case examples . Toxicol. Sci. 158 , 23 – 35 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Li M. , Gehring R. , Riviere J. E. ( 2015 ). Development and application of a multiroute physiologically based pharmacokinetic model for oxytetracycline in dogs and humans . J. Pharm. Sci. 104 , 233 – 243 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Vahl C. I. , Riviere J. E. ( 2016b ). Human food safety implications of variation in food animal drug metabolism . Sci. Rep. 6 , 27907 . Google Scholar CrossRef Search ADS Liu C. J. , Wang H. , Jiang Y. B. , Du Z. X. ( 2011 ). Rapid and simultaneous determination of amoxicillin, penicillin G, and their major metabolites in bovine milk by ultra-high-performance liquid chromatography-tandem mass spectrometry . J. Chromatogr. B 879 , 533 – 540 . Google Scholar CrossRef Search ADS Loccisano A. E. , Longnecker M. P. , Campbell J. L. , Andersen M. E. , Clewell H. J. III ( 2013 ). Development of PBPK models for PFOA and PFOS for human pregnancy and lactation life stages . J. Toxicol. Environ. Health Part A 76 , 25 – 57 . Google Scholar CrossRef Search ADS PubMed Lumen A. , McNally K. , George N. , Fisher J. W. , Loizou G. D. ( 2015 ). Quantitative global sensitivity analysis of a biologically based dose-response pregnancy model for the thyroid endocrine system . Front. Pharmacol . 6 , 107. Google Scholar CrossRef Search ADS PubMed Lyons N. A. , Kerrisk K. L. , Garcia S. C. ( 2014 ). Milking frequency management in pasture-based automatic milking systems: A review . Livest. Sci. 159 , 102 – 116 . Google Scholar CrossRef Search ADS Martin-Jimenez T. , Riviere J. E. ( 1998 ). Population pharmacokinetics in veterinary medicine: Potential use for therapeutic drug monitoring and prediction of tissue residues . J. Vet. Pharmacol. Ther. 21 , 167 – 189 . Google Scholar CrossRef Search ADS PubMed McNally K. , Cotton R. , Loizou G. D. ( 2011 ). A workflow for global sensitivity analysis of PBPK models . Front. Pharmacol . 2 , 31. Google Scholar CrossRef Search ADS PubMed Mercer H. , Geleta J. , Carter G. ( 1974 ). Absorption and excretion of penicillin G from the mastitic bovine udder . J. Am. Vet. Med. Assoc . 164 , 613 – 617 . Google Scholar PubMed Mercer H. , Geleta J. , Schultz E. , Wright W. ( 1970 ). Milk-out rates for antibiotics in intramammary infusion products used in the treatment of bovine mastitis: Relationship of somatic cell counts, milk production level, and drug vehicle . Am. J. Vet. Res. 31 , 1549 – 1560 . Google Scholar PubMed National Research Council (NRC) ( 1999 ). The use of drugs in food animals. In The Use of Drugs in Food Animals: Benefits and Risks ( National Research Council , Ed.), pp. 110 – 141 . National Academy Press , Washington, DC . Nouws J. F. , Ziv G. ( 1978 ). Tissue distribution and residues of benzylpenicillin and aminoglycoside antibiotics in emergency-slaughtered ruminants . Tijdschr. Diergeneeskd . 103 , 140 – 151 . Google Scholar PubMed Novak P. , Vokralova J. , Broucek J. ( 2009 ). Effects of the stage and number of lactation on milk yield of dairy cows kept in open barn during high temperatures in summer months . Archiv Tierzucht 52 , 574 – 586 . Paige J. C. , Tollefson L. , Miller M. ( 1997 ). Public health impact on drug residues in animal tissues . Vet. Hum. Toxicol . 39 , 162 – 169 . Google Scholar PubMed Papich M. G. ( 1987 ). The beta-lactam antibiotics: Clinical pharmacology and recent developments . Compend. Contin. Educ. Pract. Vet . 9 , 68 – 76 . Papich M. G. , Riviere J. E. ( 2009 ). Lactam antibiotics: penicillins, cephalosporins, and related drugs. In Veterinary Pharmacology and Therapeutics. J.E. Riviere, and M.G. Papich, Eds., 9th ed., pp. 866–888. Wiley, Ames, IA. Papich M. G. , Korsrud G. O. , Boison J. O. , Yates W. D. , MacNeil J. D. , Janzen E. D. , Cohen R. D. , Landry D. A. ( 1993 ). A study of the disposition of procaine penicillin G in feedlot steers following intramuscular and subcutaneous injection . J. Vet. Pharmacol. Ther. 16 , 317 – 327 . Google Scholar CrossRef Search ADS PubMed Payne M. A. , Craigmill A. , Riviere J. E. , Webb A. I. ( 2006 ). Extralabel use of penicillin in food animals . J. Am. Vet. Med. Assoc. 229 , 1401 – 1403 . Google Scholar CrossRef Search ADS PubMed Portis E. , Lindeman C. , Johansen L. , Stoltman G. ( 2012 ). A ten-year (2000-2009) study of antimicrobial susceptibility of bacteria that cause bovine respiratory disease complex–Mannheimia haemolytica, Pasteurella multocida, and Histophilus somni—In the United States and Canada . J. Vet. Diagn. Invest. 24 , 932 – 944 . Google Scholar CrossRef Search ADS PubMed Raison-Peyron N. , Messaad D. , Bousquet J. , Demoly P. ( 2001 ). Anaphylaxis to beef in penicillin-allergic patient . Allergy 56 , 796 – 797 . Google Scholar CrossRef Search ADS PubMed Randall W. , Durbin C. , Wilner J. , Collins J. ( 1954 ). Antibiotic Concentration and Duration in Animal Tissue and Body Fluids. Antibiotics Annual 1953-1954 . Medical Encyclopedia, Inc ., New York, NY . Rasmussen F. ( 1959 ). Mammary excretion of benzylpenicillin, erythromycin, and penethamate hydroiodide . Acta Pharmacol. Toxicol. 16 , 194 – 200 . Google Scholar CrossRef Search ADS Riviere J. E. , Craigmill A. L. , Sundlof S. F. ( 1986 ). Food Animal Residue Avoidance Databank (FARAD): An automated pharmacologic databank for drug and chemical residue avoidance . J. Food Prot. 49 , 826 – 830 . Google Scholar CrossRef Search ADS Riviere J. E. , Tell L. A. , Baynes R. E. , Vickroy T. W. , Gehring R. ( 2017 ). Guide to FARAD resources: Historical and future perspectives . J. Am. Vet. Med. Assoc. 250 , 1131 – 1139 . Google Scholar CrossRef Search ADS PubMed Riviere J. E. , Webb A. I. , Craigmill A. L. ( 1998 ). Primer on estimating withdrawal times after extralabel drug use . J. Am. Vet. Med. Assoc. 213 , 966 – 968 . Google Scholar PubMed Schadewinkel-Scherkl A. M. , Rasmussen F. , Merck C. C. , Nielsen P. , Frey H. H. ( 1993 ). Active transport of benzylpenicillin across the blood-milk barrier . Pharmacol. Toxicol. 73 , 14 – 19 . Google Scholar CrossRef Search ADS PubMed Shankaran H. , Adeshina F. , Teeguarden J. G. ( 2013 ). Physiologically-based pharmacokinetic model for Fentanyl in support of the development of Provisional Advisory Levels . Toxicol. Appl. Pharmacol. 273 , 464 – 476 . Google Scholar CrossRef Search ADS PubMed Smith G. W. , Gehring R. , Riviere J. E. , Yeatts J. L. , Baynes R. E. ( 2004 ). Elimination kinetics of ceftiofur hydrochloride after intramammary administration in lactating dairy cows . J. Am. Vet. Med. Assoc. 224 , 1827 – 1830 . Google Scholar CrossRef Search ADS PubMed Sterner T. R. , Ruark C. D. , Covington T. R. , Yu K. O. , Gearhart J. M. ( 2013 ). A physiologically based pharmacokinetic model for the oxime TMB-4: Simulation of rodent and human data . Arch. Toxicol. 87 , 661 – 680 . Google Scholar CrossRef Search ADS PubMed Suhren G. ( 1996 ). Influence of residues of antimicrobials in milk on commercially applied starter cultures-model trials . Kieler. Milchw. Forsch. 48 , 131 – 149 . Tan Y. M. , Liao K. H. , Conolly R. B. , Blount B. C. , Mason A. M. , Clewell H. J. III ( 2006 ). Use of a physiologically based pharmacokinetic model to identify exposures consistent with human biomonitoring data for chloroform . J. Toxicol. Environ. Health A 69 , 1727 – 1756 . Google Scholar CrossRef Search ADS PubMed United States Department of Agriculture (USDA) ( 2015 ). U.S. National Residue Program: 2014 Residue Sample Results. United States Department of Agriculture , Washington, DC . Available at: http://www.fsis.usda.gov/wps/wcm/connect/2428086b-f8ec-46ed-8531-a45d10bfef6f/2014-Red-Book.pdf?MOD=AJPERES. Last Accessed March 29, 2018. United States Department of Agriculture (USDA) ( 2017a ). U.S. National Residue Program: 2015 Residue Sample Results . United States Department of Agriculture , Washington, DC . Available at: https://www.fsis.usda.gov/wps/wcm/connect/f57333e5-9ff8-43ed-a787-6824f44bbac4/2015-red-book.pdf?MOD=AJPERES. Last Accessed March 29, 2018. United States Department of Agriculture (USDA) ( 2017b ). U.S. National Residue Program: 2016 Residue Sample Results . United States Department of Agriculture , Washington, DC . Available at: https://www.fsis.usda.gov/wps/wcm/connect/d84a5cac-5b4e-4e60-85b4-8886d0dc1660/2016-Red-Book.pdf?MOD=AJPERES. Last Accessed March 29, 2018. United States Pharmacopeia ( 2003 ). Penicillin G Veterinary—Intramammary-Local . United States Pharmacopeia , Rockville, MD . Available at: www.aavpt.org/resource/resmgr/imported/penicillinG-in.pdf. Last Accessed March 29, 2018. PubMed PubMed U.S. Food and Drug Administration (FDA) ( 2005 ). M-I-05-5: Tolerance and/or Safe Levels of Animal Drug Residues in Milk . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://wayback.archive-it.org/7993/20170113001136/http://www.fda.gov/Food/GuidanceRegulation/GuidanceDocumentsRegulatoryInformation/Milk/ucm077350.htm. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2006 ). Guidance for Industry—General Principles for Evaluating the Safety of Compounds Used in Food-Producing Animals . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/ohrms/dockets/98fr/05d-0219-gdl0001.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2010 ). Freedom of Information . NADA 065-010. Norocillin. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://wayback.archive-it.org/7993/20170406083042/https://www.fda.gov/downloads/AnimalVeterinary/Products/ApprovedAnimalDrugProducts/FOIADrugSummaries/UCM218419.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2012 ). Regulation of Carcinogenic Compounds Used in Food-Producing Animals. Code of Federal Regulations Title 21, Parts 500.82. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.accessdata.fda.gov/scripts/cdrh/cfdocs/cfCFR/CFRSearch.cfm?fr=500.82. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2013 ). Penicillin G Procaine Implantation and Injectable Dosage Forms. 21 CFR 522.1696. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.gpo.gov/fdsys/pkg/CFR-2013-title21-vol6/pdf/CFR-2013-title21-vol6-chapI-subchapE.pdf. Last Accessed March 29, 2018. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2015a ). Pasteurized Milk Ordinance 2015 Revision . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/downloads/food/guidanceregulation/guidancedocumentsregulatoryinformation/milk/ucm513508.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2015b ). Milk Drug Residue Sampling Survey . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/downloads/animalveterinary/guidancecomplianceenforcement/complianceenforcement/ucm435759.pdf. Last Accessed March 29, 2018. PubMed PubMed U.S. Food and Drug Administration (FDA) ( 2016 ). Guidance for Industry—General Principles for Evaluating the Human Food Safety of New Animal Drugs Used In Food-Producing Animals. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/downloads/AnimalVeterinary/GuidanceComplianceEnforcement/GuidanceforIndustry/ucm052180.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2017 ). Tolerances for Residues of New Animal Drugs in Food . Code of Federal Regulations Title 21, Parts 556.510. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.accessdata.fda.gov/scripts/cdrh/cfdocs/cfcfr/CFRSearch.cfm?fr=556.510. Last Accessed March 31, 2018. Van OS J. , Burtelaar J. , Gouldsward J. ( 1974 ). Intramuscular treatment of bovine mastitis with various penicillins. Penicillin concentrations in the milk . Tijdschr. Diergeneeskd . 99 , 114 – 122 . Vetharaniam I. , Davis S. R. , Soboleva T. K. , Shorten P. R. , Wake G. C. ( 2003 ). Modeling the interaction of milking frequency and nutrition on mammary gland growth and lactation . J. Dairy Sci. 86 , 1987 – 1996 . Google Scholar CrossRef Search ADS PubMed Vilim A. , Larocque L. , Macintosh A. ( 1979 ). Depletion of brilliant blue F.C.F. and penicillin G in milk from treated cows . J. Food Prot. 42 , 491 – 494 . Google Scholar CrossRef Search ADS Vilim A. , Larocque L. , Macintosh A. ( 1980 ). Depletion of brilliant blue F.C.F., penicillin G, and dihydrostreptomycin in milk treated does with experimentally induced mastitis . J. Food Prot. 43 , 356 – 359 . Google Scholar CrossRef Search ADS Vogel G. , Nicolet J. , Martig J. , Tschudi P. , Meylan M. ( 2001 ). Pneumonia in calves: Characterization of the bacterial spectrum and the resistance patterns to antimicrobial drugs . Schweiz. Arch. Tierheilkd. 143 , 341 – 350 . Google Scholar PubMed Whittem T. , Whittem J. H. , Constable P. D. ( 2012 ). Modelling the concentration-time relationship in milk from cattle administered an intramammary drug . J. Vet. Pharmacol. Ther. 35 , 460 – 471 . Google Scholar CrossRef Search ADS PubMed Wishart D. S. , Knox C. , Guo A. C. , Shrivastava S. , Hassanali M. , Stothard P. , Chang Z. , Woolsey J. ( 2006 ). DrugBank: A comprehensive resource for in silico drug discovery and exploration . Nucleic Acids Res . 34 , D668 – D672 . Google Scholar CrossRef Search ADS PubMed World Health Organization (WHO) ( 2010 ). Characterization and Application of Physiologically Based Pharmacokinetic Models in Risk Assessment . World Health Organization, International Programme on Chemical Safety , Geneva, Switzerland . Yang F. , Liu H. W. , Li M. , Ding H. Z. , Huang X. H. , Zeng Z. L. ( 2012 ). Use of a Monte Carlo analysis within a physiologically based pharmacokinetic model to predict doxycycline residue withdrawal time in edible tissues in swine . Food Addit. Contam. Part A Chem. Anal. Control Expo. Risk Assess . 29 , 73 – 84 . Google Scholar CrossRef Search ADS PubMed Yang X. , Doerge D. R. , Teeguarden J. G. , Fisher J. W. ( 2015 ). Development of a physiologically based pharmacokinetic model for assessment of human exposure to bisphenol A . Toxicol. Appl. Pharmacol. 289 , 442 – 456 . Google Scholar CrossRef Search ADS PubMed Yoon J. , Lee L. , Kim C. , Chung Y. , Kim C. ( 2004 ). Effects of milk production, season, parity and lactation period on variations of milk urea nitrogen concentration and milk components of Holstein dairy cows . Asian-Australas. J. Anim. Sci. 17 , 479 – 484 . Google Scholar CrossRef Search ADS Yoon M. , Nong A. , Clewell H. J. III , Taylor M. D. , Dorman D. C. , Andersen M. E. ( 2009 ). Evaluating placental transfer and tissue concentrations of manganese in the pregnant rat and fetuses after inhalation exposures with a PBPK model . Toxicol. Sci. 112 , 44 – 58 . Google Scholar CrossRef Search ADS PubMed Yoon M. , Schroeter J. D. , Nong A. , Taylor M. D. , Dorman D. C. , Andersen M. E. , Clewell H. J. III ( 2011 ). Physiologically based pharmacokinetic modeling of fetal and neonatal manganese exposure in humans: Describing manganese homeostasis during development . Toxicol. Sci. 122 , 297 – 316 . Google Scholar CrossRef Search ADS PubMed Zeng D. , Lin Z. , Fang B. , Li M. , Gehring R. , Riviere J. E. , Zeng Z. ( 2017 ). Pharmacokinetics of mequindox and its marker residue 1, 4-bisdesoxymequindox in swine following multiple oral gavage and intramuscular administration: An experimental study coupled with population physiologically based pharmacokinetic modeling . J. Agric. Food Chem. 65 , 5768 – 5777 . Google Scholar CrossRef Search ADS PubMed Zhu X. , Huang L. , Xu Y. , Xie S. , Pan Y. , Chen D. , Liu Z. , Yuan Z. ( 2017 ). Physiologically based pharmacokinetic model for quinocetone in pigs and extrapolation to mequindox . Food Addit. Contam. Part A 34 , 192 – 210 . Google Scholar CrossRef Search ADS Ziv G. , Rasmussen F. ( 1975 ). Distribution of labeled antibiotics in different components of milk following intramammary and intramuscular administrations . J. Dairy Sci . 58 , 938 – 946 . Google Scholar CrossRef Search ADS PubMed Ziv G. , Sulman F. G. ( 1974 ). Effects of probenecid on the distribution, elimination, and passage into milk of benzylpenicillin, ampicillin and cloxacillin . Arch. Int. Pharmacodyn. Ther. 207 , 373 – 382 . Google Scholar PubMed © The Author(s) 2018. Published by Oxford University Press on behalf of the Society of Toxicology. All rights reserved. For permissions, please e-mail: journals.permissions@oup.com This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/about_us/legal/notices) http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Toxicological Sciences Oxford University Press

Probabilistic Physiologically Based Pharmacokinetic Model for Penicillin G in Milk From Dairy Cows Following Intramammary or Intramuscular Administrations

Loading next page...
 
/lp/ou_press/probabilistic-physiologically-based-pharmacokinetic-model-for-xF0SuRy8Cp
Publisher
Oxford University Press
Copyright
© The Author(s) 2018. Published by Oxford University Press on behalf of the Society of Toxicology. All rights reserved. For permissions, please e-mail: journals.permissions@oup.com
ISSN
1096-6080
eISSN
1096-0929
D.O.I.
10.1093/toxsci/kfy067
Publisher site
See Article on Publisher Site

Abstract

Abstract Penicillin remains one of the most frequently identified violative drug residues in food-producing animals. The predominant violations of penicillin were found in cull dairy cows. In the United States, procaine penicillin G is approved to be used in dairy cows through intramuscular (IM) and intramammary (IMM) administrations. Physiologically based pharmacokinetic (PBPK) models are useful tools to predict withdrawal intervals and tissue residues of drugs in food animals to ensure food safety, especially for extralabel drug use due to the scarcity of experimental data after extralabel administrations. Currently, no PBPK model is available to predict penicillin concentrations in milk. A population PBPK model with a physiologically based compartment for the mammary gland was established for penicillin G in dairy cows. The model predicted the tissue and milk residues well based on comparison with data from previous pharmacokinetic studies. The predicted milk discard interval of procaine penicillin G administered at 10 times the label dose for 3 repeated IM administrations was 182 h, and 122 h at 4 times the label dose after 3 repeated IMM infusions. Predicted results showed that even 4 times label dose did not lead to violative tissue residues in healthy dairy cows with IMM infusions. The predominant violations found in cull dairy cows may be caused by altered pharmacokinetics due to mastitis, other diseases, and/or interactions with other drugs, which have impacts on penicillin distribution and elimination. The current PBPK model can help predict milk discard interval for penicillin following extralabel use through IM and IMM administrations. PBPK modeling, dairy cows, food safety, drug residue, Food Animal Residue Avoidance Databank (FARAD), penicillin INTRODUCTION Penicillin is one of the top 3 most commonly detected violative residues in animal-derived foods reported by the National Residue Program of United States Department of Agriculture (USDA) from 2014 to 2016 (USDA, 2015, 2017a,b). This is mainly because penicillin G is one of the most widely used antimicrobials and tested regularly (FDA, 2013; Portis et al., 2012; Vogel et al., 2001). Extralabel use of penicillin G in food-producing animals is very common (Chiesa et al., 2006) due to reduced sensitivity to penicillin in bacteria. The extralabel (also called off-label) use of veterinary drugs means to use a drug in a manner that is different from the approved labeling under the supervision of a veterinarian, and extralabel use may only occur if prescribed by a veterinarian (FDA, 2017). Consumption of beef or pork products containing violative penicillin residues can lead to anaphylactic reactions for sensitive individuals (Dayan, 1993; Gomes and Demoly, 2005; Raison-Peyron et al., 2001). In addition, violative residues of penicillin G in milk could interfere with starter cultures for fermented dairy products (Payne et al., 2006; Suhren, 1996) and may result in suspension of the producer’s permit or certification (NRC, 1999). Drug residues above the regulatory safe level in animal-derived products challenge the global food safety (Baynes and Riviere, 2014; Baynes et al., 2016; Paige et al., 1997). The U.S. Food and Drug Administration (FDA) has established a zero tolerance for penicillin in milk products (FDA, 2017). In actual practice, FDA uses the “safe level” of 5 ng/ml for penicillin residues in milk. The USDA National Residue Program reported that 135 out of 213 (63%) and 153 out of 216 (71%) violative samples in cull dairy cows were found with penicillin G violations in 2015 and 2016, respectively (USDA, 2017a,b). The predominant violations of penicillin found in cull dairy cows lead to the concerns of potential drug residues in the cow-derived food products, including milk. In the United States, procaine penicillin G is approved to be used in dairy cows through intramuscular (IM) and intramammary (IMM) administrations. The approved doses of procaine penicillin G are at a daily dose of 6600 IU/kg of body weight (6.5 mg/kg) for no more than 7 consecutive days through the IM injection (Papich et al., 1993), and 100 000 IU (99 mg) into each affected quarter of udders every 12 h for no more than 3 doses (United States Pharmacopeia, 2003). The procaine penicillin G is slowly absorbed after IM administration (Papich, 1987). The major metabolite of penicillin G, penicilloic acid, does not have the antimicrobial activity and accounts for about 16%–30% of an IM dose of penicillin G in humans (DrugBank, 2018; Wishart et al., 2006). Unlike in humans, the metabolism of penicillin G in food-producing animals is not well characterized, and the elimination of penicillin G is mainly through renal secretion (Papich and Riviere, 2009). The determination of withdrawal intervals (WDIs) for extralabel use of veterinary drugs is important to avoid the violations of tissue residues. The WDI or milk discard time is the time for drug residues in the edible tissues or milk to deplete below concentrations that are considered safe for human consumptions (FDA, 2012). The milk produced before the milk discard time must be discarded to avoid drug residues in milk (EMA, 2000). According to the Pasteurized Milk Ordinance (FDA, 2015a), milk violation tests carried out by FDA are based on the samples from bulk milk tank. In the United States, the 99th percentile tolerance limit with 95% confidence is used to determine the withdrawal period and milk discard time (FDA, 2006, 2016). In Europe, it is recommended to determine the withdrawal period with the upper one-sided 95% tolerance limit for the drug residue below the maximum residue limit with 95% confidence (EMA, 2016). Due to the lack of pharmacokinetic data after extralabel use of veterinary drugs, it is difficult to use statistical tolerance methods or classic pharmacokinetic methods to determine the WDIs after extralabel administrations. The physiologically based pharmacokinetic (PBPK) model is a useful tool to predict tissue residues and WDIs of veterinary drugs for food safety assessment (Henri et al., 2016a; Lin et al., 2016b). PBPK modeling is a mechanism-based approach to simulate the absorption, distribution, metabolism, and elimination of chemicals (Andersen, 2003; Lin et al., 2016a; WHO, 2010) and to predict the concentrations in the edible tissues and products beyond times when data are not available. By changing the related parameter values, the PBPK model can be extrapolated to different species and therapeutic scenarios. There is no whole-body PBPK model for predicting the milk discard intervals of drugs after extralabel use. Recently, we published a multiroute population PBPK model for penicillin G in beef cattle (Li et al., 2017). On the basis of cattle model, in the present study we developed a PBPK model with physiologically based compartment of mammary glands to predict milk concentrations and discard intervals after IM or IMM administration of penicillin G at extralabel doses in dairy cows. The systemic absorption from the mammary gland and the distribution of penicillin G were considered. The volume change by milk secretion and episodic milking process were also included. Probabilistic analysis using Monte Carlo simulation was applied to address the uncertainties of parameter values. This PBPK model contributes to the field of toxicological science and biological modeling by establishing a quantitative tool that allows to conduct drug residue safety assessment of cow-derived food products and creating a new physiologically based approach to simulate milk production and episodic milking processes, which can be applied to other species, including rodents, humans, and goats. MATERIALS AND METHODS Data source for model calibration The pharmacokinetic data used in the calibration and evaluation of the PBPK model were acquired from the Comparative Pharmacokinetic Database from the Food Animal Residue Avoidance Databank (FARAD). FARAD is supported by USDA National Institute of Food and Agriculture with an overarching goal to provide veterinary practitioners the most current and accurate information to facilitate the production of safe foods of animal origin through the prevention and mitigation of violative chemical residues in food animal products (Craigmill et al., 2006; Riviere et al., 1986, 2017). Pharmacokinetic data in dairy cows after IM or IMM administration of procaine penicillin G were selected. The formulation of aqueous suspension was used for IM injection of procaine penicillin G. For IMM injection, the formulation of procaine penicillin G in peanut oil, which is more relevant to the current clinical formulations, was used. Brief description and key information of selected pharmacokinetic studies are summarized in Table 1. The graphic pharmacokinetic data were extracted from selected studies using WebPlotDigitizer (version 3.10, https://automeris.io/WebPlotDigitizer/; last accessed March 29, 2018.). Table 1. Summary of Pharmacokinetic Studies of Penicillin G in Dairy Cows Used for Calibration and Evaluation of the PBPK Model PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg All data used for model calibration and evaluation were from healthy cows. The abbreviations for routes: IM, intramuscular injection; IMM, intramammary infusion. The abbreviations for matrices: P, plasma; M, milk. The abbreviations for assays: HPLC, high performance liquid chromatography; MM, microbiological methods. NA, not available or not applicable. Table 1. Summary of Pharmacokinetic Studies of Penicillin G in Dairy Cows Used for Calibration and Evaluation of the PBPK Model PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg PK Studies/Purposes Routes n Matrices Analytical Methods Treated Quarters Dose Regimen Dose per Administration Calibration  Mercer et al. 1974 IMM 6 P, M MM All One infusion 396.5 mg per quarter  Hogh and Rasmussen 1966 IM 13 P, M NA NA One injection 10 mg/kg  Van OS et al. 1974 IM 4 P, M MM NA One injection 40 mg/kg Evaluation  Vilim et al. 1979 IMM 6 M MM Two quarters 24-h interval 3 doses 100 mg per quarter  Vilim et al. 1980 IMM 6 M MM All 24-h interval 3 doses 100 mg per quarter  Knappstein et al. 2003 IMM 4 to 6 M HPLC All 24-h interval 3 doses 1898 mg per quarter  Randall et al. 1954 IM 2 P, M MM NA One injection 11 mg/kg  Edwards 1966 IM 2 P, M MM NA One injection 2.2, 4.4 or 11 mg/kg  FDA 2010 (NADA 065-010) IM 20 M HPLC NA 24-h interval 4 doses 6.5 mg/kg All data used for model calibration and evaluation were from healthy cows. The abbreviations for routes: IM, intramuscular injection; IMM, intramammary infusion. The abbreviations for matrices: P, plasma; M, milk. The abbreviations for assays: HPLC, high performance liquid chromatography; MM, microbiological methods. NA, not available or not applicable. Model structure The PBPK model for penicillin G with the mammary gland compartment was established based on our pervious generic PBPK model for penicillin G in beef cattle (Li et al., 2017). The model consisted of 8 compartments corresponding to different tissues in the body, including liver, kidney, muscle, fat, lung, the rest of body and udder connected by the circulating blood system (Figure 1A). For food safety assessment purposes, all the major edible tissues, including liver, kidney, muscle, and fat, were included. The udder was considered as a compartment, because the model considered the absorption of penicillin G from milk to the systemic circulation and the excretion of penicillin G into the milk. The lung was also included as a compartment as it is the therapeutic target for penicillin G. Each compartment was defined by a tissue weight and tissue blood flow rate. The compartments for urine and feces were established as virtual excretory compartments without volume changes. The flow-limited model, which performed well for the previous PBPK model of penicillin G, was applied in the current model. Figure 1. View largeDownload slide A schematic diagram of the physiologically based pharmacokinetic (PBPK) model for penicillin G in dairy cows and the 1- and 2-compartment udder models. A, The PBPK model for penicillin G in dairy cows. Two different administration routes including intramuscular (IM) and intramammary (IMM) administrations are presented in the model. The dash line for the compartment of udder indicates that the volume changes since the preceding milking and this could be simulated using the Hill-Langmuir equation. B, The 1-compartment model for the udder. The milk space represents the combined milk space including alveoli, ducts, and the cistern. C, The 2-compartment model for the udder. Figure 1. View largeDownload slide A schematic diagram of the physiologically based pharmacokinetic (PBPK) model for penicillin G in dairy cows and the 1- and 2-compartment udder models. A, The PBPK model for penicillin G in dairy cows. Two different administration routes including intramuscular (IM) and intramammary (IMM) administrations are presented in the model. The dash line for the compartment of udder indicates that the volume changes since the preceding milking and this could be simulated using the Hill-Langmuir equation. B, The 1-compartment model for the udder. The milk space represents the combined milk space including alveoli, ducts, and the cistern. C, The 2-compartment model for the udder. Model calibration and parameterization The current model was based on the published PBPK model for penicillin G in beef cattle, and the physiologically based mammary gland compartment and IMM infusion were added to the model. The physiological parameters related to dairy cows were acquired from previous research. The average values of blood flow to the udder (QUC), the tissue volume fractions of the udder (VUC), venous blood volume (VvenC), and arterial blood volume (VartC) were obtained or calculated based on previous experimental studies (Campbell and Marshall, 2016; Gionbelli et al., 2015). For chemical-specific parameters (eg, partition coefficients) of penicillin G, the original values used in the cattle model (Li et al., 2017) were applied in the current dairy cow model. The parameters related to the udder compartment and IMM infusions, including VmaxC, Km, and KabC were fitted with pharmacokinetic data chosen for calibration (Table 1). VmaxC and Km were optimized by fitting with calibration datasets following IM administrations (Hogh and Rasmussen, 1966; Van OS et al. 1974), and KabC was estimated by fitting to the dataset after IMM administration (Mercer et al., 1974). The values of all physiological parameters and chemical-specific parameters used in the PBPK model are shown in Table 2. Table 2. Parameter Values Used in the PBPK Model for Penicillin G in Milk From Dairy Cows Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 “Model fitted” indicates that the parameter values were determined by fitting the model with calibration pharmacokinetic data. Table 2. Parameter Values Used in the PBPK Model for Penicillin G in Milk From Dairy Cows Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 Parameter Abbreviation Mean Reference Cardiac output (l/h/kg) QCC 5.97 Li et al. 2017 Tissue volume (fraction of body weight, unitless)  Arterial blood (lactating cows) VartC 0.012 Gionbelli et al. 2015  Venous blood (lactating cows) VvenC 0.037 Gionbelli et al. 2015  Liver VLC 0.014 Li et al. 2017  Kidney VKC 0.002 Li et al. 2017  Muscle VMC 0.27 Li et al. 2017  Fat VFC 0.15 Li et al. 2017  Lung VLuC 0.008 Li et al. 2017  Udder VUC 0.008 Gionbelli et al. 2015  Rest of body VrestC 0.507 Total adds to 1 Blood flow (fraction of cardiac output, unitless)  Liver QLC 0.405 Li et al. 2017  Kidney QKC 0.09 Li et al. 2017  Muscle QMC 0.18 Li et al. 2017  Fat QFC 0.08 Li et al. 2017  Udder QUC 0.081 Campbell and Marshall, 2016  Rest of body QrestC 0.164 Total adds to 1 The maximum volume of milk space (l) VMmilksp 26.8 Whittem et al., 2012 Time to reach half of the maximum volume of milk space (h) S 23.2 Whittem et al., 2012 Milk residue ratio after milking (unitless) F 0.04 Whittem et al., 2012 Absorption rate constant (/h) Intramuscular Absorption rate constant of IM administration (/h) Kim 0.07 Li et al. 2017 Fraction of undissolved procaine penicillin G (unitless) Frac 0.6 Li et al. 2017 Dissolution rate constant into penicillin G moieties (/h) Kdiss 1.00E-05 Li et al. 2017 Intramammary Maximum velocity for penicillin G excretion in the mammary gland scaling to BW (mg/h/kg BW) VmaxC 0.0022 Model fitted Concentration of penicillin to reach half Vmax (mg/l) Km 0.7 Model fitted Rate constant of penicillin G absorbed from milk (/h) KabC 1.00E-04 Model fitted Tissue: plasma partition coefficient for the parent drug (unitless)  Liver PL 3 Li et al. 2017  Kidney PK 2.5 Li et al. 2017  Muscle PM 0.3 Li et al. 2017  Fat PF 0.04 Li et al. 2017  Lung PLu 0.18 Li et al. 2017  Udder PU 0.2 Model fitted  Rest of body Prest 0.479 Li et al. 2017 Hepatic metabolic rate constant (/h/kg) KmetC 0.0025 Li et al. 2017 Percentage of plasma protein binding (unitless) PB 0.483 Li et al. 2017 Urinary elimination rate constant (l/h/kg) KurineC 0.45 Li et al. 2017 “Model fitted” indicates that the parameter values were determined by fitting the model with calibration pharmacokinetic data. Establishment and validation of the udder compartment The compartments representing the milk and udder were established according to the physiologically based mathematical model described by Whittem et al. (2012). In brief, the Hill-Langmuir equation was used to describe the volume-time relationship for milk production. The concentration change of penicillin G due to milk production and episodic milk emptying was considered. Both the 1-compartment and 2-compartment models were created and validated using previously reported data (Davis et al., 1998). The 1-compartment model was established considering the milk space in the udder as 1 compartment with the function of milk secretion and storage (Figure 1B). The 2-compartment model divided the milk space into the alveoli, which serve the function of milk secretion and storage, and the cistern, which only has the function for milk storage (Figure 1C). For the 2-compartment model, the distribution of penicillin G between alveoli and cistern compartments was incorporated. Berkeley Madonna (Version 8.3.23.0; University of California at Berkeley, California) was used to develop the PBPK model with the mammary gland compartment and run all simulations. Codes of the model are provided in the Supplementary Data and will be available on our website (http://iccm.k-state.edu/; last accessed March 29, 2018). As the model was constructed based on the previous PBPK model, only key and new mathematical equations are described in detail in following paragraphs. For all the other equations, details have been provided in the Supplementary Data and the previously published PBPK model for penicillin G in cattle and swine (Li et al., 2017). For the 1-compartment udder model, the volume change of the milk space in udder was modeled the same as the volume change of the milk production. The milk volume in the udder compartment depends on the episodic milking and the milk production between milkings. The milk production was described using the Hill-Langmuir equation (equation 1) reported in the previous milk model (Whittem et al., 2012). According to the previous experimental data (Davis et al., 1998), the equation describes the milk production rises rapidly in a linear trend at the beginning, and the rate of milk secretion slows down with the increased hydrostatic pressure in the udder compartment. The volume of residual milk in udder can be determined by the fraction of maximum volume of the milk space (equation 2). Previous experimental research reported the residual milk fractions are in the range of 4%–14% (Carruthers et al., 1993; Isaksson and Arnarp, 1988; Knight et al., 1994), here 4% was used as the value for the residual fraction as reported in the previous model (Whittem et al., 2012). The total volume of milk space should account for the residual volume after milking and the volume of milk secreted after the preceding milking (equation 3). Vmilk=tm×VMmilksptm+s, (1) Vresidue=VMmilksp×F, (2) Vmilksp=Vresidue+Vmilk, (3) where Vmilk is the volume of milk produced since the preceding milking (l); tm is the time after the preceding milking (h); VMmilksp is the maximum volume of the milk space including spaces in alveoli, ducts, and the cistern (l); s is the time for Vmilk reaching half of VMmilksp (h); Vmilksp is the volume of the milk space (l); Vresidue is the residual volume of the milk after milking (l); F is the ratio of the residual volume to the maximum volume of the milk space (unitless). Compared with the 1-compartment model, the 2-compartment model is more physiologically based, and can be used to simulate the drug distribution and movement in different parts of udder compartment. However, the 2-compartment model is relatively complex, includes more parameters and requires robust datasets to validate. Briefly, the volume changes in the total milk space and the alveolar compartment were both simulated with the Hill-Langmuir equation (equation 4). The volume in cistern compartment was determined by subtracting the volume of alveoli from the total milk space (equation 5). The residual volume of milk in the milk space of udder was mostly retained in the alveoli (Vetharaniam et al., 2003), and the residual volume in cistern was assumed as zero immediately after milking (Knight et al., 1994). Valv= Valvresidue+tm×VMalvtm+salv, (4) Vcis=Vmilksp-Valv, (5) where Valv is the volume of milk in the alveoli (l); tm is the time after the preceding milking (h); VMalv is the maximum volume of the alveoli (l); salv is the time for the volume of milk in alveoli reaching half of VMalv (h); Valvresidue is the residual volume of the alveoli after milking (l); Vcis is the volume of milk in the cistern (l). Both the 1-compartment and 2-compartment udder models were validated (Figs. 2A and 2B) based on published experimental data (Davis et al., 1998). As available pharmacokinetic data in milk are mostly sparse data and to avoid over parameterization, the 1-compartment udder model was used in the final PBPK model. The code describing the 2-compartment udder model is provided in the code in the Supplementary Data and can be used in the future when needed and when more robust datasets are available. Figure 2. View largeDownload slide Validation and application of the milk models. Model predictions for milk volumes (solid line) and observed data (black squares, circles, and diamonds) are compared in panel A for the 1-compartment model and in panel B for the 2-compartment model. The observed data are from previous study by Davis et al. (1998). In panel C, the milk production volumes by different milking intervals were simulated using 1-compartment milking model based on the parameter values reported in Whittem et al. (2012) (VMmilksp = 26.8 l, s = 23.2 h, F = 0.04). Figure 2. View largeDownload slide Validation and application of the milk models. Model predictions for milk volumes (solid line) and observed data (black squares, circles, and diamonds) are compared in panel A for the 1-compartment model and in panel B for the 2-compartment model. The observed data are from previous study by Davis et al. (1998). In panel C, the milk production volumes by different milking intervals were simulated using 1-compartment milking model based on the parameter values reported in Whittem et al. (2012) (VMmilksp = 26.8 l, s = 23.2 h, F = 0.04). The previous milk model assumed no absorption and mass transfer of drugs between the udder compartment and systemic circulation (Whittem et al., 2012). In order to incorporate the udder compartment into a whole-body PBPK model, the absorption of penicillin G from the udder to the systemic circulation and the excretion of penicillin G from plasma to milk were considered so that the current model could be applied to both IMM and IM administrations. The absorption and excretion processes of penicillin G involve both the passive diffusion and active carrier-mediated transport (Rasmussen, 1959; Schadewinkel-Scherkl et al., 1993; Ziv and Sulman, 1974). On the basis of previously reported experimental data (Schadewinkel-Scherkl et al., 1993), the carrier-mediated transport accounts for around 25% of the absorption of penicillin G at mammary gland (Supplementary Figure 1A) and for over 75% of the penicillin G excretion (Supplementary Figure 1B). Due to very limited experimental data for both processes, the parameters related to the kinetics of the absorption and excretion need to be determined by fitting with pharmacokinetic data. In order to reduce uncertainties, the excretion process was considered as the active carrier-mediated transport only, and the absorption process was described as only through the passive diffusion (equation 6). The Michaelis-Menten equation was used for the carrier-mediated transport of penicillin excreted into milk (equation 7), and the first-order linear equation was applied for the passive diffusion. The rate of change of penicillin G in the udder compartment was described using mass balance differential equations as previously described (DeWoskin et al., 2013; Leavens et al., 2012; Lin et al., 2016b). Only the penicillin G not bound to plasma proteins was considered as freely available for excretion. The mass balance equations of penicillin G in the udder compartment are described below. RU=QU (CAfree−CVU) – Rmilkex+Kab×Amilk (6) Rmilkex =Vmax×CVUCVU+Km, (7) CVU=CU/PU, (8) where RU is the rate of change in the concentration of penicillin G in the udder (mg/h); QU represents the volume of blood flow to the udder per hour (l/h); CAfree is the arterial blood concentration of penicillin G not bound with plasma proteins (mg/l); CVU is the concentration of penicillin G in the udder venous blood (mg/l); Rmilkex is the rate for penicillin excreted into milk (mg/h); Kab is the rate constant of penicillin G absorbed from milk (/h); Amilk is the amount of penicillin G in milk (mg); Vmax represents the maximum rate for penicillin G excretion from udder tissue to milk (mg/h); Km is the concentration of penicillin G to achieve the half Vmax (mg/l); CU is the penicillin G concentration in the udder tissue (mg/l); PU is the udder tissue: plasma partition coefficient (unitless). The penicillin input through IMM administration, excretion from plasma to milk, penicillin absorption from milk, and elimination through the episodic milking process were considered for the change of penicillin concentrations in milk (equation 9). The concentration of penicillin in milk was simulated by dividing amount of penicillin in milk by the milk volume (equations 10 and 11). Rmilk=Rinputimm+Rmilkex−Kab*Amilk−Rmilking, (9) d/dt(Amilk)=Rmilk, (10) Cmilk=Amilk/Vmilk, (11) where Rmilk is the rate of change in the concentration of penicillin G in milk (mg/h); Rinputimm stands for the input rate of penicillin through IMM administration (mg/h); Rmilkex is the rate for penicillin excreted into milk (mg/h); Kab is the rate constant of penicillin G absorbed from milk (/h); Amilk is the amount of penicillin G in milk (mg); Rmilking is the rate of change in the concentration of penicillin G by milking process (mg/h); Vmilk is the volume of milk produced since the preceding milking (l); Cmilk is the penicillin G concentration in milk (mg/l). The milking intervals and frequencies have impacts on the milk production (Gehring and Smith, 2006) and lead to the change of drug concentrations in milk. Briefly, higher milk production could lead to further dilution of drug concentrations in milk. In pharmacokinetic studies for penicillin G in milk, many employed unequal intervals such as 10 + 14 h or 15 + 9 h intervals for milking instead of equal 12-h milking intervals. In order to better fit the experimental data, the array equations were used (Supplemenatry equations 1 and 2). The IM and IMM administrations were incorporated into this PBPK model because they are the label routes for procaine penicillin G in dairy cows. The IM administration was simulated using a 2-compartment dissolution model based on the approach used to simulate the absorption of long-acting oxytetracycline (Lin et al., 2015). Detailed description of the dissolution model and the IM injection was included in our previous papers (Li et al., 2017; Lin et al., 2015). The multiple-dose therapeutic scenario was achieved by using the conditional operator of “IF…THEN…ELSE…” to create a control factor (Supplementary equations 3 and 7). The volume of milk was simulated using array functions, which helped to simulate the milking processes as the discrete steps. However, the amount of penicillin G in the milk and the udder was simulated using continuous functions. The DELAY functions were also applied here to help overcome the false peaks due to the function differences. For more details, please refer to the Supplementary Data. Evaluation and sensitivity analysis The performance of the PBPK model was evaluated by comparing model simulations with experimental pharmacokinetic data of milk and plasma not used in the model calibration (Table 1). On the basis of World Health Organization (WHO) guidelines (WHO, 2010), if the results from simulations matched the measured kinetic profiles and were generally within a factor of 2 of the measured values, the model was considered reasonable and validated. These criteria are based on the considerations that calibration datasets and evaluation datasets are obtained under different conditions (eg, different experimental animals/human subjects, laboratories, and detection methods), so some level of discordance is to be expected (WHO, 2010). The goodness of fit between log-transformed values of observed and predicted plasma and milk concentrations were further analyzed with linear regression and the determination coefficient (R2) was calculated for both calibration and evaluation results. The goodness of fit was evaluated using the determination coefficient (R2) for linear regression analysis and the mean absolute percentage error (MAPE) value. The R2 is an indicator for the overall model simulation performance, and a model simulation with a R2 value equal to or higher than 0.75 is considered acceptable. The analysis for MAPE was also carried out based on the previously reported methods (Cheng et al., 2016; Lin et al., 2017). The MAPE value lower than 50% was considered as the acceptable prediction. A local sensitivity analysis was performed for a discrete time point of 24 h to determine which parameters were most influential on the 24-h area under the curve (AUC) of the plasma, liver, kidney, and milk concentrations of penicillin G. Briefly, each parameter was increased by 1% and the corresponding 24-h AUC of penicillin concentrations was computed. Normalized sensitivity coefficient (NSC) was calculated using equations reported previously (Lin et al., 2011; Yoon et al., 2009). We conduced local sensitivity analysis rather than global sensitivity analysis because local sensitivity analysis is sufficient to identify highly sensitive parameters on our selected dose metrics of interest based on the multiple relevant studies (Craigmill, 2003; Leavens et al., 2014; Zeng et al., 2017; Zhu et al., 2017). Global sensitivity analysis could be done in the future if the goal is to study interaction effects between parameters (Lumen et al., 2015; McNally et al., 2011). Monte Carlo analysis Monte Carlo simulation obtains numerical results based on the repeated random sampling of parameter values from specified distribution of each parameter. This method has been used in the applications of PBPK modeling to estimate drug tissue residues and WDIs in the area of food safety (Buur et al., 2006; Yang et al., 2012; Zeng et al., 2017). In the current PBPK model for penicillin G in milk, Monte Carlo simulation was also applied to estimate the effects of parameter uncertainties and intraspecies variability of dairy cows on milk penicillin G concentrations. For these simulations, the relatively sensitive parameters, which were identified from the local sensitivity analysis (Supplementary Table 1), distributed randomly around the mean values were specified in Table 3. Log-normal distributions of model parameters were assumed for all chemical-specific parameters, such as partition coefficients, absorption rate constants, and elimination rate constants. The lognormal transformation for parameters in Monte Carlo simulations was achieved by using Supplementary equations 8–10. Physiological parameters, including body weights, cardiac outputs, and fractions of blood flows and tissue volumes were assumed to be normally distributed (Clewell III et al., 2000; Li et al., 2017; Shankaran et al., 2013; Sterner et al., 2013; Tan et al., 2006; Yang et al., 2015). Table 3. Values and Distributions of Parameters Used in the Monte Carlo Analysis for the PBPK Model of Dairy Cows Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 A pound sign (#) indicates the CV was calculated based on previous experimental data. Default values reported by previous models were used for CVs of the other parameter values. The default values of CVs for physiological parameters were 0.3, and for chemical-specific parameters were 0.2. The 5th and 95th percentiles of each parameter were calculated as the lower and upper bounds to represent 95% confidence interval. Please refer to Table 2 for the description of each parameter. Table 3. Values and Distributions of Parameters Used in the Monte Carlo Analysis for the PBPK Model of Dairy Cows Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 Parameter Distribution Mean SD CV Lower Bound Upper Bound BW Normal 299.96 46.180 0.15# 209.45 390.464 QCC Normal 5.97 1.990 0.33# 2.07 9.87 QKC Normal 0.09 0.027 0.30 0.037 0.143 VMmilksp Lognormal 26.8 8.040 0.30 11.042 42.558 S Lognormal 23.2 6.960 0.30 9.559 36.841 F Lognormal 0.04 0.012 0.30 0.016 0.064 Kim Lognormal 0.07 0.021 0.30 0.029 0.111 Frac Lognormal 0.6 0.012 0.30 0.576 0.624 VmaxC Lognormal 0.0022 0.0003 0.30 0.0014 0.0026 Km Lognormal 0.7 0.006 0.30 0.688 0.712 KabC Lognormal 1.00E-04 2.40E-05 0.30 5.30E-05 1.47E-04 PL Lognormal 3 0.600 0.20 1.824 4.176 PK Lognormal 2.5 0.500 0.20 1.52 3.48 KmetC Lognormal 0.0025 7.50E-4 0.30 0.001 0.004 KurineC Lognormal 0.45 0.135 0.30 0.185 0.715 A pound sign (#) indicates the CV was calculated based on previous experimental data. Default values reported by previous models were used for CVs of the other parameter values. The default values of CVs for physiological parameters were 0.3, and for chemical-specific parameters were 0.2. The 5th and 95th percentiles of each parameter were calculated as the lower and upper bounds to represent 95% confidence interval. Please refer to Table 2 for the description of each parameter. The detail of Monte Carlo simulation in Berkeley Madonna was introduced in our previous generic PBPK model for penicillin G in beef cattle (Li et al., 2017), and the code of the current model is available in the Supplementary Data. Briefly, model parameters were varied randomly around the values (central tendencies) used or estimated in model calibration by their probabilistic distributions (variability). The values of sensitive physiological and chemical-specific parameters were randomly selected based on their distributions, it is necessary to use adjustment factors to ensure the plausibility and mass balance for the PBPK model (Covington et al., 2007; Gelman et al., 1996; Sterner et al., 2013). The Monte Carlo simulation was set up to batch run for 1000 times with model parameters randomly selected from the defined distributions. Application of the model for milk production and tissue concentration prediction The model of milk production can be used as a stand-alone program to predict the milk production based on the different milking intervals. The 12, 12-h interval, 8, 16-h interval, and 8, 8, 8-h interval were used to simulate and compare the milk productions with different milking intervals. The PBPK model for penicillin G with the mammary gland compartment can be used to predict the systemic absorption and distribution of penicillin G after the IMM injection and to predict the milk concentration of penicillin G following IM injections. For IMM infusions, the label dose 100 000 IU (99 mg per quarter of udders) and 2 commonly used extralabel doses (2× label dose, 198 mg per quarter; 4× label dose, 396 mg per quarter) to 2 of the 4 quarters were simulated for 3 repeated doses with 12-h intervals. For IM administrations, the label dose 6000 IU/kg (6.5 mg/kg) and 2 commonly used extralabel doses (5× label dose, 32.5 mg/kg; 10× label dose, 65 mg/kg) were simulated for 3 repeated doses with 24-h intervals. A schematic illustrating the dosing regimens is provided in the Supplementary Figure 2. The commonly used extralabel doses for IMM infusions were based on the phone call records from the FARAD call center, and for IM administrations were obtained from the previous report (Payne et al., 2006). Infusions of penicillin G for all 4 quarters are only for the worst scenario of mastitis, and infusions for 2 quarters are common (Vilim et al., 1979). Here, IMM infusions for 2 quarters were simulated. These doses were applied as input for the current PBPK model to predict the concentrations of penicillin G in the milk, plasma, liver, and kidney. The Monte Carlo analysis was also used to simulate the distribution of penicillin concentrations based on the distribution of parameter values. Different therapeutic scenarios were analyzed using Monte Carlo simulations. Each simulation was run 1000 times. The mean value, 1st and 99th percentile values of simulated results were calculated and plotted. The extended meat WDI and the milk discard interval were determined when the values of 99th percentiles of the target tissue and milk concentrations, respectively, fell below the tolerance or the safe level. The tolerance of penicillin G for edible tissues in cattle is 0.05 µg/g (FDA, 2013). The safe level of penicillin G for milk is 5 ng/ml (FDA, 2005, 2015b) and the tolerance is 0 by U.S. FDA (Brynes, 2005; FDA, 2017). The limits of detection (LODs) for penicillin G in the milk from current methods are lower than or equal to 1 ng/ml, such as 1 ng/ml (Holstege et al., 2002), 0.25 ng/ml (Liu et al., 2011), and 0.23 ng/ml (Huang et al., 2013). In addition, we compared the model-predicted milk discard intervals from Monte Carlo simulations between 1000 and 10 000 iterations based on a 3 daily IM administration paradigm at the labeled dose level. RESULTS Validation and Application of the Milk Model The 1- and 2-compartment milk models were used to simulate the volume of milk produced. Due to the lack of published data for the volumes of milk produced by different milking intervals, the model was only calibrated with experimental data from a study that determined the pattern of milk accumulation in the milk space of udder over time (Davis et al., 1998), and the parameter values were based on the parameter values from the previous milk model (VMmilksp = 26.8 l, s = 23.2 h, F = 0.04) (Whittem et al., 2012). There were no other data available for the model evaluation. The 1-compartment model predicted the total milk volume in the milk space well (Figure 2A). The simulated results from 2-compartment model were also in good agreement with milk volumes measured in alveoli, cistern, and total milk space in udder (Figure 2B). Both the 1- and 2-compartment milk models could be used as the udder compartment for the PBPK model. The milk model could also be applied to simulate the impact of different milking intervals on the daily milk production (Figure 2C represents the results from the 1-compartment model). The results showed that the 12-h milking intervals produced slightly higher milk volumes compared with 6–18 h alternative intervals. Milking 4 times per day achieved the highest daily milk yields compared with milking 2 or 3 times a day, which is in accordance with previous reported research (Allen et al., 1986; Erdman and Varner, 1995). However, milking more than twice a day does not fit the circadian rhythm and the real-life working shifts well and requires more labor (Lyons et al., 2014). Model Calibration The PBPK model with the 1-compartment milk model was used to simulate penicillin G concentrations in plasma and milk after different therapeutic regimens. Model predictions were compared with observed concentrations in dairy cows exposed to penicillin G through IMM infusion or IM injection (results are shown in Figs. 3 and 4, respectively). Overall, the model well simulated the kinetic profiles of penicillin G in the milk and plasma in dairy cows through both IM and IMM administrations. In particular, the model can accurately predict the penicillin G residues in the milk and plasma even at the later time points, which are important for residue monitoring and the determination of times when the concentrations of residues in the milk fall below or near the tolerance of penicillin G. Figure 3. View largeDownload slide Calibration of the physiologically based pharmacokinetic model using pharmacokinetic data through the intramammary (IMM) infusion. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from dairy cows exposed to procaine penicillin G via single IMM infusion (396 mg per quarter for all quarters). Experimental data are from the previous study (Mercer et al. 1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 3. View largeDownload slide Calibration of the physiologically based pharmacokinetic model using pharmacokinetic data through the intramammary (IMM) infusion. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from dairy cows exposed to procaine penicillin G via single IMM infusion (396 mg per quarter for all quarters). Experimental data are from the previous study (Mercer et al. 1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 4. View largeDownload slide Calibration of the model using the pharmacokinetic data through the intramuscular (IM) injection. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from daily cows exposed to procaine penicillin G via single IM injection (10 mg/kg, A, B; 40 mg/kg, C, D). Experimental data of panels A and B are from the study carried out by Hogh and Rasmussen (1966), and panels C and D are from the pharmacokinetic data reported by Van OS et al. (1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 4. View largeDownload slide Calibration of the model using the pharmacokinetic data through the intramuscular (IM) injection. Comparison of model predictions (solid line) and observed data (red circles) for concentrations of penicillin G in the plasma and milk from daily cows exposed to procaine penicillin G via single IM injection (10 mg/kg, A, B; 40 mg/kg, C, D). Experimental data of panels A and B are from the study carried out by Hogh and Rasmussen (1966), and panels C and D are from the pharmacokinetic data reported by Van OS et al. (1974). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Model Evaluation The pharmacokinetic datasets not used for the model calibration were applied to evaluate the model performance. Measured concentrations of penicillin G in the milk and plasma of dairy cows after IMM and IM administrations for single or multiple doses were compared with the model predictions (Figs. 5 and 6). Overall, the model predicted the concentrations of penicillin G adequately at different time points. For the IMM administration, the evaluation results shown in Figure 5 indicate the model well simulated the penicillin G concentrations in milk following different milking intervals. The R2 values of the regression analyses between log-transformed values of measured and simulated concentrations of penicillin G in the milk and plasma from both calibration and evaluation datasets were all higher than 0.75 (Supplementary Figure 3). For the evaluation datasets, the value of overall R2 for the IMM administration was 0.99 (Supplementary Figure 3C), and for the IM administration was 0.92 (Supplementary Figure 3D). The R2 for data in milk was 0.81, and for plasma was 0.97 for the IM administration. All the MAPE values for calibration and evaluation results were lower than 40% (Supplementary Figure 4), which indicates that the model predictions are acceptable. In particular, the MAPE values for plasma for different simulations were generally lower than or near 20% indicating good predictions of plasma concentrations following either IM or IMM administrations. Overall, the PBPK model adequately simulates the measured results in the milk and plasma from the independent studies of dairy cows via both IMM and IM administrations. Figure 5. View largeDownload slide Evaluation of the model using intramammary (IMM) data with different milking intervals. Comparison of model predictions (solid line) and observed data (blue circles) for concentrations of penicillin G in the milk of dairy cows exposed to procaine penicillin G via IMM infusions for repeated 3 doses (1898 mg per quarter for all quarters, A, B; 100 mg per quarter for 2 quarters, C; 100 mg per quarter for all quarters, D) is shown. Experimental data for A, B (mean) are from the study reported by Knappstein et al. (2003), for panel C are from the previous pharmacokinetic study (Vilim et al. 1979), and for panel D are from the study performed by Vilim et al. (1980). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 5. View largeDownload slide Evaluation of the model using intramammary (IMM) data with different milking intervals. Comparison of model predictions (solid line) and observed data (blue circles) for concentrations of penicillin G in the milk of dairy cows exposed to procaine penicillin G via IMM infusions for repeated 3 doses (1898 mg per quarter for all quarters, A, B; 100 mg per quarter for 2 quarters, C; 100 mg per quarter for all quarters, D) is shown. Experimental data for A, B (mean) are from the study reported by Knappstein et al. (2003), for panel C are from the previous pharmacokinetic study (Vilim et al. 1979), and for panel D are from the study performed by Vilim et al. (1980). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 6. View largeDownload slide Evaluation of the model using intramuscular (IM) data with different milking intervals. Comparison of model predictions (solid line and/or dotted line) and observed data (blue squares or circles) for concentrations of penicillin G in the plasma and milk of dairy cows exposed to procaine penicillin G via 1 single IM injection (11 mg/kg, A, B; 2.2 mg/kg, D; 4.4 mg/kg, E; 11 mg/kg, F), and IM injection for repeated 4 doses (6.5 mg/kg, C) is shown. Experimental data for A, B (mean) are from the study reported by Randall et al. (1954), for panel C are from reported pharmacokinetic data (FDA, 2010), and for panels D–F are from the previous study by Edwards et al. (1966). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Figure 6. View largeDownload slide Evaluation of the model using intramuscular (IM) data with different milking intervals. Comparison of model predictions (solid line and/or dotted line) and observed data (blue squares or circles) for concentrations of penicillin G in the plasma and milk of dairy cows exposed to procaine penicillin G via 1 single IM injection (11 mg/kg, A, B; 2.2 mg/kg, D; 4.4 mg/kg, E; 11 mg/kg, F), and IM injection for repeated 4 doses (6.5 mg/kg, C) is shown. Experimental data for A, B (mean) are from the study reported by Randall et al. (1954), for panel C are from reported pharmacokinetic data (FDA, 2010), and for panels D–F are from the previous study by Edwards et al. (1966). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article). Sensitivity Analysis The sensitivity analysis was applied to screen the sensitive parameters in the PBPK model. Thirty-four model parameters in the PBPK model for dairy cows were analyzed using the method of local sensitivity analysis. Results of the local sensitivity analysis based on 1% variation of the parameter values are shown in Supplementary Table 1. Only parameters with at least 1 absolute value of NSC > 0.15 are shown in the table. For physiological parameters, only the body weight (BW), the cardiac output (QCC), and the blood flow to kidney (QKC) were the relatively sensitive parameters for the model. The results indicate that the PBPK model was not that sensitive to physiological parameters. The rate constant of penicillin G absorbed from milk back into the udder tissue (KabC) was the sensitive parameter for AUCs of the plasma, liver and kidney, and the fraction of undissolved procaine penicillin G (Frac) was sensitive for all AUCs. The AUC of liver was also highly sensitive to the liver partition coefficient (PL) with NSC value of 1.00. The AUC of kidney was highly sensitive to the urine elimination rate constant and the kidney partition coefficient with NSC values of −0.90 and 1.00, respectively. For AUC of milk, it was sensitive to the body weight (BW), the maximum volume of the milk space in udder (VMmilksp), the time for Vmilk reaching half of VMmilksp (s), the ratio of the residual volume to the maximum volume of the milk space (F), and the maximum rate of penicillin G excretion in the mammary gland (VmaxC; body weight-scaled) with NSC values of 0.99, −0.99, 0.61, −0.24, and 1.00, respectively. Probabilistic Analysis The relatively sensitive parameters listed in Supplementary Table 1 were subsequently used for probabilistic analysis based on the current PBPK model for penicillin G in dairy cows. The Monte Carlo sampling technique was applied for the probabilistic analysis. From the simulation results of the PBPK model after label and extralabel use of procaine penicillin G through IM and IMM administrations, the concentrations of penicillin G in the liver are higher than other edible tissues and plasma (Supplementary Figs. 5B and 5D). Therefore, the liver was chosen as the tissue to determine the WDIs for extralabel use of penicillin G in dairy cows. The Monte Carlo simulations showed that the WDIs after 3 repeated IMM infusions with label dose and 4× label dose in dairy cows were 98 and 122 h for milk, and were both 0 h for edible tissues, respectively (Figs. 7A–D). The predicted WDIs after 3 doses via IM injections with label dose and 10× label dose were 125 and 182 h for milk, and were 4 and 5 days for edible tissues (Figs. 7E–H). The exact model-predicted milk discard intervals are reported in this article, but in practice, the recommended milk discard intervals should be adjusted accordingly to fit within a milking schedule for a particular dairy farm. If the estimated WDI of edible tissues was a fraction of a day, it was rounded up to the next whole day. The results of the probabilistic analysis of the PBPK model for 2× label dose through IMM infusions and 5× label dose following IM administrations are available in the Supplementary Figure 6, and the results for probabilistic analysis in the plasma and kidney are provided in the Supplementary Figure 7. No considerable differences were seen between the Monte Carlo simulation results between 1000 and 10 000 iterations (Supplementary Figure 8). The label withdrawal periods were obtained from the Veterinarian’s Guide to Residue Avoidance Management (VetGRAM) of FARAD (Riviere et al., 2017). The labeled milk discard time via IMM infusions is 84 h, and for edible tissues the label withdrawal period is 4 days. The label withdrawal period of milk via IM injections is 48 h, and of edible tissues is 14 days in dairy cows. Figure 7. View largeDownload slide Probabilistic analysis for penicillin G concentrations in milk and liver via intramammary (IMM) and intramuscular (IM) administrations. IMM infusions with label dose (99 mg per quarter for 2 quarters) and 4× label dose (396 mg per quarter for 2 quarters) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (A–D). IM injections with label dose (6.5 mg/kg) and 10× label dose (65 mg/kg) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (E–H). The milking intervals were 12 h. Though zero tolerance was established for penicillin G in milk, the safe level (SL) of 5 µg/kg was used by FDA in actual practice. Tolerances (TOL) for penicillin G in edible tissues are 50 µg/kg. Figure 7. View largeDownload slide Probabilistic analysis for penicillin G concentrations in milk and liver via intramammary (IMM) and intramuscular (IM) administrations. IMM infusions with label dose (99 mg per quarter for 2 quarters) and 4× label dose (396 mg per quarter for 2 quarters) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (A–D). IM injections with label dose (6.5 mg/kg) and 10× label dose (65 mg/kg) for 3 repeated doses were simulated as the therapeutic scenario for dairy cows (E–H). The milking intervals were 12 h. Though zero tolerance was established for penicillin G in milk, the safe level (SL) of 5 µg/kg was used by FDA in actual practice. Tolerances (TOL) for penicillin G in edible tissues are 50 µg/kg. DISCUSSION In this study, a PBPK model for penicillin G in dairy cows was established for labeled routes of IM and IMM administrations. On the basis of model, probabilistic analysis was performed to predict the milk discard intervals and the extended WDIs for edible tissues after IM or IMM administration of procaine penicillin G at extralabel doses. The PBPK model for dairy cows was established and calibrated with pharmacokinetic data of penicillin G through IM and IMM administrations. As the systemic absorption and the drug excretion were both taken into account, the model can be used to predict the concentrations of penicillin G in milk and edible tissues through either the systemic or the local administration. Application of this PBPK model can help predict the extended milk discard intervals after extralabel use of penicillin G. In the field of toxicology, several PBPK models with a milk compartment are available, such as the lactational PBPK models for atrazine in rodents (Lin et al., 2013), perfluorooctanoic acid in humans (Loccisano et al., 2013), manganese in humans (Yoon et al., 2011), and flunixin in dairy cows (Leavens et al., 2014). The milk compartments in all these models were described as virtual compartments and were not physiologically based, which limits the models’ capability to extrapolate to other species. In the present model, a physiologically based approach simulating the mammary gland compartment was established based on the milk production and episodic milking processes to facilitate realistic prediction of drug concentrations in the milk and systemic distribution through the IMM infusion. The Hill-Langmuir equation was adapted from previous milk model (Whittem et al., 2012) to simulate the decrease of milk production due to the increased hydrostatic pressure (Mercer et al., 1970; Smith et al., 2004; Whittem et al., 2012). The 1- and 2-compartment milk models were both adapted to the current PBPK model for penicillin G. The 1-compartment is good enough for the sparse data available from existing pharmacokinetics studies. The 2-compartment model could help predict the drug movements within the 2 compartments of mammary glands, which is useful for the estimation of the drug concentration at the diseased state. The 2-compartment model would help design therapeutic regimens of penicillin G and other veterinary drugs for cows. Both 1- and 2-compartment models provide a basis to simulate pharmacokinetics of penicillin G in mastitis cows. The systemic and local routes of drug administrations, the milk production process including the milk producing rate and milking intervals, and the movement of penicillin G such as the episodic milking process and transfer of penicillin G between the plasma, udder tissue and milk were all considered in the current PBPK model. This model greatly improves our understanding of the pharmacokinetic property of penicillin G in dairy cows, especially in the mammary glands after both IM and IMM administrations. This physiologically based approach simulating the milk compartment can be extrapolated to other species, including rats, humans, and goats. Currently, different statistical methods, such as the “Time To Safe Concentration,” the “Safe Concentration based on Linear Regression” and the “Safe Concentration Per Milking” methods, were used to calculate the withdrawal period for milk (Henri et al., 2017b). However, these methods have various limitations, including the assumption of linearity of the log-transformed tissue/milk depletion data, which may not be true for all drugs (Chevance et al., 2017; Henri et al., 2017b). The nonlinear mixed-effects pharmacokinetic modeling has many advantages (Martin-Jimenez and Riviere, 1998) as a pharmaco-statistical approach to estimate the withdrawal period for milk (Chevance et al., 2017). If the pharmacokinetic data are not available for extralabel use of veterinary drugs, the method of half-life multipliers could be applied to calculate the extended WDIs (Riviere et al., 1998; Smith et al., 2004). However, the nonlinear mixed-effects pharmacokinetic model and the half-life multipliers method are not based on the mechanisms of drug action and the physiology, and cannot be extrapolated to different therapeutic scenarios and animal species. The current PBPK model for milk is a physiologically- and mechanistic-based model and could be extrapolated to the different species, drugs and therapeutic scenarios to predict tissue residues and milk concentrations. The current model is based on the assumption that the 4 quarters of the udder tissues are identical in dairy cows, and no milk production changes were considered for different lactating stages. As the model simulations are about 2 weeks, the potential impact of lactational stages on the milk production was not considered in this model. On the basis of available research (Novak et al., 2009; Yoon et al., 2004), the lactational stages could be incorporated into the model if needed. The parameters of the absorption and excretion processes were determined through model fitting with pharmacokinetics data. These parameters are related to the organic cation transporters and organic anion transporters expressed in the mammary glands of dairy cows (Al-Bataineh et al., 2009). However, the information on the expression of transporters related to the drug transport in mammary glands of dairy cows is limited. In the future, with more research on transporters in food-producing animals, the strategy of in vitro to in vivo extrapolation can be used to calculate the relevant parameter values for drug distribution through mammary glands. For some sensitive parameters related to the milk secretion, such as VMmilksp, s, and F, although their values were validated in an earlier study (Whittem et al., 2012), the experimentally measured values are not available and not commonly reported by pharmacokinetic studies for drugs in milk. These values, whenever available experimentally, would help reduce uncertainties for the current milk model. The binding of penicillin G to milk proteins is not considered in the current model, due to the limited information on protein binding of penicillin G in milk. The current model could be improved when more experimental data on protein binding of penicillin G in milk become available in the future. Also note that the milk is a suspension of fat droplets in an aqueous phase containing proteins, carbohydrates, and electrolytes, in which proteins can extensively bind some antibiotics, resulting in uneven distribution of the drug in different components of the milk (Ziv and Rasmussen, 1975). It has been shown that drug distribution in different fractions of milk can be very different for a systemic route of administration versus a local infusion depending on the drug lipophilicity (Ziv and Rasmussen, 1975). Thus, drug concentrations can be different between different milk samples (whole milk vs skim milk; foremilk vs hindmilk) depending on the drug lipophilicity and sampling methods. At this stage, there are no sufficient data that are granular enough to predict drug concentrations in different components of milk. Therefore, the present model may need some refinements depending on the investigated drugs (lipophilic vs hydrophilic), the dairy products of interest (whole milk vs skim milk), and the sampling methods (cisternal milk vs total milk). Additional studies are needed in order to establish a milk model considering different milk components, which will help model extrapolation to different animal species and different drugs. This is important because the milk compositions are different among different species (eg, 3.3%–5.4% fat in cow milk vs 5.3%–9.0% in buffalo vs 10.2%–21.5% in deer milk) (Claeys et al., 2014) and different drugs have different interactions with milk components. The milk discard interval for penicillin G following labeled IMM infusion predicted by the current model was in accordance with the labeled milk discard time from FDA (98 vs 84 h). The model-predicted milk discard interval following labeled IM injection was longer than the labeled milk discard time by FDA (125 vs 48 h), suggesting that the model-predicted result after labeled IM administration is relatively conservative. The inconsistence between the model-predicted results and the labeled milk discard time after IM administration may be due to different pharmacokinetic studies used in our model compared with the FDA method. The current model was based on a number of datasets from different studies using different drug formulations in different animal populations, thus the results consider the variability of a large diverse population of animals. On the other hand, the labeled milk discard time is typically obtained based on the data from a small set of animals. Overall, the present model results support the recommendations for the appropriate time to harvest milk following the use of penicillin G. On the basis of simulation results, the model-predicted WDIs in the target tissue liver after labeled IMM or IM exposure were consistently much lower than the labeled withdrawal periods. Also, the IMM infusion did not lead to the violations in edible tissues even at 4× label dose. However, major violations of penicillin G were identified in cull dairy cows by USDA. Also, a prolonged elimination of penicillin G in emergency-slaughter dairy cows after IM administration of penicillin G has been reported (Nouws and Ziv, 1978). These inconsistencies may be due to the pharmacokinetic data used for the model calibration of the current model were all from healthy dairy cows. Mastitis and other diseased conditions, as well as the interactions with other drugs in diseased animals may lead to the different distribution and elimination patterns of penicillin G. Due to the fact that pathophysiological data of mastitis in dairy cows are limited, the diseased model for dairy cows with mastitis is not available now. The diseased model is needed to more accurately simulate tissue distribution of penicillin G and to avoid tissue residue violations in cull dairy cows. Thus, the present model cannot be used to predict extended WDIs in target tissue after extralabel use of penicillin G in cull dairy cows. The PBPK model in dairy cows with mastitis would be a future direction for our research when additional data become available. In conclusion, the PBPK model for penicillin G in dairy cows successfully predicts the drug concentrations in the milk and plasma following IM or IMM administration of procaine penicillin G. This is the first reported PBPK model with a physiologically based compartment of mammary glands. By utilizing probabilistic Monte Carlo analysis, the model could help estimate a conservative milk discard time after extralabel use of penicillin G, thereby facilitating food safety assessment of cow-derived milk products. This PBPK framework can be applied to different antibiotics used for dairy cows to help avoid residue violations in the milk and reduce the amount of waste milk. This modeling strategy could also be extrapolated to other species, such as rodents, humans, and goats. SUPPLEMENTARY DATA Supplementary data are available at Toxicological Sciences online. ACKNOWLEDGMENTS The authors would like to acknowledge Dr Pritam Sidhu, Dr Yi-Hsien Cheng, Dr Dongping Zeng, and Shiqiang Jin in the Institute of Computational Comparative Medicine at Kansas State University for helpful discussions. The authors declare no conflict of interest. FUNDING United States Department of Agriculture (USDA) National Institute of Food and Agriculture (NIFA) for the Food Animal Residue Avoidance Databank (FARAD) Program (Award No. 2016-41480-25729 and 2017-41480-27310); The Kansas Bioscience Authority funds to the Institute of Computational Comparative Medicine (ICCM) at Kansas State University. REFERENCES Al-Bataineh M. M. , Van Der Merwe D. , Schultz B. D. , Gehring R. ( 2009 ). Cultured mammary epithelial monolayers (BME-UV) express functional organic anion and cation transporters . J. Vet. Pharmacol. Ther. 32 , 422 – 428 . Google Scholar CrossRef Search ADS PubMed Allen D. B. , DePeters E. J. , Laben R. C. ( 1986 ). Three times a day milking: Effects on milk production, reproductive efficiency, and udder health . J. Dairy Sci. 69 , 1441 – 1446 . Google Scholar CrossRef Search ADS PubMed Andersen M. E. ( 2003 ). Toxicokinetic modeling and its applications in chemical risk assessment . Toxicol. Lett. 138 , 9 – 27 . Google Scholar CrossRef Search ADS PubMed Baynes R. , Riviere J. E. ( 2014 ). Importance of veterinary drug residues. In Strategies for Reducing Drug and Chemical Residues in Food Animals: International Approaches to Residue Avoidance, Management, and Testing ( Baynes R. , Riviere J.E. , Eds.), pp. 1 – 8 . Wiley , New York, NY . Baynes R. E. , Dedonder K. , Kissell L. , Mzyk D. , Marmulak T. , Smith G. , Tell L. , Gehring R. , Davis J. , Riviere J. E. ( 2016 ). Health concerns and management of select veterinary drug residues . Food Chem. Toxicol. 88 , 112 – 122 . Google Scholar CrossRef Search ADS PubMed Brynes S. D. ( 2005 ). Demystifying 21 CFR Part 556–tolerances for residues of new animal drugs in food . Regul. Toxicol. Pharmacol. 42 , 324 – 327 . Google Scholar CrossRef Search ADS PubMed Buur J. , Baynes R. , Smith G. , Riviere J. ( 2006 ). Use of probabilistic modeling within a physiologically based pharmacokinetic model to predict sulfamethazine residue withdrawal times in edible tissues in swine . Antimicrob. Agents Chemother. 50 , 2344 – 2351 . Google Scholar CrossRef Search ADS PubMed Campbell J. R. , Marshall R. T. ( 2016 ). Physiology of Lactation, Dairy Production and Processing: The Science of Milk and Milk Products . Waveland Press , Long Grove, IL . Carruthers V. R. , Davis S. R. , Bryant A. M. , Henderson H. V. , Morris C. A. , Copeman P. J. A. ( 1993 ). Response of Jersey and Friesian cows to once a day milking and prediction of response based on udder characteristics and milk composition . J. Dairy Res. 60 , 1 – 11 . Google Scholar CrossRef Search ADS Cheng Y.-H. , Lin Y.-J. , You S.-H. , Yang Y.-F. , How C. M. , Tseng Y.-T. , Chen W.-Y. , Liao C.-M. ( 2016 ). Assessing exposure risks for freshwater tilapia species posed by mercury and methylmercury . Ecotoxicology 25 , 1181 – 1193 . Google Scholar CrossRef Search ADS PubMed Chevance A. , Jacques A. M. , Laurentie M. , Sanders P. , Henri J. ( 2017 ). The present and future of withdrawal period calculations for milk in the European Union: Focus on heterogeneous, nonmonotonic data . J. Vet. Pharmacol. Ther. 40 , 218 – 230 . Google Scholar CrossRef Search ADS PubMed Chiesa O. A. , Von Bredow J. , Smith M. , Heller D. , Condon R. , Thomas M. H. ( 2006 ). Bovine kidney tissue/biological fluid correlation for penicillin . J. Vet. Pharmacol. Ther. 29 , 299 – 306 . Google Scholar CrossRef Search ADS PubMed Claeys W. L. , Verraes C. , Cardoen S. , De Block J. , Huyghebaert A. , Raes K. , Dewettinck K. , Herman L. ( 2014 ). Consumption of raw or heated milk from different species: An evaluation of the nutritional and potential health benefits . Food Control 42 , 188 – 201 . Google Scholar CrossRef Search ADS Clewell H. J. III , Gentry P. R. , Covington T. R. , Gearhart J. M. ( 2000 ). Development of a physiologically based pharmacokinetic model of trichloroethylene and its metabolites for use in risk assessment . Environ. Health Perspect . 108 , 283 – 305 . Google Scholar CrossRef Search ADS PubMed Covington T. R. , Robinan Gentry P. , Van Landingham C. B. , Andersen M. E. , Kester J. E. , Clewell H. J. III ( 2007 ). The use of Markov chain Monte Carlo uncertainty analysis to support a Public Health Goal for perchloroethylene . Regul. Toxicol. Pharmacol. 47 , 1 – 18 . Google Scholar CrossRef Search ADS PubMed Craigmill A. L. ( 2003 ). A physiologically based pharmacokinetic model for oxytetracycline residues in sheep . J. Vet. Pharmacol. Ther . 26 , 55 – 63 . Google Scholar CrossRef Search ADS PubMed Craigmill A. L. , Riviere J. E. , Webb A. I. ( 2006 ). Tabulation of FARAD Comparative and Veterinary Pharmacokinetic Data . Wiley-Blackwell Publishing , Ames, IA . Davis S. R. , Farr V. C. , Copeman P. J. , Carruthers V. R. , Knight C. H. , Stelwagen K. ( 1998 ). Partitioning of milk accumulation between cisternal and alveolar compartments of the bovine udder: Relationship to production loss during once daily milking . J. Dairy Res. 65 , 1 – 8 . Google Scholar CrossRef Search ADS PubMed Dayan A. D. ( 1993 ). Allergy to antimicrobial residues in food: Assessment of the risk to man . Vet. Microbiol. 35 , 213 – 226 . Google Scholar CrossRef Search ADS PubMed DeWoskin R. S. , Sweeney L. M. , Teeguarden J. G. , Sams R. , Vandenberg J. ( 2013 ). Comparison of PBTK model and biomarker based estimates of the internal dosimetry of acrylamide . Food Chem. Toxicol. 58 , 506 – 521 . Google Scholar CrossRef Search ADS PubMed DrugBank ( 2018 ). Benzylpenicillin. Available at: https://www.drugbank.ca/drugs/DB01053. Last Accessed March 29, 2018. Edwards S. ( 1966 ). Penicillin levels in the milk following intramuscular injection . Vet. Record 78 , 583 – 585 . Google Scholar CrossRef Search ADS Erdman R. A. , Varner M. ( 1995 ). Fixed yield responses to increased milking frequency . J. Dairy Sci. 78 , 1199 – 1203 . Google Scholar CrossRef Search ADS PubMed European Medicines Agency (EMA) ( 2000 ). Guideline on Withdrawal Periods for Milk. European Medicines Agency , London, UK . Available at: http://www.ema.europa.eu/docs/en_GB/document_library/Scientific_guideline/2009/10/WC500004496.pdf. Last Accessed March 29, 2018. European Medicines Agency (EMA) ( 2016 ). Guideline on Approach Towards Harmonisation of Withdrawal Periods. European Medicines Agency , London, UK . Available at: http://www.ema.europa.eu/docs/en_GB/document_library/Scientific_guideline/2016/07/WC500210929.pdf. Last Accessed March 29, 2018. Gehring R. , Smith G. W. ( 2006 ). An overview of factors affecting the disposition of intramammary preparations used to treat bovine mastitis . J. Vet. Pharmacol. Ther. 29 , 237 – 241 . Google Scholar CrossRef Search ADS PubMed Gelman A. , Bois F. , Jiang J. ( 1996 ). Physiological pharmacokinetic analysis using population modeling and informative prior distributions . J. Am. Stat. Assoc. 91 , 1400 – 1412 . Google Scholar CrossRef Search ADS Gionbelli M. P. , Duarte M. S. , Valadares Filho S. C. , Detmann E. , Chizzotti M. L. , Rodrigues F. C. , Zanetti D. , Gionbelli T. R. S. , Machado M. G. ( 2015 ). Achieving body weight adjustments for feeding status and pregnant or non-pregnant condition in beef cows . PLoS One 10 , e0112111. Google Scholar CrossRef Search ADS PubMed Gomes E. R. , Demoly P. ( 2005 ). Epidemiology of hypersensitivity drug reactions . Curr. Opin. Allergy Clin. Immunol. 5 , 309 – 316 . Google Scholar CrossRef Search ADS PubMed Henri J. , Carrez R. , Méda B. , Laurentie M. , Sanders P. ( 2017a ). A physiologically based pharmacokinetic model for chickens exposed to feed supplemented with monensin during their lifetime . J. Vet. Pharmacol. Ther. 40 , 370 – 382 . Google Scholar CrossRef Search ADS Henri J. , Jacques A. M. , Sanders P. , Chevance A. , Laurentie M. ( 2017b ). The present and future of withdrawal period calculations for milk in the European Union: Dealing with data below the limit of quantification . J. Vet. Pharmacol. Ther. 40 , 116 – 122 . Google Scholar CrossRef Search ADS Hogh P. , Rasmussen F. ( 1966 ). Mammary excretion of penicillin fellowing intramuscular injection of Leocillin_1tnR (penethamate hydriodide BAN) and penicillin procaine in cows . Nord. Vet. Med . 18 , 545 – 554 . Holstege D. M. , Puschner B. , Whitehead G. , Galey F. D. ( 2002 ). Screening and mass spectral confirmation of β-lactam antibiotic residues in milk using LC-MS/MS . J. Agric. Food Chem. 50 , 406 – 411 . Google Scholar CrossRef Search ADS PubMed Huang X. , Chen L. , Chen M. , Yuan D. , Nong S. ( 2013 ). Sensitive monitoring of penicillin antibiotics in milk and honey treated by stir bar sorptive extraction based on monolith and LC-electrospray MS detection . J. Sep. Sci. 36 , 907 – 915 . Google Scholar CrossRef Search ADS PubMed Isaksson A. , Arnarp L. ( 1988 ). Quantitative estimation of residual milk in bovine udders—A methodological study . Acta Vet. Scand. 29 , 259 – 262 . Google Scholar PubMed Knappstein K. , Suhren G. , Walte H. ( 2003 ). Influence of milking frequency on withdrawal period after application of beta-lactam antibiotics-based drugs . Anal. Chim. Acta 483 , 241 – 249 . Google Scholar CrossRef Search ADS Knight C. H. , Hirst D. , Dewhurst R. J. ( 1994 ). Milk accumulation and distribution in the bovine udder during the interval between milkings . J. Dairy Res. 61 , 167 – 177 . Google Scholar CrossRef Search ADS PubMed Leavens T. L. , Tell L. A. , Clothier K. A. , Griffith R. W. , Baynes R. E. , Riviere J. E. ( 2012 ). Development of a physiologically based pharmacokinetic model to predict tulathromycin distribution in goats . J. Vet. Pharmacol. Ther. 35 , 121 – 131 . Google Scholar CrossRef Search ADS PubMed Leavens T. L. , Tell L. A. , Kissell L. W. , Smith G. W. , Smith D. J. , Wagner S. A. , Shelver W. L. , Wu H. , Baynes R. E. , Riviere J. E. ( 2014 ). Development of a physiologically based pharmacokinetic model for flunixin in cattle (Bos taurus) . Food Addit. Contam. Part A Chem. Anal. Control Expo. Risk Assess. 31 , 1506 – 1521 . Google Scholar CrossRef Search ADS PubMed Li M. , Gehring R. , Riviere J. E. , Lin Z. ( 2017 ). Development and application of a population physiologically based pharmacokinetic model for penicillin G in swine and cattle for food safety assessment . Food Chem. Toxicol. 107 , 74 – 87 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Fisher J. W. , Ross M. K. , Filipov N. M. ( 2011 ). A physiologically based pharmacokinetic model for atrazine and its main metabolites in the adult male C57BL/6 mouse . Toxicol. Appl. Pharmacol. 251 , 16 – 31 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Fisher J. W. , Wang R. , Ross M. K. , Filipov N. M. ( 2013 ). Estimation of placental and lactational transfer and tissue distribution of atrazine and its main metabolites in rodent dams, fetuses, and neonates with physiologically based pharmacokinetic modeling . Toxicol. Appl. Pharmacol. 273 , 140 – 158 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Gehring R. , Mochel J. P. , Lave T. , Riviere J. E. ( 2016a ). Mathematical modeling and simulation in animal health—Part II: Principles, methods, applications, and value of physiologically based pharmacokinetic modeling in veterinary medicine and food safety assessment . J. Vet. Pharmacol. Ther. 39 , 421 – 438 . Google Scholar CrossRef Search ADS Lin Z. , Jaberi-Douraki M. , He C. , Jin S. , Yang R. S. , Fisher J. W. , Riviere J. E. ( 2017 ). Performance assessment and translation of physiologically based pharmacokinetic models from acslX to Berkeley Madonna, MATLAB, and R language: Oxytetracycline and gold nanoparticles as case examples . Toxicol. Sci. 158 , 23 – 35 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Li M. , Gehring R. , Riviere J. E. ( 2015 ). Development and application of a multiroute physiologically based pharmacokinetic model for oxytetracycline in dogs and humans . J. Pharm. Sci. 104 , 233 – 243 . Google Scholar CrossRef Search ADS PubMed Lin Z. , Vahl C. I. , Riviere J. E. ( 2016b ). Human food safety implications of variation in food animal drug metabolism . Sci. Rep. 6 , 27907 . Google Scholar CrossRef Search ADS Liu C. J. , Wang H. , Jiang Y. B. , Du Z. X. ( 2011 ). Rapid and simultaneous determination of amoxicillin, penicillin G, and their major metabolites in bovine milk by ultra-high-performance liquid chromatography-tandem mass spectrometry . J. Chromatogr. B 879 , 533 – 540 . Google Scholar CrossRef Search ADS Loccisano A. E. , Longnecker M. P. , Campbell J. L. , Andersen M. E. , Clewell H. J. III ( 2013 ). Development of PBPK models for PFOA and PFOS for human pregnancy and lactation life stages . J. Toxicol. Environ. Health Part A 76 , 25 – 57 . Google Scholar CrossRef Search ADS PubMed Lumen A. , McNally K. , George N. , Fisher J. W. , Loizou G. D. ( 2015 ). Quantitative global sensitivity analysis of a biologically based dose-response pregnancy model for the thyroid endocrine system . Front. Pharmacol . 6 , 107. Google Scholar CrossRef Search ADS PubMed Lyons N. A. , Kerrisk K. L. , Garcia S. C. ( 2014 ). Milking frequency management in pasture-based automatic milking systems: A review . Livest. Sci. 159 , 102 – 116 . Google Scholar CrossRef Search ADS Martin-Jimenez T. , Riviere J. E. ( 1998 ). Population pharmacokinetics in veterinary medicine: Potential use for therapeutic drug monitoring and prediction of tissue residues . J. Vet. Pharmacol. Ther. 21 , 167 – 189 . Google Scholar CrossRef Search ADS PubMed McNally K. , Cotton R. , Loizou G. D. ( 2011 ). A workflow for global sensitivity analysis of PBPK models . Front. Pharmacol . 2 , 31. Google Scholar CrossRef Search ADS PubMed Mercer H. , Geleta J. , Carter G. ( 1974 ). Absorption and excretion of penicillin G from the mastitic bovine udder . J. Am. Vet. Med. Assoc . 164 , 613 – 617 . Google Scholar PubMed Mercer H. , Geleta J. , Schultz E. , Wright W. ( 1970 ). Milk-out rates for antibiotics in intramammary infusion products used in the treatment of bovine mastitis: Relationship of somatic cell counts, milk production level, and drug vehicle . Am. J. Vet. Res. 31 , 1549 – 1560 . Google Scholar PubMed National Research Council (NRC) ( 1999 ). The use of drugs in food animals. In The Use of Drugs in Food Animals: Benefits and Risks ( National Research Council , Ed.), pp. 110 – 141 . National Academy Press , Washington, DC . Nouws J. F. , Ziv G. ( 1978 ). Tissue distribution and residues of benzylpenicillin and aminoglycoside antibiotics in emergency-slaughtered ruminants . Tijdschr. Diergeneeskd . 103 , 140 – 151 . Google Scholar PubMed Novak P. , Vokralova J. , Broucek J. ( 2009 ). Effects of the stage and number of lactation on milk yield of dairy cows kept in open barn during high temperatures in summer months . Archiv Tierzucht 52 , 574 – 586 . Paige J. C. , Tollefson L. , Miller M. ( 1997 ). Public health impact on drug residues in animal tissues . Vet. Hum. Toxicol . 39 , 162 – 169 . Google Scholar PubMed Papich M. G. ( 1987 ). The beta-lactam antibiotics: Clinical pharmacology and recent developments . Compend. Contin. Educ. Pract. Vet . 9 , 68 – 76 . Papich M. G. , Riviere J. E. ( 2009 ). Lactam antibiotics: penicillins, cephalosporins, and related drugs. In Veterinary Pharmacology and Therapeutics. J.E. Riviere, and M.G. Papich, Eds., 9th ed., pp. 866–888. Wiley, Ames, IA. Papich M. G. , Korsrud G. O. , Boison J. O. , Yates W. D. , MacNeil J. D. , Janzen E. D. , Cohen R. D. , Landry D. A. ( 1993 ). A study of the disposition of procaine penicillin G in feedlot steers following intramuscular and subcutaneous injection . J. Vet. Pharmacol. Ther. 16 , 317 – 327 . Google Scholar CrossRef Search ADS PubMed Payne M. A. , Craigmill A. , Riviere J. E. , Webb A. I. ( 2006 ). Extralabel use of penicillin in food animals . J. Am. Vet. Med. Assoc. 229 , 1401 – 1403 . Google Scholar CrossRef Search ADS PubMed Portis E. , Lindeman C. , Johansen L. , Stoltman G. ( 2012 ). A ten-year (2000-2009) study of antimicrobial susceptibility of bacteria that cause bovine respiratory disease complex–Mannheimia haemolytica, Pasteurella multocida, and Histophilus somni—In the United States and Canada . J. Vet. Diagn. Invest. 24 , 932 – 944 . Google Scholar CrossRef Search ADS PubMed Raison-Peyron N. , Messaad D. , Bousquet J. , Demoly P. ( 2001 ). Anaphylaxis to beef in penicillin-allergic patient . Allergy 56 , 796 – 797 . Google Scholar CrossRef Search ADS PubMed Randall W. , Durbin C. , Wilner J. , Collins J. ( 1954 ). Antibiotic Concentration and Duration in Animal Tissue and Body Fluids. Antibiotics Annual 1953-1954 . Medical Encyclopedia, Inc ., New York, NY . Rasmussen F. ( 1959 ). Mammary excretion of benzylpenicillin, erythromycin, and penethamate hydroiodide . Acta Pharmacol. Toxicol. 16 , 194 – 200 . Google Scholar CrossRef Search ADS Riviere J. E. , Craigmill A. L. , Sundlof S. F. ( 1986 ). Food Animal Residue Avoidance Databank (FARAD): An automated pharmacologic databank for drug and chemical residue avoidance . J. Food Prot. 49 , 826 – 830 . Google Scholar CrossRef Search ADS Riviere J. E. , Tell L. A. , Baynes R. E. , Vickroy T. W. , Gehring R. ( 2017 ). Guide to FARAD resources: Historical and future perspectives . J. Am. Vet. Med. Assoc. 250 , 1131 – 1139 . Google Scholar CrossRef Search ADS PubMed Riviere J. E. , Webb A. I. , Craigmill A. L. ( 1998 ). Primer on estimating withdrawal times after extralabel drug use . J. Am. Vet. Med. Assoc. 213 , 966 – 968 . Google Scholar PubMed Schadewinkel-Scherkl A. M. , Rasmussen F. , Merck C. C. , Nielsen P. , Frey H. H. ( 1993 ). Active transport of benzylpenicillin across the blood-milk barrier . Pharmacol. Toxicol. 73 , 14 – 19 . Google Scholar CrossRef Search ADS PubMed Shankaran H. , Adeshina F. , Teeguarden J. G. ( 2013 ). Physiologically-based pharmacokinetic model for Fentanyl in support of the development of Provisional Advisory Levels . Toxicol. Appl. Pharmacol. 273 , 464 – 476 . Google Scholar CrossRef Search ADS PubMed Smith G. W. , Gehring R. , Riviere J. E. , Yeatts J. L. , Baynes R. E. ( 2004 ). Elimination kinetics of ceftiofur hydrochloride after intramammary administration in lactating dairy cows . J. Am. Vet. Med. Assoc. 224 , 1827 – 1830 . Google Scholar CrossRef Search ADS PubMed Sterner T. R. , Ruark C. D. , Covington T. R. , Yu K. O. , Gearhart J. M. ( 2013 ). A physiologically based pharmacokinetic model for the oxime TMB-4: Simulation of rodent and human data . Arch. Toxicol. 87 , 661 – 680 . Google Scholar CrossRef Search ADS PubMed Suhren G. ( 1996 ). Influence of residues of antimicrobials in milk on commercially applied starter cultures-model trials . Kieler. Milchw. Forsch. 48 , 131 – 149 . Tan Y. M. , Liao K. H. , Conolly R. B. , Blount B. C. , Mason A. M. , Clewell H. J. III ( 2006 ). Use of a physiologically based pharmacokinetic model to identify exposures consistent with human biomonitoring data for chloroform . J. Toxicol. Environ. Health A 69 , 1727 – 1756 . Google Scholar CrossRef Search ADS PubMed United States Department of Agriculture (USDA) ( 2015 ). U.S. National Residue Program: 2014 Residue Sample Results. United States Department of Agriculture , Washington, DC . Available at: http://www.fsis.usda.gov/wps/wcm/connect/2428086b-f8ec-46ed-8531-a45d10bfef6f/2014-Red-Book.pdf?MOD=AJPERES. Last Accessed March 29, 2018. United States Department of Agriculture (USDA) ( 2017a ). U.S. National Residue Program: 2015 Residue Sample Results . United States Department of Agriculture , Washington, DC . Available at: https://www.fsis.usda.gov/wps/wcm/connect/f57333e5-9ff8-43ed-a787-6824f44bbac4/2015-red-book.pdf?MOD=AJPERES. Last Accessed March 29, 2018. United States Department of Agriculture (USDA) ( 2017b ). U.S. National Residue Program: 2016 Residue Sample Results . United States Department of Agriculture , Washington, DC . Available at: https://www.fsis.usda.gov/wps/wcm/connect/d84a5cac-5b4e-4e60-85b4-8886d0dc1660/2016-Red-Book.pdf?MOD=AJPERES. Last Accessed March 29, 2018. United States Pharmacopeia ( 2003 ). Penicillin G Veterinary—Intramammary-Local . United States Pharmacopeia , Rockville, MD . Available at: www.aavpt.org/resource/resmgr/imported/penicillinG-in.pdf. Last Accessed March 29, 2018. PubMed PubMed U.S. Food and Drug Administration (FDA) ( 2005 ). M-I-05-5: Tolerance and/or Safe Levels of Animal Drug Residues in Milk . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://wayback.archive-it.org/7993/20170113001136/http://www.fda.gov/Food/GuidanceRegulation/GuidanceDocumentsRegulatoryInformation/Milk/ucm077350.htm. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2006 ). Guidance for Industry—General Principles for Evaluating the Safety of Compounds Used in Food-Producing Animals . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/ohrms/dockets/98fr/05d-0219-gdl0001.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2010 ). Freedom of Information . NADA 065-010. Norocillin. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://wayback.archive-it.org/7993/20170406083042/https://www.fda.gov/downloads/AnimalVeterinary/Products/ApprovedAnimalDrugProducts/FOIADrugSummaries/UCM218419.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2012 ). Regulation of Carcinogenic Compounds Used in Food-Producing Animals. Code of Federal Regulations Title 21, Parts 500.82. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.accessdata.fda.gov/scripts/cdrh/cfdocs/cfCFR/CFRSearch.cfm?fr=500.82. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2013 ). Penicillin G Procaine Implantation and Injectable Dosage Forms. 21 CFR 522.1696. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.gpo.gov/fdsys/pkg/CFR-2013-title21-vol6/pdf/CFR-2013-title21-vol6-chapI-subchapE.pdf. Last Accessed March 29, 2018. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2015a ). Pasteurized Milk Ordinance 2015 Revision . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/downloads/food/guidanceregulation/guidancedocumentsregulatoryinformation/milk/ucm513508.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2015b ). Milk Drug Residue Sampling Survey . U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/downloads/animalveterinary/guidancecomplianceenforcement/complianceenforcement/ucm435759.pdf. Last Accessed March 29, 2018. PubMed PubMed U.S. Food and Drug Administration (FDA) ( 2016 ). Guidance for Industry—General Principles for Evaluating the Human Food Safety of New Animal Drugs Used In Food-Producing Animals. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.fda.gov/downloads/AnimalVeterinary/GuidanceComplianceEnforcement/GuidanceforIndustry/ucm052180.pdf. Last Accessed March 29, 2018. U.S. Food and Drug Administration (FDA) ( 2017 ). Tolerances for Residues of New Animal Drugs in Food . Code of Federal Regulations Title 21, Parts 556.510. U.S. Food and Drug Administration , Silver Spring, MD . Available at: https://www.accessdata.fda.gov/scripts/cdrh/cfdocs/cfcfr/CFRSearch.cfm?fr=556.510. Last Accessed March 31, 2018. Van OS J. , Burtelaar J. , Gouldsward J. ( 1974 ). Intramuscular treatment of bovine mastitis with various penicillins. Penicillin concentrations in the milk . Tijdschr. Diergeneeskd . 99 , 114 – 122 . Vetharaniam I. , Davis S. R. , Soboleva T. K. , Shorten P. R. , Wake G. C. ( 2003 ). Modeling the interaction of milking frequency and nutrition on mammary gland growth and lactation . J. Dairy Sci. 86 , 1987 – 1996 . Google Scholar CrossRef Search ADS PubMed Vilim A. , Larocque L. , Macintosh A. ( 1979 ). Depletion of brilliant blue F.C.F. and penicillin G in milk from treated cows . J. Food Prot. 42 , 491 – 494 . Google Scholar CrossRef Search ADS Vilim A. , Larocque L. , Macintosh A. ( 1980 ). Depletion of brilliant blue F.C.F., penicillin G, and dihydrostreptomycin in milk treated does with experimentally induced mastitis . J. Food Prot. 43 , 356 – 359 . Google Scholar CrossRef Search ADS Vogel G. , Nicolet J. , Martig J. , Tschudi P. , Meylan M. ( 2001 ). Pneumonia in calves: Characterization of the bacterial spectrum and the resistance patterns to antimicrobial drugs . Schweiz. Arch. Tierheilkd. 143 , 341 – 350 . Google Scholar PubMed Whittem T. , Whittem J. H. , Constable P. D. ( 2012 ). Modelling the concentration-time relationship in milk from cattle administered an intramammary drug . J. Vet. Pharmacol. Ther. 35 , 460 – 471 . Google Scholar CrossRef Search ADS PubMed Wishart D. S. , Knox C. , Guo A. C. , Shrivastava S. , Hassanali M. , Stothard P. , Chang Z. , Woolsey J. ( 2006 ). DrugBank: A comprehensive resource for in silico drug discovery and exploration . Nucleic Acids Res . 34 , D668 – D672 . Google Scholar CrossRef Search ADS PubMed World Health Organization (WHO) ( 2010 ). Characterization and Application of Physiologically Based Pharmacokinetic Models in Risk Assessment . World Health Organization, International Programme on Chemical Safety , Geneva, Switzerland . Yang F. , Liu H. W. , Li M. , Ding H. Z. , Huang X. H. , Zeng Z. L. ( 2012 ). Use of a Monte Carlo analysis within a physiologically based pharmacokinetic model to predict doxycycline residue withdrawal time in edible tissues in swine . Food Addit. Contam. Part A Chem. Anal. Control Expo. Risk Assess . 29 , 73 – 84 . Google Scholar CrossRef Search ADS PubMed Yang X. , Doerge D. R. , Teeguarden J. G. , Fisher J. W. ( 2015 ). Development of a physiologically based pharmacokinetic model for assessment of human exposure to bisphenol A . Toxicol. Appl. Pharmacol. 289 , 442 – 456 . Google Scholar CrossRef Search ADS PubMed Yoon J. , Lee L. , Kim C. , Chung Y. , Kim C. ( 2004 ). Effects of milk production, season, parity and lactation period on variations of milk urea nitrogen concentration and milk components of Holstein dairy cows . Asian-Australas. J. Anim. Sci. 17 , 479 – 484 . Google Scholar CrossRef Search ADS Yoon M. , Nong A. , Clewell H. J. III , Taylor M. D. , Dorman D. C. , Andersen M. E. ( 2009 ). Evaluating placental transfer and tissue concentrations of manganese in the pregnant rat and fetuses after inhalation exposures with a PBPK model . Toxicol. Sci. 112 , 44 – 58 . Google Scholar CrossRef Search ADS PubMed Yoon M. , Schroeter J. D. , Nong A. , Taylor M. D. , Dorman D. C. , Andersen M. E. , Clewell H. J. III ( 2011 ). Physiologically based pharmacokinetic modeling of fetal and neonatal manganese exposure in humans: Describing manganese homeostasis during development . Toxicol. Sci. 122 , 297 – 316 . Google Scholar CrossRef Search ADS PubMed Zeng D. , Lin Z. , Fang B. , Li M. , Gehring R. , Riviere J. E. , Zeng Z. ( 2017 ). Pharmacokinetics of mequindox and its marker residue 1, 4-bisdesoxymequindox in swine following multiple oral gavage and intramuscular administration: An experimental study coupled with population physiologically based pharmacokinetic modeling . J. Agric. Food Chem. 65 , 5768 – 5777 . Google Scholar CrossRef Search ADS PubMed Zhu X. , Huang L. , Xu Y. , Xie S. , Pan Y. , Chen D. , Liu Z. , Yuan Z. ( 2017 ). Physiologically based pharmacokinetic model for quinocetone in pigs and extrapolation to mequindox . Food Addit. Contam. Part A 34 , 192 – 210 . Google Scholar CrossRef Search ADS Ziv G. , Rasmussen F. ( 1975 ). Distribution of labeled antibiotics in different components of milk following intramammary and intramuscular administrations . J. Dairy Sci . 58 , 938 – 946 . Google Scholar CrossRef Search ADS PubMed Ziv G. , Sulman F. G. ( 1974 ). Effects of probenecid on the distribution, elimination, and passage into milk of benzylpenicillin, ampicillin and cloxacillin . Arch. Int. Pharmacodyn. Ther. 207 , 373 – 382 . Google Scholar PubMed © The Author(s) 2018. Published by Oxford University Press on behalf of the Society of Toxicology. All rights reserved. For permissions, please e-mail: journals.permissions@oup.com This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/about_us/legal/notices)

Journal

Toxicological SciencesOxford University Press

Published: Mar 21, 2018

There are no references for this article.

You’re reading a free preview. Subscribe to read the entire article.


DeepDyve is your
personal research library

It’s your single place to instantly
discover and read the research
that matters to you.

Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month

Explore the DeepDyve Library

Search

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Organize

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Access

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Your journals are on DeepDyve

Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.

All the latest content is available, no embargo periods.

See the journals in your area

DeepDyve

Freelancer

DeepDyve

Pro

Price

FREE

$49/month
$360/year

Save searches from
Google Scholar,
PubMed

Create lists to
organize your research

Export lists, citations

Read DeepDyve articles

Abstract access only

Unlimited access to over
18 million full-text articles

Print

20 pages / month

PDF Discount

20% off