Marije Niemeijer (1), Witold Wiecek (2), Suzanna Huppelschoten (1), Audrey Baze (3), Celine Parmentier (3), Richard S. Paules (4), Lysiane Richert (3), Bob van de Water (1), Frederic Y. Bois (2)
(1) Division of Drug Discovery and Safety, LACDR, Leiden University, The Netherlands. (2) CERTARA, Simcyp Division, Sheffield, UK. (3) KaLy-Cell, Plobsheim, France. (4) National Toxicology Program, NIEHS, NIH, USA
Objectives: The activation of adaptive stress response pathways is a key-event in the development of chemical-induced liver injury. Between-patient differences in activation of those pathways could explain differences in susceptibility to drug-induced liver injury and allow mitigating them. Most of the time, those differences are estimated from small panels of primary human hepatocytes (PHHs) exposed in vitro. Using extensive transcriptome data on a large panel of PHHs, we aimed to assess, through modelling and simulation, the impact of panel size on inter-individual variance estimates. We modelled statistically the inter-individual variance for activation of four important stress-signalling pathways: oxidative stress, unfolded protein response, DNA damage response and NFκB signalling.
Methods: The BioSpyder TempO-seq technology (Yeakley et al. 2017) was used to interrogate systematically the transcriptome of a panel of 51 cryo-preserved PHHs from different donors. Diethyl maleate in vitro exposure was used to probe the activation of the oxidative stress pathway, tunicamycin was used for the unfolded protein response pathway, cisplatin for DNA damage response, and TNFα treatment for NFκB signalling. In each case, six concentrations encompassing a broad range were tested and responses were recorded after 8 and 24 hours. Individual concentration-gene activation modelling was performed with software BMDExpress 2 to obtain benchmark concentration levels (BMC) and maximum fold-change (maxFC) for the top activated genes for each individual (Phillips et al. 2019). As a next step, the joint (approximately bivariate lognormal) distribution of BMCs and maxFCs of these genes for each subject and each chemical at each time point was modelled using a population mixed-effect framework in Stan and GNU MCSim (Carpenter et al. 2017, Bois 1997). The model was then used to run Monte Carlo simulations of experiments with various PHH panel sizes. The results of those virtual experiments were analysed to create forest plots of between-donors variability estimates of gene activation, as they would be in routine assessments.
Results: There was evidence of a systematic burst effect (bimodal gene responses, see Judson et al. 2016) to various degrees. Tunicamycin (and the unfolded protein response pathway) was the most affected with a first cluster of genes responding at lower doses (geometric mean, GM, for BMCs at 0.01 μM, with a geometric SD, GSD, of about 10 fold) and a second cluster responding only at concentrations higher than 2 μM (GM 10 μM, GSD at 3). Diethyl maleate (and oxidative stress) responses were all grouped within only one cluster (GM 300 μM, GSD 4 fold). Within the early responding cluster, differences in gene between donors were significant (population GSDs between 2.6 and 6 fold for BMCs and between 1.5 and 2.2 fold for maxFCs). Responses were more homogeneous after 24 hours compared to 8 hours, indicating higher variability at the initial phase of stress. In general, low numbers of donors systematically under-estimated BMCs or maxFC CVs (by 25% to a factor 2 with only 3 donors, for example) and are imprecise in capturing the population variance of stress response activation (only 5% to 10% chance of falling between the 5th and 95th percentiles of the best CVs estimate).
Conclusions: Our results confirm that stress pathway activation analyses should focus on specific responses by excluding burst effect responses. Typical panel study designs have low power to identify the correct inter-individual variance. Their size should be increased, or “safety” factors derived from our analyses could be used to obtain conservative estimates of variability and DILI susceptibilities. We plan to refine our analyses to identify sub-groups of genes contributing most to susceptibility in responses for specific donors.
References:
[1] Yeakley JM et al. 2017. A trichostatin a expression signature identified by tempo-seq targeted whole transcriptome profiling. PLoS One, 12(5), e0178302.
[2] Phillips JR et al. 2019. BMDExpress 2: enhanced transcriptomic dose-response analysis workflow. Bioinformatics, 35(10), pp.1780-1782.
Reference: PAGE () Abstr 9331 [www.page-meeting.org/?abstract=9331]
Poster: Drug/Disease Modelling - Safety