Sekvenování, sestavení a anotace genomu

Sekvenovali jsme O. sinensis z okresu Nyingchi v Tibetu, Čína. Provedli jsme analýzu WGS pomocí platforem pro sekvenování nové generace Roche 454 a Illumina HiSeq 2000. Tím jsme získali čisté soubory sekvenačních dat o velikosti ~5,4 Gb, čímž jsme získali přibližně 45,1násobné pokrytí genomu (doplňková tabulka S1). Velikost genomu jsme odhadli na ~124,08 Mb a ~119,8 Mb na základě průtokové cytometrie a rozložení hloubky 17 merů sekvenovaných čtení (doplňkové obrázky S1-2 a doplňková tabulka S2). Genom O. sinensis byl nejprve sestaven z dlouhých čtení Roche 454 pomocí programu Newbler16 a následně byl sestaven scaffolding předsestavených kontigů se sekvenačními čteními Illumina mate pair pomocí programu SSPACE17. To nakonec poskytlo sestavu genomu o velikosti ~116,4 Mb, která pokrývá ~97 % odhadované velikosti genomu a obsahuje 156 scaffoldů (>2 Kb) s hodnotou ScafN50 ~3 Mb a 9 141 kontigů (N50 = 21 423 bp) (tabulka 1 a doplňkové tabulky S3-4). Pro ověření kvality sestavy genomu jsme nejprve zarovnali všechny DNA a exprimované sekvenční značky (EST) O. sinensis dostupné ve veřejných databázích a získali jsme míru mapování 98,85 %, resp. 95,33 % (Doplňková tabulka S5). Zadruhé jsme mapovali všechna čistá dlouhá čtení Roche 454 (~1,84 Gb) k sestaveným sekvencím genomu, což ukázalo téměř dokonalé zarovnání s mírou mapování 99,01 % (Doplňková tabulka S5). Zatřetí, námi sestavené transkripty vykazovaly dobré zarovnání k sestavenému genomu; z 11 742 transkriptů bylo mapováno 91,29 % (pokrytí transkriptů ≥ 80 % a identita ≥ 90 %; doplňková tabulka S5). Nakonec jsme vyhodnotili úplnost naší sestavy O. sinensis pomocí BUSCO18; 94,0 %, 4,0 % a 1,8 % z 1 315 očekávaných BUSCO konzervovaných genů Ascomycota bylo v naší sestavě O. sinensis identifikováno jako úplné, fragmentované a chybějící (Doplňková tabulka S5).

Vytvořili jsme ~15,05 Gb dat sekvenování RNA (RNA-Seq) získaných z celkem šesti knihoven reprezentujících tři hlavní vývojová stadia, která nám pomohla při predikci genů (Doplňkový obrázek S4 a Doplňkové tabulky S6,7). V kombinaci s ab initio predikcí, zarovnáním proteinů a EST, EvidenceModeler combing a dalším filtrováním jsme definovali 7 939 genů kódujících proteiny (tabulka 1 a doplňková tabulka S8). Z těchto předpovězených genů bylo možné přibližně 97,0 % a 71,51 % funkčně klasifikovat a podpořit daty RNA-Seq (Doplňkové tabulky S9-11). Pomocí BUSCO z linie Ascomycota jsme dále zjistili, že 94,4 %, 3,6 %, 1,8 % a 0,2 % genů bylo kompletních, fragmentovaných, chybějících a duplikovaných, což ukazuje na dobrou kvalitu naší genové anotace (Doplňková tabulka S11). Provedli jsme také vyhledávání homologů a anotovali nekódující RNA (ncRNA) geny, čímž jsme získali 146 genů přenosové RNA (tRNA), 33 genů ribozomální RNA (rRNA), 70 genů malé nukleární RNA (snoRNA) a 15 genů malé jaderné RNA (snRNA) (doplňkový obrázek S6 a doplňková tabulka S12). Anotace opakovaných sekvencí prezentovala, že transponovatelné elementy (TE) tvořily přibližně 74,67 % sestaveného genomu a 80,07 % nezpracovaných čtení, což znamená, že ~5,45 % nesestaveného genomu tvoří TE (Doplňkové tabulky S13-14). Obsah GC byl 43,09 % v celém genomu a 61,49 % v kódujících sekvencích (doplňkový obrázek S3; doplňkové tabulky S4 a S8). Anotovali jsme 8 918 jednoduchých sekvenčních opakování, která poskytnou cenné genetické markery pro pomoc budoucím šlechtitelským programům housenice čínské (Supplemental Tables S15-16 a Supplemental Figure S7).

