小麥miRNA啟動子的基因組分析

任治鵬，王多佳，田 宇，婁貴成，李 暢，王政委，張 達，蒼 晶

(東北農業(yè)大學生命科學學院，黑龍江哈爾濱 150030)

MicroRNA(miRNA)是一類長度為21～24 nt的非編碼RNA，廣泛存在于植物中，通過負調控其靶基因，參與調控植物的生長發(fā)育和逆境脅迫響應[1-2]。miRNA的生物合成過程主要包括miRNA基因的轉錄、初始轉錄本加工為成熟miRNA以及成熟miRNA裝載形成RNA誘導的沉默復合體(RNA-induced silencing complex，RISC)[3-7]。RISC通過酶切降解靶基因mRNA或者抑制靶基因mRNA的翻譯，從而對靶基因進行轉錄后水平上的調控[8]。

在植物中，能夠轉錄形成miRNA的miRNA基因大部分位于基因間隔區(qū)，作為獨立的轉錄單位，只有部分miRNA基因位于蛋白質編碼基因內，能與宿主基因共同轉錄[9]。研究表明，miRNA基因由RNA聚合酶Ⅱ(RNA polymerase Ⅱ，Pol Ⅱ)轉錄[10]。Pol Ⅱ型啟動子包括核心啟動子區(qū)和上游作用元件，核心啟動子區(qū)主要由TATA-box、轉錄起始位點(transcription start site，TSS)等構成[11]。了解miRNA基因的位置、啟動子的TSS、特定順式作用元件等上游序列特征，對于研究miRNA的表達模式及miRNA介導的調控網絡具有重要意義[12]。近年來，通過生物信息學分析結合高通量測序，對植物miRNA基因的啟動子開展了一定的研究。如在擬南芥中，Megraw等[13]和Xie等[14]通過5′-RACE的方法，發(fā)現大部分擬南芥miRNA啟動子包含TATA-box。Zhou等[15]通過CoVote的方法，在擬南芥、水稻等植物中鑒定了基因間miRNA基因的啟動子，結果表明，miRNA基因與蛋白質編碼基因均由Pol Ⅱ型啟動子啟動，并具有特定的上游元件。Zhao等[16]利用cDNA數據對水稻和擬南芥兩種植物miRNA啟動子元件進行比較，同時通過ChIP方法對擬南芥miRNA基因的TSS進行了預測[17]。隨著植物基因組研究的發(fā)展，促進了miRNA啟動子的鑒定和研究。如Cui等[18]通過基因組數據，定位了水稻miRNA前體(miRNA precursor，pre-miRNA)在染色體上的位置，并通過TSSP軟件預測了miRNA基因的TSS、TATA-box等核心啟動子區(qū)。Liu等[19]和Han等[20]利用大豆基因組數據對miRNA基因的啟動子特征進行了相關分析。Kanjanawattanawong[21]等發(fā)現，橡膠樹中對乙烯響應的miRNA啟動子具有多種植物激素相關作用元件。Zhou等[22]利用TSSP-TCM軟件對擬南芥、毛果楊、水稻、高粱4種植物的miRNA啟動子進行生物信息學分析，發(fā)現基因間和基因內以及保守和非保守miRNA的啟動子具有不同的基因組分布特征及特異性作用元件。此外，研究者在擬南芥[23]、水稻[24]miRNA啟動子中也發(fā)現具有與脅迫相關的特異性轉錄因子結合元件。

六倍體(2n=6x=42，AABBDD)普通小麥(TriticumaestivumL.)是全球種植最廣泛的農作物之一，為人類提供了20%的消耗能量[25]。目前，對小麥miRNA的研究主要集中于克隆鑒定、表達特征分析以及通過預測靶基因進行功能研究等方面[25-26]，然而關于小麥miRNA啟動子的研究報道較少。近年來，國際小麥基因組測序聯盟(International Wheat Genome Sequencing Consortium，IWGSC)對中國春小麥基因組的組裝工作已經完成，其公布的小麥全基因組序列信息對于小麥miRNA啟動子的分析研究具有極大的促進作用。本研究通過生物信息學方法對miRNA基因組位置分布、miRNA啟動子預測以及順式作用元件的富集和特異性進行研究，以期在基因組水平對小麥miRNA啟動子有一個較為全面的了解，為小麥miRNA的轉錄調控探究以及新miRNA的預測提供依據。

