日本落葉松SNP和InDel位點及其所在抗生物脅迫基因挖掘

摘要：[目的]挖掘日本落葉松SNP和InDel位點及其所在抗生物脅迫基因，為日本落葉松分子育種提供分子標記和候選基因。[方法]使用158份來自中國、日本和英國地區的日本落葉松轉錄組數據，首先進行分子標記位點的鑒定和分類，隨后比較了活動期和休眠期基因的表達水平，最后將帶有可靠分子標記的差異表達基因進行功能注釋。[結果]本研究共鑒定到515 935個SNP位點和1 056個InDel位點，它們分布在35 827個基因上。通過比較不同地區的位點數量，推測日本地區的日本落葉松遺傳多樣性較為豐富。至少在50份轉錄組中出現的非同義突變SNP位點有6 444個，InDel位點為10個，它們分布在3 742個基因上，可以作為可靠的分子標記進行利用。活動期和休眠期的轉錄組比較后，發現帶有可靠分子標記的2 569個基因差異表達；GO注釋后發現其中101個基因與植物對真菌、細菌、卵菌、病毒、昆蟲和線蟲的抗性反應有關。[結論]這些結果不僅為日本落葉松全基因組關聯分析和全基因組選擇育種提供了分子標記，也為利用轉基因和基因編輯手段進行遺傳育種提供了候選基因。

關鍵詞：日本落葉松；分子標記；抗性基因；轉錄組

中圖分類號：S722;S791.223 文獻標識碼：A 文章編號：1001-1498（2024）05-0054-11

單核苷酸多態性標記（single nucleotidepolymorphism，SNP）是由單堿基的轉換及顛倒等突變而產生的不同個體之間堿基變異，是物種可遺傳變異中最常見的一種，目前已被應用于植物抗性育種研究中，是全基因組關聯分析的基礎。Calic等人在山毛櫸（Fagus grandifolia Ehrh.）中檢測出4個與抗病性狀極顯著相關的SNP。插入缺失性標記（insertion and deletion，InDel）是指同一物種不同個體之間基因組同一位點的序列發生不同大小核苷酸片段的插入或缺失，目前在植物抗性育種研究中報道較少。

日本落葉松（Larix kaempferi （Lamb.） Carr.）原產日本，被引種到世界各地作為重要造林樹種。目前，關于日本落葉松林染病的案例越來越多，例如因松材線蟲（Bursaphelenchus xylophilus）引起的松材線蟲病，因疫霉（Phytophthoraramorum）引發的疫霉病。為了適應生存環境，在一個生長周期中日本落葉松需要經歷兩個時期：活動期和休眠期。在活動期，落葉松生長旺盛，能夠為病原物提供豐富的資源和營養，且生長環境也更適宜病原物生長和繁殖；而在休眠期，落4p/YAJ1nS+ij3rSQpcwiEw==葉松停止生長，為病原物提供的資源和營養有限，再加上環境的影響，病原物的活動和繁殖受到限制。在這兩個時期中有很多基因差異表達，構成了日本落葉松適應生存環境的分子基礎，因此“在活動期和休眠期差異表達的部分基因是否在日本落葉松防御病原物攻擊的過程中發揮重要作用”以及“這些基因中存在的序列變異是否可能造成不同基因型擁有不同抗性”值得探討。

日本落葉松基因組序列[9-10]以及大量轉錄組的測定，為鑒定序列變異提供了豐富的序列資源，也為推動分子標記輔助育種技術在日本落葉松抗性育種的應用奠定基礎。本研究利用日本落葉松158份轉錄組數據對SNP和InDel位點進行鑒定和分析，并比較了這些位點所在基因在活動期和休眠期的表達，以期為日本落葉松抗性育種提供分子標記和候選基因。

1 材料和方法

本研究的全流程如圖1所示。

1.1 轉錄組數據收集

78份轉錄組測序材料采自中國林業科學研究院，液氮冷凍后-80℃超低溫冰箱保存，干冰運至北京諾禾致源科技股份有限公司進行轉錄組測序，轉錄組數據已提交到國家生物信息中心（https：//www.cncb.ac.cn）數據庫（登錄號：PRJCA022796）。并且在公共數據庫NCBI（www.ncbi.nlm.nih.gov）上，篩選和下載測序平臺盡量一致、Strategy為RNA-Seq、Layout為PAIRED的日本落葉松轉錄組數據（表1），共同用于鑒定分子標記。參考轉錄本來自日本落葉松的50 690條編碼序列（Codingsequence， CDS）（http：//btg.kazusa.or.jp/blast.html）。

1.2 轉錄組數據處理

利用fastp軟件（v0.32）進行數據過濾，去掉測序接頭，保留長度在20 bp以上的序列；將過濾后的數據按照長度進行篩選，去掉長度小于30bp，或者只有一端的read。將質控后得到的cleanreads，使用BWA軟件比對到參考轉錄本上，獲取reads在參考轉錄本上的定位信息。使用Salmon軟件統計比對到每個轉錄本上的read數量，并計算每個轉錄本的表達量，使用Transcripts perkilobase million（TPM）表示。

