宋逸釩，倪 挺，魏 剛

(復(fù)旦大學(xué) 生命科學(xué)學(xué)院，上海 200438)

不孕不育日漸成為現(xiàn)代社會生殖健康面臨的嚴(yán)重問題，其中男性不育約占50%。從生物學(xué)的角度來看，男性生殖問題的關(guān)鍵是如何產(chǎn)生正常的生殖細(xì)胞，而精子發(fā)生又是生殖細(xì)胞產(chǎn)生的核心環(huán)節(jié)。深入研究并認(rèn)知精子發(fā)生的分子機制，對人類的生殖健康、男性的生育調(diào)節(jié)、男性不育的診治均有十分重要的理論和實踐意義。臨床樣品上發(fā)現(xiàn)的與精子發(fā)生障礙相關(guān)的遺傳突變通常需要在小鼠等重要模式動物中進行功能驗證和機制探索，因此對于小鼠精子發(fā)生不同時期細(xì)胞的基因表達調(diào)控研究可為人類精子發(fā)生調(diào)節(jié)機制的闡明提供重要線索。除了信使RNA，新近研究還發(fā)現(xiàn)了不少lncRNA(long non-coding RNA，長非編碼RNA)在生殖細(xì)胞增殖、分化過程中起到重要作用[1]。環(huán)狀RNA(circular RNA,circRNA)作為lncRNA家族的一個新成員，其是否在精子發(fā)生過程中發(fā)揮特定的生物學(xué)功能也引起了廣泛的關(guān)注。

circRNA和過去研究的許多種類RNA不同，其序列并不以基因組外顯子正常的順序排列組成，而是RNA的5’末端和3’末端頭尾相連，形成一個像質(zhì)粒一樣共價閉合的單鏈環(huán)狀分子[2]。按形成circRNA的序列在原始母基因中的類別，circRNA可分為外顯子、內(nèi)含子和基因間區(qū)域來源，其中，目前發(fā)現(xiàn)的circRNA絕大多數(shù)來源于外顯子[2]。研究表明外顯子側(cè)翼的內(nèi)含子中的反向重復(fù)序列(比如人中的Alu序列)可促進兩個線性外顯子的反向成環(huán)而形成circRNA[3]。

近年來的研究也不斷發(fā)現(xiàn)：circRNA具有多種基因調(diào)控功能，如高豐度的circRNA可以競爭所在母基因正常的剪接或可變剪接，從而使得正常剪接的RNA水平下降，引起下游的一系列反應(yīng)[4]；circRNA可作為microRNA(miRNA)海綿吸附特定miRNA，繼而影響該miRNA對其下游靶基因的調(diào)控[5]。此外，circRNA上的m6A(腺嘌呤6位甲基化)修飾可以促進翻譯[6]。circRNA還可以通過結(jié)合RNA結(jié)合蛋白(RBP)[7]后與Pol II互作調(diào)控宿主基因的表達[8]等方式體現(xiàn)其功能。但是在精子發(fā)生過程中，circRNA的生成機制和下游功能尚不清楚。

精子發(fā)生是一個高度有序的過程，它的每一個階段都被精細(xì)調(diào)控。目前已知有一些RNA結(jié)合蛋白如ELAVL1/HuR、hnRNP G-T、NANOS2等在精子發(fā)生過程中起著重要調(diào)控作用[7,9-10]。同時，多種miRNA如miR-21[11]、miR-17-92簇[12]、mir-184[13]等也被證實在精子發(fā)生過程中起重要調(diào)節(jié)作用。有研究發(fā)現(xiàn)circRNA在睪丸中高度富集[14]，與合作者的研究也發(fā)現(xiàn)小鼠精子發(fā)生過程中l(wèi)ncRNA和circRNA均表現(xiàn)出了顯著的表達變化[15]，提示circRNA在精子發(fā)生過程中也可能起著重要作用，其上游的形成機制值得深入研究。利用之前發(fā)表的數(shù)據(jù)，系統(tǒng)分析了小鼠精子發(fā)生五個時期(小鼠精原干細(xì)胞、原始精原細(xì)胞、前細(xì)線期精母細(xì)胞，粗線期精母細(xì)胞及圓形精子細(xì)胞)中的circRNA及其可能成環(huán)機制，并分析探討了隨精子發(fā)生不同時期表達顯著變化的circRNA行使功能的可能機制。

