牙鲆感染遲鈍愛德華氏菌的miRNA轉(zhuǎn)錄組分析

山亞男，鄭津輝，高 虹

（1.天津師范大學(xué)生命科學(xué)學(xué)院，天津 300387；2.天津師范大學(xué)天津市動植物抗性重點實驗室，天津 300387）

牙鲆（Paralichthys olivaceus）屬鰈形目（Pleuronectiformes）牙鲆科（Paralichthyidae）牙鲆屬（Paralichthys），是重要的水產(chǎn)養(yǎng)殖魚類.牙鲆生長過程中易感染細(xì)菌、病毒和真菌等，從而引發(fā)疾病，影響牙鲆的產(chǎn)量和品質(zhì).由于細(xì)菌感染導(dǎo)致的發(fā)病率較高且傳染性強(qiáng)，比病毒和真菌的覆蓋面更廣泛[1-2]，因此細(xì)菌感染致病是牙鲆養(yǎng)殖研究中需要重點解決的問題.遲鈍愛德華氏菌（Edwardsiella tarda）是對水產(chǎn)動物有嚴(yán)重危害的致病菌，主要通過魚類的消化道、鰓或受傷表皮侵入魚體[3]，目前發(fā)現(xiàn)該種細(xì)菌能夠在二十多種魚類中引發(fā)病害，如牙鲆、中華鱉、鯉魚等.遲鈍愛德華氏菌是一種人、魚共患病原菌[4-5]，對人類健康存在潛在威脅.全面深入研究牙鲆對遲鈍愛德華氏菌免疫應(yīng)答的分子機(jī)理，能夠為建立一種有效的可預(yù)防和控制該病菌感染的免疫防治新方法提供理論基礎(chǔ).

MicroRNAs（miRNAs）是一類小的非編碼RNA，通過識別同源序列和干擾轉(zhuǎn)錄、翻譯或表觀遺傳過程來調(diào)控基因表達(dá).目前對魚類miRNA的研究相對較少，大部分miRNA的研究還是集中在疾病、常見動物或者植物方面[6-7].mi RNA在魚類的免疫系統(tǒng)中有重要作用，很多差異表達(dá)miRNA參與入侵感染所產(chǎn)生的炎癥反應(yīng)，表達(dá)量被上調(diào)或下調(diào).如Sha等[8]研究發(fā)現(xiàn)，半滑舌鰨被鰻弧菌侵染后，其肝臟、頭腎、脾和腸中有10個差異顯著表達(dá)的miRNA；Xia等[9]研究發(fā)現(xiàn)，尖吻鱸感染弧菌后，體內(nèi)檢測到了63種miRNA，其中34種表達(dá)量上升，12種表達(dá)量下降；Ordas等[10]使用沙門氏傷寒桿菌對斑馬魚胚胎進(jìn)行攻毒，發(fā)現(xiàn)Toll樣受體信號通路激活并誘導(dǎo)miR-146a/b表達(dá)上調(diào).另外，聚肌苷酸胞苷酸誘導(dǎo)的大黃魚、柱狀黃桿菌侵染的鯉魚、嗜水氣單胞菌侵染的草魚和虹彩病毒屬腫大細(xì)胞病毒侵染的牙鲆中，也存在相似的miRNA差異性表達(dá)結(jié)果[11-13].因此，找到魚類參與免疫應(yīng)答的關(guān)鍵miRNA并鑒定其功能，將有助于深入理解魚類感染病菌的生理病理活動以及深度開發(fā)和應(yīng)用mi RNA.本研究以牙鲆的頭腎組織為實驗材料，用遲鈍愛德華氏菌侵染牙鲆，利用高通量測序技術(shù)分析牙鲆感染遲鈍愛德華氏菌后miRNA的差異表達(dá)情況，為進(jìn)一步評估m(xù)i RNA在牙鲆免疫調(diào)控系統(tǒng)中的作用提供重要的理論依據(jù).

1 材料與方法

1.1 材料及處理

實驗用牙鲆購自天津市西青區(qū)某水產(chǎn)批發(fā)市場，從中選取9尾生長發(fā)育狀態(tài)良好的個體，體長約35 cm，大小均勻.將牙鲆分為感染組和對照組，每組3尾，設(shè)置3個平行對照.選擇消毒滅菌后的可控溫、供氧和水過濾的水循環(huán)魚缸，溫度控制在20℃左右，鹽的質(zhì)量分?jǐn)?shù)為1.7%～1.8%，充氧處理7 d，使牙鲆充分適應(yīng)環(huán)境.牙鲆營底棲生活，不適合強(qiáng)光照射，因此需要對魚缸進(jìn)行遮光處理.實驗開始前將牙鲆放在上述水環(huán)境中暫養(yǎng)1 d使其適應(yīng)實驗環(huán)境，每天喂養(yǎng)1次.

