Limitations and possibilities of low cell number ChIP-seq

Chromatin immunoprecipitation (ChIP) is used to determine the genomic interaction sites between nuclear proteins and nucleic acids. Standard ChIP protocols used for genome-wide studies typically require large quantities of starting material, in the range of 10⁷ cells. The amount of material immunoprecipitated varies depending on the target protein and antibody employed, but is usually in the range of a few hundred picograms to tens of nanograms.

Over the past decade, genome-wide analysis of ChIP material has been possible by employing whole genome amplification techniques to produce microgram quantities of DNA for hybridisation to microarrays, commonly referred to as ChIP-chip. Recently, the application of high-throughput sequencing (HTS) to analyse the immunoprecipitates, commonly referred to as ChIP-seq, has replaced ChIP-chip as the preferred means of harvesting genome-wide data, and confers a number of advantages including higher resolution, improved quantification range, greater genome coverage, fewer artifacts, and lower cost[1]. However, the library preparation methods required to render immunoprecipitated DNA ready for HTS sequencing involve inefficient enzymatic steps and multiple purifications, each resulting in sample loss. As a result, ChIP-seq requires a similar starting amount of immunoprecipitated DNA (1–10 ng) to ChIP-chip, and nonetheless involves many cycles of PCR (typically 15–18 cycles). These limitations mean that it is still challenging to apply ChIP-seq to low numbers of cells.

There have been several recent papers detailing ChIP protocols applicable to low cell numbers (down to 100 cells) based on the inclusion of carrier DNA[2] or improvements in the efficiency of the immunoprecipitations themselves[3–7]. However, these techniques have been limited to locus-specific analysis of the precipitates by quantitative PCR. By incorporating whole genome amplification (WGA) methods, these techniques have been extended to allow genome-wide analysis by ChIP-chip from as few as 10,000 cells[8] or 1,000 cells[9].

Three recent publications detail the development of ChIP-seq protocols for use with low cell numbers, all based on formaldehyde crosslinked chromatin: Using standard Illumina library preparation procedures, Hitchler & Rice demonstrated ChIP-seq from 1 × 10⁶ human stem cells and 5 ng immunoprecipitated DNA[10]. Greater reductions have been possible using alternative library preparation methods, where cell requirements have been lowered to 10,000[11] (the method is also presented in more detail elsewhere[12]), and even 5000[13], with input DNA amounts below 50 pg. However, these methods employ either lengthy linear amplification procedures, or primer extension (4 cycles) and 15 cycles of PCR – all prior to a standard Illumina library prep entailing a further 17–18 cycles of PCR.

An alternative method of performing ChIP avoids the use of formaldehyde cross-linking and is thus known as native ChIP (N-ChIP)[14–17]. The advantages of N-ChIP over cross-linked ChIP (X-ChIP) include higher resolution, and lack of unspecific interaction caused by formaldehyde crosslinking. It has also been suggested that the N-ChIP method is more sensitive than X-ChIP, as epitopes may not be masked by cross-linked proteins or be themselves denatured by formaldehyde[15] making N-ChIP ideally suited to studies aiming to examine small cell numbers. However, N-ChIP is generally considered only applicable to histone proteins, although successful N-ChIP of MeCP2[18] and a handful of transcription factors including Runx2, Dlx6 and Sin3a have been reported[19].

We present here an N-ChIP method for genome-wide analysis by ChIP-seq, optimised for use with low starting cell numbers (here 200,000, divided into two immunoprecipitations of 100,000 each). This demonstrates ChIP-seq with 200-times fewer cells than a previously published method used as a benchmark for comparison. The performance of the optimised method was evaluated for read mapping, sensitivity and specificity at a range of starting cell numbers covering three orders of magnitude, starting with the published amount of 2 × 10⁷ cells / IP and reduced to a point where sensitivity was compromised, to determine the limits of the technique.

The ChIP method described here was developed using the N-ChIP method of Zhao and colleagues[20, 21] as a starting point. We therefore set the existing technique, using the published amount of 2 × 10⁷ cells per immunoprecipitation, as a benchmark against which to compare the performance of our method at a series of decreasing cell numbers. The new method presented here significantly shortens the procedure by eliminating the need for dialysis, and incorporates modifications optimised for low cell numbers.