1 材料和方法

1.1 數(shù)據(jù)獲取

所用數(shù)據(jù)來自于本實驗室之前發(fā)表的數(shù)據(jù)[15]，其生物學(xué)樣本均來自于DBA/2和C57BL6兩種品系的老鼠產(chǎn)生的雜交后代，樣本獲取、RNA提取及RNA文庫構(gòu)建流程參見文獻[15]。原始數(shù)據(jù)也可從NCBI SRA數(shù)據(jù)庫下載(SRP067167)。每個時期樣品(處理)均有兩份生物學(xué)重復(fù)，各樣本及它們對應(yīng)的原始reads數(shù)據(jù)見表1。為了讓后續(xù)分析更加精準(zhǔn)，利用NGSQCToolkit[16]去除了原始測序數(shù)據(jù)中的低質(zhì)量reads，并用fastuniq對質(zhì)量過濾后的reads進行去冗余。

表1 參試樣本信息及測序reads數(shù)Table 1 Reference sample information and sequencing reads

1.2 小鼠精子發(fā)生過程中差異表達circRNA分析

CIRI[17]是一款常用的從RNA-seq測序數(shù)據(jù)中預(yù)測鑒定circRNA的軟件，使用CIRI的默認(rèn)參數(shù)，預(yù)測各參試樣本中的circRNA，并將得到的circRNA和CIRCpedia數(shù)據(jù)庫[18]中的circRNA進行比對分析。采用sailfish-cir軟件[19](默認(rèn)參數(shù))對circRNA在各時期的表達進行定量分析。同時，借助sailfish_cir計算出的circRNA read數(shù)目，分析獲得到各個時期的環(huán)形-線性比。最后，為了確保circRNA差異表達分析的可靠性，僅對Ribozero文庫來源的circRNA使用Next maSigPro[20]分析五個時期中差異表達的circRNA，并采用R中的clusterProfiler包[21]對差異表達circRNA的來源基因進行了Gene Ontoloty(GO)功能類分析。

1.3 小鼠精子發(fā)生過程中circRNA上游生成機制分析

為了分析在小鼠精子發(fā)生過程中哪些因素可能影響了外顯子型circRNA的生成，首先分析了circRNA側(cè)翼內(nèi)含子的長度，并將它和隨機抽取的內(nèi)含子的長度進行比較。借助RepeatMasker[22]分析獲得基因組上重復(fù)序列，并歸類到幾種常見的反向重復(fù)序列類型；從它們在基因組中的位置信息判定哪些circRNA的側(cè)翼內(nèi)含子含有反向重復(fù)序列，并計算每種重復(fù)序列在circRNA側(cè)翼內(nèi)含子中的反向互補配對情況。

1.4 小鼠精子發(fā)生過程中circRNA功能機制預(yù)測

使用CIRI-FULL[23]和CIRI-AS[24]從基因組信息中獲取circRNA的序列信息。隨后用miRanda-3.3a預(yù)測circRNA中miRNA的結(jié)合位點(設(shè)置miRanda的參數(shù)為-sc 170 -en 25，相比較miRanda的默認(rèn)參數(shù)-sc 140和-en 1，提升了score和energy閾值，減少了預(yù)測結(jié)果的假陽性)。結(jié)合circRNA序列信息和m6A基序(RRm5ACH)，利用IRESfinder[25]和自己編寫的程序搜索circRNA上的m6A潛在位點預(yù)測circRNA作為翻譯模板的可能性。結(jié)合POSTAR2數(shù)據(jù)庫[26]收集的RBP結(jié)合基序，搜索circRNA中潛在的RBP結(jié)合位點。

2 結(jié)果分析

2.1 小鼠精子發(fā)生5個時期中的環(huán)形RNA

