Population Life-course exposure to health effects model (PLETHEM): An R package for PBPK modeling

Journal Pre-proofs Editorial Population Life-course Exposure to Health Effects Model (PLETHEM): An R package for PBPK modeling Salil N. Pendse, Alina Efremenko, C. Eric Hack, Marjory Moreau, Pankajini Mallick, Michael Dzierlenga, Chantel I. Nicolas, Miyoung Yoon, Harvey J. Clewell, Patrick D. McMullen PII: DOI: Reference:

S2468-1113(19)30065-9 https://doi.org/10.1016/j.comtox.2019.100115 COMTOX 100115

To appear in:

Computational Toxicology

Received Date: Accepted Date:

26 November 2019 9 December 2019

Please cite this article as: S.N. Pendse, A. Efremenko, C. Eric Hack, M. Moreau, P. Mallick, M. Dzierlenga, C.I. Nicolas, M. Yoon, H.J. Clewell, P.D. McMullen, Population Life-course Exposure to Health Effects Model (PLETHEM): An R package for PBPK modeling, Computational Toxicology (2019), doi: https://doi.org/10.1016/ j.comtox.2019.100115

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier B.V.

Population Life-course Exposure to Health Effects Model (PLETHEM): An R package for PBPK modeling Salil N. Pendse1 Alina Efremenko1 C. Eric Hack1 Marjory Moreau1 Pankajini Mallick1 Michael Dzierlenga1 Chantel I. Nicolas1 Miyoung Yoon1,2 Harvey J. Clewell1,3 Patrick D. McMullen1,*

1ScitoVation,

Durham, NC, USA

2ToxStrategies, 3Ramboll,

*To

Cary, NC, USA

Research Triangle Park, NC, USA

whom correspondence should be addressed

Patrick D. McMullen 100 Capitola Drive, Suite 106, Durham, NC, USA [email protected] 847-477-5938

For submission to: Computational Toxicology Running title: The PLETHEM PBPK modeling package Keywords: IVIVE; Monte Carlo; Source to Outcome Continuum; Pharmacokinetics

1 Abstract An outstanding challenge in the acceptance of alternatives to animal testing is the systematic incorporation of computational models into decision making pipelines. Fifteen years ago, the US EPA Office of Research and Development's framework for computational toxicology emphasized the need for computational methods to bridge the source-to-outcome continuum. This can be achieved by linking exposure estimation methods, physiologically based pharmacokinetic (PBPK) modeling, and computational systems biology pathway modeling tools into a standardized framework. To that end, we have developed the Population Life-course Exposure to Health Effects Model (PLETHEM) suite, a modular open source modeling platform that provides users the ability to create, run, share, and audit PBPK models. The platform consists of a database of chemicals, QSAR models, life-stage specific physiological and metabolic parameters needed to parameterize PBPK models, an R-based engine to perform model simulations, and an interactive user interface to define and select parameter sets for the models. PLETHEM implements easy to use interfaces for a generic PBPK model and a high-throughput IVIVE model. These model interfaces along with the included databases provide capabilities necessary for rapid analysis of chemicals using PBPK modeling. PLETHEM includes the ability to run Monte Carlo analyses to investigate population variance and a set of life-stage equations to investigate life-stage-based sensitivities. The PLETHEM database also incorporates ontogeny profiles for key metabolic enzymes that can be used to calculate in vivo metabolic clearance using measured in vitro clearance. In addition, PLETHEM has an ability to link to a number of EPA and OECD exposure estimation programs. These models, which estimate exposures in the workplace and the general populations, can be used to drive PBPK model-based estimates of resulting internal exposures to support risk assessments. PLETHEM is now freely available as an R package through the Bitbucket and GitHub open source repositories.

2 Introduction While estimates vary, the number of distinct chemical entities in global commerce total in the tens of thousands . Natural and anthropogenic chemicals in the environment can interact with human biology—as well as the biology of other species—in a myriad of ways. Because chemicals are widely used in diverse applications (e.g., consumer products, pesticides, and industrial applications) and are encountered in a variety of exposure contexts, the nature of these interactions and the means by which the chemicals enter the body is varied. Understanding the potential impact of these chemicals on biological systems requires the consideration of two factors. First is the qualitative nature of the impact on the biological system. Unlike pharmaceuticals, these environmental compounds are not typically synthesized to act on a specific cellular pathway. In fact, many times the potential mode of action is not understood or is too broad to be fully characterized . However, these interactions at the molecular level in some cases contribute to a disease state, referred to as the chemical’s adverse effect. The second consideration for evaluating potential chemical hazard is the amount of chemical required to cause the adverse effect. This is commonly expressed in terms of external dose: the amount of chemical entering the body via an exposure route (e.g., ingestion, inhalation, dermal absorption). While these metrics are necessary in order to make comparisons with exposure information, there can be significant variation across chemicals, between individuals, and between species in pharmacokinetics: the way the chemical agent is absorbed, distributed throughout the body, metabolized, and ultimately eliminated. Ultimately, the concentration of chemical in a target tissue may be more directly related to potential health effects than the external exposure. Physiologically based pharmacokinetic (PBPK) modeling has been used extensively to establish internal dose metrics in animals and humans given an exposure . PBPK modeling can be used to further ascertain the mode of action and refine the point of departure for a chemical of interest as well as to extrapolate across species and to different exposure routes [6]. Life-stage modeling can be used to identify susceptible subpopulations [8]. Reverse-dosimetry workflows can be used to estimate external exposures corresponding to internal biomarkers [3]. PBPK

modeling also allows for reduced animal use by incorporating in vitro and in silico experimental data into model development . The use of PBPK modeling to overcome dosimetry challenges and inform chemical safety decisions is increasing with the rise of new approach methodologies in the global regulatory frameworks. However, PBPK models are expensive to develop and the majority of PBPK models are typically characterized for a limited domain of application (e.g., specific chemicals and/or specific exposure scenarios). There is a parallel need for a pharmacokinetic modeling framework that could be used to model the balance of the chemical universe for which detailed PBPK models have been developed [10]. Generalized, chemical-agnostic PBPK models could serve as the basis for screening-level estimates of the relationship between external and internal dose, informing the design of downstream experiments and the development of additional models to refine understanding of pharmacokinetics for target compounds. The National Research Council in their 2007 report stated that there is a need to develop more in vitro and in silico approaches to reduce the dependence on animal testing [15]. Simultaneously, the United States Environmental Protection Agency (EPA), developed the concept of the source-to-outcome continuum as a way to elucidate the steps involved in tracking the fate of the compound from environmental release all the way to the adverse outcome in a biological system [16] (Figure 1). They also expressed the need for a tool that can track the progression of a compound under consideration along this continuum. Each of the steps along the continuum account for a different type of uncertainty involved in assessing chemical risk. More importantly, this specifies exposure estimation as integral to assessing risk for chemicals in a high-throughput manner. Exposure estimates can be used in conjunction with an external point of departure for an adverse effect for risk-based decision making. This concept—often encoded as a ratio of hazard and exposure—is powerful in a number of decision-making contexts, such as when trying to prioritize a large list of chemicals for further expensive testing. Exposure across certain sections of the population can be used by regulators to better understand the risk associated with a new chemical. Here we describe a modeling tool—the Population Life-course Exposure to Health Effects

Model (PLETHEM) framework—that integrates various aspects of the source to outcome continuum. At its core, PLETHEM is a general-use PBPK model called rapidPBPK. PLETHEM is bundled as an R package and is freely available via the Comprehensive R Archive Network (PLETHEM CRAN Page) and GitHub (PLETHEM GitHub Page). We have also setup a helpdesk at [email protected] to assist end users with any questions regarding PLETHEM.

Figure 1. The source-to-outcome continuum as described previously [16]. Reproduced with permission.

3 Methods PLETHEM is written as a package for the statistical language R . R is a freely available opensource software that is used extensively in computational toxicology. PLETHEM was designed to reduce technical challenges involved with PBPK modeling by providing an easy-to-use interface to parameterize, save, simulate models, and to view and export the results of those simulations. Broadly, the PLETHEM package can be divided into three distinct parts: the user interface, the database, and the underlying logic. The package contains multiple workflows used to evaluate the toxicity of chemicals. In the following sections, we explain the different workflows within PLETHEM and their utility in risk assessment.

3.1 PBPK Model: The rapidPBPK model (Figure 2) is the default PBPK model implemented in PLETHEM. The model is derived from the IndusChemFate model [11] and has 11 flow- or diffusion-limited

compartments . The standard approach to PBPK modeling encourages users to define the compartments as flow limited and only moving to diffusion limited definitions if needed. By default, the permeability area coefficients between the exchange sub-compartment and the tissue sub-compartment for all the compartments of the rapidPBPK are set to 1000, forcing flow-limited behavior. If required, users can set the permeability area coefficient of the tissue to simulate diffusion-limited distribution. The model also allows for four routes of exposure (oral, inhalation, intravenous, and dermal). Parent chemical clearance in the model can be defined via hepatic metabolism, urinary excretion, fecal excretion, exhalation, and plasma clearance. The list of parameters for the model can be found in the supplementary materials (section S1).