實驗用遲鈍愛德華氏菌由天津市水產(chǎn)養(yǎng)殖病害防治中心提供，保存在本研究實驗室.將其接種于pH值為7.2的LB液體培養(yǎng)基中，放入28℃恒溫振蕩培養(yǎng)箱（200 r/min）中培養(yǎng)24 h.當(dāng)菌液OD600等于1.0時，在4℃、12 000 r/min條件下離心10 min，收集菌體.感染前用PBS離心洗滌菌體3次，用血球計數(shù)板計數(shù).從6尾牙鲆中任意選取3尾作為感染組，用PBS將菌液稀釋到1×107CFU/mL，立即注射進(jìn)感染組牙鲆體內(nèi).分別在感染3 h和6 h時對牙鲆進(jìn)行解剖，取出頭腎組織，將其剪成小塊放入冷凍管中，立即放入液氮中冷凍，之后保存于-80℃冰箱中，用于RNA提取、轉(zhuǎn)錄組測序和生物信息學(xué)分析.委托上海元莘生物醫(yī)藥科技有限公司完成轉(zhuǎn)錄組測序.

1.2 miRNA測序結(jié)果的生物信息學(xué)分析

通過高通量測序平臺（Illumina HiSeq）對牙鲆的miRNA進(jìn)行測序，得到原始的FASTQ數(shù)據(jù)，里面含有帶接頭、低質(zhì)量的序列（raw reads）.為了保證后續(xù)分析的準(zhǔn)確性，運用Fastx-Toolkit軟件對raw reads進(jìn)行質(zhì)量控制，得到高質(zhì)量的測序結(jié)果（clean reads），基于clean reads進(jìn)行生物信息學(xué)分析.為了確定已知的miRNA，將得到的高質(zhì)量測序結(jié)果與miRBase數(shù)據(jù)庫中斑馬魚（Danio rerio）的miRNA前體及成熟體序列進(jìn)行比對.

1.3 差異表達(dá)miRNA的篩選及其靶基因功能分析

對各樣本的已知miRNA進(jìn)行表達(dá)量統(tǒng)計分析，利用TPM對表達(dá)量進(jìn)行均一化處理.應(yīng)用DEGseq2軟件對表達(dá)數(shù)據(jù)進(jìn)行統(tǒng)計學(xué)分析，篩選不同樣本間差異顯著的基因.應(yīng)用miRanda靶基因預(yù)測軟件預(yù)測差異表達(dá)miRNA的靶基因.應(yīng)用Gene Ontology、KEGG pathway軟件對靶基因的代謝通路或信號轉(zhuǎn)導(dǎo)途徑進(jìn)行富集分析.

1.4 差異表達(dá)miRNA的熒光定量PCR驗證

為了驗證本次篩選差異表達(dá)miRNA的準(zhǔn)確性，選擇6個差異表達(dá)的miRNA（miR-730、miR-727-5p、miR-127-3p、miR-205、miR-34a-5p、miR-202-3P），采用熒光定量PCR方法，分析它們在牙鲆注射遲鈍愛德華氏菌后的相對表達(dá)量變化.使用HiScript 1st Strand cDNA Synthesis Kit（cDNA一鏈合成試劑盒）對所提取的RNA進(jìn)行逆轉(zhuǎn)錄處理，使用Primer5設(shè)計差異表達(dá)miRNA的PCR引物，引物序列如表1所示.使用AceQ qPCR SYBR Green Master Mix（without ROX）試劑盒進(jìn)行實時熒光定量PCR反應(yīng)，每個樣本設(shè)置3個平行管，取擴(kuò)增后所得3個平行管Ct值（即每個反應(yīng)管內(nèi)的熒光信號到達(dá)設(shè)定的域值時所經(jīng)歷的反應(yīng)循環(huán)數(shù)）的平均值.

表1 6個差異表達(dá)miRNA的實時定量表達(dá)引物序列Tab.1 Real time quantitative expression primer sequences of six differentially expressed miRNAs

2 結(jié)果與分析

2.1 牙鲆對照組和感染組的miRNA測序數(shù)據(jù)

