An atlas of the tomato epigenome reveals that KRYPTONITE shapes TAD-like boundaries through the control of H3K9ac distribution

In order to set the basis for an extensive understanding of genome organization and its interplay with gene activity in tomato, we undertook complementary cell biology and biochemical approaches. First, we investigated tomato nuclei organization using transmission electron microscopy. This analysis revealed that heterochromatin is strongly associated with the nuclear envelope (SI Appendix, Fig. S1 A and B). Second, we performed immuno-localization of a subset of 26 histone marks and RNA polymerase II (RNA Pol II, Fig. 1A). The localization of RNA Pol II across the nucleus coincided with 17 euchromatic histone marks, as well as with the facultative heterochromatin (FH)-related marks H3K27me2 and H3K27me3. As expected, the two constitutive heterochromatin related marks H3K9me2 and H3K27me1 exhibited a peripheric distribution, in proximity to the nuclear envelope. Interestingly, our results also revealed the localization of the histone mark H3K23me3 in the center of the nucleus and four other marks: H3K14me1, H3K23me1, H3K56me1, and H3K56me2, localized in the constitutive heterochromatic area. This suggested that these four histone marks are characteristic heterochromatic features (Fig. 1A and SI Appendix, Fig. S1C).

To go further, we mapped by Chromatin immunoprecipitation followed by sequencing (ChIP-seq) the distribution of these histone marks and RNA Pol II across the tomato genome. This analysis revealed a genome-wide view of the distribution of histone marks in heterochromatic and euchromatic genome compartments (see Fig. 1B for a global analysis and Fig. 1C for a detailed plot of chromosome 3). As expected, we observed a correlation between the immuno-localization analysis and the positioning of histone marks. The deposition of RNA Pol II revealed that gene expression occurs mostly at chromosome arms, in association with the distribution of permissive and facultative heterochromatin histone marks. Notably, we confirmed that the histone mark H3K23me3, and the histone marks H3K14me1, H3K23me1, H3K56me1, and H3K56me2, serve as euchromatin and constitutive heterochromatin signatures, respectively (Fig. 1C).

To identify the relationships between histone mark signatures and gene transcriptional activity, we conducted RNA sequencing (RNA-seq) analysis. Protein coding genes were classified into 16 groups according to their expression level (SI Appendix, Figs. S2 and S3). This classification of gene expression revealed that repressed and active genes were distinct in terms of epigenetic profiles and that a large repertoire of histone modifications is associated with each state (Fig. 1D and SI Appendix, Fig. S3). We found that six constitutive heterochromatinmarks were almost exclusively on genes belonging to the lowest quantile of expression, and this was also the case for the facultative heterochromatin marks H3K27me2 and H3K27me3. By contrast, different euchromatin-related marks were found to be associated with highly and lowly expressed genes (Fig. 1D). Interestingly, genes exhibiting low expression levels were enriched in H3K4me1 and H3K4me2, and their enrichment decreased along with higher transcript levels (Fig. 1D). Moreover, H3K4me1 enrichment was notably reduced in nonexpressed genes (SI Appendix, Fig. S3R). On the other hand, a high accumulation of gene transcripts is highly correlated with the enrichment of the other 18 marks which are signatures of euchromatin. With this in mind, we conducted a comprehensive correlation of the main constitutive heterochromatin signatures, including the constitutive histone marks (H3K14me1, H3K23me1, H3K56me1, and H3K56me2), as well as facultative heterochromatin and euchromatin signatures over genes belonging to the lowest quantile (i.e., not expressed). We found distinct epigenetic profiles that are associated with heterochromatin formation and gene silencing (Fig. 1E). Five different chromatin signatures were identified from our correlation analysis (SI Appendix, Fig. S4). The first, representing 21% of genes, was mainly enriched with heterochromatin-related marks, with a small proportion of genes showing association with the facultative heterochromatin mark H3K27me2. The second, accounting for 9% of genes, was highly enriched for H3K27me2 and to a lesser extent for constitutive heterochromatin marks. Notably, the third signature (8%) was associated with H3K27me3, indicating that they belong to polycomb repressive complex II-controlled genes and represent only a small subset of transcriptionally repressed genes. Remarkably, the fourth signature corresponded to a significant proportion of repressed genes (47%) that are associated with the histone mark H3K4me2. Finally, the fifth group represents 15% of repressed and has a less clear signature, appearing to be enriched for the constitutive heterochromatin-related histone marks H3K14me1, H3K23me1, and H3K27me1. These chromatin signatures revealed by ChIP-seq were in agreement with the immunolocalization using specific antibodies (Fig. 1A), revealing that the central and peripheral distribution of chromatin marks defines two distinct transcriptionally inactive compartments.