Figure 2. Model schematic for the rapidPBPK model in PLETHEM. The model is based on the IndusChemFate model and has 11 diffusion limited compartments with 3 routes of exposure. Clearance can be modeled as fecal clearance (fraction absorbed in the gut), plasma clearance (first order non-saturable), urinary clearance (first order non-saturable), and hepatic clearance (saturable or linear).

3.1.1 Exposure The PLETHEM rapidPBPK model described above allows users to define exposure through inhalation, oral, intravenous or dermal routes. Only one route of exposure can be active for a single simulation, and dosing through multiple routes of exposure simultaneously is not currently possible. Inhalation exposure scenarios are flexible and can simulate environmental exposure, consumer exposures, or occupational exposures. They allow the user to set continuous exposures or exposures that

are intermittent, as might be appropriate for occupational exposures (e.g., active for a set period per day, repeated for a specified number of days per week). Intravenous (IV) exposure in the model is simulated as a continuous infusion that for a specific length of time. The IV exposure has the units of mg/h exposed for a length of time expressed in h/day which can be setup to repeat daily. Bolus IV infusion can be modeled by setting the length to be 0.01h and multiplying the exposure value by 100. Oral

exposure scenarios can

be

setup

as

repeated bolus dose with exposure in mg/day/kg body weight, or as a continuous exposure through drinking water with dissolved chemical concentration in mg/L. PLETHEM has 2 options for the oral absorption model. The simplest model uses a parameter to Figure 4. Simple Oral Exposure model in PLETHEM. Fa

describe the fraction absorbed from the gut (fa, unitless denotes the fraction of the chemical available for between 0 and 1), with 1 ― 𝑓𝑎 being the fraction not

absorption into the gut tissue from the lumen. Ka defines the rate at which the available fraction is absorbed per hour.

absorbed and excreted in feces (Figure 3). The second oral absorption model can simulate delayed gavage vehicle release and fecal excretion using first-order rate constants, and enterocyte metabolism (Figure 4). For this model, kVtoL (1/h) defines the rate at which the chemical leaves the vehicle, thus making it available for absorption, and kfec (1/h) is the rate at Figure 3. Oral exposure model for simulating oral gavage dissolved in vehicle as used in oral in vivo exposures. PLETHEM currently does not have QSAR models to estimate any of the parameters for this exposure model

which chemical is excreted in feces. Thus, the Fa in the first model is a factor that simulates the combined effect of kfec and kVtoL parameters. In both

absorption models, ka (/h) provides the rate at which the chemical is absorbed from the lumen into portal

circulation and kent (L/h) is used to denote metabolism in the gut as it is absorbed and passes through the enterocytes. Dermal routes of exposure expand the utility of PLETHEM models to occupational exposures and exposures from consumer products. We have adapted the dermal model described previously

[18].

Reference

Error!

source

not

found. shows the schematic for the dermal exposure added

to

the

Figure 5. Dermal exposure model implemented in PLETHEM. The exposure can either be

skin defined as a continuous deposition or bolus dermal dose deposited on the skin surface.

compartment in the model.

QSAR models are used to estimate permeation of the chemical through different layers of the skin.

QSAR models previously described are used to predict diffusion of the chemical through the different layers of the skin [18]. The units for the exposure values for different routes of exposure were selected to reflect the ways in which these are most commonly reported in the literature. The model however simulates distribution in the units of µmoles. PLETHEM uses the molecular weight of the chemical being studies expressed in mol/g to convert mass units to the one required by the model.

3.1.2 Metabolism The model allows for metabolic clearance in the liver and a first order clearance in the blood. The liver metabolic clearance is modeled as Michaelis-Menten kinetics of the form:

𝑉𝑚𝑎𝑥 ∗ 𝐶𝑜𝑛𝑐𝑝𝑎𝑟𝑒𝑛𝑡 𝑑 𝑀𝑒𝑡𝑎𝑏𝑜𝑙𝑖𝑡𝑒 = 𝑑𝑡 𝐾𝑀 + 𝐶𝑜𝑛𝑐𝑝𝑎𝑟𝑒𝑛𝑡 Where Metabolite is the amount of metabolite formed, Vmax is the maximum metabolism rate (µmol/h), KM is the Michaelis-Menten constant for the chemical (µmol/L), and Concparent is the concentration of the parent compound in the liver (µmol/L). The user can enter Vmax directly while setting up the model or use in vitro-to-in vivo extrapolation (IVIVE) to estimate Vmax using measured in vitro intrinsic clearance. Metabolism in the liver can also be modelled as linear non-saturable kinetics. 𝑑 𝑀𝑒𝑡𝑎𝑏𝑜𝑙𝑖𝑡𝑒 = 𝑉𝐾𝑀1 ∗ 𝐶𝑜𝑛𝑐𝑝𝑎𝑟𝑒𝑛𝑡 𝑑𝑡 Where VKM1 is the first order metabolism rate in the liver expressed in L/h. The values for the metabolic parameters (Vmax, KM,. and VKM1) can be defined by the user. If there exists in vitro hepatic intrinsic clearance data for the chemical of interest, then the IVIVE module described below can be used to estimate the hepatic clearance in vivo. Clearance in the blood is modeled as a first order rate dependent on the concentration of the parent in the plasma. 𝑑 𝑀𝑒𝑡𝑎𝑏𝑜𝑙𝑖𝑡𝑒𝑃𝑙𝑎𝑠𝑚𝑎 = 𝐾𝑏𝑙𝑑 ∗ 𝐶𝑜𝑛𝑐𝑝𝑎𝑟𝑒𝑛𝑡 𝑑𝑡 Where 𝐾𝑏𝑙𝑑 is a first order metabolism rate in blood expressed in L/h.

3.1.3 Excretion PLETHEM models Fecal and Urinary routes of excretion. Fecal excretion is modeled as described above. Urinary excretion is modeled as follows: 𝑑 𝐴𝑚𝑜𝑢𝑛𝑡 𝑖𝑛 𝑈𝑟𝑖𝑛𝑒 = 𝐺𝐹𝑅 ∗ 𝐹𝑅𝑊𝑆𝑂𝐿 ∗ (1 ― 𝐹𝑅𝐸𝑆𝑅𝑃) ∗ 𝐶𝑜𝑛𝑐𝑝𝑎𝑟𝑒𝑛𝑡 𝑑𝑡

Where GFR is the Glomerular Filtration Rate in L/h, FRWSOL is the fraction of chemical dissolved in the liquid phase of blood, FRESRP is the fraction of chemical resorbed by the kidney, and Concparent.is the concentration of chemical in the kidney.

3.2 Life-course Equations: Changes in physiology—both macroscopically (e.g., body weight) and microscopically (e.g., protein concentration)—over the life of an organism can have a dramatic impact on the predictions from a PBPK model. Accounting for these age-dependent differences allows appropriate consideration of dosimetry for potentially exposed and susceptible subpopulations (e.g., children, pregnant women, and/or elderly). PLETHEM contains equations that describe the changes in many physiological parameters over the lifespan of an individual (see Supplemental Material S2). These relationships between physiological parameters and age are calibrated via data from the United States Centers for Disease Control and Prevention’s National Health and Nutrition Examination Survey (NHANES) conducted in 2013 [19]. These relationships govern age-dependent parameters by default. The values of these parameters however can be set by the users directly if required.

3.3 Population and Exposure Variability: The model by default simulates an average human for the age and gender selected. Similarly, the life-course equations return the values of physiological parameters for a representative person (defined by the arithmetic mean of each physiological attribute) from the NHANES population. For some applications, it is important to understand the chemical kinetics on the population of interest rather than just the average person from that population to better establish chemical dosimetry. PLETHEM provides the user the ability to define distribution parameters for all population-relevant physiological parameters in the model. One of the important goals of PLETHEM is to provide a unified platform for tracking a chemical along the source-to-outcome continuum. Exposure variability is an important factor that needs to be considered along this continuum. The variability in exposure and physiological parameters is helpful in obtaining a complete picture of population sensitivity for any given chemical. Normal, lognormal, or uniform distributions can be selected to simulate variability for relevant parameters as identified in Supplementary Table 1.

3.4 In Vitro to In Vivo Extrapolation (IVIVE): IVIVE workflows allow the users to estimate hepatic metabolic clearance in vivo from measured in vitro hepatic clearance. Error! Reference source not found. shows the different in vitro systems that are frequently used for estimating hepatic metabolism parameters such as Vmax and KM. Data obtained from any of these systems can be extrapolated using PLETHEM by applying the following extrapolation with values for

the

scaling constants dependen t on the type

of

invitro system used.

