Larvesvampen Ophiocordyceps sinensis’ genom giver indsigt i højlandstilpasning af svampes patogenicitet
Genomsekventering, samling og annotation
Vi har sekventeret O. sinensis fra Nyingchi-distriktet i Tibet, Kina. Vi udførte en WGS-analyse med næste generations sekventeringsplatforme Roche 454 og Illumina HiSeq 2000. Dette genererede rene sekvensdatasæt på henholdsvis ~5,4 Gb og gav således en ca. 45,1-dobbelt genomdækning (Supplerende tabel S1). Vi anslog, at genomets størrelse er ~124,08 Mb og ~119,8 Mb baseret på henholdsvis flowcytometri og 17-mer-dybdefordeling af sekventerede læsninger (Supplerende figurer S1-2 og Supplerende tabel S2). O. sinensis-genomet blev først samlet fra Roche 454 long reads ved hjælp af Newbler16, efterfulgt af scaffolding af præ-assemblerede contigs med Illumina mate pair sekventeringslæsninger ved hjælp af SSPACE17. Dette gav i sidste ende en genomsamling på ~116,4 Mb, der dækker ~97 % af den anslåede genomstørrelse og indeholder 156 scaffolds (>2 Kb) med en ScafN50-værdi på ~3 Mb og 9 141 contigs (N50 = 21 423 bp) (Tabel 1 og Supplerende tabeller S3-4). For at validere kvaliteten af genomsamlingen, justerede vi først alt DNA og udtrykte sekvensmærker (ESTs) af O. sinensis, der er tilgængelige i de offentlige databaser, og opnåede kortlægningsrater på henholdsvis 98,85% og 95,33% (Supplerende tabel S5). For det andet kortlagde vi alle rene Roche 454 long reads (~1,84 Gb) til de sammensatte genomsekvenser, hvilket viste en næsten perfekt tilpasning med en kortlægningsrate på 99,01 % (Supplerende tabel S5). For det tredje viste de transkripter, vi samlede, en god tilpasning til det samlede genom; af 11,742 transkripter blev 91,29% kortlagt (transkriptdækning ≥80% og identitet ≥90%; Supplerende tabel S5). Endelig evaluerede vi fuldstændigheden af vores O. sinensis-samling ved hjælp af BUSCO18; 94.0%, 4.0% og 1.8% af de 1,315 forventede Ascomycota BUSCO-konserverede gener fra Ascomycota blev identificeret som henholdsvis komplette, fragmenterede og manglende i vores O. sinensis-samling (Supplerende tabel S5).
Vi genererede ~15.05 Gb RNA-sekventering (RNA-Seq) data opnået fra i alt seks biblioteker, der repræsenterer de tre vigtigste udviklingsstadier, til hjælp for genprædiktion (Supplerende figur S4 og Supplerende tabeller S6,7). I kombination med ab initio-prædiktion, protein- og EST-aligneringer, EvidenceModeler-kæmning og yderligere filtrering definerede vi 7 939 proteinkodende gener (tabel 1 og supplerende tabel S8). Af disse forudsagte gener kunne ca. 97,0 % og 71,51 % klassificeres funktionelt og understøttes af henholdsvis RNA-Seq-data (Supplerende tabeller S9-11). Ved hjælp af BUSCO fra Ascomycota-lineage fandt vi endvidere, at 94,4 %, 3,6 %, 1,8 % og 0,2 % af generne var henholdsvis komplette, fragmenterede, manglende og duplikerede, hvilket indikerer en god kvalitet af vores genannotation (Supplerende tabel S11). Vi udførte også homologsøgninger og annoterede ikke-kodende RNA (ncRNA) gener, hvilket gav 146 transfer RNA (tRNA) gener, 33 ribosomale RNA (rRNA) gener, 70 små nukleare RNA (snoRNAs) gener og 15 små nukleare RNA (snRNA) gener (Supplerende figur S6 og Supplerende tabel S12). Annotationen af gentagelsessekvenser viste, at transposable elementer (TE’er) tegnede sig for ca. 74,67 % af det samlede genom og 80,07 % af de rå læsninger, hvilket indikerer, at ~5,45 % af det ikke-samlede genom består af TE’er (Supplerende tabeller S13-14). GC-indholdet var 43,09 % i hele genomet og 61,49 % i de kodende sekvenser (Supplerende figur S3; Supplerende tabeller S4 og S8). Vi annoterede 8.918 simple sekvensgennemgange, som vil give værdifulde genetiske markører til hjælp for fremtidige avlsprogrammer for kinesisk larvesvamp (Supplerende tabeller S15-16 og Supplerende figur S7).
Retrotransposon-drevet genomudvidelse og massiv fjernelse af ikke-kollineære gener
Sammenligning af genomstørrelser viste, at O. sinensis-genomet var næsten 3,4 gange større end andre entomopatogene svampe i Hypocreales-familien (Supplerende tabel S17 og Supplerende figur S8A). Gentagelsessekvensanalyse viste, at denne udvidelse primært skyldtes en hurtig spredning af transposable elementer. Ca. 74,67 % af O. sinensis-genomprogrammet bestod af gentagne sekvenser (Supplerende tabeller S13-14), hvilket er usædvanligt større end dem, der er rapporteret i Metarhizium anisopliae (~0,98 %)19, Metarhizium acridum (~1,52 %)19, Cordyceps militaris (~3,04 %)20 og Beauveria bassiana (~2,03 %)21 (P < 4,822e-07) (Supplerende figur S8B). MULE-elementerne var især de mest hyppige og tegnede sig for ~1,6% (~1,9 Mb) af O. sinensis-genomet og mere end 59% af DNA-transposonerne i O. sinensis. Retrotransposoner, hovedsagelig lang-terminale gentagelser (LTR) retrotransposoner, omfattede ~59,76% af O. sinensis-genomet, og udbredelse i stor skala af disse fandt sted for ca. ~38 millioner år siden (MYA) (Supplerende figur S9).
I modsætning til den hurtige amplifikation af LTR retrotransposoner, der driver ekspansionen af O. sinensis-genomet, er et andet bemærkelsesværdigt træk det dramatiske tab af protein-kodende gener i O. sinensis-linjen i sammenligning med andre entomopatogene svampe. Sammenlignet med i alt 7 939 proteinkodende gener i O. sinensis var der i gennemsnit mere end 10 095 gener i andre entomopatogene svampe, f.eks. Metarhizium anisopliae (10 582)19, Metarhizium acridum (9 849)19, Cordyceps militaris (9 684)20, Beauveria bassiana (10 366)21 og Tolypocladium inflatum (9 998)22 (tabel 1). En sådan reduktion af genantallet blev yderligere påvist ved identifikation af ikke-kollineære gener og sammenlignende analyse af syntenyblokke mellem genomerne af O. sinensis og C. militaris. Vi identificerede i alt 308 synteniske blokke, der dækker næsten 72,7 % (~23,4 Mb i C. militaris vs. ~43,5 Mb i O. sinensis) af C. militaris-genomet (Fig. 1A; Supplerende figur S10 og Supplerende tabeller S18-19). Af disse synteniske genomiske regioner var der et fald i antallet af ikke-kollineære gener i O. sinensis (2.127) sammenlignet med C. militaris (3.259), men en stigning i antallet af gentagne sekvenser (23,8 Mb i O. sinensis vs. 0,40 Mb i C. militaris) (Fig. 1B og Supplerende tabel S19). Funktionel annotation af de 2.468 gener, der gik tabt i O. sinensis, viste, at de hovedsageligt var involveret i aminosyremetabolisme, såsom biosyntese af aminosyrer (ko01230), arginin- og prolinmetabolisme (ko00330) og tyrosinmetabolisme (ko00350) (Supplerende figur S11 og Supplerende tabel S20). Det er bemærkelsesværdigt, at næsten 81% af gentagelsessekvenserne i disse 308 synteniske blokke var LTR retrotransposoner, hvoraf 40,4% var Gypsy retroelementer (Fig. 1B og Supplerende tabel S19). Molekylær datering anslog, at denne særlige klasse af LTR retrotransposoner blev amplificeret ~38 Mya, hvilket er i overensstemmelse med opløftningen af Qinghai-Tibetan Plateau (Fig. 1C).
Hurtig evolution af genfamilier relateret til svampes patogenicitet
Et af de mest slående kendetegn ved O. sinensis-genomet er manglen på meget homologe genpar. Af de forudsagte 7.939 proteinkodende gener delte ingen par >90% aminosyreidentitet i de kodende sekvenser, og der var kun ét par, der delte >80% aminosyreidentitet (Fig. 2A og Supplemental Table S21). Dette træk blev også observeret i den nært beslægtede C. militaris og ektomykorrhizasvampen Tuber melanosporum 23. Sammenlignet med andre entomopatogene svampe som B. bassiana og C. militaris var multigenfamilier i O. sinensis begrænsede i antal og udgjorde kun 8,7 % af det forudsagte proteom; de fleste genfamilier havde kun to medlemmer (Supplerende figur S12). Gentilvæksten var markant lavere end gentabet, og blandt de 7.800 genfamilier, der blev fundet i den seneste fælles forfader (MRCA) af Hypocreales, blev 1.756 tilsyneladende tabt i O. sinensis (Fig. 2B). Et så kompakt genkodningsrum i O. sinensis-genomet tyder på karakteren af denne højt specialiserede svamp med en lav kapacitet til at tilpasse sig flere miljømæssige signaler.
For at forstå udviklingen af genfamilier, der er relateret til svampes patogenicitet og højlandstilpasning til barske miljøer, undersøgte vi de funktionelle egenskaber af genfamilier, der har undergået udvidelser eller sammentrækninger i O. sinensis. O. sinensis-genomet viste en betydelig ekspansion af genfamilier, der hovedsagelig er involveret i svampepatogenitet, herunder peroxidaseaktivitet (PF01328; P < 0.01), serinhydrolase (PF03959; P < 0,01), deuterolysinmetalloprotease (M35) peptidase (PF02102; P < 0,01) og cytokrom P450 (PF00067; P < 0,01) (Supplerende tabel S23). Interessant nok fandt vi, at de udvidede genfamilier også er funktionelt beriget i Pfam-kategorien af glukose-methanol-cholin (GMC) oxidoreduktase, der er involveret i ecdysteroidmetabolismen i moltning hos insekter (Supplerende tabel S23). I sammenligninger med andre entomopatogene svampe blev genfamilieudvidelsen i O. sinensis-linjen også observeret med overrepræsentation af Pfam-termer, der putativt er relateret til tilpasningen til lav temperatur (PF06772; P < 0,01) (Supplerende tabel S23) (Supplemental Table S23).
I modsætning hertil var genfamilier, der udviste kontraktionsstatus, hovedsagelig involveret i transportprocessen og energistofskiftet, såsom ABC-transportører (PF00005; P < 0,01), aminosyrepermease (PF00324; P < 0,01) og ATP-syntase (PF00306; P < 0,05) (Supplerende tabel S24) (Supplemental Table S24). Bortset fra dynamisk evolution af disse genfamilier påviste vi yderligere 1,077 (~13.57%) artsspecifikke gener i O. sinensis (Fig. 2C). Af dem kunne 318 (~29.53%) gener annoteres funktionelt og var signifikant beriget i GO-kategorier forbundet med stivelsesbinding (GO: 2001070; P < 0.01), patogenese (GO: 0009405; P < 0.01), patogenese (GO: 0009405; P < 0.01) og cellevægmakromolekyl-katabolisk proces (GO: 0016998; P < 0.01) (Supplerende tabel S25).
For at undgå infektion af svampepatogener producerer insektværter ofte hurtigt masser af reaktive oxygenarter (ROS) for direkte at dræbe patogener. Som et svar udviklede patogenerne ROS-antioxidantforsvarssystemet i løbet af evolutionen, hvoraf peroxidaser, der fungerer som ROS-scavenging-enzymer, betragtes som en af de mest fremtrædende og integrerede komponenter24, 25. Blandt de udvidede gener i O. sinensis var peroxidaseaktivitet en af de stærkt berigede funktionelle kategorier (Supplerende tabel S23). Søgninger efter skjulte Markov-modeller (HMM) afslørede 42 (~0,53%) peroxidase gener i O. sinensis, hvis antal var bemærkelsesværdigt større end i C. militaris (28) og gær (21), hvilket tyder på, at en dobbelt udvidelse af peroxidase gener potentielt kan resultere i en stærk kapacitet til at hjælpe med ROS-afgiftning i O. sinensis (Fig. 3A og Supplerende tabel S26). Blandt disse 42 peroxidasegener er haloperoxidase (haem) den mest hyppige og udgør ~16,67 % af de samlede peroxidaser, der blev påvist (Fig. 3B). I modsætning til andre nært beslægtede svampearter, der helt mangler det typiske 2-Cystein-peroxiredoxin, bevarer O. sinensis stadig en kopi (Fig. 3B). 2-Cys peroxiredoxin blev tidligere vist at spille en rolle i responsen på forskellige niveauer af oxidativ stress i Vibrio vulnificus 26. En sammenlignende analyse viste, at det bevarede gen i O. sinensis tilhører Prx1-typen, som blev rapporteret at være funktionelt bevaret27 og kun udtrykkes, når cellerne udsættes for høje niveauer af H2O2 genereret eksogent26.
I modsætning til infektionsmekanismen hos plantepatogener (PP’er), som kræver kulhydrataktive enzymer (CAZymes) til nedbrydning af den komplekse plantecellevæg28, inficerer insektpatogener (IP’er) typisk deres værter ved at trænge igennem kutikula29. For at teste dette sammenlignede vi O. sinensis og fire andre insektpatogener (M. anisopliae, M. acridum, C. militaris og B. bassiana) med de fire plantepatogener (Fusarium graminearum, Magnaporthe grisea, Grosmannia clavigera og Botrytis cinerea) (Supplemental Table S17) (Supplerende tabel S17). Vores resultater viste, at insektpatogener havde flere proteaser (i gennemsnit 362 i IP’er vs. 342 i PP’er; P < 0,43) og proteinkinaser (i gennemsnit 151 i IP’er vs. 119 i PP’er; P < 0,0014) til nedbrydning af insektkutikula sammenlignet med plantepatogener (Supplerende tabeller S27-29). I modsætning hertil havde plantepatogener flere CAZymer end insektpatogener til nedbrydning af plantecellevæggen (i gennemsnit 161 i IP’er vs. 231 i PP’er) (Supplerende tabeller S30-32) (Supplerende tabeller S30-32). Hvis man ser bort fra plantepatogener, havde O. sinensis bemærkelsesværdigt nok færre gener, der kodede for proteaser (260) end andre insektpatogener, såsom M. anisopliae (437), M. acridum (361), C. militaris (355) og B. bassiana (396). Imidlertid indeholder ~35% af disse proteaser, der er karakteriseret i O. sinensis, et signalpeptid, som med større sandsynlighed har været involveret i interaktioner mellem patogen og vært (Supplerende tabeller S10 og S34), hvilket er større end i andre entomopatogene svampe (i gennemsnit 20%). I lighed med de andre insektpatogener faldt flere cellulasefamilier, herunder GH7, GH45 og GH51, også eller var fraværende i O. sinensis (Supplerende tabel S30).
Vi undersøgte også genekspressionsprofiler på tværs af de tre udviklingsstadier af O. sinensis med længdeforhold af svamp vs. insekt, der nåede ~1.20×, ~1.75× og ~2.20×. Resultaterne viser, at i alt 411 gener blev differentielt udtrykt (DEG) mellem de tre udviklingsstadier (Supplerende figur S14). Funktionel annotation af disse 411 DEG’er viste, at de hovedsagelig var involveret i svampepatogenitet, såsom glycosylhydrolaser familie 16 (PF00722; FDR < 0,01), cytokrom P450 (PF00067; FDR < 0,01) og major facilitator superfamily (PF07690; FDR < 0,05). Desuden var gener, der kodede for enzymer forbundet med mitokondriel respiratorisk kæde, også funktionelt beriget, såsom NAD-afhængig epimerase/dehydratase-familie (PF01370; FDR < 0,01) og BCS1 N-terminalt domæne (PF08740; FDR < 0,01) (Supplerende tabel S33).
Positiv darwinistisk udvælgelse tjener som drivkræfter for svampes patogenicitet
Positiv udvælgelse har utvivlsomt spillet en afgørende rolle i udviklingen af forskellige organismer, der lever i højhøjdemiljøer på Qinghai-Tibetan Plateau, og mange af de fænotypiske træk viser sandsynligvis sådanne udvælgelsessignaturer3,4,5. Af de 1,499 højkonfidens single-copy ortologer, der deles mellem O. sinensis og de andre 12 svampearter, blev 163 positivt udvalgte gener (PSG’er) identificeret i O. sinensis ved hjælp af branch-site likelihood ratio test (LRT; P < 0,05) (Supplerende tabel S35) (Supplemental Table S35). Af dem er et gen (OSIN3929; her kaldet OsPRX1) blevet funktionelt impliceret i peroxidaseaktivitet (Fig. 3C). Dette gen er et medlem af peroxiredoxin-familien med 1-cystein og er meget homologt med PRX1 (YBL064C) i S. cerevisiae 30. PRX1-gener i S. cerevisiae og to humanpatogene svampe, A. fumigatus og C. albicans, er funktionelt bevarede og er nødvendige for detoxificering af den oxidative burst i værtscellerne31, 32. Navnlig resulterede sletning af PRX1 i det velkendte rispatogen, Magnaporthe oryzae, i et næsten fuldstændigt tab af patogenicitet, hvilket tyder på, at denne peroxidase er nøglen til interaktioner mellem vært og patogen27. Påfaldende nok blev flere gener, der er involveret i vært-patogen-interaktioner, herunder peroxisomal biogenese, proteinkinase og metallopeptidaser, også påvist at være under positiv selektion (Fig. 3C).
Evolution af parringssystemet
I ascomycetøse svampe kontrolleres parringssystemet normalt af parringstype (MAT)-lokus33. Vores genomsekventeringsanalyse viste, at O. sinensis ikke kun besad MAT1-2-1 parringstypegenet inden for MAT1-2 idiomorfen, men også havde tre parringstypegener (dvs. MAT1-1-1, MAT1-1-1, MAT1-1-2 og MAT1-1-3) inden for MAT1-1 idiomorfen (Supplemental Figure S15B). Dette træk blev verificeret ved hjælp af gensekventering af hele genomet af 31 naturlige populationer over næsten hele det geografiske område, hvilket indikerer, at O. sinensis faktisk er homothalsk (Supplerende figur S15A og Supplerende tabel S36). Denne egenskab er ekstremt forskellig fra dens nært beslægtede svampepatogener, såsom Tolypocladium inflatum (MAT1-2)22, C. militaris (MAT1-1)20, B. bassiana (MAT1-1)21, M. anisopliae (MAT1-1)19 og M. acridum (MAT1-2)19, som er heterotalliske og kun besidder et enkelt parringstypelokus. I lighed med det velkendte homotalliske plantepatogen Fusarium graminearum 34 afslørede organiseringen af disse to MAT-loci i O. sinensis fusionsstatus inden for den idiomorfe genomiske region, som var særlig beriget med LTR retrotransponsoner. Splittelsen mellem den homotalliske O. sinensis og den heterotalliske C. militaris blev anslået til at finde sted næsten 174,2 MYA (Supplerende figur S13C) og blev udsat for flere konverteringer af parringssystemet fra selvinkompatibel til selvkompatibel i løbet af dens udviklingshistorie (Supplerende figur S15C), hvilket ligner den filamentøse ascomycete-slægt Neurospora 35.
Populationsdivergens baseret på breddegrader på Qinghai-Tibetan Plateau
For at undersøge de genom-dækkende relationer og populationsdivergens indsamlede og resekventerede vi 31 accessioner af O. sinensis på tværs af dens kendte udbredelsesområde, herunder Qinghai, Sichuan, Yunnan og Gansu-provinserne og den Tibet Autonome Region på Qinghai-Tibetan Plateau (Supplerende figur S16 og Supplerende tabel S37). Vi genererede i alt 183 millioner parvis afsluttede læsninger (~36,68 Gb sekvenser) med en gennemsnitlig dybde på ~10,1× (rå data) (Supplerende tabel S38) (Supplemental Table S38). Ud fra disse data genererede vi et sæt af 816 960 enkeltnukleotidpolymorfismer (SNP’er) og 48 092 strenge indels (indsættelser og sletninger) for at vurdere slægtskab mellem populationer af O. sinensis (Supplerende figurer S18-19 og Supplerende tabel S39). Størstedelen af genomiske varianter (71,1 %) blev kortlagt til intergeniske regioner med en delmængde kortlagt til de kodende regioner (23,3 % bestående af 101.997 synonyme og 88.224 nonsynonyme SNP’er med en substitutionsratio på 0,86) (Supplerende figur S19 og Supplerende tabel S39). Det fylogenetiske træ, der blev konstrueret ved hjælp af SNP-datasættene, opdelte de 31 accessioner i tre geografisk adskilte grupper, der spænder fra regioner med lav breddegrad til regioner med høj breddegrad (Fig. 4A) – et resultat, der blev forstærket af PCA ved hjælp af den første og anden egenvektor (Fig. 4B og Supplerende figur S20A). Variation af antallet af formodede forfødte populationer (K) viste, at når K = 3, er de tre forskellige grupper i overensstemmelse med PCA og den fylogenetiske rekonstruktion (Fig. 4C og Supplerende figur S20B). Nogle accessioner fra gruppen på lave breddegrader udviser stærke tegn på blanding og er mere spredt sammenlignet med de to andre grupper, hvilket indikerer større genetisk diversitet, muligvis på grund af fælles forfædres polymorfismer og/eller nylige introgressionshændelser (Fig. 4D,E). Den estimerede befolkningsdifferentieringsstatistik (F ST ) blandt disse tre breddegrader-baserede grupper afslørede yderligere den basale karakter af det lave breddeområde, især populationer fra Nyingchi-distriktet i Tibet, hvilket yderligere blev bevist af dens væsentligt forhøjede nukleotiddiversitet (π) inden for gruppen og lavere befolkningsdifferentiering med de to andre grupper på høje breddegrader. (Fig. 4C-F).
Vi undersøgte yderligere de gener, der blev påvirket af forskellige niveauer af SNP-indhold og ikke-synonyme mutationer (Supplerende tabel S40). Funktionel berigelsesanalyse af de 100 øverste gener med det højeste SNP-indhold og/eller ikke-synonyme mutationer viser, at de hovedsageligt er involveret i metabolismen af svampes sekundære metabolitter, såsom polyketidsyntase dehydratase (PF14765; FDR < 0,01), KR-domæne (PF08659; FDR < 0,01) og kondensationsdomæne (PF00668; FDR < 0,01). De funktionelle kategorier, der var forbundet med fedtsyrebiosyntese, blev også beriget, såsom acyltransferase-domæne (PF00698; FDR < 0,01) og beta-ketoacylsyntase (PF00109 og PF02801; FDR < 0,01) (Supplerende tabeller S41-42).