基于基因組重測序的栗屬種質資源遺傳特性分析

中圖分類號：S664.2 文獻標志碼：A 文章編號：1002-2910（2025）02-0011-08

Analysis of the genetic characteristics of Castanea germplasms based on genome resequencing

WANG Jinping1， TIAN Shoule1， CAI Gongzhan2， SUN Xiaoli1， TAO Diyang3， SHEN Guangning1* （1.Shandong Institute of Pomology，Tai'an，Shandong 2710oo， China; 2.Shandong Shangnong Agricultural Technology Co.，Ltd.，Heze，Shandong，China;3.Culaishan Forest Farm ofTai'an City，Tai'an， Shandong ， China）

Abstract： In this study， genome resequencing technology was employed to conduct an analysis of the genetic characteristics of a total of 180 Castanea germplasms，which included Chinese chestnut （Castanea molissima） cultivars from 16 provinces in China， as wellas European chestnut （Castanea sativa） and Japanese chestnut （Castanea crenata）. As a result， 48 374 108 genetic variation sites were successfully identified. According to the analysis of these variation sites，the l80 cultivars could be classified into eight major groups. There was a phenomenon of gene introgresson among different groups. The population dynamic changes of Chinese chestnut， European chestnut， and Japanese chestnut during the evolutionary process were significantly different. The study also discovered the gene flow of the Chinese chestnut populations from Beijing to Shandong. However， the study has limitations such as uneven distribution of the sample size and limited analysis dimensions. In the future， it is necessary to expand the sample scope and conduct in-depth research by integrating multi-omics approaches.

Key words： chestnut; germplasms; genome re-sequencing; genetic characteristics; population structure; gene flow

板栗（CastaneamollissimaBlume）在中國農林產業與生態建設中占據重要地位[1]。深入剖析其種質資源的遺傳特性，是推動品種優化、培育高產優質且適應性強新品種的關鍵。以往關于板栗的遺傳分析多采用形態學標記，通過觀察植株的形態特征，如樹形、葉形、果實形狀等，以及基于簡單序列重復（SSR）等分子標記技術。這些方法雖能初步揭示板栗種群間的遺傳差異，但由于易受環境因素干擾，且所能檢測到的遺傳變異位點有限，無法全面、精確地展現板栗復雜的遺傳全貌，滿足現代遺傳育種研究的需求[2]。

基因組重測序能夠對已有參考基因組的物種進行深度測序，通過與參考基因組細致比對，能夠高效、精準地識別出單核昔酸多態性（SNP）、短片段插入缺失（InDeI）以及大片段結構變異（SV）等豐富的遺傳變異信息[3]，助力解析群體遺傳結構、挖掘關鍵基因，極大地推動了植物遺傳育種進程[4]。

盡管基因組重測序技術在植物研究中成果豐碩，但栗屬種質資源遺傳特性分析方面的應用研究報道較少。本研究運用基因組重測序技術，對廣泛收集的栗屬種質資源展開全面深入的遺傳分析，旨在明晰栗屬種質資源的遺傳多樣性、明確群體結構劃分以及梳理親緣關系脈絡，為栗屬種質資源的保護利用和遺傳改良提供理論與數據支撐。

1材料與方法

1.1 試驗材料

本研究收集了180份來自不同地區的栗屬種質資源（表1）。這些種質資源在形態特征、生態適應性等方面具有一定的差異，能夠較好地代表栗屬種質資源的多樣性。

1.2 試驗方法

1.2.1DNA提取與質量檢測采用改良的CTAB法提取板栗葉片的基因組DNA，通過 0 . 7 5 % 瓊脂糖凝膠電泳檢測DNA片段大小和DNA降解程度；

NanoDropOne分光光度計檢測DNA純度，檢測OD260/280比值在 1 . 8 ～ 2 . 2 ，無蛋白質和肉眼可見雜質污染；Qubit3.0熒光儀（LifeTechnologies，Carlsbad，CA，USA）檢測DNA濃度，檢測濃度大于 ，總量大于 2 μ g 。

1.2.2基因組重測序文庫構建及庫檢DNA樣品檢測合格后，使用Covaris超聲波破碎儀將DNA隨機打斷成 3 0 0 ～ 5 0 0 b p 的片段，打斷后的樣品磁珠進行片段選擇，使得樣品條帶集中在 2 0 0 ～ 4 0 0 b p 左右。再經末端修復、加A尾、加測序接頭、純化、PCR擴增、PCR產物環化等步驟完成整個文庫制備工作。

文庫構建完成后，先使用Qubit2.0熒光定量儀進行初步定量，隨后使用Agilent2100生物分析儀對文庫的插入片段進行檢測，插入片段大小符合預期后，使用Q-PCR方法對文庫的有效濃度進行準確定量，以保證文庫質量。測序時，采用雙端150bp的測序策略，測序深度達到10X以上，以保證數據的準確性和可靠性。