𝐶𝑙𝑒𝑎𝑟𝑎𝑛𝑐𝑒𝑖𝑛𝑣𝑖𝑣𝑜 = 𝐶𝑙𝑒𝑎𝑟𝑎𝑛𝑐𝑒𝑖𝑛𝑣𝑖𝑡𝑟𝑜 × 𝑉𝑆𝐹 × 𝑇𝑆𝐹 × 𝐸𝐶𝑆𝐹 × 𝐿𝑖𝑣𝑒𝑟 𝑊𝑒𝑖𝑔ℎ𝑡 Where: Clearanceinvitro is the clearance measured in the invitro system (e.g., µL/min/mg protein) VSF is the factor used to scale the in vitro volume units to liters (e.g., if volume is measured in µL Figure 6. Different in vitro systems used to measure metabolism. Data measured in these systems can be scaled in PLETHEM using IVIVE to estimate the in vivo metabolic clearance. Each system has its own drawbacks and advantages that can make it suitable for certain application. Tissue based systems make use of hepatocytes or isolated microsomes and cytosol from hepatocytes to measure clearance. Recombinant enzymes on the other hand use cell free methods to estimate clearance. Clearance measured from tissue based systems, cell free systems and enzyme ontogeny data can be used to calculate metabolic clearance at different life stages.

then VSF has

a

value of 10-6) TSF is the factor used to scale the in vitro time units to per hour (e.g., if time is measured in minutes then TSF has a value of 60) ECSF is the scaling factor used to scale in vitro enzyme content or hepatocellularity to per gram liver. Liver Weight is the weight of the liver at the life stage for which the ECSF is defined. It is expressed in grams. Table 1 provides the values for ECSF in rats and humans [20-22]. Table 1 Enzyme content scaling factors (ECSF) for different measurement systems.

ECSF In vitro system Microsomes Cystosol

Units 𝑚𝑔 𝑀𝑖𝑐𝑟𝑜𝑠𝑜𝑚𝑎𝑙 𝑃𝑟𝑜𝑡𝑒𝑖𝑛𝑠 𝑚𝑔 𝐶𝑦𝑡𝑜𝑠𝑜𝑙𝑖𝑐 𝑃𝑟𝑜𝑡𝑒𝑖𝑛𝑠

𝑔 𝐿𝑖𝑣𝑒𝑟 𝑤𝑒𝑖𝑔ℎ𝑡 𝑔 𝐿𝑖𝑣𝑒𝑟 𝑤𝑒𝑖𝑔ℎ𝑡

S9 fraction

𝑚𝑔 𝑃𝑟𝑜𝑡𝑒𝑖𝑛𝑠 𝑖𝑛 𝑆9 𝐹𝑟𝑎𝑐𝑡𝑖𝑜𝑛

Hepatocytes

# 𝑜𝑓 𝑚𝑖𝑙𝑙𝑖𝑜𝑛 ℎ𝑒𝑝𝑎𝑡𝑜𝑐𝑦𝑡𝑒𝑠

𝑔 𝐿𝑖𝑣𝑒𝑟 𝑤𝑒𝑖𝑔ℎ𝑡 𝑔 𝐿𝑖𝑣𝑒𝑟 𝑤𝑒𝑖𝑔ℎ𝑡

Human

Rat

40

45

80

91

120

136

99

110

The extrapolated clearance value is then used in the PBPK model. The IVIVE workflow extrapolates intrinsic clearance as both saturable and linear hepatic clearance. The type of metabolism to use is dependent on the chemical to be modeled. Only Vmax and VKM1 values can be extrapolated using PLETHEM. KM values measured in in vitro studies can be used directly. The complement of metabolic enzymes and their expressions changes in all animals during their lifespan. This change goes beyond allometric scaling, the normal increase in capacity for metabolism with

increase in body weight. Different metabolic enzymes reach their peak abundance at different times during early life. These differences can sometimes be crucial to understand the toxicity of chemicals, especially for potentially exposed and susceptible subpopulations. This is especially true when the metabolism of a specific chemical is driven largely by a small number of enzymes. The PLETHEM database contains information on abundance of 9 metabolic enzymes that are commonly implicated in metabolism of environmental chemicals and pharmaceuticals (Table 2). Abundance information for cytochrome P450 enzymes [23] and for carboxylesterase enzymes [24] in adult liver were described previously [23]. The abundance information for carboxylesterase enzymes was obtained from[24]. CES2C abundance was calculated based on the relative abundance of this enzyme, compared to the microsomal CES2M using the expression level ratio of CES2m/CES2c from Hines et al. (2016). Table 2. Abundance information on 9 enzymes in PLETHEM database. ISEF[25-27] is the inter-system extrapolation factor and fu(mic) is the unbound fraction in the media.

Name

Abundance (at age 25

ISEF

fu(mic)

Location

Y) CYP1A2

39

0.43

1

Microsomes

CYP2B6

16

0.43

1

Microsomes

CYP3A4

93

0.24

1

Microsomes

CYP2C19

11

0.25

1

Microsomes

CYP3A5

17

0.24

1

Microsomes

CES1M

1664

1

1

Microsomes

CES1C

556

1

1

Cytosol

CES2M

174

1

1

Microsomes

CES2C

98

1

1

Cytosol

In addition to that, PLETHEM also contains equations for calculating the fractional abundance of these enzymes at different ages relative to a 25-year-old adult (Supplementary Excel file). An additional IVIVE model is included in PLETHEM to scale in vitro clearance measured in recombinant enzyme

systems. This scaled value can be used in conjunction with scaled hepatic intrinsic clearance from cellbased systems and the included enzyme ontogeny to define age-dependent changes in metabolism. These age-dependent metabolism values can then be used for model-specific life stages using the rapidPBPK model.

3.5 High-Throughput In Vitro to In Vivo Extrapolation (HT-IVIVE) As established earlier, in vitro measures of metabolism can be used to estimate in vivo hepatic clearance of compounds. In vitro assays are also useful systems to measure a point of departure (POD) for specific modes of action. A challenge with these in vitro assays, however, lies within the interpretation of the results they generate. Often it is difficult to relate in vitro POD values like AC50 or EC50 to real world exposure values. HT-IVIVE was developed to address this specific challenge [9, 12, 28, 29]. PLETHEM implements HT-IVIVE algorithms described previously [9, 29] in conjunction with IVIVE for scaling hepatic metabolism – as described in the section above – to estimate oral or inhalation equivalent dose for volatile or non-volatile compounds from measured in vitro PODs. The external equivalent dose for each individual POD can then be used to establish hazard. Table 3 details the options available to the user for performing HT-IVIVE within PLETHEM. Table 3. Options and in vitro data types accepted for HT-IVIVE modeling.

Options

Types Oral Exposure Non-Volatile Chemicals, Oral Exposure Volatile

HT-IVIVE Type Chemicals, Inhalation Exposure Volatile Chemicals Intrinsic Clearance

Hepatic Clearance, Renal Clearance, Blood Clearance

Hepatic Clearance Scaling

Rowland Equation, Restrictive Clearance, Non-restrictive Clearance

Hepatic Clearance

Whole Hepatocytes, Sub Cellular Fraction, S9 Fraction, Recombinant

Measuring Systems

Enzymes

3.5.1 Exposure Estimation Tools A central tenet of chemical safety assessment is that human health risk is a function of two factors: the potential of a chemical to impact health at a dose and the probability of human exposure at that dose. As such, estimating exposure assessment and including it in risk characterization is essential. Exposure estimation tools are being developed by the toxicology community for this purpose [30, 31]. PLETHEM aims to unify these tools under a singular platform. We have incorporated multiple exposure tools within this version of the package. Specifically, the current version of the package provides integration to Stochastic Human Exposure and Dose Simulation (SHEDS-HT) , ExpoCast Systematic Empirical Evaluation of Models (SEEM) [33], and ECETOC Targeted Risk Assessment (TRA) [34]. The degree of integration of these tools is primarily determined by two factors: accessibility and complexity. First, integration at the programmatic level requires availability of model code and data, either in a raw format or through an application programming interface. Open source tools are typically available for this. Second, if the exposure estimating tool requires substantial input from the user, it is less suited for integration into PLETHEM due to its complexity. We defined three levels of integration between PLETHEM and third-party exposure estimation tools (Figure 7Figure 7). Level 1 tools are an inherent part of the platform, and the user interacts with them directly through the PLETHEM interface. ExpoCast SEEM is an example of a Level 1 tool because it is open source and implemented in R. Level 2 tools cannot easily directly interface with PLETHEM but produce output files that have a well-defined format and are suitable for importing programmatically. For these, the user interacts with the third-party tool and imports the results using a PLETHEM-side application programming interface. Examples are SHEDS-HT and TRA. Level 3 tools are completely independent of PLETHEM; the user would produce an exposure estimate directly using the third-party tool and enter this value into the PLETHEM interface. An example of an exposure tool that would be Level 3 is the Consumer Exposure Model , which relies on a proprietary database technology and outputs a human-readable format that varies depending on usage.

Figure 7. Different levels of integration of exposure tools within PLETHEM. Level 1 exposure tools are integrated completely within the package. Level 2 exposure tools need to run independently of PLETHEM, though we have developed user interfaces to read in in estimates obtained from these tools. Level 3 exposures tools do not integrate with PLETHEM in any way. The user is expected to run the tools interpedently of PLETHEM and then manually transfer exposure estimates from these tools into PLETHEM projects

3.5.2 High-Throughput Toxicokinetic Package USEPA’s National Center for Computational Toxicology (NCCT) has created general-purpose pharmacokinetics models to aid in the interpretation of bioactivities determined via high-throughput screening efforts and other in vitro approaches, as well as an accompanying open-source package to aid in their use [10, 36]. The computational toxicology workflows and datasets included in that package are complementary to the workflows and models included in PLETHEM. The package however is meant to run using the R console command line interface. This creates a technical barrier for many users attempting to use the package for risk assessment. Pursuant to our Memorandum of Understanding with the USEPA, we have developed a user interface in R to run the HTTK package from within PLETHEM. This will allow the user to run most of the functions within the HTTK package without having to program in R.

