Environmental and Genetic Factors Affecting Apospory Expressivity in Diploid Paspalum rufum

In angiosperms, gametophytic apomixis (clonal reproduction through seeds) is strongly associated with polyploidy and hybridization. The trait is facultative and its expressivity is highly variable between genotypes. Here, we used an F1 progeny derived from diploid apomictic (aposporic) genotypes of Paspalum rufum and two F2 families, derived from F1 hybrids with different apospory expressivity (%AES), to analyze the influence of the environment and the transgenerational transmission of the trait. In addition, AFLP markers were developed in the F1 population to identify genomic regions associated with the %AES. Cytoembryological analyses showed that the %AES was significantly influenced by different environments, but remained stable across the years. F1 and F2 progenies showed a wide range of %AES variation, but most hybrids were not significantly different from the parental genotypes. Maternal and paternal genetic linkage maps were built covering the ten expected linkage groups (LG). A single-marker analysis detected at least one region of 5.7 cM on LG3 that was significantly associated with apospory expressivity. Our results underline the importance of environmental influence in modulating apospory expressivity and identified a genomic region associated with apospory expressivity at the diploid level.


Introduction
A small proportion of seed plants are able to avoid meiosis and fertilization, to produce seeds identical to their maternal parent, by apomixis [1,2]. This natural reproductive system owns great agronomic importance as it allows hybrid combinations to be perpetuated without the need to repeat the original crosses every season. It also accelerates breeding programs and enables seed propagation of crops reproduced vegetatively [1,2]. In this sense, apomictic plants are potentially a constant source of renewable seeds, which would be of great advantage for crop production, especially for developing countries, if their uses became freely available to the public and private sectors [3,4].
Gametophytic apomixis in nature is strongly associated with polyploidy and almost all gametophytic apomictic are polyploid while their sexual relatives are typically diploid [5,6]. Gametophytic apomixis involves the formation of unreduced embryo sacs from the megaspore mother cell itself after avoiding meiosis (diplospory) or from a nucellar cell of the ovule (apospory) [7,8]. In both cases, the embryo originates by parthenogenesis of the non-reduced egg cell. In addition, the endosperm formation may require the fertilization of polar nuclei (pseudogamy) or develop autonomously [9]. Although apomixis is widely distributed among angiosperms, it is rare in crop gene pools. Despite its presence in wild relatives of some cereal crops, introgression approaches had failed to transfer the character [10,11]. For this reason, the  Norrmann et al. 1989 In the present work, 49 additional hybrids of the population (hereafter named F 1-Z ), were characterized, together with the parental genotypes, in Zavalla, in the central region of Argentina ( Figure 1 and Table 2). Reproductive classification (%AES) was performed by cytoembryological observation of cleared ovaries at anthesis. Parental genotypes, R6#45 and R5#49, characterized in Zavalla, showed 3.6% and 6.4% of AES, respectively, while in the F 1-Z apospory expressivity ranged from 0% to 17.8% (mean: 4.56% ± 3.82) ( Figure 1 and Table 2). The comparative analysis between each hybrid of the F 1-Z and the parental genotypes revealed that the %AES of most of the hybrids was included within the range of both parental genotypes values ( Figure 1). However, higher apospory expressivity was detected in two plants (#72 and #86) which showed a significant increase (13.75% and 17.71%, respectively) with respect to both parental genotypes. On the other hand, three individuals: #67, #74, and #89 did not show any AES ( Figure 1 and Table 2). Considering all F 1 plants ((F 1-C and F 1-z ), 84 hybrids (95.45%), out of 88 showed AES and only four hybrids (4.54%) had exclusively MES.

