A total of six libraries of inter-specific CSSLs were constructed in two O. sativa (Indica and Japonica) backgrounds using three wild (ORSC) donor accessions. To help select wild donors for CSSL development, genetic relationships were examined among 286 ORSC accessions and 54 O. sativa varieties initially using 50 Simple Sequence Repeat (SSR) markers and subsequently using a set of 113,739 GBS SNPs (Kim et al., 2016). A rooted Neighbor Joining (NJ) dendrogram using an Oryza officinalis, Wall ex Watt accession set as the outgroup was constructed with 100 bootstrap replicates (Tamura and Nei, 1993) in Geneious v10.0.9. ORSC donors were selected to represent geographical, phenotypic and genomic diversity, and donor selection was also influenced by cross compatibility with the cultivated recurrent parents (RPs). Three donor plants per wild accession were grown out for crossing with multiple individuals from the two RPs in the Guterman Greenhouse at Cornell University (Kim et al., 2016). Efforts to develop purified donor lines from the ORSC individuals used for crossing were not successful, and therefore additional plants were grown from the original donor seed stocks and used as proxies in some of the subsequent analyses.

Three criteria were used to select lines for generation advance over the course of development: (1) presence of a wild introgression in a target region, i.e., a region that contributed to the goal of providing coverage of the entire donor genome with independent, overlapping substituted segments, (2) minimum number of wild introgressions in non-target regions, contributing to the goal of creating a library consisting of a set of near isogenic lines, each carrying a single wild introgression in an RP genetic background and (3) recovery of a “lost introgression,” in cases where a heterozygous donor segment present in one generation was not recovered in the subsequent generation, requiring the retrieval of a line from an earlier generation. The breeding scheme is presented in Supplementary Figure S1.

Marker-assisted selection (MAS) was carried out at every generation. Foreground selection was performed in the first three generations using 384 Oligo Pool Assays (OPAs) to establish sets of lines that covered the donor genome (Thomson et al., 2012). Subsequent generations utilized 6 K or 7 K SNP arrays (the C6AIR and C7AIR; Thomson et al., 2017; Morales et al., 2020), to eliminate non-target introgressions (negative selection) and increase the proportion of recurrent parent background. Length of each donor segment was determined based on the physical position of the marker genotypes. For a given introgression line, the preferred range for a target donor introgression was 4–5 Mb (13–17 cM), with maximum recurrent parent background. Each introgression line was selected such that there was a 1–2 Mb overlap within its neighboring introgression line on both sides of the introgressed segment. Donor segment selection was given preference over background recovery. Therefore, if donor segments were lost during backcrossing, plants from an earlier generation with less recurrent parent recovery were selected for the final library. For each population, the set of lines representing the final library were analyzed using CSSL Finder v.0.9.722¹. The target segment size was calculated in Microsoft Excel, whereas the number of donor and background segments, and the percent recurrent parent genome, were calculated using CSSL finder. The libraries were visualized using the CSSL finder or GGT (van Berloo, 2008) software during population development.

Genomic DNA was extracted from leaf tissue from individual plants in each backcross generation using a modified CTAB protocol (Dellaporta et al., 1983; Fulton et al., 1995). As new technology became available, various genotyping methods were employed for MAS during library development. These platforms were SSRs, 384-OPA, the Cornell_6K_Array_Infinium_Rice (C6AIR) and the C7AIR (McCouch et al., 2002; Orjuela et al., 2010; Thomson et al., 2012, 2017; Morales et al., 2020). Details on each of these genotyping methods are provided in Supplementary Data 1. Genotyping of the final libraries was carried out using the C7AIR SNP array (Morales et al., 2020). The “ACGT” nucleotide file from Genome Studio was exported in the PLINK format and uploaded on TASSEL GUI². The taxa (samples) were divided into six groups, each representing one of the six CSSL populations. In each library, SNPs were called based on major vs. minor allele, where it was assumed that the minor allele was inherited from the wild donor parent. SNPs with ≥39% missing data or heterozygous calls, monomorphic SNPs, and singleton SNPs were removed from each dataset. Specific SNPs were removed when it was known that one of the parents could carry either the major or the minor allele at that variant due to heterogeneity of the parental accessions. In addition to genome-wide assays, all five parental accessions in addition to CSSLs that had red or red/brown pericarp were analyzed using functional InDel markers for Rc developed by Sweeney et al. (2006) and sequencing primers for Rd described by Furukawa et al. (2007). Likewise, donors and RPs were sequenced at black hull genes, Ph and BH4 using primers described in Yu et al. (2008) and Zhu et al. (2011).