通過高通量測序平臺得到牙鲆感染組和對照組的轉(zhuǎn)錄組數(shù)據(jù).原始數(shù)據(jù)統(tǒng)計：對照組有11 071 405條raw reads，感染組（3 h）有15 770 243條raw reads，感染組（6 h）有15 122 498條raw reads.對raw reads進(jìn)行質(zhì)量控制，得到clean reads：對照組有9 571 465條clean reads，GC含量為68.59%；感染組（3 h）有13 277 647條clean reads，GC含量為68.08%；感染組（6 h）有12 868 934條clean reads，GC含量為70.44%.不同處理組中miRNA的長度分布如圖1所示.由圖1可以看出，本研究中clean reads的序列長度大部分在21～23 nt，其中22 nt的比例最大，這與miRNA的普遍分布長度一致.

圖1 不同處理組中miRNA長度分布Fig.1 Length distribution of miRNA in different treatment groups

2.2 牙鲆對照組和感染組的miRNA堿基分析

分析對照組和感染組中牙鲆miRNA的堿基組成：對照組中，腺嘌呤（A）、尿嘧啶（U）、胞嘧啶（C）、鳥嘌呤（G）所占比例分別為25.81%、28.42%、19.38%、26.39%；感染組（3 h）中，4種堿基所占比例分別為26.15%、28.04%、19.51%、26.30%；感染組（6 h）中，4種堿基所占比例分別為26.17%、28.12%、19.26%、26.45%.胞嘧啶在3組轉(zhuǎn)錄組中的比例均最低，這與房路京等[14]對大黃魚的研究結(jié)果一致，說明不同物種之間的miRNA堿基種類所占比例存在著相似之處.

2.3 牙鲆差異表達(dá)miRNA分析

根據(jù)篩選標(biāo)準(zhǔn)|log2（foldchange）|＞1和Q值＜0.01，對差異表達(dá)miRNA進(jìn)行評估.從感染組（3 h）中篩選得到12個差異表達(dá)miRNA，其中2個上調(diào)，10個下調(diào)；從感染組（6 h）中篩選得到32個差異表達(dá)miRNA，其中20個上調(diào)，12個下調(diào).感染組（3 h）與感染組（6 h）的差異表達(dá)miRNA對比交集如圖2所示.由圖2可以看出，相較于感染3h，感染6h有更多miRNAs的表達(dá)發(fā)生了變化，說明不同的miRNAs作用時間不同，隨著感染時間的延長，越來越多的miRNAs參與到免疫反應(yīng)中.

圖2 感染組（3 h）和感染組（6 h）的差異表達(dá)miRNA對比交集Fig.2 Comparison of miRNA between the infection group（3 h）and the infection group（6 h）

2.4 差異表達(dá)miRNA的靶基因GO富集分析

對靶基因進(jìn)行GO通路富集分析，結(jié)果顯示，感染組（3 h）中有2 793條差異顯著的富集通路（P＜0.05），其中包括2 021個生物過程（biological process，BP）、315個細(xì)胞組分（cellular component，CC）和457個分子功能（molecular function，MF）；感染組（6 h）有3 048條差異顯著的富集通路（P＜0.05），其中包括2 235個生物過程、342個細(xì)胞組分和471個分子功能.對比分析2個感染組，結(jié)果顯示，差異顯著的富集通路有2 803條（P＜0.05），其中包括1 984個生物過程、341個細(xì)胞組分和478個分子功能.這些差異顯著的富集通路大部分跟生物過程有關(guān)，如生長發(fā)育、信號轉(zhuǎn)導(dǎo)、細(xì)胞通訊等.每個部分選取20個最富集的通路，做GO富集柱狀圖，結(jié)果如圖3所示.由圖3（a）可以看出，感染組（3 h）生物過程中候選靶基因富集程度最高的GO條目是生物調(diào)節(jié)（biological regulation），其次是生物過程調(diào)節(jié)（regulation of biological process）；細(xì)胞組分中富集程度最高的GO條目是質(zhì)膜（plasma membrane），其次是細(xì)胞外圍（cell periphery）；分子功能中富集程度最高的GO條目是結(jié)合（binding），其次是細(xì)胞骨架蛋白結(jié)合（cytoskeletal protein binding）.由圖3（b）可以看出，感染組（6 h）生物過程中候選靶基因富集程度最高的GO條目是生物調(diào)節(jié)，其次是多細(xì)胞生物過程（multicellular organismal process）；細(xì)胞組分中富集程度最高的GO條目是細(xì)胞外圍，其次是質(zhì)膜；分子功能中富集程度最高的GO條目是結(jié)合，其次是蛋白質(zhì)的結(jié)合（protein binding）.由圖3（c）可以看出，感染組（3 h vs 6 h）生物過程中候選靶基因富集程度最高的GO條目是發(fā)育過程（developmental process），其次是身體結(jié)構(gòu)發(fā)育（anatomical structure development）；細(xì)胞組分中富集程度最高的GO條目是細(xì)胞外圍，其次是質(zhì)膜；分子功能中富集程度最高的GO條目是結(jié)合，其次是蛋白質(zhì)的結(jié)合.

