Secvențierea, asamblarea și adnotarea genomului

Am secvențiat O. sinensis din districtul Nyingchi din Tibet, China. Am efectuat o analiză WGS cu ajutorul platformelor de secvențiere de generație următoare Roche 454 și Illumina HiSeq 2000. Aceasta a generat seturi de date de secvențe curate de ~5,4 Gb, obținându-se astfel o acoperire a genomului de aproximativ 45,1 ori, respectiv (tabelul suplimentar S1). Am estimat că dimensiunea genomului este de ~124,08 Mb și ~119,8 Mb pe baza citometriei de flux și, respectiv, a distribuției adâncimii de 17 meri a citirilor secvențiate (figurile suplimentare S1-2 și tabelul suplimentar S2). Genomul O. sinensis a fost mai întâi asamblat din citirile lungi Roche 454 cu ajutorul Newbler16, urmat de schelărirea contigurilor preasamblate cu citirile de secvențiere a perechilor mate Illumina cu ajutorul SSPACE17. Acest lucru a dus în cele din urmă la un ansamblu genomic de ~116,4 Mb care acoperă ~97% din dimensiunea estimată a genomului și conține 156 de schele (>2 Kb) cu o valoare ScafN50 de ~3 Mb și 9 141 de contigi (N50 = 21 423 pb) (tabelul 1 și tabelele suplimentare S3-4). Pentru a valida calitatea ansamblului genomului, am aliniat mai întâi tot ADN-ul și etichetele de secvență exprimată (EST) din O. sinensis disponibile în bazele de date publice și am obținut rate de cartografiere de 98,85 % și, respectiv, 95,33 % (tabelul suplimentar S5). În al doilea rând, am cartografiat toate citirile lungi Roche 454 (~1,84 Gb) cu secvențele genomului asamblat, prezentând o aliniere aproape perfectă cu o rată de cartografiere de 99,01% (tabelul suplimentar S5). În al treilea rând, transcriptele pe care le-am asamblat au prezentat o aliniere bună la genomul asamblat; din 11 742 de transcripte, 91,29% au fost cartografiate (acoperire a transcriptelor ≥80% și identitate ≥90%; tabelul suplimentar S5). În cele din urmă, am evaluat caracterul complet al ansamblului nostru O. sinensis utilizând BUSCO18; 94,0%, 4,0% și 1,8% din cele 1.315 gene conservate Ascomycota BUSCO așteptate au fost identificate ca fiind complete, fragmentate și, respectiv, lipsă în ansamblul nostru O. sinensis (tabelul suplimentar S5).

Am generat ~15,05 Gb de date de secvențiere ARN (RNA-Seq) obținute dintr-un total de șase biblioteci reprezentând cele trei stadii majore de dezvoltare pentru a ajuta la predicția genelor (figura suplimentară S4 și tabelele suplimentare S6,7). În combinație cu predicția ab initio, alinierea proteinelor și EST, pieptănarea EvidenceModeler și filtrarea ulterioară, am definit 7 939 de gene codificatoare de proteine (tabelul 1 și tabelul suplimentar S8). Dintre aceste gene prezise, aproximativ 97,0% și 71,51% au putut fi clasificate din punct de vedere funcțional și, respectiv, susținute de datele RNA-Seq (tabelele suplimentare S9-11). Utilizând BUSCO din neamul Ascomycota, am constatat, de asemenea, că 94,4 %, 3,6 %, 1,8 % și 0,2 % dintre gene au fost complete, fragmentate, lipsă și, respectiv, duplicate, ceea ce indică o bună calitate a adnotării noastre genetice (tabelul suplimentar S11). Am efectuat, de asemenea, căutări de omologi și am adnotat gene de ARN necodificatoare (ncRNA), obținând 146 de gene de ARN de transfer (ARNt), 33 de gene de ARN ribozomal (ARNr), 70 de gene de ARN nucleolar mic (snoRNAs) și 15 gene de ARN nuclear mic (snRNA) (Figura suplimentară S6 și Tabelul suplimentar S12). Adnotarea secvențelor repetate a prezentat faptul că elementele transpozabile (TE) au reprezentat aproximativ 74,67 % din genomul asamblat și 80,07 % din citirile brute, ceea ce indică faptul că ~5,45 % din genomul neasamblat este format din TE (tabelele suplimentare S13-14). Conținutul de GC a fost de 43,09 % în întregul genom și de 61,49 % în secvențele codificatoare (Figura suplimentară S3; Tabelele suplimentare S4 și S8). Am adnotat 8 918 repetări de secvențe simple care vor oferi markeri genetici valoroși pentru a sprijini viitoarele programe de ameliorare a ciupercii chinezești a omizilor (tabelele suplimentare S15-16 și figura suplimentară S7).