DNBSEQ上機測序。檢測合格文庫安排上機測序（DNBSEQ）：單鏈環狀DNA分子通過滾環復制，形成一個包含300多個拷貝的DNA納米球（DNB）。將得到的DNBs采用高密度DNA納米芯片技術，加到芯片上的網狀小孔內，通過聯合探針錨定聚合技術（cPAS）進行測序。

1.2.3測序數據篩選測序平臺得到的原始圖像數據文件，經過堿基識別分析即可轉化為原始測序序列，稱為RawData或RawReads，對Rawreads過濾，得到Cleanreads。使用fastp[5]對原始數據進行質控與過濾，具體過濾標準如下：剔除序列中的接頭序列，即adapter；去除reads尾部的polyG和polyX；以滑窗的方式統計窗口內的堿基計算平均質量值，將低質量的滑窗剪裁掉；剔除N個數大于5的reads；剔除質量低于15的堿基占比高于 40 % 的reads；剔除過濾后長度低于 1 5 b p 的reads。

1.2.4變異鑒定利用生物信息學軟件將測序數據與板栗參考基因組[（表2）進行比對。首先，使用BWA[7軟件將測序reads比對到參考基因組上（參數：mem-R，其余參數采用軟件默認參數），然后利用SAMtools[8軟件對比對結果進行排序（參數：sort）、去重（參數：markdup-r）等處理，將比對結果從sam（SequenceAlignment/MAP）文件轉為排序后的bam文件（binaryAlignment/Map），統計比對情況。為了更具體更細致的評估各樣本比對參考基因組上的情況，采用滑窗的方式將基因組劃分成大小相等的子區間，統計各個子區間內的平均測序深度（比對上該區間總的堿基數除以區間大小）以及覆蓋度（區間內有reads覆蓋的區域占區間大小的比例），從而更全面的了解各個樣本的比對情況。

由于二代測序數據reads長度較短且基因組上存在著重復序列從而導致reads的錯誤比對，同時測序數據的分布不均勻。這些原因可能會導致變異鑒定結果的假陽性存在。因此為了降低這一比例，獲得高質量的SNP和INDEL，使用GATK[9]（版本：4.2.5.0；參數：VariantFiltration），根據GATK官方推薦的hard-filtering標準對鑒定到的SNP和Indel分別進行過濾，具體過濾標準如下：

基于GATKhard-filtering過濾后的剩下的變異結果文件，使用vcftools（版本：0.1.16；參數：-maf，-max-missing），剔除掉maf（次要等位基因頻率）小于0.03，基因組型缺失比例大于 20 % 的位點。最后使用剩余的變異位點進行群體結構分析。

1.2.5群體結構分析使用PHYLIP[10]（版本：3.696;參數：neighbor）中的鄰接法（neighbor－joiningmethods，簡稱NJ），構建進化樹。后續基于樹文件（newick格式）使用ggtree進行可視化。

使用fastStructure[11，2]（版本：1.0；參數：-K- seed - c v= 1 0 ）進行群體遺傳結構分析，K值取2～10，基于不同的seed值針對每個K值獨立重復10次（最佳K值根據 fastStructure 提供的chooseK.py程序以及具體的群體背景進行確定）。

使用PSMC[13] （pairwise sequentially Markoviancoalescent）（版本：latest；參數： - N2 5 - t1 5 - r5 -b-p ，堿基突變速率： 每代年數：15），利用個體重測序數據推測該個體所屬的種群在歷史上各個時期的有效群體大小。

利用TreeMix[14]軟件通過從多個種群中獲得等位基因頻率，返回該種群的最大似然（ML）樹，并推斷基因流。

2 結果與分析

2.1栗屬種質資源統計

根據來源對重測序栗屬資源進行了分析，180份栗屬種質資源分別來源于中國16個省（市）以及歐洲和日本國家，其中山東49份、河北13份、安徽11份、湖南15份、湖北16份、江蘇24份、浙江14份、河南9份、陜西8份、北京6份、云南3份、廣東3份、廣西3份、福建2份、江西1份、貴州1份、歐洲栗1份、日本栗1份（表1）。

2.2 測序數據質控

通過對180份栗屬資源的測序，每個樣本得到的Rawreads數量 6 3 0 1 0 8 6 6 ～ 1 8 8 5 7 2 2 0 4 ，通過進行質控與過濾獲得Cleanreads的數量范圍為： ratio在 9 7 . 2 3 % 以上，Q30ratio在 9 0 . 2 0 % 以上，GC含量 34 . 3 9 %～3 6 . 8 0 % 。以2020年發表的板栗全基因組組裝結果為參考基因組，基因組大小為6 8 8 . 9 9 M b 。

通過對比對結果分析，180個樣本的比對率在9 7 . 6 2 % ～ 9 9 . 4 6 % ，對參考基因組的覆蓋度在

，平均測序深度在13.11X～3 2 . 4 3 X 。

2.3等位基因類型及統計

基于堿基類型信息，對48374108個變異位點進行了分類，其中二等位變異、三等位變異、四等位變異、五等位變異、六等位變異、七等位變異數量分別為43564584、4015291、405120、130080、91839、167194。對每個品種具體的變異情況（INDEL變異、堿基轉換、顛換、顛換轉換比、純合基因型、雜合基因型以及雜合比等）進行了統計。整體變異中各類型顛換以及轉換變異的比例進行了可視化（圖2），針對變異在基因組上的分布進行可視化（圖3）。

2.4群體結構分析

群體進化樹構建。根據所測180份資源的重測序結果，基于鄰接法構建了系統進化樹，180份資源可以分成八大類群，大的聚類結果與品種地理來源沒有明顯相關。

如來自于湖南的z-305（沅陵2號）、z-365（沅陵1號），江蘇的z-280（句容短刺）、z-385（青毛軟扎）、z-373（炮車7號）、Z-101（青扎）、z-103（焦扎），湖北的 z - 0 9 3 （鄂栗1號），山東的z-251（輻3）、z-239（青毛軟刺）、z-284（魯35）、z-341（沭河7號），河南的z-353（磚橋處暑紅），安徽的z-042（大紅袍）、z-094（舒城大紅袍），河北的 z - 2 9 0 （西溝2號）等16個品種聚為一類。有部分小規模聚類值得關注，如來自廣東的z-271（韶栗18號）、Z-363（早香1號）聚為一類；來自江蘇的z-019（優選處暑紅）與來自安徽的z-055（處暑紅）聚到一起，來自湖北的z-046（紅毛早）和z-374（淺刺紅毛早）聚到一起；河北的 z - 0 0 7 （東密塢無花）與z-201（貴州野毛栗）聚到一起；z-198（廣西14-5）、z-011（高店10號）、z-324（合肥大紅袍）z-320（白毛栗）與Shenyinou（沈引歐）聚為一類。

群體Structure分析。基于檢測過濾后的品種SNP數據，使用fastStructure進行群體遺傳結構分析，當 K=8 時Cv-Error較小，與系統發育樹分群結果相一致，各來源地表現出基因結構相互滲透的現象（圖4）。因此認為收集的180份樣本分為8個類群比較合適。

2.5 有效群體大小分析

本研究從山東、河北、河南、江蘇、安徽、日本栗、歐洲栗等群體中分別選擇1個樣本的重測序數據推測個體所屬的種群在歷史上各個時期的有效群體大小。如圖5所示，除廣東外，中國各個產區的有效群體變化趨勢相近，而日本栗有效群體變化趨勢則與中國栗差異較大。由圖可知，大約在1 Ma.B.P.（Million-anniversary Before Present）中國栗、日本栗、歐洲栗三個群體大小大致重合在一起。中國栗在 0 . 4 ～ 0 . 5 M a . B . P 有效群體達到最大值，日本栗在0.2～0.3Ma.B.P.有效群體達到最大值，且有效群體數量明顯高于中國栗。此后經歷了長期的衰退。大約在 2 0 ～ 2 5 K a .B.P.（kilion -anniversaryBeforePresent）又開始逐漸擴張。

2.6 基因流分析

根據各群體中每個SNP類型變異位點等位基因數量統計，以每省（市）的品種為一個種群，使用Treemix進行基因流預測。由于廣西、福建、廣東等省份的品種數量較少，不進行分析。對湖北、湖南、山東、陜西、河北、北京、河南、安徽、浙江、江蘇等10個省（市）的品種進行基因流分析。如圖6所示，河南、安徽、浙江、江蘇、湖北、湖南漂移參數較低，群體間遺傳差異小；檢測出從北京板栗種群向山東板栗種群的基因流，反映出山東早期從北京引種交流活動比較頻繁。

3小結與討論

3.1群體結構分析在育種中的利用

本研究首次運用基因組重測序技術對來自中國、歐洲和日本共180份栗屬種質資源展開了全面的遺傳特性分析，成功獲取了高質量的測序數據，并精準鑒定出豐富的遺傳變異位點。

通過系統進化樹和Structure等群體結構分析表明，180份資源可劃分為八大類群，且各地區品種間存在明顯的基因滲透現象，這充分證實了中國板栗資源遺傳背景的多樣性以及不同地區品種間頻繁的基因交流。有效群體大小分析揭示了中國栗、日本栗和歐洲栗在進化歷程中的群體動態變化，為深入理解板栗的進化歷史提供了關鍵線索。基因流分析進一步解釋了不同省份板栗群體間的遺傳關系，為板栗品種的遺傳改良提供了信息。

研究鑒定出的大量遺傳變異位點，為后續開展板栗功能基因挖掘、分子標記開發、進化生物學研究以及深入解析板栗重要農藝性狀的遺傳機制提供了豐富的素材。

板栗在相互引種和選育過程中出現了同名異物或者同物異名的現象，如本研究通過分析發現優選處暑紅與處暑紅，紅毛早與淺刺紅毛早等很可能為同物異名。明確的群體結構和遺傳關系，能夠指導種質資源的合理收集、保存和管理，避免遺傳資源的重復收集和流失。育種工作者可依據遺傳距離和基因流信息，精準選擇具有較大遺傳差異的親本進行雜交，從而提高育種效率，培育出更具優良性狀的新品種，有力推動板栗產業的高質量發展。

3.2 研究展望

盡管本研究取得了一系列重要成果，但仍存在一定的局限性。在試驗材料方面，雖然收集了來自多個地區的栗屬種質資源，但歐洲栗、日本栗以及中國的廣西、云南等部分地區的樣本數量相對較少，可能對群體遺傳結構分析的準確性產生一定影響。未來研究可進一步擴大樣本收集范圍，增加樣本數量，以更全面地揭示栗屬種質資源的遺傳多樣性。

在數據分析方面，本研究主要聚焦于群體結構、有效群體大小和基因流等常規分析，后續研究將結合表型數據，運用全基因組關聯分析（GWAS）等方法，深入挖掘與板栗產量、品質、抗病性等重要性狀緊密相關的遺傳標記和功能基因，為分子標記輔助育種提供更精準的靶點。

此外，本研究僅涉及基因組層面的分析，而植物的生長發育和性狀表現是一個復雜的調控網絡，受到轉錄組、蛋白質組和代謝組等多個層面的協同調控。未來將結合開展多組學聯合分析，多個維度解析板栗的遺傳調控機制，為栗屬種質資源的創新利用和遺傳改良拓寬路徑。

參考文獻：

[1］闞黎娜，李倩，謝爽爽，等．我國板栗種質資源分布及營養成分比較［J]．食品工業科技，2016，37（20）：396-400.

［2］李沛，何治霖，談月霞，等．基于重測序數據與表型性狀的寬皮柑橘遺傳多樣性分析與優異種質篩選［J].中國農業科學，2024，57（23）：4761-4795.

［3］姚遠，鄧利君，胡娟，等．‘脆紅李’及其早熟芽變全基因組重測序分析［J]．園藝學報，2024，51（10）：2255-2266.

[4]VARSHNEYRK，GRANERA，SORRELLSME.Genomicsequencing and discovery of markers for crop improvement[J]．Trends in biotechnology，2009，27（7）：442-450.

[5]CHEN S， ZHOU Y，CHEN Y，et al. Fastp： an ultra- fast all- in-oneFASTQ preprocessor[J].Bioinformatics，2018，34（17）：1884 -1890.

[6]WANG J， TIAN S，SUN X，et al. Construction of Pseudomoleculesfor the Chinese Chestnut （Castanea mollissima） Genome[J].G3-Genes GenomesGenetics，2020，10（10）：3565-3574.

[7]HENG L， RICHARD D. Fast and accurate short read alignment withBurrows-Wheeler transform[J].Bioinformatics，2010，14（25）：1754-1760.

[8]LI H，HANDSAKER B，WYSOKER A，et al.The SequenceAlignment/Map format and SAMtools[J].Bioinformatics，2009，25（16）：2078-2079.

[9]DANECEK P，AUTONA，ABECASIS G，et al.The variant allformatandVCFtools[J].Bioinformatics，201l，27（15）：2156-2158.

[10］RAJA，STEPHENS M，PRITCHARD JK. Fast STRUCTURE：Variational Inference of Population Structurein Large SNPDataSets[J].Genetics，2014，197（2）：573-589.

[11]LAM HM， XU X， LIU X，et al. Addendum： Resequencing of 31wild and cultivated soybean genomes identifies paterns of geneticdiversity and selection[J].Nature Genetics，2011，43（4）：387.

[12］LAYER R M，CHIANGC，QUINLAN AR，et al. LUMPY：aprobabilisticframeworkforstructuralvariantdiscovery[J]．GENOMEBIOLOGY，2014，15（6）：1-19.

[13]SHENGLIN，LIU，MICHAEL，et al. PSMC （pairwise sequentiallyMarkovian coalescent） analysis of RAD （restriction site associatedDNA） sequencing data[J]. Molecular Ecology Resources，2017，（17）：631-641.

[14]FITAK R R. Opt M： estimating the optimal number of migrationedges on population trees using Treemix[J].Biology Methods andProtocols，2021，6（1）：1-10.