Retrotranspozony řízená expanze genomu a masivní odstranění nekolineárních genů

Srovnání velikostí genomu ukázalo, že genom O. sinensis je téměř 3,4krát větší než u ostatních entomopatogenních hub z čeledi Hypocreales (Doplňková tabulka S17 a Doplňkový obrázek S8A). Analýza opakovaných sekvencí odhalila, že toto rozšíření bylo způsobeno především rychlým rozšířením transponovatelných elementů. Přibližně 74,67 % sestavy genomu O. sinensis bylo tvořeno opakovanými sekvencemi (Doplňkové tabulky S13-14), což je výjimečně více než u Metarhizium anisopliae (~0,98 %)19, Metarhizium acridum (~1,52 %)19, Cordyceps militaris (~3,04 %)20 a Beauveria bassiana (~2,03 %)21 (P < 4,822e-07) (Doplňkový obrázek S8B). Nejhojnější byly zejména elementy MULE, které tvořily ~1,6 % (~1,9 Mb) genomu O. sinensis a více než 59 % DNA transpozonů u O. sinensis. Retrotranspozony, většinou retrotranspozony s dlouhými koncovými repeticemi (LTR), tvořily ~59,76 % genomu O. sinensis a k jejich rozsáhlému rozšíření došlo přibližně před ~38 miliony let (MYA) (doplňkový obrázek S9).

V kontrastu s rychlou amplifikací retrotranspozonů LTR, které jsou hnací silou expanze genomu O. sinensis, je dalším pozoruhodným rysem dramatický úbytek protein-kódujících genů v linii O. sinensis ve srovnání s jinými entomopatogenními houbami. Ve srovnání s celkovým počtem 7 939 protein-kódujících genů u O. sinensis bylo u jiných entomopatogenních hub, např. u Metarhizium anisopliae (10 582)19, Metarhizium acridum (9 849)19, Cordyceps militaris (9 684)20, Beauveria bassiana (10 366)21 a Tolypocladium inflatum (9 998)22 , v průměru více než 1095 genů (tab. 1). Takové snížení počtu genů bylo dále prokázáno identifikací nekolineárních genů a srovnávací analýzou bloků syntézy mezi genomy O. sinensis a C. militaris. Identifikovali jsme celkem 308 syntenických bloků, které pokrývají téměř 72,7 % (~23,4 Mb u C. militaris vs. ~43,5 Mb u O. sinensis) genomu C. militaris (obr. 1A; doplňkový obrázek S10 a doplňkové tabulky S18-19). Z těchto syntenických genomových oblastí došlo u O. sinensis ke snížení počtu nekolineárních genů (2 127) ve srovnání s C. militaris (3 259), ale ke zvýšení počtu opakovaných sekvencí (23,8 Mb u O. sinensis vs. 0,40 Mb u C. militaris) (obr. 1B a doplňková tabulka S19). Funkční anotace 2 468 genů, které byly u O. sinensis ztraceny, ukázala, že se podílejí především na metabolismu aminokyselin, jako je biosyntéza aminokyselin (ko01230), metabolismus argininu a prolinu (ko00330) a metabolismus tyrosinu (ko00350) (doplňkový obrázek S11 a doplňková tabulka S20). Pozoruhodné je, že téměř 81 % opakovaných sekvencí v těchto 308 syntenických blocích tvořily LTR retrotranspozony, z nichž 40,4 % byly cikánské retroelementy (obr. 1B a doplňková tabulka S19). Molekulární datování odhadlo, že tato konkrétní třída retrotranspozonů LTR se amplifikovala ~38 Mya, což je v souladu s vyzdvižením Qinghai-Tibetské plošiny (obr. 1C).