圖3 GO富集柱狀圖Fig.3 Histogram of GO enrichment

2.5 靶基因KEGG富集分析

對牙鲆感染組差異顯著表達(dá)的miRNAs進(jìn)行KEGG富集分析.感染組（3 h）差異顯著表達(dá)的miRNAs共獲得71個KEGG富集通路（P＜0.05），其中包括31個有機(jī)系統(tǒng)、17個環(huán)境信息處理、6個細(xì)胞過程、13個人類疾病和4個代謝通路.選取差異顯著的30條通路繪制散點圖，結(jié)果如圖4（a）所示，富集程度最大的KEGG通路是軸突導(dǎo)向（axon guidance），其次是谷氨酸突觸（glutamatergic synapse），富集靶基因數(shù)目最多的是絲裂原活化蛋白激酶信號通路（MAPKsignaling pathway）.感染組（6 h）差異顯著表達(dá)的miRNAs共獲得69個KEGG富集通路（P＜0.05），其中包括33個有機(jī)系統(tǒng)、18個環(huán)境信息處理、6個細(xì)胞過程、10個人類疾病和2個代謝通路.選取差異顯著的30條通路繪制散點圖，結(jié)果如圖4（b）所示，富集程度最大的KEGG通路是軸突導(dǎo)向，其次是Rap1信號通路（Rap1 signaling pathway），富集靶基因數(shù)目最多的是癌癥通路（pathways in cancer）.將感染組（3 h）和感染組（6 h）兩組進(jìn)行對比分析，共得到79個KEGG富集通路（P＜0.05），其中包括33個有機(jī)系統(tǒng)、19個環(huán)境信息處理、6個細(xì)胞過程、20個人類疾病和1個代謝通路.選取差異顯著的30條繪制散點圖，結(jié)果如圖4（c）所示，富集程度最大的通路是軸突導(dǎo)向，其次是鈣信號通路（calcium signaling pathway），富集靶基因數(shù)目最多的是絲裂原活化蛋白激酶信號通路.這些信號通路共同參與牙鲆細(xì)胞的免疫調(diào)節(jié)以及分裂分化等生物學(xué)過程.

圖4 KEGG富集散點圖Fig.4 Catter plot of KEGG enrichment

2.6 差異表達(dá)miRNA的實時熒光定量PCR驗證

為了驗證高通量測序結(jié)果的可靠性，選取靶基因，對與免疫相關(guān)的6個差異表達(dá)miRNA（miR-730、miR-727-5p、miR-127-3p、miR-205、miR-34a-5p、miR-202-3p，其中3個表達(dá)上調(diào)，3個下調(diào)）進(jìn)行實時熒光定量PCR驗證，結(jié)果如表2所示.miR-730預(yù)測的靶基因為甘露糖受體（mannose receptor，MR），包括牙鲆陽離子非依賴性甘露糖-6-磷酸受體樣、牙鲆巨噬細(xì)胞甘露糖受體1樣（X1、2、3）、牙鲆甘露糖受體C型1和牙鲆胰島素樣生長因子2受體；miR-727-5p預(yù)測到的靶基因是牙鲆C型甘露糖受體2樣；miR-127-3p預(yù)測到的靶基因是牙鲆C型甘露糖受體2樣；miR-205預(yù)測的靶基因之一是TIRAP，這是Toll樣受體（Toll-like receptor，TLR）信號通路中的接頭分子，上接TLR，下接MyD88；miR-34a-5p預(yù)測的靶基因是TLR13；miR-202-3p預(yù)測的靶基因之一是TLR7.

表2 實時定量PCR和轉(zhuǎn)錄組分析結(jié)果Tab.2 Results of qPCR and RNA-seq

由表2可知，6個miRNA的實時熒光定量PCR相對表達(dá)量結(jié)果和高通量測序結(jié)果趨勢一致，miR-730、miR-727-5p、miR-127-3p表達(dá)下調(diào)，miR-205、miR-34a-5p、miR-202-3p表達(dá)上調(diào)，由此表明轉(zhuǎn)錄組測序結(jié)果可信.

3 討論與結(jié)論