利用CIRI2系統(tǒng)分析了小鼠精子發(fā)生各參試樣本的RNA-seq數(shù)據(jù)，共發(fā)現(xiàn)30 960個circRNA，其中有14 920個(48.2%)circRNA可被CIRCpedia數(shù)據(jù)庫中的circRNA交互驗證(見圖1a)。由于CIRCpedia數(shù)據(jù)庫中只收集了小鼠睪丸的circRNA，并未收集完整小鼠精子發(fā)生各個時期的所有數(shù)據(jù)[18]，因此推測余下的 circRNA(約50%)可能是小鼠精子發(fā)生特有的circRNA。RNaseR處理可以大量去除線性RNA從而富集circRNA，分析結(jié)果顯示，ribozero文庫來源的RNA-seq數(shù)據(jù)鑒定出的circRNA的數(shù)量在5個不同細(xì)胞類型中變化幅度不大；而在精原干細(xì)胞中，RNaseR處理文庫中發(fā)現(xiàn)的circRNA數(shù)要大大多于從ribozero文庫中發(fā)現(xiàn)的circRNA數(shù)目(見圖1b)，這也提示了CIRI發(fā)現(xiàn)circRNA的可靠性。

對精子發(fā)生各個時期的circRNA進行序列來源分析發(fā)現(xiàn)，它們絕大多數(shù)主要來源于外顯子區(qū)域，即外顯型circRNA(見圖1c)。由于circRNA可能通過競爭性剪接影響所在母基因的正常線性mRNA的生成，計算了各時期環(huán)形-線性比，發(fā)現(xiàn)十個樣本的環(huán)形線性比均小于3%(見圖1d)。但有趣的是，圓形精子細(xì)胞中環(huán)形線性比要明顯高于其他精子發(fā)生時期(見圖1d)。因為circRNA不易降解，所以推測圓形精子細(xì)胞中的circRNA增多更可能是circRNA上游生成機制的改變引起的。

2.2 精子發(fā)生過程中差異表達的circRNA及其所在基因的功能富集分析

分析小鼠精子發(fā)生不同時期circRNA本身的豐度差異和所屬母基因的功能情況可為circRNA在精子發(fā)生中的功能推測提供重要線索。判別差異表達的基礎(chǔ)是對circRNA豐度的定量，之前circRNA定量都是直接計算circRNA反向剪切連接處的讀長(Junction reads)數(shù)目，Li等開發(fā)了sailfish-cir軟件[19]，該軟件應(yīng)用改進過的EM(Expectation-Maximization)算法對circRNA進行定量，能校正多種已知的系統(tǒng)性偏差，也克服了之前定量方法存在依賴測序深度、數(shù)據(jù)離散的缺陷，提高了circRNA定量的準(zhǔn)確度。采用sailfish-cir算法對小鼠精子發(fā)生不同時期circRNA的穩(wěn)態(tài)表達量進行了定量，隨后我們使用maSigPro基于線性回歸的方法分析了在不同時期樣本中表達量顯著不同(P≤0.001，軟件默認(rèn)差異判斷標(biāo)準(zhǔn))的circRNA，共得到409個差異表達circRNA。從這409個差異circRNA的熱圖中可以看出，circRNA表達量上升和下降的趨勢表現(xiàn)出了極高的一致性，并且組間差異較小(見圖2a)。為了探索這些差異表達的circRNA與小鼠精子發(fā)生過程的聯(lián)系，對這些circRNA來源的母基因做了基因本體(GO)富集分析，結(jié)果發(fā)現(xiàn)它們富集的生物學(xué)過程包括精子發(fā)生、纖毛運動和形態(tài)、微管束的形成和運動以及染色質(zhì)修飾和組蛋白修飾等(見圖2b)，暗示相應(yīng)circRNA可能在精子發(fā)生中起作用。上述結(jié)果暗示circRNA在精子發(fā)生過程中呈動態(tài)變化且所在基因與精子發(fā)生功能具有相關(guān)性。

圖1 小鼠精子發(fā)生各時期細(xì)胞中circRNA的數(shù)目、類型及環(huán)形線性比
Fig.1 Number,type,and circular-linear ratio of circRNA in mouse spermatogenesis

注：(b)中1和2指同一細(xì)胞類型的兩個生物學(xué)重復(fù)。RNaseR指對SSC細(xì)胞進行RNaseR酶處理去除線性RNA的樣品。(c)中外顯子來源即exon、內(nèi)含子來源即intron、基因間區(qū)域即intergenic。

圖2 小鼠精子發(fā)生不同時期circRNA的差異表達及所在基因功能富集分析Fig.2 Differential expression of circRNA and functional enrichment analysis of the gene in mouse spermatogenesis at different stages