Rychlá evoluce genových rodin souvisejících s patogenitou hub

Jednou z nejnápadnějších charakteristik genomu O. sinensis je nedostatek vysoce homologních genových párů. Z předpokládaných 7939 genů kódujících proteiny nesdílel žádný pár >90% identitu aminokyselin v kódujících sekvencích a pouze jeden pár měl >80% identitu aminokyselin (obr. 2A a doplňková tabulka S21). Tento rys byl pozorován také u blízce příbuzné C. militaris a ektomykorhizní houby Tuber melanosporum 23. Ve srovnání s jinými entomopatogenními houbami, jako jsou B. bassiana a C. militaris, byl počet multigenních rodin u O. sinensis omezený a tvořily pouze 8,7 % předpokládaného proteomu; většina genových rodin měla pouze dva členy (Supplemental Figure S12). Míra přírůstku genů byla nápadně nižší než míra jejich úbytku a ze 7 800 genových rodin nalezených u nejnovějšího společného předka (MRCA) Hypocreales jich bylo u O. sinensis zřejmě ztraceno 1 756 (obr. 2B). Takto kompaktní prostor genového kódování genomu O. sinensis naznačuje povahu této vysoce specializované houby s nízkou schopností adaptace na mnohočetné podněty prostředí.

Abychom porozuměli evoluci genových rodin, které souvisejí s patogenitou hub a adaptací vysokohorských hub na drsné prostředí, zkoumali jsme funkční vlastnosti genových rodin, které u O. sinensis prošly expanzí nebo kontrakcí. V genomu O. sinensis došlo k výraznému rozšíření genových rodin, které se podílejí především na patogenitě hub, včetně peroxidázové aktivity (PF01328; P < 0.01), serinové hydrolázy (PF03959; P < 0,01), deuterolysinové metaloproteázy (M35) peptidázy (PF02102; P < 0,01) a cytochromu P450 (PF00067; P < 0,01) (Doplňková tabulka S23). Zajímavé je, že jsme zjistili, že rozšířené rodiny genů jsou také funkčně obohaceny v kategorii Pfam glukóza-metanol-cholin (GMC) oxidoreduktázy, která se podílí na metabolismu ekdysteroidů při línání hmyzu (Doplňková tabulka S23). Při srovnání s ostatními entomopatogenními houbami bylo u linie O. sinensis rovněž pozorováno rozšíření genových rodin s nadměrným zastoupením termínů Pfam domněle souvisejících s adaptací na nízkou teplotu (PF06772; P < 0,01) (Doplňková tabulka S23).

Naproti tomu rodiny genů vykazující stav kontrakce byly zapojeny především do transportních procesů a energetického metabolismu, jako jsou ABC transportéry (PF00005; P < 0,01), permeáza aminokyselin (PF00324; P < 0,01) a ATP syntáza (PF00306; P < 0,05) (Doplňková tabulka S24). Kromě dynamické evoluce těchto genových rodin jsme u O. sinensis dále zjistili 1077 (~13,57 %) druhově specifických genů (obr. 2C). Z nich 318 (~29,53 %) genů bylo možné funkčně anotovat a byly významně obohaceny v kategoriích GO spojených s vazbou škrobu (GO: 2001070; P < 0,01), patogenezí (GO: 0009405; P < 0.01) a katabolickým procesem buněčné stěny (GO: 0016998; P < 0,01) (Doplňková tabulka S25).