3.6 User Interface One of the biggest challenges in utilizing PBPK modeling and its workflows on a large scale is a high degree of technical knowledge needed to implement these workflows. This prevents more biologically-minded scientists from effectively implementing and exercising these models. One of the

main goals of PLETHEM was to provide users the ability to perform computational modeling workflows without placing the burden of programming on them. We achieved this by creating user interfaces for each of the workflows mentioned above. Figure 5 shows a screenshot of the interface for the rapidPBPK model.

Figure 8. Screenshot for the UI used to define physiological parameter sets in PLETHEM. In this example the user has decided to not include Skin, Muscle, and Bone compartments in the model. The volume and blood flows for these compartments are added to the slowly perfused tissue compartment and rapidly perfused tissue compartment as appropriate. The partitions can be calculated by selecting the QSAR model to use and the chemicals for which they need to be calculated.

The user interfaces allow users to select workflows, parametrize models, simulate models, and compare the results using elements everyone is familiar with, such as text boxes, sliders, buttons, and interactive charts. Additionally, all these workflows can be invoked using RStudio Addins (Figure 9). RStudio is available as an open source desktop graphical integrated development environment for working with R (https://www.rstudio.com/).

Figure 9. RStudio Addins for PLETHEM. All the workflows within PLETHEM can be invoked by selecting the appropriate Addin. This removes the need for the user to use the R console in any way to run commands. The same functions can still be run from the command lime in case the user does not have access to the RStudio IDE for running PLETHEM.

3.7 PLETHEM Databases The second biggest challenge preventing widespread use of these modeling workflows is data. Datasets needed to parameterize the model are not always readily available. We have addressed that issue by providing all the datasets and functions needed to parametrize PBPK models within the package itself. PLETHEM database contains default parameter set for adult human, rat and mouse. Chemical parameters from the set of chemicals found in the original IndusChemFate paper are also available [11]. These parameters can be used alongside quantitative structure-activity relationship (QSAR) models found in PLETHEM to estimate partition coefficients for the model. PLETHEM consists of three distinct internal SQLite databases. The main database is static and contains all the information needed to parameterize the models that are implemented within PLETHEM. The table “ParamNames” contains data for all the parameters that are a part of any model in PLETHEM. This table also identifies if variability can be defined for the given parameter. The main database also stores information related to the current project the user is running. This allows PLETHEM to correctly update the PLETHEM project file when needed. Additionally, this

database also stores the location of the user database. This value is overwritten every time the user sets a new user database. The “user database” allows ad hoc creation of parameter sets that are not a part of PLETHEM by default. These parameter sets can be imported into PLETHEM projects just like a parameter set from the main database. The “project database” is a temporary database that is created when the user is doing modeling in PLETHEM. It has the same structure as the project .Rdata file.

3.8 Modeling Projects in PLETHEM PLETHEM uses a project management system that allows different users to share their modeling projects easily. This ensures reproducibility of results obtained using PLETHEM. When a user is creating different models for any given project, the data associated with those models is stored in the project database described above (Section 1.7) on the user’s machine. When they quit the project, the entire temporary database is converted to a .Rdata file. If the same user wants to open the project again, they can load the project file using a RStudio Addin. Similarly, they can send the file to other users. The recipient can then load the project on their local machine in a similar manner. All the datasets created by the first user will still be available for them to test or expand upon. These files can also be distributed by the author of any publication that uses PLETHEM for modeling.

4 Results To demonstrate the IVIVE and QSAR components of PLETHEM, we performed test simulations in PLETHEM using three compounds with different physicochemical properties. The model structure was kept unchanged for all the simulations and built-in QSAR models were used to estimate partition coefficients. PLETHEM currently has no tools to estimate the free (unbound) fraction of chemical in plasma (fupls), fraction absorbed in the gut (fa), and the rate of absorption from the gut per hour (ka). The values for these parameters were kept the same for all chemicals and were set to the values given in Table

4. We selected chemicals with in vitro clearance values measured using human hepatocytes or subcellular fractions, and with time course plasma concentration data following oral or inhalation exposure (Table 5). Table 4. Absorption parameters used in forward dosimetry example.

Parameter

Value

Units

fupls

1

None

Fa

1

None

Ka

5

h-1

Table 5. In vitro metabolism data used for performing IVIVE. The clearance was measured in different in vitro systems and different organisms.

Compound

Species

Parameter

Value

Assay

Source

Trichloroethylene

Human

Lipscomb et al, 1998 [38]

Human

microsomes

Clewell et al, 1997 [39]

Coumarin

Human

28.9 34.45 2.89 13 171.29 0.67

microsomes

ATRA

Vmax KM Vmax KM Vmax KM

microsomes

Draper et al, 1997

Table 6. Corresponding in-life kinetic data for three case compounds. The data used is for the same organism for which in vitro metabolism values were measured (see Table 5).

Compound

Species

Endpoint

Source

Trichloroethylene

Human

[Plasma] timecourse

Muller et al, 1974 [41]

ATRA

Human

[Plasma] timecourse

Clewell et al, 1997 [39]

Coumarin

Human

[Plasma] timecourse

Ritschel et al, 1979 [42]

The results from the simulation and the corresponding experimental data are as shown in Figure 10.

Figure 10. Comparing the simulation results obtained from PLETHEM to measured in vivo kinetics. These simulations were run using the same model structure for all chemicals. PLETHEM predictions agree with most simulation results. These results can be improved by entering chemical specific values or curve fitting for fraction absorbed in the gut, fraction unbound in plasma, and the rate of absorption in the gut.

Note that the purpose of this example is to demonstrate the ability to combine in vitro clearance, physiology database, and measured in vivo PK data in PLETHEM to generate a simulation of in vivo data. The aim was not to “fit” the experimental data, or evaluate the accuracy of the simulations, though this could certainly be done by fitting parameters or replacing Figure 11. Running the ATRA simulation using chemical specific parameters for absorption improves the predictions. None of the data associated with the

defaults

with

chemical-specific compartments have been changed between the two simulations.

information. Figure 11 shows the results from running the ATRA simulation with both default parameters and chemical specific parameters for fupls, Ka and Fa (Table 7). No other parameters were changed for the new prediction. The model predictions are markedly improved by using chemical specific values for these parameters. Table 7. Chemical specific absorption parameters for the ATRA model

Parameter

Value

Units

Source

fupls

0.05

None

Muindi et al 1992b [43]

Fa

0.48

None

Jing et al 2017 [44]

Ka

0.9

h-1

Jing et al 2017 [44]

5 Discussion The chemical risk assessment community is collectively moving toward a paradigm based on human biology rather than based on testing in live rodent models. New approach methodologies (NAMs) including in vitro assays and computational models have started gaining ground as a viable alternative to in vivo testing of chemical toxicity. Similarly, there is an understanding that potential for adverse health effects should be investigated early in the development pipeline for new chemicals. PBPK modeling is uniquely positioned to play a major role in this transition as it can bridge the gap between in vivo studies and data obtained from these new techniques. Traditional PBPK modeling is very exact and does not lend itself readily to rapid prioritization of existing chemicals or lead selection for new chemicals. PLETHEM was designed as a tool that could fill this role. The PBPK model within PLETHEM enables users to interpret information obtained from NAMs within an in vivo context by removing the significant challenges involved in parameterizing a PBPK modeling. It does this through the QSAR models, life course equations and workflows for performing IVIVE for metabolism. Most stakeholders trying to answer one of the two questions raised above are not technical experts in PBPK modeling and would like to arrive at an understanding of their set of chemicals with minimal effort spent on the actual modeling. The interfaces within the tool allow users with limited technical knowledge of PBPK modeling to still use it for making safety and prioritization decisions. Because the goal of PLETHEM is to provide an easy-to-use generic PBPK model, it does not lend itself readily to performing chemical specific safety evaluations. These evaluations are performed when large amount of information regarding the chemical class and its mode of action is already known. Specifically, the lack of pharmaco-dynamic compartments in the model and full body PBPK models for metabolites limit its utility for traditional PBPK modeling applications like reducing uncertainty factors.

However, the use of interfaces for performing PBPK modeling is highly desired even in those applications and we have successfully used customized versions of PLETHEM to create models used in support of risk assessment [45]. Multiple workflows within PLETHEM can be used for prioritizing a large set of compounds using any of the built-in workflows. The datasets, QSAR models, and IVIVE workflows built into PLETHEM can be used to simulate chemical kinetics using chemical properties and metabolism measured in vitro. An additional roadblock preventing widespread adoption of these in silico techniques for risk assessment is the relatively high technical expertise required to run these models. While there is certainly a large community of modelers versed in the art of PBPK, there is a disparate set of stakeholders involved in generating, interpreting, and making policy decisions based on data on PBPK models. PLETHEM lowers the barrier of entry for interacting with PBPK models, building confidence among stakeholders in this increasingly relevant technology.

6 References 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11.

12. 13. 14. 15. 16. 17. 18. 19.

Schmidt, C.W., TSCA 2.0: A New Era in Chemical Risk Management. Environ Health Perspect, 2016. 124(10): p. A182-A186. Fay, K.A., et al., Differentiating Pathway-Specific From Nonspecific Effects in HighThroughput Toxicity Data: A Foundation for Prioritizing Adverse Outcome Pathway Development. Toxicol Sci, 2018. 163(2): p. 500-515. Clewell, H.J., et al., Quantitative interpretation of human biomonitoring data. Toxicol Appl Pharmacol, 2008. 231(1): p. 122-33. Clewell, R.A. and H.J. Clewell, 3rd, Development and specification of physiologically based pharmacokinetic models for use in risk assessment. Regul Toxicol Pharmacol, 2008. 50(1): p. 129-43. Lipscomb, J.C., et al., Physiologically-based pharmacokinetic (PBPK) models in toxicity testing and risk assessment. Adv Exp Med Biol, 2012. 745: p. 76-95. Clewell, H.J., et al., Comparison of cancer risk estimates for vinyl chloride using animal and human data with a PBPK model. Sci Total Environ, 2001. 274(1-3): p. 37-66. Clewell, H.J., et al., Evaluation of the potential impact of age- and gender-specific pharmacokinetic differences on tissue dosimetry. Toxicol Sci, 2004. 79(2): p. 381-93. Yoon, M. and H.J. Clewell, 3rd, Addressing Early Life Sensitivity Using Physiologically Based Pharmacokinetic Modeling and In Vitro to In Vivo Extrapolation. Toxicol Res, 2016. 32(1): p. 15-20. Yoon, M., et al., Quantitative in vitro to in vivo extrapolation of cell-based toxicity assay results. Crit Rev Toxicol, 2012. 42(8): p. 633-52. Wambaugh, J.F., et al., Toxicokinetic Triage for Environmental Chemicals. Toxicol Sci, 2015. 147(1): p. 55-67. Jongeneelen, F.J. and W.F. Berge, A generic, cross-chemical predictive PBTK model with multiple entry routes running as application in MS Excel; design of the model and comparison of predictions with experimental results. Ann Occup Hyg, 2011. 55(8): p. 841-64. Wetmore, B.A., et al., Incorporating High-Throughput Exposure Predictions With Dosimetry-Adjusted In Vitro Bioactivity to Inform Chemical Toxicity Testing. Toxicol Sci, 2015. 148(1): p. 121-36. Niederalt, C., et al., A generic whole body physiologically based pharmacokinetic model for therapeutic proteins in PK-Sim. J Pharmacokinet Pharmacodyn, 2018. 45(2): p. 235257. Loizou, G. and A. Hogg, MEGen: A Physiologically Based Pharmacokinetic Model Generator. Front Pharmacol, 2011. 2: p. 56. National Research Council, Toxicity Testing in the 21st Century: A Vision and a Strategy. 2007, The National Academies Press: Washington, DC. Kavlock, R., et al., Computational toxicology: framework, partnerships, and program development. September 29-30, 2003, Research Triangle Park, North Carolina. Reprod Toxicol, 2005. 19(3): p. 265-80. R Core Team, R: A Language and Environment for Statistical Computing. 2018. ten Berge, W., A simple dermal absorption model: derivation and application. Chemosphere, 2009. 75(11): p. 1440-5. Centers for Disease Control and Prevention. National Health and Nutrition Examination Survey Data. 2014 [cited 2017 May 8]; Available from: https://wwwn.cdc.gov/nchs/nhanes/ContinuousNhanes/Default.aspx?BeginYear=2013.

20. 21. 22. 23. 24.

25. 26.

27. 28. 29. 30. 31. 32. 33. 34. 35. 36. 37. 38. 39. 40.

Yoon, M., M.C. Madden, and H.A. Barton, Developmental expression of aldehyde dehydrogenase in rat: a comparison of liver and lung development. Toxicol Sci, 2006. 89(2): p. 386-98. Houston, J.B. and A. Galetin, Methods for predicting in vivo pharmacokinetics using data from in vitro assays. Curr Drug Metab, 2008. 9(9): p. 940-51. Barter, Z.E., et al., Covariation of human microsomal protein per gram of liver with age: absence of influence of operator and sample storage may justify interlaboratory data pooling. Drug Metab Dispos, 2008. 36(12): p. 2405-9. Achour, B., J. Barber, and A. Rostami-Hodjegan, Expression of hepatic drugmetabolizing cytochrome p450 enzymes and their intercorrelations: a meta-analysis. Drug Metab Dispos, 2014. 42(8): p. 1349-56. Boberg, M., et al., Age-Dependent Absolute Abundance of Hepatic Carboxylesterases (CES1 and CES2) by LC-MS/MS Proteomics: Application to PBPK Modeling of Oseltamivir In Vivo Pharmacokinetics in Infants. Drug Metab Dispos, 2017. 45(2): p. 216-223. Proctor, N.J., G.T. Tucker, and A. Rostami-Hodjegan, Predicting drug clearance from recombinantly expressed CYPs: intersystem extrapolation factors. Xenobiotica, 2004. 34(2): p. 151-78. Crewe, H.K., et al., Are there differences in the catalytic activity per unit enzyme of recombinantly expressed and human liver microsomal cytochrome P450 2C9? A systematic investigation into inter-system extrapolation factors. Biopharm Drug Dispos, 2011. 32(6): p. 303-18. Wetmore, B.A., et al., Incorporating population variability and susceptible subpopulations into dosimetry for high-throughput toxicity testing. Toxicol Sci, 2014. 142(1): p. 210-24. Wetmore, B.A., et al., Relative impact of incorporating pharmacokinetics on predicting in vivo hazard and mode of action from high-throughput in vitro toxicity assays. Toxicol Sci, 2013. 132(2): p. 327-46. Yoon, M., et al., Evaluation of simple in vitro to in vivo extrapolation approaches for environmental compounds. Toxicol In Vitro, 2014. 28(2): p. 164-70. European Centre for Ecotoxicology and Toxicology of Chemicals, heatDB. European Centre for Ecotoxicology and Toxicology of Chemicals, Guidance for Effective Use of Human Exposure Data in Risk Assessment of Chemicals. 2016. Isaacs, K.K., et al., SHEDS-HT: an integrated probabilistic exposure model for prioritizing exposures to chemicals with near-field and dietary sources. Environ Sci Technol, 2014. 48(21): p. 12750-9. Wambaugh, J.F., et al., High-throughput models for exposure-based chemical prioritization in the ExpoCast project. Environ Sci Technol, 2013. 47(15): p. 8479-88. European Centre for Ecotoxicology and Toxicology of Chemicals, ECETOC TRA version 3: Background and Rationale for the Improvements. 2012. USEPA Office of Pollution Prevention and Toxics, Consumer Exposure Model. 2017. Pearce, R.G., et al., httk: R Package for High-Throughput Toxicokinetics. J Stat Softw, 2017. 79(4): p. 1-26. RStudio Team, RStudio: Integrated Development for R. 2016 Boston, MA: RStudio, Inc. Lipscomb, J.C., et al., In vitro to in vivo extrapolation for trichloroethylene metabolism in humans. Toxicol Appl Pharmacol, 1998. 152(2): p. 376-87. Clewell, H.J., 3rd, et al., A physiologically based pharmacokinetic model for retinoic acid and its metabolites. J Am Acad Dermatol, 1997. 36(3 Pt 2): p. S77-85. Draper, A.J., A. Madan, and A. Parkinson, Inhibition of coumarin 7-hydroxylase activity in human liver microsomes. Arch Biochem Biophys, 1997. 341(1): p. 47-61.

41. 42. 43.

44. 45.

7

Muller, G., M. Spassovski, and D. Henschler, Metabolism of trichloroethylene in man. II. Pharmacokinetics of metabolites. Arch Toxicol, 1974. 32(4): p. 283-95. Ritschel, W.A., M.E. Brady, and H.S. Tan, First-pass effect of coumarin in man. Int J Clin Pharmacol Biopharm, 1979. 17(3): p. 99-103. Muindi, J., et al., Continuous treatment with all-trans retinoic acid causes a progressive reduction in plasma drug concentrations: implications for relapse and retinoid "resistance" in patients with acute promyelocytic leukemia. Blood, 1992. 79(2): p. 299303. Jing, J., et al., Physiologically Based Pharmacokinetic Model of All-trans-Retinoic Acid with Application to Cancer Populations and Drug Interactions. J Pharmacol Exp Ther, 2017. 361(2): p. 246-258. Mallick, P., et al., Development and Application of a Life-Stage Physiologically-Based Pharmacokinetic (PBPK) Model to the Assessment of Internal Dose of Pyrethroids in Humans. Toxicol Sci, 2019.

8 Supplementary Material S1

PBPK model

Supplemental Table 1. Table of parameters for the PBPK model. The package only allows users to define distributions for parameters which have a “Yes” in the variability column during a Monte-Carlo (MC) simulation

Name

Units

ParamSet

Can define variability

Daily oral dose Total length of dosing Number of doses Concentration in drinking water Number of drinking water doses per day Inhalation dose Length of inhalation dose Days of dosing in a week Intravenous Infusion Length of intravenous dose Volume of drinking water Total Bolus Dose Replicates Inhaled Concentration Dermal deposition rate Exposed Skin Area Daily oral dose Total length of dosing Number of doses Total Bolus Dose Replicates Molecular Weight Water Solubility Fraction resorbed in urine Fraction unbound in plasma Fraction Dissolved in Water Phase of Blood Bone tissue to total bone volume ratio Body Weight Cardiac Output Hematocrit Factor Brain tissue to total brain volume ratio Fat tissue to total fat volume ratio Respiration rate GI tissue to total GI

mg/kg BW h/day None mg/L

Exposure Exposure Exposure Exposure

Yes No No Yes

None

Exposure

No

ppm h days mg/h h/day

Exposure Exposure Exposure Exposure Exposure

Yes No No Yes No

L None

Exposure Exposure

Yes No

um/L mg/cm2/h cm2 mg/kg BW h/day None None

Exposure Exposure Exposure Exposure Exposure Exposure Exposure

No Yes Yes Yes No No No

g/mol mg/L None None

Chemical Chemical Chemical Chemical

No No Yes Yes

None

Chemical

Yes

None

Physiological

No

kg L/h None None

Physiological Physiological Physiological Physiological

Yes Yes Yes No

None

Physiological

No

L/h None

Physiological Physiological

Yes No

volume ratio Dead space Urinary flow rate Glomerular filtration rate Heart tissue to total heart volume ratio Kidney tissue to total kidney volume ratio Fraction fat perfusion Liver tissue to total liver volume ratio Lung tissue to total lung volume ratio Fraction skin perfusion Skin partition coefficient Muscle tissue to total muscle volume ratio Fractional perfused tissue volume Fraction muscle perfusion Muscle partition coefficient Rapidly perfused tissue to total rapidly perfused tissue volume ratio Skin tissue to total skin volume ratio Fraction bone perfusion Slowly perfused tissue to total slowly perfused tissue volume ratio Tidal volume Fraction brain perfusion Fractional Blood Compartment Volume Fractional bone volume Fraction lung perfusion Fractional brain volume Fraction heart perfusion Fractional fat volume Fractional GI volume Fraction GI perfusion Fractional heart volume Fractional kidney volume Fractional arterial liver perfusion Fractional venous liver perfusion Fractional liver volume

None L/kg/day None None

Physiological Physiological Physiological Physiological

Yes Yes Yes No

None

Physiological

No

None None

Physiological Physiological

Yes No

None

Physiological

No

None None None

Physiological Physiological Physiological

Yes Yes No

None

Physiological

Yes

None None

Physiological Physiological

Yes Yes

None

Physiological

No

None

Physiological

No

None None

Physiological Physiological

Yes No

L None None

Physiological Physiological Physiological

Yes Yes Yes

None None None None None None None None None None

Physiological Physiological Physiological Physiological Physiological Physiological Physiological Physiological Physiological Physiological

Yes Yes Yes Yes Yes Yes Yes Yes Yes Yes

None

Physiological

Yes

None

Physiological

Yes

Fractional lung volume Fraction kidney perfusion Fractional muscle volume Fractional rapidly perfused tissue volume Fraction rapidly perfused tissue perfusion Fractional skin volume Fractional slowly perfused tissue volume Fraction slowly perfused tissue perfusion First order metabolism in blood Permeability area coefficient for fat tissue Fat partition coefficient Permeability area coefficient for skin tissue Permeability area coefficient for muscle tissue Permeability area coefficient for bone tissue Bone partition coefficient Permeability area coefficient for brain tissue Brain partition coefficient Permeability area coefficient for lung tissue Lung partition coefficient Permeability area coefficient for heart tissue Heart partition coefficient Permeability area coefficient for GI tissue GI partition coefficient Fraction absorbed in gut lumen Rate of absorption in gut lumen Permeability area coefficient for liver tissue Liver partition coefficient Permeability area coefficient for kidney tissue Kidney partition coefficient Permeability area coefficient for rapidly perfused tissue

None None None None

Physiological Physiological Physiological Physiological

Yes Yes Yes Yes

None

Physiological

Yes

None None

Physiological Physiological

Yes Yes

None

Physiological

Yes

L/h

ADME

Yes

None

ADME

Yes

None None

ADME ADME

Yes Yes

None

ADME

Yes

None

ADME

Yes

None None

ADME ADME

Yes Yes

None None

ADME ADME

Yes Yes

None None

ADME ADME

Yes Yes

None None

ADME ADME

Yes Yes

None None

ADME ADME

Yes Yes

None

ADME

Yes

None

ADME

Yes

None None

ADME ADME

Yes Yes

None

ADME

Yes

None

ADME

Yes

Rapidly perfused tissue partition coefficient Permeability area coefficient for slowly perfused tissue Slowly perfused tissue partition coefficient Plasma-Air Partition Coefficient Michaelis Menten Constant Scaled saturable metabolism rate Scaled linear metabolism Total Stratum Corneum permeation coefficient Evaporation rate from Stratum Corneum Maximum capacity of Stratum Corneum Transfer Rate from vehicle to gut lumen Rate of metabolism in the gut lumen Transfer rate to elimination via fecal excretion Total Days Number of Monte-Carlo Runs Duration

S2

None

ADME

Yes

None

ADME

Yes

None

ADME

Yes

None

ADME

Yes

µMolar

ADME

Yes

µmol/h

ADME

No

L/h

ADME

No

cm/h

ADME

No

cm/h

ADME

No

mg/cm2

ADME

No

/h

ADME

Yes

/h

ADME

Yes

/h

ADME

Yes

days none

Simulation Simulation

No No

h

Simulation

No

Life-course equations

Life course equations encoded in PLETHEM are described below. These relationships were first described across multiple publications [1-3]that are referenced as appropriate in the following pages.

S2.1

Body Weight (BW)

Body weight data from the 2015 NHANES survey was fitted using least squares to obtain an age dependent equation . S.2.1.1

Females

𝑩𝑾

{

𝟒𝟐.𝟖𝟓𝟒 𝟐𝟐 ∗ 𝒂𝒈𝒆 + 𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟐𝟔 𝟑 + 𝒂𝒈𝒆 𝟏 + 𝟏𝟒𝟐.𝟏𝟐 ∗ 𝒆𝟎.𝟒𝟒𝟎𝟓 ∗ 𝒂𝒈𝒆 = ( 𝟔.𝟏𝟕 ∗ 𝟏𝟎 ―𝟗 ∗ (𝒂𝒈𝒆 ∗ 𝟏𝟐)𝟑 ― 𝟗.𝟔𝟖 ∗ 𝟏𝟎 ―𝟓 ∗ (𝒂𝒈𝒆 ∗ 𝟏𝟐)𝟐 + 𝟎.𝟏𝟎𝟗 ∗ (𝒂𝒈𝒆 ∗ 𝟏𝟐) + 𝟒𝟑.𝟖𝟐𝟒) ∗ 𝟒𝟐.𝟖𝟓𝟒 𝟑.𝟒 +

𝟔𝟔.𝟏𝟕

𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟐𝟔

S.2.1.2

Males

𝑩𝑾

{

𝟐𝟖 ∗ 𝒂𝒈𝒆 𝟓𝟔 + 𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟐𝟔 𝟔 + 𝒂𝒈𝒆 𝟏 + 𝟕𝟎 ∗ 𝒆 ―𝟎.𝟑𝟏𝟕𝟑 ∗ 𝒂𝒈𝒆 = ( 𝟐.𝟏𝟗 ∗ 𝟏𝟎 ―𝟖 ∗ (𝒂𝒈𝒆 ∗ 𝟏𝟐)𝟑 ― 𝟏.𝟐𝟐 ∗ 𝟏𝟎 ―𝟒 ∗ (𝒂𝒈𝒆 ∗ 𝟏𝟐)𝟐 + 𝟎.𝟏𝟐𝟎 ∗ (𝒂𝒈𝒆 ∗ 𝟏𝟐) + 𝟓𝟎.𝟖𝟐𝟒) ∗ 𝟖𝟏.𝟐𝟎𝟖

S2.2

𝟒+

𝟕𝟔.𝟔

𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟐𝟔

Body Height (BH)

The equations defining body height were fitted to data from the 2015 NHANES survey . S.2.2.1

Females

𝑩𝑯 =

S.2.2.2

{

(

𝟏𝟔𝟐.𝟏𝟎𝟕 ∗ 𝟏 ―

𝟏𝟔𝟐.𝟏𝟓 ― 𝟐 ∗

(

𝒂𝒈𝒆 + 𝟎.𝟕𝟓 𝟏+ 𝟐.𝟑 𝟏𝟏.𝟏𝟓

) ) 𝟎.𝟕

𝒇𝒐𝒓 𝒂𝒈𝒆 ≤ 𝟑

𝒆𝟎.𝟏𝟑𝟓 ∗ (𝒂𝒈𝒆 ― 𝟏𝟏.𝟐𝟓𝟑𝟔) + 𝒆𝟏.𝟐𝟕 ∗ (𝒂𝒈𝒆 ― 𝟏𝟏.𝟐𝟓𝟑𝟔)

)

𝒇𝒐𝒓 𝒂𝒈𝒆 > 𝟑

Males

𝑩𝑯 =

S2.3

(

𝟏

{

(

𝟏𝟖𝟐.𝟎𝟗𝟕𝟒 ∗ 𝟏 ―

𝟏𝟕𝟕.𝟎 ― 𝟐 ∗

((

𝟏

(

𝒂𝒈𝒆 + 𝟎.𝟕𝟓 𝟏+ 𝟎.𝟓 𝟖.𝟐

) ) 𝟎.𝟗

𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟑

)

𝒆𝟎.𝟏𝟑𝟓 ∗ (𝒂𝒈𝒆 ― 𝟏𝟒.𝟖𝟓) + 𝒆𝟏.𝟓 ∗ (𝒂𝒈𝒆 ― 𝟏𝟒.𝟖𝟓))

Body mass index (BMI)

BMI is described by the following equation adapted from Keys et al [4] :

𝑩𝑴𝑰 =

𝑩𝑾 𝑩𝑯𝟐

𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟑

S2.4

Body surface area (BSA)

BSA is described by the following equation adapted from Gehan and George [5]:

𝑩𝑺𝑨 = 𝐞𝐱𝐩 ( ―𝟑.𝟕𝟓 + 𝟎.𝟒𝟐 ∗ 𝐥𝐧 (𝑩𝑯) + 𝟎.𝟓𝟐 ∗ 𝐥𝐧 (𝑩𝑾)) Where 𝑩𝑺𝑨 𝒊𝒔 𝑩𝒐𝒅𝒚 𝑺𝒖𝒓𝒇𝒂𝒄𝒆 𝑨𝒓𝒆𝒂 𝒊𝒏 𝒎𝟐 𝑩𝑯 𝒊𝒔 𝒃𝒐𝒅𝒚 𝒉𝒆𝒊𝒈𝒉𝒕 𝒊𝒏 𝒎𝒆𝒕𝒆𝒓𝒔 𝑩𝑾𝒊𝒔 𝒃𝒐𝒅𝒚 𝒘𝒆𝒊𝒈𝒉𝒕 𝒊𝒏 𝒌𝒊𝒍𝒐𝒈𝒓𝒂𝒎𝒔

S2.5

Hematocrit Factor

Hematocrit factor at different ages was calculated using the following equations [6]:

{

𝟎.𝟑𝟓𝟗 𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟐 𝒉𝒄𝒕 = (𝟏.𝟏𝟐𝟖𝟏𝟓𝟖 ∗ 𝟏𝟎 ―𝟔 ∗ 𝒂𝒈𝒆𝟑) ― (𝟏.𝟕𝟐𝟑𝟔𝟐 ∗ 𝟏𝟎 ―𝟒 ∗ 𝒂𝒈𝒆𝟐) + (𝟖.𝟏𝟓𝟐𝟔𝟒 ∗ 𝟏𝟎 ―𝟑 ∗ 𝒂𝒈𝒆) + 𝟎.𝟑𝟐𝟕𝟑𝟔𝟑 𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟐

S2.6

Fat volume

The percentage of body fat has been linked to body measurements such as BMI in many studies. NHANES collected data on body measurements including fat volume BW, BH, and BMI in the general population. We obtained these data from the website of the Centers for Disease Control and Prevention (CDC, 2014). Models for predicting fat volume were fitted to the NHANES data by performing multiple linear regression with age and BMI. The relationship between the fraction of body fat, age, and BMI was best described by the following equation :

8.1.1.1

Females (𝟏.𝟓𝟑𝟑𝟒 ∗ 𝒆 ―𝟎.𝟏𝟎𝟑 ∗ 𝒂𝒈𝒆 + 𝟎.𝟔𝟕) ∗ 𝑩𝑴𝑰 + 𝟎.𝟔𝟐𝟕𝟔 ∗ 𝒂𝒈𝒆 + 𝟏.𝟎𝟑𝟎𝟏 𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟐𝟓 𝒗𝒇𝒂𝒕𝒄 = 𝟏.𝟗𝟐𝟐𝟒 ∗ 𝑩𝑴𝑰 ― 𝟎.𝟎𝟏𝟖𝟓𝟏𝟕 ∗ 𝑩𝑴𝑰𝟐 + 𝟎.𝟎𝟓𝟓𝟑𝟕 ∗ 𝒂𝒈𝒆 ― 𝟎.𝟕𝟗𝟒𝟖𝟗𝟒 𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟐𝟓

{

8.1.1.2

Males

𝒗𝒇𝒂𝒕𝒄 =

(𝟐.𝟖𝟗𝟕𝟓 ∗ 𝒆 ―𝟎.𝟏𝟐𝟗 ∗ 𝒂𝒈𝒆 + 𝟎.𝟔𝟕) ∗ 𝑩𝑴𝑰 + 𝟎.𝟐𝟔𝟑𝟓 ∗ 𝒂𝒈𝒆 ― 𝟒.𝟖𝟒𝟑 𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟐𝟎

{ ―𝟓.𝟑𝟑𝟕𝟗𝟖 ∗ 𝑩𝑴𝑰 + 𝟎.𝟏𝟏𝟏𝟒𝟗 ∗ 𝑩𝑴𝑰 + 𝟎.𝟎𝟗𝟕𝟗𝟓 ∗ 𝒂𝒈𝒆 + 𝟖𝟓.𝟐𝟒𝟓𝟐𝟏 𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟐𝟎 𝟐

The fraction of body fat 𝒗𝒇𝒂𝒕𝒄 defined above is scaled by the body weight (𝑩𝑾) to get the volume of fat in the body, which is expressed as 𝒗𝒇𝒂𝒕.

8.1.2 Plasma Volume The plasma volume is given by the following relation based on the study by Cropp et al [7]. The same equation is used for males and females. 𝑽𝒑𝒍𝒔 = 𝟏𝟎𝟏.𝟐𝟎𝟖𝟐 ∗ 𝒍𝒐𝒈𝟏𝟎(𝑩𝑺𝑨) + 𝟑.𝟐𝟖𝟔𝟗 ∗ (𝟏 ― 𝒉𝒄𝒕)

S2.7

Liver Volume

Age dependence of liver volume previously defined by Noda et al [8] was implemented in the model S.2.7.1

Females

{

𝟎.𝟎𝟓𝟎𝟏𝟐 ∗ 𝑩𝑾𝟎.𝟕𝟖 𝒇𝒐𝒓 𝒂𝒈𝒆 < 𝟐𝟑 𝑽𝒍𝒊𝒗 = (𝟏.𝟎𝟕𝟐𝟖 ∗ 𝑩𝑺𝑨 ― 𝟎.𝟑𝟒𝟓𝟕) ∗ 𝒍𝒊𝒗𝒘𝒕𝟐𝟔𝟖 𝒇𝒐𝒓 𝒂𝒈𝒆 ≥ 𝟐𝟑 𝟏.𝟎𝟕𝟐𝟖 ∗ 𝑩𝑺𝑨𝟐𝟔𝟖 ― 𝟎.𝟑𝟒𝟓𝟕 Where, 𝒍𝒊𝒗𝒘𝒕𝟐𝟔𝟖 is liver weight 268 weeks after birth. This value defaults to 1.36 kg. And 𝑩𝑺𝑨𝟐𝟔𝟖 is the body surface area 268 weeks after birth. This value defaults to 1.75 m2.

S.2.7.2

Males

The equation defined above is also used for estimating liver volume in males. The default value for livwt268 in

males is 1.58 kg and the value for BSA268 is 2.01 m2

S2.8

Brain Volume

Brain volume in liters is calculated using the relation defined in Willman et al. [9]. For both males and females, the following equation defines the volume of brain

𝑽𝒃𝒓𝒏 = 𝟏𝟎 ∗

S2.9

𝟎.𝟑𝟏𝟓 + 𝒂𝒈𝒆 𝟗 + 𝟔.𝟗𝟐 ∗ 𝒂𝒈𝒆

Gut Volume

The total GI volume in the body is calculated using the equation defined by Price et al. [10]. This volume includes all individual sections of the GI tract. The model does not currently support a complex GI model. S.2.9.1

Females. 𝑽𝒈𝒖𝒕 = 𝟎.𝟎𝟐𝟕 ∗ 𝑳𝑩𝑴 𝑳𝑩𝑴 = 𝑩𝑾 ― 𝑽𝒇𝒂𝒕

S.2.9.2

Males 𝑽𝒈𝒖𝒕 = 𝟎.𝟎𝟐𝟏 ∗ 𝑳𝑩𝑴 𝑳𝑩𝑴 = 𝑩𝑾 ― 𝑽𝒇𝒂𝒕

Where, 𝑳𝑩𝑴 is the lean body weight in kilograms. Lean body weight is calculated by subtracting the body fat weight from the total body weight. Though Vfat is in L and BW is in kg, we assume a tissue density to be equal to that of water (𝟏𝟎𝟎𝟎 𝒌𝒈/𝒎𝟑).

S2.10 Rapidly Perfused Tissue Volume The relation between rapidly perfused tissue volume with age is given by the following equations for females and males .

S.2.10.1 Females 𝑽𝒓𝒑𝒇 = 𝟐.𝟒𝟔𝟒 ∗ 𝑽𝒈𝒖𝒕 S.2.10.2 Males 𝑽𝒓𝒑𝒇 = 𝟐.𝟓𝟗𝟔 ∗ 𝑽𝒈𝒖𝒕

In both cases 𝑽𝒈𝒖𝒕 is the GI volume at the given age as defined earlier

S2.11 Kidney Volume The age dependence of kidney volume reported by Bosgra et al. [3] is used in the model 𝑩𝑯 ) ((𝟏𝟎𝟎

𝟏.𝟗𝟑

―𝟐.𝟑𝟎𝟔 ∗

𝑽𝒌𝒊𝒅𝒏𝒆𝒚 = 𝒆

)

S2.12 Lung Volume The age dependence of Lung Volume reported by Bosgra et al. [3] is used in the model

―𝟐.𝟎𝟗𝟐 ∗

(( ) ) 𝑩𝑯 𝟏𝟎𝟎

𝟐.𝟏

𝑽𝒍𝒖𝒏𝒈 = 𝒆

S2.13 Skin Volume Equations defined by Clewell et al. are used to simulate age dependence of skin volume

𝑽𝒔𝒌𝒊𝒏 =

𝑩𝑺𝑨 ∗ 𝟏𝟎𝟎𝟎𝟎 ∗ 𝒅𝒆𝒑𝒕𝒉𝒔𝒌𝒊𝒏 𝟏𝟎𝟎𝟎

where 𝒅𝒆𝒑𝒕𝒉𝒔𝒌𝒊𝒏 is the average skin depth in cm. We use a default value of 0.15 cm.

S2.14 Muscle and Bone volume Bone and muscle volumes at various ages have been defined in the ICRP reference publication [12]. For these tissues, we linearly interpolate values between ages reported in the publication.

S2.15 Cardiac output Cardiac output was calculated for male and female based on BSA using the following equation [2] used to estimate the resting cardiac output in 𝑳/𝒉 𝑸𝒄 = 𝑩𝑺𝑨 ∗ 𝟑.𝟓 ∗ 𝟔𝟎

S2.16 Pulmonary Parameters Based on the data from IRCP reference publication , a nonlinear regression analysis was performed to describe the dead space in the lung (𝑫𝑺), the tidal volume (𝑻𝑽) and the breathing rate (𝑸𝒂𝒍𝒗) at resting. Lognormal, hyperbola, allosteric-sigmoidal, dose-response, exponential and Gompertz growth curves were used for fitting. Based on the degree of freedom and visual inspection, the best fit curve was chosen, and the equation derived from this curve was used to describe age-dependent changes in pulmonary parameters. S.2.16.1 Females 𝐥𝐧

𝑫𝑺 = 𝟎.𝟏𝟐𝟕𝟒 ∗ 𝒆

∗𝒆 (𝟎.𝟏𝟑𝟒𝟓 𝟎.𝟐𝟕𝟒 )

―𝟎.𝟏𝟔𝟎𝟗 ∗ 𝒂𝒈𝒆

𝟎.93278

𝑻𝑽 = 𝟎.𝟎𝟐𝟕𝟏𝟐 + (𝒂𝒈𝒆𝟎.𝟗𝟏𝟕𝟑) ∗

𝑎𝑔𝑒0.9173 + 20.56542 𝑄𝑎𝑙𝑣 = 0.62 ∗ 𝑄𝑐 + 82

S.2.16.2 Males ln

𝐷𝑆 = 0.17 ∗ 𝑒

∗𝑒 (0.1371 0.17 )

―0.13 ∗ 𝑎𝑔𝑒

―0.1357 ∗ 𝑎𝑔𝑒

𝑇𝑉 = 0.6842 ∗ 𝑒 ―2.765 ∗ 𝑒 𝑄𝑎𝑙𝑣 = 0.96 ∗ 𝑄𝑐 + 58

S2.17 Brain Tissue Perfusion Brain tissue perfusion was estimated using the equation described by Yoon et al. . S.2.17.1 Females

ln

(

―0.5 ∗

( )

3.064 𝑄𝑏𝑟𝑛 = ∗𝑒 𝑎𝑔𝑒

𝑎𝑔𝑒 (47.57 )

)

2

1.94

∗ 𝑄𝑐

S.2.17.2 Males ln

―0.5 ∗

𝑄𝑏𝑟𝑛 =

( )

3 ∗𝑒 𝑎𝑔𝑒

(

2

( ) 𝑎𝑔𝑒 41.92

1.876

)

∗ 𝑄𝑐

S2.18 Remaining Tissue Perfusion The relation between tissue volumes and perfusion defined by Clewell et al. has been used to determine perfusion for all other tissues.

𝑄𝑡𝑖𝑠 = 𝑄𝑡𝑖𝑠𝐶𝐴𝑑𝑢𝑙𝑡 ∗

𝑉𝑡𝑖𝑠𝐶 (𝑉𝑡𝑖𝑠𝐶𝐴𝑑𝑢𝑙𝑡 ) ∗ 𝑄𝑐

Where, 𝑄𝑡𝑖𝑠𝐶𝐴𝑑𝑢𝑙𝑡 is perfusion of the tissue as a fraction of the cardiac output in adults (25Y). This was obtained from the ICRP reference publication . 𝑉𝑡𝑖𝑠𝐶 is the volume of the tissue as a fraction of body weight. 𝑉𝑡𝑖𝑠𝐶𝐴𝑑𝑢𝑙𝑡 is the volume of the tissue as a fraction of body weight in adults (25Y).

9 Supplemental References 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12.

13.

Ruark, C.D., et al., Quantitative bias analysis for epidemiological associations of perfluoroalkyl substance serum concentrations and early onset of menopause. Environ Int, 2017. 99: p. 245-254. Wu, H., et al., Can the observed association between serum perfluoroalkyl substances and delayed menarche be explained on the basis of puberty-related changes in physiology and pharmacokinetics? Environ Int, 2015. 82: p. 61-8. Bosgra, S., et al., An improved model to predict physiologically based model parameters and their inter-individual variability from anthropometry. Crit Rev Toxicol, 2012. 42(9): p. 751-67. Keys, A., et al., Indices of relative weight and obesity. J Chronic Dis, 1972. 25(6): p. 32943. Gehan, E.A. and S.L. George, Estimation of human body surface area from height and weight. Cancer Chemother Rep, 1970. 54(4): p. 225-35. Yip, R., C. Johnson, and P.R. Dallman, Age-related changes in laboratory values used in the diagnosis of anemia and iron deficiency. The American journal of clinical nutrition, 1984. 39(3): p. 427-436. Gerd Cropp, J.A., Changes in blood and plasma volumes during growth. The Journal of Pediatrics, 1971. 78(2): p. 220-229. Noda, T., et al., Liver volume in children measured by computed tomography. Pediatr Radiol, 1997. 27(3): p. 250-2. Willmann, S., et al., Development of a physiology-based whole-body population model for assessing the influence of individual variability on the pharmacokinetics of drugs. J Pharmacokinet Pharmacodyn, 2007. 34(3): p. 401-31. Price, P.S., et al., Modeling interindividual variation in physiological factors used in PBPK models of humans. Crit Rev Toxicol, 2003. 33(5): p. 469-503. Clewell, H.J., et al., Evaluation of the potential impact of age- and gender-specific pharmacokinetic differences on tissue dosimetry. Toxicol Sci, 2004. 79(2): p. 381-93. Basic anatomical and physiological data for use in radiological protection: reference values. A report of age- and gender-related differences in the anatomical and physiological characteristics of reference individuals. ICRP Publication 89. Ann ICRP, 2002. 32(3-4): p. 5-265. Yoon, M., et al., Physiologically based pharmacokinetic modeling of fetal and neonatal manganese exposure in humans: describing manganese homeostasis during development. Toxicol Sci, 2011. 122(2): p. 297-316.

Abstract An outstanding challenge in the acceptance of alternatives to animal testing is the systematic incorporation of computational models into decision making pipelines. Fifteen years ago, the US EPA Office of Research and Development's framework for computational toxicology emphasized the need for computational methods to bridge the source-to-outcome continuum. This can be achieved by linking exposure estimation methods, physiologically based pharmacokinetic (PBPK) modeling, and

computational systems biology pathway modeling tools into a standardized framework. To that end, we have developed the Population Life-course Exposure to Health Effects Model (PLETHEM) suite, a modular open source modeling platform that provides users the ability to create, run, share, and audit PBPK models. The platform consists of a database of chemicals, QSAR models, life-stage specific physiological and metabolic parameters needed to parameterize PBPK models, an R-based engine to perform model simulations, and an interactive user interface to define and select parameter sets for the models. PLETHEM implements easy to use interfaces for a generic PBPK model and a high-throughput IVIVE model. These model interfaces along with the included databases provide capabilities necessary for rapid analysis of chemicals using PBPK modeling. PLETHEM includes the ability to run Monte Carlo analyses to investigate population variance and a set of life-stage equations to investigate life-stage-based sensitivities. The PLETHEM database also incorporates ontogeny profiles for key metabolic enzymes that can be used to calculate in vivo metabolic clearance using measured in vitro clearance. In addition, PLETHEM has an ability to link to a number of EPA and OECD exposure estimation programs. These models, which estimate exposures in the workplace and the general populations, can be used to drive PBPK model-based estimates of resulting internal exposures to support risk assessments. PLETHEM is now freely available as an R package through the Bitbucket and GitHub open source repositories.