基于轉錄組數據的印度谷螟微衛星位點分析

唐培安, 陶冶心, 薛 昊, 袁明龍

(1. 南京財經大學食品科學與工程學院/江蘇省現代糧食流通與安全協同創新中心, 南京 210023; 2. 草地農業生態系統國家重點實驗室/蘭州大學草地農業科技學院, 蘭州 730020)

研究報告
ResearchReports

基于轉錄組數據的印度谷螟微衛星位點分析

唐培安1, 陶冶心1, 薛 昊1, 袁明龍2

(1. 南京財經大學食品科學與工程學院/江蘇省現代糧食流通與安全協同創新中心, 南京 210023; 2. 草地農業生態系統國家重點實驗室/蘭州大學草地農業科技學院, 蘭州 730020)

基于印度谷螟轉錄組測序數據篩選其功能微衛星(EST-SSR)分子標記,設計該蟲的微衛星(SSR)引物并驗證其適用性。用MicroSAtellite微衛星搜索軟件查找印度谷螟轉錄組中微衛星的數量、重復次數以及所有微衛星的位置信息,采用Primer Premier 5軟件設計SSR引物,并進行PCR驗證。結果表明含EST-SSR的序列3 173個,占總搜索序列的8.52%,平均每13 kb中就出現1個EST-SSR。共發現192種微衛星類型,其中重復單元為單堿基、二堿基和三堿基的比例較高,依次是28.21%、39.84%、25.43%。除單堿基重復單元外,出現頻率最高的是AT/AT,有593次。在66對印度谷螟EST-SSR引物中,有28對PCR擴增成功。這說明基于印度谷螟轉錄組數據開發微衛星標記是可行的,本研究獲得的印度谷螟EST-SSR位點為今后開展該蟲的種群遺傳學研究提供了基礎數據。

轉錄組; 高通量測序; 螟蛾科; 微衛星; 分子標記

微衛星(microsatellite)又稱SSR(simple sequence repeat, 簡單重復序列),是指以1～6個核苷酸為重復單元串聯組成的重復序列[1],1981年,微衛星首次被發現,它們大都以2～3個堿基重復類型為主,主要分布于原核生物和真核生物的基因組中[2-3]。微衛星具有高多態性、多等位性、共顯性、高可重復性、數量豐富和對基因組有很好覆蓋性等特征[4-5]。因此,該標記已被廣泛應用于生物遺傳研究中[6]。

印度谷螟Plodiainterpunctella屬鱗翅目Lepidoptera螟蛾科Pyralidae,是世界性分布的重大儲糧害蟲。該蟲在我國南北方廣泛分布[7],依據2004-2005年第6次全國儲糧昆蟲調查,印度谷螟被列為危害嚴重且分布廣泛的主要鱗翅目儲糧害蟲之一[8]。目前,關于印度谷螟SSR的研究還未見報道,在GenBank中也沒有發現登錄的印度谷螟SSR序列,不能滿足其種群遺傳結構研究的需要。因此,開發印度谷螟的微衛星位點,是有效開展其種群遺傳結構和基因交流、系統發育與分子進化等方面研究的重要前提和首要工作,同時對于揭示其危害規律和綜合防治等都具有重要意義。傳統開發SSR標記的方法需要構建文庫、篩選SSR克隆、測序等復雜且工作量大的步驟[9-11]。隨著新一代高通量測序技術的成長,尤其是轉錄組測序的飛速發展,降低了微衛星分子標記的開發成本,提高了開發效率,為非模式物種SSR位點的大批量開發提供了一種經濟、高效的途徑[12-13]。從轉錄組數據中獲得的與功能基因相關SSR位點又叫表達序列標簽SSR(expressed sequence tag SSR,EST-SSR)。目前,從昆蟲轉錄組數據庫中篩選功能EST-SSR標記已有許多成功的報道[14-17]。本研究基于已測印度谷螟轉錄組數據,挖掘其EST-SSR位點,分析這些微衛星位點的類型和分布特點,探討鑒定到的EST-SSR位點的適用性,以期獲得可供該蟲種群遺傳學研究所需的核分子標記,也為印度谷螟近緣物種微衛星位點的發掘提供可借鑒的方法。

1 材料與方法

1.1 材料及轉錄組數據來源