Aby se hmyzí hostitelé vyhnuli infekci houbovými patogeny, často rychle produkují velké množství reaktivních forem kyslíku (ROS), které patogeny přímo ničí. Jako odpověď se u patogenů během evoluce vyvinul antioxidační obranný systém ROS, jehož peroxidasy, působící jako enzymy zachraňující ROS, jsou považovány za jednu z nejvýznamnějších a nedílných složek24, 25 . Mezi rozšířenými geny u O. sinensis byla peroxidasová aktivita jednou z vysoce obohacených funkčních kategorií (Supplemental Table S23). Vyhledávání pomocí skrytého Markovova modelu (HMM) odhalilo u O. sinensis 42 (~0,53 %) peroxidázových genů, jejichž počet byl pozoruhodně vyšší než u C. militaris (28) a kvasinek (21), což naznačuje, že dvojnásobné rozšíření peroxidázových genů může potenciálně vést k silné schopnosti napomáhat detoxikaci ROS u O. sinensis (obr. 3A a doplňková tabulka S26). Mezi těmito 42 peroxidázovými geny je nejhojněji zastoupena haloperoxidáza (hem), která tvoří ~16,67 % všech zjištěných peroxidáz (obr. 3B). Na rozdíl od jiných blízce příbuzných druhů hub, které zcela postrádají typický 2-cysteinový peroxiredoxin, si O. sinensis stále zachovává jednu kopii (obr. 3B). Již dříve bylo prokázáno, že 2-Cys peroxiredoxin hraje u Vibrio vulnificus 26 roli v reakci na různé úrovně oxidačního stresu. Srovnávací analýza odhalila, že zachovaný gen u O. sinensis patří k typu Prx1, o kterém se uvádí, že je funkčně konzervovaný27 a exprimuje se pouze tehdy, když jsou buňky vystaveny vysokým hladinám H2O2 generovaným exogenně26.

Na rozdíl od mechanismu infekce rostlinných patogenů (PP), který vyžaduje sacharidově aktivní enzymy (CAZymes) k degradaci složité buněčné stěny rostlin28, hmyzí patogeny (IP) obvykle infikují své hostitele pronikáním kutikulou29. Abychom to ověřili, porovnali jsme O. sinensis a další čtyři hmyzí patogeny (M. anisopliae, M. acridum, C. militaris a B. bassiana) se čtyřmi rostlinnými patogeny (Fusarium graminearum, Magnaporthe grisea, Grosmannia clavigera a Botrytis cinerea) (Doplňková tabulka S17). Naše výsledky ukázaly, že hmyzí patogeny mají více proteáz (v průměru 362 u IP vs. 342 u PP; P < 0,43) a proteinkináz (v průměru 151 u IP vs. 119 u PP; P < 0,0014) k rozkladu hmyzí kutikuly ve srovnání s rostlinnými patogeny (Doplňkové tabulky S27-29). Naproti tomu rostlinné patogeny obsahovaly více CAZymů než hmyzí patogeny pro degradaci rostlinné buněčné stěny (v průměru 161 v IPs vs. 231 v PPs) (Doplňkové tabulky S30-32). S vyloučením rostlinných patogenů měl O. sinensis pozoruhodně méně genů kódujících proteázy (260) než ostatní hmyzí patogeny, například M. anisopliae (437), M. acridum (361), C. militaris (355) a B. bassiana (396). Nicméně ~35 % těchto proteáz charakterizovaných u O. sinensis obsahuje signální peptid, který se s větší pravděpodobností podílel na interakcích mezi patogenem a hostitelem (doplňkové tabulky S10 a S34), což je více než u jiných entomopatogenních hub (v průměru 20 %). Podobně jako u ostatních hmyzích patogenů i u O. sinensis došlo k poklesu nebo absenci několika celulázových rodin, včetně GH7, GH45 a GH51 (Doplňková tabulka S30).

