meeting the following two criteria was considered passed: 1) a band of the correct size was detected in the wildtype sample, accounting for more than 50% of the total signal in the lane; 2) the intensity of the specific band decreases to less than 50% in mutants or knockdown samples. If immunoblots were unsuccess- ful, immunofluorescence was used as a characterization measure. For immunofluorescence tests, an antibody meeting the following two criteria were considered passed: 1) the immunofluorescence pattern must conform to expectations (for example, nuclear staining for a chromatin protein); 2) no immunofluorescence signal is detected in mutants. Some antibodies were considered validated if their ChIP profiles were consistent with those of a second, validated antibody to the same protein or to a known complex member. Antibody characterization data are part of the metadata provided with each dataset; they are available at www. modENCODE.org. Chromatin immunoprecipitation and microarray processing Protocols for the preparation of formaldehyde cross-linked chromatin from cultured cells, ChIP conditions, and array hybridization conditions are described in a recent article by Kharchenko and colleagues [14]. For all analyses, heterochroma- tin/euchromatin border positions previously defined by H3K9me2 enrichment were used [15]. Data analysis Processing of ChIP–chip data. The M-value (log2 ratio of signal intensities between ChIP and input) was calculated for each array dataset. Data normalization and identification of regions (or peaks) with significant enrichment were performed as described in Kharchenko et al [14]. At least two biological replicates were performed for each ChIP profile included in the analysis. The independent biological replicates were considered consistent if their target lists overlapped more than 75% or if the top 40% of the targets in each replicate had more than 80% in common. For correlation and other analyses, 500 bp bins were used to average the enrichment levels. For the heatmap visualization, metagene profiles were obtained with a scaled gene body of 3 kb. The extended regions of +/22 kb from the TSS and TES were included. Chromatin states model. The five-state chromatin annota- tions (based on the K-means algorithm) for the heterochromatin regions in BG3 and S2 cells were obtained from [15]. The number of states in this model was derived by combining states with similar enrichment patterns after starting with a higher number of states. Gene expression analysis. RNA-seq data for BG3 cells, S2 cells, and third instar larvae from Cherbas and colleagues were used for this analysis [64]. Reads Per Kilobase of exon model per Million mapped reads (RPKM) was calculated for each gene. Based on the distributions of expression levels, genes with an log10(RPKM+1).0.6 were considered expressed for BG3 and S2 cells; the threshold was log10(RPKM+1).0.4 for third instar larval data. This definition of expressed genes was used throughout the paper, e.g. in metagene analyses. Polymerase pausing analysis. For GRO-seq analysis, the pausing index (PI) was calculated according to the method developed by Larschan and colleagues [26], defining the PI of a gene as the ratio of signal at the 59end (first 500 bp) to the first 25% of the remaining gene body. For our analysis, genes were divided into three groups: euchromatin, pericentric heterochro- matin, and chromosome 4 using the border positions for S2 cells defined by Riddle et al. [15]. A PI threshold value of 10 was used. Drosophila Chromosome 4 Chromatin Structure To confirm these findings, we also estimated occurrence of pausing from ChIP-chip data using the PI proposed by Zeitlinger and colleagues [31], defined as the ratio between the maximum enrichment value around TSS (+/2300 bp) and the medium enrichment values over the gene body (600 bp downstream of TSS to the end of the gene) (Table S3B). Genes with a PI.4 are considered paused by this method. Genes shorter than 500 bp as well as overlapping genes were excluded from the analysis. Significance was determined using a permutation test. To compare method 1 and method 2, the per transcript data from method 1 was converted to per gene data, and overlap was estimated by comparing the top scoring 1,000, 1,500, and 2,000 genes. The overlap is 41% for 1,000 genes, 46% for 1,500 genes, and 50% for 2,000 genes. For GO analysis, we used the GOToolBox (http://genome.crg. es/GOToolBox/) to compare the 76 chromosome 4 genes with GO annotation to the whole genome reference [65]. A hypergemetrical test with Bonferroni correction was used to determine significance. For motif analysis, the fraction of promoters including a pause button motif (KCGRWCG) [29] was determined using CisGenome [66] with a window of +/ 260 bp around the TSS. Significance was determined using a permutation test. For the TRL (GAGA), inverted TRL, and Inr motifs, a window of +/2200 bp around the TSS was used. The Tm analysis was carried out as described by Nechaev and colleagues [33]. Mutant analysis. Differential gene expression comparing HP1a depletion to wildtype was performed using RNA-seq data (GEO accession GSE39083). RNA was prepared from third instar larvae using Trizol according the manufacturer’s recommenda- tions. rRNA-depleted cDNA libraries suitable for Illumina sequencing were prepared and sequencing was carried out by the Genome Technology Access Center (GTAC) at Washington University. FPKM (Fragment per Per Kilobase of exon model per Million mapped reads) values for each gene were calculated using the output from CuffLinks [67]. The cutoff value between expressed/silent genes was log2(FPKM+1) = 1.4. To obtain statistical significances of gene expression changes, we used several independent tools: BaySeq [68], EdgeR [69], and DESeq [70]. The set of genes detected in all analyses as significant was used for analysis. Pausing index for HP1a and POF mutants was calculated as described in [31]. To compare enrichment levels of H3K9me2, H3K9me3, POF and HP1a in various mutants with those in wildtype, we normalized the profiles using noise level signals as proposed in [71]. For each profile, the scaling factor was calculated as the ratio of the median absolute deviation to the lagged differences between mutant and WT. This is defined as median|dWTi - median(dWTi)|/ median|dMutanti - median(dMutanti)|, where dWTi = xWT i+1 - xWT i, , and xWT i, is the log-ratio of the ith probe in WT data. Similarly, dMutanti = xMutant i+1 - xMutant i, xMutant i, is the log-ratio of the ith probe in mutant data. The M-value profiles were then normalized by the factors. For RNA pol II profiles, we performed quantile normalization. Supporting Information Figure S1 Enrichment of histone marks and chromosomal proteins in S2 cells. A. Enrichment levels for the novel chromatin marks reported here are mapped onto the five main combinatorial chromatin states for heterochromatin as defined in Riddle et al 2011 [15]. Histone marks are shown in panel 1, chromosomal proteins in panel 2. Repeat enrichment and expression status for each state are shown in panel 3. Panel 4 illustrates the relationship PLOS Genetics | www.plosgenetics.org 16 September 2012 | Volume 8 | Issue 9 | e1002954