Expansiunea genomului determinată de retrotranspozoni și eliminarea masivă a genelor necoliniarizate

Compararea dimensiunilor genomului a arătat că genomul O. sinensis era de aproape 3,4 ori mai mare decât alte ciuperci entomopatogene din familia Hypocreales (tabelul suplimentar S17 și figura suplimentară S8A). Analiza secvențelor repetate a arătat că această expansiune s-a datorat în principal unei proliferări rapide a elementelor transpozabile. Aproximativ 74,67% din ansamblul genomului O. sinensis a fost compus din secvențe repetate (tabelele suplimentare S13-14), excepțional de mari decât cele raportate la Metarhizium anisopliae (~0,98%)19, Metarhizium acridum (~1,52%)19, Cordyceps militaris (~3,04%)20 și Beauveria bassiana (~2,03%)21 (P < 4,822e-07) (figura suplimentară S8B). Elementele MULE au fost în mod notabil cele mai abundente, reprezentând ~1,6 % (~1,9 Mb) din genomul O. sinensis și mai mult de 59 % din transpozonii ADN din O. sinensis. Retrotranspozonii, în cea mai mare parte retrotranspozoni cu repetiție terminală lungă (LTR), au cuprins ~59,76% din genomul O. sinensis și a căror proliferare pe scară largă a avut loc aproximativ la ~38 milioane de ani în urmă (MYA) (Figura suplimentară S9).

În contrast cu amplificarea rapidă a retrotranspozonilor LTR care determină expansiunea genomului O. sinensis, o altă caracteristică remarcabilă este pierderea dramatică a genelor codificatoare de proteine în descendența O. sinensis în comparație cu alte ciuperci entomopatogene. În comparație cu un total de 7 939 de gene codificatoare de proteine la O. sinensis, existau mai mult de 10 095 de gene în medie la alte ciuperci entomopatogene, de exemplu, Metarhizium anisopliae (10 582)19, Metarhizium acridum (9 849)19, Cordyceps militaris (9 684)20, Beauveria bassiana (10 366)21 și Tolypocladium inflatum (9 998)22 (tabelul 1). O astfel de reducere a numărului de gene a fost evidențiată și de identificarea genelor necoliniare și de analiza comparativă a blocurilor de sintenie între genomurile O. sinensis și C. militaris. Am identificat un total de 308 blocuri de sintenie care acoperă aproape 72,7 % (~23,4 Mb în C. militaris față de ~43,5 Mb în O. sinensis) din genomul C. militaris (Fig. 1A; Figura suplimentară S10 și tabelele suplimentare S18-19). Dintre aceste regiuni genomice sintenice, a existat o scădere a genelor necoliniare în O. sinensis (2 127) în comparație cu C. militaris (3 259), dar o creștere a secvențelor repetate (23,8 Mb în O. sinensis vs. 0,40 Mb în C. militaris) (Fig. 1B și Tabelul suplimentar S19). Adnotarea funcțională a celor 2 468 de gene care au fost pierdute în O. sinensis a arătat că acestea erau implicate în principal în metabolismul aminoacizilor, cum ar fi biosinteza aminoacizilor (ko01230), metabolismul argininei și prolinei (ko00330) și metabolismul tirozinei (ko00350) (figura suplimentară S11 și tabelul suplimentar S20). În mod special, aproape 81% din secvențele repetate în aceste 308 blocuri sintenice erau retrotranspozoni LTR, dintre care 40,4% erau retroelemente țigănești (Fig. 1B și tabelul suplimentar S19). Datarea moleculară a estimat că această clasă particulară de retrotranspozoni LTR s-a amplificat ~38 Mya, ceea ce este în concordanță cu ridicarea platoului Qinghai-Tibetan (Fig. 1C).

