{'fraction_non_alphanumeric': 0.03557405186804919, 'fraction_numerical': 0.020949425125202793, 'mean_word_length': 5.055239179954442, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Reorganized 3D Genome Structures Support Transcriptional Regulation in Mouse Spermatogenesis

Summary Three-dimensional chromatin structures undergo dynamic reorganization during mammalian spermatogenesis; however, their impacts on gene regulation remain unclear. Here, we focused on understanding the structure-function regulation of meiotic chromosomes by Hi-C and other omics techniques in mouse spermatogenesis across five stages. Beyond confirming recent reports regarding changes in compartmentalization and reorganization of topologically associating domains (TADs), we further demonstrated that chromatin loops are present prior to and after, but not at, the pachytene stage. By integrating Hi-C and RNA-seq data, we showed that the switching of A/B compartments between spermatogenic stages is tightly associated with meiosis-specific mRNAs and piRNAs expression. Moreover, our ATAC-seq data indicated that chromatin accessibility per se is not responsible for the TAD and loop diminishment at pachytene. Additionally, our ChIP-seq data demonstrated that CTCF and cohesin remain bound at TAD boundary regions throughout meiosis, suggesting that dynamic reorganization of TADs does not require CTCF and cohesin clearance.


INTRODUCTION
The process of chromatin organization in three-dimensions (3D) has been largely mysterious since the first recognition of chromatin structures (Dounce et al., 1972). Hi-C technology (high-throughput genome-wide chromatin conformation capture, Lieberman-Aiden et al., 2009) is developed from chromosome conformation capture (3C) (Dekker et al., 2002), which enables the global characterization of the 3D chromatin architectures. Results from Hi-C (original version and its derived versions such as in situ Hi-C) revealed that chromatin architecture comprises a hierarchy of structures in mammals, with chromatin A/B compartments at multi-megabase scale, topologically associating domains (TADs) at hundreds of kilobases scale, and chromatin loops at kilobases to hundreds of kilobases scale (Dixon et al., 2012;Lieberman-Aiden et al., 2009;. Previous studies have indicated that TADs (Flavahan et al., 2016;Hnisz et al., 2016;Lupiá ñ ez et al., 2015) and/or enhancer-promoter loops impact gene expression and cellular physiology (Bonev et al., 2017;Isoda et al., 2017;Weintraub et al., 2017) and transcription elongation also affects the 3D genome organization in specific regions (Heinz et al., 2018). However, the elimination of TADs and chromatin loops by rapid degradation of the cohesin complex only modestly affects gene transcription programs (Rao et al., 2017), and inhibiting transcription with a-Amanitin does not substantially affect the establishment of 3D chromatin organization during early development (Du et al., 2017;Ke et al., 2017). Thus, much remains to be learned about how differentially arranged 3D chromatin architectures control key biological processes.
Mammalian spermatogenesis is known to involve extensive chromatin re-organization, of which analyses using imaging techniques have revealed these dynamic changes (Hao et al., 2019;Sassone-Corsi, 2002). During spermatogenesis, a small number of primitive type A spermatogonia (priSG-A) differentiate into more developed spermatogonia including type A spermatogonia (SG-A) and type B spermatogonia (SG-B) in sequence via mitosis (de Rooij, 2001). Subsequently, meiosis is initiated and marked by doublestrand break (DSB) formation and is followed by full synapses of homologous chromosomes in pachytene spermatocytes (pacSC) (Tong and Lin, 2018). After meiosis I, spermatocytes rapidly divide in the ensuing meiosis II stage to form haploid round spermatids (rST) and spermatozoa (SZ) via spermiogenesis. Recent

Higher-Order Chromatin Structures Are Reorganized during Mouse Spermatogenesis
To investigate the relationship between 3D genome structures and transcriptional regulation during mammalian spermatogenesis, we isolated spermatogenic cells from five different stages of mouse spermatogenesis (priSG-A, SG-A, pacSC, rST, and SZ) ( Figure 1A) using the unit gravity sedimentation procedure (STA-PUT method) (Bellvé et al., 1977;Bryant et al., 2013;Gan et al., 2013;Hur et al., 2016;Korhonen et al., 2015;Liu et al., 2015;Luense et al., 2016;Wang et al., 2019) except for the SZ (see Transparent Methods). After validating the purities of these isolated cells ( Figure S1A) by two methods, morphological characterization ( Figure 1A) and immunofluorescence staining with cell-specific markers ( Figure S1B), we performed in situ Hi-C with these different spermatogenic cells, respectively. The correlation analyses of the Hi-C data supported the high quality of our datasets among the biological replicates (Table S1 and Figure S2A), and we randomly sampled the same amount of Hi-C valid pairs (172,252,595) from spermatogenic cells of each stage to do the following analysis.
The chromatin interaction heatmaps revealed dramatic architectural reorganization throughout mouse spermatogenesis, as represented by chromosome 1 (Figure 1B), which confirmed recent reports using sub-stages (Vara et al., 2019;Wang et al., 2019). The A/B compartment patterns on autosomes were largely unaltered during mouse spermatogenesis ( Figure 1C), whereas the strength, shown by saddle plots (Imakaev et al., 2012) ( Figure S2B) and compartment strength (A-A and B-B compartments interaction strength relative to A-B compartments interaction strength) ( Figure S2C) varied throughout the process. The strength of compartmentalization decreased from priSG-A to SG-A and dropped to the lowest level in pacSC before again increasing in rST and finally reaching its highest detected level in SZ. We also iScience 23, 101034, April 24, 2020 observed that the distinct TAD structures in priSG-A and SG-A largely disappeared in pacSC and rST. Interestingly, they were largely restored in the SZ stage ( Figures 1D and 2A). This trend of TAD reorganization was most clearly manifested by global analysis of averaged observed/expected interaction frequencies centered on these boundaries ( Figure 1D) and was also very apparent as seen on a random representative genomic region such as Chromosome 1, 85-90 Mb region ( Figure 2A). We further adopted an ''insulation score (IS)'' approach to calculate the strength of TAD boundaries, which analyzes the relative chromatin interaction frequency across a boundary (Crane et al., 2015). We found that the IS values decreased gradually from priSG-A to pacSC, then increased from rST to SZ again, and the IS value of SZ almost exactly matched that of priSG-A ( Figure 2B).

priSG-A and SZ Share Similar 3D Chromatin Interactions and Folding Patterns
With these IS values, we performed hierarchical clustering and found that TAD structures that once existed in priSG-A, and disappeared in pacSC, had reappeared in SZ, making these two stages (priSG-A and SZ) clustered together ( Figure 2C). These data led us to examine the 3D genome structures at the different spermatogenic stages in more detail. We found that the global 3D chromatin interactions were similar between priSG-A, SG-A, and SZ ( Figures 3A-3C), whereas pacSC showed stronger short-range chromosome interactions (<5 Mb) but weaker long-range chromosome interactions (>5 Mb) as compared with SG-A (Figure 3A). Specifically, the contact frequency P(s) with genomic distances demonstrated that 3D genome structures were reorganized during spermatogenesis ( Figure 3C), showing that priSG-A, SG-A, and SZ had similar chromatin folding patterns and that pacSC represented an extreme case, consistent with the above IS hierarchical clustering data.
Recent studies showed that mitotic chromatin exhibits a folding pattern distinct from that of interphase chromatin (Gibcus et al., 2018;Naumova et al., 2013). Because both meiotic and mitotic cells undergo chromosome condensation, we asked what is the difference between meiosis (pacSC) and mitosis in terms of 3D chromatin structures. The P(s) curves showed that at a scale less than 1 Mb, mitotic and meiotic chromatin shared similar folding patterns with a P(s)$s slope at about À0.5 ( Figure 3D). However, they displayed a striking difference at longer distances, where the meiotic P(s) curve had a steep drop at $8 Mb, whereas the mitotic curve dropped at $30 Mb ( Figure 3D). This result suggested that meiotic pacSC chromatin is condensed via distinctive mechanisms as compared with that of mitotic chromatin.

Active Compartments for Genes and piRNA Clusters with Meiosis-Related Functions Are Enriched in pacSC
We further examined A/B compartmentalization in these spermatogenic stages in detail. The A compartments are defined as regions with a positive first eigenvector value when eigenvector decomposition is applied on observed/expected intra-chromatin interaction matrices (50 kb bin size in this study). These computed A compartments correspond to the active, euchromatic nuclear compartments. By contrast, the B compartments largely correspond to those inactive, heterochromatic genomic regions (Lieberman-Aiden et al., 2009). We investigated if there is A/B compartment switching between priSG-A and pacSC. We found A compartments accounted for 43.2% in priSG-A and 13.9% of them changed to B compartments in pacSC. Of the other 56.8% B compartments in priSG-A, 15.3% changed to A compartments in pacSC ( Figure 4A). Notably, the number of genome regions in the A (active) compartments was higher in pacSC, suggesting that pacSC chromatin is in a more transcriptionally active state.
We then identified the genes harbored by compartments that showed switching behaviors and examined their expression using RNA-seq data at the corresponding priSG-A and pacSC stages . Although there were genes that switched from A to B or B to A compartments with both up-and downregulation of their expression, we found that genes that were originally located in compartment B regions in priSG-A but switched to compartment A regions in pacSC, and at same time increased their expression, were the most abundant ( Figure 4B, 1,037 genes). This result implies that compartment switching from B to A is correlated with up-regulation of gene expression for their proper biological functions. Indeed, a number of these genes function in cilium formation and DNA DSBs repair ( Figure S3A), which are critical for normal sperm physiology and function ( Figure S3B). One example for both compartments switching and increased gene expression in pachytene cells is the DAZ family protein, Boll ( Figure 4C), which is required for the meiotic G2/M transition and germ cell development (VanGompel and Xu, 2010).
PIWI proteins and PIWI-interacting RNAs (piRNAs) direct the silencing of target nucleic acids in animal germline and soma (Ozata et al., 2019). In early pachytene, piRNAs are produced from large genomic loci referred to as piRNA clusters; these RNAs function in controlling transposons and have been shown to contribute to paternal (B) Heatmaps of log 2 ratio comparisons of the observed/expected interaction frequencies (100 kb bins, chromosome 1) for pacSC versus priSG-A, SZ versus pacSC, SZ versus priSG-A cells, respectively. (C) The P(s) curves (relationship between interaction probability and genomic distance) for each tested spermatogenic stage. (D) The P(s) curves (relationship between interaction probability and genomic distance) for pacSC and mitosis, and slopes (k) À1 (blue) and À0.5 (orange) were shown by the dotted lines.
iScience 23, 101034, April 24, 2020 imprinting (Aravin et al., 2007;Watanabe et al., 2011) and are required for the subsequent meiosis process (Fu and Wang, 2014). In addition to identifying genes in the compartments switched regions, we also analyzed the expression of piRNA clusters that locate to compartments showing stage-specific switching (piRNA clusters were identified from a published database, Rosenkranz, 2016). We found that 79.9% of all piRNA clusters had compartment A status in priSG-A, of which 88% remained in A compartments in pacSC; in contrast, 80% of the piRNA clusters with compartment B status in priSG-A switched to compartment A status in pacSC, correlating with their increased expression changes ( Figure 4D).
To explore the difference between piRNA clusters that had compartment-A status in both priSG-A and pacSC cells (not switched) and that switched compartments from B to A between the priSG-A and pacSC stages, we examined the conservation status of these two categories of piRNA clusters using phastCons score, which measures the probability that each nucleotide belongs to a conserved element (Siepel et al., 2005). We averaged phastCons scores of every nucleotide in each piRNA cluster and compared the average scores of not switched and switched piRNA clusters between the priSG-A and pacSC stages. The former (not switched) was more conserved than the switched between these two stages in a phylogenetic analysis of 60 vertebrate species ( Figure 4E). This result indicated that the piRNA clusters with compartment A status in both priSG-A and pacSC cells exhibit signs of selection pressure and are conserved among different species, suggesting functional constraints on them and thus their potential functions in these species. These data suggested that the functional activity of piRNAs might be strongly impacted by their locations in compartment A or B in the 3D genome.

Chromatin Accessibility and CTCF/Cohesin Binding Are Retained throughout Mouse Spermatogenesis
Since TADs dramatically changed throughout spermatogenesis and disappeared in pacSC cells, we asked whether the 3D architecture of pacSC cells correlates with other genomic/epigenomic features. We conducted ATAC-seq to measure chromatin accessibility (Buenrostro et al., 2015) in four spermatogenic cell types (priSG-A, SG-A, pacSC, and rST) ( Figure S4A). We first found that the accessible chromatin regions were mostly located at gene promoters or enhancers and the intensity of the accessible chromatin signal decreased as spermatogenesis proceeds ( Figures S4B and S4C). When we compared the chromatin accessibility in priSG-A and pacSC, we found minimal differences between them ( Figure 5A). This similarity of chromatin accessibility was rather surprising because the 3D genome architecture in these two cell types was largely different ( Figure 5A). These data suggested that 3D genome reorganization might not be functionally linked to chromatin opening, at least at the level of TADs.
Chromatin accessibility is closely related to transcriptional activity. We thus examined transcriptional activity using ChIP-seq with antibodies against the elongating form of RNA Polymerase II (Pol II S2P). The enrichment of Pol II S2P on gene transcriptional start sites (TSSs) showed that these genes were actively iScience 23, 101034, April 24, 2020 transcribing in both pacSC and rST ( Figure S5A). However, the chromatin accessibility showed distinct patterns in pacSC and rST samples ( Figure S5B). The chromatin accessibility in these two cell types was also very different from patterns in the other two cell types: priSG-A and SG-A ( Figure S5B). Specifically, between pacSC and priSG-A, we identified 6,580 (out of 15,815) differentially accessible chromatin regions, whereas between pacSC and rST, there were 1,288 (of 15,815) differentially accessible regions. The chromatin regions displaying altered accessibilities between priSG-A and pacSC were analyzed for GO terms, and we found that they enriched for genes related to meiosis, and especially to the key events in meiosis I, e.g., meiotic recombination or synapsis ( Figure S5C). A similar comparison of differentially accessible chromatin regions between pacSC and rST also revealed enrichment of meiotic genes, but notably, there was no longer any enrichment for gene groups specific to meiotic recombination or synapsis ( Figure S5D). These results provided a molecular explanation of the biological knowledge that pacSC are at the prophase of the first meiotic division and that rST are produced after meiosis completion. Moreover, these data highlighted the utility of assessing chromatin accessibility states for inferring the functional status of transcriptional activity.
We then further investigated whether the differential chromatin accessibility may provide a molecular basis to the observed disappearance of TADs at the pacSC stage. As CTCF is enriched at TAD boundaries and the depletion of CTCF results in loss of TADs (Dixon et al., 2012;Nora et al., 2017), we checked if the chromatin accessibility of CTCF motifs may be affected at the TAD boundaries. Intriguingly, the chromatin accessibilities at CTCF motifs appeared highly similar between priSG-A and pacSC ( Figure 5B), indicating that CTCF could bind to chromatin, which does not form conspicuous TADs (i.e., at the pacSC and rST stages). We experimentally examined CTCF binding in pacSC and rST cells via CTCF ChIP-seq. We found that, at both the pacSC stage and the rST stage, CTCF was still enriched at the TAD boundaries that were defined from the priSG-A stage ( Figure 5C). In addition to CTCF, we also checked the binding of cohesin in pacSC and rST by Rad21 ChIP-seq. The result showed that in both pacSC and rST, Rad21 was also enriched at the TAD boundaries that were defined from the priSG-A stage ( Figure 5C). These results suggested that the TAD disappearance in the meiotic stages of spermatogenesis could not be explained by the loss of CTCF or cohesin binding. Our data also renewed the understanding of CTCF and cohesin function in 3D genome structure, i.e., although CTCF and cohesin may be the molecules required for TAD maintenance in mammals (Nora et al., 2017;Rao et al., 2017), our data suggest that CTCF or cohesin binding at TAD boundaries is insufficient to maintain TADs at the pacSC stage.
During mammalian spermiogenesis, a highly condensed chromatin structure of sperm nuclei is in part formed by histone-protamine exchange and 1%-8% of histones in mouse are retained in sperm chromatin (Jung et al., 2019). Sperm histones were found to be mainly retained in distal intergenic regions, as shown by ChIP-seq using nucleoplasmin (NPM)-treated sperm . However, whether TAD boundaries play a role in retaining these histones in the sperm is unknown. We analyzed published protamine 1 targeted ChIP-seq data from sperms  and averaged the ChIP-seq signals around TAD boundaries in the sperm. Although protamine 1 was shown to occupy TSS regions, there was no enrichment of protamine 1 at the TAD boundaries ( Figure S5E). We also analyzed the published histone H3 ChIP-seq data  and found no enriched localization of histone H3 around TAD boundaries, either ( Figure S5F). These results implied that the TAD structure in the sperm may not be related to the histone or protamine constitution. Another question we asked was whether sperm-retained histones play substantial roles in the next generation. We found that, although histone H3 was depleted from the TSS regions, 171 genes had at least one exon in histone H3 retaining regions, and 16 of these 171 genes had been previously reported to be expressed in the early embryonic stage (Theiler stage 1-5, Table S2) (Smith et al., 2019). Further work examining the distribution of different modifications of retained histones and their potential roles in subsequent embryogenesis will provide more accurate answers to the function of retained histones.
Meiotic DSB Sites Are Pre-opened at the priSG-A Stage and piRNA Clusters Are Uniquely Open at the pacSC Stage DNA double-strand break (DSB) is essential for the recombination in meiosis, and it begins at the leptotene stage and ends at the zygotene stage (Baudat et al., 2013;. We examined the accessibility of chromatin regions surrounding previously reported DSB hotspots midpoints by analyzing our ATAC-seq data. We first identified meiotic DSB hotspots midpoints by the published Dmc1 ChIP-seq data obtained from DSB-positive cells in Hop2 knockout mice, whose meiotic DSBs are not repaired and meiotic progression is blocked at the pachytene-like stage . We then compared the chromatin accessibility around these meiotic DSB hotspots midpoints at the priSG-A and pacSC stages. Remarkably, this analysis revealed that the numbers of the accessible DSB hotspots midpoints were largely unchanged between these two stages ( Figure 5D). These results raised an intriguing possibility that meiotic DSB hotspots might be pre-opened (i.e., "primed") prior to their starting point, perhaps at a time even earlier than the priSG-A stage, and remained open after their repairing at the pacSC stage.
Since we observed that there were more piRNA clusters located in active compartments at the pacSC stage than at the priSG-A stage ( Figure 4D), we further compared the numbers of accessible pacSC piRNA clusters between priSG-A and pacSC. A larger number of accessible piRNA clusters were found in pacSC but not in priSG-A ( Figure 5E). These data were consistent with the above observation that more piRNA clusters were located in active compartments at the pacSC stage compared with that at the priSG-A stage, and it was also in line with the fact that pacSC-specific piRNAs are transcribed at the pacSC stage. These data together supported a refined coordinated functional contribution of DSBs and piRNA clusters to the meiosis process of murine spermatogenesis.

Reorganized Chromatin Loops in Sperm Contribute to Early Embryo Development
To examine the finer-scale architecture of the chromatin, we identified chromatin loops from the Hi-C data in each spermatogenic stage using cLoops (Cao et al., 2020). We observed that these loop interactions were reprogrammed throughout the spermatogenesis process ( Figures 6A and 6B), similar to the dynamics at TAD or compartment strength levels. The number of loops was significantly decreased at the pacSC and rST stages, followed by an increase to a maximum value in SZ ( Figure 6C). We aggregated the interaction values around loops defined in priSG-A and found that the strength of loops, as shown by the APA values , showed the same trend as the loop numbers ( Figure S6A).
To test the functional roles of the most abundant loops formed in SZ, we used ChIP-seq data of histone marks and RNA Pol II binding in testis from previous reports (Shen et al., 2012) to annotate potential enhancers. This analysis revealed that 1,572 genes were putatively involved in the formation of promoterenhancer loops at the SZ stage ( Figure 6C). Of these 1,572 genes, only a small percentage (317, or 20%) were also anchored in loops in priSG-A samples ( Figure S6B), suggesting that, although chromatin architecture in SZ samples restored the pattern in priSG-A at the larger scale (of compartments and TADs), at the smaller scale, the reestablished loops in SZ samples were distinct from those in the priSG-A samples. Genes that were involved in loops in SZ might tend to be transcriptionally active in subsequent embryonic stages, and previous studies have demonstrated the expression of a majority (1,174, or 75%) of these 1,572 genes in embryos, with 217 genes specifically detected in the developmental window of Theiler stages 1-5, although many of them repeatedly appeared in different stages of these early stages (Smith et al., 2019) ( Figure 6D). Specifically, 63 genes could be detected from Theiler stage 2 and 26 genes could be detected from Theiler stage 1 until Theiler stage 4 ( Figure S6C). These 217 genes were found to be involved in biological processes like regulation of gene expression and cell differentiation ( Figure S6D), and their reported knockout phenotypes in mammals include embryonic growth abnormality, arrest, and lethality (Figure 6E). As examples, Actr3Gt ðA009F03ÞVauA =Actr3Gt ðA009F03ÞVauA mice showed embryonic growth arrest (Vauti et al., 2007) and Adar tm2Knk =Adar + mice showed abnormal embryonic erythropoiesis and abnormal liver development (Wang et al., 2000). These results implied that the reorganized loops we detected in SZ might facilitate the expression of genes required for cell differentiation during the subsequent embryogenesis.

DISCUSSION
Recent studies have shown the dynamic changes in TADs and chromatin compartmentalization during mammalian spermatogenesis (Alavattam et al., 2019;Patel et al., 2019;Vara et al., 2019;Wang et al., 2019). Here, we further explored dynamic higher-order chromatin structures and their functional impacts on transcriptional regulation during mouse spermatogenesis combining Hi-C, ATAC-seq, ChIP-seq, and RNA-seq data. Our results complement recent reports by demonstrating that chromatin loops are dynamically reorganized during spermatogenesis in a pattern similar to that of TADs: they are present at stages prior to pacSC, are ''erased'' in pacSC, and get re-established after it. We also showed that the switching of A/B compartments between spermatogenic stages is tightly associated with meiosis-specific mRNAs and piRNAs expression. For mechanism, we showed that chromatin accessibility and the binding of CTCF and cohesin at TAD boundaries are not responsible for the TAD and loop diminishment at the pachytene stage.
It was reported that A/B compartment switching is correlated with the expression of specific genes during spermatogenesis (Vara et al., 2019). We also found some meiosis-related genes switch A/B compartment states and change expression between priSG-A and pacSC. In addition, we focused on piRNA clusters and demonstrated that most piRNA clusters are located in active A compartments and their chromatin states are quite stable, yet a majority of the piRNA clusters that are located in inactive B compartments in priSG-A are switched to active A compartments in pacSC. This is in line with the fact that pacSC cells are known to transcribe a characteristic set of pacSC-specific piRNAs (Aravin et al., 2006;Girard et al., 2006;Ozata et al., 2019); thus, our findings indicate that the chromatin compartment switch likely affects the activity of piRNAs.
Our results showed that, although the 3D genome organization changes greatly at the pacSC stage, chromatin at regions containing meiotic genes are still accessible, indicating that chromatin accessibility per se is not functionally associated with TADs or chromatin loops at the pacSC stage. In addition, enrichment of meiotic DSB sites in active A compartments has been reported during the pachytene stage (Patel et al., 2019). Chromatin accessibility was suggested to correlate with the occurrence of meiotic DSBs by some authors (Patel et al., 2019), but others have shown that accessible chromatin regions in pacSC have limited overlap with DSB sites (Maezawa et al., 2018). In the present study, we showed that many meiotic DSB sites are in accessible chromatin regions, and our data unexpectedly revealed that these sites are already "primed" in the priSG-A stage, which may play a functional role in generating DSB sites.
Chromatin tracing by super-resolution imaging at the single-cell level revealed TAD-like structures and showed that TAD boundaries tend to be located near CTCF and cohesin binding sites (Bintu et al., 2018). These TAD-like structures were still found in individual cells after removal of cohesin from chromatin, although the signal for the TAD-like structures was buried when the cells were analyzed as a population (Bintu et al., 2018). Upon specific chemical degradation of cohesin, the boundaries of the TAD-like structures in single cells appeared to be randomly distributed, suggesting that cohesin might not be required for TAD maintenance but may rather function to establish developmentally appropriate TAD boundaries (Bintu et al., 2018). Moreover, cohesin-mediated transcription in genomic regions out of the chromosomal axes in primary spermatocytes can provide an environment conducive to both gene expression and the formation of DSBs (Vara et al., 2019). We speculate that the population-averaged signal for TADs in the CTCF-depleted cells of the previous report (Nora et al., 2017) might also mask the occupancy of TAD structures in individual cells. However, this is not the case for meiotic chromatin: our ChIP-seq data showed that specific binding of CTCF and cohesin (at TAD boundaries that were defined at the priSG-A stage) is retained at the pacSC stage when TADs disappear, suggesting that there are additional factors that regulate the formation or maintenance of TAD structures, at least during meiosis. Multiple methodological strategies that combine 3D genome detection and immunoprecipitation or chromatin isolation and mass spectrometry and experiments better control for variables could be employed to identify candidates of such factors. In addition, loss of TADs and CTCF binding was reported to occur in mitotic chromatin (Gibcus et al., 2018;Naumova et al., 2013;Oomen et al., 2019), which, when taken together with our findings, strongly implies that TADs maintenance may be controlled via distinct mechanisms in mitosis versus meiosis.
In addition to lost TADs and chromatin loops, the global chromatin organization between mitotic and meiotic chromatin is also similar when compared with that in interphase, as revealed by our P(s) curve analysis. Nevertheless, meiotic chromatin exhibits A/B compartments ( Figure 1C), whereas mitotic chromatin does not (Gibcus et al., 2018). Thus, a comparison of factors modulating 3D genome organization at the compartmentalization level in meiotic and mitotic chromatin may discover unknown factors involved in higher-order chromatin organization. In addition to CTCF binding, transcriptional activity represents another major difference between meiosis and mitosis and is known to be affected by higher-order chromatin organization. Although both meiotic and mitotic chromatin are compacted, meiotic chromatin remains transcriptionally active ( Figure S5A), whereas mitotic chromatin is largely transcriptionally inert (Gottesfeld and Forbes, 1997;Wang et al., 2019). We propose that the distinct chromatin interaction patterns may be one of the causes underlying this difference.
Previous studies have reported that TADs and chromatin loops are associated with gene expression (Bonev and Cavalli, 2016;Flavahan et al., 2016;Hnisz et al., 2016;Lupiá ñ ez et al., 2015), whereas there are also some findings that show that loss of TADs and chromatin loops by rapid degradation of cohesin or CTCF have only modest effects on gene expression (Nora et al., 2017;Rao et al., 2017). Thus, the impacts of chromatin organization at TAD and chromatin loop scales on gene expression are still not very clear. A recent study by Micro-C that explores the chromatin organization from single nucleosomes to whole genome-scale found that there are microTADs and gene-level folding within TADs that form in a gene-dependent manner and are highly associated with chromatin accessibility and transcriptionally active chromatin (Hsieh et al., 2020). Here, our results together with several recent reports show that TADs and chromatin loops are almost lost at the pacSC stage during mouse spermatogenesis (Patel et al., 2019;Vara et al., 2019;Wang et al., 2019), whereas we also show that CTCF/cohesin binding on chromatin, chromatin accessibility at promoter/ enhancer regions, and Pol II binding remain at the pacSC stage. Thus, although TADs and chromatin loops are largely lost at the pacSC stage, the microTADs and gene-level folding, which cannot be detected by in situ Hi-C, might remain on pacSC chromatin and contribute to transcription regulation and gene expression at the pacSC stage. Our chromatin accessibility and Pol II S2P ChIP-seq data support this possibility, and future studies with nucleosome-resolution Micro-C and single-cell methods will further deepen our understanding of the relationship between chromatin organization and transcription regulation during mouse spermatogenesis. Altogether, beyond confirming several major conclusions about TAD reorganization and the dynamics of A/B compartments strength, our study also reveals mechanistic insights about the functional relationships between 3D genome organization, CTCF or cohesin binding, and TAD maintenance and demonstrates the roles of A/B compartments switching in controlling meiotic gene transcription programs as well as potential function of chromatin loops in SZ in preparing for subsequent embryogenesis.

Limitations of the Study
In this study, we explored the 3D genome organization during mouse spermatogenesis with isolated spermatogenic cells from five different stages. We showed that chromatin loops and TADs are reorganized during mouse spermatogenesis, which almost disappear at the pacSC stage, whereas the transcriptional activity, CTCF binding on chromatin, and chromatin accessibility of promoters and enhancers are largely maintained at the pacSC stage. Thus, it is unclear how the transcriptional regulation is going on at the pacSC stage when the pacSC chromatin has no chromatin loops and TADs. A recent study about transcription-linked mammalian 3D chromatin folding with Micro-C uncovered a finer-scale genome organization that is associated with transcription and chromatin accessibility more tightly, compared with chromatin loops and TADs (Hsieh et al., 2020). Future studies that focus on the finer-scale genome organization at the pacSC stage might reveal the relationship between 3D genome structures and transcriptional regulation, as whether the finer-scale genome organization remains at the pacSC stage while chromatin loops and TADs are almost lost is still unclear.
In addition, we and others only explored the 3D genome organization during mouse and rhesus monkey spermatogenesis with limited isolated cell types (Alavattam et al., 2019;Patel et al., 2019;Vara et al., 2019;Wang et al., 2019), which restrict us to get a more comprehensive understanding about the 3D genome organization during mouse spermatogenesis. Future studies with single-cell Hi-C technologies and related methods might give us a more comprehensive view about the 3D genome organization during mouse spermatogenesis and uncover the underlying mechanisms about the reorganized chromatin structures.

METHODS
All methods can be found in the accompanying Transparent Methods supplemental file.

SUPPLEMENTAL INFORMATION
Supplemental Information can be found online at https://doi.org/10.1016/j.isci.2020.101034. Aravin     ChIP-Seq as well as their input control were obtained from  and .

Purification of male germ cells during spermatogenesis
C57BL/6 mice were originally purchased from Vital River Laboratories in Beijing,

Immunofluorescence
Immunofluorescence analysis was performed to identify the cell type and purity of the isolated cells. Since the isolation methods were routine for us with high purity and success (Gan et al., 2013), only part of the isolated cell types was used to perform Immunofluorescence. Specifically, the cell pellet was

In situ Hi-C
In situ Hi-C was used to explore the higher order chromatin organization with isolated cells during mouse spermatogenesis . Briefly, about 5 x 10 5 to 5 x 10 6 cells were crosslinked using formaldehyde ( sequenced on the HiSeq X ten platform.

ATAC-Seq
ATAC-Seq was performed following a previously published method (Buenrostro et al., 2013). Briefly, 5 × 10 4 cells were spun at 500 g for 5 min at 4 °C, and the cell pellet was washed with 50 μl of cold 1× PBS and collected by centrifugation at 500 g for 5 min at 4 °C.

Native ChIP-Seq
Native ChIP was performed following a previously published report (Gilfillan et al., 2012). Briefly, the cells were washed once with 1x PBS and then pelleted at 800 g for 5 min.

Hi-C: alignment
Hi-C dataset from cells at the SZ stage was obtained from published work (Jung et al., 2017). Paired-end .fastq files of Hi-C libraries were mapped using HiC-Pro (v 2.8.1) (Servant et al., 2015). Matrix binning and balancing was performed with the cooler software package (0.8.6.post0). In brief, paired reads were iteratively mapped to the Mus musculus mm10 genome via bowtie2 (v 2.2.5) (Langmead and Salzberg, 2012); the unmapped reads contain the MboI ligation sites were then trimmed and aligned back to the mm10 genome again. After combing the two-round mapping results, uniquely aligned paired reads were kept and assigned to MboI restriction fragments. The paired reads with dangling ends and or predicted self-circle structures pairs were excluded from valid pairs. This alignment and filtering process resulted in sequencing reads described in (Table S1). We merged replicates and sampled all valid pairs of each sample to the same 172,252,595 reads, which were further converted to Hi-C contact matrices binned into 100kb, 50kb, 25kb, and 5kb sizes respectively (Abdennur and Mirny, 2019). The binned contact matrices were further normalized by using an iterative correction method (Servant et al., 2015). The Hi-C matrices were further transferred into the .hic format to facilitate visualization with juicebox (Durand et al., 2016).

Hi-C: interaction frequency curves
The interaction frequency curves (P(s)) were calculated with normalized interaction matrices in 100kb resolution following methods used in (Naumova et al., 2013). We first divided distances into logarithmically spaced bins and counted the number of interactions at corresponding distances at each bin. To obtain the probability P(s), we divided the number of interactions in each bin by the total number of possible region pairs.

Hi-C: A/B Compartment, TADs, and chromatin loops analysis
The compartment profiles of Hi-C matrices at 50kb resolution were identified using cooltools (https://github.com/mirnylab/cooltools) eigdecomp.py. To do this, we first calculate the intra-chromatin observed/expected interaction matrix. We used the interaction matrix at 500kb resolution and calculated the mean number of interactions at a given genomic distance, d, using a sliding window with linearly increasing size, and get the expected contact probability. Then the observed/expected matrix was used to generate a Pearson correlation matrix, which was subject to principal component analysis (PCA). The eigenvector of the first principal component defines the A/B compartment profile (Lieberman-Aiden et al., 2009). With convention, genomic compartments were assigned to 'A' (active/euchromatic compartments) and 'B' (inactive/heterochromatic compartments) primarily based on associated mRNA transcription.
Compartment saddle analysis was performed by using cooltools (Gibcus et al., 2018;Schwarzer et al., 2017). Briefly, we used the observed/expected Hi-C maps, obtained by dividing each diagonal of a matrix by its chromosome-wide average value from 50kb binned iteratively corrected interaction maps of cis interactions. In each observed/expected map, we rearranged the rows and the columns in the order of increasing eigenvector value. Finally, we aggregated the rows and the columns of the resulting matrix into 30 equally sized aggregated bins, thus obtaining a compartmentalization plot ("saddle plot"). The The topologically associating domains (TADs) were identified following an insulation square analysis (Crane et al., 2015). Briefly, to calculate the 'insulation' score of each bin in the 25kb binned Hi-C data, we sided a sequential size of squares (from 100kb to 1000kb, 50kb as steps) along the diagonal of the interaction matrix for each chromosome, and aggregated all signals within the square. A 125kb window was used for calculation of the delta vector. Only boundaries for which more than half of the differentially sized squares had a 'boundary strength' value larger than 0.25 were retained.
Insulation scores were plotted around all priSG-A TADs as well as their nearby regions (+/-0,5 TAD length). Using all of the priSG-A TADs as control TADs, we then plotted the composite interaction by averaging interactions near TAD boundaries (range) for each of all the sampled stages of spermatogenesis.
The chromatin loops were identified using cLoops (Cao et al., 2018) with parameters -w -j -s -m 3. The APA (Aggregation Peak Analysis) was performed following the methods developed by . Briefly, to measure the enrichment of loops over the local background, the KR normalized contact frequency of pixels of loops, as well as the surrounding pixels up to 10 bins away in both x and y directions, i. e., 50kb*50kb local contact matrices, were collected and plotted. APA scores (P2LL, the ratio of the central pixel to the mean of the pixels in the lower-left corner, representing the strength of chromatin loops) were determined by dividing the center pixel value by the mean value of the 25 (5*5) pixels in the lower right section of the APA plot.
We downloaded lists of active promoters (marked by H3K4me3 or PoI II binding signals) and enhancers (marked by H3K4me1 or H3K27ac) in the sperm (Shen, et al., 2012). We used their genomic loci to annotate our loop anchors in each stage of spermatogenesis. We then queried the genes with promoter-enhancer loops to MGI-Mouse Gene Expression Database (Smith et al., 2018) and found genes with detectable expression in Theiler Stage 1-5.

ATAC-Seq data analysis
All ATAC-Seq sequencing data were mapped to the mm10 reference genome using snap-aligner (v 1.0) (Zaharia et al., 2011). Samtools was then used to further sorted, index, and remove duplicates from the aligned reads (Li et al., 2009). The ATAC-Seq peaks of all samples were called using MACS2 (v 2.1.0) with the default parameters . Pearson correlation coefficients were calculated based on the normalized reads number for merged peaks from all samples between different samples. The ATAC-Seq reads at identified peaks were quantified using bedtools (Quinlan and Hall, 2010) and normalized using the quantile norm function in R. The differential accessible chromatin regions were identified using DESeq2 . The functional enrichment analysis of the differentially accessible chromatin regions was performed using GREAT (McLean et al., 2010). The motif enrichment analysis of differentially accessible chromatin regions was performed using the HOMER findMotif function (Heinz et al., 2010). The list of meiotic DNA double strand break sites was generated from previously published DMC1 ChIP-Seq . The list of piRNA cluster sites from 14.5dpp mouse was obtained from a piRNA cluster database (Rosenkranz, 2015). phastCons scores were downloaded from ucsc genome browser