Dynamics of transcriptome and chromatin accessibility revealed sequential regulation of potential transcription factors during the brown adipose tissue whitening in rabbits

Abstract

Brown adipose tissue (BAT) represents a valuable target for treating obesity in humans. BAT losses of thermogenic capacity and gains a "white adipose tissue-like (WAT-like)" phenotype (BAT whitening) under thermoneutral environments, which could lead to potential low therapy responsiveness in BAT-based obesity treatments. However, the epigenetic mechanisms of BAT whitening remain largely unknown. In this study, BATs were collected from rabbits at day0 (D0), D15, D85, and 2 years (Y2). RNA-sequencing (RNA-seq) and the assay for transposase-accessible chromatin with high-throughput sequencing (ATAC-seq) were performed to investigate transcriptome and chromatin accessibility of BATs at the four whitening stages, respectively. Our data showed that many genes and chromatin accessible regions (refer to as "peaks") were identified as significantly changed during BAT whitening in rabbits. The BAT-selective genes downregulated while WAT-selective genes upregulated from D0 to Y2, and the de novo lipogenesis-related genes reached the highest expression levels at D85. Both the highly expressed genes and accessible regions in Y2 were significantly enriched in immune response-related signal pathways. Analysis of different relationships between peaks and their nearby genes found an increased proportion of the synchronous changes between chromatin accessibility and gene expression during BAT whitening. The synergistic changes between the chromatin accessibility of promoter and the gene expression were found in the key adipose genes. The upregulated genes which contained increased peaks were significantly enriched in the PI3K-Akt signaling pathway, steroid biosynthesis, TGF-beta signaling pathway, osteoclast differentiation, and dilated cardiomyopathy. Moreover, the footprinting analysis suggested that sequential regulation of potential transcription factors (TFs) mediated the loss of thermogenic phenotype and the gain of a WAT-like phenotype of BAT. In conclusion, our study provided the transcriptional and epigenetic frameworks for understanding BAT whitening in rabbits for the first time and might facilitate potential insights into BAT-based obesity treatments.

Keywords: ATAC-seq; BAT; epigenetics; rabbits; whitening.

FIGURE 1

Transcriptome analysis of BAT whitening in rabbits. (A) The appearance of BATs at four stages. (B) Differential analysis of the transcriptome. The scatter plot show the genes with log2(FC) > 1 or log2(FC) < - 1. The red and orange points show the genes with FDR <0.01. The grey points show the genes with FDR ≥0.01. The top 10 up- and downregulated DEGs that had been assigned official gene symbols were shown in corresponding directions. (C) qRT-PCR validation of selected genes in RNA-seq. The Ct-values of qRT-PCR are normalized to D0 and RN18S. Three independent samples are set per group and two technical replicates are set for one individual experimental replicate in qRT-PCR. (D) Heatmap analysis of differentially expressed genes (DEGs) based on the K-means clustering method. TPM values of genes are Z-scaled by row. The expression patterns of DEGs in different clusters are shown in corresponding box-plots. The numbers of genes in different clusters are presented in the left annotation of the heatmap. (E,F) Bar plots show the top 10 enriched Gene Ontology enrichment in biological process category (GO-BP) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways of genes in RC6 and RC2, respectively. Landscape and dynamics of chromatin accessibility during BAT whitening.

FIGURE 2

The landscape and dynamics of chromatin accessibility during rabbit BAT whitening. (A) Fragment lengths detected in a representative ATAC-seq library. (B) Genome-wide chromatin accessibility of samples in D0, D15, D85, and Y2. (C) Genomic annotation of peaks. (D) Enrichment of ATAC-seq reads around the gene bodies of a representative sample. The red and blue colors show the higher and lower ATAC-seq signals, respectively. (E) Hierarchical clustering analysis based on the signal intensity of genome-wide chromatin accessible regions (peaks). The signal intensities of peaks were normalized using CPM. (F) Heatmap analysis of 25,464 differential peaks (DPs) using K-means method. The color scale shows the relative signal intensities of peaks. SMT, the summit of peaks. The chromatin accessible patterns of DPs in different clusters of each replicate are shown in corresponding box-plots. (G) TF binding motif enrichments of peaks in cluster K4, K6, and K3. Integrated analysis of ATAC-seq and RNA-seq data revealed chromatin accessibility regulating gene expression.

FIGURE 3

Combined analysis of ATAC-seq data and RNA-seq data. (A) Analysis of four types of relationship between chromatin accessibility and gene expression based on different genomic elements. (B) Tracks of key adipose genes. The blue and green tracks showing the reads coverage of RNA-seq and ATAC-seq. Gene exons are marked in black. Gene information is marked by “gene name (chromosome strand: start - end)”. The locations of DPs were shadowed in grey. (C) KEGG pathway analysis of upregulated genes with increased peaks. Genome-wide footprinting analysis revealed sequential regulation of potential TF groups.

FIGURE 4

Genome-wide TF footprinting analysis during BAT whitening in rabbits. (A,B) The TF footprints at the loci of UCP1 and LEP. The arrows showed the gene transcription directions. The blue and orange tracks showed read coverages of RNA-seq and footprint binding score of ATAC-seq. The bound TFs in highlighted regions were shown in the zoomed axes. (C) Clustering of TF based on footprint score. Each row represents 1 TF, and each column represents a BAT whitening stage. TF footprint binding scores are Z-scaled by row. The blue and red color indicates the low and the high footprint score, respectively. Each cluster is associated with a mean trend line (left to right) and timepoint-specific boxplots respectively. (D) Integrative analysis of gene expression and footprint of TFs in SC2 and SC7. The bar plot showing the mean TPM of TF genes and the black bar showing the DEGs. The left and right panels of heatmap showing the expression profiles and footprint score of TFs across the four stages, respectively. The MA plots showing the chromatin accessibility changes (from D0 to D15) of target regions of corresponding TFs, with annotation of the number of increased or decreased peaks. (E) Integrative analysis of gene expression and footprint of TFs in SC1, SC4, SC5, SC6, and SC3. The bar plot showing the mean TPM of TF genes and the black bar showing the DEGs. The left and right panel of heatmap showing the expression profiles and footprint score of TFs across the four stages.