1.3 基因組變異檢測和分類

使用變異檢測軟件bcftools，從1 58份轉錄組數據中分別檢測SNP和InDel。隨后對原始結果進行過濾。其中mpileup設定參數：-C50-d 80 000 -L 10 000-q 20；call設定參數-vMm；filter設定參數-i＇%QUAL＞10＆（DP4[2]+DP4[3]）＞2‘。最后使用SnpE仟工具對檢測到的變異進行功能分類。使用TBtools軟件，對來自不同地區（中國黑龍江，中國遼寧，日本和英國）轉錄組中的SNP和InDel進行韋恩分析。

1.4 差異基因分析

158份轉錄組數據包含本課題組之前測定的12份來自活動期和休眠期的日本落葉松轉錄組數據。本研究使用TBtools軟件對這12份轉錄組數據進行了差異基因分析（p＜0.05），隨后利用eggNOG和SwissProt數據庫對得到的差異基因進行注釋。

2 結果

2.1 轉錄組數據統計

158份轉錄組數據統計顯示，每份數據大小在17 292 864～159 873 888 bp之間；質控后得到的數據大小在16 053 498～155 706 686 bp之間，Q30在80.76%～97.3g%之間，高質量數據比例在80.89%～99.74%之間。將質控后得到clean reads Lk對到參考轉錄本上，發現共有48 094個轉錄本得到比對，占總轉錄本數量的94.88%，每份轉錄組中有55.72%～88.00%的clean readsLC對到參考轉錄本上。

2.2 SNP、InDel鑒定及分類

每份轉錄組中SNP的數量在50 520～154 343個之間，InDel的數量在3 045～9 457個之間。158份轉錄組共有SNP位點515 935個（表2），共有InDel位點1 056個（表3），這些位點分布在35 827個轉錄本上，占總轉錄本數量的70.68%。在SNP中，注釋為錯義突變的有302 087個，占比最大，為58.55%，注釋為同義突變的有203 529個，占比39.45%；注釋為終止密碼子獲得的有9 091個，注釋為起始密碼子丟失的有318個；同時注釋為兩種類型的有861個。在InDel中，注釋為移碼突變的最多，有742個，占比70.27%；同時注釋為兩種類型的有34個，且其中一種類型都為終止密碼子獲得。

分別從日本、英國、中國遼寧和中國黑龍江來源的轉錄組中隨機選擇6份進行SNP和InDel位點分析。其中，日本獨有的位點有52 798個，英國獨有的有25 785個，中國黑龍江獨有的有46 1 70個，中國遼寧獨有的有24 056個，這4個地區共有的有12 879個（圖2）。

2.3 位點篩選

CDS內的非同義突變可能會導致基因功能的改變。本研究篩選出注釋為非同義突變的SNP位點有312 399個，InDel位點全部為非同義突變，有1 056個，這些位點所對應的基因共有33 101個。在158份轉錄組中都出現的位點極少，SNP位點僅有755個（0.24%），InDel位點僅有1個（0.09%）。在158份轉錄組中，至少在25份里出現的SNP和InDel分別為9 642個（3.09%）和16個（1.52%）；至少在50份里出現的SNP和InDel分別為6 444個（2.06%）和10個（0.95%）；至少在75份里出現的SNP和InDel分別為5 281個（1.69%）和6個（0.57%）（圖3）。

基于此，本研究選擇至少在50份里出現的SNP和InDel可作為可靠分子標記，它們所對應的基因共有3 742個。

2.4 抗性相關位點與基因分析

本研究利用日本落葉松活動期和休眠期的轉錄組進行差異基因分析，鑒定到1 7 828個差異表達基因（p＜0.05）（圖1）。將這些差異表達基因與帶有可靠分子標記的3 742個基因進行韋恩分析（圖1），發現有2 569個基因既在活動期和休眠期差異表達，又攜帶可靠分子標記（圖1、4）。

對這2 569個基因進行GO注釋：注釋到“對真菌的反應（G0：0009620）”，“對真菌的防御反應（GO：0050832）”，“對共生真菌的反應（GO：0009610）”和“對真菌防御反應的調控（GO：1900150）”的基因有26個，其中在活動期上調的有17個（表4）；注釋到“對細菌的反應（GO：0009617）”，“對細菌的防御反應（GO：0042742）”和“對細菌防御反應的調控（GO：1900424）”的基因有50個，其中在活動期上調的有31個（表5）；注釋到“對卵菌的反應（GO：0002239）”，“對卵菌的防御反應（GO：0002229）”和“對卵菌防御反應的調控（GO：1902288）”的基因有5個（表6）；注釋到“對線蟲的反應（GO：0009624）”，“對線蟲的防御反應（GO：0002215）”，“對昆蟲的反應（GO：0009625）”和“對昆蟲的防御反應（GO：0002213）”的基因有8個（表7）；注釋到“對病毒的反應（GO：0009615）”，“細胞對病毒的反應（GO：0098586）”，“病毒的傳播（GO：0046794）”，“對病毒的防御反應（GO：0051607）”，“病毒誘導的基因沉默（GO：0009616）”和“病毒在多細胞宿主中的轉運（GO：0046739）”的基因有11個（表8）。