2.3 小鼠精子發(fā)生中circRNA的上游生成機制分析

已有研究表明，在人細(xì)胞中側(cè)翼內(nèi)含子中反向重復(fù)序列的配對可促進外顯子來源的circRNA的生成，外顯子型circRNA的兩個側(cè)翼內(nèi)含子中Alu重復(fù)序列的反向配對(Inverted repeated across，即IRacross)數(shù)目與同一內(nèi)含子內(nèi)部的Alu序列的配對(Inverted repeated within，即IRwithin)數(shù)目之間的競爭是circRNA形成的重要因素(見圖3)，IRacross配對數(shù)目越大，越能促進circRNA的生成[3]。

圖3 circRNA形成及側(cè)翼內(nèi)含子中重復(fù)序列反向互補配對示意圖Fig.3 CircRNA formation and reverse complementary repeated sequences in flanking introns

注：一個circRNA由綠色和藍色兩個exon組成，位于兩個側(cè)翼內(nèi)含子中的重復(fù)序列反向互補配對(IRacross)可以促進circRNA的生成，而位于同一內(nèi)含子中的重復(fù)序列反向互補配對(IRwithin)則會抑制circRNA的生成。

小鼠中除了SINE/Alu序列，還有SINE序列(SINE/B2和SINE/B4)、LINE/L1、ERVL-MaLR和ERVK等內(nèi)源逆轉(zhuǎn)錄病毒來源的重復(fù)序列。那么小鼠精子發(fā)生過程中這些不同類型的重復(fù)序列間的互補配對是否也對circRNA的形成有貢獻呢？為了回答這個問題，對小鼠精子發(fā)生中外顯子來源的circRNA的側(cè)翼內(nèi)含子進行了系統(tǒng)分析。首先將circRNA側(cè)翼內(nèi)含子的長度和隨機抽取的內(nèi)含子進行了比較，結(jié)果發(fā)現(xiàn)，小鼠精子發(fā)生相關(guān)細(xì)胞中所存在的circRNA的側(cè)翼內(nèi)含子的長度明顯比隨機抽取的內(nèi)含子更長(見圖4a)(Wilcoxon Signed Rank Test，P<0.001)。接著評估外顯子來源的circRNA兩個側(cè)翼的內(nèi)含子之間所形成的反向互補配對的數(shù)目，結(jié)果顯示，小鼠精子發(fā)生過程中circRNA側(cè)翼內(nèi)含子中各種類型的IRacross數(shù)目均顯著多于隨機獲取的內(nèi)含子中的IRacross數(shù)目(見圖4b)，提示小鼠精子細(xì)胞中多種重復(fù)序列對circRNA的形成均有潛在貢獻。比較了circRNA的側(cè)翼內(nèi)含子和隨機抽取的內(nèi)含子中IRacross-IRwithin差值，發(fā)現(xiàn)成環(huán)的外顯子中該值也顯著大于對照(見圖4c)，進一步提示側(cè)翼內(nèi)含子中反向重復(fù)序列的配對可能是小鼠精子發(fā)生過程中細(xì)胞內(nèi)circRNA生成的重要促進因素。

2.4 circRNA下游功能預(yù)測

在小鼠精子發(fā)生過程中發(fā)現(xiàn)了許多circRNA，這些circRNA在精子發(fā)生進程中的潛在生物學(xué)功能是值得探討的問題。circRNA可以通過多種方式來發(fā)揮其生物學(xué)功能，如促進所在基因的轉(zhuǎn)錄、競爭所在基因成熟mRNA的產(chǎn)生、作為microRNA(miRNA)的分子海綿、作為多個蛋白質(zhì)結(jié)合的分子海綿、作為翻譯模板等[5-6]。主要從兩個層面來分析這些circRNA潛在的功能，即作為miRNA分子海綿的circRNA和具有翻譯潛能的circRNA，同時對circRNA和RBP的互作進行了初步分析。

2.4.1 circRNA與miRNA結(jié)合作為miRNA“海綿”的功能預(yù)測

