Chromatin particle spectrum analysis: a method for comparative chromatin structure analysis using paired-end mode next-generation DNA sequencing

Abstract

Microarray and next-generation sequencing techniques which allow whole genome analysis of chromatin structure and sequence-specific protein binding are revolutionizing our view of chromosome architecture and function. However, many current methods in this field rely on biochemical purification of highly specific fractions of DNA prepared from chromatin digested with either micrococcal nuclease or DNaseI and are restricted in the parameters they can measure. Here, we show that a broad size-range of genomic DNA species, produced by partial micrococcal nuclease digestion of chromatin, can be sequenced using paired-end mode next-generation technology. The paired sequence reads, rather than DNA molecules, can then be size-selected and mapped as particle classes to the target genome. Using budding yeast as a model, we show that this approach reveals position and structural information for a spectrum of nuclease resistant complexes ranging from transcription factor-bound DNA elements up to mono- and poly-nucleosomes. We illustrate the utility of this approach in visualizing the MNase digestion landscape of protein-coding gene transcriptional start sites, and demonstrate a comparative analysis which probes the function of the chromatin-remodelling transcription factor Cbf1p.

Figure 1.

(A) DNA from MNase digested chromatin fractions purified from wild-type reference yeast strain BY4742 and an isogenic Δcbf1 mutant and separated by agarose gel electrophoresis. The lanes marked: ‘Total’ and ‘Input’ show the DNA species purified from MNase digested cells before and after a brief centrifugation step (to remove high-molecular weight material; see ‘Materials and Methods’ section) respectively. The ‘Input’ fractions were used for sequencing. Marker sizes are shown in kilobases and bands corresponding to DNA from mono-, di- and tri-nucleosomes are indicated to the right of the gel. (B) End-to-end distances of paired sequence reads reflect the distribution of chromatin particle input DNA. Graph of the number of Bowtie-aligned (9) paired-end reads obtained by Illumina GAIIx sequencing of material shown in Figure 1A ‘Input’ lanes versus paired-read end-to-end distance (Bowtie SAM format insert size value). The values shown on the x-axis indicate aligned sequences with the initial read mapping to the R/+/Crick strand of the reference genome. F/–/Watson strand reads show an identical distribution (Supplementary Figure S1). Peaks at ∼150, 300 and 450 bp are marked and correspond to mono-, di- and tri-nucleosome DNA fractions respectively. (C) Paired read size class dyad frequencies describe a landscape of MNase-protected species surrounding protein-coding gene TSSs. Using wild-type cell data, cumulative dyad frequencies for 15-bp bins, centred on and surrounding the protein coding gene TSSs mapped by (12), were plotted for paired read size classes 36 and 50– 700 bp in 25-bp intervals as a surface graph. Dyad frequencies within each size class were normalized to the mean number of reads per bin for that size class in order to plot each data set on the same y-axis. Normalized cumulative frequency values >1 are coloured yellow; values <1 are coloured blue. (D) Paired read dyad frequencies reflect distinct chromatin particle distributions surrounding TSSs. Cumulative frequency distributions of dyads from paired read size classes 50, 150 and 300 bp, were normalized and plotted relative to protein-coding gene TSSs as described above. The 150-bp dyad graph, which should yield peaks arising from nucleosomes, is shown relative to a cartoon diagram of the previously characterized TSS nucleosome landscape (4) with nucleosomes –1 and +1 flanking the nucleosome-free region (NFR) marked. (E) Nucleosome distributions derived by microarray analysis are identical to those determined with 150-bp dyad frequency distributions. Frequency distributions of 150-bp dyads were plotted with respect to budding yeast Chr I sequence as shown on the x-axis. The frequency distributions were plotted in 15-bp bins, with peaks smoothed to a three bin moving average. The y-axis indicates the number of aligned reads. The track marked ‘MNase array’ shows intensity values for MNase digested mono-nucleosome DNA hybridized to a tiling array taken from the data set described in ref. (18). Positive peaks in the array track are indicative of the presence of positioned nucleosomes and are co-incident with 150-bp dyad peaks. Open reading frames are shown as black boxes with +/F/Watson strand ORFs on top and –/R/Crick strand ORFs on bottom. The positions of nucleosome-free intergenic-regions are marked: ‘NFR’.

Figure 2.