Chromatin prepared using the two methods was prepared from cultured CD4⁺ lymphocytes, and immunoprecipitated with anti-H3K4me3 antibody. Enrichment at positive and negative control loci (see methods section) was measured by quantitative PCR prior to generation of Illumina sequencing libraries. Each ChIP-seq library was then sequenced on a single lane of an Illumina GAIIx sequencer, and generated chromatin profiles typical of H3K4 trimethylation, from which peaks were called. At this point it was determined that additional sequencing was required to saturate peak calling (see below) in the lowest cell number employed, so additional sequencing was performed for this sample using a single lane of an Illumina HiSeq 2000 machine.

The total number of reads generated for each library and the results of aligning these libraries to the human genome are summarized in Figure 1. As cell numbers are reduced, the number of unmapped reads increases. In addition, the percentage of mapped reads derived from duplicate reads increases. A sample of unmapped reads was aligned against the GenBank nucleotide database[22] using blast[23]. The results showed that a small proportion of unmapped reads in all cases represent those containing sequencing errors that fail to map to the human genome (the BWA algorithm does not tolerate more than 2 mismatches in its 32 bp seed). The remainder of unmapped reads, which increase with decreasing input cell number, fail to map with high confidence to any sequence in the GenBank database and are apparently PCR amplification artifacts. The increased level of duplicate reads seen at lower cell numbers is also assumed to be introduced during the PCR amplification (18 cycles) required as part of Illumina library preparation. The same number of PCR cycles were applied to each sample, to facilitate inter-sample comparisons. A disadvantage of this approach is that a greater number of cycles than necessary to generate sufficient library to sequence, were applied to the higher cell number samples. Despite this, only in the low cell number samples did the decreased amount and complexity of the input material, lead to high proportions of duplication during amplification.

ChIP-seq profiles from each experiment can be seen in Figure 2a and2b, and by visual inspection, peaks of H3K4me3 can be seen to occupy promoter regions as expected from previous studies[24]. To demonstrate the increased sensitivity of the new method relative to the benchmark, an additional profile generated using the benchmark at a low cell number is shown (Figure 2a and2b). The resolution afforded by the MNase digestion allows identification of individual nucleosome positions (Figure 2b). The relationship between H3K4me3 and transcription start sites (TSS) was confirmed by plotting H3K4me3 levels relative to transcription start sites (Figure 2c), which reproduces earlier findings showing depletion of histones at the TSS itself, and a series of clearly positioned nucleosomes upstream and downstream[20]. Furthermore, immunoprecipitation with H3K4me3 is strongest at highly expressed genes, and diminishes with decreasing expression level (Figure 2c).

To define genomic regions of H3K4 trimethylation, peak calling was performed using MACS[26], using only uniquely mapping, non-duplicate reads. Performing peak calling whilst including duplicate reads led to the appearance of high numbers of non-specific peaks, particularly in the lowest cell number sample (data not shown). Calling peaks when correcting for background control (sequencing libraries prepared from each sample’s ChIP input DNA) made negligible difference to the number of peaks recognised; 0.03 – 1.43% of peaks were no longer called when using a control dataset. The total number of peaks called for each cell number is summarized in Table 1. The number of peaks called falls appreciably only at the lowest cell input number tested, to below 75% of the number called in the benchmark. Crucially, despite the lower numbers of uniquely mapped reads recovered from the lowest cell number sample, peak calling was saturated for all samples: When decreasing proportions of total available reads were used to call peaks (Figure 3a), the absence of a reciprocal relationship between read count and peaks called indicates that all samples are approaching saturation (i.e. where adding more reads will not increase the number of peaks called). Therefore, the lower number of peaks called using only 20,000 cells / IP is due to the reduced number of useful reads (non duplicated and uniquely mapping) available as cell numbers fall, and not simply because more sequencing is required.

Protocol & cell number / ChIP	Bench’	New	New	New	New	New
	2 x 107	2 x 107	2.5 x 106	5 x 105	1 x 105	2 x 104