本研究用遲鈍愛德華氏菌感染牙鲆，利用高通量測序平臺對牙鲆的頭腎組織進(jìn)行轉(zhuǎn)錄組測序，并對數(shù)據(jù)進(jìn)行質(zhì)量控制，得到clean reads.對照組中得到9 571 465條clean reads，GC含量為68.59%；感染組（3 h）中得到13 277 647條clean reads，GC含量為68.08%；感染組（6 h）中得到12 868 934條clean reads，GC含量為70.44%.分析牙鲆對照組與感染組中miRNA的堿基組成發(fā)現(xiàn)，胞嘧啶所占比例均最低.利用RNA-seq測序技術(shù)及生物信息學(xué)分析，得到牙鲆感染遲鈍愛德華氏菌過程中頭腎組織的差異表達(dá)miRNA結(jié)果為：感染組（3 h）有12個差異表達(dá)miRNA，其中2個上調(diào)，10個下調(diào)；感染組（6 h）有32個差異表達(dá)miRNA，其中20個上調(diào)，12個下調(diào).牙鲆被遲鈍愛德華氏菌感染的時間不同，miRNA的反應(yīng)程度也不同，感染6 h時應(yīng)激反應(yīng)處于高峰，預(yù)測miRNA可能參與了牙鲆的穩(wěn)態(tài)調(diào)節(jié).選取與模式識別受體（TLR、MR）信號通路相關(guān)的6個差異表達(dá)miRNA進(jìn)行實時熒光PCR驗證，所得結(jié)果與轉(zhuǎn)錄組測序結(jié)果一致，表明轉(zhuǎn)錄組測序結(jié)果可信.感染組（6 h）的3個miRNA均發(fā)生上調(diào)，這3個miRNA的靶基因與TLR信號通路有關(guān)，可能與牙鲆的免疫負(fù)調(diào)節(jié)相關(guān)，miRNA的表達(dá)上調(diào)，它對應(yīng)的靶基因TLR的表達(dá)應(yīng)下調(diào).本課題組前期研究結(jié)果顯示，牙鲆感染遲鈍愛德華氏菌過程中，感染6 h時多數(shù)TLR基因處于表達(dá)峰值.

miRNA在抗細(xì)菌的先天免疫反應(yīng)中起著非常重要的負(fù)反饋調(diào)控作用，這種調(diào)控可以確保免疫系統(tǒng)發(fā)生免疫反應(yīng)時不會出現(xiàn)過激反應(yīng).TLR信號必須受到嚴(yán)格控制以避免產(chǎn)生過度炎癥，并且在組織感染或損傷后允許組織修復(fù)恢復(fù)穩(wěn)態(tài).miRNA是TLR信號的重要控制因子，部分miRNA由固有免疫細(xì)胞中的TLR激活誘導(dǎo)，同其他miRNA一起通過靶向TLR通路中的信號分子miRNA起作用[15].TLR信號通路中存在著一些具有重要調(diào)控作用的調(diào)節(jié)因子，其調(diào)控可以是正向也可以是負(fù)向的，miRNA還可以通過作用于這些調(diào)節(jié)因子間接發(fā)揮對TLR信號通路的調(diào)控作用.miRNA被證明是固有免疫和適應(yīng)性免疫之間的重要連接，它們的失調(diào)可能在炎癥疾病的發(fā)病機(jī)制中起到重要作用[16-22].本研究發(fā)現(xiàn)，差異表達(dá)miRNA對應(yīng)的靶基因所參與的信號通路大部分跟生物過程以及細(xì)胞組分中的膜成分有關(guān)，如生長發(fā)育、信號轉(zhuǎn)導(dǎo)、細(xì)胞通訊、質(zhì)膜等.生物過程中候選靶基因富集程度最高的GO條目是生物調(diào)節(jié)，富集程度高的通路大多跟免疫調(diào)節(jié)以及分裂分化等生物學(xué)過程相關(guān)聯(lián)，質(zhì)膜具有識別和傳遞信息的功能，以上結(jié)果可以為牙鲆免疫基因的篩選和免疫系統(tǒng)的研究提供理論依據(jù).這些差異表達(dá)miRNA及其靶基因可能與牙鲆抗遲鈍愛德華氏菌感染過程密切相關(guān)，對研究牙鲆的免疫調(diào)控機(jī)制具有重要指導(dǎo)意義，其所提供的差異表達(dá)量變化的信息可以為后續(xù)研究牙鲆應(yīng)對遲鈍愛德華氏菌刺激時靶基因的篩選以及靶基因所起作用提供理論依據(jù).