The distribution of the 26 different histone marks and RNA Pol II, allowed us to define 16 different chromatin states in the tomato genome (Fig. 1F and SI Appendix, Fig. S5). Among those, one corresponded to Polycomb-repressed genes, and two to H3K4me2 repressed genes, further highlighting that H3K27me3 deposition does not appear to be the only histone mark associated with gene transcriptional repression. In addition, each chromatin state is linked to positional features such as distal enhancers, proximal enhancers intergenic regions among others (SI Appendix, Fig. S6). This classification allows the definition of chromatin features associated with the regulation of gene expression and reveals the distinctive euchromatic and heterochromatic features that characterize the tomato genome.

To investigate the link between histone marks, chromatin states, and genome topology, we performed Hi-C assays and calculated the insulation index of regions enriched for RNA Pol II or each of the assessed histone marks (Fig. 2A) and chromatin states (SI Appendix, Fig. S7). The median of the frequency of interaction with the neighboring regions was lower for regions enriched with euchromatin-related marks than for constitutive heterochromatin, whereas facultative heterochromatin associated regions exhibited an intermediate contact behavior (Fig. 2B). Hence, we further analyzed a distinctive euchromatin mark (H3K9ac), as well as facultative and constitutive heterochromatin histone marks (H3K27me3 and H3K9me2, respectively). A Visual inspection of the interaction matrix revealed the widespread presence of interaction hotspots between genomic bins associated with the same histone modifications (Fig. 2C). Therefore, to evaluate the potential of a genomic region to establish a distinct compartmental domain characterized by particular chromatin marks, we produced a Pearson correlation map and organized it based on various chromatin features. In this process, we restructured the rows and columns within the correlation matrix. Instead of arranging them based on their linear sequence position, we sorted the bins in ascending order of the selected feature’s signal strength (Fig. 2D and SI Appendix, Fig. S8). Sorting by H3K9ac signal value, we observed a strong compartmentalization, as seen with chromosome 3 (Fig. 2D), and on all chromosomes (SI Appendix, Fig. S8). The matrix displayed a clear segregation which is characterized by an anticorrelated enrichment of H3K9ac and H3K9me2 histone marks. Reciprocally, sorting by the H3K9me2 signal value also resulted in a clear compartmentalization due to the strong anticorrelated enrichment between these two histone marks. By contrast, when sorting the chromosome by H3K27me3 signal value, we also observed compartmental segregation, but less pronounced as compared to H3K9ac sorting. Notably, this patterning occurred in all the tomato chromosomes (SI Appendix, Fig. S8), suggesting that H3K9ac and H3K9me2 are important driving elements for genome compartmentalization. Finally, to reveal the most frequent interactions among the 26 histone marks and RNA Pol II, we calculated the odds ratios (ORs) of their interactions (Materials and Methods). Analysis of all feature combinations of interacting loci across the genome revealed that chromatin interactions occur more frequently between regions with concordant epigenetic status (Fig. 2E and SI Appendix, Fig. S9). Altogether, our results hint at an intricate relationship between histone modifications and genome topology organization, and subsequently gene activity.

The thorough characterization of tomato chromatin states and their correlation with genome topology revealed a distinct and exclusive pattern of H3K9ac and H3K9me2 distribution that are closely tied to the spatial organization of euchromatin and heterochromatin compartments. To test the role of these two chromatin marks on genome topology we use the mutant kyp, which is known to affect H3K9me2 deposition (25). We mapped H3K9me2 and H3K9ac by ChIP-seq in WT and kyp plants (Fig. 3A and SI Appendix, Fig. S10). This analysis revealed that the levels of H3K9me2 in kyp plants diminished in discrete pericentromeric regions. Notably, the regions where H3K9me2 was lost in kyp showed a sharp deposition of H3K9ac. To investigate why H3K9ac deposition occurs in specific regions, we analyzed the associated motifs and found a significant enrichment for the TATA binding protein (TBP) motif bound by the TBP. This protein is known to interact with the transcription factor IID and SAGA two protein complexes that display a histone acetyl-transferase activity (26) (SI Appendix, Fig. S11). Further analysis revealed that the global reduction of H3K9me2 in kyp correlated with a strong deposition of H3K9ac (Fig. 3 A and B and SI Appendix, Fig. S12).