Total number reads	13 402 262	13 798 839	14 136 243	16 693 519	15 601 563	56 904 707
Number unique, non-duplicate reads	6 011 891	7 0330 709	5 794 519	3 463 886	2 423 126	661 591
Number of peaks called (fraction relative to benchmark)	16 545	16 244	17 054	15 636	14 771	12 296
	(1)	(0.98)	(1.03)	(0.95)	(0.89)	(0.74)
Sensitivity relative to benchmark	1	0.93	0.96	0.89	0.85	0.69
Specificity relative to benchmark	1	0.94	0.93	0.95	0.97	0.98

Peaks were called using MACS, allowing no ambiguously mapping or duplicate reads.

Using the method of Zhao and colleagues as a benchmark, the overlap with peaks called from the new method was evaluated as a measure of sensitivity (Table 1). Sensitivity was well maintained down to 1 × 10⁵ cells / IP, where 85% of peaks could still be detected. As expected from the reduced number of peaks called, sensitivity fell in the lowest cell number sample to 70%. Peaks were not lost randomly with reduced cell numbers. Rather, the same peaks were affected in each sample, with the preferential loss of those with lowest significance (fewest reads) evident as cell numbers were reduced.

Importantly, peak position was not adversely affected by lowering cell number. As can be seen in Figure 3b, peak overlaps between the benchmark and other datasets were close to 100% (i.e. not partially overlapping, or within an arbitrarily chosen window size). However, peak width was reduced in the lowest cell number samples (Figure 3c). Whilst every effort was made to ensure similar MNase digestion between the samples, we cannot exclude that the narrower peaks seen with lower cell numbers are due to increased digestion in these samples. However, the effect may be entirely explained by the lower number of reads available for peak calling in these samples. We are currently unable to separate the two possibilities.

In addition, extra peaks not present in the benchmark dataset were used to calculate a measure of specificity (Table 1). Specificity was not affected by scaling down cell numbers, with all datasets having greater than 90% of called peaks “on target” when comparing to the benchmark. It should be noted that this assumes that the dataset of Zhao and colleagues represents a gold standard and that additional peaks are false positives, which is not necessarily the case. Comparing the signal intensity of peaks at each location revealed a strong correlation, which deteriorated at the lowest cell number of 2 × 10⁴ cells per IP (Figure 3d).

To demonstrate the application of the new method to other histone modifications, the transcriptionally repressive H3K27me3 mark was also examined. Figure 4a shows the ChIP-seq profile of H3K4me3 and H3K27me3 generated from 100,000 cells / IP at the active STAT4 and inactive MYO1B loci (a comparable image can be found in the paper of Barski et al.,[20]). The mutually exclusive nature of H3K4 and H3K27 trimethylation is clearly visible in these profiles.

In this study we have developed and employed a rapid N-ChIP technique applicable to small cell numbers, which functions well down to 100,000 cells / IP. Whilst this limit is higher than the requirements reported using alternative HTS library preparation methods[11–13], it minimizes the use of nucleic acid amplification and associated risk of bias in the data. As the first N-ChIP protocol tailored to low cell numbers, it therefore offers an attractive alternative method to map the genome-wide distribution of histone modifications. The success of the method using standard HTS library preparation techniques may reflect the reported higher efficiency of N-chip relative to X-ChIP[15]. We would also like to stress the importance of selecting specific antisera, such as using a peptide array as employed here, to the success of this or any ChIP protocol using low input cell numbers. We have demonstrated the application of this method to immuno-purified CD4⁺ and CD8⁺ primary lymphocytes, thus avoiding the need for cell culture, which risks altering epigenetic modifications. Nonetheless, the study highlights the need for careful monitoring of sequence read mapping in the analysis of ChIP-seq data from limited cell numbers, to identify sources of wasted reads and ensure sufficient coverage for reliable peak calling.