精子發(fā)生中miRNA可通過轉(zhuǎn)錄后調(diào)控調(diào)節(jié)其靶基因的表達進而影響精子發(fā)生進程，如miR-19a、miR-19b通過調(diào)節(jié)PTEN的表達來影響原始生殖細(xì)胞的增殖[27-28]， miR-122a通過結(jié)合TNP2對精子生成的后期階段發(fā)揮調(diào)控作用[29]等。那么研究發(fā)現(xiàn)的409個精子發(fā)生中差異表達的circRNA中有多少可能通過miRNA分子海綿的方式起作用呢？利用miRanda-3.3a[30]對這409個差異表達的circRNA進行了miRNA結(jié)合位點預(yù)測，結(jié)果發(fā)現(xiàn)137個circRNA具有miRNA結(jié)合位點，涉及124個miRNA。值得提及的是，研究發(fā)現(xiàn)具有miRNA結(jié)合位點的circRNA中都只有一個miRNA結(jié)合位點(見圖5a)。雖然早期的研究曾報導(dǎo)circRNA可結(jié)合很多的miRNA，但研究表明，只有一個miRNA結(jié)合位點的circRNA也可調(diào)控相應(yīng)miRNA的效應(yīng)濃度并對表型產(chǎn)生影響[31]。發(fā)現(xiàn)MMU_CIRCpedia_39694(chr4:45987462|45990230)這個circRNA來源于精子發(fā)生相關(guān)基因TDRD7，預(yù)測分析發(fā)現(xiàn)它可以結(jié)合mmu-miR-7042-5p，從而有可能減少該miRNA對有功能的TDRD7線性mRNA的影響。精子發(fā)生過程中這些潛在的circRNA-miRNA互作在小鼠精子發(fā)生中的調(diào)控作用有待后續(xù)實驗進一步驗證。

圖4 circRNA側(cè)翼內(nèi)含子長度及其中的反向互補配對的重復(fù)序列分析Fig.4 Length of circRNA flanking introns and repeated sequence analysis of their reverse complemen tary pairs

2.4.2 circRNA上的m6A基序促進翻譯的功能預(yù)測

有研究表明，一些circRNA可作為翻譯模板產(chǎn)生蛋白質(zhì)，這些可被翻譯的circRNA具有m6A修飾和特定的基序(RRm5ACH)，這一特征可招募翻譯起始復(fù)合物并促進核糖體組裝，使circRNA最終可以翻譯出蛋白質(zhì)[6]。為了考察精子發(fā)生過程中是否存在具有潛在翻譯功能的circRNA，以普通的mRNA序列作為對照組，對全體circRNA、差異表達的circRNA和來自精子發(fā)生相關(guān)基因的差異表達circRNA分別進行m6A基序預(yù)測分析，發(fā)現(xiàn)全體circRNA中有85%的序列含有RRm5ACH基序，而作為對照的全部mRNA的序列中該基序的比例僅為69%，兩者間存在顯著差異(見圖5b左，p值 < 0.001，Wilcoxon Signed Rank Test)。而使用IRESfinder[25]發(fā)現(xiàn)所有circRNA中有46%的序列含有IRES(Internal ribosome entry site，內(nèi)部核糖體進入位點)位點，而作為對照的全部mRNA的序列中僅有32%，兩者間亦存在顯著差異(見圖5b右，p值 < 0.001，Wilcoxon Signed Rank Test)。這一結(jié)果暗示精子發(fā)生過程中的部分circRNA有潛在的蛋白編碼功能。對相應(yīng)circRNA進行更為深入的實驗和功能驗證或許有可能為深入理解精子發(fā)生的分子調(diào)控機制開辟新的視角。

2.4.3 circRNA與RNA結(jié)合蛋白(RBP)的互作分析

