It Is Just a Matter of Time: Balancing Homologous Recombination and Non-homologous End Joining at the rDNA Locus During Meiosis

Ribosomal RNA genes (rDNAs) are located in large domains of hundreds of rDNA units organized in a head-to-tail manner. The proper and stable inheritance of rDNA clusters is of paramount importance for survival. Yet, these highly repetitive elements pose a potential risk to the genome since they can undergo non-allelic exchanges. Here, we review the current knowledge of the organization of the rDNA clusters in Arabidopsis thaliana and their stability during meiosis. Recent findings suggest that during meiosis, all rDNA loci are embedded within the nucleolus favoring non-homologous end joining (NHEJ) as a repair mechanism, while DNA repair via homologous recombination (HR) appears to be a rare event. We propose a model where (1) frequent meiotic NHEJ events generate abundant single nucleotide polymorphisms and insertions/deletions within the rDNA, resulting in a heterogeneous population of rDNA units and (2) rare HR events dynamically change rDNA unit numbers, only to be observed in large populations over many generations. Based on the latest efforts to delineate the entire rDNA sequence in A. thaliana, we discuss evidence supporting this model. The results compiled so far draw a surprising picture of rDNA sequence heterogeneity between individual units. Furthermore, rDNA cluster sizes have been recognized as relatively stable when observing less than 10 generations, yet emerged as major determinant of genome size variation between different A. thaliana ecotypes. The sequencing efforts also revealed that transcripts from the diverse rDNA units yield heterogenous ribosome populations with potential functional implications. These findings strongly motivate further research to understand the mechanisms that maintain the metastable state of rDNA loci.

Ribosomal RNA genes (rDNAs) are located in large domains of hundreds of rDNA units organized in a head-to-tail manner. The proper and stable inheritance of rDNA clusters is of paramount importance for survival. Yet, these highly repetitive elements pose a potential risk to the genome since they can undergo non-allelic exchanges. Here, we review the current knowledge of the organization of the rDNA clusters in Arabidopsis thaliana and their stability during meiosis. Recent findings suggest that during meiosis, all rDNA loci are embedded within the nucleolus favoring non-homologous end joining (NHEJ) as a repair mechanism, while DNA repair via homologous recombination (HR) appears to be a rare event. We propose a model where (1) frequent meiotic NHEJ events generate abundant single nucleotide polymorphisms and insertions/deletions within the rDNA, resulting in a heterogeneous population of rDNA units and (2) rare HR events dynamically change rDNA unit numbers, only to be observed in large populations over many generations. Based on the latest efforts to delineate the entire rDNA sequence in A. thaliana, we discuss evidence supporting this model. The results compiled so far draw a surprising picture of rDNA sequence heterogeneity between individual units. Furthermore, rDNA cluster sizes have been recognized as relatively stable when observing less than 10 generations, yet emerged as major determinant of genome size variation between different A. thaliana ecotypes. The sequencing efforts also revealed that transcripts from the diverse rDNA units yield heterogenous ribosome populations with potential functional implications. These findings strongly motivate further research to understand the mechanisms that maintain the metastable state of rDNA loci.

