基于RNA-seq數(shù)據(jù)大規(guī)模挖掘蜜蜂球囊菌的SSR分子標記

李汶東, 熊翠玲, 王鴻權(quán), 侯志賢, 童新宇, 張 璐, 付中民, 鄭燕珍, 陳大福, 郭 睿

(福建農(nóng)林大學蜂學學院，福建 福州 350002)

?

基于RNA-seq數(shù)據(jù)大規(guī)模挖掘蜜蜂球囊菌的SSR分子標記

李汶東, 熊翠玲, 王鴻權(quán), 侯志賢, 童新宇, 張 璐, 付中民, 鄭燕珍, 陳大福, 郭 睿

(福建農(nóng)林大學蜂學學院，福建 福州 350002)

基于已獲得的球囊菌轉(zhuǎn)錄組數(shù)據(jù)預測微衛(wèi)星標記(SSR)，并進行SSR位點的信息分析和SSR引物的挖掘.利用軟件MISA對球囊菌轉(zhuǎn)錄組中42610條unigenes進行搜索，共預測出7968個SSRs，分布于5233條unigenes中，其中最主要的重復類型為三核苷酸重復(53.15%)，其次為二核苷酸重復(32.32%)和四核苷酸重復(8.46%).二核苷酸重復中的基序主要是AG/CT(占總量的15.8%).針對所有的SSR位點，利用Primer Premier 5軟件設計出6956對SSR特異性引物，隨機選取20對引物對兩個不同來源的球囊菌樣品進行SSR位點擴增，共有6對引物成功擴增出符合預期的目的片段.研究結(jié)果表明，利用球囊菌轉(zhuǎn)錄組數(shù)據(jù)開發(fā)SSR引物是可行的.

蜜蜂球囊菌; 微衛(wèi)星標記; 轉(zhuǎn)錄組; RNA-seq

微衛(wèi)星(simple sequence repeat, SSR)是以1～6個核苷酸堿基為重復單元(motif)組成的簡單串聯(lián)重復序列.與同類分子標記技術(shù)相比，SSR具有多態(tài)性高、呈共顯性遺傳、重復性好、實驗操作容易等特點[1]，在動植物及微生物等學科領(lǐng)域得到了廣泛應用，但是SSR的特異性一定程度上限制了其引物的通用性[2].對于基因序列未知的物種而言，SSR的開發(fā)較為困難.傳統(tǒng)SSR開發(fā)以基因文庫構(gòu)建法(包括SSR富集文庫)為主[3-4]，不僅試驗過程繁雜，而且效率低下.SSR還可以利用公共基因數(shù)據(jù)庫(NCBI, EMBL, DDBJ)中的共享序列來發(fā)掘，然而對于模式生物或新物種，有限的基因序列信息嚴重制約SSR的開發(fā).2005年以來，第二代高通量測序技術(shù)發(fā)展迅猛，為規(guī)?；z傳變異檢測[5]和標記位點開發(fā)[6-8]帶來了新機遇.如何利用高通量測序數(shù)據(jù)高效、快速地挖掘SSR位點，是當前分子遺傳學領(lǐng)域研究的熱點之一.

蜜蜂球囊菌(Ascosphaeraapis)簡稱球囊菌，是一種特異性侵染蜜蜂幼蟲的致死性真菌病原.白堊病的發(fā)病時期一般為春季和初夏，雖然不對蜂群造成毀滅性打擊，但卻能造成成年蜜蜂數(shù)量及蜂群群勢的大幅下降，從而嚴重影響蜂蜜等產(chǎn)品的產(chǎn)量[9-10]，據(jù)報道，白堊病可造成蜂蜜產(chǎn)量下降5%～37%[11].目前尚無抗真菌藥物被批準用于養(yǎng)蜂生產(chǎn)，一些替代的防治方法，如飼養(yǎng)管理[12-13]等，雖取得一定效果，但總體并不理想.已鑒定的球囊菌SSRs極少，嚴重阻礙球囊菌的分子遺傳學研究.前期研究中，本課題組已成功組裝并注釋了球囊菌的參考轉(zhuǎn)錄組(未發(fā)表數(shù)據(jù))，可為球囊菌的分子研究提供依據(jù)，也為球囊菌的SSR位點挖掘提供寶貴資源.本研究基于已獲得的球囊菌轉(zhuǎn)錄組數(shù)據(jù)，利用MISA軟件進行大規(guī)模的SSR位點預測，并設計特異性SSR引物對不同來源的球囊菌樣本進行SSR位點擴增.

