芒果炭疽病抗感品種全基因組重測序分析

刁興旺， 吳莉君， 何 紅， 姚全勝， 柳 鳳

(1.廣東海洋大學濱海農業學院，廣東湛江 524088； 2.中國熱帶農業科學院南亞熱帶作物研究所/海南省熱帶園藝產品采后生理與保鮮重點實驗室/農業部熱帶果樹生物學重點實驗室，廣東湛江 524091)

芒果(MangiferaindicaL.)屬漆樹科(Anacradiaceae)芒果屬(Mangifera)常綠喬本，原產于印度和東南亞，目前主要分布在熱帶和亞熱帶地區。芒果樹具有速生易長、抗逆性強、結果早、產量高、易栽培管理、經濟壽命長等優點，在我國熱帶農業經濟中占有重要地位[1]，但芒果炭疽病對產業危害嚴重。芒果炭疽病是芒果產業中普遍且危害性最大的一種真菌病害，其病原菌主要是膠孢炭疽菌[Colletotrichumgloeosporioides(Penz.) Sacc.]。該病害主要危害芒果的新梢、嫩葉、花和果實，也是芒果最嚴重的采后病害之一，在貯運期會引起果實的腐爛，占果實總病害量的70%以上，嚴重影響芒果的品質和經濟價值[2]。選育抗病品種是防控芒果炭疽病經濟有效的方法。與傳統選育技術相比，分子標記輔助選擇育種以及基因組輔助選擇育種等技術是基于全基因組范圍內發掘變異位點，在果樹遺傳育種中不僅選擇準確性高而且不受環境影響，可以極大提高選擇效率并縮短育種年限[3]。

全基因重測序可為分子標記等育種技術提供大量突變位點信息，隨著測序技術的發展及成本降低，已有至少20個果樹全基因組測序并公布基因草圖[4]，包括蘋果(Malusdomestica)[5]、梨(Pyrusspp)[6]、棗(ZiziphusjujubaMill.)[7]、龍眼(DimocarpuslonganLour.)[8]、柑橘(CitrusreticulataBlanco)[9-10]、芒果[11]等。以此為基礎開展的全基因組重測序分析也在果樹重要經濟性狀研究領域取得重要成果。2020年高立志教授團隊選用芒果栽培種紅象牙(雜合度約為1.8%)，結合單分子實時測序(SMRT)、Illumina測序技術和遺傳圖譜，生成了一個包含 34 529 個預測蛋白質編碼基因 371.6 MB 大小的芒果基因組，并公開發布了染色體級別的芒果參考基因組序列，為基于全基因組重測序技術發掘芒果全基因組范圍內重要性狀相關的單核苷酸多態性(SNP)、插入/缺失(InDel)位點和進一步開發相關分子標記進行輔助育種奠定了可行的基礎。

目前，關于芒果炭疽病抗感品種遺傳群體構建以及全基因組重測序相關分析尚未見報道。本研究選取芒果炭疽病高感品種愛文和炭疽病高抗品種金煌，通過全基因組重測序，并以品種“紅象牙”的基因組為參考比對后分析檢測SNP、InDel位點，對多態性進行差異分析及功能注釋。一方面可用于后期開發多態性分子標記，芒果抗炭疽病品種快速鑒定、種質資源評價及遺傳作圖選擇利用；另一方面可為全基因組關聯分析、功能基因定位及分子標記輔助育種提供研究基礎。

1 材料與方法

1.1 試驗材料

2019年秋季，在中國熱帶農業科學院南亞熱帶作物研究所(110.28°E，21.17°N)將采自我國芒果主產區廣東省、廣西壯族自治區、云南省、四川省的芒果炭疽病高感品種愛文和高抗品種金煌的葉片進行初步處理。選取健康嫩葉送至廣州基迪奧生物科技有限公司進行全基因組重測序。

1.2 愛文和金煌全基因組重測序

經質檢合格的樣品用于構建DNA文庫，用FASTP (版本0.18.0)[12]對Illumina平臺的原始數據進行過濾，最終得到質量較好的有效數據用于后續分析。

1.3 基因組變異檢測與注釋