To determine the role of the ectopic H3K9ac deposition on genome topology, we performed Hi-C on WT and kyp plants. In the WT, TAD-like structures rich in H3K9me2 exhibit borders with high levels of H3K9ac (Fig. 3C and SI Appendix, Figs. S13–S16). Notably, the emergence of H3K9ac deposition within the same TAD-like structure kyp resulted in a reorganization of chromatin 3D conformation, determining different TAD-like borders and splitting the chromatin into two sub-self-interacting domains. The determination of the insulation index over genomic loci that display a switch from H3K9me2 to H3K9ac in kyp showed that this chromatin 3D reorganization is a global feature across the tomato genome (Fig. 3D and SI Appendix, Fig. S17). To understand the relationship between the level of H3K9ac and the insulation index, we divided the sites that gain H3K9ac in kyp into three categories (low, medium, and high levels) and plotted the insulation associated with each category. Interestingly, we observed a correlation between the H3K9ac level and the insulation (SI Appendix, Fig. S18). In addition, to determine to what extent this gain of H3K9ac and insulation correlate with chromatin accessibility, and gene expression, we performed ATAC-seq and RNA-seq. We found a positive correlation between these four features, suggesting that the local chromatin structure may play a role in 3D genome organization (SI Appendix, Fig. S18). Remarkably, a pile-up of chromatin interactions between genomic loci that display a switch from H3K9me2 to H3K9ac revealed that local associations (within TAD-like structures) are depleted (Fig. 3E and SI Appendix, Fig. S19), whereas distal contact points emerge (Fig. 3F and SI Appendix, Fig. S20), supporting the view that H3K9ac spots act as anchor points for newly formed TAD-like structures. Collectively, our results uncover the role of H3K9ac as a major determinant of chromatin organization and KYP as a modulator of chromatin topology in tomato (Fig. 4).

Solanum lycopersicum cv M82 was used for this study. The kyp mutant in M82 was donated by David C. Baulcombe (25). All seeds were sown in a commercial mix substrate under controlled conditions in growth chambers at 24 °C with a 16/8 h light/dark period. Genotyping of kyp plants was performed according to ref. 23. The fourth true leaf was used for all experiments in this study.

Leaves were cut into thin pieces in a petri dish and then fixed for 10 min using vacuum infiltration with a fixation solution containing 3% glutaraldehyde and 1% PFA in 0.1 M cacodylate and 250 µM CaCl₂. Following several washes with 0.1 M sodium cacodylate buffer, the plant tissue was postfixed with 1% OsO₄ and 1.5% potassium ferricyanide in 0.1 M sodium cacodylate for 1 h at room temperature. Subsequently, the plant tissues were washed with H₂O. Afterward, the tissues were dehydrated using a series of ethanol baths (10%, 20%, 30%, 50%, 70%, 90%, 100%), each for 30 min. The plant tissues were subsequently embedded in a mixture of resin and propylene oxide, then in pure epoxy resin, for a total duration of 72 h. Samples were polymerized for 24 h at 60 °C. Then, 80-nm ultrathin sections were prepared with an ultramicrotome (UC6, Leica Microsystems), counterstained with 2% uranyl acetate and lead citrate. Finally, the prepared plant material was observed with a transmission electron microscope operating at 80 kV (JEOL 1400).

Immunofluorescence labeling assays were conducted according to ref. 39. Four-week-old leaves were fixed, and nuclei were isolated, placed on a poly-lysine slide, and incubated overnight at 4 °C with primary antibodies, 400× diluted, listed in SI Appendix, Table S1. Slides were washed and incubated during 1 h at room temperature in dark with Goat anti-Rabbit Alexa Fluor Plus 488 (A11034 Invitrogen) and Goat anti-Mouse Alexa Fluor Plus 555 (A32727 Invitrogen) or with Goat anti-Rabbit Alexa Fluor Plus 555 (A32732 Invitrogen) and Goat anti-Mouse Alexa Fluor Plus 488 (A32723 Invitrogen) secondary antibodies (400× diluted). DNA was counterstained with DAPI in SlowFade Diamond Antifade mounting media. Slides were directly imaged on a confocal microscope (Zeiss Microsystems).

