Nonparametric Tests for Differential Histone Enrichment with ChIP-Seq Data

Abstract

Chromatin immunoprecipitation sequencing (ChIP-seq) is a powerful method for analyzing protein interactions with DNA. It can be applied to identify the binding sites of transcription factors (TFs) and genomic landscape of histone modification marks (HMs). Previous research has largely focused on developing peak-calling procedures to detect the binding sites for TFs. However, these procedures may fail when applied to ChIP-seq data of HMs, which have diffuse signals and multiple local peaks. In addition, it is important to identify genes with differential histone enrichment regions between two experimental conditions, such as different cellular states or different time points. Parametric methods based on Poisson/negative binomial distribution have been proposed to address this differential enrichment problem and most of these methods require biological replications. However, many ChIP-seq data usually have a few or even no replicates. We propose a nonparametric method to identify the genes with differential histone enrichment regions even without replicates. Our method is based on nonparametric hypothesis testing and kernel smoothing in order to capture the spatial differences in histone-enriched profiles. We demonstrate the method using ChIP-seq data on a comparative epigenomic profiling of adipogenesis of murine adipose stromal cells and the Encyclopedia of DNA Elements (ENCODE) ChIP-seq data. Our method identifies many genes with differential H3K27ac histone enrichment profiles at gene promoter regions between proliferating preadipocytes and mature adipocytes in murine 3T3-L1 cells. The test statistics also correlate with the gene expression changes well and are predictive to gene expression changes, indicating that the identified differentially enriched regions are indeed biologically meaningful.

Keywords: kernel smoothing; nonparametric testing; normalization; spatial histone profiles.

Figure 1

Histograms of two test statistics for the mouse adipogenesis ChIP-seq data, (A) Z_0λ and (B) Z_{0λ,H W}, for 9,874 genes with the maximum number of read counts in both day −2 and day 7 fewer than 5. The red curve in each plot represents the standard normal density.

Figure 2

Observed mouse adipogenesis ChIP-seq bin-counts for top 12 genes ranked by the test statistics Z_{0λ,W H} over the promoter region for day −2 (red) and day 7 (black). Vertical line represents the transcription starting site.

Figure 3

Observed ChIP-seq bin-counts over the promoter region for day −2 (red) and day 7 (black) for 12 genes with large Z_{0λ,W H} but small fold changes. Vertical line represents the transcription starting site.

Figure 4

Plots of gene expression fold changes as a function of two different test statistics. Top: proposed smoothing-kernel test statistics; bottom: fold changes. Left panel: genes with enriched H3K27ac binding at day −2; right panel: genes with enriched H3K27ac binding at day 7.

Figure 5

Plots of proportions of up/down-regulated genes in different intervals of the test statistics for the mouse adipogenesis ChIP-seq data, (A)–(B): proposed smoothing-kernel test statistics; (C)–(D): fold change statistics; (A), (C): genes with enriched H3 K27ac at day −2; (B), (D): genes with enriched H3 K27ac at day 7.

Figure 6

Model fit (left panel) and prediction (right panel) for log of the gene expression fold changes using the proposed statistics Z_{0λ,H W} (top panel) and fold changes (bottom panel) of six histone-modification ChIP-seq data in promoter, gene body, and downstream region.

Figure 7

Histogram of the test statistics $Z_{{\lambda }_i,WH}$ with the different bandwidths: (A) λ₁ = 5/280, (B) λ₂ = 20/280, (C) λ₃ = 60/280, (D) λ₄ = 90/280 for 9,874 genes with the maximum number of read count in both day −2 and day 7 fewer than 5 in mouse adipogenesis ChIP-seq data.

Figure 8

Top: Histogram of differential enrichment test statistics Z_new between two biological replicates of the ENCODE data for all 23,807 genes. Bottom: Histogram of differential enrichment test statistics Z_new between two cell types (B-lymphoblastoid cell vs HeLa-S3 cervical carcinoma cells) of the ENCODE data for all 23,807 genes. The red curve represents the standard normal density.