印度谷螟采自江蘇吳江,以大米粉為飼料,在溫度為(30±1)℃、濕度為75%±5%、24 h無光照條件的恒溫培養箱中培養10代以上,選取4齡幼蟲作為試驗材料。采用德國Qiagen公司的RNeasy Plus Universal Mini Kit試劑盒提取RNA,由華大科技進行高通量轉錄組測序(轉錄組數據的NCBI登錄號:SRP060836)。通過Illumina Hiseq2000平臺測序,總計產出5.7 GB數據,共組裝Unigene 37 246個,總長49 682 542 bp,平均長度1 334 bp,N50達到2 307 bp,過濾后質量不低于20的堿基比例(Q20,測序錯誤率<1%)為98.31%,表明測序結果良好。

1.2 微衛星序列篩選

以組裝出來的Unigene作為參考序列,使用Perl操作平臺下的SSR分析軟件MicroSAtellite (MISA) (http:∥pgrc.ipkgatersleben.de/misa/)在Unigene中搜索SSR位點。篩選條件為:單堿基重復次數≥12,二堿基重復次數≥6,三堿基重復次數≥5,四、五、六堿基重復次數≥4,并且兩個相鄰微衛星重復單元之間的最小長度設為100 bp。

1.3 印度谷螟SSR引物設計及PCR擴增

基于SSR重復單元兩翼序列設計引物,采用Primer 5.0對印度谷螟EST-SSR設計了66對引物。引物篩選條件如下:①引物不含有重復序列;②將獲得的引物比對到Unigene序列,引物5′端允許有3個堿基的錯配,而3′端允許有1個堿基的錯配;③去除比對到非目標Unigene的引物,保證引物與Unigene唯一匹配。④使用ssr_finder校驗SSR,使用產物序列來尋找SSR,檢驗結果是否與MISA結果相同,并篩選出相同的SSR產物。設計好的引物由金斯瑞生物技術公司(江蘇南京)合成。

采用TIANamp Genomic DNA Kit(北京天根)試劑盒提取單頭印度谷螟4齡幼蟲的基因組DNA。PCR擴增體系:基因組DNA 1 μL,上游引物2 μL,下游引物2 μL,滅菌水7.5 μL,ExTaqMix 12.5 μL。PCR反應條件:94℃預變性3 min,再進行30個循環,每個循環包括94℃變性30 s,Tm復性30 s, 72 ℃延伸1 min;最后72 ℃延伸10 min。PCR產物用2%的瓊脂糖凝膠電泳檢測。

2 結果與分析

從構建的印度谷螟轉錄組數據庫中,獲得EST-SSR序列3 173個,占總搜索基因的8.52%,平均每13 380 bp中就出現1個SSR位點。

2.1 微衛星數量及分布特征

單堿基、二堿基、三堿基、四堿基、五堿基和六堿基重復的EST-SSR數分別為895、1 264、807、100、28、79,依次占總SSR個數的28.21%、39.84%、25.43%、3.15%、0.88%、2.49%(表1)。可見,在印度谷螟基因數據庫中二堿基重復序列占據優勢,其次是單堿基和三堿基重復序列。此外,復合型SSR有72個,占SSR總數的2.27%。在單堿基、二堿基、三堿基、四堿基、五堿基、六堿基SSR中出現頻率最多的分別是A/T、AT/AT、CCG/CGG、AAAT/ATTT、ACG-CGC/GCGCGT(圖1),在各自重復單元類型中所占比例依次為98.55%、46.91%、36.31%、33.00%、17.86%、12.66%。在五堿基和六堿基中各類重復基元分布較為平均,這是由于此類重復較少的原因。

單堿基重復次數集中在12～15次,二堿基重復次數集中在6～8次,三堿基重復次數集中在5～8次,四堿基重復次數集中于5～6次,五堿基和六堿基則以重復4次為主。且隨著堿基重復單元的重復次數增加,其數量均呈下降趨勢。

表1 印度谷螟轉錄組中不同微衛星類型分布

圖1 印度谷螟轉錄組不同堿基類型微衛星位點分布Fig.1 Distribution of the microsatellites with different repeat motifs in the Plodia interpunctella transcriptome

2.2 微衛星長度分布特征