Evoluția rapidă a familiilor de gene legate de patogenitatea fungică

Una dintre cele mai izbitoare caracteristici ale genomului O. sinensis este lipsa perechilor de gene foarte omoloage. Din cele 7.939 de gene codificatoare de proteine prezise, nicio pereche nu a împărtășit >90% identitate de aminoacizi în secvențele codificatoare și a existat doar o singură pereche care a împărtășit >80% identitate de aminoacizi (Fig. 2A și tabelul suplimentar S21). Această caracteristică a fost observată, de asemenea, la C. militaris și la ciuperca ectomicorizală Tuber melanosporum 23, strâns înrudite. În comparație cu alte ciuperci entomopatogene, cum ar fi B. bassiana și C. militaris, familiile multigene din O. sinensis au fost limitate ca număr și au reprezentat doar 8,7 % din proteomul prezis; majoritatea familiilor de gene au avut doar doi membri (figura suplimentară S12). Rata câștigului de gene a fost izbitor de mică decât cea a pierderii de gene, iar dintre cele 7 800 de familii de gene găsite în cel mai recent strămoș comun (MRCA) al Hypocreales, 1 756 au fost aparent pierdute în O. sinensis (Fig. 2B). Un astfel de spațiu compact de codificare a genelor din genomul O. sinensis sugerează natura acestei ciuperci extrem de specializate, cu o capacitate redusă de adaptare la multiplele indicii de mediu.

Pentru a înțelege evoluția familiilor de gene care au legătură cu patogenitatea fungică și adaptarea în zonele înalte la medii dure, am investigat proprietățile funcționale ale familiilor de gene care au suferit expansiuni sau contracții la O. sinensis. Genomul O. sinensis a prezentat o expansiune considerabilă a familiilor de gene care sunt implicate în principal în patogenitatea fungică, inclusiv activitatea peroxidazei (PF01328; P < 0.01), serin hidrolază (PF03959; P < 0,01), deuterolizină metaloprotează (M35) peptidază (PF02102; P < 0,01) și citocrom P450 (PF00067; P < 0,01) (tabelul suplimentar S23). În mod interesant, am constatat că familiile de gene extinse sunt, de asemenea, îmbogățite din punct de vedere funcțional în categoria Pfam a oxidoreductazei de glucoză-metanol-colină (GMC) implicată în metabolismul ecdysteroid al mucegaiului la insecte (Tabelul suplimentar S23). În comparațiile cu alte ciuperci entomopatogene, expansiunea familiilor de gene în descendența O. sinensis a fost, de asemenea, observată cu o suprareprezentare a termenilor Pfam putativ legați de adaptarea la temperaturi scăzute (PF06772; P < 0,01) (tabelul suplimentar S23).

În schimb, familiile de gene care prezentau o stare de contracție au fost implicate în principal în procesul de transport și în metabolismul energetic, cum ar fi transportatorii ABC (PF00005; P < 0,01), aminoacizii permează (PF00324; P < 0,01) și ATP sinteza (PF00306; P < 0,05) (Tabelul suplimentar S24). În afară de evoluția dinamică a acestor familii de gene, am mai detectat 1.077 (~13,57%) gene specifice speciilor în O. sinensis (Fig. 2C). Dintre acestea, 318 (~29,53%) gene au putut fi adnotate din punct de vedere funcțional și au fost îmbogățite semnificativ în categoriile GO asociate cu legarea amidonului (GO: 2001070; P < 0,01), patogeneza (GO: 0009405; P < 0,01), patogeneza (GO: 0009405; P < 0.01), și procesul catabolic al macromoleculelor peretelui celular (GO: 0016998; P < 0,01) (Tabelul suplimentar S25).

