Characterization of cell specific transcriptional responses to hepatotoxicants is lost in the averages of bulk RNA-sequencing (RNA-seq). Single-cell/nuclei RNA-seq technologies enable the transcriptomes of individual cell (sub)types to be assessed within the context of in vivo models. Single-nuclei RNA-sequencing (snSeq) of frozen liver samples from male C57BL/6 mice gavaged with sesame oil vehicle or 30 μg/kg 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) every 4 days for 28 days was used to demonstrate the application of snSeq for the evaluation of xenobiotics. A total of 19,907 genes were detected across 16,015 nuclei from control and TCDD-treated livers. Eleven cell (sub)types reflected the expected cell diversity of the liver including distinct pericentral, midzonal, and periportal hepatocyte subpopulations. TCDD altered relative proportions of cell types and elicited cell-specific gene expression profiles. For example, macrophages increased from 0.5% to 24.7%, while neutrophils were only present in treated samples, consistent with histological evaluation. The number of differentially expressed genes (DEGs) in each cell type ranged from 122 (cholangiocytes) to 7625 (midcentral hepatocytes), and loosely correlated with the basal expression level of Ahr, the canonical mediator of TCDD and related compounds. In addition to the expected functions within each cell (sub)types, RAS signaling and related pathways were specifically enriched in nonparenchymal cells while metabolic process enrichment occurred primarily in hepatocytes. snSeq also identified the expansion of a Kupffer cell subtype highly expressing Gpnmb, as reported in a dietary NASH model. We show that snSeq of frozen liver samples can be used to assess cell-specific transcriptional changes and population shifts in models of hepatotoxicity when examining freshly isolated cells is not feasible.

The environmental contaminant 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) elicits dose-dependent hepatotoxicity that includes fat accumulation, inflammation, and fibrosis that may progress to hepatocellular carcinoma. To further investigate these effects, RNA-Seq data were integrated with computationally identified putative dioxin response elements, and complementary targeted metabolomic and aryl hydrocarbon receptor (AhR) ChIP-Seq data from female C57BL/6 mice gavaged with TCDD every 4 days for 28 days. Data integration using CytoKEGG with manual curation identified dose-dependent alterations in central carbon and amino acid metabolism. More specifically, TCDD increased pyruvate kinase isoform M2 (PKM2) gene and protein expression. PKM2 has lower catalytic activity resulting in decreased glycolytic flux and the accumulation of upstream intermediates that were redirected to the pentose phosphate pathway and serine/folate biosynthesis, 2 important NADPH producing pathways stemming from glycolysis. In addition, the GAC:KGA glutaminase (GLS1) protein isoform ratio was increased, consistent with increases in glutaminolysis which serves an anaplerotic role for the TCA cycle and compensates for the reduced glycolytic flux. Collectively, gene expression, protein, and metabolite changes were indicative of increases in NADPH production in support of cytochrome P450 activity and ROS defenses. This AhR-mediated metabolic reprogramming is similar to the Warburg effect and represents a novel advantageous defense mechanism to increase anti-oxidant capacity in normal differentiated hepatocytes.

Dose-dependent differential gene expression provides critical information required for regulatory decision-making. The lower costs associated with RNA-Seq have made it the preferred technology for transcriptomic analysis. However, concordance between RNA-Seq and microarray analyses in dose response studies has not been adequately vetted. We compared the hepatic transcriptome of C57BL/6 mice following gavage with sesame oil vehicle, 0.01, 0.03, 0.1, 0.3, 1, 3, 10, or 30 μg/kg TCDD every 4 days for 28 days using Illumina HiSeq RNA-Sequencing (RNA-Seq) and Agilent 4×44 K microarrays using the same normalization and analysis approach. RNA-Seq and microarray analysis identified a total of 18,063 and 16,403 genes, respectively, that were expressed in the liver. RNA-Seq analysis for differentially expressed genes (DEGs) varied dramatically depending on the P1(t) cut-off while microarray results varied more based on the fold change criteria, although responses strongly correlated. Verification by WaferGen SmartChip QRTPCR revealed that RNA-Seq had a false discovery rate of 24% compared to 54% for microarray analysis. Dose–response modeling of RNA-Seq and microarray data demonstrated similar point of departure (POD) and ED50 estimates for common DEGs. There was a strong correspondence between RNA-Seq and Agilent array transcriptome profiling when using the same samples and analysis strategy. However, RNA-Seq provided superior quantitative data, identifying more genes and DEGs, as well as qualitative information regarding identity and annotation for dose response modeling in support of regulatory decision-making.

The aryl hydrocarbon receptor (AhR) is a highly conserved transcription factor that mediates a broad spectrum of species-, strain-, sex-, age-, tissue-, and cell-specific responses elicited by structurally diverse ligands including 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD). Dose-dependent effects on liver-specific and sexually dimorphic gene expression were examined in male and female mice gavaged with TCDD every 4 days for 28 or 92 days. RNA-seq data revealed the coordinated repression of 181 genes predominately expressed in the liver including albumin (3.7-fold), α-fibrinogen (14.5-fold), and β-fibrinogen (17.4-fold) in males with corresponding AhR enrichment at 2 hr. Liver-specific genes exhibiting sexually dimorphic expression also demonstrated diminished divergence between sexes. For example, male-biased Gstp1 was repressed 3.0-fold in males and induced 4.5-fold in females, which were confirmed at the protein level. Disrupted regulation is consistent with impaired GHR-JAK2-STAT5 signaling and inhibition of female specific CUX2-mediated transcription as well as the repression of other key transcriptional regulators including Ghr, Stat5b, Bcl6, Hnf4a, Hnf6, Foxa1/2/3, and Zhx2. Attenuated liver-specific and sexually dimorphic gene expression was concurrent with the induction of fetal genes such as alpha-fetoprotein. The results suggest AhR activation causes the loss of liver-specific and sexually dimorphic gene expression producing a functionally “de-differentiated” hepatic phenotype.