Integrated Epigenetic Mapping of Human and Mouse Salivary Gene Regulation

Abstract

Significant effort has been applied to identify the genome-wide gene expression profiles associated with salivary gland development and pathophysiology. However, relatively little is known about the regulators that control salivary gland gene expression. We integrated data from DNase1 digital genomic footprinting, RNA-seq, and gene expression microarrays to comprehensively characterize the cis- and trans-regulatory components controlling gene expression of the healthy submandibular salivary gland. Analysis of 32 human tissues and 87 mouse tissues was performed to identify the highly expressed and tissue-enriched transcription factors driving salivary gland gene expression. Following RNA analysis, protein expression levels and subcellular localization of 39 salivary transcription factors were confirmed by immunohistochemistry. These expression analyses revealed that the salivary gland highly expresses transcription factors associated with endoplasmic reticulum stress, human T-cell lymphotrophic virus 1 expression, and Epstein-Barr virus reactivation. DNase1 digital genomic footprinting to a depth of 333,426,353 reads was performed and utilized to generate a salivary gland gene regulatory network describing the genome-wide chromatin accessibility and transcription factor binding of the salivary gland at a single-nucleotide resolution. Analysis of the DNase1 gene regulatory network identified dense interconnectivity among PLAG1, MYB, and 13 other transcription factors associated with balanced chromosomal translocations and salivary gland tumors. Collectively, these analyses provide a comprehensive atlas of the cis- and trans-regulators of the salivary gland and highlight known aberrantly regulated pathways of diseases affecting the salivary glands.

Keywords: exocrine glands; gene expression regulation; gene regulatory networks; messenger RNA; salivary glands; systems biology.

Figure 1.

Integrated genomic analysis approach. (A) Identification of highly expressed and organ-specific transcription factors via analysis of adult mouse and human RNA expression atlases. (B) Utilization of DNase1 digital genomic footprinting for whole genome network inference and transcription factor activity profiling. (C) Immunohistochemical protein expression confirmation and tissue-type scoring. GRN, gene regulatory network; IHC, immunohistochemistry; SG, salivary gland; TF, transcription factor.

Figure 2.

Analysis of gene expression atlases reveals the architecture of transcription factor (TF) expression across multiple organ types and species. (A) Heat map of TF Spearman rank correlations across 37 human tissues. The salivary gland TF expression profile closely matches that of the pancreas. (B) Ortholog expression rank conservation among all expressed human and mouse genes. TF expression ranks for both species were compared in bins of increasing size between species. When all genes are analyzed, highly expressed human genes are unlikely to be equally highly expressed in the mouse. (C) Ortholog rank conservation between mouse and human for TFs only. Forty percent of the top 10 human TFs are also in the top 10 expression ranks within the mouse. The pancreas is a notable exception and exhibits a high degree of rank conservation for all genes.

Figure 3.

Whole genome network construction and transcription factor (TF) activity inference via DNase1 digital genomic footprinting. (A) Bivariate genomic footprinting (BaGfoot) TF activity profiling of the salivary gland (SG) versus the mouse heart DNase1 profile. PLAG1 and ZFP281 activity levels are clearly increased as compared with the heart. (B) BaGfoot TF activity profile of the SG relative to the mouse lung. Increased activity levels of PLAG1, EGR1, and ZFP281 are observed in both profiles. (C) The complexity of the cis-regulatory architecture for the Aqp5 gene. Ninety-one unique TF footprints are observed within 5 kb of the Aqp5 transcription start site. (D) The core SG gene regulatory network constructed from the DNase1 footprinting. Blue nodes: General TFs. Green nodes: TFs associated with SG development. Orange nodes: TFs associated with salivary pathology. Purple nodes: TFs associated with SG oncogenesis. The network connectivity of the highly expressed, specific, and active TFs was extracted from the DNase1 gene regulatory network and visualized. Increased edge density represents increasing numbers of TF footprints within a gene’s promoter. Sixty-six ZFP281 footprints were detected within 5 kb of the Nfib promoter. (E) Pathway enrichment testing identifies overrepresented connections between TFs and gene pathways. A total of 105 TFs exhibited statistically significant enrichments in connections to specific gene pathways. The connectivity of 8 highly active, expressed, or salivary-specific factors is presented here. ZFP281 regulates 16 gene pathways. PLAG1 was detected as enriched to regulate 10 pathways, including the endoplasmic reticulum protein–processing pathway and the Itpr2/Itpr3 genes. PWM, position weight matrix.

Figure 4.

Immunohistochemical staining of selected transcription factors (TFs) reveals 5 distinct TF expression clusters within the salivary gland. (A) Representative images for PLAG1, XBP1, ETV1, and IgG control. (B) Hierarchical clustering of TF subcellular localization. Five clusters of TF expression were noted within the 39 assessed TFs corresponding to acinar-specific, ductal-specific TFs expressed within both cell types and TFs with a primarily cytoplasmic localization.