使用比對軟件BWA(版本 0.7.12)[13]將2個樣本全基因組重測序數據對比到參考基因組上(https：//bigd.big.ac.cn/search/？dbId=gwh&q=PRJCA002248)，結果使用軟件Picard (版本1.129)(Http：//sourceforge.net/projects/picard/)標記重復讀長(reads)，統計標記后的reads深度及覆蓋度。使用BEDTools(版本 2.25.0)[14]軟件進行覆蓋度統計。使用變異檢測軟件GATK(版本 3.4-46)[15]進行群體SNP檢測，并使用ANNOVAR(版本 2)[16]對檢測出的變異(variant)進行功能注釋。對檢測、過濾得到最終的SNP位點集進行統計分析。使用BreakDancer(版本 1.1.2)[17]進行大片段結構變異的檢測。使用CNVnator(版本0.3.2)[18]來進行拷貝數變異的檢測。

1.4 變異基因分析

在獲得的基因組變異檢測分析基礎上，查找愛文和金煌與參考基因組之間可能存在功能差異的基因序列，進而通過Blastx將基因序列比對到蛋白數據庫NR、SwissProt、KEGG和KOG(E值<0.000 01)，得到給定基因具有最高序列相似性的蛋白，從而得到該基因的蛋白功能注釋信息，進行統計分析。

2 結果與分析

2.1 全基因組重測序數據分析

提取芒果愛文、金煌樣本的DNA經質檢合格后進行全基因組重測序建庫，分別獲得4.73、5.58 Gb原始數據通過測序得到的原始圖像數據轉化序列數據后進行測序評估，芒果愛文、金煌重測序數據經質量控制后分別得到37 207 272、31 501 486個過濾后的數據，正確識別率大于Q30 (測序正確率為99.90%)的堿基分別占92.27%、92.76%，基因組GC含量愛文略低于金煌，分別為36.18%、37.80% (表1)。測序深度為15.01×和13.73×，基因組覆蓋度分別為92.23%、92.27%。測序質量良好，所得數據量和序列長度可準確反映基因組信息。

2.2 SNP檢測與注釋

通過與參考基因組比對(https：//ngdc.cncb.ac.cn/search/？dbId=gwh&q=PRJCA002248)，分別檢測到愛文SNP位點 3 316 161個，金煌SNP位點 3 075 003 個，其中愛文的轉換數是2 288 222個，顛換1 027 939個；金煌的轉換數是2 129 307個，顛換945 696個，2個芒果品種的SNP位點類型多為轉換。金煌的雜合率為57.12%，略高于愛文的雜合率(表2)。

表1 芒果愛文、金煌的重測序數據質量對比

表2 芒果愛文、金煌SNP位點數量統計

對芒果愛文、金煌2個樣本基因組SNP進行注釋統計(表3)，經分析發現，芒果愛文、金煌的SNP變異主要分布在基因間區，分別檢測到1 928 078、1 775 226 個；2個品種在編碼區的非同義突變分別為126 159、123 930個；分別因SNP突變獲得終止子 1 927、1 866個，丟失終止子470、460個；分別在可變剪切位點2 bp以內的基因組區域內檢測到SNP位點 1 608、1 573個。

表3 芒果愛文、金煌基因組中SNP位點注釋

2.3 芒果愛文、金煌基因組 InDel檢測與注釋

結合芒果愛文、金煌2個樣本全基因組重測序數據與參考基因組比對結果，檢測發現2個樣本分別獲得244 686、201 520個InDel位點，其中2個樣本純合基因型的InDel數量分別為154 530、117 125個，雜合率分別為36.85%、41.88%。分別統計芒果愛文、金煌2個樣本基因組InDel位點進行注釋并對基因組位置信息和編碼信息進行統計發現，芒果愛文、金煌2個樣本InDel位點主要分布在基因間區分別有118 467、94 428個。2個樣本InDel位點主要為插入類型分別有4 048、3 861個；2個芒果品種在外顯子區域InDel位點存在顯著差異分別有 8 312、7 605個(表4、表5)。

表4 芒果愛文、金煌基因組中InDel編碼信息統計

表5 芒果愛文、金煌基因組中InDel位點注釋

2.4 芒果愛文、金煌在DNA水平變異基因分析

將芒果愛文、金煌間DNA水平變異基因進行KEGG、KOG、NR、SwissProt數據庫注釋，發現愛文、金煌基因組間存在多種類型的變異基因，與功能數據庫對比得出2個樣本間對應變異基因數量分別為28 981、17 151、30 013、22 910個(圖1)。

為進一步研究含有突變位點的基因功能，向GO數據庫各術語標簽(term)進行映射(圖2)。對芒果愛文、金煌基因組重測序結果進行GO 注釋統計，二者差異基因主要存在于參與生物過程基因中，其中參與水楊酸(SA)介導的信號通路(GO：0009863)的有354個(1.93%)，參與水楊酸代謝過程(GO：0009696)的有269個(1.47%)，富集到對細菌的反應(GO：0009617)相關基因為656個(3.58%)，對真菌的反應(GO：0009620)相關基因有523個(2.85%)，對線蟲的反應(GO：0009624)相關基因9個(0.05%)。

蛋白相鄰類的聚簇(KOG) 變異基因分析結果顯示，芒果愛文、金煌2個樣本間重測序數據在以下幾類功能基因存在變異，其中信號轉導途徑類基因2 529(12%)個，能量生產和轉化類基因778(4%)個，氨基酸轉運與代謝類基因798(4%)個，防御機制類基因181(1%)個(圖3)。

根據芒果愛文、金煌基因組重測序結果進行KEGG注釋(圖4)可知，含有差異位點的基因可注釋到139條代謝途徑，其中與植物抗病性緊密相關的代謝途徑：植物激素信號轉導代謝途徑(ko04075)有差異基因數395 (5.93%)，植物-病原互作代謝途徑(ko04626)差異基因數為303(4.55%)，淀粉和蔗糖代謝途徑(ko00500)基因差異數為200個(3.00%)。

2.5 芒果愛文、金煌抗病相關基因差異性分析

通過對芒果愛文、金煌2個樣本全基因組重測序數據進行樣本間差異基因統計共獲得33 466個差異基因，對其中功能性變異基因進行初步分類和篩選，得到2個樣本間與E3泛素蛋白連接酶相關差異基因數為349個，與WRKY轉錄因子相關差異基因數為100個，與植物抗性相關差異性基因數為383個，其中包括抗病蛋白、耐藥性、抗線蟲蛋白等基因。在所獲得的抗病蛋白差異基因中包括抗病蛋白RGA2差異基因48個、葉銹病10號抗病基因座受體差異基因31個、抗銹激酶基因24個，抗病蛋白RPM1差異基因24個、抗病蛋白DSC1差異基因19個、增強抗病性蛋白基因19個、NBS-LRR抗病蛋白基因18個、含NB-ARC結構域的蛋白質變異基因16個、抗病蛋白At4g27190差異基因15個，抗霜霉病蛋白基因7個、抗稻瘟病蛋白基因9個。含BTB/POZ結構域和錨蛋白重復序列的NPR1蛋白差異基因5個。這些抗病蛋白與植物抗病性密切相關，2個樣本中抗性相關基因的變異可能是導致二者間抗病性差異的原因。

3 討論與結論

植物對病原菌的抗性由多方面原因組成，其中最重要的原因是遺傳因素中的抗病基因，決定抗病基因的產物的結構、功能、信號傳遞途徑的因素是多方面的[19]。Li等公開發表了全球首個栽培芒果染色體級別的參考基因組序列，使開展芒果種質資源的全基因組重測序分析成為可能[20]。芒果品種愛文肉質膩滑，纖維少，味甜，品質較好，但是易感炭疽病；芒果品種金煌樹勢強，果實耐貯藏，中熟，對炭疽病表現為高抗。對芒果炭疽病感病品種愛文和芒果炭疽病抗病品種金煌的全基因組重測序，通過全基因組變異分析技術可以深度挖掘二者存在抗感差異的原因，在基因水平闡明芒果抗炭疽病機制，并且開展芒果優異資源挖掘和基因功能分析。

近年來，基于全基因重測序的SNP檢測分析已經成為植物學科尤其是果樹學遺傳研究的高效便捷的手段[21]。陳璇等對我國野生型大麻和栽培型大麻進行測序深度為10×全基因組重測序分析，共檢測到2 264 150個單核苷酸多態性位點，2個樣本的雜合度分別為0.12%、0.11%，在一定程度上反映我國野生大麻和栽培大麻在基因組水平上的差異性[12]。邵勤等在盤菜品種溫盤2號測序深度為46×全基因組重測序中共檢測到1 424 852個SNP位點，298 244 個Small InDel位點，3 068個結構變異，共導致了42 334個基因變異，并發現了大量參與肉質根形成發育的關鍵基因發生突變[22]。非同義SNP位點變異會使氨基酸序列發生改變進而影響蛋白質序列。洪森榮等對上饒梨2個品種全基因組重測序分析發現，在六月雪、黃皮消2個品種中分別檢測到6 171 357、6 140 603個SNP位點，共有 2 282 個SNP 關聯了2 067 個基因，2個樣本InDel位點分別為800 388、799 603個，共5 115個差異InDel關聯到了3 682個基因，并檢測到耐藥蛋白和抗病蛋白等重要基因發生變異[23]。通過單核苷酸多態性以及結構變異分析可以發掘同一物種間存在表型差異的變異本質。楊赟等在紅葉杜仲、小葉杜仲的10×深度全基因組中發現嚴重影響蛋白質功能的SNP位點有1 516、1 640個，中度影響功能的SNP位點有 41 328、47 192個，篩選了12個特異性的SNP位點，并從SNP位點的差異出發探討了導致杜仲顏色差異的原因[24]。劉冰浩等通過對桂柚一號、沙田柚進行全基因組重測序，分別檢測到了 2 985 308、2 997 229 個SNP差異位點，并且發現結構變異(SV)在兩者變異中起到了很大作用[25]。

本研究通過對芒果易感炭疽病品種愛文、炭疽病高抗品種金煌進行全基因組重測序，以“紅象牙”的基因組為參考，分別檢測到3 316 161、3 075 003個單核苷酸多態性位點，通過對比到數據庫檢測到二者間SNP突變共導致了33 466個基因變異，這些變異基因導致了2個品種間植物-病原互作代謝途徑、激素介導信號通路、抗病蛋白表達等多種生物過程的差異。植物激素不僅可以作為內源信號分子參與植物逆境脅迫應答，也可以作為外援信號分子誘導植物抗病反應[26]。在2個樣本間檢測到大量與E3泛素蛋白連接酶相關基因變異，E3泛素蛋白連接酶不僅可以將泛素與底物蛋白鏈接參與細胞功能調節[27]，還可以作為正負調節因子參與植物對病原菌防御中的效應器誘導免疫反應[28]。蛋白質泛素化修飾廣泛地參與植物免疫應答的調控，尤其廣泛地參與了茉莉酸、水楊酸和乙烯等植物激素防御途徑中的信號傳導和生物合成等環節。泛素介導的蛋白酶水解在水楊酸的信號傳導和生物合成等多個過程中發揮作用，水楊酸是系統獲得性免疫中的關鍵激素，在植物抵御真菌病害侵染中發揮重要作用[29]，在愛文、金煌2個芒果品種間SNP突變而引起的差異基因中，發現有354個參與水楊酸介導的信號通路，有269個參與水楊酸代謝過程，水楊酸以及其類似物可激活抗病相關基因表達并對病原菌的侵害產生明顯抗性[30]。張慶雨等研究發現，水楊酸能夠提高草莓對炭疽病侵染的抗性，水楊酸可能參與了表達水平和炭疽病抗性水平相關的FaNBS20基因的重要響應途徑[31]。水楊酸還可以作為誘導信號激活NPR1表達，NPR1的過量表達轉基因植物株系增強了對細菌和真菌等病原菌的抗性[32]，Silva 等通過草莓過表達NPR1發現其對炭疽病的抗性明顯增加[33]。NPR1包含的BTB結構域增加了其自身參與泛素化的可能性，另外NPR1包含的錨定蛋白重復序列可以介導蛋白之間互作，SA很可能參與調節此過程[34-35]。同時NPR1也是水楊酸途徑中關鍵的輔基活因子，研究表明，NPR1位點突變會導致植物不能激活系統性獲得免疫而易感病[36]。本研究中2個樣本間檢測到的5個NPR1基因變異位點分布在芒果5條不同染色體上，可能是導致2個樣本表現對炭疽病抗性差異的關鍵。

炭疽病在不同植物中的抗性遺傳表現不同，一般是由單基因控制的，大多數受顯性基因控制，少數情況下是由隱性抗病基因或者不完全顯性抗病基因控制，不同抗病基因間相互獨立或存在互作關系[37]。通過對芒果炭疽病抗感品種金煌、高感品種愛文全基因組重測序數據分析，深度挖掘與炭疽病抗性相關的差異基因并解析其所參與的代謝途徑，可為進一步闡明其遺傳特點，進行抗病基因定位、克隆、轉化及特異分子標記開發并且選育芒果抗病品種提供理論研究基礎和良好的數據支撐。