1 材料與方法

1.1 球囊菌轉(zhuǎn)錄組數(shù)據(jù)

球囊菌轉(zhuǎn)錄組數(shù)據(jù)委托廣州基迪奧生物科技有限公司完成，測序平臺為Illumina HiSeq 2500，采用雙端(PE125)測序，RNA-seq數(shù)據(jù)共組裝得到42610條unigenes.前期獲得的球囊菌轉(zhuǎn)錄組數(shù)據(jù)已上傳美國國家生物技術(shù)信息中心(NCBI)SRA數(shù)據(jù)庫(SRA464366).

1.2 供試球囊菌樣品

本研究中用于驗證SSR引物的球囊菌樣品來源于福建福州和浙江蒼山的白堊病幼蟲尸體.從上述蟲尸上刮取部分菌絲和孢子，放入研缽，加入適量無菌水充分研磨，將研磨液轉(zhuǎn)移至干凈的EP管，利用真菌基因組DNA提取試劑盒(TIANamp Genomic，中國)抽提總DNA，作為PCR模板，-20 ℃保存?zhèn)溆?

1.3 SSR鑒定

利用軟件MISA(http://pgrc.ipk-gatersleben.de/misa/)搜索unigenes的微衛(wèi)星標記，按照以下標準從unigenes中查找SSR位點：二核苷酸重復≥6次，三核苷酸重復≥5次，四核苷酸重復≥5次，五核苷酸重復≥5次和六核苷酸重復≥5次.

1.4 SSR引物設計與驗證

根據(jù)MISA的輸出結(jié)果，利用Primer Premier 5(PREMIER Biosofe Int., Palo Alto, CA)對每個含有16 bp堿基重復的SSR設計引物.目標擴增片段設置為必須包含SSR起始-3 bp，終止+6 bp，擴增片段80～300 bp.引物長度設置為18～25 bp，最適長度為22 bp，引物最大允許有1個不能識別的堿基；引物的退火溫度(Tm)設置為55～65 ℃，最適Tm為58 ℃，上下游引物間的Tm差異最大允許3 ℃，引物末端穩(wěn)定性最大為250.引物由上海生工生物工程有限公司合成.為了驗證所設計的引物能否穩(wěn)定擴增，從批量設計出的引物中隨機挑選20對引物進行驗證.PCR反應體系為20 μL：上述DNA模板1 μL，Mixture 10 μL，上下游引物(10 μmol·L-1)各1 μL，無菌水補足20 μL.反應程序：94 ℃預變性5 min；94 ℃變性50 s，55 ℃退火30 s，72 ℃延伸30 s，共33個循環(huán)；72 ℃延伸10 min，4 ℃保存.PCR產(chǎn)物用4%瓊脂糖凝膠電泳檢測.

2 結(jié)果與分析

2.1 球囊菌轉(zhuǎn)錄組中SSR位點的分布特點

利用軟件MISA對球囊菌轉(zhuǎn)錄組42 610條unigenes的數(shù)據(jù)進行搜索，共找到7 968個SSR位點，分布于5 233條unigenes，發(fā)生頻率(含有SSR的unigene數(shù)量與總unigene數(shù)量之比)為12.28%.含1個以上SSR位點的unigene序列有1 615條，以復合型形式存在的SSR序列數(shù)目1 012條，SSR的分布頻率(SSR個數(shù)與總unigene數(shù)量比)為18.70%.這些SSR基序包含1～5 bp的串聯(lián)重復序列.在球囊菌轉(zhuǎn)錄組數(shù)據(jù)中，SSR基序的重復類型占SSR總數(shù)最多的是三核苷酸重復，達53.15%；其次是二核苷酸重復與四核苷酸重復，分別占32.32%和8.46%；五核苷酸和六核苷酸重復的含量相對較少，分別占SSR總數(shù)的3.31%和2.76%(表1).球囊菌轉(zhuǎn)錄組SSR重復單元的重復次數(shù)分布在4～30次之間.二核苷酸重復次數(shù)介于6～30次，三核苷酸重復次數(shù)介于5～30次，四、五和六核苷酸重復次數(shù)最多分布在4次.其中重復次數(shù)最多為5次，占30.40%.

表1 基于重復單元數(shù)目中SSRs在蜜蜂球囊菌轉(zhuǎn)錄組中的出現(xiàn)頻率

根據(jù)不同SSR基序出現(xiàn)頻率分析(圖1)，二核苷酸重復中主要是AG/CT基序，占總量的15.8%，其它基序百分比超過2%的有AC/GT、AT/AT、AAC/GTT、AAG/CTT、AAT/ATT、ACC/GGT、ACG/CGT、AGC/CTG、AGG/CCT、ATC/ATG和CCG/CGG，其比例分別為4.1%、11.9%、5.1%、7.3%、3.3%、4.6%、5.6%、11.6%、6.4%、6.1%和2.3%，三核苷酸中ACT/AGT基序最少，占SSR總數(shù)的0.8%.

圖1 基于基序類型中SSRs在蜜蜂球囊菌轉(zhuǎn)錄組中的出現(xiàn)頻率

2.2 球囊菌SSR引物設計

基于篩選的SSRs，運用Primer Premier 5軟件進行引物的批量設計，按照已設置好的參數(shù)，共對6956個unigenes設計出15963對SSR引物.隨機挑選其中的20對引物對來源于福建福州與浙江蒼山的白堊病幼蟲樣本DNA進行SSR位點擴增(表2)，電泳結(jié)果顯示共有6對引物可成功擴增，且擴增片段符合預期(圖2)，有效擴增率達到30%，說明基于RNA-seq數(shù)據(jù)開發(fā)SSR引物的方法可行，這些SSR位點多態(tài)性的豐富度有待于進一步的試驗評估.

表2 SSR引物的特征

M：DNA marker；泳道1、3、5、7、9、11：來源于福建福州的球囊菌樣品；泳道2、4、6、8、10、12：來源于浙江蒼山的球囊菌樣品；泳道1和2使用的是SSR21引物；泳道3和4使用的是SSR22引物；泳道5和6使用的是SSR35引物；泳道7和8使用的是SSR36引物；泳道9和10使用的是SSR37引物；泳道11和12使用的是SSR38引物.

3 討論

此前，本課題組利用RNA-seq技術(shù)對球囊菌進行深度測序，成功組裝并注釋其參考轉(zhuǎn)錄組.在此基礎(chǔ)上，本研究利用生物信息學方法從球囊菌的轉(zhuǎn)錄組數(shù)據(jù)中預測出多達7968個SSR位點，這些SSR位點種類豐富，以三核苷酸重復類型為主，其次為二核苷酸重復.目前，有關(guān)球囊菌的微衛(wèi)星標記的研究報道幾乎沒有，其它近緣真菌微衛(wèi)星標記的研究報道也較少.Qu et at[14]基于糙皮側(cè)耳(Pleurotusostreatus)基因組預測了2114個SSRs并對其進行了功能富集，發(fā)現(xiàn)外顯子中的三核苷酸重復類型SSR最為豐富，與本研究結(jié)果一致.Catalina et al[15]基于單軸霉(Plasmoparaobducens)的基因組開發(fā)了37個SSR分子標記.Liu et al[16]基于金針菇(Flammulinavelutipes)的基因組預測出1321個SSRs，從中挑取115個SSRs對124個F.velutipes品種進行PCR擴增，有109個可以擴增出目的片段，其中有25個SSR位點具有多態(tài)性.

在一些昆蟲如黃粉蟲(Tenebriomolitor)[17]、扶桑綿粉蚧(Phenacoccussolenopsis)[18]和沙蔥螢葉甲(Galerucadaurica)[19]中，基于轉(zhuǎn)錄組數(shù)據(jù)挖掘SSR位點的研究結(jié)果顯示，單核苷酸重復是最主要的重復類型，其次為三核苷酸重復.在本課題組的另一研究中，基于中華蜜蜂(Apisceranacerana)幼蟲腸道轉(zhuǎn)錄組數(shù)據(jù)開發(fā)SSR分子標記的研究結(jié)果顯示，SSR位點種類以二核苷酸重復類型為主，其次為三核苷酸重復.在表達序列標簽(ESTs)中，多數(shù)情況下單核苷酸重復最多，其次為三核苷酸重復.對于昆蟲，三核苷酸重復是EST-SSR中優(yōu)勢的核心重復類型，因為三核苷酸核心基元在編碼區(qū)較其它重復基元類型更加穩(wěn)定，極少產(chǎn)生編碼框滑動突變現(xiàn)象[20].動物和植物的轉(zhuǎn)錄組或基因組中，二堿基重復的SSR中GC/CG是非常稀有的重復基元.本研究發(fā)現(xiàn)球囊菌轉(zhuǎn)錄組中二核苷酸重復的SSR雖然占比達到32.32%，但未發(fā)現(xiàn)GC/CG重復單元存在，這與動植物轉(zhuǎn)錄組或基因組中的情況相似.

本研究中，為了驗證MISA軟件預測出的球囊菌的SSR引物，隨機挑選20對SSR引物，對國內(nèi)兩個不同來源的球囊菌樣品SSR位點進行PCR擴增，結(jié)果顯示有6對引物可有效擴增出SSR位點，有效擴增率達30%，說明基于RNA-seq數(shù)據(jù)預測出的SSR位點只有部分真實存在，然而，通過此法開發(fā)出的可用于多態(tài)性分析的SSR位點數(shù)目仍非?？捎^.本研究只針對兩個來源的球囊菌樣本進行驗證，樣本數(shù)量較少，下一步將繼續(xù)擴大采樣區(qū)域并采集更多的球囊菌樣本，以獲得不同地域的球囊菌在遺傳分類學上的依據(jù).本研究采用的是瓊脂糖凝膠電泳法，雖然此法較聚丙烯酰胺凝膠電泳(PAGE)法分辨率低，但研究結(jié)果表明部分SSR位點真實存在，這些SSR位點多態(tài)性的確定有待于進一步研究.未來的工作重點是進一步利用已得到的球囊菌轉(zhuǎn)錄組數(shù)據(jù)挖掘更多的SSR位點，并采用分辨率更高的PAGE法和熒光標記的毛細管電泳法對已初步鑒定的SSR位點進行多態(tài)性驗證.

開發(fā)SSR分子標記的傳統(tǒng)試驗方法繁瑣且效率低下，以RNA-seq為代表的二代測序技術(shù)解決了傳統(tǒng)微衛(wèi)星標記開發(fā)的瓶頸問題，使微衛(wèi)星標記的大規(guī)模開發(fā)成為現(xiàn)實.球囊菌是非模式生物，其基因組序列信息早在2006年就已公布[21]，但其基因功能注釋信息直到2016年才公布.本研究說明基于高通量測序數(shù)據(jù)開發(fā)非模式生物的SSR分子標記是一種快速、高效的途徑.目前，以單分子測序技術(shù)為代表的三代測序技術(shù)逐漸興起，有望使SSR分子標記開發(fā)及功能基因定位變得更為方便和快捷.

[1] JARNE P, LAGODA P J L. Microsatellites, from molecules to populations and back[J]. Trends in Ecology and Evolution, 1996,11(10):424-429.

[2] GLENN T C, SCHABLE N A. Isolating microsatellite DNA loci[J]. Methods in Enzymology, 2005,395:202-222.

[3] HAMILTON M B, PINCUS E L, DI-FIORE A, et al. Universal linker and ligation procedures for construction of genomic DNA libraries enriched for microsatellites[J]. Biotechniques, 1999,27(3):500-507.

[4] ZANE L, BARGELLONI L, PATARNELLO T. Strategies for microsatellite isolation: A review[J]. Molecular Ecology, 2002,11(1):1-16.

[5] O′ NEILL E M, SCHWARTZ R, BULLOCK C T, et al. Parallel tagged amplicon sequencing reveals major lineages and phylogenetic structure in the North American tiger salamander (Ambystomatigrinum) species complex[J]. Molecular Ecology, 2013,22(1):111-129.

[6] DAVEY J W, HOHENLOHE P A, ETTER P D, et al. Genome-wide genetic marker discovery and genotyping using next-generation sequencing[J]. Nature Reviews Genetics, 2011,12(7):499-510.

[7] YU J N, WON C, JUN J, et al. Fast and cost-effective mining of microsatellite markers using NGS technology: an example of a Korean water deerHydropotesinermisargyropus[J]. PLoS One, 2011,6(11):e26933.

[8] PANDEY G, MISRA G, KUMARI K, et al. Genome-wide development and use of microsatellite markers for large-scale genotyping applications in foxtail millet [Setariaitalica(L.)][J]. DNA Research, 2013,20(2):197-207.

[9] Bailey L, Ball B V. Honey Bee Pathology[M]. San Diego: Academic Press Inc, 1991:13.

[10] WOOD M. Microbes help bees battle chalkbrood[J]. Agricultural Research Washington D.C, 1998,46(8):16-17.

[11] ZAGHLOUL O A, MOURAD A K, EL K M B, et al. Assessment of losses in honey yield due to the chalkbrood disease, with reference to the determination of its economic injury levels in Egypt[J]. Communications in Agricultural & Applied Biological Sciences, 2005,70(4):703-714.

[12] ANDERSON D L, GIACON H, GIBSON N. Detection and thermal destruction of the chalkbrood fungus (Ascosphaeraapis) in honey[J]. Journal of Apicultural Research, 1997,36(3-4):163-168.

[13] BAGGIO A, GALLINA A, DAINESE N, et al. Gamma radiation: a sanitating treatment of AFB-contaminated beekeeping equipment[J]. Apiacta, 2005,40:22-27.

[14] QU J, HUANG C, ZHANG J. Genome-wide functional analysis of SSR for an edible mushroomPleurotusostreatus[J]. Gene, 2016, 575(2): 524-530.

[15] CATALINA S S, YAZMIN R, DANIEL V, et al. Polymorphic SSR markers forPlasmoparaobducens(Peronosporaceae), the newly emergent downy mildew pathogen of impatiens (Balsaminaceae)[J]. Applications in Plant Sciences, 2015,3(11):1500073.

[16] LIU X B, FENG B, LI J, et al. Genetic diversity and breeding history of Winter Mushroom (Flammulinavelutipes) in China uncovered by genomic SSR markers[J]. Gene, 2016,591:227-235.

[17] 朱家穎,吳國星,楊斌.基于轉(zhuǎn)錄組數(shù)據(jù)高通量發(fā)掘黃粉甲微衛(wèi)星引物[J].昆蟲學報,2013,56(7):724-728.

[18] 羅梅,張鶴,賓淑英,等.基于轉(zhuǎn)錄組數(shù)據(jù)高通量發(fā)掘扶桑綿粉蚧微衛(wèi)星引物[J].昆蟲學報,2014,57(4):395-400.

[19] 張鵬飛,周曉榕,龐保平,等.基于轉(zhuǎn)錄組數(shù)據(jù)高通量發(fā)掘沙蔥螢葉甲微衛(wèi)星引物[J].應用昆蟲學報,2016,53(5):1 058-1 064.

[20] WANG B, EKBLOM R, CASTOE T A, et al. Transcriptome sequencing of black grouse (Tetraotetrix) for immune gene discovery and microsatellite development[J]. Open Biology, 2012,2(4):120 054.

[21] QIN X, EVANS J D, ARONSTEIN K A, et al. Genome sequences of the honey bee pathogensPaenibacilluslarvae, andAscosphaeraapis[J]. Insect Biochemistry and Molecular biology, 2006,15(5):715-718.

(責任編輯:吳顯達)

Large scale development of SSR molecular markers ofAscosphaeraapisbased on RNA-seq data

LI Wendong, XIONG Cuiling, WANG Hongquan, HOU Zhixian, TONG Xinyu, ZHANG Lu,FU Zhongmin, ZHENG Yanzhen, CHEN Dafu, GUO Rui

(College of Bee Science, Fujian Agriculture and Forestry University, Fuzhou, Fujian 350002, China)

To illustrate the molecular genetics and functional genomics ofAscosphaeraapis, simple sequence markers (SSRs) were predicted based on previous transcriptome data ofA.apis, and SSR loci were analyzed and the specific SSR primers were developed. Firstly, 42610 unigenes assembled from RNA-seq data were interpreted by MISA software, and 7968 SSR loci were detected, they were distributed in 5233 unigenes. Among these SSR loci, tri-nucleotide repeats were the dominant repeat type, accounting for 53.15% of the total, which was followed by di-nucleotide repeats (32.32%) and tetra-nucleotide repeats (8.46%). The AG/CT motif was the majority (15.8%) in di-nucleotide repeats. In total, 6956 pairs of primers were designed from all SSRs. Among them, 20 SSR primer pairs were randomly selected and used to amplify the SSR loci inA.apissamples from 2 regions in China. It turned out that 6 pairs of primers were able to amplify the target fragments with expected sizes, indicating that large-scale development ofA.apisSSR primers on the basis of transcriptome data is feasible.

Ascosphaeraapis; microsatellite; transcriptome, RNA-seq

2017-01-03

2017-02-20

現(xiàn)代農(nóng)業(yè)產(chǎn)業(yè)技術(shù)體系建設專項資金(CARS-45-KXJ7)；福建農(nóng)林大學科技發(fā)展資金(KF2015123)；國家自然科學基金(30800806).

李汶東(1996-)，男.研究方向：蜜蜂分子生物學.Email：wendongli22@163.com.通訊作者郭睿(1987-)，男，博士，講師.研究方向：蜜蜂分子生物學.Email：ruiguo@fafu.edu.cn.

S895.1+3

A

1671-5470(2017)04-0434-05

10.13323/j.cnki.j.fafu(nat.sci.).2017.04.013