印度谷螟轉錄組中發現的3 173個EST-SSR長度從12～192個堿基不等,平均長度為16.24個堿基。序列長度12～15個堿基的有2 175個,占總SSR個數的68.5%,其中單堿基、二堿基和三堿基重復單元所占比例分別為35.54%、39.72%、24.74%,分布較為平均。序列長度16～20個堿基有560個,所占比例為17.6%,其中二堿基重復單元比例最高(48.04%),其次是三堿基重復單元(27.50%)。序列長度21～25個堿基的有297個,占總EST-SSR的9.3%,其中三堿基重復單元占29.97%,其次是二堿基重復單元(27.27%),再次是六堿基重復單元(24.58%),其他重復單元數量較少。序列長度大于25的數量較少,僅占總EST-SSR的4.4%,且全部是復合型(圖2)。

圖2 印度谷螟轉錄組中EST-SSR的長度分布Fig.2 Length distribution of EST-SSRs in the Plodia interpunctella transcriptome

2.3 微衛星在轉錄組編碼區分布特征

印度谷螟轉錄組數據中,有23 310個Unigene序列功能注釋成功。在23 310個注釋成功的Unigene序列中發現EST-SSR位點2 226個,總長度13 182 617 bp,其中編碼區858個EST-SSR位點。編碼區EST-SSR中二堿基重復319個,三堿基重復378個,在非編碼區二堿基重復517個,三堿基重復291個。說明在編碼區以三堿基重復為主,非編碼區以二堿基重復為主。

2.4 微衛星引物設計與PCR擴增

根據印度谷螟可用于引物設計的EST-SSR序列共隨機設計了66對引物,部分引物信息見表2。以單頭印度谷螟4齡幼蟲DNA為模板,進行PCR擴增,結果表明66對引物中的31對引物可以擴增得到條帶(圖3)。其中兩對引物擴增產物與預期相差較大,一對引物產生非特異性條帶,其余28對引物與預期相符。

表2 印度谷螟31對EST-SSR引物信息

圖3 印度谷螟31個EST-SSR位點的擴增電泳圖Fig.3 Agarose gel showing the amplification of 31 EST-SSR loci in the Plodia interpunctella transcriptome

3 討論

本研究通過生物信息學方法從印度谷螟轉錄組數據庫37 246條Unigenes中共獲得EST-SSR位點3 173個,其出現頻率(8.52%)高于大多數的昆蟲,如扶桑綿粉蚧Phenacoccussolenopsis(6.33%)[18],煙粉虱Bemisiatabaci(5.07%)[19],橘小實蠅Bactroceradorsalis(4.23%)[11],黃粉蟲Tenebriomolitor(1.67%)[20],云南切梢小蠹Tomicusyunnanensis(1.29%)[21]和蔥地種蠅Deliaantiqua(1.12%)[22]等;而目前僅有少數昆蟲的EST-SSR位點出現頻率高于印度谷螟,如黑翅土白蟻Odontotermesformosanus(9.98%)[23]。導致這種差異的原因,可能與EST-SSR篩選的標準不同以及用于測序的RNA質量相關。當然,最本質的原因應是物種本身的差異所致。

本研究表明,印度谷螟轉錄組中微衛星的種類較為豐富,1～6 核苷酸重復類型都有出現,主要重復類型為1～3核苷酸重復,占SSR總數的93.48%。在單堿基重復基元中, A/T是占優勢的重復基元,這與已研究的大多數昆蟲類似,如褐飛虱Nilaparvatalugens[24]、黏蟲Mythimnaseparata[15]和扶桑綿粉蚧P.solenopsis[18]等。研究發現,基于轉錄組數據庫發掘的昆蟲微衛星以三堿基重復SSR出現的頻率最高,如二點委夜蛾Athetislepigone[25]、灰飛虱Laodelphaxstriatellus[26]、橘小實蠅B.dorsalis[11]、黃粉蟲T.molitor[20]和溫帶臭蟲Cimexlectularius[16]等。普遍的觀點認為,三堿基重復序列占據主導地位,是由于三聯體密碼子的選擇作用所致[27]。因為密碼子的基本單位就是三堿基,其他類型的重復單元會引起閱讀框的改變,增加錯配或者移碼突變的機會。本研究發現,印度谷螟編碼區中三堿基重復也為優勢重復單元,支持該論點。然而,在印度谷螟的所有EST-SSR中最優勢重復單元為二堿基重復,這與麗蠅蛹集金小蜂Nasoniavitripennis、N.giraulti和N.longicornis[28]、黑腹果蠅Drosophilamelanogaster[29]和細梢小卷蛾Rhyacionialeptotubula[17]等類似。以往對昆蟲EST-SSR的研究中,以GC/GC為核心基元的兩堿基重復SSR的存在非常稀少,幾乎為零。但是在印度谷螟二堿基重復單元中GC/GC含量占了1/3左右,僅僅低于重復單元AT/AT,這與同屬于鱗翅目的黏蟲[15]、二點委夜蛾[25]和細梢小卷蛾[17]相似。因此,我們推斷GC/GC重復基元可能與鱗翅目的某些特殊功能有關系。此外,研究表明,SSR的多態性與其長度以及核心基元的重復次數的多少具有相關性,即長度越長、核心基元的重復次數越多,其多態性越高[30]。印度谷螟的EST-SSR長度集中于12～15個堿基占據了70%,而且還有12%(多為二、三堿基重復單元)長度大于20個堿基。因此,印度谷螟可能具有較高的遺傳多態性,這也是該蟲高度適應環境而成為世界性害蟲的遺傳基礎。

