Concentration-response gene expression analysis in zebrafish reveals phenotypically-anchored transcriptional...

Polycyclic aromatic hydrocarbons (PAHs) are ubiquitous environmental contaminants and are associated with human disease. Canonically, many PAHs induce toxicity via activation of the aryl hydrocarbon receptor (AHR) pathway. While the interaction between PAHs and the AHR is well-established, understanding which AHR-regulated transcriptional effects directly result in observable phenotypes and which are adaptive or benign is important to better understand PAH toxicity. Retene is a frequently detected PAH in environmental sampling and has been associated with AHR2-dependent developmental toxicity in zebrafish, though its mechanism of toxicity has not been fully elucidated. To interrogate transcriptional changes causally associated with retene toxicity, we conducted whole-animal RNA sequencing at 48 h post-fertilization after exposure to eight retene concentrations. We aimed to identify the most sensitive transcriptomic responses and to determine whether this approach could uncover gene sets uniquely differentially expressed at concentrations which induce a phenotype. We identified a concentration-response relationship for differential gene expression in both number of differentially expressed genes (DEGs) and magnitude of expression change. Elevated expression of cyp1a at retene concentrations below the threshold for teratogenicity suggested that while cyp1a expression is a sensitive biomarker of AHR activation, it may be too sensitive to serve as a biomarker of teratogenicity. Genes differentially expressed at only non-teratogenic concentrations were enriched for transforming growth factor-β (TGF-β) signaling pathway disruption while DEGs identified at only teratogenic concentrations were significantly enriched for response to xenobiotic stimulus and reduction-oxidation reaction activity. DEGs which spanned both non-teratogenic and teratogenic concentrations showed similar disrupted biological processes to those unique to teratogenic concentrations, indicating these processes were disrupted at low exposure concentrations. Gene co-expression network analysis identified several gene modules, including those associated with PAHs and AHR2 activation. One, Module 7, was strongly enriched for AHR2-associated genes and contained the strongest responses to retene. Benchmark concentration (BMC) of Module seven genes identified a median BMC of 7.5 µM, nearly the highest retene concentration with no associated teratogenicity, supporting the hypothesis that Module seven genes are largely responsible for retene toxicity.