circRNA可與RBP結(jié)合從而影響特定的生物過程。Du等發(fā)現(xiàn)circ-FOXO3與CDK2蛋白結(jié)合可抑制細(xì)胞周期進程[32]，Abdelmohsen等發(fā)現(xiàn)circ-PABPN1可競爭性結(jié)合HuR并進而抑制HuR與PABPN1的mRNA結(jié)合，從而降低PABPN1翻譯效率[33]。NF90/NF110能和成熟的circRNA直接結(jié)合形成circRNA-蛋白復(fù)合體(circRNP)，并在抗病毒過程中發(fā)揮重要的免疫功能[34]。小鼠精子發(fā)生過程中的circRNA是否也有結(jié)合RBP的潛力？為了回答這個問題，利用公共數(shù)據(jù)庫中的RBP數(shù)據(jù)庫POSTAR2對精子發(fā)生過程中的circRNA進行了系統(tǒng)分析。結(jié)果發(fā)現(xiàn)，共有10 517個circRNA具有RBP結(jié)合位點，其中除了約3 547個circRNA只有1個RBP結(jié)合位點外，其余circRNA都有2個以上的潛在RBP結(jié)合位點(見圖5c)，甚至部分circRNA有10個以上的潛在RBP結(jié)合位點(見圖5c)，如chr19:5800494|5800738含有31個，chr17:39845133|39845215含有30個。另外，除了一個circRNA可以結(jié)合多個RBP，一個RBP也可能被多個circRNA所吸附。比如我們發(fā)現(xiàn)總共有558個circRNA含有MSI2蛋白結(jié)合位點，而MSI2蛋白被報導(dǎo)在精子發(fā)生過程中起重要作用[35]。上述結(jié)果強烈暗示，circRNA可能通過“RBP海綿”的作用來調(diào)節(jié)細(xì)胞中相應(yīng)RBP的效應(yīng)濃度。circRNA-RBP互作在精子發(fā)生中的作用很值得后續(xù)深入研究。

圖5 小鼠精子發(fā)生相關(guān)細(xì)胞中circRNA的m6A基序及miRNA、RBP結(jié)合位點分析Fig.5 Analysis of m6

3 討 論

circRNA是一種特殊形式的內(nèi)源性RNA，其閉環(huán)結(jié)構(gòu)使之可能逃脫細(xì)胞內(nèi)核酸外切酶的作用，并可能通過多種方式發(fā)揮其生物學(xué)功能。circRNA的上游生成機制雖已在人類細(xì)胞中進行了分析和部分功能證明，但小鼠精子發(fā)生過程中大量產(chǎn)生的circRNA是否通過類似的機制生成仍不清楚。通過系統(tǒng)的生物信息學(xué)分析發(fā)現(xiàn)各種類型的反向重復(fù)序列(SINE/L1、SINE/B2、SINE/B4、ERVL-MaLR和ERVK)在小鼠精子發(fā)生過程中產(chǎn)生的circRNA側(cè)翼內(nèi)含子中均有富集(見圖4b、4c)，提示它們也可能具有類似人circRNA側(cè)翼內(nèi)含子中Alu序列促進circRNA生成的功能。本研究結(jié)果拓展了外顯子型circRNA生成機制的順式作用元件類別，為理解circRNA生成的分子調(diào)控機制提供了新的線索，但具體的作用機制還需進一步的實驗驗證(如在circRNA表達載體中引入這些不同類別的重復(fù)序列)。

目前尚未有circRNA在小鼠精子發(fā)生中起作用的報道。雖然我們發(fā)現(xiàn)了精子發(fā)生過程中有大量的circRNA產(chǎn)生，并初步推測了其形成機制，但它們中哪些circRNA對精子發(fā)生有調(diào)控作用，以何種方式起作用仍不清楚。功能驗證需要尋找有潛能的候選circRNA，而本研究預(yù)測分析發(fā)現(xiàn)的部分circRNA的miRNA海綿功能、翻譯模板潛能預(yù)測和RBP分子海綿等可為篩選功能性circRNA提供有價值的參考線索。雖然生信分析不能直接證明這些circRNA在精子發(fā)生中的作用，但其分析結(jié)果暗示部分circRNA有可能通過多種方式來參與精子發(fā)生的調(diào)節(jié)，后續(xù)的研究可在此基礎(chǔ)上篩選候選circRNA，通過設(shè)計反向引物并結(jié)合定量PCR(qRT-PCR)驗證相應(yīng)circRNA的存在，繼而構(gòu)建circRNA過表達載體并注射小鼠睪丸，使得相應(yīng)細(xì)胞中過表達該circRNA；或者在小鼠中通過基因編輯刪除配對的反向重復(fù)序列來下調(diào)circRNA的產(chǎn)生，從而驗證其是否具有調(diào)控精子發(fā)生表型的功能。