Mechanistic and Quantitative Understanding of Pharmacokinetics in Zebrafish Larvae through Nanoscale Blood Sampling and Metabolite Modeling of Paracetamol ========================================================================================================================================================== * Rob C. Van Wijk * Elke H.J. Krekels * Vasudev Kantae * Anita Ordas * Thijs Kreling * Amy C. Harms * Thomas Hankemeier * Herman P. Spaink * Piet H. van der Graaf ## Abstract Zebrafish larvae are increasingly used for pharmacological research, but internal drug exposure is often not measured. Understanding pharmacokinetics is necessary for reliable translation of pharmacological results to higher vertebrates, including humans. Quantification of drug clearance and distribution requires measurements of blood concentrations. Additionally, measuring drug metabolites is of importance to understand clearance in this model organism mechanistically. We therefore mechanistically studied and quantified pharmacokinetics in zebrafish larvae, and compared this to higher vertebrates, using paracetamol (acetaminophen) as a paradigm compound. A method was developed to sample blood from zebrafish larvae 5 days post fertilization. Blood concentrations of paracetamol and its major metabolites, paracetamol-glucuronide and paracetamol-sulfate, were measured. Blood concentration data were combined with measured amounts in larval homogenates and excreted amounts and simultaneously analyzed through nonlinear mixed-effects modeling, quantifying absolute clearance and distribution volume. Blood sampling from zebrafish larvae was most successful from the posterior cardinal vein, with a median volume (interquartile range) of 1.12 nl (0.676–1.66 nl) per blood sample. Samples were pooled (*n* = 15–35) to reach measurable levels. Paracetamol blood concentrations at steady state were only 10% of the external paracetamol concentration. Paracetamol-sulfate was the major metabolite, and its formation was quantified using a time-dependent metabolic formation rate. Absolute clearance and distribution volume correlated well with reported values in higher vertebrates, including humans. Based on blood concentrations and advanced data analysis, the mechanistic and quantitative understanding of paracetamol pharmacokinetics in zebrafish larvae has been established. This will improve the translational value of this vertebrate model organism in drug discovery and development. **SIGNIFICANCE STATEMENT** In early phases of drug development, new compounds are increasingly screened in zebrafish larvae, but the internal drug exposure is often not taken into consideration. We developed innovative experimental and computational methods, including a blood-sampling technique, to measure the paradigm drug paracetamol (acetaminophen) and its major metabolites and quantify pharmacokinetics (absorption, distribution, elimination) in zebrafish larvae of 5 days post fertilization with a total volume of only 300 nl. These parameter values were scaled to higher vertebrates, including humans. ## Introduction In drug discovery and development, the zebrafish (*Danio rerio*) larva is a promising vertebrate model organism (Zon and Peterson, 2005; MacRae and Peterson, 2015). It has many advantages, including 70% genetic homology to humans (Howe et al., 2013), the possibility of high-throughput experimentation, easy genetic modification, and imaging of internal processes due to its transparency in early life ((Van Wijk et al., 2016; Schulthess et al., 2018)). In experiments with zebrafish larvae, it is common practice to expose the larvae to study drugs by waterborne treatment, dissolving the drug of interest in the medium in which the larvae swim. Ignoring the internal exposure, as is done in this approach, potentially leads to false positives or negatives (Diekmann and Hill, 2013). More importantly, translation of pharmacological findings to higher vertebrates is very limited without an internal exposure-response relationship ((Morgan et al., 2012; Van Wijk et al., 2016)). To be able to use zebrafish larvae to their full potential in pharmacological research, a mechanistic understanding of the pharmacokinetic processes in zebrafish larvae is needed together with a framework to quantitatively scale this to higher vertebrates. Methods to measure internal drug amounts in larval homogenates over time have been established only recently (Kühnert et al., 2013, 2017; Vogs et al., 2015; Kantae et al., 2016). In the case of our paradigm compound paracetamol (acetaminophen), internal paracetamol amounts were subsequently used to quantify the pharmacokinetic processes of absorption and elimination through quantification of the rate of absorption and a clearance relative to total larval volume. Although the clearance derived from relative values in zebrafish larvae scaled quite reasonably to reported absolute clearance values in higher vertebrates ((Kantae et al., 2016; Van Wijk et al., 2019)), accurate scaling of clearance requires absolute clearance values in zebrafish larvae. However, quantification of absolute clearance requires blood concentrations, which will also allow for the estimation of a distribution volume, but blood-sampling methods for zebrafish larvae that are only a few hundred nanoliters in total volume (Guo et al., 2017) have not yet been described. In addition, mechanistic characterization of metabolic clearance in zebrafish larvae is of importance, as the utility of zebrafish larvae in drug development depends on their drug metabolism being similar to higher vertebrates in a qualitative (similar metabolites) and quantitative (similar rate and extent) manner. This is because drug metabolites can be pharmacologically active or toxic. A first step in assessing mechanistic similarities in drug metabolism between vertebrate species is genetic confirmation, through DNA sequencing and gene expression of the enzymes and possible cofactors, that similar enzymatic pathways are present in the zebrafish. For example, our paradigm compound paracetamol is a substrate for two important metabolic pathways, glucuronidation and sulfation (Janus et al., 2003; Neirinckx et al., 2010; Lee et al., 2012; Owens et al., 2012). In zebrafish, the enzyme systems for these pathways have been documented based on genetic homology with humans ((Van den Boom et al., 2012; Van Wijk et al., 2016)). However, functional and quantitative confirmation that drug metabolites are actually formed *in vivo* and at what rate they are formed and eliminated is crucial. This is currently lacking for most enzymatic pathways that metabolize drugs in zebrafish larvae. Here, the objective is to mechanistically and quantitatively study pharmacokinetics in zebrafish larvae. To that aim, a method is developed to take nanoscale blood samples to measure paracetamol and its major metabolites. A mechanistic metabolite model is developed to quantify paracetamol absolute clearance and distribution volume, which are compared with those from higher vertebrates, including humans. ## Materials and Methods #### Experimental Design. Experiments were performed in zebrafish (*D. rerio*) larvae of 5 days post fertilization (dpf), as this age is within the ethically acceptable time frame for larval experiments, and it offers the largest capacity to metabolize paracetamol within that time frame (Van Wijk et al., 2019). The investigation was divided into two experiments. In the first experiment, larvae were continuously treated with 1 mM waterborne paracetamol concentration in embryo medium for up to 200 minutes. In the second experiment, larvae were treated with 1 mM waterborne paracetamol concentration for 60 minutes, then washed and transferred to clean washout embryo medium, and elimination was studied over 240 minutes. A new method was developed to sample blood from zebrafish larvae. Blood samples were taken during the first experiment. Observations between 10 and 125 minutes of waterborne treatment consisted of minimal three replicates of 15–35 pooled blood samples. These data were combined with data from a previous study, which included six replicates of measured paracetamol amounts in whole larval homogenates of five zebrafish larvae at 10–180 minutes in the first experiment and 60–300 minutes in the second experiment (Van Wijk et al., 2019). Additionally, six replicates of washout medium of the second experiment were sampled at 60–300 minutes to measure excreted paracetamol and metabolites. Paracetamol and its two major metabolites, paracetamol-glucuronide and paracetamol-sulfate, were measured in blood, medium, and homogenate samples by liquid chromatography–mass spectrometry (LC-MS). Nonlinear mixed-effects modeling was performed on all data simultaneously to quantify absolute clearance, including metabolic formation and elimination rates and distribution volume. Finally, obtained parameter values for absolute clearance and distribution volume were compared with published values in higher vertebrates. #### Zebrafish Larvae Husbandry. Maintenance and handling of zebrafish followed international consensus protocols (Westerfield, 2000), and the planning and execution of all experiments were compliant with European regulation (European Union, 2010). Adult wild-type AB/TL zebrafish were used to fertilize eggs and kept in glass aquaria (max 6/l, volume 10 l, 120 × 220 × 490 mm; Fleuren & Nooijen BV, Nederweert, The Netherlands) with circulating water at 27.7 ± 0.1°C on a 14-hour/10-hour light/dark cycle (lights on at 08:00) and twice daily feeding with artemia or feed particles (Gemma Micro/Diamond, Skretting; Nutreco NV, Amersfoort, The Netherlands). Water quality was controlled by JUMO Acquis touch S (JUMO GmbH & Co., Weesp, The Netherlands). Fertilized eggs were collected within 20 minutes of fertilization. Eggs and larvae were kept at 28°C in embryo medium, which was refreshed daily. Exposure experiments were performed at room temperature. #### Blood Sampling. To develop a method for blood sampling from zebrafish larvae at 5 dpf, larvae were washed with embryo medium using Netwell insert filters (Corning Life Sciences B.V., Amsterdam, The Netherlands) and transferred to a microscope slide coated with agarose. Larvae were not anesthetized to prevent decreased blood flow; instead, superficial drying with soft lens paper prevented movement. Sampling was performed using a needle pulled by a micropipette puller (Sutter Instruments, Novato, CA) from a 0.75-mm borosilicate glass capillary (Sutter Instruments) without filament, positioned in a micromanipulator (World Precision Instruments, Berlin, Germany), connected to a manual CellTram Vario oil pump (Eppendorf Nederland B.V., Nijmegen, The Netherlands), under 20× magnification (Leica, Amsterdam, The Netherlands). After decapitation or tail cut was found to result in limited yields, direct sampling from the circulation was explored. To identify the most suitable location, sampling from different anatomic locations was tested (Isogai et al., 2001), including the heart, the dorsal aorta, the caudal vein, and the posterior cardinal vein. An image was taken from each blood sample in the needle to determine the sample volume. The image analysis program Fiji (Schindelin et al., 2012) was used to calculate the volume of the blood sample within the needle with the volume formula of a truncated cone (eq. 1):![Formula][1](1)where *r*1 and *r*2 are the upper and lower radii of the blood sample in the pulled needle, and *h* is the length of the blood sample within the needle taking into account the 45° angle of the needle with respect to the microscope. After the blood sample was imaged, it was injected into a 2-*µ*l droplet of heparin solution (5 International Units/ml) under microscope both to prevent coagulation and for sample handling. Sample replicates consisted of blood samples from 15–35 larvae to reach measurable levels, of which the total number depended on experimental and time constraints, and were pooled into a 0.5-ml Eppendorf tube and stored at −80°C. #### Measurements of Paracetamol and Its Metabolites. On the day of measurement, blood samples were thawed and 200 *µ*l of methanol with paracetamol-D4 internal standard (1.4 pg/*µ*l) was added. Samples were centrifuged (16,000*g*, 10 minutes) and 180 *µ*l of supernatant was transferred to a 0.5-ml Eppendorf tube to be evaporated by vacuum centrifuge (Beun-de Ronde, Abcoude, The Netherlands) until dry. Samples were reconstituted in 10 *µ*l of 80/20 (v/v) purified water/methanol and transferred to an LC-MS vial for randomized injection of 7 *µ*l into the ultraperformance liquid chromatography system (Acquity; Waters Chromatography B.V., Etten-Leur, The Netherlands) linked to a quadrupole ion trap MS/MS (QTRAP-6500; AB Sciex B.V., Nieuwekerk aan den IJssel, The Netherlands). An electrospray ionization source in positive (paracetamol) and negative (paracetamol-glucuronide and paracetamol-sulfate) mode was used as published before ((Kantae et al., 2016; Van Wijk et al., 2019)). Method criteria were 90%–100% accuracy and a precision corresponding to a relative S.D. less than 10%. Lower limit of quantification (LLOQ) in blood samples was 0.05 pg/*μ*l for paracetamol and paracetamol-sulfate and 5 pg/*μ*l for paracetamol-glucuronide. The homogenate samples were measured as described before ((Kantae et al., 2016; Van Wijk et al., 2019)). LLOQ in homogenates was 0.09 pg/*μ*l for paracetamol and paracetamol-sulfate and 9.0 pg/*μ*l for paracetamol-glucuronide. Of the washout medium of the second experiment, samples of 1850 *µ*l were collected and stored at −80°C. On the day of measurement, samples were thawed, 10 *µ*l of 125-pg/*µ*l paracetamol-D4 internal standard was added, and samples were evaporated by vacuum centrifuge until dry. Samples were reconstituted in 50 *µ*l of purified water and centrifuged (16,000*g*, 15 minutes), and 25 *µ*l of supernatant was transferred to an LC-MS vial for randomized injection of 5 *µ*l into the LC-MS system as described earlier. Background measurement of the washout medium at *t* = 60 was subtracted from the sample measurements. A calibration curve ranging from 0.05 to 100 pg/*µ*l was prepared in 50/50 (v/v) methanol/purified water and was used to calculate compound-excreted amounts in picomoles per larva. LLOQ in washout medium was 0.05 pg/*µ*l for paracetamol and paracetamol-sulfate and 5 pg/*µ*l for paracetamol-glucuronide. #### Pharmacokinetic Data Analysis. A pharmacokinetic model was developed for paracetamol and its metabolites using nonlinear mixed-effects modeling. NONMEM (version 7.3) (Beal et al.) was used through interfaces Pirana (version 2.9.6) (Keizer et al., 2011) and PsN (version 4.7.0) (Lindbom et al., 2005), and graphical output was created using R (version 3.4.2) (R Core Team., 2014) through the RStudio interface (version 1.1.383; RStudio Inc., Boston, MA). The first-order conditional estimation algorithm was used. A zero-order absorption rate constant was estimated to quantify paracetamol absorption from the surrounding medium, as we have previously shown that paracetamol concentrations in the treatment medium remain constant during the experiments (Van Wijk et al., 2019). One- and two-compartment models were tested for the distribution of each compound. Linear and nonlinear metabolic clearance for paracetamol was tested, including one- and two-substrate Michaelis-Menten kinetics (Ingalls, 2012), the latter of which was used to account for possible saturation of the sulfation pathway caused by depletion of the sulfate group donor 3′-phosphoadenosine 5′-phosphosulfate (PAPS) (Reddyhoff et al., 2015). Last, time-dependent metabolism using both exponential and sigmoidal relationships was tested (Brochot et al., 2011). First-order rates of excretion into medium with or without recovery fractions were estimated per compound. Because of the destructive sampling, interindividual variability could not be distinguished from residual variability; therefore, only the latter was estimated, using an additive, proportional, or a combination of additive and proportional error per compound per sample type. Level 2 covariance between the residual error for compounds measured in the same sample was tested (Karlsson et al., 1995). All compounds had less than 10% of data points below the limit of quantification, with the sole exception of paracetamol-glucuronide in blood samples (71%), warranting the M2 method of ignoring these values in the analysis (Beal, 2001). Model selection was based on the likelihood ratio test between nested models, using a drop in objective function value of 6.63 between nested models to indicate statistical significance, corresponding to *P* < 0.01 assuming a *χ*2 distribution. Additionally, physiologic plausibility of parameter estimates and goodness-of-fit plots were assessed (Nguyen et al., 2017). Structural parameter precision was considered acceptable when relative standard errors (RSE) were below 50%. #### Comparison of Absolute Clearance and Distribution Volume of Paracetamol to Higher Vertebrates. The estimate of absolute paracetamol clearance was compared with values reported in higher vertebrates as previously published (Kantae et al., 2016), as was the estimate of the distribution volume. For the latter, similar to the comparison of clearance values, a literature search was performed in PubMed using “distribution volume OR volume of distribution AND paracetamol OR acetaminophen.” Only values published in the last 20 years were included. Values obtained from specific disease models, from combination treatment, or from obese patients, and from models with more than one compartment for paracetamol distribution were excluded. A regression through the log-transformed parameter values based on log-transformed bodyweight was calculated including 95% confidence interval using R (v.3.4.2) (R Core Team., 2014) through user interface RStudio (v.1.1.383; RStudio Inc.). This regression was based on data from mature individuals only, comparable to that of paracetamol clearance. The bodyweight of the larvae was calculated from their previously published volume (Guo et al., 2017) and an assumed density of 0.997 g/ml, corresponding to that of water (Kantae et al., 2016). #### Materials. Paracetamol, paracetamol-glucuronide, paracetamol-sulfate, and paracetamol-D4 internal standard were acquired from Sigma-Aldrich Chemie B.V. (Zwijndrecht, The Netherlands). Heparin was acquired from Academic Hospital Pharmacy Leiden (Leiden University Medical Center, Leiden, The Netherlands). Embryo medium consisted of demineralized water containing 60 *µ*g/ml Instant Ocean sea salts (Sera, Heinsberg, Germany). Agarose was acquired from Sphaero (Gorinchem, The Netherlands). Ultraperformance liquid chromatography/mass spectrometry–grade MeOH was acquired from Biosolve B.V. (Valkenswaard, The Netherlands). PURELAB (Veolia Water Technologies B.V., Ede, The Netherlands) was used to purify water. ## Results #### Blood Sampling. From the different anatomic locations, the posterior cardinal vein proved most efficient for blood sampling and resulted in a median volume of 1.12 nl per larva with an interquartile range of 0.676–1.66 nl. Figure 1 and Supplementary Video showthis procedure. The sample volume was determined based on a separately captured image of the blood sample within the needle. [Supplemental Table 1](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1) lists the number blood samples that were pooled for each sampling time point. ![Fig. 1.](http://jpet.aspetjournals.org/http://jpet.aspetjournals.org/content/jpet/371/1/15/F1.medium.gif) [Fig. 1.](http://jpet.aspetjournals.org/content/371/1/15/F1) Fig. 1. Blood sampling from posterior cardinal vein in a zebrafish larva at 5 days post fertilization using a pulled needle. Blood sample is shown in needle tip. #### Measurement of Paracetamol and Metabolites. Paracetamol, paracetamol-glucuronide, and paracetamol-sulfate could be measured in the blood samples, of which representative chromatograms are shown in [Supplemental Fig. 1](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1). The symbols in Fig. 2 show paracetamol and metabolite concentrations and total amounts in blood samples and in homogenates and washout medium, respectively. Although the data are variable, clear trends can be observed. ![Fig. 2.](http://jpet.aspetjournals.org/http://jpet.aspetjournals.org/content/jpet/371/1/15/F2.medium.gif) [Fig. 2.](http://jpet.aspetjournals.org/content/371/1/15/F2) Fig. 2. Paracetamol (blue circles, left column), paracetamol-glucuronide (magenta triangles, middle column), and paracetamol-sulfate (orange squares, right column) observations in blood samples (top row, concentrations) as well as in homogenates (middle row, amounts) and excreted into medium (bottom row, amounts) after constant waterborne treatment (solid symbols) and washout (open symbols) experiments. Blood concentrations are shown in picomoles per nanoliters (millimolar). The model prediction is shown as a solid or dotted line for the constant waterborne treatment or washout experiment, respectively. The lower number of data points of paracetamol-glucuronide is explained by its lower limit of quantification being higher than that of the other compounds. Paracetamol observations in homogenates have been reported previously (Van Wijk et al., 2019). The paracetamol and metabolite concentrations (Fig. 2, top panel) could be measured in the nanoscale blood samples and show that the sulfate metabolite reaches similar concentrations as the parent, while concentrations of the glucuronide metabolite are about 10-fold lower. This is similar to the observations in homogenates (Fig. 2, middle panel). Paracetamol concentrations reach a steady-state value around 0.1 mM or 15 mg/l, which is 10% of the paracetamol concentration in the treatment medium. The paracetamol amounts in larval homogenate from the washout experiment (Fig. 2, open symbols) show a monoexponential decline of paracetamol. For the metabolites, the washout experiment shows an increase in internal amounts between 60 and 120 minutes, after which paracetamol-glucuronide decreases to 1 pmol/larva and paracetamol-sulfate stabilizes between 10 and 20 pmol/larva at 300 minutes. The amounts of paracetamol, paracetamol-glucuronide, and paracetamol-sulfate excreted into the medium during the washout experiment are shown in Fig. 2 (bottom panel). These data show low values with high variability and a minor increase over the washout period. #### Pharmacokinetic Data Analysis. Figure 3 shows a schematic representation of the developed parent-metabolite model for paracetamol. A model with zero-order absorption and one compartment per compound best fit the data. The formation of paracetamol-glucuronide was best described with a first-order metabolic formation rate, whereas the paracetamol sulfation rate could not be captured with a zero- or first-order metabolic formation rate, as both methods led to underpredictions at earlier time points and overprediction at later time points. Conventional Michaelis-Menten kinetics depending on parent compound concentrations could not be identified for paracetamol-sulfate, nor could a two-substrate Michaelis-Menten kinetic model, which represents the potential rate-limiting influence of depletion of the sulfate group donor as well. Regarding functions based on time-related changes in sulfation rate, a sigmoidal function as described by eq. 2 was statistically significantly better than an exponentially declining time-dependent sulfation rate and provided the best fit:![Formula][2](2)where *k*PS,f is the metabolic formation rate for sulfation in min−1, *k*PS,f,0 is the metabolic formation rate for sulfation at time point 0 in min−1, *t* is the time in minutes, and *t*50 is the time in minutes at which the formation rate for the sulfate metabolite is at 50% of its value at time 0. ![Fig. 3.](http://jpet.aspetjournals.org/http://jpet.aspetjournals.org/content/jpet/371/1/15/F3.medium.gif) [Fig. 3.](http://jpet.aspetjournals.org/content/371/1/15/F3) Fig. 3. Schematic representation of the structural pharmacokinetic model of paracetamol and its metabolites in zebrafish larvae. Compartments depicted with solid lines represent compartments in which both amounts in homogenized larvae or blood concentrations in larvae were available. Compartments indicated with dashed lines represent compartments in which excreted amounts in the medium were available. *k*a, absorption rate constant; *k*P,e, paracetamol excretion rate; *k*PS,f, metabolic formation rate for sulfation; *k*PG,f, metabolic formation rate for glucuronidation; *k*G,e, paracetamol-glucuronide excretion rate; *k*S,e, paracetamol-sulfate excretion rate; P, paracetamol; G, paracetamol-glucuronide; S, paracetamol-sulfate. The excretion of paracetamol and its metabolites was best captured with a first-order rate constant. A recovery fraction on the excreted paracetamol and metabolite amounts was included to prevent overestimation of the excreted compounds, which resulted in a significantly improved fit both graphically and statistically. Figure 4 shows the paracetamol clearance over time as a sum of the constant metabolic formation rate for glucuronidation, the excretion rate, and the time-dependent metabolic formation rate for sulfation according to the model. ![Fig. 4.](http://jpet.aspetjournals.org/http://jpet.aspetjournals.org/content/jpet/371/1/15/F4.medium.gif) [Fig. 4.](http://jpet.aspetjournals.org/content/371/1/15/F4) Fig. 4. Paracetamol clearance over time. Total clearance (transparent green solid line) is the sum of the time-dependent metabolic formation rate for sulfation (orange dashed line), the metabolic formation rate for glucuronidation (magenta dotted line), and the unchanged excretion rate (blue solid line). Residual variability for homogenate samples was found to be best described with a combination of additive and proportional error for paracetamol and paracetamol-glucuronide and with a proportional error for paracetamol-sulfate, including covariance between the proportional errors of paracetamol and paracetamol-glucuronide and between the proportional errors of paracetamol-glucuronide and paracetamol-sulfate. Residual variability for blood samples was best described with an additive error for paracetamol and proportional errors for paracetamol-glucuronide and for paracetamol-sulfate. The residual variability for medium samples was best described by an additive error for paracetamol, a proportional error for paracetamol-glucuronide, and a combination of additive and proportional for paracetamol-sulfate. The model predictions (lines) relative to the observed values (symbols) are shown in Fig. 2 and show a good description of the data. Table 1 contains the parameter estimates including RSE as a measure of their precision of the final model. RSEs for structural parameters were acceptable. Additional goodness-of-fit-plots, showing good model accuracy, can be found in [Supplemental Figs. 2](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1)–[4](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1). View this table: [TABLE 1](http://jpet.aspetjournals.org/content/371/1/15/T1) TABLE 1 Parameter estimates for the final model The final model code and data set are available through the Drug Disease Model Resources (DDMoRe) Repository, Model ID DDMODEL00000300 ([http://repository.ddmore.foundation/model/DDMODEL00000300](http://repository.ddmore.foundation/model/DDMODEL00000300)). #### Comparison of Absolute Clearance and Distribution Volume of Paracetamol to Higher Vertebrates. A correlation between paracetamol clearance and bodyweight in 12 higher vertebrates, including humans, has been reported before (calculated exponent = 0.78, 95% confidence interval: 0.65–0.91, *R*2 = 0.90) ((Kantae et al., 2016; Van Wijk et al., 2019)). Paracetamol clearance estimated from homogenate data alone, assuming homogeneous distribution of the drug throughout the whole larvae, fell outside of the 95% confidence interval of this correlation (Van Wijk et al., 2019). Here, absolute paracetamol clearance is estimated by simultaneously modeling observed concentrations from blood samples, which could be quantified as the result of the blood-sampling method developed here, and observed internal amounts from homogenates. Figure 5 shows the estimated absolute paracetamol clearance of the zebrafish larva in relation to that of 12 higher vertebrates. Because of the time-dependent metabolic formation rate for sulfation, the full range from paracetamol clearance at *t* = 0 (solid square) to paracetamol clearance at *t* = ∞ (open square) is shown. This range is partly within the 95% confidence interval. Noticeably, clearance values from higher immature vertebrates are also overpredicted by the correlation based on the values reported for mature organisms. ![Fig. 5.](http://jpet.aspetjournals.org/http://jpet.aspetjournals.org/content/jpet/371/1/15/F5.medium.gif) [Fig. 5.](http://jpet.aspetjournals.org/content/371/1/15/F5) Fig. 5. Relationship between paracetamol clearance and bodyweight of 13 vertebrate species. Clearance in zebrafish larvae at 5 days post fertilization is depicted as the clearance value at *t* = 0 (closed square) and *t* = ∞ (open square), with the vertical line depicting the range in clearance over time. Mature individuals are shown in blue, immature individuals in red. The correlation (dashed line; shaded area, 95% confidence interval) is based on values from mature individuals of the higher vertebrates only. Adapted with permission from Kantae et al. (2016). The literature search resulted in reported values for the distribution volume of paracetamol in 13 higher vertebrates, which are provided with their references in [Supplemental Table 2](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1). A correlation between the distribution volume of paracetamol and bodyweight based on data obtained from mature individuals of the higher vertebrates only is shown in Fig. 6 (dashed line including 95% confidence interval). The exponent in this relationship was calculated to be 0.94 (95% confidence interval: 0.67–1.2, *R*2 = 0.68). Although 5 orders of magnitude removed from the smallest higher vertebrate, the distribution volume of paracetamol in zebrafish larvae lies well within the 95% confidence interval of the correlation between distribution volume and bodyweight. Even though the values from immature organisms are not taken into account for the correlation between the distribution volume of paracetamol and bodyweight, the data points are in agreement with that correlation. ![Fig. 6.](http://jpet.aspetjournals.org/http://jpet.aspetjournals.org/content/jpet/371/1/15/F6.medium.gif) [Fig. 6.](http://jpet.aspetjournals.org/content/371/1/15/F6) Fig. 6. Relationship between distribution volume of paracetamol and bodyweight of 14 vertebrate species, including zebrafish larvae at 5 days post fertilization. Mature individuals are shown in blue, immature individuals in red. The correlation (dashed line; shaded area, 95% confidence interval) is based on values from mature individuals in higher vertebrates only. Raw data can be found in [Supplemental Table 2](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1). ## Discussion Mechanistic and quantitative characterization of pharmacokinetics in zebrafish larvae is crucial for this vertebrate model organism to deliver on its promising role in pharmacological research. To date, a method to take blood samples, which are essential to characterize pharmacokinetics, was lacking. The quantification of absolute clearance and distribution volume in zebrafish larvae, needed for reliable extrapolation of pharmacokinetics between species, required measurements of these blood concentrations in addition to knowing drug amounts in larval homogenates as a determinant of the bioavailable dose. This is analogous to the well known requirement of knowing both blood concentrations and absolute administered drug amounts by intravenous drug administration to quantify absolute pharmacokinetic parameters of extravascular-administered drugs. Additionally, confirmation of functional drug metabolism is of importance, because drug metabolites can be pharmacologically active or toxic. We report, for the first time, a method for nanoscale blood sampling from zebrafish larvae at five dpf, in which not only the paradigm compound paracetamol but also its two major metabolites were measured. Absolute clearance of paracetamol in zebrafish larvae correlated reasonably well with that reported in higher vertebrates. Results from our study suggest that, in addition to the presence of genetic homology for the enzymatic systems established before (Van Wijk et al., 2016), these metabolic pathways also behave functionally similarly to higher vertebrates, as both metabolites are present. Quantitatively, the ratio between the two metabolites was different from those reported in mature humans (Reith et al., 2009). Paracetamol-sulfate proved to be the major metabolite at this larval stage with respect to paracetamol-glucuronide. A similar profile, however, has been shown in human newborns and infants, resulting from limited glucuronidation capacity, which increases with increasing age (Anderson et al., 2002; Krekels et al., 2015). This is not unexpected, as zebrafish larvae have also not reached maturity. Consequently, the impact of enzymatic maturation needs to be taken into account when translating results from immature zebrafish larvae to mature individuals of higher vertebrates. Metabolic formation rate for sulfation could not be accurately captured using time-constant parameters. A conventional Michaelis-Menten model could not accurately describe the observed pharmacokinetic profiles. Therefore, a two-substrate Michaelis-Menten kinetic model was tested, accounting for the concentrations of both the substrate paracetamol and a hypothetical donor group, representing the sulfate group donor PAPS, which may limit the biotransformation (Slattery et al., 1987; Liu and Klaassen, 1996; Reddyhoff et al., 2015). The data, however, did not allow for reliable and stable estimation of the relevant model parameters. Therefore, an empirical time-dependent function was tested, with a sigmoidal time-dependent metabolic formation rate resulting in the best fit but with relatively high RSE values of the estimated parameters. Normalization of the metabolic formation rate for sulfation at time point 0, *k*PS,f,0, at a different time point reduced its RSE value below 50%, indicating no overparameterization (Goulooze et al., 2019). [Supplemental Fig. 5](http://jpet.aspetjournals.org/lookup/suppl/10.1124/jpet.119.260299/-/DC1) shows the paracetamol-sulfate fit of the metabolite model, in which the time-dependent metabolic formation rate for sulfation was replaced by a time-constant first-order rate. The parameters of that model are estimated with acceptable precision (RSE < 50%), but a clear misfit of the paracetamol-sulfate data can be observed. The low *t*50 estimate resulted in a steep decline of sulfation in the first minutes of the experiment (Fig. 4). We believe this function represents the impact of depletion of the sulfate group donor PAPS on the metabolic formation rate, as the shape of this curve is similar to the shape of simulated PAPS levels from a reported physiologic mathematical model at the cellular level (Reddyhoff et al., 2015). As time approached infinity, sulfation nears the asymptote of zero metabolic formation rate. This is unlikely, as, for example, the sulfate group donor PAPS is expected to be produced by the organism at a constant rate, but this production rate could not be estimated with the current data and time scale (Reddyhoff et al., 2015). Consequently, interpretation of this model is limited to the time course presented here. The availability of blood concentrations in addition to absolute internal amounts allowed for the estimation of absolute values for distribution volume for both the parent and the major metabolites, which is unique. Figure 2 shows variable data yielded from extremely small samples and low concentrations; however, this did not limit accurate estimation of distribution volumes, as indicated by RSE values being 25% or lower. The estimated distribution volume of paracetamol in zebrafish larvae correlated well with those reported in higher vertebrates. The estimated distribution volumes for the glucuronide and sulfate metabolites are 0.6× and 1.2× the distribution volume of paracetamol, respectively. Reported values show large variability and range 0.2–1.6× that of paracetamol for paracetamol glucuronide and 0.09–0.9× that of paracetamol for paracetamol-sulfate in humans ((Lowenthal et al., 1976; Owens et al., 2012, 2014; Van Rongen et al., 2016)), and 1.2–2.7× for paracetamol glucuronide and 0.25–1.5× for paracetamol-sulfate in rats (Yamasaki et al., 2011; Lee et al., 2012). These values are in line with our estimates, but care must be taken when interpreting these data. Determining distribution volume of drug metabolites had to be done in patients who all had specific conditions which also impacted normal physiology. In the zebrafish, the total amounts necessary were readily available for that determination. The developed method of blood sampling results in blood samples at the nanoliter scale. Multiple blood samples have to be pooled into a single analytical replicate to quantify drug and metabolite concentrations by conventional LC-MS. To improve throughput and efficiency, nanoelectrospray ionization mass spectrometry could be considered (Fujii et al., 2015). This technique directly injects a sample of sub-nanoliter volumes together with ionization solvent into an electron spray ionization source and has, for instance, been used to sample and quantify concentrations of the anticancer drug tamoxifen in individual cancer cells (Ali et al., 2019). Currently, the high-throughput potential in zebrafish larvae for the purpose of obtaining absolute pharmacokinetic parameter values required for interspecies scaling is limited by the practicalities of the blood sampling and by the need for the development of highly sensitive analytical methods for these samples. We, however, envision that this step is only performed for promising drug candidates that have been selected for further development. Additionally, microinjection of zebrafish zygotes and embryos has been automated, resulting in higher throughput (Spaink et al., 2013; Cordero-Maldonado et al., 2019). In the future, this technique might also be applicable to sampling. We have confirmed that the paracetamol metabolites are formed by metabolizing enzymes of the zebrafish larvae and not by metabolizing enzymes of the microbiota by repeating the experiment with germ-free larvae that lack microbiota (data not shown). To establish whether paracetamol-glucuronide and paracetamol-sulfate were the major metabolites in zebrafish, other metabolites were identified and measured by LC-MS. Two minor metabolites, 3-methoxy-paracetamol and n-acetylcysteine-paracetamol, were detected around or below the limit of quantification, resulting in a negligible impact on the mass balance. Low fractions of paracetamol and its metabolites excreted according to the model were recovered in the washout experiment, for which an empirical recovery fraction was estimated. This might be due to adhesion of these compounds to the wells. Moreover, the accuracy of the measurements of the excreted amounts may be negatively impacted by measurements being close to the limit of quantification. A sensitivity analysis excluding these measurements showed similar parameter estimates, suggesting negligible impact of these data on the final model outcome. Limitations in estimating distribution volume and thereby absolute clearance values are overcome by having measurements of both drug concentrations in the blood as well as total amounts in the larvae. However, it is not clear whether this approach will also be applicable to studying the pharmacokinetics of lipophilic drugs, as these drugs may be more difficult to wash off the larvae due to skin adhesion, which complicates measurements of total internal amounts. However, alternative methods for oral dosing using nanoparticles might provide a solution for this limitation. In conclusion, an improved understanding of pharmacokinetics in the vertebrate model organism zebrafish larva is presented. Based on the method developed to sample blood from these larvae, absolute clearance resulting from different metabolic routes, as well as distribution could be quantified for the paradigm compound paracetamol. A comparison of the estimated absolute clearance and distribution volume with those reported for higher vertebrates shows a good correlation. This improved confidence in the translational value of the zebrafish can contribute to the role of this small vertebrate in drug discovery and development. ## Acknowledgments The authors thank Sebastiaan Goulooze for code review of the R and NONMEM scripts, Bjørn Koch for assistance in the germ-free experiments, and Parth Upadhyay for proofreading the manuscript. ## Data Availability Statement The full data set and model file are publicly available through the DDMoRe Repository, Model ID DDMODEL00000300 ([http://repository.ddmore.foundation/model/DDMODEL00000300](http://repository.ddmore.foundation/model/DDMODEL00000300)). ## Authorship Contributions *Participated in research design:* Van Wijk, Krekels, Spaink, Van der Graaf. *Conducted experiments:* Van Wijk, Ordas, Kreling. *Contributed new reagents or analytic tools:* Kantae, Harms, Hankemeier. *Performed data analysis:* Van Wijk, Krekels. *Wrote or contributed to the writing of the manuscript:* Van Wijk, Krekels, Kantae, Ordas, Kreling, Harms, Hankemeier, Spaink, Van der Graaf. ***Note Added in Proof***—A supplemental video of the blood sampling procedure mentioned in Figure 1 was added after the Fast Forward version appeared online August 15, 2019. [http://jpet.aspetjournals.org/content/early/2019/08/15/jpet.119.260299/tab-figures-data](http://jpet.aspetjournals.org/content/early/2019/08/15/jpet.119.260299/tab-figures-data) ## Footnotes * Received May 22, 2019. * Accepted July 31, 2019. * 1 Current affiliation: Mechanistic Biology and Profiling, Discovery Sciences, BioPharmaceuticals R&D, AstraZeneca, Cambridge, United Kingdom. * This work was supported by the Leiden University Fund (A.O. and H.P.S.). * P.H.v.d.G. is an employee of Certara. All authors declare no conflict of interest. * [https://doi.org/10.1124/jpet.119.260299](https://doi.org/10.1124/jpet.119.260299). * ![Graphic][3]This article has supplemental material available at [jpet.aspetjournals.org](http://jpet.aspetjournals.org). ## Abbreviations dpf : days post fertilization LC-MS : liquid chromatography–mass spectrometry LLOQ : lower limit of quantification PAPS : 3′-phosphoadenosine 5′-phosphosulfate RSE : relative standard error * Copyright © 2019 by The Author(s) This is an open access article distributed under the [CC BY-NC Attribution 4.0 International license](http://creativecommons.org/licenses/by-nc/4.0/). ## References 1. Ali A, Abouleila Y, Shimizu Y, Hiyama E, Watanabe TM, Yanagida T, and Germond A (2019) Single-cell screening of tamoxifen abundance and effect using mass spectrometry and Raman-spectroscopy. Anal Chem 91:2710–2718. 2. Anderson BJ, van Lingen RA, Hansen TG, Lin Y-C, and Holford NHG (2002) Acetaminophen developmental pharmacokinetics in premature neonates and infants: a pooled population analysis. Anesthesiology 96:1336–1345. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1097/00000542-200206000-00012&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=12170045&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000175965500011&link_type=ISI) 3. 1. Beal S, 2. Sheiner L, 3. Boeckmann A, and 4. Bauer RJ , editors *NONMEM 7.3.0 Users Guides. (1989-2013)*, ICON Development Solutions, Hanover, MD. 4. Beal SL (2001) Ways to fit a PK model with some data below the quantification limit. J Pharmacokinet Pharmacodyn 28:481–504. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1023/A:1012299115260&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=11768292&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000172744200004&link_type=ISI) 5. Brochot A, Dunne A, Poggesi I, and Vermeulen A (2011) Specifying models with time-dependent pharmacokinetic parameters in NONMEM. PAGE 20 Abstr 2176. 6. Cordero-Maldonado ML, Perathoner S, van der Kolk KJ, Boland R, Heins-Marroquin U, Spaink HP, Meijer AH, Crawford AD, and de Sonneville J (2019) Deep learning image recognition enables efficient genome editing in zebrafish by automated injections. PLoS One 14:e0202377. 7. Diekmann H and Hill A (2013) Zebrafish as a platform for in vivo drug discovery ADMETox in zebrafish. Drug Discov Today Dis Models 10:e31–e35. 8. European Union (2010) Council directive 2010/63/EU on the protection of animals used for scientific purposes. Off J Eur Union L276/33. 9. Fujii T, Matsuda S, Tejedor ML, Esaki T, Sakane I, Mizuno H, Tsuyama N, and Masujima T (2015) Direct metabolomics for plant cells by live single-cell mass spectrometry. Nat Protoc 10:1445–1456. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1038/nprot.2015.084&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=26313480&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 10. Goulooze SC, Völler S, Välitalo PAJ, Calvier EAM, Aarons L, Krekels EHJ, and Knibbe CAJ (2019) The influence of normalization weight in population pharmacokinetic covariate models. Clin Pharmacokinet 58:131–138. 11. Guo Y, Veneman WJ, Spaink HP, and Verbeek FJ (2017) Three-dimensional reconstruction and measurements of zebrafish larvae from high-throughput axial-view *in vivo* imaging. Biomed Opt Express 8:2611–2634. 12. Howe K, Clark MD, Torroja CF, Torrance J, Berthelot C, Muffato M, Collins JE, Humphray S, McLaren K, Matthews L, et al. (2013) The zebrafish reference genome sequence and its relationship to the human genome [published correction appears in *Nature* (2014) 505:248]. Nature 496:498–503. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1038/nature12111&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=23594743&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000317984400040&link_type=ISI) 13. Ingalls B (2012) *Mathematical Modelling in Systems Biology: An Introduction*, The MIT Press, Boston MA. 14. Isogai S, Horiguchi M, and Weinstein BM (2001) The vascular anatomy of the developing zebrafish: an atlas of embryonic and early larval development. Dev Biol 230:278–301. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1006/dbio.2000.9995&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=11161578&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000167132900014&link_type=ISI) 15. Janus K, Grochowina B, Antoszek J, Suszycki S, and Muszczynski Z (2003) The effect of food or water deprivation on paracetamol pharmacokinetics in calves. J Vet Pharmacol Ther 26:291–296. [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=12887612&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 16. Kantae V, Krekels EHJ, Ordas A, González O, Van Wijk RC, Harms AC, Racz PI, van der Graaf PH, Spaink HP, and Hankemeier T (2016) Pharmacokinetic modeling of paracetamol uptake and clearance in zebrafish larvae: expanding the allometric scale in vertebrates with five orders of magnitude. Zebrafish 13:504–510. 17. Karlsson MO, Beal SL, and Sheiner LB (1995) Three new residual error models for population PK/PD analyses. J Pharmacokinet Biopharm 23:651–672. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1007/BF02353466&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=8733951&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 18. Keizer RJ, van Benten M, Beijnen JH, Schellens JH, and Huitema AD (2011) Piraña and PCluster: a modeling environment and cluster infrastructure for NONMEM. Comput Methods Programs Biomed 101:72–79. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1016/j.cmpb.2010.04.018&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=20627442&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000286955500006&link_type=ISI) 19. Krekels EHJ, van Ham S, Allegaert K, de Hoon J, Tibboel D, Danhof M, and Knibbe CAJ (2015) Developmental changes rather than repeated administration drive paracetamol glucuronidation in neonates and infants. Eur J Clin Pharmacol 71:1075–1082. 20. Kühnert A, Vogs C, Altenburger R, and Küster E (2013) The internal concentration of organic substances in fish embryos--a toxicokinetic approach. Environ Toxicol Chem 32:1819–1827. 21. Kühnert A, Vogs C, Seiwert B, Aulhorn S, Altenburger R, Hollert H, Küster E, and Busch W (2017) Biotransformation in the zebrafish embryo -temporal gene transcription changes of cytochrome P450 enzymes and internal exposure dynamics of the AhR binding xenobiotic benz[a]anthracene. Environ Pollut 230:1–11. 22. Lee SH, An JH, Lee HJ, and Jung BH (2012) Evaluation of pharmacokinetic differences of acetaminophen in pseudo germ-free rats. Biopharm Drug Dispos 33:292–303. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1002/bdd.1799&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=22806334&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 23. Lindbom L, Pihlgren P, and Jonsson EN (2005) PsN-Toolkit--a collection of computer intensive statistical methods for non-linear mixed effect modeling using NONMEM [published correction appears in Comput Methods Programs Biomed (2005) 80:277]. Comput Methods Programs Biomed 79:241–257. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1016/j.cmpb.2005.04.005&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=16023764&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000231971600005&link_type=ISI) 24. Liu L and Klaassen CD (1996) Different mechanism of saturation of acetaminophen sulfate conjugation in mice and rats. Toxicol Appl Pharmacol 139:128–134. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1006/taap.1996.0151&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=8685895&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=A1996UW60700016&link_type=ISI) 25. Lowenthal DT, Oie S, Van Stone JC, Briggs WA, and Levy G (1976) Pharmacokinetics of acetaminophen elimination by anephric patients. J Pharmacol Exp Ther 196:570–578. [Abstract/FREE Full Text](http://jpet.aspetjournals.org/lookup/ijlink/YTozOntzOjQ6InBhdGgiO3M6MTQ6Ii9sb29rdXAvaWpsaW5rIjtzOjU6InF1ZXJ5IjthOjQ6e3M6ODoibGlua1R5cGUiO3M6NDoiQUJTVCI7czoxMToiam91cm5hbENvZGUiO3M6NDoianBldCI7czo1OiJyZXNpZCI7czo5OiIxOTYvMy81NzAiO3M6NDoiYXRvbSI7czoxOToiL2pwZXQvMzcxLzEvMTUuYXRvbSI7fXM6ODoiZnJhZ21lbnQiO3M6MDoiIjt9) 26. MacRae CA and Peterson RT (2015) Zebrafish as tools for drug discovery. Nat Rev Drug Discov 14:721–731. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1038/nrd4627&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=26361349&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 27. Morgan P, Van Der Graaf PH, Arrowsmith J, Feltner DE, Drummond KS, Wegner CD, and Street SDA (2012) Can the flow of medicines be improved? Fundamental pharmacokinetic and pharmacological principles toward improving phase II survival. Drug Discov Today 17:419–424. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1016/j.drudis.2011.12.020&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=22227532&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000304297100004&link_type=ISI) 28. Neirinckx E, Vervaet C, De Boever S, Remon JP, Gommeren K, Daminet S, De Backer P, and Croubels S (2010) Species comparison of oral bioavailability, first-pass metabolism and pharmacokinetics of acetaminophen. Res Vet Sci 89:113–119. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1016/j.rvsc.2010.02.002&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=20211479&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 29. Nguyen THT, Mouksassi M-S, Holford N, Al-Huniti N, Freedman I, Hooker AC, John J, Karlsson MO, Mould DR, Pérez Ruixo JJ, et al., and Model Evaluation Group of the International Society of Pharmacometrics (ISoP) Best Practice Committee (2017) Model evaluation of continuous data pharmacometric models: metrics and graphics. CPT Pharmacometrics Syst Pharmacol 6:87–109. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1002/psp4.12161&link_type=DOI) 30. Owens KH, Medlicott NJ, Zacharias M, Curran N, Chary S, Thompson-Fawcett M, and Reith DM (2012) The pharmacokinetic profile of intravenous paracetamol in adult patients undergoing major abdominal surgery. Ther Drug Monit 34:713–721. [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=23149443&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 31. Owens KH, Murphy PG, Medlicott NJ, Kennedy J, Zacharias M, Curran N, Sreebhavan S, Thompson-Fawcett M, and Reith DM (2014) Population pharmacokinetics of intravenous acetaminophen and its metabolites in major surgical patients. J Pharmacokinet Pharmacodyn 41:211–221. 32. R Core Team (2014) *R: A Language and Environment for Statistical Computing*, R Foundation for Statistical Computing, Vienna, Austria. 33. Reddyhoff D, Ward J, Williams D, Regan S, and Webb S (2015) Timescale analysis of a mathematical model of acetaminophen metabolism and toxicity. J Theor Biol 386:132–146. 34. Reith D, Medlicott NJ, Kumara De Silva R, Yang L, Hickling J, and Zacharias M (2009) Simultaneous modelling of the Michaelis-Menten kinetics of paracetamol sulphation and glucuronidation. Clin Exp Pharmacol Physiol 36:35–42. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1111/j.1440-1681.2008.05029.x&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=18759860&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 35. Schindelin J, Arganda-Carreras I, Frise E, Kaynig V, Longair M, Pietzsch T, Preibisch S, Rueden C, Saalfeld S, Schmid B, et al. (2012) Fiji: an open-source platform for biological-image analysis. Nat Methods 9:676–682. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1038/nmeth.2019&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=22743772&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000305942200021&link_type=ISI) 36. Schulthess P, Van Wijk RC, Krekels EHJ, Yates JWT, Spaink HP, and van der Graaf PH (2018) Outside-in systems pharmacology combines innovative computational methods with high-throughput whole vertebrate studies. CPT Pharmacometrics Syst Pharmacol 7:285–287. 37. Slattery JT, Wilson JM, Kalhorn TF, and Nelson SD (1987) Dose-dependent pharmacokinetics of acetaminophen: evidence of glutathione depletion in humans. Clin Pharmacol Ther 41:413–418. [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=3829578&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=A1987G942400010&link_type=ISI) 38. Spaink HP, Cui C, Wiweger MI, Jansen HJ, Veneman WJ, Marín-Juez R, de Sonneville J, Ordas A, Torraca V, van der Ent W, et al. (2013) Robotic injection of zebrafish embryos for high-throughput screening in disease models. Methods 62:246–254. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1016/j.ymeth.2013.06.002&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=23769806&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000324454000009&link_type=ISI) 39. Van den Boom J, Heider D, Martin SR, Pastore A, and Mueller JW (2012) 3′-Phosphoadenosine 5′-phosphosulfate (PAPS) synthases, naturally fragile enzymes specifically stabilized by nucleotide binding. J Biol Chem 287:17645–17655. [Abstract/FREE Full Text](http://jpet.aspetjournals.org/lookup/ijlink/YTozOntzOjQ6InBhdGgiO3M6MTQ6Ii9sb29rdXAvaWpsaW5rIjtzOjU6InF1ZXJ5IjthOjQ6e3M6ODoibGlua1R5cGUiO3M6NDoiQUJTVCI7czoxMToiam91cm5hbENvZGUiO3M6MzoiamJjIjtzOjU6InJlc2lkIjtzOjEyOiIyODcvMjEvMTc2NDUiO3M6NDoiYXRvbSI7czoxOToiL2pwZXQvMzcxLzEvMTUuYXRvbSI7fXM6ODoiZnJhZ21lbnQiO3M6MDoiIjt9) 40. Van Rongen A, Välitalo PAJ, Peeters MYM, Boerma D, Huisman FW, van Ramshorst B, van Dongen EPA, van den Anker JN, and Knibbe CAJ (2016) Morbidly obese patients exhibit increased CYP2E1-mediated oxidation of acetaminophen. Clin Pharmacokinet 55:833–847. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1007/s40262-015-0357-0&link_type=DOI) 41. Van Wijk RC, Krekels EHJ, Hankemeier T, Spaink HP, and Van der Graaf PH (2016) Systems pharmacology of hepatic metabolism in zebrafish larvae. Drug Discov Today Dis Models 22:27–34. 42. Van Wijk RC, Krekels EHJ, Kantae V, Harms AC, Hankemeier T, Van der Graaf PH, and Spaink HP (2019) Impact of post-hatching maturation on the pharmacokinetics of paracetamol in zebrafish larvae. Sci Rep 9:2149. 43. Vogs C, Kühnert A, Hug C, Küster E, and Altenburger R (2015) A toxicokinetic study of specifically acting and reactive organic chemicals for the prediction of internal effect concentrations in Scenedesmus vacuolatus. Environ Toxicol Chem 34:100–111. 44. Westerfield M (2000) *The Zebrafish Book. A Guide for the Laboratory Use of Zebrafish (Danio rerio)*, 4th ed, University of Oregon Press, Eugene, OR. 45. Yamasaki I, Uotsu N, Yamaguchi K, Takayanagi R, and Yamada Y (2011) Effects of kale ingestion on pharmacokinetics of acetaminophen in rats. Biomed Res 32:357–362. [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=22199125&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) 46. Zon LI and Peterson RT (2005) In vivo drug discovery in the zebrafish. Nat Rev Drug Discov 4:35–44. [CrossRef](http://jpet.aspetjournals.org/lookup/external-ref?access_num=10.1038/nrd1606&link_type=DOI) [PubMed](http://jpet.aspetjournals.org/lookup/external-ref?access_num=15688071&link_type=MED&atom=%2Fjpet%2F371%2F1%2F15.atom) [Web of Science](http://jpet.aspetjournals.org/lookup/external-ref?access_num=000226153200016&link_type=ISI) [1]: /embed/graphic-1.gif [2]: /embed/graphic-4.gif [3]: /embed/inline-graphic-2.gif