Gene expression profiling For a single microarray experiment, we isolated 20 central brains from 3rd instar larvae for both w1118 (i.e. wild type) and In(3LR)sep fly strains. Tissue samples from each fly strain were collected in 1 ml TriZol (Invitrogen) reagent, and total RNA was isolated according to the manufacturer’s protocol. For each fly strain we collected 4 independent RNA isolates. We amplified 1 mg of total RNA of each isolate using the Superscript RNA amplification system (Invitrogen). We labeled 1 mg of amplified RNA with the ULS system (Kreatech), using either Cy3 or Cy5 dyes (Amersham). We paired independent aRNA samples from w1118 (i.e. wild type) and In(3LR)sep fly strains and performed dual channel microarray hybridizations, including technical replicates (i.e. dye swaps). Our microarrays contain a 15k INDAC oligomer set that probes 93% of all protein coding genes in the D.melanogaster genome. We used the LIMMA package (R package for statistical computing; http://www.r-project.org) to obtain loess normalized log2-transformed In(3LR)sep over w1118 ratios. In addition, we determined log2-transformed signal intensities from all w1118 channels to obtain wild type expression data. Data processing Each independent replicate 4C experiment was normalized using a loess algorithm (R package for statistical computing; http://www.r-project.org). For each bait two replicates were combined and analyzed. We designed a ‘‘background filter’’ to correct for interactions that are a consequence of short nucleotide distance to the bait sequence. We only applied this filter to data that comes from the bait containing chromosome (e.g. chromosome 3R for Abd-B), because we consider other chromosome arms to be independent entities. First, we ranked data based on distance of the probe to the bait sequence (Figure S5). Next, we applied a loess algorithm to this ranked data (for each probe in the data a quadratic polynomial is fitted using a subset of the data), yielding smoothened representation of the data as a loess curve. The curve shows high values for probes located at short physical distance from the bait and gradually decreases as physical distance increases. As such, the loess curve represents a hypothetical linear conformation of the chromosome. Any long-range interaction will be detectable as a region with 4C signals that exceed the loess curve in that same position. To obtain these regions we subtracted the loess curve from the normalized 4C signals. Importantly, the degree of smoothening in a loess algorithm is determined by the smoothening parameter alpha (a), which determines the size of the subset of the data that is used to fit each local polynomial. We arbitrarily chose an a of 2, because it yields the most gradual degree of smoothening (Figure S5). However, we also tested smaller a values, which all gave similar results (Table S6), indicating that our results are independent of the parameters that we used. We used this corrected data as input for the domainogram and DID detection algorithms [24]. Originally, the DID algorithm discretizes the domainograms into a set of multi-layered and partially overlapping discrete domains of enrichment, referred to as BRICKs [24]. However, we extracted only the DID with the smallest window size of each set of partially overlapping discrete domains, because this will represent the center of gravity of a long-range chromatin interaction. Stringency of the DID algorithm is provided by the user defined bias factor c [24]. We empirically determined determined c, such that no DIDs were identified in randomly permutated data. To test the significance of overlap between DIDs and PcDs, we randomly permuted the DIDs’ coordinates across the bait containing chromosome a thousand times, while maintaining their size aspects. We used Pearson’s chi-square (x2) goodness of fit test to evaluate the difference between the observed coincidence of Polycomb Domains Interact in Nuclear Space DIDs and PcDs/K27Ds and the expected/random coincidence (average of 1000 random shufflings). To test whether interacting PcDs have distinctive genomic features, we divided all PcDs from the third chromosome into two categories: 1) interacting with at least one of our 4 PcG baits (Abd- B, Antp, Ptx1, trh); and 2) non-interacting with our baits. We used the non-parametric Mann-Whitney U test to test for significant differences between the two groups. We analyzed the following genomic features (based on the April 2004 release of the UCSC genome browser): size of domains, PRE sequences frequency, gene density, repetitive DNA sequences frequency and non-coding conserved sequences frequency. Our data demonstrated that interacting PcDs are significantly larger. In order to avoid size bias, we corrected the frequencies of PREs, genes, and non-coding conserved sequences by calculating the feature frequencies per 1 Kbp domain. All further statistical analyses were performed using the R package for statistical computing (http://www.r-project.org); scripts are available upon request. We deposited all our microarray data sets (4C, DamID, ChIP- chip, and gene expression profiles) to the GEO website under accession number: GSE26694 (http://www.ncbi.nlm.nih.gov/ geo/query/acc.cgi?acc = GSE26694). Supporting Information Figure S1 PcG proteins repress expression of Homeotic genes in central brain tissue. (A–B) In wild type larvae, Antp antibody stains nuclei in the ventral nerve cord (VNC) cells, while it is absent from the central brain cells (CB). (C–D) Reduction of the PcG protein PH by RNA interference results in detection of Antp in subsets of central brain cells. (E–F) In wild type larvae, Abd-B antibody only stains nuclei of cells in the most posterior compartment of ventral nerve cord (arrow). (G–H) Reduction of the PcG protein PH by RNA interference results in detection of Abd-B in a subset of in central brain cells. Left panels show projections of consecutive confocal sections through the CNS recorded at 4 mm intervals. Right panels depict a single confocal section through the right central brain. Scale bar is 50 mm. Anterior is always up. Images are false colored for better contrast. Found at: doi:10.1371/journal.pgen.1001343.s001 (3.80 MB TIF) Figure S2 Validation of PC binding and H3K27me3 chromatin maps in central brain tissue. (A–D) Chromatin maps of specific gene loci: Bithorax complex (A), Antennapedia complex (B), Ptx1 locus (C), and trachealess (trh) locus (D). Polycomb binding (top) and H3K27me3 levels (middle) are highly similar in these loci, yielding approximately the same boundaries for PcDs (blue) and K27Ds (green). Gene positions are given at the bottom, with target genes highlighted in cyan, while nontargets are grey. Some gene symbols are indicated, and the gene that contains the bait sequence is highlighted in bold. Scale bar is 50 Kbp. (E) On a genome-wide scale, a large proportion of all PcDs and K27Ds co-occur/ (partially) overlap in central brain tissue (observed), such overlap is not seen in a random overlap (expected). Random overlap was based on 1000 random permutations of domain positions, while maintaining their size aspects. Expected is the average random overlap after 1000 samplings, error bars indicate standard deviation. (F–G). The borders of all domains were aligned (relative position = 0) and the average H3K27me3 levels (F, green line) or PC binding (G, blue line) is plotted. Positions outside domains are indicated by negative values, and positive values are inside domains (x-axes). The y-axes represent average fold-change in signal. (F) Inside PcDs (blue area) the H3K27me3 levels are ,2 fold higher than outside (white area). (G) Likewise, Pc binding is PLoS Genetics | www.plosgenetics.org 12 March 2011 | Volume 7 | Issue 3 | e1001343