INTRODUCTION
The central importance of the rRNA genes for the biology of any organism is evident, as they are essential for survival and for all cellular processes. They are among the evolutionary oldest and also most highly transcribed genomic regions forming the RNA building blocks of ribosomes. Most eukaryotic genomes contain clusters with hundreds to thousands of rRNA gene copies arranged in tandem which are transcribed and processed within the nucleolus.
In eukaryotes, the 18S, 5.8S, and 25S rRNAs form the scaffold for the small and large ribosomal subunits and all three are encoded together in functional units and transcribed as a single polycistronic 45S precursor transcript by RNA polymerase I (Wallace and Birnstiel, 1966;Moss et al., 2007;Layat et al., 2012). The 45S rRNA gene units (also termed rDNA units) are arranged in a head-to-tail manner in large clusters known as nucleolus organizing regions (NORs; Ritossa and Spiegelman, 1965;Wallace and Birnstiel, 1966). In the Arabidopsis thaliana reference ecotype Col-0, the 45S rDNA units are approximately 10 kb long and arranged in two clusters, each with ~400 repeats, at the top of chromosomes 2 and 4 (Copenhaver and Pikaard, 1996;Sims et al., 2021). A further component of the large ribosomal subunit, the 5S rRNA, is located on chromosomes 3, 4, and 5 in the A. thaliana Col-0 ecotype, also arranged in clusters and transcribed by RNA polymerase III (Murata et al., 1997;Layat et al., 2012). One 45S rDNA unit is also found in proximity of the 5S rDNA located on chromosome 3 (Abou-Ellail et al., 2011). Although rRNA transcripts account for approximately 50% of all transcribed RNAs in a cell, only a fraction of the rRNA genes is transcribed at a given time (Warner, 1999;Grummt and Pikaard, 2003;Pontvianne et al., 2010Pontvianne et al., , 2012. Recent studies have shown that individual 45S and also 5S rDNA units are not identical. Instead, they display a substantial amount of variability, not only within the intergenic regions but also in the genic regions transcribing the conserved ribosomal RNA subunits (Chandrasekhara et al., 2016;Havlová et al., 2016;Rabanal et al., 2017b). These variants have been exploited as molecular markers to study rDNA cluster-specific expression (see below).
The high level of transcriptional activity leads to torsional stress in rDNA and requires the activity of topoisomerases to relieve the positive and negative torsions (French et al., 2011). During replication, highly transcribed regions of the genome, such as the rDNA loci, may encounter frequent collisions between transcription and replication machineries, which need to be resolved (Castel et al., 2014;García-Muse and Aguilera, 2016;Sims et al., 2021). Both processes mentioned above are sources of DNA damage and genome instability in general. The unique nature of the rDNA loci not only makes them especially vulnerable to various types of DNA damage, it also demands special attention during DNA repair. The highly repetitive rDNA loci, with their hundreds of nearly identical rDNA units arranged head-to-tail, may undergo dramatic re-arrangements during homologous recombination (HR) DNA repair. HR may ultimately lead to lengthening or shortening of the rDNA arrays and in general to copy number instability (Warmerdam et al., 2016). Furthermore, the presence of rDNA clusters on multiple chromosomes adds the additional risk of inter-chromosomal recombination. As outlined in more detail below, during meiosis, a developmental program essential for the recombination of genetic traits, numerous DNA doublestrand breaks (DSBs) are introduced. In this context, the rDNA loci are sequestered away from the canonical HR pathway and a different DNA repair pathway is employed, termed non-homologous end joining (NHEJ; Sims et al., 2019). NHEJ will less likely lead to genome re-arrangements and rDNA copy number loss, but may lead to the introduction of single nucleotide polymorphisms (SNPs) and short-range insertions/ deletions (InDels; Chang et al., 2017;Wright et al., 2018;Xu and Xu, 2020).
In this review, we summarize the recent findings concerning the stability of the rDNA loci and their inheritance from a perspective of meiosis. We also provide a model, in agreement with the current data, that defines HR and NHEJ as the major determinants of rDNA cluster size and rDNA unit sequence variability.

DNA DOUBLE-STRAND BREAK FORMATION AND REPAIR
DNA damage, if not appropriately repaired, leads to loss of genetic material, genome re-arrangements, and cell cycle arrest. One of the most deleterious DNA insults are DNA DSBs which can for instance be generated by genotoxic agents or molecular tools like homing endo-nucleases, TALE nucleases, and CRISPR/ Cas9 (Wu et al., 2014;Lopez et al., 2021) by endogenous processes like the re-establishment of collapsed replication forks or by a dedicated machinery during meiosis. As mentioned above, DSBs can be repaired by different repair pathways, the most prominent being NHEJ and HR. NHEJ has been found to be active during all cell cycle stages, whereas HR is the dominant repair pathway during S and G2 and is briefly introduced below.

NON-HOMOLOGOUS END JOINING
Non-homologous end joining was first described in mammals, where it is the predominant mechanism for DSB repair in non-cycling, somatic cells. It is differentiated in c-NHEJ (canonical) and a-NHEJ (alternative) pathways, the latter including microhomology-mediated end joining (MMEJ), all with the direct ligation of processed DNA ends as common denominator. Re-joining of blunt ends or ends with a few overlapping bases occurs without regard for preserving the sequence or context integrity (Hays, 2002;McVey and Lee, 2008;Chang et al., 2017). c-NHEJ is initiated with the recognition and the juxtaposition of the broken ends. In mammals, this step is promoted by the DNA-dependent protein kinase (DNA-PK), a complex composed of the KU heterodimer (Xrcc5/6) and Frontiers in Plant Science | www.frontiersin.org 3 October 2021 | Volume 12 | Article 773052 the kinase DNA-PK catalytic subunit (DNA-PKcs; Blackford and Jackson, 2017). It is important to note that DNA-PKcs have not been found to be encoded in plant genomes, indicating that in plants NHEJ is orchestrated differently (Templeton and Moorhead, 2005;Yoshiyama et al., 2013). The Artemis protein and the Xrcc4/DNA ligase IV heterodimer are subsequently recruited, with Artemis involved in the maturation of the DSB ends and the Xrcc4/DNA ligase IV complex catalyzing the resealing of the ends (Lees-Miller and Meek, 2003;Meek et al., 2004;Bleuyard et al., 2006). The KU heterodimer is composed of Ku70 and Ku80 and is involved in recognition, protection, and juxtaposition of the ends of a DSB. DNA-PKcs proteins are recruited to the DSB sites via interactions with the Ku/DNA complex and by phosphorylating various substrates (e.g.: Ku70, Ku80, Artemis, Xrcc4; Fell and Schild-Poulter, 2015). Artemis possesses both exo-and endo-nuclease activities and performs phospho-regulated maturation of the DSB ends as it cleaves DNA hairpins and other DNA structures (Lobrich and Jeggo, 2017). The final step, consisting of the ligation of broken ends, is carried out by the Xrcc4/DNA ligase IV heterodimer, which is recruited by DNA-PK. The MRN complex, composed of the proteins Mre11, Rad50, and Nbs1, stimulates this ligase activity in vitro and is also implicated in the juxtaposition of the ends of the break (Grawunder et al., 1997;Durdikova and Chovanec, 2017). The alternative NHEJ pathway MMEJ is promoted in the absence of c-NHEJ factors and involves the alignment of microhomologies at the DSB site (Seol et al., 2018). DNA ends are bound by PARP1 (potentially competing with Ku proteins; Wang et al., 2006;Cheng et al., 2011). Following DNA binding, PARP1 gets activated and poly-ADP-ribosylates itself and various targets in the vicinity leading to more accessible chromatin (Polo and Jackson, 2011;Beck et al., 2014). Subsequently, the MRE11-complex is recruited to process the DNA and prepares them for ligation via Ligase I or Ligase III.
The counterparts of most NHEJ proteins have been identified and characterized in plants (Bleuyard et al., 2006;Charbonnel et al., 2011). For instance, LIGASE4 is a well-conserved hallmark factor, also in plants, in the c-NHEJ DNA repair pathway (Friesner and Britt, 2003). MRE11 and its complex partners have also been identified and characterized in plants, and they together are required for both HR and MMEJ. Importantly, no homologs of DNA-PK and Ligase III and some further factors have been identified in plant genomes (Manova and Gruszka, 2015;Yoshiyama, 2016), highlighting some fundamental differences in the DNA damage response in plants and other organisms.

HOMOLOGOUS RECOMBINATION
In contrast, DNA DSB repair via the HR pathway preserves sequence integrity. Following DSB formation (see above), initiation of HR depends on the localization of the Mre11-Rad50-Nbs1/Xrs2 (MRN/X) complex and its partner CtIP/ Sae2/Com1 to the DSB sites (Wright et al., 2018). The MRN complex bridges the two ends, is involved in DNA end processing, and recruits further processing proteins (e.g., a 5' to 3' exo-nuclease). The nucleolytic activities yield a 3' ssDNA overhang, competent to invade dsDNA to probe for a homologous repair template. In addition, the MRN/X complex recruits the DNA damage kinase ATM/Tel1 which phosphorylates a large number of downstream targets (including Rad9, Rad17, Rad53, Rpa1, Xrs1 Com1/Sae2, and Exo1) involved in DNA repair and checkpoint control (Clerici et al., 2005;Roitinger et al., 2015). The ssDNA ends are coated with the replication protein A (RPA), thereby stimulating the recruitment of recombinases [in yeast via Rad52; in higher eukaryotes via BRCA2 (Krogh and Symington, 2004)]. The recombinase Rad51 (and in meiosis its relative Dmc1; see below) mediates subsequent strand invasion to probe for homologue sequences, assisted and stimulated by a battery of accessory proteins (Sung et al., 2003;Chan et al., 2019). In S/G2, the cell cycle stage during which HR is promoted, the sister chromatid and the chromatids of the homologous chromosome are available as repair templates. Following invasion and successful homology check, the invading strand is elongated and the displaced strand captured by the ssDNA overhang at the DSB site. Subsequently, the elongated strands are ligated yielding a double holiday junction (dHJ) that can lead, after resolution, to restoration of the original chromosome or to a cross-over and therefore a mutual exchange of chromosome arms. In case the sister chromatid has been used as a repair template, such an exchange is genetically neutral; in case a chromatid of the homologous chromosome has been used, such an exchange yields a chimeric chromosome. The latter is the desired repair product during meiosis to support meiotic chromosome disjunction and increase genetic diversity (Ohkura, 2015).
Alternatively, prior to second-end capture, the recombination intermediate can be dismantled by helicases and the invading, now elongated, strand anneals to the DSB site it originated from (also known as SDSA -synthesis dependent strand annealing). Subsequent DNA synthesis and ligation repairs the lesion, with the potential of some genetic information transfer (gene conversion) in case the template strand contained sequence polymorphisms, but without exchange of chromosome arms.
In this sense, in canonical non-repetitive regions of the genome, HR delivers a more faithful repair outcome with a high likelihood to re-establish the original DNA sequence, while NHEJ leads mostly to short-range deletions and to some extent to insertions and SNPs (Betermier et al., 2014;Liu and Huang, 2014;Ceccaldi et al., 2016).

MEIOSIS
Meiosis is a specific developmental process required for the formation of gametes, carrying the genetic information for the next generation. Meiosis is characterized by two consecutive cell divisions that reduce the genome size by half and by recombination of the paternal and maternal genomes. Novel allelic combinations are created by the mutual exchange of genetic information between parental chromosomes. This depends on meiotic DNA DSBs which are enzymatically induced by the conserved SPO11 protein (together with less conserved partners; Hunter, 2015;Mercier et al., 2015;Robert et al., 2016). About 250-300 DSBs are introduced in each individual meiocyte in Arabidopsis (Edlinger et al., 2011), and they all have to be reliably repaired for successful completion of meiosis. As mentioned, meiotic DSBs are introduced following DNA replication; therefore, cells are in G2-phase with HR being the predominant DNA repair pathway. Meiotic HR is specifically tuned to generate genetic diversity, preferentially using a chromatid of the homologue, and not the sister chromatid, as a repair template [inter-homolog (IH) bias]. Multiple such events along a chromosome ensure that homologous chromosomes recognize each other. At least one IH interaction per chromosome pair has to mature into a cross-over to ensure correct segregation of homologs during the first meiotic division (Gray and Cohen, 2016). In non-repetitive regions, recognition of the homologous partners works very reliably and non-allelic recombination events are not observed. This process is also aided by a meiosis-specific chromosome organization ("bouquet"), clustering telomers (and often also centromeres) to reduce the search space for the ssDNA nucleoprotein filaments (Harper et al., 2004). Genomic loci that are comprised of repetitive sequences, like the rDNA clusters, create a liability during recombination since they can undergo non-allelic exchanges and are a potential source of deletions, duplications, inversions, or translocations (Sasaki et al., 2010).

DSB FORMATION AND REPAIR AT THE rDNA LOCUS
Most of the studies concerning DSB repair at the rDNA region involve the use of induced DSBs by exogenous factors or the use of mutants that perturb the stability of the rDNA (Harding et al., 2015;Sluis et al., 2015;Warmerdam et al., 2016). In plants, a recent study employed CRISPR-Cas9 to induce DSBs at the rDNA locus. This led to a large population of plants each containing a varying number of rDNA repeats ranging from about 20 to 200% of the wild-type copy number (Lopez et al., 2021). While these plants represent a powerful resource to study rDNA dynamics in the future, the actual response to the Cas9-mediated DNA lesions has not yet been studied. In mammalian cells, it has been established that the DNA damage response at the rDNA and within the nucleolus depends on a critical threshold: low levels of DSB formation activate NHEJ, excessive DSB formation within the rDNA is repaired via HR, concomitant with transcriptional downregulation and nucleolus re-organization (van Sluis and McStay, 2017).
Studying rDNA repair in a meiotic environment is advantageous since a relatively defined number of endogenous DSBs are formed in a tightly regulated fashion. This allows monitoring DSB repair at the rDNA loci under physiological conditions (Sims et al., 2019). In plants, only a handful of factors are known to be involved in the repair process and stability of the rDNA in somatic and meiotic tissues after DSB formation. The RECQ/TOP3/RMI1 complex partner RMI2, the DNA helicases RTEL1, and FANCJ have been shown to be independently needed for maintaining the stability of the 45S rDNA loci in somatic tissues of some plants (Rohrig et al., 2016;Dorn et al., 2019). Furthermore, several additional studies have shown the importance of the chromatin assembly complex CAF-1in preventing DSB formation at the rDNA loci and maintaining rDNA copy numbers (Mozgová et al., 2010;Pavlistova et al., 2016). In addition, low amounts of 45S rDNA copies have shown to promote genomic instability in a genome-wide manner by generating large genomic re-arrangements (Picart-Picolo et al., 2020;Lopez et al., 2021). In meiosis, c/a-NHEJ factors, such as LIG4 and MRE11, have been shown to be important for DNA repair within the rDNA region, whereas HDA6 and NUC2, which are involved in regulating rDNA transcription and nucleolus integrity, are essential for limiting HR at the rDNA (Sims et al., 2019).

A BALANCE BETWEEN HR AND NHEJ
Studies in human cells, employing artificially induced DSBs, have described a re-organization of the nucleolus and a shift from NHEJ repair to HR upon reaching a certain threshold of DNA damage (van Sluis and McStay, 2017). This is concomitant with the formation of the nucleolar caps (Reynolds et al., 1964) and a shutdown of rRNA transcription while breaks in the rDNA persist. Nucleolar caps have not yet been described in other organisms other than humans and mice. In yeast, sites of DSBs within the rDNA re-localize to an extra nucleolar site for repair by HR (Horigome et al., 2019).
Work performed in A. thaliana shows that in physiological conditions, such as meiosis, the DNA lesions in the rDNA are preferentially repaired by NHEJ. The nucleolus creates a HR-refractory zone with strongly reduced numbers of HR events at the NORs (Sims et al., 2019). It is anticipated that sporadic events of HR can still occur, and they may leave noticeable traces, like rDNA unit duplication/ amplification/loss (copy number variation) and variable numbers of sequence repeats within the 45S rDNA units. Maintaining this unique HR-refractory domain depends on specific chromatin modifications which are distinct from the meiotic nucleus (Sims et al., 2019). In general, an increase in HR at the rDNA locus leads to the loss of units and reduced cell fitness. This is, for instance, well described in the FAS1 mutant background, in which an HR-dependent shortening of the NORs has been reported (Mozgová et al., 2010;Muchova et al., 2015). Nevertheless, HR events likely Frontiers in Plant Science | www.frontiersin.org 5 October 2021 | Volume 12 | Article 773052 also occur in wild-type plants, since a dramatic divergence in rDNA copy numbers is detectable within tens of generations, suggesting the presence of low recurring HR events that lead to a change in NOR length (Rabanal et al., 2017a). Furthermore, the presence of long segments of identical rDNA units is an indication of homogenization events at the rDNA locus, likely mediated by HR repair (Copenhaver and Pikaard, 1996;Sims et al., 2021). We suggest that DSBs within the rDNA, occurring at physiological levels (e.g., as a results of transcription/replication collisions or generated during the meiotic program), are repaired via NHEJ, preserving rDNA unit numbers and unit-internal repeat structures, but at the cost of producing errors. Currently, it is unclear whether there is a preference toward canonical or alternative NHEJ pathways (Sims et al., 2019). However, there are indications that both pathways are necessary for repairing lesions within the rDNA since LIG4 and MRE11 have an equal impact on rDNA stability. It is important to mention that DSBs generated within the rDNA by transcription/ replication conflicts or during meiosis are still rare events (Sims et al., 2019).
In general, the errors produced by the NHEJ pathways have the potential to generate sequence diversity between the rDNA units, and one would expect for them to accumulate at the transcriptional start and termination sites of each unit, due to selection against mutations in the portions of the rDNA that yield rRNA integrated into ribosomes. In fact, the highest number of SNPs and InDels is found in the external transcribed sequences (ETS), particularly close to the promoter and terminator regions (Chandrasekhara et al., 2016;Rabanal et al., 2017b). In contrast, very few SNPs/InDels are found within the portions transcribing the ribosomal RNA subunits (18S, 5.8S and 25S). This correlates well with the suggested high levels of transcriptional stress in the rDNA, with purifying selection acting on the regions transcribing rRNA subunits and with DNA lesions being repaired via NHEJ.
Imbalanced accumulation of polymorphisms is also apparent between the two NORs of A. thaliana. A recent study combining long-and short-read sequencing technologies to define the nucleotide composition and organization of 405 individual rDNA units of NOR2 of ecotype Col-0 identified less SNPs/InDels on the transcriptionally less active NOR2, than on NOR4 (Sims et al., 2021). To display the sequence diversity of these 405 rDNA units, we utilized their data and generated a phylogenetic network. For the analysis, we excluded the highly repetitive region of the SalI repeat boxes from each unit. The TCS network was inferred (Supplementary Figure S1) using the integrated method of the TCS approach (Templeton et al., 1992;Clement et al., 2002), which is based on the concept of statistical parsimony in PopArt (Leigh and Bryant, 2015). The network shows a lack of phylogenetic structure in the data indicating that a lot of parallel and reverse mutations obscure the relations between the units and that the conservative nature of the rDNA units in general may possibly mask local phylogenetic information. To address this latter point, we repeated the TCS analysis for short stretches of rDNA units (represented in the 59 BACs as published in Sims et al., 2021). Indeed, the majority of the BACs show a clear tree-like structure, with only very little reticulation. Thus, locally, the evolutionary process follows a classical tree-like pattern. Moreover, directly adjacent units on a BAC tend to be next to each other in the tree (data and visualization available upon request). The contigs identified in (Sims et al., 2021) provide additional evidence of tree-like evolution of the NOR2 region (Figure 1). Though the tree-like relation breaks down, the more rDNA units are analyzed due to multiple identical units occurring along the NOR2 region.
A plausible explanation for the higher abundance of SNPs/ InDels on NOR4 could be derived from the fact that in the ecotype Col-0, NOR4 is transcriptionally active in all analyzed tissues, while NOR2 is selectively silenced during development, and it is only active in certain tissue types (Chandrasekhara et al., 2016;Rabanal et al., 2017a). Transcriptional stress per se is a prime source of DSBs, and the rDNA is considered a hotspot of transcription and replication stress (Takeuchi et al., 2003). Since the pattern of NOR expression varies greatly among Arabidopsis ecotypes with some expressing predominantly one and some the other NOR (and some both), it would be interesting to analyze whether rDNA polymorphisms are positively correlated with transcriptional activity in different ecotypes.
It is interesting to speculate that the nucleolus represents an HR-refractory sub-compartment within the nucleus during meiosis (and after pre-meiotic DNA replication). As stated above, both NORs are transcriptionally activated in order to be recruited to the nucleolus and embedded in its HR-refractory zone (Sims et al., 2019). Perturbing the rDNA transcriptional activity or the nucleolar architecture generates an imbalance in the rDNA protective mechanism. In this sense, rDNA transcriptional activation, and subsequent recruitment into the nucleolus, could be a key regulatory mechanism to determine the mode of rDNA repair after DNA damage. The recruitment into the nucleolus following transcription is a conserved feature of rDNA (Pontvianne et al., 2013;Sims et al., 2019).
The protective mechanisms surrounding the 45 rDNA regions could not be limited to the nucleolus itself, since in certain tissues, the majority of 45S rDNA genes are not transcribed and excluded from the nucleolus. Inactive NOR4 rDNA genes are generally located at the nucleolar periphery, whereas NOR2 rDNA genes are completely excluded from the nucleolus area (Pontvianne et al., 2013).
It remains unknown whether the nucleolus plays a protective role in other plant tissues or in other organisms. In human cells, massive DNA damage of the rDNA leads to the formation of nucleolar caps. It has been shown that these caps contain broken rDNA which then becomes available to the HR machinery of the nucleus (Sluis et al., 2015), lending support to the idea that the nucleolus represents a general and conserved HR-refractory sub-compartment. Hence, the nucleolus might have the intrinsic property of excluding HR-related proteins. In line with this idea, the nucleolar proteomes of Arabidopsis and of humans showed no evidence of the presence of HR proteins (Andersen et al., 2005;Montacié et al., 2017).

CONTROLLING SEQUENCE HOMOGENEITY AND HETEROGENEITY
The repetitive rDNA loci are considered intrinsically unstable genomic regions since they are prone to various types of DNA damage and repair events. The sequence variations identified in individual rDNA units (Chandrasekhara et al., 2016;Havlová et al., 2016;Rabanal et al., 2017b;Sims et al., 2021) may represent past DNA repair events following an error-prone pathway (NHEJ). Taking into consideration rDNA copy numbers, it is possible to evaluate the history of DNA repair events following an error-free pathway (HR). While sequence variations of rDNA units can readily be analyzed in individual plants, the evaluation of rDNA unit copy number variations demands the analysis of large populations or multiple successive generations (Rabanal et al., 2017b;Sims et al., 2021). The rDNA copy number can also be considered as a genetic trait and studied in pedigrees. Indeed, when analyzing the trait of "rDNA copy number" over a few generations (two generations in F2s, about eight in recombinant inbred lines -RILs), it appears stable enough that it can be mapped FIGURE 1 | The TCS networks were inferred from rDNA units of the contigs F2N4-F1B23-F2G13, F1F17-F2C3-F2J17, F1F11-F1N27 and F2G18-F19A6 identified in (Sims et al., 2021). The first unit of BAC F2N4 and the fourth unit of BAC F1F11 were excluded from the analysis. Furthermore, the highly repetitive SalI boxes were not taken into account for this data analysis. In the network, each node represents a unique sequence and its size is proportional to its frequency within the data. Short vertical bars on the lines connecting similar sequences represent the number of variations between them. Visualizations of the analyses of the 59 individual BACs are available upon request.
to either NOR in segregating populations. Moreover, the apparent lack of F1-like rDNA copy number phenotypes after several generations of inbreeding in a RIL population further strengthens the notion that the NORs of homologous chromosomes rarely recombine, in agreement with the idea that the nucleolus is a HR-refractory sub-compartment of the nucleus. Importantly, analyzing a wider generational time window, a progressive divergence in the number of rDNA units in single seed descent A. thaliana plants was apparent within tens of generations. As a consequence of this unstable inheritance, and in spite of the fact that rDNA unit numbers vary considerably in natural A. thaliana populations (Davison et al., 2007;Long et al., 2013), genome-wide association studies failed to map the source of the variation to either of the NORs (Long et al., 2013). This means that rare events of HR might take place, only evident in large populations or when observing multiple successive generations, which lead to dramatic rDNA unit number variations.
In contrast, within plants containing a small amount of 45S rDNA units, the rDNA gene copy numbers can be quickly restored and amplified to wild-type levels. This indicated that there is mechanism in place to restore the 45S rDNA copy numbers within individuals with low amount of rRNA genes (Pavlistova et al., 2016).

FUNCTIONAL AND EVOLUTIONARY IMPACT OF rDNA HETEROGENEITY
Different studies on different organisms (including humans, flies, worms, and plants) have shown that the rDNA genes are not identical either within or among individuals of the same species (Gonzalez et al., 1985;Keller et al., 2006;Stage and Eickbush, 2007;Pillet et al., 2012;Bik et al., 2013). Nevertheless, certain SNPs/Indels are stable and abundant enough in either of the two NORs in A. thaliana that they qualify to serve as reporters of NOR-specific expression (Chandrasekhara et al., 2016;Rabanal et al., 2017a). There is unequivocal evidence of selective silencing of one of the two NORs during vegetative development in A. thaliana, with the majority of all rRNAs being generated just from one locus. Nevertheless, there is also compelling evidence that (1) there is selective transcriptional activation of certain rDNA units from the otherwise silenced NOR locus in some tissues and (2) that not all rDNA units at the active NOR locus are transcribed at the same time (Pontvianne et al., 2013). rDNA unit variants are not randomly distributed along the NORs [at least established for NOR2 (Sims et al., 2021)], but rather in variant sub-clusters that share certain SNPs/Indels combinations. In some instances, these blocks of corresponding rDNA units are disrupted by rDNA units of a different subtype, but still are transcriptionally co-regulated (Sims et al., 2021). These findings provide a solid base for the future dissection of the fine-tuned regulation of expression of rDNA variant units within a NOR.
It is tempting to speculate that the heterogeneous population of rDNA units and their regulated expression has an important impact on protein translation. The presence of expressed rRNA variants has been shown in various organisms by analyzing total RNA (Kuo et al., 1996;Carranza et al., 1999;Tseng et al., 2008;Rabanal, Mandáková, et al., 2017;Simon et al., 2018). Some of the identified SNPs/InDels were located within the genic regions that encode the 25S and 18S rRNA subunits which are integrated into ribosomes. Furthermore, several studies demonstrated that variant rRNAs are incorporated into polysomes, the ribosomal fraction actively committed to protein translation (Gonzalez et al., 1988;Cloix et al., 2002;Mentewab et al., 2011;Dimarco et al., 2012;Kurylo et al., 2018;Parks et al., 2018;Sims et al., 2021). Interestingly, various of these rRNA gene variants are differentially expressed in a tissue-specific manner. Furthermore, some sequence variations are located in regions that could have a functional impact on the biology of ribosomes. Most of the genic rRNA sequence variations are located in ribosomal expansion segments, that vary greatly between species, but could have an important impact on interacting proteins. A few SNPs/InDels occur in the rRNA core domains. For instance, one G to T transition present in A. thaliana is located between the H74 and H88 ribosomal domains at the peptidyl transferase site and thus has the potential to impact ribosomal translation directly. In the parasite Plasmodium, two structurally distinct 18S rRNAs are differentially expressed during its life cycle (Gunderson et al., 1987;Waters et al., 1989). And more recently, the expansion segment 9S has been shown to selectively recruit Hox9 mRNA via its 5' UTR stem-loop (Leppek et al., 2020).
In addition, it has been shown that in the bacteria Vibrio vulnificus, from a heterogeneous population of ribosomes, it primarily uses ribosomes containing a particular ribosomal RNA variant to translate stress-related mRNA (Song et al., 2019;Leppek et al., 2020). Similarly, in Escherichia coli, a specific branch of the stress response utilizes a truncated rRNA to selectively bias translation of stress response proteins (Vesper et al., 2011).

CONCLUSIONS/PERSPECTIVES
The most current sequencing technologies, in combination with detailed and large-scale population studies and in-depth analyses of ribosomal RNA variants, have generated novel and exciting insights.
Without any doubt, the NORs cannot be regarded as stable, rigid domains comprised of -nearly -identical rDNA units anymore, but rather as dynamic chromosomal loci with high variation in rDNA unit copy numbers and sequences. We consider a delicate balance of the HR and NHEJ DNA repair mechanisms to be responsible for the dynamic nature of the NORs. We suggest that frequent (meiotic) NHEJ events generate abundant SNPs and InDels within the rDNA, resulting in a heterogeneous population of rDNA units. We also propose that rare HR events dynamically change rDNA unit numbers. The latter may only be observed in large populations and/or over many generations (Figure 2). Furthermore, the ribosomes are no longer seen as invariant machines that translate proteins from available mRNAs but rather as a heterogeneous population of ribonuclear complexes, differing in rRNA and protein composition, with defined functions controlling protein translation (Figure 2).
In the future, it will be interesting to generate the detailed sequence information of NORs from various organisms, ecotypes, and individuals. Knowledge of the precise rDNA unit sequences will allow detailed analyses of the dynamic changes of the NORs, their (potentially context dependent) differential transcriptional regulation, and the integration of rRNA variants into actively translating ribosomes (with the potential to impact protein translation).

AUTHOR CONTRIBUTIONS
JS, FR, and PS: conceptualization. JS and CE: data analysis and visualization. PS, AH, and FR: funding acquisition. All authors contributed to the writing of the article and approved the submitted version.

ACKNOWLEDGMENTS
We thank the members of the PS and AH labs for critical discussions and constructive feedback.

SUPPLEMENTARY MATERIAL
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2021.773052/ full#supplementary-material Supplementary Figure S1 | TCS network was inferred from all sequences of the rDNA units of 59 BACs without the region of the highly repetitive SalI boxes. To build the network, the integrated method of the TCS approach, which is based on the concept of statistical parsimony, in PopArt was used. In the network, each node represents a unique sequence, its size is proportional to its frequency within the data, and its color indicates from which BAC a rDNA originated from. In other words, a big, multicolored node is a collection of different rDNA units which are identical and come from different BACs (and therefore different locations within the NOR2). Short vertical bars on the lines connecting similar sequences represent the number of variations between them. Furthermore, sequences, not present in the data, were inferred and represented as small black dots. The TCS network contains 39 nodes which comprise more than one rDNA unit. In total, 238 rDNA units occur in the 39 nodes. The biggest node includes 62 rDNA units.

A B
FIGURE 2 | (A) Illustration of transcription of variant rRNAs from non-identical 45 rDNA units and their integration into translating ribosomes. The concept of heterogeneous ribosomes has been introduced considering different protein compositions. Here, this concept is extended, also considering different rRNA variants.
(B) Diagram illustrating the occurrence of non-homologous end joining (NHEJ) and homologous recombination (HR) as DNA repair modes in the highly repetitive nucleolus organizing regions (NORs) during meiosis. NHEJ is considered to be the commonly deployed repair pathway, leading to short-range repair scars in the affected rDNA units, contributing to sequence heterogeneity and preserving the integrity of the NOR. Meiotic DNA repair events via HR are considered rare events and will only be evident in large populations, over multiple generations. HR may contribute to NOR size variability and rDNA unit homogenization.