Morphological variation for 16 traits was documented for donor accessions and RPs based on observations of two to five plants per accession. Parent plants were grown in 4 inch pots in the Guterman greenhouse at Cornell University (Ithaca, NY, United States) during summer 2015 and evaluated as described in Supplementary Table S1. The phenol reaction of the grain was determined on the parents following the method of Chen et al. (2014). CSSLs were also evaluated for pericarp color, hull color, seed shattering, seed set, and delayed flowering (noted as very late flowering under long (∼15 h) days.

Custom R scripts³ were written to calculate per-library and per-line statistics and to make figures of marker distribution and graphical genotypes. Percent donor genome, percent heterozygosity, and percent missingness were computed based on the total number of informative markers per library. To map QTL in the interspecific CSSL population, stepwise regression within IciMapping v4.0 software⁴ (RSTEP-LRT option) can be implemented as described in Wang et al. (2017) and Balakrishnan et al. (2020).

Seeds from the CSSLs and the RPs, IR64 and Cybonnet, are available from the Genetics Stocks-Oryza (GSOR) center in Stuttgart, AR, United States⁵. As is standard practice, five seeds per genetic stock will be provided for researchers to amplify. For partially sterile and low-yielding CSSLs, seeds from sibling or progeny lines may be provided in lieu of the original CSSL seed. Seeds from genebank accessions of the three ORSC donors are available as follows: OrA (W1944) from the National Institute of Genetics, Japan; OrB (IRGC106148) and OrC (IRGC105567) from the International Rice Germplasm Collection in the Philippines.

Genetic relationships among a collection of 286 accessions from the ORSC and a set of 54 O. sativa control samples were examined using a rooted NJ dendrogram and a previously generated set of 113,739 GBS SNPs (Kim et al., 2016; Figure 1). Based on geographical, phenotypic and genomic diversity, as well as cross-compatibility, three ORSC accessions were selected for crossing with both Indica and Japonica elite varieties. They were (1) NIAS W1944 (hereafter referred to as OrA) from China, classified in the Institute of Genetics in Japan as O. rufipogon Griff.; (2) IRGC106148 (hereafter referred to as OrB) from Laos, identified as Oryza nivara Sharma and Shastry; and (3) IRGC105567 (hereafter referred to as OrC) from Indonesia, classified as O. rufipogon Griff (Table 1). These accessions represented wild subpopulations W6, W4/W1 admix, and W1, respectively, as reported by Kim et al. (2016). Genetically, this suggests that OrA is more closely related to the Japonica varietal group, OrB shares ancestry with Aus and Indica, and OrC represents a group of ancestors that do not share close ancestry with any O. sativa subpopulation.

Recurrent parents were selected to represent both the Japonica and the Indica varietal groups, with emphasis on elite, publicly available breeding material commonly used in public-sector breeding programs and widely grown in different parts of the world. Cultivars Cybonnet (CV-122, PI636726, GSOR301380) and IR64 (IRGC117268, GSOR301401) were selected as the RP (Table 1). Cybonnet is an elite tropical japonica from the United States adapted to the sub-tropics, and was released by the Arkansas Agricultural Experiment Station in 2004. It is early maturing, high-yielding, and semi-dwarf, with long-grain, improved milling yield and good blast resistance (Gibbons et al., 2006). IR64 is an elite Indica developed by the International Rice Research Institute (IRRI) in the Philippines and released in 1985. It is a semi-dwarf variety with high yield potential, short growth duration, enhanced resistance to diseases and insect pests and is well-known for its superior eating quality (Khush, 2005; Mackill and Khush, 2018). IR64 is adapted to the tropics and is widely grown throughout South and SE Asia. It is still commonly used as a parent for breeding of both inbred and hybrid varieties (Virmani and Kumar, 2004; Toriyama and Kazama, 2016).

The wild donors and the elite RPs were evaluated in the greenhouse for 16 characters, with emphasis on plant and seed morphological traits, flowering time, and seed set (Table 2 and Supplementary Table S1). OrA from China was readily differentiated from the other two donors with its tall, upright plant type, compact panicles, lack of shattering, lack of stolons and elbows, as well as its early flowering. OrB from Laos had short stature, abundant stolons and elbows, spreading culm habit, open panicle type and abundant shattering, with intermediate flowering time (Figure 1). OrC from Indonesia was also short, flowered later than the other donors, had open panicles that shattered readily, and very poor seed set. All donors had black hull, long awns and positive phenol reaction. Conversely, neither recurrent parent displayed awns, elbows or a shattering habit during any greenhouse or field observations, and both had a high rate of seed production. Cybonnet had erect culm angle, compact panicles and straw colored grains. IR64 had erect culm angle with intermediate panicle type and golden hull color. Cybonnet seeds showed no reaction to phenol while IR64 seeds had a positive phenol reaction whereby seeds turned a dark purple-black color when exposed to phenol.

Table 3 provides a summary of the six CSSL libraries in terms of number of lines, donor contribution, marker coverage, and percent heterozygosity. Individual libraries are comprised of 63–81 lines each, with an estimated percent donor contribution per line ranging from 2.9 to 6.9%, based on the total number of informative markers per library. Libraries constructed using OrA gave rise to populations with the largest number of spurious background introgressions per line compared to libraries constructed using OrB or OrC, though the reason for this is not clear (Figures 2, 3).

The number of informative SNPs per population ranged from 1,496 to 2,962, and was correlated with the degree of genetic relatedness between the parents. Thus, in the Cybonnet × OrA library, the number of SNPs detected was lower than in the Cybonnet × OrB or OrC libraries, and the reverse was true for the IR64 libraries. This was consistent with the ORSC subpopulation relationships reported by Kim et al. (2016) and with the phylogenetic analysis summarized in Figure 1. The average inter-marker distance is negatively correlated with marker number, but does not adequately portray the patterns of marker distribution observed in the different libraries; these are summarized in Figure 4. The large gaps in marker coverage are most pronounced in the IR64 populations (with maximum sizes of 8.3, 9.3, and 8.8 Mb in the OrA, OrB, and OrC populations, respectively), but are also apparent in the Cybonnet populations (with maximum sizes of 4.8, 3.3, and 1.9 Mb in the OrA, OrB, and OrC populations, respectively) (regions highlighted in green in Figure 4). On chromosome 4, a marker gap at 12,954,360–21,271,163 bp in the IR64 × OrA population was found to be partially shared with a gap at 14,482,427–23,294,656 bp in the IR64 × OrC and also with a gap from 20,352,551 bp until the end of the chromosome at 35,502,694 bp in the Cybonnet × OrA population. Other examples include the gaps on chromosome 8 observed from 9,289,473 to 17,516,403 bp in IR64 × OrA, 8,897,653–18,185,078 bp in IR64 × OrB, and 15,439,243–18,108,496 bp in IR64 × OrC. These regions that are devoid of informative markers are referred to as “SNP deserts,” and are not believed to be the result of major bias in the SNP assays, based on the fact that marker coverage in the Cybonnet × OrB and Cybonnet × OrC populations is fairly uniform and that apparent gaps in one population are invariably well covered regions in another. They are indicated as white spaces in the CSSL maps in Figures 2, 3. Overall, the amount of missing data after SNP filtering was very low, ranging from 0.03 to 0.32%, and the level of heterozygosity is <1% in all libraries.

We hypothesize that the “SNP deserts” may be due to one or more of the following: (1) these regions are common by state (cbs) and inferred to be common by descent (cbd) in the different cross combinations, thereby resulting in swaths of monomorphic SNPs; (2) sterility challenges led to irrecoverable introgressions, whereby target donor segments in one generation could not be propagated into future generations, resulting in an apparent lack of polymorphism in these regions in the final library dataset; and/or (3) heterogeneity in recurrent parent accessions used for crossing gave rise to polymorphic regions in final libraries that were difficult to attribute to either donor or recurrent parent and thus markers in these regions were dropped during SNP filtering (see section “Materials and Methods”). It is difficult to explicitly attribute marker gaps to cbd (hypothesis one) because the exact individuals used for crossing may not have been the same ones used for genotyping; however, given our current understanding of the genetic relationships between wild and cultivated rice subpopulations, we postulate that it is unlikely the gap on chromosome 4 shared between IR64 × OrA, IR64 × OrC, and Cybonnet × OrA libraries is a result of cbd. As an example of hypothesis two, we report efforts to maintain a line containing an introgression covering the region on chromosome 11 from 14 to 25.7 Mb in the Cybonnet × OrA population. The line was eventually lost due to sterility, and its absence from the final library led to apparent monomorphism in the genotyping dataset, explaining the marker gap in this region for that library (Figure 4). As for hypothesis three, within-line genetic heterogeneity is common in many improved varieties of rice, and is an essential characteristic of landrace and wild accessions found in genebanks (Olufowote et al., 1997). In this study, three different individuals per wild accession were used as pollen donors in crosses with multiple individuals from the two RPs. Regions in the genome where heterogenity is prevalent among individuals of an accession pose difficulties downstream in determining whether donor introgressions are present. Future attempts to carry out deep sequencing on these lines may likewise run into similar challenges in attributing alleles to donor or recurrent parent.

Depending on the specific donor segments carried by each of the CSSLs, wild characteristics are sometimes observed in the final lines. Classic domestication traits that help differentiate wild and cultivated rice include pericarp color, hull color, degree of seed shattering and flowering time. Some of these characteristics make the CSSLs difficult to grow and/or subject them to regulatory oversight, and thus knowing which lines carry “wild” traits is helpful in designing experiments to further evaluate the lines.

During development of the CSSL libraries, we observed instances of low seed set (less than 200 seeds per plant under greenhouse or field conditions) in several lines per population (Supplementary Table S2). Overall, 28 lines were noted to have low levels of fertility (four in the Cyb × OrA, six in the Cyb × OrB, four in the Cyb × OrC, four in the IR64 × OrA, six in the IR64 × OrB, and four lines in the IR64 × OrC populations). Six lines in the Cybonnet × OrB and OrC libraries had low or no seed set in the field (Lafayette, LA, United States, 2019) but it was possible to obtain some seed set from the same plants grown under greenhouse conditions in Ithaca, NY, United States. Over the course of development, we also noticed that some lines were fertile in the heterozygous condition, but we were unable to harvest seed from homozygous offspring. In some of these cases, we were able to select sib lines that contained segment(s) of interest and produced seed when selfed, but over time, it became more difficult to find lines containing specific regions of introgression; in some cases, we lost segments entirely due to sterility (Supplementary Table S3). For example, in both the Cyb × OrB and IR64 × OrB libraries, a ∼650 kb region on chromosome 1 (41,985,829–42,630,129 bp) was only maintained as a heterozygous region and the seed set was always very low in the CSSLs carrying this introgression, indicative of the presence of a sterility gene(s). This region was implicated only in the OrB libraries, as sterility was not a problem in the corresponding lines developed from the OrA or the OrC donors. At this time, we have chosen to maintain some lines as backcross populations as an insurance policy to avoid losing additional segments. A case in point is Line Cyb_OrB_4 which contains six introgressions in chromosomes 1, 3, and 11, and is currently propagated by backcrossing to Cybonnet.

There are hundreds of loci known to contribute to varying levels of gametic and zygotic incompatibilities distributed throughout the 12 chromosomes of rice as noted in Oryzabase⁷ (Kurata and Yamazaki, 2006) and the OGRO database⁶; (Yamamoto et al., 2012). Most are population-specific, and this makes it difficult to associate the low fertility phenotypes with particular genes or introgressions in the CSSLs. Our objective in documenting CSSLs with low levels of seed set (as well as very late flowering) is to alert researchers which lines may be difficult to propagate and therefore require special attention. In some cases, we also alert users of these genetic resources that seeds of heterozygous backcross progeny or bulked seeds harvested from multiple generations of pure-line sibs may be distributed for particular lines.

The C7AIR includes eight functional SNPs in its design that target important agronomic or physiological traits (Morales et al., 2020). They tag causal variants in genes conferring blast disease resistance (Pi-ta; Bryan et al., 2000), grain size (GS3; Fan et al., 2006; Takano-Kai et al., 2009), and eating and cooking quality (ALK, WAXY; Okagaki and Wessler, 1988; Wang et al., 1995; Gao et al., 2003). The parents of all six libraries were polymorphic for these variants, and CSSLs carrying donor SNPs at Pi-ta, GS3, ALK, and WAXY were documented (Supplementary Table S2). At GS3, both Cybonnet and IR64 carried the allele for long grain, while the three wild parents carried the wild type allele conferring short-medium grain length. At Pi-ta, the two RPs carried the allele for resistance while the wild parents carried the ancestral pi-ta allele conferring blast susceptibility (Huang et al., 2008). The evolution of resistance at the Pi-ta locus is believed to have been the result of a mutation that coincided with domestication, so it is not entirely surprising that none of the wild donors used in this study were found to carry the resistant allele (Lee et al., 2009; Yoshida and Miyashita, 2009). We were unable to predict grain quality phenotypes in the CSSLs because the traits determined by ALK and WAXY are quantitative, epistatic, and post-transcriptionally regulated (Bligh et al., 1998; Tian et al., 2009; Huggins et al., 2019). Evaluating these lines for novel grain and nutritional quality attributes is an exciting area for future research.

To facilitate downstream use and evaluation of our six libraries by the greater rice genetics and breeding community, we developed a simple, user-friendly R Shiny Application that enables searching for lines of interest across user-selected libraries by inputting a genomic region of interest. The web application returns a comprehensive list of individual lines across the selected libraries that harbor introgressions either fully or partially across the region of interest from the donor accessions so that users may shortlist lines that they may potentially want to evaluate. Users may access the C7Air genotype data in Supplementary Dataset 2 or may also download C7AIR genotype data directly from the application for their filtered list of lines. This application is available at http://cssl-identifier.rcac.purdue.edu/ and is hosted by Purdue University Research Computing.

The design and development of these CSSLs represents only the first step in an ambitious conceptualization of a public pre-breeding program that serves to develop, test, and deploy novel variation, with the long-term goal of enabling breeders to more readily tap into the hidden potential of wild and exotic plant genetic resources. The materials developed here are available for distribution and offer researchers new opportunities to undertake complementary multi-year, multi-location trials for yield and agronomic performance, response to abiotic and biotic stresses, and quality traits important to the rice community. To maximize the potential for wide utility and interpretation, phenotypes must be collected using controlled vocabularies, with accompanying climatic metadata and management details to provide environmental context. Further, experimental data describing phenotypic, genotypic and environmental variation should be deposited into a common database or data management system to make those data findable, accessible, interoperable and re-usable (Wilkinson et al., 2016) so that benefits can be collectively realized.

Pre-breeding is inevitably a long-term value proposition. Benefits will continue to accrue over many years if efforts are carried out in partnership with applied breeding programs dedicated to finished variety development. A partnership model that leverages public and private sector resources would enable the large scale experiments that are needed to drive recombination via rapid-cycle recurrent selection in order to overcome linkage drag and give rise to novel sources of useful variation. Favorable new materials could then be extracted by breeders periodically for crossing with elite lines, while the essential cauldron of diverse alleles were continuously managed, simmered and stirred. The potential of a well-designed, thoughtfully implemented pre-breeding program that is fueled by complementary partnerships, creative benefit-sharing agreements, and strategic deployment of resources has far-reaching consequences for plant breeders and for agricultural communities world-wide. Such a pre-breeding program can drive innovation, motivate new thinking, and unleash strategies that iteratively create value by reinvesting in the deployment of natural variation across the landscape, bringing rewards for decades to come.