Evaluation of the Environmental Influence on Apospory Expressivity
Following the observations described above, we analyze the stability of apospory expressivity in the two environments assayed. Therefore, five genotypes with different %AES were selected from the F 1-C population to be evaluated at both Zavalla and Corrientes. Two diploid hybrids with low and no apospory expressivity (F 1-C #9 and F 1-C #12, respectively) and three with relatively high expressivity (F 1-C #31, F 1-C #15, and F 1-C #39) were selected as well. Moreover, tetraploid genotypes with different apospory expressivity were also included: one highly apomictic natural tetraploid Q3756 and two artificially generated tetraploids LD1 and LD3 ( Table 1).

Apospory Expressivity Variation in F 2 Progenies
To analyze the effect of hybridization on apospory expressivity, two F 2 progenies obtained by crossing two pairs of plants from the F 1-C population with similar (low or high) apospory expressivity were generated. Thus, F 1-C #9 × F 1-C #12 (mean 1.40% ± 1.98, Table 4) and F 1-C #15 × F 1-C #39 (mean 8.79% ± 1.81, Table 3) were crossed to obtain the families F 2L and F 2H, respectively. The apospory proportion of both populations was determined for two consecutive years at Zavalla. The characterization revealed no significant difference between both periods for each individual, only one hybrid (F 2H #8) showed a significant variation in the %AES between both years (Tables 3 and 4).
Individuals of each F 2 progeny were compared to their parental genotypes. Samples fixed at Zavalla in the same year (2017) were used for this analysis, except for F 1-C #9 which was only analyzed in 2018. In F 2L , the mean %AES of hybrids along both years showed a range of variation from 0% to 13.3% AES (Table 4), and most individuals were not significantly different from the parental genotypes ( Figure 4A). However, the plant F 2L #10 showed a significant increase in apospory expressivity with respect to its progenitors ( Figure 4A, Table 3). In F 2H , all individuals produced AES, the mean %AES, over both years, ranged from 1.82% to 15.61% (Table 3). The comparison of the %AES, between the parental genotypes, and each hybrid during 2017, revealed that most hybrids were not different from the parental genotypes ( Figure 4B), only F 2H #7 showed a significant reduction in apospory expressivity (0.9%) with respect to parental genotypes ( Figure 4B and Table 3). Overall, these results verify the stability of the trait across different years and showed that, in the two F 2 progenies, as well as in the F 1 , the parental phenotypes seem to delimit the range of %AES variation of the progeny. Notably, the F 2H family, which descend from both aposporous parents with higher proportions of AES, showed, in both years, mean values of %AES (mean: 10.05 and 7.19, median: 10.8 and 6.5 respectively) ( Table 3, Figure 4C), significantly higher than the means of the F 2L population (mean: 3.33 and 2.75, median: 1.0 and 2.1) ( Table 3, Figure 4C), which derive from parents with only one aposporous parental genotype, with low %AES ( Figure 4C). Moreover, the F 2H median value was also significantly higher than the mean proportion of apospory expressivity of the F 1-Z, (mean: 4.60, median: 3.26), which come from both parents (R6#45 and R5#49) with intermediate %AES (Table 1, Figure 4C).  Figure 4C), significantly higher than the means of the F2L population (mean: 3.33 and 2.75, median: 1.0 and 2.1) ( Table 4, Figure 4C), which derive from parents with only one aposporous parental genotype, with low %AES ( Figure 4C). Moreover, the F2H median value was also significantly higher than the mean proportion of apospory expressivity of the F1-Z, (mean: 4.60, median: 3.26), which come from both parents (R6#45 and R5#49) with intermediate %AES (Table 1, Figure 4C).

F1 genetic Linkage Map
To develop molecular tools to study apospory inheritance and identify genomic regions associated with apospory, we decided to build a genetic linkage map diploid P. rufum. At first, to select informative AFLP markers, 105 selective primer combinations were tested on parental genotypes. This analysis identified 33 highly polymorphic combinations, each producing more than 15 polymorphic fragments (Supplementary Materials Table S2) that were selected for linkage analysis in the F1 progeny (Supplementary Materials Figure S2). Good-quality AFLP libraries were obtained for 87 out of 88 F1 offspring (F1-Z#90 was discarded). The hybrid origin of F1 plants was confirmed for all 87 individuals included in the analysis. A binary data matrix was built for each parental genotype containing 278 maternal, 239 paternal, and 110 bi-parental (segregating from both parental genotypes) markers (Table 5). Markers fitting the expected segregations, 1:1 for SDAF (single dose amplification fragment) and 3:1 for BSDF (bi-parental single dose fragment), were used for mapping. Segregation analysis confirmed the expected values (0 < χ 2 < 3.8) for 266 maternal, 221 paternal, and 94 BSDF markers (Table 5). About 7.4% of the markers showed a variable degree of distorted segregation ratios in the progeny (4.1 < χ 2 < 104.3) ( Table 5). Two linkage maps corresponding to the maternal and paternal parents were constructed.

F 1 Genetic Linkage Map
To develop molecular tools to study apospory inheritance and identify genomic regions associated with apospory, we decided to build a genetic linkage map diploid P. rufum. At first, to select informative AFLP markers, 105 selective primer combinations were tested on parental genotypes. This analysis identified 33 highly polymorphic combinations, each producing more than 15 polymorphic fragments (Supplementary Materials Table S2) that were selected for linkage analysis in the F 1 progeny (Supplementary Materials Figure S2). Good-quality AFLP libraries were obtained for 87 out of 88 F 1 offspring (F 1-Z #90 was discarded). The hybrid origin of F 1 plants was confirmed for all 87 individuals included in the analysis. A binary data matrix was built for each parental genotype containing 278 maternal, 239 paternal, and 110 bi-parental (segregating from both parental genotypes) markers (Table 5). Markers fitting the expected segregations, 1:1 for SDAF (single dose amplification fragment) and 3:1 for BSDF (bi-parental single dose fragment), were used for mapping. Segregation analysis confirmed the expected values (0 < χ 2 < 3.8) for 266 maternal, 221 paternal, and 94 BSDF markers (Table 5). About 7.4% of the markers showed a variable degree of distorted segregation ratios in the progeny (4.1 < χ 2 < 104.3) ( Table 5). Two linkage maps corresponding to the maternal and paternal parents were constructed. The data matrix of the maternal progenitor (R6#45), including a total of 360 markers (266 SDAF and 94 BSDF), was first filtered and thus 75 markers were discarded due to similarity and two individuals were removed due to 96.5% of genetic similarity. As described in the Section 4, an initial grouping was assayed at LOD > 6.0 which fixed the order of 219 markers, with a mean distance of 7.5 cM (±6.4), in 20 LG in the coupling phase. A second test at LOD > 2.0, using the above-fixed order, added 35 additional markers, extending the total distance from 1478.4 cM to 1672.0 cM and shortening the average inter-markers distance to 7.14 cM (±6.5). In the next step, maternal homologous LGs were combined, including markers that segregated in the repulsion phase, as described in the Materials and Methods. Finally, the maternal map incorporated 234 markers distributed in 10 LGs, covering a total distance of 1071.8 cM with a mean inter-marker distance of 4.78 cM (±4.8 cM) and a maximum distance of 28.74 cM. Each LG contained an average of 23.4 (±11.7) markers with a maximum of 44 and a minimum of six markers. The largest LG, LG4, covered 151.9 cM and contained 33 markers. The shortest was LG10, which covers 62.7 cM and included 10 markers ( Figure 5 and Supplementary Materials Table S3). Based on bi-parental markers segregation, it was possible to identify eight out of the ten female and male homologous LGs (Figure 7, Supplementary Materials Table S4). As there were no bi-parental markers on the LG6 and LG10 numbers, they were arbitrarily assigned ( Figures 5 and 6). The relative order between bi-parental markers in both maps showed that most of them are equally ordered in maternal and paternal maps (Figure 7). In the paternal (R5#49) map, the initial matrix included a total of 315 markers (221 SDAF and 94 BSDF) ( Table 5) and after filtering 74 markers were discarded due to their similarities. Following the same procedure as for the maternal map, the first approach ordered 185 markers in 20 LGs, covering 965.8 cM and a mean inter-marker distance of 5.85 cM (±5.3). The next step at LOD > 2 added 29 additional markers, extending the total distance to 1193 cM and the mean distance between markers to 6.15 cM (±6.3). Finally, after combining homologous LGs, the paternal map included 212 markers, distributed in 10 LGs, covering a total distance of 914 cM. The average distance between markers was 4.53 cM (±5.4) reaching a maximum distance of 34.9 cM. Each LG contained an average of 21.2 (±11.8) markers, LG1 was the largest with 50 markers, covering 117.2 cM, and LG10 which was the shortest containing 14 markers distributed over 46.7 cM (Figure 6, Supplementary Materials Table S3).
Based on bi-parental markers segregation, it was possible to identify eight out of the ten female and male homologous LGs (Figure 7, Supplementary Materials Table S4). As there were no bi-parental markers on the LG6 and LG10 numbers, they were arbitrarily assigned (Figures 5 and 6). The relative order between bi-parental markers in both maps showed that most of them are equally ordered in maternal and paternal maps (Figure 7).

Identification of Genomic Regions Associated with Apospory Expressivity
A linear regression analysis was performed combining phenotypic information of F1

Identification of Genomic Regions Associated with Apospory Expressivity
A linear regression analysis was performed combining phenotypic information of F 1 hybrids, as the dependent variable, and genotypic data derived from the mapping procedure, as the independent variable [70,71]. Due to the variation of the %AES between locations, F 1-C and F 1-Z populations were analyzed separately. This single-marker analysis detected a total of 40 significant associated markers (p < 0.05) with an R 2 from 0.1 to 0.34. The evaluation of the F 1-C showed 26 markers associated with apospory expressivity (0.1 < R 2 < 0.34), 23 of which were included in the map. As three of them had been previously eliminated by similarity, their position specifications correspond to those similar mapped markers (B16_164.6, C29_332.2, C32_244.8), the other three remaining unlinked ( Table 6). As shown in Table 6, the LGs with the highest number of associated markers were: LG3, with a total of seven markers, six maternal (M), and one paternal (P), and LG4 with a total of 10 markers, seven M, one P and two bi-parental (BP). The other LGs containing associated markers were: LG2 (one marker), LG5 (one marker), and LG10 (three markers). The markers with the highest R 2 (>0.15) were found in LG2, LG3, and LG4 (Table 6). Regarding the F 1-Z population, 13 markers showed a significant association (p<0.05) with R 2 ranging from 0.10 to 0.17, 12 markers were mapped, one of which was assigned by similarity as mentioned above (C12_83), ( Table 7). As for the F 1-C , LG3 showed the highest number of associated markers, with a total of 5 markers (two M, two P, and one BP), while LG6 showed two paternal markers, LG7 two maternal markers, and LG1, LG8, and LG10 showed one marker each (P, M, and M, respectively) ( Table 7). The markers with the highest R 2 (>0.15) were found in LG3 and LG7. Interestingly, both F 1 populations showed markers associated with the trait on LG3 and also on LG10, but the latter has an arbitrary number in both parental maps as it was not possible to identify female and male homologous. Then, focusing on LG3, both maternal LG3 and paternal LG3 were combined using the "Combine Groups for Map Integration" function (JoinMap 4.0) to characterize the relative order between the two groups of markers. The integration generated a new LG3 of 131.2 cM with a total of 39 markers. Most of the markers associated with apospory expressivity, identified from each population, clustered in different regions of the LG3, however, one marker from the F 1-C , B4_490.5 (R 2 = 0.13; p = 0.03) and another from the F 1-Z , C4_173.5 (R 2 = 0.12; p = 0.02) enclosed a region of 5.7 cM (Figure 8).  Interestingly, both F1 populations showed markers associated with the trait on LG3 and also on LG10, but the latter has an arbitrary number in both parental maps as it was not possible to identify female and male homologous. Then, focusing on LG3, both maternal LG3 and paternal LG3 were combined using the "Combine Groups for Map Integration" function (JoinMap 4.0) to characterize the relative order between the two groups of markers. The integration generated a new LG3 of 131.2 cM with a total of 39 markers. Most of the markers associated with apospory expressivity, identified from each population, clustered in different regions of the LG3, however, one marker from the F1-C, B4_490.5 (R 2 = 0.13; p = 0.03) and another from the F1-Z, C4_173.5 (R 2 = 0.12; p = 0.02) enclosed a region of 5.7 cM (Figure 8).

Discussion
The inheritance and genetic analysis of gametophytic apomixis have long been carried out mainly from a qualitative perspective, on natural polyploid systems, leading to the inherent difficulties of the ploidy level and also limiting the understanding of the variation of apomixis expressivity. Paspalum rufum offers a diploid system to study the inheritance of apomixis that would avoid the problems related to polyploidy. In the current work, we studied apomixis from a quantitative perspective by analyzing the variation of apospory expressivity in two consecutive generations of P. rufum diploid hybrids, and in different environments. In addition, a genetic linkage map was built and used for searching genetic regions that could be modulating the extent of the trait.

Apospory Inheritance
Although the main purpose of our work was to analyze the expressivity of apospory, it is necessary to discuss how apospory capacity is transmitted to the offspring at the diploid level. Our results showed a stable inheritance of the capacity to produce AES over two generations. In the F 1 population, 84 individuals out of 88 were able to produce AES, and also in both F 2 progeny, almost all hybrids were able to produce AES. As in F 1 , both parents are natural aposporic genotypes, and their genetic origin is unknown, we could not rule out the hypothesis of a single dominant gene, in homozygous or multiallelic configuration, which would confer aposporic ability to all progeny. In this sense, individuals considered as sexual would also carry aposporic potential but in a low proportion that is difficult to detect. This could be explained by the fact that inter-annual analyses revealed that some genotypes classified as completely sexual one year were able to produce a very low proportion of ovules with AES another year. In this context, our observation would be consistent with previous genetic analyses on the inheritance of apomixis at the tetraploid level in the genus Paspalum, which mainly point to a single dominant locus involved [7,12,13]. Then, the genetic determinant(s) of apospory would be constitutively present in the F 1 offspring of P. rufum, but the expressivity of the trait would be modulated by other unknown internal and/or external factors.

Apospory Expressivity Variation
Based on the results described above, two main questions arise: how the environment influences apospory expressivity and whether genetic/epigenetic background determines the degree of apospory expressivity. To start unravelling these issues concerning the first question, we analyzed the behavior of the trait in different environmental conditions. Our analyses have shown that when the same genotypes were grown in different locations, they showed significant variation in the degree of apospory expressivity. Several studies have reported variation in apospory under different circumstances. It was previously found that in facultative apomictic tetraploid genotypes of P. notatum, the frequency of apomixis was found to vary during the flowering period, being higher at peak flowering time and decreasing until the end of the flowering season [65]. Similar variations had been reported for highly apomictic wild P. notatum and P. cromyorrhizon, where apomictic reproduction varies throughout the flowering period [64,72]. Furthermore, tetraploid genotypes of P. cromyorrhizon were exposed to different day length treatments and the increment of light exposure from 12 to 14 h increased SES and decreased AES production [64]. Accordingly, a similar analysis in facultative apomictic Ranunculus auricomus reported that an extended photoperiod resulted in a higher proportion of sexual ovules [73]. Therefore, the authors concluded that environmental conditions, which affect flowering, may also influence the expression of apospory [64] and that residual sexuality could occur under less favorable environmental conditions [65]. A recent study on the sexual versus apomixis expression in P. intermedium reported that sexual reproduction increased at lower temperatures in facultative apomictic plants [74]. Regarding our results, we must point out that the main changes between Corrientes and Zavalla locations are the latitude and average temperature. These differences are associated with extended photoperiods [75] and lower temperatures [76], throughout the flowering months (from October to December) in Zavalla with respect to Corrientes. In addition, during sample fixation, along the years of our study, we noticed that there was a delay in flowering time at Zavalla with respect to Corrientes (data not shown) from the end of September to late October/early November, which also implies longer days. Thus, our observations reinforce previous results showing that apospory expressivity, in diploid and tetraploid genotypes of P. rufum, is significantly affected by environmental variations. As P. rufum plants are native to Corrientes [58], reallocation to Zavalla could result in a stressful condition that induces sexuality, which is consistent with the previously reported effect of stress on reproductive development in other apomictic systems [77,78]. Despite the location (latitude) effects, our work also showed that apospory expressivity in P. rufum is stable across the years when plants remain at the same location.
To deepen our understanding of the genetic/epigenetic influence on apospory expressivity, we analyzed the quantitative variation of the trait across successive generations. To rule out environmental influences on the trait, all plants were maintained in the same location. Then, two crosses were performed one between F 1 hybrids with low apospory and the other between hybrids with high apospory, generating two F 2 progenies, F 2L and F 2H , respectively. We observed that the differences detected between the parents were maintained between the two F 2 progenies so that the expressivity of the parental apospory was transmitted to the progeny and was also stable over different years. Accordingly, apospory expressivity in the F 1-Z population also showed low ranges of variation, similar to those of F 2L , in line with the expressivity of its parents. A similar study was carried out in Ranunculus where apospory expressivity was analyzed in both F 1 and F 2 progenies. These results revealed that apospory expressivity was significantly increased in F 2 progeny when both F 1 parents, rather than just one, were aposporous [79]. In another study in Hieracium, on the assessment of the autonomous endosperm formation, it was revealed that the trait is conferred by a single dominant locus but additional genetic factors modulated its expressivity [44]. These outcomes are also in agreement with our previous results, where we found that both autopolyploidy and hybridisation of diploid P. rufum genotypes, induced higher expressivity of the trait [69]. Consequently, we propose that the degree of apospory expressivity, in the diploid system of P. rufum, is transmitted to successive generations in a dose-dependent manner and could be modulated by genetic/epigenetic factors. Regarding the latter, a recent comparative study of epigenetic patterns between the F 1-C hybrids used as parental plants for F 2L versus those used for F 2H revealed that the differential apospory expressivity of these genotypes was associated with different methylation patterns [80]. Thus, our observations of the transgenerational transmission of apospory expressivity could be associated with the fact that in plants, epiallelic variation can be stably propagated over several generations [81,82].
From an evolutionary perspective, our results are in agreement with a recent analysis of P. intermedium populations suggesting that environmental modulation of sexuality and apomixis expression would provide genetic variability to facultative apomictic linages allowing them to adapt to ecological challenges [74]. Furthermore, the effect we observed of hybridization on apospory expressivity would favor polyploidization within natural diploid sexual populations through the triploid bridge. This would indirectly help to establish an apomictic cytotype, as recently proposed by Hojsgaard et al. [63].

Linkage Map
To better understand the genetic basis modulating the expressivity of apospory at the diploid level, we built a genetic linkage map with AFLP markers. Two genetic linkage maps were generated for each parental genotype. The size of maternal and paternal maps were 1071.8 cM and 914 cM, respectively, which were similar to the maps of related diploid species such as Paspalum notatum, (991 cM) [83], Setaria italic (964 cM) [84], and Sorghum bicolor (1095 cM), [85]. Both linkage maps had ten LG which corresponds to the basic chromosome number of P. rufum (n = x = 10). The bi-parental markers identified eight of the ten homologous LGs between the two maps and also confirmed the relative order of the markers.
The maternal and paternal maps presented a mean distance between markers of 4.8 cM and 5.8 cM, respectively. These values are comparable and even lower, than those obtained in the previous map of Paspalum notatum [22,83] indicating a good relative density. The DNA content of P. rufum is 0.75 pg/Cx [57] which is equivalent to 733.5 Mbp, according to the 1pg/978 Mbp ratio previously reported [86]. Consequently, considering the total length (cM) of both recombination maps, 1 cM would be equivalent to a mean of 743.4 kpb. These are expected values for plant genomes and comparable to those estimated for the diploid P. notatum map (0.57 pg/Cx) in which 1cM is equivalent to 559 kpb [83]. Additionally, we found that the recombination frequency was similar for the maternal and paternal genotypes in contrast to the higher recombination frequency found in the genetic linkage map of paternal apomictic tetraploid P. notatum [22].
The map generated will be useful in future works to map markers linked to apomixis that have been previously found in other Paspalum species. It can also be useful to search for markers associated with parthenogenesis, as one of the parental genotypes of the F 1 is capable of completing apomixis [68].

Quantitative Approach
Apospory expressivity regulation is in itself a very important aspect to take into account as in P. notatum breeding programs, sexual × apomictic offspring show a high range of variation in apospory expressivity [66,87]. Several works have dealt with the problem of obtaining a low proportion of highly apomictic hybrids after crosses between sexual and apomictic progenitors [45,46,66,87]. Therefore, although molecular markers linked to apospory are very useful in saving time and resources in breeding programmes, the identification of obligate apomictic hybrids with high agronomic performance is still limited [66]. These led us to analyze apospory expressivity from a quantitative perspective. Some research works have already considered apospory as a quantitatively modulated trait [29,[48][49][50]. For P. rufum we observed that the F 1 generation shows a wide range of apospory variation [69], in addition, some progeny genotypes showed a significant increase in apospory expressivity compared to the parental genotypes, which is in agreement with transgressive segregation [88]. This behavior could be attributed to recombination between parental genotypes, which have quantitative trait loci (QTL) with antagonistic effects [89]. A similar variation of apospory expressivity was previously reported in hybrids of Paspalum notatum [46] and Hieracium [29], leading the authors to suggest that unknown factors segregate in the genetic background modifying apospory expressivity. In line with these observations, recent comparative gene expression work between apomictic and sexual genotypes of Ranunculus auricomus reported that the pattern of gene expression was reflective of transgressive and genome dosage effects, supporting the hypothesis of a dominant factor controlling apomixis and modulated by secondary modifiers [90].
In this context, we looked for genetic factors that would be modifying apospory expressivity in the mapping population of P. rufum. We performed a single-marker analysis between markers and apospory expressivity. As we saw that different locations affect apospory expressivity, F 1 sub-populations were analyzed separately. This initial analysis showed that the markers significantly associated with apospory expressivity were scattered across several LG. In F 1-C , associated markers were able to explain 10-34% of the phenotypic variation, and LG3 and 4 had both the highest number of markers and the highest R 2 . In F 1-Z , the number of significantly associated markers was about half of those found in F 1-C , with a lower range of R 2 from 10-17%. The wide dispersion of markers along the genome was in line with recent work in Eragrostis curvula where diplospory has, for the first time, been assessed as a quantitative trait. Although qualitative analysis confirmed that diplospory is determined by a single dominant locus, QTL analysis revealed that diplospory expressivity could be regulated by five different regions, two very closely linked to diplospory locus and three additional regions distributed along two different LG [50].
We observed that each of the F 1 sub-populations showed different sets of associated markers, which could be explained by the strong influence of the environment on apospory expressivity. However, some of them coincide in the LG3, which had a total of 14 markers in both F 1-C (8) and F 1-Z (6). The integration of maternal and paternal LG3s showed that the markers found in each location clustered in different genomic regions. Notwithstanding, a marker detected in Zavalla was close (5.7 cM) to another marker detected in Corrientes, so this region could be considered in future analyses as a putative constitutive QTL.
The low proportion of phenotypic variation explained by the associated markers could be attributed to the low range of variation in apospory expressivity. This is evidenced by the difference observed between the two sub-populations as F 1-C , which ranges from 0 to 35%, had both a higher number of associated markers and higher R 2 values than F 1-Z which ranges from 0-17%. On the other hand, the small size of each sub-population could also have a negative impact on the results. Both aspects should be taken into account in future analyses. The molecular tools developed here will be the starting point for future molecular analyses of apomixis in the diploid system of P. rufum.

Plant Material
The starting plant material was the F 1 progeny previously developed from the cross between two natural diploid genotypes: R6#45 used as pistillate parent and R5#49 used as pollen donor [69] ( Table 1). The 39 individuals of the F 1 offspring (named F 1-C in advance) had been previously established in the IBONE experimental field (Latitude: −27.469213. Longitude: −58.830635) and their reproductive mode had already been determined by cytoembryological observation [69]. In the present work, 49 additional individuals were added to the original F 1-C population, increasing its size to 88 individuals, which were established at the Agronomy Collage of the National University of Rosario, located in Zavalla (FCA-UNR, Latitude: −33.018538. Longitude: −60.878858) (named F 1-Z in advance). F 2 progenies were developed by crossing F 1-C individuals, selected by their apospory expressivity. Two individuals with low apospory expressivity: F 1-C #9 (1.2% AES) and F 1-C #12 (0% AES) and two with high apospory expressivity: F 1-C #15 (32.8% AES) and F 1-C #39 (35.8% AES) were crossed as previously described [69], to generate F 2L and F 2H progenies respectively. Both F 1-C #9 and F 1-C #15 were used as pistillate parents, while the other two individuals as pollen donors, respectively (Table 1). Filled seeds were germinated in sterilized soil, and seedlings were planted in small pots in a greenhouse. Subsequently, the plants were transferred to a field nursery and grown under natural conditions at the FCA-UNR. To verify the hybrid origin of F 1 and F 2 plants and avoid individuals generated by self-fertilization, apomixis, or un-intended cross-contamination marker segregation analysis was performed. F 1 and F 2 progenies were evaluated by AFLP (see below) and RAPD markers respectively. A total of 13 primers of the UBC series were tested and five (UBC301, UBC310, UBC329, UBC344, and UBC349) showed polymorphism between the F 1-C genotypes used as parental genotypes of the F 2 progenies (Supplementary Materials Figure S1, Tables S5 and S6). Three of the five polymorphic markers detected were used to assess the hybridity of F 2L and four for F 2H . Each plant was scored for the presence of paternal and maternal bands and was considered hybrid if it showed at least two pollen donors specific bands (Supplementary Materials Figure S1 and Table S6). The data confirmed the existence of non-maternal offspring and parental combinations in the F 2 generations (Supplementary Materials Table S6). In addition, previously characterized tetraploid genotypes were used, one natural (Q3756) and two artificially obtained by colchicine treatment, LD1 and LD3 [69], (Table 1).

Cytoembriological Analysis
The %AES of each genotype was determined by counting the number of ovaries containing AES, in anthesis inflorescences. Samples were treated according to Soliman et al., (2019) [91]. Briefly, they were initially fixed in FAA (70% ethanol: glacial acetic acid: formaldehyde in the proportion 90:5:5) for 24 h and then transferred to 70% ethanol for at least 24 h, then spikelets were dissected and pistils were cleared following the protocol described by Young et al. (1979) [92]. The ovaries were observed with a Leica DM2500 light transmission microscope equipped with a differential interference contrast (DIC) system and a digital camera (Leica Microsystems DFC 295. Wetzlar, Germany). At least 47 ovaries per plant (mean number of 73) were analyzed for the presence of single meiotic embryo sacs (MES), single or multiple AES, or one MES plus one or more AES (MES + AES), aborted ES (AbES) were also recorded. In individuals with a low number of ovaries containing AES (less than 0.05%); the number of ovaries examined was increased up to at least 134. Apospory expressivity was estimated as the percentage of ovaries containing at least one AES (alone or in combination with MES), out of the total number of ovaries scored. Embryo sac types were classified according to the method described by Norrmann et al. (1989) [55]. In particular, embryo sacs showing the egg apparatus, two polar nuclei, and a cluster of antipodal cells were classified as meiotic of the Polygonum type. Embryo sacs showing the egg apparatus and two polar nuclei, but lacking antipodal cells, were considered AESs of the Paspalum type [67]. The environmental effects of the different locations were evaluated by quantifying the %AES in anthesis ovaries fixed in Corrientes and Zavalla. The temporal stability of the trait was estimated by fixing the ovaries in different years in Zavalla.

AFLP Markers
DNA was extracted by using the CTAB method according to Martínez et al. (2003) [19]. Library construction was performed following the protocol described by Vos et al. (1995) [93] with modifications reported in Cnops et al. (1996) [94]. Briefly, genomic DNA (300 ng) was digested and ligated at 37 • C for 4 h using 5 U Eco-RI, . Cycling conditions to ensure optimal primer selectivity were: four cycles of 45" at 94 • C. 30 at 65 • C, 1 at 72 • C followed by 12 cycles decreasing annealing temperature by 0.7 • C each cycle and 22 cycles as mentioned before. Each product (1 µL) was mixed with 10 µL of formamide and 0.1 µL of standard size (Genescan LIZ 500, Thermo Scientific, Waltham, MA, USA), then were denatured at 94 • C for 5 . Samples were analyzed with the ABI 3130xl Genetic Analyzer (Thermo Scientific, Waltham, MA, USA) using Genemapper 4.0 software. The F 1 hybrid condition was verified by AFLP during the development of the linkage map. Eighty-seven individuals, out of a total of 88, from which good quality AFLP libraries were obtained, were found to have pollen donor-specific bands. Each plant had at least 98 bands of the paternal genotypes, and no plants of apomictic or self-fertilized origin, with only maternal bands, were detected (Supplementary Materials Table S7).

Genetic Linkage Analysis
Segregation data from each parental genotype were analyzed independently. As the diploid F 1 mapping populations come from non-inbred parents, different allelic configurations were expected for each locus [95]. For single dose amplification fragments (SDAFs), which are polymorphic between parents (i.e., present in R6#45 and absent R5#49 and vice versa), the expected segregation ratio was 1:1. Markers that were monomorphic between parents, but segregating in the population, were termed Bi-parental Single Dose Fragments (BSDF) and their expected segregation ratio was 3:1 for presence:absence respectively. A χ2 test was used to determine the goodness of fit between the observed and the expected number of genotypes for each segregation class. For linkage analysis, each parental data file included both SDAF and BSDF markers. Two genetic linkage maps were constructed, one for each parental genotype, using JoinMap 4.1 software [96] and selecting the crosspollinator full-sib population (CP) option. Bi-parental markers were used as allelic bridges to identify homologous linkage groups (LG) between the parental genotypes [95]. Markers with identical segregation in different individuals (>98%) were removed using the JoinMap 4.0 command "similarity of individuals". Clustering analysis was carried out using a LOD score threshold of 6.0 or higher and a default linkage value of 0.4 for recombination frequency, the maximum detectable recombination fraction (maxR). Marker order and genetic distance within each linkage group were calculated using a regression mapping algorithm and the mapping function of Kosambi [97] with default options. Subsequently, a second set of markers was added to the originally established groups with a LOD score < 2.0 using the "fixed order" function. In addition, the segregation data for the presence/absence of those markers linked in the "repulsion phase" were re-coded (inverted) to include them in one linkage group per homologous pair according to al-Janabi et al. (1993) [98].

Linear Regression Analysis
To identify the markers associated with the trait, a linear regression analysis was performed using molecular markers as an independent variable and the %AES as a dependent variable [71]. Each F 1 sub-population, 37 individuals from the F 1-C and 47 individuals from F 1-Z , was analyzed separately due to the variation of the trait, verified in the present work, between the two locations. As apospory expressivity did not have a normal distribution (Shapiro-Wilk's test showed a < 0.05), the data were transformed (Y = Log (%AES + 1)) reaching normal distribution (p = 0.95 and 0.34 for Corrientes and Zavalla, respectively). Markers were considered to be associated when p < 0.05, the R 2 was understood as the proportion of phenotypic variation explained by each QTL [71].

Statistical Analysis
Quantitative comparison between %AES of different genotypes or between the same genotype grown under different conditions was performed according to Delgado et al. (2016) [69]. In brief, confidence intervals (CIs), around the AES frequency of each plant, were calculated using the online resource http://vassarstats.net/ (accessed on 17 August 2021). [99], following the method described by Wilson [100], without continuity correction. Fisher s exact test [101] was used to determine the differences between the AES frequency of the experimental samples; calculations were performed using the online resource MeasuringU (https://measuringu.com/calculators/ab-cal/) (accessed on 17 August 2021). [102] to compare independent proportions, while the non-parametric Mood median test, using STATGRAPHIC Centurion XVI (Version 16.1.03), was used to compare the median of the frequency of apospory between F 1 and F 2 progenies.

Conclusions
Here, we report that apospory expressivity is influenced by the environment in both diploid and tetraploid genotypes of P. rufum. This is an important aspect to be considered for future application of apomixis technology and its introduction into sexual species. Our results show that it is possible to transfer the degree of apospory expressivity from parents to the offspring. The marker analysis allowed us to generate the first linkage map for the diploid system of P. rufum and to identify at least one 5.7 cM genomic region, in an LG, significantly associated with apospory expressivity. The molecular tools developed will be essential in the search for markers associated with other apomictic components from a qualitative and quantitative perspective in the diploid system of P. rufum.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/plants10102100/s1. Figure S1: Amplification profiles of RAPD markers; Figure S2: Amplification profiles of AFLP markers obtained using capillary electrophoresis with the ABI 3130xl sequencer; Table S1: Cytoembryological analysis of P. rufum diploid hybrids and tetraploid genotypes in Zavalla; Table S2: AFLP primers sequences and combinations used to build the genetic linkage maps; Table  S3: Specifications of maternal and paternal map; Table S4: List of bi-parental markers shared between R6#45 and R5#49 maps; Table S5: UBC primer used to verify hybridity in F 2 progenies; Table S6: Segregations of RAPD markers in F 2 offsprings and Table S7: Single dose AFLP markers detected in each F 1 offspring.