Pentru a evita infectarea agenților patogeni fungici, gazdele insectelor produc adesea rapid o mulțime de specii reactive de oxigen (ROS) pentru a ucide direct agenții patogeni. Ca răspuns, agenții patogeni au evoluat în timpul evoluției sistemul de apărare antioxidant ROS, dintre care peroxidazele, care acționează ca enzime de captare a ROS, sunt considerate una dintre cele mai proeminente și integrale componente24, 25. Printre genele extinse la O. sinensis, activitatea peroxidazei a fost una dintre categoriile funcționale foarte îmbogățite (tabelul suplimentar S23). Căutările modelului Hidden Markov (HMM) au evidențiat 42 (~0,53 %) de gene peroxidază în O. sinensis, numărul acestora fiind remarcabil mai mare decât cel din C. militaris (28) și drojdie (21), ceea ce sugerează că o expansiune de două ori a genelor peroxidază ar putea avea ca rezultat potențial o capacitate puternică de a ajuta la detoxifierea ROS în O. sinensis (Fig. 3A și tabelul suplimentar S26). Dintre aceste 42 de gene de peroxidază, haloperoxidază (haem) este cea mai abundentă, reprezentând ~16,67 % din totalul peroxidazelor detectate (Fig. 3B). Spre deosebire de alte specii fungice înrudite îndeaproape care sunt complet lipsite de peroxiredoxina tipică de 2-cisteină, O. sinensis păstrează încă o copie (Fig. 3B). S-a demonstrat anterior că peroxiredoxina 2-Cys joacă un rol în răspunsul la diferite niveluri de stres oxidativ la Vibrio vulnificus 26. O analiză comparativă a arătat că gena reținută în O. sinensis aparține tipului Prx1, care a fost raportată ca fiind conservată funcțional27 și exprimată numai atunci când celulele sunt expuse la niveluri ridicate de H2O2 generate exogen26.

În comparație cu mecanismul de infectare a agenților patogeni ai plantelor (PP), care necesită enzime active cu carbohidrați (CAZime) pentru a degrada peretele celular complex al plantelor28, agenții patogeni ai insectelor (IP) își infectează de obicei gazdele prin penetrarea cuticulei29. Pentru a testa acest lucru, am comparat O. sinensis și alte patru insecte patogene (M. anisopliae, M. acridum, C. militaris și B. bassiana) cu cei patru agenți patogeni pentru plante (Fusarium graminearum, Magnaporthe grisea, Grosmannia clavigera și Botrytis cinerea) (tabelul suplimentar S17). Rezultatele noastre au demonstrat că agenții patogeni ai insectelor au avut mai multe proteaze (în medie, 362 în PI față de 342 în PP; P < 0,43) și proteine kinaze (în medie, 151 în PI față de 119 în PP; P < 0,0014) pentru a degrada cuticula insectelor în comparație cu agenții patogeni ai plantelor (tabelele suplimentare S27-29). În schimb, agenții patogeni ai plantelor au adăpostit mai multe CAZime decât agenții patogeni ai insectelor pentru degradarea peretelui celular al plantelor (în medie, 161 în IPs vs. 231 în PPs) (tabelele suplimentare S30-32) (tabelele suplimentare S30-32). Excluzând agenții patogeni ai plantelor, O. sinensis a avut în mod remarcabil mai puține gene care codificau proteaze (260) decât alți agenți patogeni ai insectelor, cum ar fi M. anisopliae (437), M. acridum (361), C. militaris (355) și B. bassiana (396). Cu toate acestea, ~35% din aceste proteaze caracterizate la O. sinensis conțin o peptidă semnal care este mai probabil să fi fost implicată în interacțiunile agent patogen-gazdă (tabelele suplimentare S10 și S34), care este mai mare decât cea din alte ciuperci entomopatogene (în medie, 20%). La fel ca în cazul altor insecte patogene, mai multe familii de celuloză, inclusiv GH7, GH45 și GH51, au scăzut sau au fost absente la O. sinensis (tabelul suplimentar S30).