ChIP-seq assays were performed on the fourth leaves of 4-wk-old tomato plants using the antibodies listed in SI Appendix, Table S1, following the procedure described in ref. 40. ChIP-seq libraries were prepared from 10 ng of DNA using the NEBNext Ultra II DNA Library Prep Kit for Illumina (NEB) according to the manufacturer’s instructions. Two independent biological replicates were generated for each condition. The quality of DNA libraries was assessed with the Agilent 2100 Bioanalyzer. Single-end 75 bp reads were produced by the Illumina NextSeq 500 platform.

Adapter removal from raw sequencing data was performed using Trimmomatic v0.38 (41). Parameters for read quality filtering were set as follows: minimum length of 30 bp; mean Phred quality score greater than 30; and leading and trailing base removal with base quality below 5. Cleaned reads were aligned to tomato M82 genome v1 (42) using Bowtie2 v2.3.5 using mode “--very-sensitive” (43). Unmapped reads, PCR duplicates, and low mapping quality alignments (MAPQ < 30) were discarded for further analysis. MACS2 v 2.2.7.1 was used for peak detection (44). For broad peaks, parameters were set as follows: --broad --nomodel --extsize 147 -g 829069930 -p 0.05 --extsize 150 --bw 500. For sharp peaks, parameters were set as “-g 829069930 -p 0.05 --extsize 150 --bw 500.” Visualization of histone modification profiles was conducted with JBrowse (45) and r package EnrichedHeatmap (46).

To define genome-wide chromatin states for all histone modifications, ChromHMM v 1.22 was used with default parameters (47, 48). Peaks for each histone mark were selected individually using a P-value cutoff of 0.05, pileup > 36, and fold change > 1.5. The peaks were binarized with ChroHMM’s BinarizeBed method setting parameters “-b 150 -peaks -center.” Chromatin state models were learned jointly on all histone markers ranging from 15 to 20 states. Finally, a model of 16 chromatin states was selected for further analysis.

Total RNA was extracted using the Nucleospin RNA kit (Macherey-Nagel), according to the manufacturer’s instructions. Complementary DNA libraries were produced using 2 μg of total RNA with the NEBNext Ultra II Directional RNA library Preparation Kit (NEB) according to the manufacturer’s instructions. Quantification and quality control of libraries were performed using an Agilent 2100 Bioanalyzer. Libraries were sequenced in pair-end on the Illumina NovaSeq 6000 platform. The average targeted length obtained per read was 150 bp. Three biological replicates were used for each genotype.

Adapters removal from raw sequencing data was performed using Trimmomatic v0.38 (41). Cleaned reads were mapped to the tomato M82 genome v1 using Bowtie2 algorithm v2.3.5 (43) with “very sensitive” mode. After removing unmapped reads, PCR duplicates, and low mapping quality alignments (MAPQ < 30), aligned bam files were used for transcript expression estimation. Read density was assessed using the bamCoverage function of deeptools v3.5.0 (49) software with default parameters. Sixteen quantiles of gene expression values for the whole genome were calculated on the basis of the coverage bigwig generated in the previous step, considering only wild-type tomato RNAseq samples, using deeptools v3.5.0 (49) software with the multiBigwigSummary function.

Two biological replicates for both the wild type and kyp mutant were processed according to ref. 50 using the DpnII enzyme (New England Biolabs) and the DdeI enzyme (New England Biolabs). DNA library preparation was performed with the NEBNext Ultra II DNA library preparation kit (NEB). For the PCR amplification step, 10 cycles were used for libraries. Libraries were quantified and quality-checked using the Agilent 2100 Bioanalyzer, then sequenced on Illumina NextSeq 500 and NovaSeq 6000 instruments.

M.B. designed research; J.A., R.B.C., H.Z., Y.H., C.D.L., C. Bergounioux, C.D., C.G., and D.L. performed research; J.A., Q.W., X.H., J.A.-S., C. Bergounioux, J.M.P.-P., O.M., N.B., S.F., Y.Z., S.Z., M.C., M.M.M., J.G.-M., C.R., D.L., and M.B. analyzed data; C. Boulogne, C.D., and C.G. performed the transmission electron microscopy experiments; and J.A., L.I.P.-B., F.A., and M.B. wrote the paper.

The authors declare no competing interest.