眾所周知,鱗翅目昆蟲SSR開發難度較大,由于其SSR的多態性較低、存在微衛星家族現象并有大量的無效等位基因[31]。為得到質量較高且具有潛在多態性的引物,在篩選引物時應盡量選擇側翼長度大于30 bp的位點,且相似度低于70%,核心重復單元重復次數二堿基在10次以上,三堿基應在8次以上,可降低家族現象的影響并且避免非特異性擴增[11]。本文利用轉錄組數據,設計得到的66對EST-SSR引物,28對擴增得到預期大小的條帶,占合成引物的42.42%。盡管這些位點的多態性還有待進一步檢驗,但其數量足夠開展印度谷螟種群遺傳學等相關研究。導致EST-SSR位點擴增失敗有諸多原因:擴增得到條帶遠大于預期,可能由于基因組DNA中含有內含子而在轉錄后未出現,干擾了擴增過程;而其中一對引物產生非特異性條帶可能是由于EST-SSR位點位于同源基因序列上;擴增失敗無條帶,可能是因為EST-SSR在基因組中的表達豐度較低而導致產物濃度很低,未能檢測到。

目前,SSR研究主要應用于遺傳學分析。對于非模式昆蟲的研究其難點在于SSR引物的篩選,因為在昆蟲等節肢動物基因組內SSR的含量較少,且存在家族現象,這加大了研究難度。傳統微衛星篩選方法費時費力,且費用高昂。轉錄組測序技術的發展解決了這一難題,從轉錄組測序結果中可以直接篩選出質量合格的SSR,經過驗證其在相近物種間具有較高的通用性。隨著轉錄組測序技術的進一步發展,測序通量加大,成本降低,這為昆蟲SSR的大規模開發提供了可能,且SSR所在基因的表達水平也可以監測到,為深入研究SSR在基因調控和修飾中所產生的作用提供了技術保障[32]。但基于轉錄組測序獲得的EST-SSR位點,其核心重復次數相較于基因組的SSR少,這是因為在組裝過程中重復次數多的SSR會發生斷裂或者丟失,這影響了EST-SSR位點的多態性。另外,cDNA中缺少內含子或者內含子位于SSR引物之間,這就導致部分EST-SSR位點的擴增失敗。另一方面,基于轉錄組測序數據發掘EST-SSR位點具有快速、方便、低廉以及數據量大等優點,已成為當前發掘昆蟲微衛星位點的重要策略。

綜上所述,本研究基于印度谷螟轉錄組數據篩選獲得了大量EST-SSR位點,并對其數量分布、相關特性及適用性進行了分析。結果表明,基于轉錄組數據發掘印度谷螟微衛星位點是高效可行的,也是發掘非模式生物EST-SSR位點的一種經濟、快速、有效的方法。本研究獲得的印度谷螟EST-SSR位點,為今后開展該蟲的遺傳多樣性及種群遺傳結構研究,闡明其種群遺傳分化及抗藥性基因流動規律,以及制定有效防控策略等提供了基礎數據。