3 討論

本研究在CDS區共鑒定到515 935個SNP位點。注釋后發現，39.45%的SNP位點為同義突變，即這些突變不引起氨基酸改變，這與芒果（Mangifera odorata （Huani））和木瓜（Caricapapaya L．）中同義突變SNP的分布頻率一致。同義突變SNP本身不會改變蛋白質序列，因此編碼蛋白的穩定性仍可以維持，這可能構成了日本落葉松、芒果和木瓜等植物適應復雜生存環境的遺傳基礎。

本研究共鑒定到1 056個InDel位點。與SNP相比，InDel發生在CDS區更可能會導致蛋白質結構和功能的改變，從而影響生物體的性狀。所以為了維護蛋白質序列和結構的穩定，InDel發生在CDS區上的數量比SNP少，這與在哺乳動物和果蠅（Drosophila melanogaster）中的研究結果一致。

分子標記位點可以揭示林木種群間及種群內的遺傳多樣性。李培等人利用SRAP標記對來自中國的29個紅椿（Toona ciliata Roem.）種源及1個澳大利亞種源進行遺傳多樣性分析，發現在總的遺傳變異中，種源間分化占79.26%，種源內分化僅占20.74%。使用來自不同地區的日本落葉松的轉錄組數據，本研究發現日本地區的分子標記位點最多（78 219個，41 .42%），說明日本種群的遺傳多樣性較為豐富。日本是日本落葉松的原生地，隨后被引種到不同地區，包括中國、英國。引種后，和原生地存在時空上的隔離，使得不同地區日本落葉松之間的基因交換頻率降低。因此，這四個地區共有的位點相對較少，僅有12 879個（6.82%）（圖2）。

本研究通過比較活動期與休眠期轉錄組，發現在3 742個可靠位點所在基因中有2 569個差異表達，其中有101個與真菌、細菌、卵菌、病毒、昆蟲和線蟲抗性反應相關。這些基因不僅構成了日本落葉松適應生存環境的分子基礎，可能也在日本落葉松防御病原物攻擊的過程中發揮重要作用；另外，這些基因中存在的SNP和InDel可以作為分子標記在日本落葉松抗性育種中進一步研究和利用。

在真菌感染植物細胞時，植物細胞會通過激活特定的防御機制。本研究共鑒定到17個與真菌相關且在活動期表達量高的基因（表4），其中PUX2存在著錯義突變SNP。擬南芥（Arabidopsisthaliana（L．）Heynh．）PUX2的一種突變會降低白粉菌（Golovinomyces orontii）在擬南芥上的繁殖。因此，日本落葉松PUX2的SNP突變也可能影響日本落葉松對真菌的抗性。

本研究共鑒定到31個與細菌相關且在日本落葉松活動期表達量高的基因（表5）。有研究發現，EIN3和EIL1突變后擬南芥增強了對丁香假單胞菌（Pseudomonas syringae）的抵抗力。具有類似情況的基因還有NFXL1，PUB13等。這些結果表明在活動期表達量高的基因在日本落葉松對細菌的抗性反應中可能發揮負調控作用。針對這些負調控抗性的基因，可以通過人工合成miRNA的方式干擾它們的表達以增強日本落葉松的抗性，這種方法在擬南芥、煙草（Nicotianatabacum L．）及番茄（Solanum lycopersicum L．）等植物上的抗病毒研究中已有應用。

本研究鑒定到8個與昆蟲和線蟲相關且在活動期和休眠期差異表達的基因（表7），其中包括和線蟲相關的兩個基因ZIF/和MTPc2。已有研究發現，接種根結線蟲（Meloidogyne spp.）后擬南芥ZIF/的表達量降低，而MTPc2的表達量在接種后第一周表達先升高，隨后降低。這些數據表明，ZIF1和MTPc2在日本落葉松對線蟲的防御中發揮一定的作用。

4 結論

通過分析來自不同地區的158份日本落葉松轉錄組，共獲取了515 935個SNP位點和1 056個InDel位點，它們分布在35 827個基因上。根據不同地區轉錄組中的位點數量，推測日本地區的日本落葉松遺傳多樣性較為豐富。至少在50份樣品中出現的非同義突變SNP位點有6 444個，InDel位點為10個，它們分布在3 742個基因上，可以作為可靠的分子標記進行利用。活動期和休眠期的轉錄組比較后，發現帶有可靠分子標記的2 569基因差異表達；GO注釋后發現其中101個基因與植物對真菌、細菌、卵菌、病毒、昆蟲和線蟲的抗性反應有關。這些結果不僅為日本落葉松全基因組關聯分析和全基因組選擇育種提供了分子標記，也為利用轉基因和基因編輯手段進行遺傳育種提供了候選基因。

（責任編輯：張研）

基金項目：科技創新2030 -重大項目（2022ZD0401602; 2022ZD0401705）