Zkoumali jsme také profily genové exprese napříč třemi vývojovými stadii O. sinensis, přičemž poměr délky houby a hmyzu dosahoval ~1,20×, ~1,75× a ~2,20×. Výsledky ukazují, že mezi třemi vývojovými stadii bylo diferenciálně exprimováno celkem 411 genů (DEG) (doplňkový obrázek S14). Funkční anotace těchto 411 DEG zjistila, že se podílejí především na patogenitě hub, jako jsou glykosylhydrolázy rodiny 16 (PF00722; FDR < 0,01), cytochrom P450 (PF00067; FDR < 0,01) a superrodina hlavních facilitátorů (PF07690; FDR < 0,05). Kromě toho byly funkčně obohaceny také geny kódující enzymy spojené s mitochondriálním dýchacím řetězcem, jako je rodina NAD dependentních epimeráz/dehydratáz (PF01370; FDR < 0,01) a BCS1 N-terminální doména (PF08740; FDR < 0,01) (doplňková tabulka S33).

Pozitivní darwinovský výběr slouží jako hnací síla pro patogenitu hub

Pozitivní výběr nepochybně sehrál rozhodující roli v evoluci různých organismů žijících ve vysokohorském prostředí Čchinghajsko-tibetské náhorní plošiny a mnoho fenotypových znaků pravděpodobně vykazuje takové selekční znaky3,4,5 . Z 1499 vysoce spolehlivých jednokopírovacích ortologů sdílených mezi O. sinensis a ostatními 12 druhy hub bylo u O. sinensis pomocí testu poměru pravděpodobnosti větví (LRT; P < 0,05) identifikováno 163 pozitivně selektovaných genů (PSG) (doplňková tabulka S35). Z nich byl jeden gen (OSIN3929; zde pojmenovaný OsPRX1) funkčně zapojen do peroxidázové aktivity (obr. 3C). Tento gen je členem rodiny peroxiredoxinů s 1-cysteinem a je vysoce homologní s PRX1 (YBL064C) v S. cerevisiae 30. Geny PRX1 u S. cerevisiae a dvou lidských patogenních hub, A. fumigatus a C. albicans, jsou funkčně konzervovány a jsou vyžadovány pro detoxikaci oxidačního výbuchu v hostitelských buňkách31, 32. Zejména delece PRX1 u známého patogenu rýže, Magnaporthe oryzae, vedla k téměř úplné ztrátě patogenity, což naznačuje, že tato peroxidasa je klíčová pro interakce mezi hostitelem a patogenem27. Pozoruhodné je, že několik genů zapojených do interakcí mezi hostitelem a patogenem, včetně peroxisomální biogeneze, proteinkinasy a metalopeptidasy, bylo rovněž zjištěno, že jsou pod pozitivní selekcí (obr. 3C).

Evoluce párovacího systému

Párovací systém je u askomycetních hub obvykle řízen lokusem párovacího typu (MAT)33 . Naše analýza sekvenování genomu zjistila, že O. sinensis má nejen gen pro typ páření MAT1-2-1 v rámci idiomorfu MAT1-2, ale také tři geny pro typ páření (tj. MAT1-1-1, MAT1-1-2 a MAT1-1-3) v rámci idiomorfu MAT1-1 (Supplemental Figure S15B). Tato vlastnost byla ověřena pomocí resekvenování celého genomu 31 přírodních populací v téměř celém geografickém areálu, což naznačuje, že O. sinensis je skutečně homothalický (Doplňkový obrázek S15A a Doplňková tabulka S36). Tato vlastnost se extrémně liší od jeho blízce příbuzných houbových patogenů, jako jsou Tolypocladium inflatum (MAT1-2)22, C. militaris (MAT1-1)20, B. bassiana (MAT1-1)21, M. anisopliae (MAT1-1)19 a M. acridum (MAT1-2)19, které jsou heterotální a mají pouze jeden lokus páření. Podobně jako u známého homothalického rostlinného patogenu Fusarium graminearum 34 odhalila organizace těchto dvou MAT lokusů u O. sinensis stav fúze v idiomorfní genomové oblasti, která byla obohacena zejména o LTR retrotransponony. Odhaduje se, že k rozdělení homothalické O. sinensis a heterothalické C. militaris došlo téměř 174,2 MYA (doplňkový obrázek S13C) a v průběhu evoluční historie došlo k několikanásobné přeměně párovacího systému ze samokompatibilního na samokompatibilní (doplňkový obrázek S15C), což připomíná vláknité askomycety rodu Neurospora 35 .