[1] 孫荊濤,楊現明,葛成,等. 微衛星分子標記在昆蟲分子生態學研究上的應用[J]. 南京農業大學學報, 2012, 35(5): 103-112.

[2] Spritz R A.Duplication/deletion polymorphism 5′-to the human β globin gene [J]. Nucleic Acids Research,1981,9(19): 5037-5048.

[3] Selkoe K A, Toonen R J. Microsatellites for ecologists: a practical guide to using and evaluating microsatellite markers[J]. Ecology Letters, 2006,9(5): 615-629.

[4] Sandhu K S, Karaoglu H, Zhang P, et al. Simple sequence repeat markers support the presence of a single genotype ofPucciniapsidiiin Australia[J]. Plant Pathology, 2016,65(7): 1084-1094.

[5] Madesis P, Ganopoulos I, Tsaftaris A.Microsatellites: evolution and contribution [M]∥Walker J M.Methods in molecular biology, USA: Humana Press, 2013: 1-13.

[6] Hadonou A M, Sargent D J, Wilson F, et al. Development of microsatellite markers inFragaria, their use in genetic diversity analysis, and their potential for genetic linkage mapping [J]. Genome, 2004, 47(3): 429-438.

[7] 簡富明. 印度谷螟生物學特性初步研究[J]. 西南農業學報, 1993, 6(3): 80-84.

[8] 嚴曉平,宋永成,沈兆鵬,等.中國儲糧昆蟲2005年最新名錄[J]. 糧食儲藏, 2006,35(2):3-9.

[9] Andrés J A, Bogdanowicz S M. Isolating microsatellite loci: looking back, looking ahead [J]. Molecular Methods for Evolutionary Genetics, 2011, 772: 211-232.

[10]Zane L, Bargelloni L, Patarnello T.Strategies for microsatellite isolation: a review [J]. Molecular Ecology, 2002, 11(1): 1-16.

[11]魏丹丹,石俊霞,張夏瑄,等.基于轉錄組數據的桔小實蠅微衛星位點信息分析[J].應用生態學報,2014,25(6):1799-1805.

[12]Schoebel C N, Brodbeck S, Buehler D, et al. Lessons learned from microsatellite development for nonmodel organisms using 454 pyrosequencing[J]. Journal of Evolutionary Biology, 2013, 26 (3): 600-611.

[13]Fernandez-Silva I, Whitney J, Wainwright B, et al. Microsatellites for next-generation ecologists: a post-sequencing bioinformatics pipeline [J]. PLoS ONE, 2013, 8 (2): e55990.

[14]Comeault A A, Sommers M, Schwander T, et al.Denovocharacterization of theTimemacristinaetranscriptome facilitates marker discovery and inference of genetic divergence [J]. Molecular Ecology Resources, 2012, 12 (3): 549-561.

[15]胡艷華, 李敏, 張虎芳, 等. 粘蟲轉錄組中SSR位點的信息分析[J]. 山西農業大學學報: 自然科學版, 2015, 35 (5): 484-489.

[16]Bai X, Mamidala P, Rajarapu S P, et al. Transcriptomics of the bed bug (Cimexlectularius)[J]. PLoS ONE, 2011, 6 (1): e16336.

[17]Zhu J Y, Li Y H, Yang S, et al.Denovoassembly and characterization of the global transcriptome forRhyacionialeptotubulausing illumina paired-end sequencing[J]. PLoS ONE, 2013, 8 (11): e81096.

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

[19]Xie W, Meng Q S, Wu Q J, et al. Pyrosequencing theBemisiatabacitranscriptome reveals a highly diverse bacterial community and a robust system for insecticide resistance[J]. PLoS ONE, 2012, 7 (4): e35181.

[20]朱家穎, 吳國星, 楊斌. 基于轉錄組數據高通量發掘黃粉甲微衛星引物[J]. 昆蟲學報, 2013, 56 (7): 724-728.

[21]袁遠, 張麗芳, 吳國星,等. 云南切梢小蠹微衛星的高通量發掘[J]. 環境昆蟲學報, 2014, 36(2): 166-170.

[22]Zhang Yujuan, Hao Youjin, Si Fengling, et al. Thedenovotranscriptome and its analysis in the worldwide vegetable pest,Deliaantiqua(Diptera: Anthomyiidae)[J]. G3: Genes| Genomes| Genetics, 2014, 4 (5): 851-859.