(A) 50-bp dyads distribute as distinct peaks often within intergenic DNA and in the vicinity of DNA-binding protein DNase I footprints and/or known sites of bound transcription factors. The frequency distribution of 50-bp dyads was plotted relative to a 50-kb stretch of Chr I in 5-bp bins. The track marked ‘UW DNase’ marks the positions of the University of Washington digital DNase I footprints identified by (6). The track labelled ‘ChIP’ marks the positions of conserved transcription factor-binding motifs identified by chromatin immunoprecipitation experiments (14). The map of Chr I is as described for Figure 1E. (B) 50-bp dyad read numbers form a Cbf1p-dependent peak over known Cbf1p-binding motifs. Mean dyad peak frequencies from wild-type (black line) and a cbf1 mutant (red line) were plotted in 15-bp bins centred on and surrounding 128 Cbf1p-binding sites defined by ChIP and sequence conservation (14; p005_c3 data set). The P-value refers to the result of a Wilcoxon Ranked Sum test comparing the median 50-bp peak read numbers associated with the bins containing the Cbf1p-binding motif between wild-type and cbf1 mutant data sets. (C) 50-bp dyad read numbers form a peak over Abf1p-binding sites. Mean 50-bp dyad frequencies were plotted, and a P-value calculated, as for Figure 2B but centred on the Abf1p-binding sites identified by (14). (D) 50-bp dyad reads form a peak over Bas1p-binding sites identified by (14). (E) 50-bp dyad reads form a peak over Mbp1p-binding sites (14). (F) 50-bp dyad reads form a peak over Mcm1p-binding sites (14). (G) 50-bp dyad reads form a peak over Rap1p-binding sites (14). (H) 50-bp dyad reads form a peak over Reb1p-binding sites (14). (I) The 50-bp dyad read peak associated with Reb1p-binding sites is lost when the position of the binding sites is randomized relative to its normal position within an intergenic region.

Figure 3.

(A) Changes in Cbf1p binding and chromatin structure can be detected in dyad frequency data at a known Cbf1p-dependent gene-regulatory region. Fifty and one hundred and 50-bp dyad frequency distributions from wild-type and cbf1 mutant cells were plotted relative to the DRS2/MAK16 locus on Chr I. Dyad frequencies from the cbf1 mutant data set were plotted as negative numbers to aid comparison with the wild-type data. The DRS2 upstream region contains two potential Cbf1p-binding CACGTG sequences marked CACGTG motifs 1 and 2 on the ORF map. Only motif 1 (the position of which is marked with a dotted line) is bound by Cbf1p according to ChIP analysis (14,15). A Cbf1p-dependent peak in 50-bp dyad frequency centred over motif 1 is labelled, and 150-bp dyad peaks (nucleosome positions) which shift relative to wild-type in the cbf1 mutant are marked with asterisks. (B) Plot of mean 50-bp dyad frequencies comparing data centred on, and surrounding, 80 genomic CACGTG motifs which show similar Cbf1p-dependent changes in 50-bp peak presence to DRS2/MAK16 (Supplementary Data S2). Wild-type values are shown with a black line and cbf1 mutant values in red. The P-value calculation is as described in Figure 2. (C) Plot of mean 150-bp dyad frequencies, showing the trend in changes in nucleosome positioning, comparing data centred on and surrounding the 80 Cbf1p peak CACGTG motifs described above in wild-type and cbf1 mutant cells. (D) Cbf1p-dependent changes in paired read dyad frequency distributions occur across the spectrum of dyad sizes surrounding Cbf1p-associated CACGTG motifs. Normalized cumulative dyad frequency distribution values were calculated for 15-bp bins centred on and surrounding either the Cbf1p peak CACGTG motifs described above or the Reb1p-binding motifs identified by (14). Normalized values for wild-type and cbf1 mutant cells were calculated for dyad size classes ranging from 50 to 400 bp in 50-bp intervals, and the differences between wild-type and cbf1 mutant values plotted in surface graph form similar to that shown in Figure 1. Difference values <0.25 are coloured blue and >0.25 in yellow. (E) Plot of mean 150-bp dyad frequencies in wild-type and cbf1 mutant cells, comparing data centred on and surrounding Reb1p-binding motifs described in (14), showing no change in nucleosome positioning. (F) Plot of mean 150-bp dyad frequencies (nucleosome positions) comparing data from wild-type and cbf1 mutant cells centred on and surrounding 80 CACGTG motifs in the yeast genome which do not bind Cbf1p (14,15). A 150-bp dyad peak corresponding to a nucleosome positioned, on average, over such non-bound CACGTG motifs is indicated with an asterisk.

Figure 4.

Peaks in the frequencies of sub-nucleosome-sized (100 bp) paired read dyads locate to genomic regions suggested to contain MNase accessible, unstable nucleosomes. (A) Frequency distribution plots of 100- and 150-bp dyads over the GAL1-10 locus in wild-type yeast. A peak, which specifically occurs in both size classes within the divergent promoter region and corresponds to an MNase-protected species described by (7), is marked with an asterisk. (B) Plots of the cumulative frequency distributions of 100-bp dyads in wild-type cells centred over and surrounding all protein coding region TSSs (black line) or TSSs associated with highly expressed genes (red line). A peak in the 100-bp dyad frequencies at −100 bp in the highly expressed gene group, which corresponds to a region of highly MNase accessible chromatin described by (21) is indicated with an asterisk.