Am examinat, de asemenea, profilurile de expresie a genelor în cele trei stadii de dezvoltare ale O. sinensis, cu rapoarte de lungime ale ciupercii față de insecte care au ajuns la ~1,20×, ~1,75× și ~2,20×. Rezultatele arată că un total de 411 gene au fost exprimate diferențiat (DEG) între cele trei stadii de dezvoltare (Figura suplimentară S14). Adnotarea funcțională a acestor 411 DEG-uri a constatat că acestea au fost implicate în principal în patogenitatea fungică, cum ar fi glicozilhidrolazele din familia 16 (PF00722; FDR < 0,01), citocromul P450 (PF00067; FDR < 0,01) și superfamilia de facilitatori majori (PF07690; FDR < 0,05). În plus, genele care codifică enzimele asociate cu lanțul respirator mitocondrial au fost, de asemenea, îmbogățite din punct de vedere funcțional, cum ar fi familia de epimeraze/dehidrataze dependente de NAD (PF01370; FDR < 0,01) și domeniul N-terminal al BCS1 (PF08740; FDR < 0,01) (Tabelul suplimentar S33).

Selecția darwinistă pozitivă servește ca forță motrice pentru patogenitatea fungică

Selecția pozitivă a jucat, fără îndoială, un rol critic în evoluția diverselor organisme care trăiesc în mediile de mare altitudine de pe platoul Qinghai-Tibetan, iar multe dintre trăsăturile fenotipice sunt susceptibile de a prezenta astfel de semnături de selecție3,4,5. Dintre cei 1 499 de orthologi monocopie cu grad ridicat de încredere împărtășiți între O. sinensis și celelalte 12 specii de ciuperci, 163 de gene selectate pozitiv (PSG) au fost identificate în O. sinensis prin utilizarea testului raportului de verosimilitate la nivelul ramurii (LRT; P < 0,05) (tabelul suplimentar S35). Dintre acestea, o genă (OSIN3929; numită aici OsPRX1) a fost implicată din punct de vedere funcțional în activitatea peroxidazei (Fig. 3C). Această genă este un membru al familiei peroxiredoxinelor cu 1-cisteină și este foarte omoloagă cu PRX1 (YBL064C) în S. cerevisiae 30. Genele PRX1 din S. cerevisiae și două ciuperci patogene pentru om, A. fumigatus și C. albicans, sunt conservate din punct de vedere funcțional și sunt necesare pentru detoxifierea exploziei oxidative în cadrul celulelor gazdă31, 32. În special, suprimarea PRX1 în bine-cunoscutul agent patogen al orezului, Magnaporthe oryzae, a dus la o pierdere aproape completă a patogenității, sugerând că această peroxidază este esențială pentru interacțiunile gazdă-patogen27. În mod remarcabil, s-a detectat că mai multe gene implicate în interacțiunile gazdă-patogen, inclusiv biogeneza peroxizomală, protein kinaza și metalopeptidazele, au fost, de asemenea, supuse unei selecții pozitive (Fig. 3C).

Evoluția sistemului de împerechere

La ciupercile ascomicete, sistemul de împerechere este, de obicei, controlat de locusul tipului de împerechere (MAT)33. Analiza noastră de secvențiere a genomului a constatat că O. sinensis nu numai că poseda gena tipului de împerechere MAT1-2-1 în cadrul idiomorfului MAT1-2, dar avea, de asemenea, trei gene tip de împerechere (și anume, MAT1-1-1-1, MAT1-1-2 și MAT1-1-3) în cadrul idiomorfului MAT1-1 (Figura suplimentară S15B). Această caracteristică a fost verificată prin resecvențierea întregului genom din 31 de populații naturale din aproape întreaga arie geografică, indicând că O. sinensis este într-adevăr homothallic (figura suplimentară S15A și tabelul suplimentar S36). Această caracteristică este extrem de diferită de cea a ciupercilor patogene strâns înrudite, cum ar fi Tolypocladium inflatum (MAT1-2)22, C. militaris (MAT1-1)20, B. bassiana (MAT1-1)21, M. anisopliae (MAT1-1)19 și M. acridum (MAT1-2)19, care sunt heterotalice și posedă doar un singur locus de împerechere. În mod similar cu bine-cunoscutul agent patogen vegetal homotalic Fusarium graminearum 34, organizarea acestor doi loci MAT la O. sinensis a dezvăluit statutul de fuziune în cadrul regiunii genomice idiomorfe, care a fost deosebit de îmbogățită în retrotranspononi LTR. Despărțirea dintre O. sinensis homotalic și C. militaris heterotalic a fost estimată la aproape 174,2 MYA (Figura suplimentară S13C) și a fost supusă unor multiple conversii ale sistemului de împerechere de la autoincompatibil la autocompatibil pe parcursul istoriei sale evolutive (Figura suplimentară S15C), semănând cu genul de ascomicete filamentoase Neurospora 35.

