Concentration-response gene expression analysis in zebrafish reveals phenotypically-anchored transcriptional responses to retene
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 hours post-fertilization after exposure to eight retene concentrations. The concentrations were selected to produce effects ranging from no phenotype to mortality and malformations in 100% of animals at 5 days post-fertilization. We identified a concentration-response relationship between retene teratogenicity and differential gene expression in both number of 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 AHR-dependent 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 7 genes identified a median BMC of 7.5 µM, nearly the highest retene concentration with no associated teratogenicity, supporting the hypothesis that Module 7 genes are largely responsible for retene toxicity.