As demonstrated here, PCR amplification from limited ChIP input material led to a reduction of mappable and unique reads through losses to amplification artifacts and duplicate molecules, which must be removed for reliable peak calling. The accumulation of such undesirable amplification artifacts in alternative techniques such as LinDA and nano-ChIP-seq[11–13] has not been adequately addressed to date, and requires further examination. Notably, the method presented here entails fewer amplification cycles than these other techniques. The introduction of bias by genome-amplification techniques such as linker-mediated PCR and WGA have been documented[29]. It is therefore desirable to minimise the number of amplification cycles whenever possible. It is likely that the levels of PCR artifacts could be reduced by employing alternative amplification conditions, or through the use of alternative polymerases.

It has been shown that both MNase digestion and sonication demonstrate some DNA cleavage sequence preference, which necessitates the use of control digestion / sonication of naked DNA when interpreting HTS data with the aim of identifying nucleosome positions or regions of high accessibility such as promoters[30, 31]. The existence of favourable cleavage sites will inevitably lead to a number of “duplicate” reads that are in fact biological in origin rather than amplification artifacts. This effect is likely to be more pronounced in MNase digested chromatin, as the preferential digestion of MNase in inter-nucleosomal linker regions further restricts genomic cleavage sites. It is therefore possible that the levels of duplication seen here will not be so high when using formaldehyde cross-linked chromatin, although this remains to be tested. Using the data presented here, inclusion of duplicate reads led to the appearance of many non-specific peaks, so it remains challenging to differentiate genuine biological duplicate reads from those arising as the result of the PCR process. The percentage of duplicate reads derived from PCR could be further reduced by employing paired-end reads rather than the single reads typically used, as here, for ChIP-seq analysis. Only in the situation where molecules have exactly the same length and genomic position, would paired reads be counted as a duplicate. Alternatively, utilizing multiple combinations of indexed adapters, termed digital sequencing[32], would allow the differentiation of biological and PCR-derived duplicate reads.

In addition to the problems of amplification artifacts and duplicate reads, ChIP-seq from limiting cell numbers is challenging due to the reduced complexity of immunoprecipitated material recovered. In the case of CD4⁺ cells presented here, sensitivity and reproducibility were not badly affected down to 100,000 cells / IP, but below this were reduced, although useful results were nonetheless obtained (70% of peaks detected by benchmark method). Further improvements are likely possible by increasing the efficiency of immunoprecipitation, DNA purification, and sequencing library generation methods.

The future adoption of epigenetic analyses in research and diagnostic procedures will require techniques that allow analysis of specific cell types, sub-populations and small biopsy samples. For example, the N-ChIP technique has been optimised for analysis of the Schistosoma parasite[33], and techniques have been developed for ChIP from particular structures such as the hippocampus[34]. We have demonstrated here the application of the enhanced method to lymphocytes isolated from human twins. The cells used for study were enriched using antibody-driven magnetic cell sorting, derived from a starting material of 10–20 mls blood. It would not have been possible to study such small, easily obtained biological samples using previously published N-ChIP procedures. Further improvements in both ChIP protocols and HTS library preparation methods, such as single molecule sequencing[35] promise to increase the possibilities for epigenetic studies from these and other challenging samples.

Using an enhanced native ChIP-seq method, we have detailed known but hitherto uncharacterised problems of performing ChIP-seq from limited cell numbers. Using our method, high quality results were obtained from 200,000 cells starting material (using 100,000 cells per IP), increasing possibilities for the study of rare cell populations and biopsy samples without the need for cell culture. It was possible to generate results from as low as 20,000 cells per IP, but at a cost to sensitivity, where only 70% of known peaks could be detected, so we declare the limit of this method to currently require 100,000 cells per IP. This method offers an alternative to other previously published methods for low cell number ChIP-seq that entails fewer cycles of amplification with associated risk of data bias. Nonetheless, as cell numbers are reduced, the amount and complexity of immunoprecipitated material is reduced. Amplification of this material leads to a corresponding rise in PCR duplicates and unmapped reads, which may necessitate costly additional sequencing, an effect that requires further characterisation in comparable methods. Additional improvements in the efficiency of immunoprecipitation and HTS library generation techniques are desirable to bring us closer to the goal of single cell analysis.