Forbedret Brassica rapa-referencegenom ved hjælp af enkeltmolekyl-sekventering og kromosomkonformationsopsamlingsteknologier

Genomsamling

For at styre genomsamlingen estimerede vi størrelsen af B. rapa-genometret ved hjælp af flowcytometri med ris som reference. Vi anslog oprindeligt, at B. rapa har en genomstørrelse på 455 Mb (Supplerende tabel S1). Yderligere undersøgelser med beregninger af den samlede længde af det konsensuskort, der er genereret på grundlag af BioNano-data, viste en genomstørrelse på 442,9 Mb (Supplerende tabel S2). Begge estimater var mindre end den tidligere rapporterede størrelse på 52915 eller 485 Mb2.

Vi samlede B. rapa-genomet ved hjælp af ~57-dobbelt dækning af PacBio-sekventeringssubstrækninger (~25.88 Gb), ~456-dobbelt dækning af BioNano-data (~207.70 Gb) og ~164-dobbelt dækning af Hi-C-reads (~74.64 Gb) Den resulterende samling bestod af 1476 contigs med en contig N50 på 1.45 Mb og en samlet længde på 351.06 Mb (tabel 1). Efterfølgende påviste vi diskrepanser inden for 22 contigs ved hjælp af Hi-C-reads (Supplerende tabel S4). I stedet for at fjerne disse contigs, opdelte vi disse ved konfliktområderne; dataene for Contig01464 er vist som et eksempel (Supplerende figur S1).