Divergence populací na základě zeměpisných šířek na Čchinghajsko-tibetské náhorní plošině

Pro zkoumání příbuzenských vztahů a divergence populací jsme shromáždili a resekvenovali 31 vzorků O. sinensis v celém jeho známém areálu rozšíření, včetně provincií Čching-chaj, S‘-čchuan, Jün-nan a Kan-su a Tibetské autonomní oblasti na Čchinghajsko-tibetské náhorní plošině (Doplňkový obrázek S16 a Doplňková tabulka S37). Celkem jsme vygenerovali 183 milionů párových čtení (~36,68 Gb sekvencí) s průměrnou hloubkou ~10,1× (nezpracovaná data) (Doplňková tabulka S38). Z těchto dat jsme vytvořili soubor 816 960 jednonukleotidových polymorfismů (SNP) a 48 092 striktních indelů (inzerce a delece) pro posouzení příbuznosti mezi populacemi O. sinensis (doplňkové obrázky S18-19 a doplňková tabulka S39). Většina genomových variant (71,1 %) byla mapována do intergenových oblastí s podskupinou mapovanou do kódujících oblastí (23,3 % sestávající z 101 997 synonymních a 88 224 nesynonymních SNP s poměrem substitucí 0,86) (doplňkový obrázek S19 a doplňková tabulka S39). Fylogenetický strom sestavený na základě souborů dat SNP rozdělil 31 vzorků do tří geograficky oddělených skupin od oblastí s nízkou nadmořskou výškou po oblasti s vysokou nadmořskou výškou (obr. 4A) – toto zjištění bylo posíleno pomocí PCA s použitím prvního a druhého vlastního vektoru (obr. 4B a doplňkový obrázek S20A). Změna počtu předpokládaných ancestrálních populací (K) ukázala, že při K = 3 jsou tři odlišné skupiny v souladu s PCA a fylogenetickou rekonstrukcí (obr. 4C a doplňkový obrázek S20B). Některé přírůstky ze skupiny z nízkých zeměpisných šířek vykazují silné známky příměsi a jsou více rozptýlené ve srovnání s ostatními dvěma skupinami, což naznačuje větší genetickou rozmanitost pravděpodobně v důsledku sdílených polymorfismů předků a/nebo nedávných introgresních událostí (obr. 4D,E). Odhad statistiky populační diferenciace (F ST ) mezi těmito třemi skupinami na základě zeměpisné šířky dále odhalil bazální charakter oblasti nízkých zeměpisných šířek, zejména populací z tibetského okresu Nyingchi, což dále dokládá jejich výrazně zvýšená nukleotidová diverzita (π) v rámci skupiny a snížená populační diferenciace s ostatními dvěma skupinami z vysokých zeměpisných šířek. (Obr. 4C-F).

Dále jsme zkoumali geny ovlivněné různými úrovněmi obsahu SNP a nesynonymních mutací (Doplňková tabulka S40). Analýza funkčního obohacení 100 genů s nejvyšším obsahem SNP a/nebo nesynonymních mutací ukazuje, že se podílejí především na metabolismu houbových sekundárních metabolitů, jako je dehydratáza polyketidové syntázy (PF14765; FDR < 0,01), doména KR (PF08659; FDR < 0,01) a kondenzační doména (PF00668; FDR < 0,01). Obohaceny byly také funkční kategorie spojené s biosyntézou mastných kyselin, jako je acyltransferázová doména (PF00698; FDR < 0,01) a beta-ketoacyl syntáza (PF00109 a PF02801; FDR < 0,01) (Doplňkové tabulky S41-42).

.