1 材料與方法

1.1 小麥miRNA基因位置的預測

所有的小麥miRNA序列來源于miRBase數據庫(Release 22.1，http://www.mirbase.org/)[27]。從Ensembl Plants(ftp://ftp.ensemblgenomes.org/pub/plants/release-48/fasta/triticum_aestivum/dna/)下載中國春小麥基因組序列信息。使用URGI BLAST(https://urgi.versailles.inrae.fr/blast/?dbgroup=wheat_iwgsc_refseq_v2_chromosomes&program=blastn)[28]進行小麥pre-miRNA的基因組定位，選擇identities=100%的blast結果作為miRNA基因的位置，對于identities≠100%的miRNA則將identities≥97%且mismatches≤2的結果作為miRNA基因的位置[19]。預測的miRNA基因通過Mapchart 2.30[29]軟件進行小麥染色體圖譜的繪制。所有能夠定位于小麥基因組上的miRNA基因根據兩種方法進行分類，第一種分類方法是根據miRNA保守性分為保守和非保守miRNA基因，鑒定方法如下：首先利用miRBase提供的所有物種pre-miRNA序列建立本地blast庫，然后將所有小麥pre-miRNA進行本地blast比對。如果其他植物中存在identities>85%且alignment length>90%的相似序列[22]，則該基因為小麥保守miRNA基因，否則為非保守性miRNA基因。第二種分類方法則根據miRNA基因在染色體上的位置進行分類，通過JBrowse(https://urgi.versailles.inra.fr/jbrowseiwgsc/gmod_jbrowse/)[30]判斷miRNA基因的分布情況，將miRNA基因分為基因間和基因內兩種類型。基因間miRNA位于蛋白質編碼基因之間，而基因內miRNA序列位置則與蛋白質編碼基因重疊[15]。判斷miRNA基因染色體位置參考的編碼蛋白質基因數據為IWGSC中國春Annotation v1.1數據庫[31]，包括可高信度(HC)和低信度(LC)蛋白質編碼基因座。

1.2 小麥miRNA基因啟動子的預測

首先通過Zhou等[15]的方法獲得pre-miRNA的基因間5′端上游序列，當pre-miRNA與上游蛋白質編碼基因轉錄方向相同時，如果它們之間的距離大于2 400 bp，則檢索pre-miRNA上游2 000 bp序列；如果距離小于2 400 bp，則檢索上游蛋白質編碼基因下游400 bp與pre-miRNA之間的序列。當pre-miRNA及其上游蛋白質編碼基因轉錄方向相反時，如果它們之間的距離大于4 000 bp，則獲取pre-miRNA上游的2 000 bp序列，如果距離小于4 000 bp，則檢索從pre-miRNA到中間點(上游蛋白質編碼基因與pre-miRNA之間)的序列。將以上方法獲得的序列作為潛在的啟動子預測區(qū)域，利用TSSP(http://www.softberry.com)進行小麥miRNA啟動子及TSS的預測。

1.3 小麥miRNA基因啟動子上游順式作用元件的分析

利用PlantCARE數據庫(http://bioinformatics.psb.ugent.be/webtools/plantcare/html/)[32]對miRNA啟動子TSS到上游2 000 bp序列中的順式作用元件進行分析。對于有多個啟動子的miRNA基因，為獲得盡可能多的順式作用元件信息，選擇距離pre-miRNA起始位點最近的TSS進行分析。為了進一步研究miRNA啟動子區(qū)域基序的特異性，通過MEME(https://meme-suite.org/meme//tools/meme)[33]對miRNA啟動子上游序列中長度為10 bp的基序進行鑒定，選擇結果中前20個基序進行分析，其他設定為默認值。利用全基因組蒙特卡羅模擬方法獲得基序的Z-score，從而判斷各基序在小麥miRNA啟動子的特異性[15]，具體方法如下：首先將所有獲得的miRNA啟動子序列作為目標集，然后在小麥基因組上隨機選擇長度為2 000 bp的序列作為參考集，參考集與目標集的序列數目相同；通過FIMO(https://meme-suite.org/meme//tools/fimo)統(tǒng)計特定基序在目標集和參考集miRNA序列上平均數量，分別記為Nt和Nr。Z-score的計算公式為Z=(Nt/Nr)=σ，它能測量目標集中的基序平均出現次數與參考集樣本的均值之間的歸一化差異[22]。利用CpGPlot(http://emboss.bioinformatics.nl/cgi-bin/emboss/cpgplot)對小麥miRNA TSS上游序列中的CpG島進行分析。

2 結果與分析

2.1 小麥miRNA基因的染色體定位

目前為止，miRBase數據庫(Release 22.1)共收錄122個小麥pre-miRNA序列。小麥pre-miRNA序列和中國春基因組序列blast結果表明，105個(86.1%)pre-miRNA定位于小麥染色體上的150個基因座上，而其余17個(13.9%)pre-miRNA位于未知染色體或基因組的基因座上，下文中不對此類pre-miRNA進行統(tǒng)計。pre-miRNA分布在小麥所有42條染色體上，其中A組染色體上有54個，B染色體上有56個，而D組染色體上有40個，93個(76.2%)pre-miRNA在染色體上只有1個拷貝，含有2個和2個以上拷貝的pre-miRNA分別有4和8個，共占比 9.84%，其他17個pre-miRNA的拷貝為0。

2.2 小麥miRNA基因的啟動子預測結果

150個小麥miRNA基因座中有148個能夠獲得啟動子潛在區(qū)域，對148個miRNA基因座5’上游序列進行啟動子預測，由于部分miRNA基因座能夠預測到多個啟動子，因此共獲得166個miRNA潛在啟動子。115個(77.7%)小麥pre-miRNA基因能夠預測到一個啟動子，其中，69個基因的上游序列只能預測到一個啟動子，而其他基因具有多個啟動子。

TSS是重要的啟動子核心元件，對小麥miRNA基因TSS位點與pre-miRNA距離分布進行統(tǒng)計分析,發(fā)現大部分小麥miRNA基因的TSS分布在上游0.8 kb區(qū)域內以及1.0～1.6 kb區(qū)域內，占全部啟動子TSS數的81.9%(0～0.8 kb：54.2%，1.0～1.6 kb：27.7%)。在所有上游區(qū)域中，小麥miRNA基因的TSS在上游0.2 kb區(qū)域內分布最多(24.1%)，而在上游0.8～1.0 kb區(qū)域分布較少(5.4%)。

根據miRNA在基因組的位置不同，可分為基因間miRNA和基因內miRNA，從圖1A可以看出，兩種miRNA的TSS均在基因上游0.2 kb區(qū)域內分布較多，不同的是基因內miRNA在上游0.2～0.4 kb、0.6～0.8 kb、1.4～1.6 kb間也具有較多的TSS分布，而基因間miRNA在這幾個區(qū)域內無明顯的分布特殊性。根據miRNA的保守性，可分為保守性miRNA和非保守性miRNA，從圖1B可以看出，兩種miRNA的TSS均在基因上游0.2 kb區(qū)域內分布最多，而與非保守miRNA相比，保守miRNA在上游1.4～1.6 kb區(qū)域內也具有較多分布。

A：基因間和基因內miRNA TSS的分布百分比；B：非保守和保守miRNA TSS的分布百分比。

2.3 小麥miRNA基因啟動子上游的特異性順式作用元件

利用PlantCARE對所有miRNA基因TSS上游2 000 bp序列進行順式作用元件分析，結果(圖2)表明，miRNA啟動子區(qū)域中含有的三種順式作用元件較多，分別為CAAT-box、TATA-box和Unnamed_4。此外與ABA響應相關的元件(ABRE)、與MeJA響應相關的元件(TGACG-motif、CGTCA-motif和MYC)、與光響應相關的元件(G-box)以及與多種脅迫和代謝調控相關的元件(MYB)在小麥miRNA基因上游的占比也較高。

圖中數據為啟動子上游順式作用元件所占百分比。

為進一步鑒定小麥miRNA基因啟動子上的基序特異性，通過MEME獲得在TSS上游序列出現頻率較高的且長度為10 bp的基序，然后利用全基因組的蒙特卡羅模擬計算獲得基序的Z-score。Z-score的大小在一定程度上能反應基序在miRNA啟動子上的特異性，Z-score大于2的基序具有miRNA基因啟動子特異性，與miRNA的轉錄調控有關的可能性較高[22]；而Z-score小于2的基序在其他基因組區(qū)域普遍存在，因此不作為miRNA啟動子重要基序進行研究。根據以上標準，獲得了3個Z-score≥2的小麥miRNA啟動子特異性基序(表1)。

表1 小麥miRNA基因啟動子特異性基序

除順式作用元件外，CpG島也是真核生物polⅡ型啟動子的重要特征之一。由于本研究中MIR9670和MIE979可能通過同一個啟動子進行轉錄，因此對114個miRNA基因啟動子的CpG島進行分析，CpGPlot預測結果表明， 61.4%的小麥miRNA基因TSS上游序列有CpG島分布，啟動子區(qū)域含有1、2、3和4個CpG島的miRNA基因分別有41、17、10、2個。

3 討 論

本研究首先將miRBase數據庫目前登錄的所有小麥pre-miRNA序列定位于小麥基因組上，在122個pre-miRNA中有部分序列無法通過blast獲得基因座，其可能原因為：(1)基因位于數據庫中的未知染色體上；(2)pre-miRNA的序列信息不完全，或所研究品種的pre-miRNA序列與參考的中國春基因組序列存在差異；(3)由于小麥基因組較大，組裝困難，目前提供的基因組版本存在部分染色體序列的缺失。本研究染色體定位結果表明，所有小麥染色體上均存在miRNA基因。前人研究表明，在三個染色體組中，B組染色體上的miRNA基因分布最多，根據IWGSC數據庫，編碼蛋白的基因也在B組染色體上的分布最多[21]。本研究選擇blast結果為100%的染色體位置為miRNA基因座，因此多拷貝的miRNA基因序列相同。具有多拷貝的miRNA基因中，只有MIR6197和MIR9774在三個染色體組上均具有拷貝，其他基因只在一個或兩個染色體組上具有拷貝，這說明大部分基因在不同染色體組上存在序列不同的情況。在動物中，miRNA基因通常聚簇并形成多順反子RNA共轉錄，而只有部分植物miRNA基因存在成簇miRNA。與非成簇miRNA不同，成簇的多個miRNA可通過同一個啟動子進行轉錄[18]。Singh等[34]研究表明，小麥中209個miRNA存在于89個多順反子基因座上。本研究中，由于使用的miRBase數據庫的注釋miRNA信息有限，只發(fā)現了1個小麥miRNA簇，該miRNA簇中的MIR9670和MIR9779位于6D染色體上，為非保守miRNA，未在Singh等[35]的研究中報道。對MIR9670和MIR9779的啟動子預測結果發(fā)現，只有一個miRNA啟動子位于上游序列，說明這兩個miRNA可能通過同一個啟動子進行轉錄。小麥成簇miRNA的啟動子特點可通過其他miRNA庫進行更深層次的研究。

150個miRNA基因中有148個具有啟動子潛在區(qū)域，MIR1133和MIR1135基因與上游蛋白質編碼基因距離過近，無法獲得啟動子區(qū)域。本研究通過TSSP對小麥miRNA基因的polⅡ型啟動子進行了預測，結果表明，大部分pre-miRNA(77.7%)上游具有潛在的啟動子序列，而少部分miRNA基因無法預測到啟動子，原因可能為：(1)原始miRNA序列較長，利用pri-RNA加工后形成pre-miRNA序列信息進行啟動子預測，其啟動子可能位于pre-miRNA的上游2 kb以外；(2)大多數啟動子預測軟件都使用同源搜索的方法，因此可能無法預測miRNA啟動子的非保守性啟動子；(3)由于基因組的重復性，部分基因組上具有多個拷貝的pre-miRNA序列為假基因，不發(fā)生轉錄，因此無法進行啟動子預測[19]。miRNA根據基因位置分為基因間miRNA和基因內miRNA，基因內miRNA通常與宿主基因共同轉錄，但Cui等[18]的研究表明，基因內miRNA也可能具有單獨的啟動子，形成獨立的轉錄本。本研究對基因內miRNA和基因間miRNA啟動子數目分別進行了統(tǒng)計，發(fā)現74.3%的基因內miRNA至少具有1個polⅡ型啟動子。以上結果表明，小麥中相當一部分的基因內miRNA也同樣具有獨立的啟動子，而未預測到啟動子的基因內miRNA則可能由宿主基因啟動子啟動轉錄。

TSS是基因的轉錄起始位點，因此TSS的預測對于miRNA基因轉錄特點的研究具有一定的意義。本研究對TSS與pre-miRNA的距離進行統(tǒng)計分析，結果表明，TSS大多數位于pre-miRNA序列上游800 bp內，尤其上游200 bp內。該結果與大豆和水稻中的miRNA基因TSS統(tǒng)計結果類似[18-19]，這說明小麥等植物的大多數miRNA核心啟動子區(qū)域與pre-miRNA序列比較接近，miRNA的PolⅡ啟動子近端區(qū)域的核心啟動子區(qū)可能對miRNA的轉錄起到更大的作用。

植物基因啟動子上的順式作用元件能夠識別結合特定的轉錄因子，從而對基因的轉錄進行相應的時空特異性調控[24]。本研究對TSS上游序列進行了順式作用元件分析，結果表明，miRNA基因啟動子區(qū)域存在多種順式作用元件。其中TATA-box最多，TATA-box是一種廣泛存在的DNA基序，擬南芥和水稻中miRNA特異性基序的功能尚不清楚，但新發(fā)現的基序對于在小麥中鑒定新的特異性miRNA以及進行miRNA的試驗分析具有一定的借鑒意義。除順式作用元件外，CpG島也是重要的啟動子序列特征[35]。本研究在小麥中部分miRNA基因上游序列鑒定出了CpG島，而前人在擬南芥miRNA啟動子中未鑒定到CpG島，水稻miRNA啟動子中鑒定出的CpG島也較少[15]，煙草中的MIR169基因家族中也未鑒定到CpG島的存在[36]。以上結果說明，小麥、水稻等單子葉植物中，CpG島的分布情況可能與雙子葉植物不同。

miRNA在植物的基因調控中起到了重要作用，啟動子是控制miRNA基因表達的重要結構，因此對于miRNA基因啟動子的分析研究具有較大的意義。本研究利用公布的小麥基因組信息對小麥miRNA啟動子進行了分析，研究了小麥miRNA的染色體定位以及TSS分布、順式作用元件特異性等啟動子特征，相關結果對于小麥miRNA表達調控的研究及新miRNA的預測具有一定的借鑒意義。在未來的研究中，隨著小RNA測序以及分子實驗技術的進步，可獲得更多的miRNA功能信息，并結合試驗方法對相關結果進行進一步的驗證。