Diferența populațiilor în funcție de latitudinile de pe Platoul Qinghai-Tibetan

Pentru a examina relațiile la nivelul genomului și divergența populațiilor, am colectat și resecvențiat 31 de accesiuni de O. sinensis în toată aria sa de distribuție cunoscută, inclusiv în provinciile Qinghai, Sichuan, Yunnan și Gansu și în Regiunea Autonomă Tibetană de pe Platoul Qinghai-Tibetan (Figura suplimentară S16 și Tabelul suplimentar S37). Am generat un total de 183 de milioane de lecturi paired-end (~36,68 Gb de secvențe) cu o adâncime medie de ~10,1× (date brute) (tabelul suplimentar S38). Din aceste date, am generat un set de 816 960 de polimorfisme cu un singur nucleotid (SNP) și 48 092 de indeluri stricte (inserții și ștergeri) pentru a evalua înrudirea între populațiile de O. sinensis (figurile suplimentare S18-19 și tabelul suplimentar S39). Majoritatea variantelor genomice (71,1 %) au fost cartografiate în regiunile intergenice, cu un subset cartografiat în regiunile de codificare (23,3 % constând în 101 997 de SNP sinonime și 88 224 de SNP nesinonime cu un raport de substituție de 0,86) (Figura suplimentară S19 și tabelul suplimentar S39). Arborele filogenetic construit cu ajutorul seturilor de date SNP a împărțit cele 31 de accesiuni în trei grupuri separate din punct de vedere geografic, de la regiunile de latitudine joasă la cele de latitudine înaltă (Fig. 4A) – o constatare care a fost întărită de PCA utilizând primul și al doilea vector propriu (Fig. 4B și Figura suplimentară S20A). Variația numărului de populații ancestrale presupuse (K) a arătat că, atunci când K = 3, cele trei grupuri distincte sunt în concordanță cu PCA și reconstrucția filogenetică (Fig. 4C și Figura suplimentară S20B). Unele accesiuni din grupul de latitudini joase prezintă dovezi puternice de amestec și sunt mai dispersate în comparație cu celelalte două grupuri, indicând o diversitate genetică mai mare, posibil datorită polimorfismelor ancestrale comune și/sau evenimentelor recente de introgresie (Fig. 4D,E). Statistica estimată de diferențiere a populației (F ST ) între aceste trei grupuri bazate pe latitudine a evidențiat și mai mult natura bazală a regiunii de latitudine joasă, în special a populațiilor din districtul Nyingchi din Tibet, care a fost evidențiată și mai mult prin diversitatea nucleotidă substanțial ridicată (π) în cadrul grupului și diferențierea redusă a populației cu celelalte două grupuri de latitudine înaltă. (Fig. 4C-F).

Am investigat în continuare genele afectate de diferite niveluri de conținut de SNP și mutații nesinonime (tabelul suplimentar S40). Analiza de îmbogățire funcțională a primelor 100 de gene cu cel mai mare conținut de SNP și/sau mutații nesinonime arată că acestea sunt implicate în principal în metabolismul metaboliților secundari fungici, cum ar fi dehidrataza polichetidei sintetazei (PF14765; FDR < 0,01), domeniul KR (PF08659; FDR < 0,01) și domeniul de condensare (PF00668; FDR < 0,01). Categoriile funcționale asociate cu biosinteza acizilor grași au fost, de asemenea, îmbogățite, cum ar fi domeniul acil-transferazei (PF00698; FDR < 0,01) și beta-cetoacil-sintetaza (PF00109 și PF02801; FDR < 0,01) (Tabelele suplimentare S41-42).