Efter stilladsering og estimering af gapstørrelser ved hjælp af BioNano-kort og mate-pair-reads (fra BRAD, http://brassicadb.org), opnåede vi 1301 stilladser med en stillads N50 på 4,44 Mb (tabel 1). For at tildele de resulterende stilladser til deres kromosomale positioner forankrede vi disse stilladser ved hjælp af Hi-C-dataene og det forbedrede genetiske kort (se Metoder). Vi forankrede 298,19 Mb sekvens på ti kromosomer, som omfattede 200 stilladser, der var grupperet ved hjælp af Hi-C-data, og 8 stilladser tildelt ved hjælp af det genetiske kort. Vores endelige samling, kaldet B. rapa-genom v3.0, omfattede i alt 353,14 Mb sekvens med 396 huller (2,08 Mb) (tabel 1). B. rapa-genomet v3.0 er længere end v1.5, men kortere end v2.5.

For at vurdere kvaliteten af B. rapa-genomet v3.0 anvendte vi forskellige datakilder. Først validerede vi fuldstændigheden af vores samling ved at søge efter centrale eukaryote gener (CEG’er) ved hjælp af CEGMA16. I alt 247 ud af 248 CEG’er var komplette, og 1 CEG var delvis, hvilket indikerer, at alle CEG’er kunne påvises i vores samling (Supplerende tabel S6). Dernæst blev genomkvaliteten testet ved at matche sekvenserne af udtrykte sekvensmærker (EST’er) af B. rapa (hentet fra dbEST på NCBI), hvilket viste, at 99,34 % af EST’erne kunne findes i det nyligt sammensatte B. rapa-genom v3.0.

Forbedring af kontiguitet

B. rapa-genom v3.0 har forbedret kontiguitet med hensyn til huller og kontigstørrelser. B. rapa-genomet v1.5 blev genereret fra Illumina-sekvenser, mens flere Illumina-reads og en relativt lille mængde PacBio-sekvensdata blev anvendt til assemblage v2.5. Disse to sammensætninger har begrænsninger på grund af deres fragmentering og lave kontiguitet (tabel 1). Ved at kombinere enkeltmolekyl-sekventering, optisk kortlægning og Hi-C-teknologi repræsenterer B. rapa-genom v3.0 en ~27-foldig (contig N50: 1446 Kb vs. 53 Kb, v2.5) og ~31-foldig (contig N50: 1446 Kb vs. 46 Kb, v1.5) forbedring af sammenhængen i forhold til de to tidligere samlinger (tabel 1). vi vurderede også størrelsen og mængden af huller i hver enkelt samling. Der var kun 396 huller i v3.0, herunder huller af kendt (122 fra BioNano og 74 fra mate-pair-data) og ukendt størrelse (190 fra Hi-C scaffold joining og 10 fra genetic map joining). Sammenlignet med de tidligere samlinger har v3.0 en ~10-foldig (5,89 Kb vs. 60,59 Kb, v2.5) og ~7-foldig (5,89 Kb vs. 40,09 Kb, v1.5) forbedring i størrelsen af huller pr. Mb i forhold til de to tidligere samlinger (tabel 1). Med hensyn til antallet af huller pr. Mb er v3.0 overlegen i forhold til henholdsvis v2.5 og v1.5 med ~23 gange (1,15 vs. 25,98, v2.5) og ~35 gange (1,15 vs. 40,09, v1.5) færre huller pr. Mb (tabel 1).

For at vurdere sammenhængen og nøjagtigheden af stilladsordneringen af de tre versioner af B. rapa-referencegenomet rekonstruerede vi først de genetiske kort baseret på de tre assemblies ved hjælp af det samme sæt resekventeringsdata af en dobbelt haploid (DH) population, der stammer fra en krydsning af to headings kinesiske kållinjer17. Derefter vurderede vi placeringen af binmarkører på de genetiske kort ved at integrere dem med de tilsvarende fysiske kort. Af de 892 binmarkører i vores samling blev 877 binmarkører (98,3 %) kortlagt i det genetiske kort. Vores samling stemte overens med det genetiske kort for 801 binmarkører (91,3 %), hvilket indikerer den høje kvalitet af v3.0 (Fig.1; Supplerende tabel S7). Vi bemærkede imidlertid, at 76 (8,7 %) binmarkører på kromosomerne A05, A08 og A09 blev kortlagt til tvetydige steder i det genetiske kort. Disse regioner indeholdt gentagne sekvenser, især ved centromeriske regioner, som beskrevet i den følgende analyse. Disse modstridende regioner blev imidlertid dækket af PacBio-reads og/eller BioNano-kort; dataene for kromosom A08 i v3.0 er vist som et eksempel (Supplerende figur S2).

Der var 1092 binmarkører på det genetiske kort fra v2.5 og 866 binmarkører på det genetiske kort fra v1.5. Vi kunne imidlertid kun kortlægge 88,7 % af binmarkørerne (969 ud af 1092) og 92,3 % af binmarkørerne (799 ud af 866) på det genetiske kort for henholdsvis v2.5 og v1.5 (Supplerende tabel S7). Vi fandt, at 15,1% af binmarkørerne (166 ud af 969) i v2.5 var diskrepante, herunder 146 binmarkører med uordnede genetiske og fysiske afstande inden for samme kromosom (intrachromosom) og 20 binmarkører med inkonsistente genetiske og fysiske afstande på forskellige kromosomer (interchromosom) (Supplerende figur S3; Supplerende tabel S7). For v1.5 var 10,0 % af binmarkørerne (80 ud af 799) diskrepante, herunder 71 binmarkører på intrakromosom og 9 binmarkører på interkromosom (Supplerende figur S4; Supplerende tabel S7). Imidlertid indeholdt v3.0 de mindst modstridende intrachromosomale binmarkører (8,7 %, 76 ud af 877 binmarkører) og ingen diskrepante interchromosomale binmarkører (Supplerende tabel S7), hvilket indikerer, at B. rapa-genom v3.0 har en højere kontiguitet end de to tidligere samlinger. Tilsammen tyder disse uafhængige valideringer på, at B. rapa-genom v3.0 har den højeste kontiguitet og den bedste rækkefølge af stilladser blandt de tre B. rapa-samlinger.

Sammenligning af genom-annotation

Vi forudsagde og annoterede genmodellerne som tidligere beskrevet6. Vi identificerede i alt 45.985 proteinkodende genmodeller i v3.0, som repræsenterede 14,74 % af genomsamlingen (tabel 1). I vores samling blev 98,75 % (45.411 af 45.985) af generne annoteret på kromosomer, og kun 1,25 % (574 af 45.985) var placeret på stilladser. De de novo annoterede gener i v3.0 blev navngivet efter standarden for genmodel-nomenklatur for Brassica-referencegener (http://www.brassica.info/info/genome_annotation.php). Antallet af genmodeller i den nye samling er højere end i v1.5 (41 020 gener), men lavere end i v2.5 (48 826 gener) (tabel 1). For yderligere at vurdere annotationens kvalitet blev der foretaget en sammenligning med annotationen af tidligere samlinger ved hjælp af BUSCO18, som er baseret på et benchmark af 1440 konserverede plantegener. Ca. 97,7 % af disse konserverede plantegener blev identificeret, og 1,7 % blev påvist som fragmenter, der blev præsenteret i v3.0 (Supplerende tabelS11).

En genomsyntenyanalyse blev udført blandt de tre assemblies ved hjælp af SynOrths19 for at identificere synteniske genpar og tandemgenarrays. Der blev i alt identificeret 2077 tandem arrays (svarende til 4963 tandemgener) i v3.0. Det samme antal tandem arrays (2077 arrays svarende til 5004 gener) blev også påvist i v1.5. En vurdering af synteny på tværs af genomet viste, at 1539 tandem arrays (svarende til 3757 gener) i v3.0 var syntene med 1494 tandem arrays (svarende til 3670 gener) i v1.5. Der blev imidlertid identificeret flere tandem arrays (3535 arrays, 8002 gener) i v2.5 (tabel 1). Vi påviste huller i regionerne med overflødige tandemgener i v2.5, mens der ikke blev fundet nogen huller i hverken 3.0 eller v1.5 (fig. 2a). Disse huller kan være resultatet af samlefejl, der er opstået ved lukning af huller ved hjælp af PacBio-reads i v2.5, hvilket igen førte til en ugyldig annotation af tandemgener. For andre tandemgener uden huller observerede vi, at enkelte gener i v3.0 og v1.5 blev annoteret som to eller flere gener i v2.5 (Fig. 2b).

Når man tager hvert tandemarray som et enkelt genlokus, var der 43 099 gener tilbage i v3.0, 44 359 gener i v2.5 og 38 093 gener i v1.5 (tabel 1). Vi udførte derefter en gensyntenyanalyse, som viste, at 39,858 gener (92,48%) i v3.0 tjente som modstykker til 40,442 (91,17%) og 35,464 gener (93,10%) i henholdsvis v2.5 og v1.5. Efter sammenligning af de annoterede gener med generne i de tidlige versioner identificerede vi 3241 versionsspecifikke gener i v3.0 sammenlignet med både v2.5 og v1.5. Af disse blev 2380 gener understøttet af beviser fra matchende mRNA-læsninger af B. rapa (fra BRAD, http://brassicadb.org/), og 2295 gener blev understøttet af proteinsekvenser af andre Brassicaceae-arter (Supplerende tabel S12). I alt blev 89,10 % (2888 af 3214) af de versionsspecifikke gener i v3.0 understøttet af mRNA-data fra B. rapa eller proteinsekvenser fra andre Brassicaceae-arter, mens kun 10,90 % (326 af 3214) af generne ikke blev understøttet.

En ny LTR-RT-udvidelsesbegivenhed identificeret i den opdaterede samling

Vi annoterede TE’er i v3.0 ved hjælp af de samme metoder som tidligere rapporteret20. I alt 235 683 TE’er blev identificeret fra 1244 familier i v3.0, og 562 unikke TE-familier blev fundet sammenlignet med v2.5 og v1.5. I v3.0 repræsenterede TE’er 37,51 % (134 Mb) af det sammensatte genom, hvilket var højere end i de tidligere sammensætninger (32,30 %, 126 Mb, v2.5; 25,44 %, 72 Mb, v1.5)2,6. I vores nye samling er de mest hyppige TE’er LTR-RT, som dækker en samlet længde på 57,64 Mb og udgør 16,32 % af det samlede genom. Ikke-LTR-RT gentagelser (LINEs og SINEs) tegner sig for 3,10 % af vores samling (Supplerende figur S5). Vi påviste DNA-transposoner svarende til 26,35 Mb, som udgør 7,46 % af den samlede genomsamling (Supplerende figur S5). En komplet liste over identificerede TE’er og gentagelser i v3.0 kan findes i Supplerende tabel S13. Desuden identificerede vi i alt 1231 miRNA’er, 1281 tRNA’er, 2865 rRNA’er og 3737 snRNA’er i B. rapa-genomet v3.0 (Supplerende tabel S19).

I vores nuværende samling annoterede vi flere LTR-RT’er (57 Mb) sammenlignet med v2.5 (44 Mb) og v1.5 (18 Mb). Vi identificerede 51.062 ikke-intakte LTR-RT’er i v3.0. Yderligere analyse afslørede, at 65,27 % (33 672 ud af 51 602) af de ikke-intakte LTR’er var placeret på de ti kromosomer, mens 34,73 % (17 922 ud af 51 602) af de ikke-intakte LTR-RT’er blev fundet på de ikke-forankrede stilladser. Ved hjælp af den samme metode6 blev i alt 13 318 intakte LTR-RT’er annoteret i v3.0. Derimod var der kun 4129 og 801 intakte LTR-RT’er i henholdsvis v2.5 og v1.56. Yderligere analyse viste, at kun 18,19 % af de intakte LTR-RT’er (2423 ud af 13 318) var placeret på de ti kromosomer, mens de fleste (81,81 %, 10 895 ud af 13 318) intakte LTR-RT’er blev fundet på de ikke-forankrede stilladser i v3.Indsættelsestiden for intakte LTR-RT’er blev beregnet som tidligere beskrevet4 , hvilket viste, at B. rapa-genomet gennemgik tre bølger af LTR-RT-ekspansion, siden det afveg fra B. oleracea (Fig. 3). Disse intakte LTR-RT’er havde en gennemsnitlig indsættelsesalder på 1,88 millioner år siden (MYA), med en medianindsættelsesalder på 1,59 MYA. Desuden fandt vi flere intakte LTR-RT’er med forskellige længder i v3.0 sammenlignet med v2.5 og v1.5 (Supplerende figur S6).

Med disse intakte LTR-RT’er blev der identificeret en ny LTR-RT-udvidelsesbegivenhed i B. rapa-genomet. Vi betegnede 3155 intakte LTR-RT-indsættelsesbegivenheder fra 0 MYA til 0.4 MYA som en “ung ekspansion” med en gennemsnitlig længde på 8135 bp og en gennemsnitlig indsættelsesdato på 0.2 MYA; 2283 intakte LTR-RT-indsættelsesbegivenheder fra 1.0 MYA til 1.4 MYA som en “medium ekspansion” med en gennemsnitlig længde på 11,902 bp og en gennemsnitlig indsættelsesdato på 1.2 MYA; og 1444 intakte LTR-RT-indsættelsesbegivenheder fra 3.0 MYA til 3,4 MYA som en “gammel ekspansion” med en gennemsnitlig længde på 9823 bp og indsættelsesdato (Fig. 3). de unge og gamle ekspansioner svarer nøje til de tidligere identificerede ekspansionsbegivenheder; den mellemstore ekspansion blev først identificeret i B. rapa-genomet og har et lignende indsættelsestidspunkt som den intakte LTR-RT ekspansionsbegivenhed i B. oleracea. Desuden blev 1778 Ty1/Copia-lignende LTR-RT’er og 4179 Ty3/Gypsy-lignende LTR-RT’er identificeret i v3.0, hvilket er meget mere end dem, der blev identificeret i de tidligere assemblies (353 Ty1/Copia og 632 Ty3/Gypsy i v2.5, 260 Ty1/Copia og 162 Ty3/Gypsy i v1.5) (Supplerende tabel S20; Supplerende figur S7, S8). Generelt var der flere Ty3/Gypsy-lignende LTR-RT’er end Ty1/Copia-lignende LTR-RT’er (Supplerende tabel S20). Sammenlignet med v2.5 og v1.5 blev Ty3/Gypsy-lignende LTR-RT’er i v3.0 tydeligvis forøget siden 5 MYA (Supplerende figur S7), mens Ty1/Copia-lignende LTR-RT’er blev forøget siden 2,2 MYA (Supplerende figur S8). Ud fra de fylogenetiske træer fandt vi, at hver gruppe af LTR-RT’er havde flere kopier i v3.0 end i v2.5 og v1.5 (Supplerende tabel S21, S22; Supplerende figur S9, S10).

Genomblokke og centromerer i B. rapa-genomet

Vi undersøgte relationerne mellem genomblokkene ved hjælp af den opdaterede samling v3.0. For at definere genomblokkene og centromererne i B. rapa-genomet v3.0 konstruerede vi først de tre subgenomer (LF, MF1 og MF2) baseret på det synteniske forhold mellem v3.0 og A. thaliana (Supplerende figur S11; Supplerende tabel S14). Vi opdagede 71 ud af de 72 (3 × 24) forventede genomiske blokke i v3.0, og de fleste af dem var arrangeret i overensstemmelse med dem, der tidligere er rapporteret i ref.10 (Fig. 4; Supplerende tabel S15). I v3.0 blev de to nye fragmenterede genomblokke F (LF) og F (MF1) identificeret på kromosomerne A01 og A05, og de blev ikke observeret i ref.10. Vi kunne ikke påvise to tidligere beskrevne, meget små genomblokke i v3.0, blok C (MF2) på kromosom A07 og blok B (MF1) på kromosom A08 (ref.10). I vores samling blev genomblokkene N/M (MF1), O/P (LF) og A/C(LF) imidlertid arrangeret på henholdsvis kromosom A01, A09 og A10, mens de var ordnet på modsatte sider i ref.10. De tre små tilstødende genomblokke (S (MF2), T(MF2) og B(MF1)) på kromosom A08 i v3.0 blev ordnet S/T/B, mens disse var ordnet som T/B/S i ref.10..

Vi sammenlignede også orienteringen af genomblokkene i v3.0 med den i ref.10. Genomblokkene W (MF2) på kromosom A02 samt G (LF) og E (LF) på kromosom A07 viste sig at være omvendte i forhold til de andre blokke, der stammer fra et enkelt ACK-kromosom. Imidlertid var orienteringen af genomblokken P (LF) på kromosom A09 og tre blokke af V i v3.0 i fremadgående retning, hvorimod disse var inverterede i ref.10. Disse resultater blev yderligere understøttet af de genetiske kort i henholdsvis v3.0 og v1.5.

Vi bestemte nøjagtigt placeringen af centromererne på alle kromosomer i v3.0. Ved at screene tidligere bestemte centromeriske gentagelsessekvenser, herunder centromeriske satellit gentagelser CentBr, CRB, TR238 og PCRBr21,22,23, identificerede vi signalerne for alle 21 palæocentromeriske regioner i v3.0, mens tre palæocentromeriske regioner ikke blev påvist i ref.10 (Fig. 4, Supplerende tabel S16). Paleocentromeranalyse viste, at de ti eksisterende B. rapa centromer alle blev arvet fra de 21 paleocentromer. I v3.0 havde centromererne af kromosomerne A01, A03, A04, A05, A06, A07 og A10 de samme tilknyttede genomblokke, der flankerede de tilsvarende centromerer, som rapporteret i ref.10 (Fig. 4). Centromeren på kromosom A02 var imidlertid placeret mellem genomblokkene P (MF2) og V (MF1), og centromeren på kromosom A09 var placeret mellem genomblokkene P (LF) og B (LF), hvorimod disse blev anset for at være paleocentromerer i ref.10 (Fig. 4). Centromeren på kromosom A08 var placeret mellem genomblokkene T (MF2) og B (MF1) i stedet for mellem genomblokkene C (MF1) og T (MF2), som rapporteret i ref.10. Desuden blev der påvist 1188 gener inden for centromeriske regioner i v3.0, mens der kun blev påvist 740 gener i ref.10 (Supplerende tabel S17).

For at vurdere vores samling med hensyn til centromererne i v3.0 analyserede vi sekvensegenskaberne i de centromeriske regioner. Vi fandt, at et betydeligt højere antal TE’er og centromer-specifikke gentagelser blev kortlagt til de centromeriske regioner end til andre dele af kromosomerne, og gentætheden og rekombinationsraten var markant lavere ved de centromeriske regioner annoteret i v3.0 (Fig. 5). Desuden blev der påvist flere centromer-specifikke gentagelser ved de centromeriske regioner i v3.0 sammenlignet med dem, der er rapporteret i ref.10 (Supplerende tabel S17).