[23]Huang Q Y, Sun P D, Zhou X G, et al. Characterization of head transcriptome and analysis of gene expression involved in caste differentiation and aggression inOdontotermesformosanus(Shiraki)[J]. PLoS ONE, 2012, 7 (11): e50383.

[24]劉玉娣, 侯茂林. 褐飛虱EST資源的微衛星信息分析[J]. 昆蟲學報, 2010, 53(3): 239-247.

[25]Li L T, Zhu Y B, Ma J F, et al. An analysis of theAthetislepigonetranscriptome from four developmental stages[J]. PLoS ONE, 2013, 8 (9): e73911.

[26]Zhang F G, Guo H Y, Zheng H J, et al. Massively parallel pyrosequencing-based transcriptome analyses of small brown planthopper (Laodelphaxstriatellus), a vector insect transmitting rice stripe virus (RSV)[J]. BMC Genomics, 2010, 11 (1): 303.

[27]潘海濤, 汪俊君, 王盈盈, 等. 小麥EST-SSR標記的開發和遺傳作圖[J]. 中國農業科學, 2010, 43 (3): 452-461.

[28]Pannebakker B A, Niehuis O, Hedley A, et al. The distribution of microsatellites in theNasoniaparasitoid wasp genome[J]. Insect Molecular Biology, 2010, 19(S1): 91-98.

[29]Schug M D, Wetterstrand K A, Gaudette M S, et al. The distribution and frequency of microsatellite loci inDrosophilamelanogaster[J].Molecular Ecology, 1998, 7(1): 57-70.

[30]Meglécz E, Nève G, Biffin E, et al. Breakdown of phylogenetic signal: a survey of microsatellite densities in 454 shotgun sequences from 154 non model eukaryote species[J]. PLoS ONE, 2012, 7 (7): e40861.

[31]Zhang D X. Lepidopteran microsatellite DNA: redundant but promising [J]. Trends in Ecology & Evolution, 2004, 19 (10): 507-509.

[32]Jhanwar S, Priya P, Garg R, et al. Transcriptome sequencing of wild chickpea as a rich resource for marker development [J]. Plant Biotechnology Journal, 2012, 10 (6): 690-702.

(責任編輯：田 喆)

Analysis of microsatellite loci inPlodiainterpunctellabased ontranscriptome dataset

Tang Peian1, Tao Yexin1, Xue Hao1, Yuan Minglong2

(1.College of Food Science and Engineering/Collaborative Innovation Center for Modern Grain Circulation and Safety,Nanjing University of Finance and Economics, Nanjing 210023, China; 2.State Key Laboratory of Grassland Agro-Ecosystems/College of Pastoral Agricultural Science and Technology, Lanzhou University, Lanzhou 730020, China)

This study is aimed to identify gene microsatellite (EST-SSR) loci from the transcriptome database ofPlodiainterpunctella. Based on these SSR sequences, the selected SSR primer pairs were validated. The number, repetition time and location information of all microsatellites obtained with the microsatellite search tool MicroSAtellite were analyzed.Primer Premier 5 was used to designP.interpunctellaSSR primers, and then these primer pairs were verified by PCR. A total of 3 173 SSRs were identified in 37 246 Unigenes, with one SSR per 13 kb. The majority of microsatellite loci consisted of mono-, di- and tri-nucleotide motifs (28.21%, 39.84% and 25.43%, respectively). However, different types of repeat SSRs had considerably different distribution. Except mono-nucleotide, AT/AT was the most frequent repeat motif (593 repeats) among all SSR motifs. Among the 66 designed primer pairs, 28 pairs were amplified successfully. Our study shows that it is feasible to develop microsatellite markers based onP.interpunctellatranscriptome. The EST-SSRs ofP.interpunctellaobtained in this study are helpful for population genetics studies of this important stored-product pest.

transcriptome; high-throughput sequencing; Pyralidae; microsatellite; molecular marker

2016-06-22

2016-08-02

糧食公益性行業科研專項(201413007-2,201513002-5-3);國家重點研發計劃專項(2016YFD0401004-04);江蘇省高校優勢學科建設工程

Q 968.1

A

10.3969/j.issn.0529-1542.2017.03.007

聯系方式 E-mail:tangpeian@163.com