陳昕迪 王騰宇 劉春霞 王文龍*

(1.內(nèi)蒙古農(nóng)業(yè)大學(xué) 獸醫(yī)學(xué)院/農(nóng)業(yè)農(nóng)村部動物疾病臨床診療技術(shù)重點實驗室，呼和浩特 010010;2.內(nèi)蒙古農(nóng)業(yè)大學(xué) 生命科學(xué)學(xué)院，呼和浩特 010010)

捻轉(zhuǎn)血矛線蟲是一種寄生于反芻動物消化道內(nèi)強致病性寄生蟲，可引起寄生宿主貧血和慢性消耗性臨床癥狀，嚴重可導(dǎo)致病畜死亡[1]。目前，捻轉(zhuǎn)血矛線蟲病的防治主要依靠伊維菌素和阿苯達唑等化學(xué)藥物，Dahourou等[2]對農(nóng)場進行流行病學(xué)調(diào)查發(fā)現(xiàn)，在布基納法索地區(qū)的傳統(tǒng)養(yǎng)殖場中，綿羊胃腸道寄生蟲的總體流行率為91%。早在2000年，國際獸醫(yī)局(International epizootic office, OIE)公布的全球調(diào)查數(shù)據(jù)顯示，不科學(xué)用藥已造成全球超過80%的國家存在捻轉(zhuǎn)血矛線蟲多重耐藥的流行情況，給畜牧業(yè)造成了巨大的經(jīng)濟損失[3]。但捻轉(zhuǎn)血矛線蟲產(chǎn)生耐藥性的調(diào)控機制還處于初探階段，為了防止耐藥性的產(chǎn)生并指導(dǎo)養(yǎng)殖戶合理用藥，挖掘耐藥性候選調(diào)控因子至關(guān)重要。

隨著近年來lncRNA在生命科學(xué)領(lǐng)域被高度關(guān)注，目前在人類和小鼠等模式生物中l(wèi)ncRNA的研究已較為深入[4]，長鏈非編碼RNA(Long non-coding RNA，lncRNA)是一類長度大于200 nt且很少具備編碼能力的RNA[5]，廣泛參與生物體的生長發(fā)育、免疫反應(yīng)以及細胞凋亡與分化等生物過程[6]。李騰等[7]在對微小隱孢子蟲的研究中發(fā)現(xiàn)，當微小隱孢子蟲感染人結(jié)直腸癌細胞后，lncRNA NEAT1的表達顯著上調(diào)，并通過調(diào)控細胞IL-8的轉(zhuǎn)錄,參與宿主細胞抗微小隱孢子蟲免疫反應(yīng)。王振勛等[8]通RNA-蛋白沉降、Western blotting等檢測手段，對弓形蟲慢性感染的小鼠腦組織樣本研究并得出結(jié)論，lncRNA 102796通過與靶基因oprd1結(jié)合抑制小膠質(zhì)細胞的增殖，影響細胞周期，進而導(dǎo)致神經(jīng)細胞的損傷。Gutschner等[9]在對lncRNA與癌癥發(fā)生的相關(guān)研究中發(fā)現(xiàn)，一些lncRNA表現(xiàn)出明顯的進化保守性和調(diào)控能力，通過表達水平上的差異和某些特異性lncRNA超表達來影響腫瘤細胞增殖，逃避生長抑制，進而在癌變過程中發(fā)揮關(guān)鍵作用，此類lncRNA可作為癌癥早期診斷及預(yù)后分子標記，并為后期藥物的開發(fā)提供新靶點。Zhou等[10]通過對乳腺癌細胞mTOR通路水平的檢測發(fā)現(xiàn)，lncRNA CCAT2可以通過激活mTOR信號通路降低乳腺癌細胞對5-氟尿嘧啶的耐藥性，進一步揭示lncRNA在解決疾病耐藥性問題中的重要作用。隨著新發(fā)現(xiàn)的lncRNA呈現(xiàn)指數(shù)型增長，人類疾病中大量耐藥性相關(guān)lncRNA被發(fā)掘出來，但在寄生蟲的耐藥性相關(guān)研究中，非編碼RNA的研究仍較為遲緩，尤其是捻轉(zhuǎn)血矛線蟲lncRNA的研究還相對滯后。

鑒于此，本研究通過對捻轉(zhuǎn)血矛線蟲伊維菌素敏感蟲株和耐藥蟲株的轉(zhuǎn)錄組測序，比較轉(zhuǎn)錄本間差異水平并根據(jù)編碼能力對lncRNA表達譜進行預(yù)測和篩選，以期從lncRNA分子層面探究耐藥蟲體在伊維菌素藥物作用下的抵抗機制，為進一步闡釋捻轉(zhuǎn)血矛線蟲伊維菌素耐藥調(diào)節(jié)機制及有效干擾靶點的開發(fā)提供新思路。

1 材料與方法

1.1 試驗蟲株

捻轉(zhuǎn)血矛線蟲伊維菌素耐藥蟲株采自內(nèi)蒙古自治區(qū)烏蘭浩特市科爾沁右翼前旗，通過對捻轉(zhuǎn)血矛線蟲感染強度EPG≥800的試驗羊群進行兩倍劑量的二次給藥后，經(jīng)糞便蟲卵檢查和對照試驗調(diào)查捻轉(zhuǎn)血矛線蟲感染情況，同時采用蟲卵減少試驗進行復(fù)檢測定感染強度后確定感染蟲株為伊維菌素耐藥蟲株[11]。敏感蟲株通過幼蟲發(fā)育抑制試驗測定敏感蟲株伊維菌素LD50敏感閾值低于9 ng/mL，且正常劑量伊維菌素體外驅(qū)蟲后的糞便蟲卵減少率大于90%，進而確定為捻轉(zhuǎn)血矛線蟲伊維菌素敏感蟲株[12]。對選取的綿羊進行剖檢，采集捻轉(zhuǎn)血矛線蟲成蟲并逐條鑒別雌雄，分裝至凍存管并儲存于液氮罐中。

1.2 成蟲cDNA文庫構(gòu)建及高通量測序

分別選取捻轉(zhuǎn)血矛線蟲伊維菌素敏感蟲株和耐藥蟲株雌蟲，雄蟲各1管，每管50條，并將雌雄蟲混合，使用Trizol(Invitrogen，美國)提取蟲體的Total RNA，并對RNA純度及濃度進行分析。采用去核糖體試劑盒(New england biolabs)對提取出的Total RNA進行過濾去除，并將RNA隨機打斷成短片段，再以片斷化后的RNA為模板，用六堿基隨機引物(Random hexamers)合成cDNA第一鏈和第二鏈。經(jīng)過Qiaquick PCR(Qiagen)試劑盒純化洗脫并進行末端修復(fù)、加堿基A，加測序接頭，降解第二條鏈后進行凝膠電泳檢測。將構(gòu)建好的文庫使用Illumina HiSeqTM4 000平臺進行雙末端測序。

1.3 轉(zhuǎn)錄組數(shù)據(jù)質(zhì)控及l(fā)ncRNA識別獲取

利用Fastp[13]軟件對Raw reads進行質(zhì)控，過濾去除低質(zhì)量Reads。使用Bowtie2[14]將較短的Reads比對GeneBank中的核糖體數(shù)據(jù)庫，去除比對上核糖體的Reads。采用HISAT2[15]軟件將保留下來的Unmapped reads比對到參考基因組(GenBank：GCA_000469685.2)，并與RNA-seq測序數(shù)據(jù)進行全局和局部的比對。

依據(jù)比對全部定位到基因組上的Reads(Total_mapped reads)的結(jié)果，計算參考基因組中的分布位置。將比對到基因組的區(qū)域分為外顯子區(qū)(Exon)、內(nèi)含子區(qū)(Intron)和基因間區(qū)(Intergenic)，以 Reads 在參考基因上的分布情況評價 mRNA 打斷的隨機程度。使用Stringtie進行重構(gòu)轉(zhuǎn)錄本，保留轉(zhuǎn)錄本長度≥200 bp且外顯子數(shù)目≥2的Reads。

1.4 lncRNA鑒定及差異表達分析

通過CPC2[16]和CNCI[17]2個軟件對新轉(zhuǎn)錄本的編碼能力進行分析，取交集作為候選lncRNAs。采用Stringtie軟件校正轉(zhuǎn)錄本長度，計算lncRNA的FPKM(Fragments per kilobase of transcript per million mapped reads)值。使用EdgeR[18]對Reads count進行標準化分析，計算假設(shè)檢驗概率(P-value)并校正FDR值，基于差異分析結(jié)果篩選FDR<0.05且|log2FC|>1的lncRNA為顯著差異lncRNA。

1.5 差異表達lncRNA的靶基因預(yù)測及功能注釋

在lncRNA的順式調(diào)控中，位于上游的lncRNA可能與啟動子或者共表達基因的其他順式調(diào)控元件有交集，進而在轉(zhuǎn)錄或者轉(zhuǎn)錄后對基因的表達進行調(diào)控，利用RNAplex[19]尋找2個長鏈RNA之間的相互作用，預(yù)測反義lncRNA與mRNA之間的互補結(jié)合，同時使用該程序中的ViennaRNA包[20]，根據(jù)熱力學(xué)結(jié)構(gòu)計算最小自由能來預(yù)測最佳堿基配對關(guān)系。

將lncRNA的靶基因應(yīng)用KOBAS軟件向GO數(shù)據(jù)庫映射，同時進行靶基因的KEGG富集分析，顯著富集的Pathway閾值為Q value≤0.05，注釋富集結(jié)果以分類二級柱狀圖和圈圖的形式展現(xiàn)。

1.6 lncRNA的ceRNA調(diào)控網(wǎng)絡(luò)構(gòu)建及連通性分析

采用Target Finder軟件將FDR<0.05且|log2FC|>1 的lncRNA和P<0.05且|log2FC|>1的miRNA依據(jù)ceRNA調(diào)控網(wǎng)絡(luò)關(guān)系，構(gòu)建lncRNA-miRNA調(diào)控網(wǎng)絡(luò)，篩選ceRNA中關(guān)鍵miRNA反應(yīng)元件(MicroRNA response elements, MRE)。

1.7 耐藥相關(guān)通路中的ceRNA調(diào)控網(wǎng)絡(luò)及連通性分析

根據(jù)Lam等[21]的研究方法對lncRNA和mRNA進行關(guān)聯(lián)分析，結(jié)合miRNA-lncRNA-mRNA三者間的調(diào)控關(guān)系使用Cytoscape構(gòu)建ceRNA調(diào)控網(wǎng)絡(luò)。依據(jù)lncRNA、mRNA、circRNA與其存在靶向調(diào)控關(guān)系的miRNA分子數(shù)量進行連通性分析，并對連通網(wǎng)絡(luò)進行可視化。根據(jù)KEGG數(shù)據(jù)庫富集結(jié)果構(gòu)建Pathway-Gene網(wǎng)絡(luò)，利用Igraph對耐藥相關(guān)通路中篩選出的lncRNA-miRNA-mRNA進行可視化分析。

1.8 差異表達lncRNA的熒光定量PCR驗證

為保證測序數(shù)據(jù)的可靠性，從篩選出的DE lncRNA中隨機挑選10個lncRNA進行qRT-PCR驗證，每組樣品設(shè)置3個重復(fù)，引物序列見表1。選取GAPDH作為內(nèi)參基因，用Primer 6.0設(shè)計lncRNA引物序列并送至Sangon合成，采用2-ΔΔCt法計算相對表達量。

表1 lncRNA熒光定量相關(guān)引物Table 1 lncRNA fluorescence quantitative primers

2 結(jié)果與分析

2.1 RNA-seq數(shù)據(jù)質(zhì)量評估

捻轉(zhuǎn)血矛線蟲伊維菌素敏感株和耐藥株提取總RNA后，得到各樣本濃度均大于450 ng/μL，總量大于3 μg且無DNase等雜質(zhì)污染和降解，OD260/OD280的比值為1.8～2.1。電泳結(jié)果可觀察到在18 S和28 S處有條帶且在凝膠樣圖譜的相應(yīng)位置正常起峰，表明RNA完整度良好，符合建庫上機和qRT-PCR的質(zhì)量要求(圖1)。敏感株和耐藥株分別獲得12.99和11.99 G的原始Raw data，兩組數(shù)據(jù)的Q20占比分別為96.54%和97.72%，Q30占比分別為90.56%和93.44%。在進行深度過濾后分別得到12.71 G和11.92 G的Clean data，其中比對上參考基因組唯一位置的Reads分別為73.62%和70.09%，Reads比對基因各部分均勻，且堿基組成平衡，質(zhì)量較高(圖2)，進一步保證后續(xù)分析的準確性。

圖1 捻轉(zhuǎn)血矛線蟲敏感(a)和耐藥(b)RNA凝膠樣譜圖和電泳圖Fig.1 H. contortus sensitive (a) and resistant (b) RNA gel-like spectra and electropherograms

圖2 捻轉(zhuǎn)血矛線蟲敏感蟲株(a)和耐藥蟲株(b)堿基分布圖Fig.2 Base distribution map of H. contortus sensitive (a) and resistant strains (b)

2.2 lncRNA及其靶基因的預(yù)測與鑒定

基于Stringtie重構(gòu)的5 374條轉(zhuǎn)錄本，使用CPC2和CNCI 2個軟件分別預(yù)測到1 141條和1 180 條lncRNA，共靶向到4 377個mRNA。兩款軟件預(yù)測到的沒有編碼潛能的轉(zhuǎn)錄本交集篩選出858個新lncRNA(圖3)，根據(jù)ceRNA競爭調(diào)控機制共靶向到474個mRNA。結(jié)合新預(yù)測到的lncRNA在基因組上的位置，統(tǒng)計到基因間長鏈非編碼RNA(Intergenic lncRNAs)400個、雙向長鏈非編碼RNA(Bidirectional lncRNAs)21個、內(nèi)含子長鏈非編碼RNA(Intronic lncRNAs)11個、反義長鏈非編碼RNA(Antisense lncRNAs)63個和正義長鏈非編碼RNA(Sense overlapping lncRNAs)245個。

圖3 捻轉(zhuǎn)血矛線蟲蟲體lncRNA預(yù)測韋恩圖Fig.3 Venn diagram of lncRNA prediction of H. contortus

2.3 顯著差異表達的lncRNA篩選

根據(jù)RNA-seq獲得的數(shù)據(jù)，將所有比對預(yù)測到的lncRNA進行FPKM值和Count值統(tǒng)計后發(fā)現(xiàn)，在敏感蟲株和耐藥蟲株中分別有65.85%和81.95%的lncRNA的FPKM值≥1(圖4)，說明lncRNA在捻轉(zhuǎn)血矛線蟲耐藥性的產(chǎn)生過程中呈現(xiàn)中高豐度表達且處于上升趨勢，且在所有鑒定到的lncRNA中MSTRG.9400.1的表達量最高。

圖4 捻轉(zhuǎn)血矛線蟲敏感株和耐藥株蟲體 lncRNA表達豐度小提琴圖Fig.4 Violin plot of lncRNA expression abundance ofH. contortus sensitive and resistant strains

伊維菌素敏感蟲株和耐藥蟲株預(yù)測到的858個新lnRNA均存在不同程度的差異表達，將FDR<0.05且|log2FC|≥1作為進一步顯著性篩選閾值，共獲得205個顯著差異表達的lncRNA，其中102個lncRNA上調(diào)，最顯著上調(diào)的3個lncRNA分別是MSTRG.13692.1、MSTRG.8350.1和MSTRG.3950.2；103個lncRNA下調(diào)，其中表達量最顯著下調(diào)3個lncRNA是MSTRG.10930.1、MSTRG.7170.1和MSTRG.3121.1(圖5)。

圖5 捻轉(zhuǎn)血矛線蟲敏感株和抗藥株之間 顯著差異表達lncRNA火山圖Fig.5 Volcano map of lncRNA differentially expressed between H. contortus sensitive and resistant strains

2.4 差異lncRNA的順式調(diào)控基因富集分析

將lncRNA與mRNA進行關(guān)聯(lián)分析，其中723個lncRNA通過順式調(diào)控作用對1 320個上游和下游基因進行潛在調(diào)控，共匹配到了1 566個調(diào)控關(guān)系。將lncRNA通過順式調(diào)控的基因向GO數(shù)據(jù)庫映射，得到每個GO功能富集情況的基因列表，應(yīng)用超幾何檢驗，獲得顯著富集的GO條目。敏感蟲株和耐藥蟲株的GO富集結(jié)果顯示：lncRNA所靶向到的基因在生物學(xué)過程中注釋到了代謝過程(GO:0008152)95個，細胞過程(GO:0009987)89個，對刺激的反應(yīng)(GO:0050896)24個共計18條GO terms，表明靶基因在受到藥物的刺激下廣泛參與到了各類化合物的代謝；在細胞組分中注釋到了細胞膜結(jié)構(gòu)(GO:0016020)47個，細胞器(GO:0043226)45個和細胞部分(GO:0005623)61個等9條GO terms(圖6)，細胞膜等膜結(jié)構(gòu)被大量注釋，說明膜組織在捻轉(zhuǎn)血矛線蟲耐藥性的產(chǎn)生過程中發(fā)揮重要作用。在分子功能中被大量注釋到催化活性(GO:0003824)106個，各類結(jié)合(GO:0005488)79個，轉(zhuǎn)運活動(GO:0005215)14個等8條GO terms；在分子功能注釋到的一些lncRNA靶基因中，很多基因具備離子通道催化活性劑、轉(zhuǎn)移酶活性等分子功能，大量的功能性物質(zhì)被轉(zhuǎn)運，廣泛參與藥物脅迫應(yīng)激和化合物代謝。

圖6 lncRNA cis調(diào)控的上下游基因GO功能注釋Fig.6 GO functional annotation of lncRNA cis-regulated upstream and downstream genes

KEGG通路富集顯示，預(yù)測得到的lncRNA所靶向的上下游基因共注釋到291條通路。其中注釋基因數(shù)最多的前10位分別為代謝途徑(KO01100)、光傳導(dǎo)(KO04745)、黏著力(KO04510)、肌動蛋白細胞骨架的調(diào)節(jié)(KO04810)、內(nèi)吞作用(KO04144)、血小板活化(KO04611)、卵母細胞減數(shù)分裂(KO04114)、甘油磷脂代謝(KO00564)、黏附連接(KO04520)、多巴胺能突觸(KO04728)和孕酮介導(dǎo)的卵母細胞成熟(KO04914)，由KEGG數(shù)據(jù)庫通路富集圈圖(圖7)可以看出代謝相關(guān)通路被大量富集，其中細胞色素P450及對異生物質(zhì)的代謝，視黃醇代謝和甘油磷脂代謝等耐藥相關(guān)代謝通路也被大量富集，由此進一步闡述在耐藥發(fā)生過程中捻轉(zhuǎn)血矛線蟲產(chǎn)生較大的代謝功能障礙。

圖7 lncRNA cis調(diào)控的上下游基因KEGG數(shù)據(jù)庫富集圈圖Fig.7 KEGG database enrichment circle map of lncRNA cis-regulated upstream and downstream genes

2.5 lncRNA的ceRNA調(diào)控網(wǎng)絡(luò)構(gòu)建

根據(jù)ceRNA競爭性調(diào)節(jié)機制，蛋白編碼基因mRNA、長鏈非編碼RNA(lncRNA)、假基因(Pseudogenes)、環(huán)狀RNA(circRNA)均是通過共同的MREs結(jié)合調(diào)控相同的miRNA，由此構(gòu)建lncRNA-miRNA調(diào)控網(wǎng)絡(luò)(圖8)共預(yù)測到205個lncRNA靶向結(jié)合371個miRNA，組成3 227組靶向調(diào)控關(guān)系組合。同一個miRNA可以靶向結(jié)合多個lncRNA，在眾多的調(diào)控組合中novel-m0251-3p能夠結(jié)合14個lncRNA且數(shù)量最多。miR-203-y、miR-2184-x、miR-219-z、miR-32-x和miR-92-y均能結(jié)合到12個lncRNA，具有較強的結(jié)合能力，但有約74.25%的miRNA靶向到的lncRNA小于6個。lncRNA中MSTRG.11468.4、MSTRG.11468.5、MSTRG.11468.8和MSTRG.11468.9可以靶向結(jié)合到的miRNA數(shù)量較多，分別為322 422和22個miRNA，其中89.54%的lncRNA能夠結(jié)合到miRNA的數(shù)量均小于10個。

圖8 捻轉(zhuǎn)血矛線蟲敏感蟲株和耐藥蟲株lncRNA-miRNA調(diào)控網(wǎng)絡(luò)Fig.8 lncRNA-miRNA regulatory network of H. contortus sensitive and resistant strains

2.6 耐藥相關(guān)調(diào)控信號通路的ceRNA調(diào)控網(wǎng)絡(luò)

通過對伊維菌素敏感蟲株和耐藥蟲株差異mRNA與miRNA，miRNA與lncRNA的靶向調(diào)控關(guān)系預(yù)測結(jié)果進行l(wèi)ncRNA-circRNA-miRNA-mRNA連通性分析(圖9)，結(jié)果顯示在連通性前5的特定連通性節(jié)點上分別有118、62、33、22和9條調(diào)控網(wǎng)絡(luò)，且MSTRG.11468、MSTRG.14365.2和MSTRG.5594.2與circRNA-miRNA-mRNA網(wǎng)絡(luò)具有較強的連通性，進一步表明其具有更強的潛在調(diào)控能力。通過KEGG富集結(jié)果對FoxO signaling pathway、cAMP signaling pathway、MAPK signaling pathway、AMPK signaling pathway和PI3K-Akt signaling pathway等耐藥相關(guān)信號通路與連通性高的節(jié)點網(wǎng)絡(luò)進行關(guān)聯(lián)分析，同時篩選出部分假定蛋白(表2)，根據(jù)上述分析結(jié)果構(gòu)建的耐藥性相關(guān)通路lncRNA-miRNA-mRNA調(diào)控網(wǎng)絡(luò)(圖10)可以看出lncRNA所靶向到的mRNA單一性較高，結(jié)合miRNA多靶效應(yīng)具有較強指向性的特點進一步篩選出ceRNA調(diào)控網(wǎng)絡(luò)中的關(guān)鍵調(diào)節(jié)因子。

圖9 連通性前10的ceRNA連通性?；鶊DFig.9 Sankey diagram of the top ten ceRNA connectivity

表2 捻轉(zhuǎn)血矛線蟲敏感蟲株與耐藥蟲株差異lncRNA序列互補靶標的功能注釋Table 2 Functional annotations of complementary targets with sequences of lncRNA differentially expressed between H. contortus sensitive and resistant strains

圖10 耐藥相關(guān)通路中顯著差異基因的lncRNA-miRNA-mRNA調(diào)控網(wǎng)絡(luò)關(guān)系Fig.10 The lncRNA-miRNA-mRNA regulatory network relationship of significantly differentially expressed genes in drug resistance-related pathways

2.7 lncRNA的表達驗證

利用qRT-PCR技術(shù)，在伊維菌素敏感蟲株和耐藥蟲株比較組中隨機挑選10個(上)下調(diào)表達量差異水平較大的DE lncRNA進行表達趨勢的驗證。如圖11所示，10個DE lncRNA在比較組中的變化情況均通過qRT-PCR結(jié)果在組內(nèi)真實表達出來，且與RNA-seq測序結(jié)果趨勢一致，進一步證實轉(zhuǎn)錄組測序結(jié)果的穩(wěn)定性和可靠性。

圖11 差異表達lncRNA qRT-PCR與RNA-seq比較分析Fig.11 Comparative analysis between qRT-PCR and RNA-seq of differentially expressed lncRNAs

3 討 論

反芻動物感染消化道線蟲特別是捻轉(zhuǎn)血矛線蟲給畜牧業(yè)生產(chǎn)造成較為嚴重的影響，由于無法及時準確篩查出已感染的羊群，大多數(shù)養(yǎng)殖戶選擇大劑量群體給藥進行驅(qū)蟲和預(yù)防[11]，致使羊群對伊維菌素和阿苯達唑等藥物產(chǎn)生耐藥性，此現(xiàn)象受到國際社會的廣泛關(guān)注[22]。相較于秀麗隱桿線蟲等模式生物，捻轉(zhuǎn)血矛線蟲的組學(xué)數(shù)據(jù)較為匱乏且研究相對滯后。目前，lncRNA已被證實在擬南芥(Arabidopsisthaliana)、蜜蜂球囊菌(Ascosphaeraapis)和羅伯茨綠僵菌(Metarhiziumrobertsii)[23-24]等真核生物中均參與到代謝和免疫等過程并調(diào)控特定蛋白的表達水平。有研究證實耐藥性的產(chǎn)生可能與藥物代謝的增加、藥物靶點的突變以及細胞凋亡和自噬存在密切聯(lián)系，同時信號通路的異常激活也可能導(dǎo)致耐藥強度的改變。非編碼 RNA 能夠極大地增加調(diào)控的復(fù)雜性并起到微調(diào)轉(zhuǎn)錄組的作用，并能快速調(diào)整蛋白質(zhì)組以響應(yīng)刺激[25]。

本研究以伊維菌素敏感蟲株和耐藥蟲株為研究對象，通過高通量測序技術(shù)將獲得的Raw data進行重新組裝，經(jīng)預(yù)測發(fā)現(xiàn)捻轉(zhuǎn)血矛線蟲中不存在已知和已存在的lncRNA，根據(jù)序列在基因組的位置信息，長度大小和編碼能力鑒定到858個新lncRNA，可以與474個mRNA靶向結(jié)合，鑒定到的lncRNA以|log2FC|≥1且FDR<0.05為顯著性閾值，存在205個顯著差異lncRNA。本研究通過表達量分析發(fā)現(xiàn)，預(yù)測到的所有新lncRNA在敏感蟲株與耐藥蟲株中均存在表達，但顯著差異的lncRNA數(shù)目相較于總數(shù)目而言占比較少，且高表達基因差異水平未呈現(xiàn)斷崖式變化，由此推測捻轉(zhuǎn)血矛線蟲在對伊維菌素產(chǎn)生抗性過程中僅通過較少數(shù)量的調(diào)控網(wǎng)絡(luò)便可對抗性水平進行調(diào)節(jié)。依據(jù)GC含量高的基因組間區(qū)存在更多的編碼基因的假說[26]，對轉(zhuǎn)錄本數(shù)據(jù)中獲得的大量lncRNA座進行預(yù)測發(fā)現(xiàn)，lncRNA中的GC含量明顯低于mRNA，本研究結(jié)果與該假說猜想一致。Azlan等[27]通過對埃及伊蚊(Aedesaegypti)lncRNA表達譜鑒定發(fā)現(xiàn)，lncRNA呈現(xiàn)低表達、低 GC 含量、長度短和低保守性等特點，在本研究所有預(yù)測到的新lncRNA中，有262個lncRNA的FPKM值小于1，其余大部分均呈現(xiàn)出中高豐度表達，但相較于mRNA的表達情況而言，lncRNA依舊處于低表達水平，與前人研究結(jié)果一致，同時也說明lncRNA的表達情況往往與鄰近的反義蛋白編碼基因相關(guān)。在lncRNA-mRNA的調(diào)控網(wǎng)絡(luò)中有356組調(diào)控關(guān)系是上調(diào)lncRNA抑制mRNA下調(diào)表達，378組下調(diào)lncRNA從而上調(diào)mRNA表達，這種反向調(diào)控模式占所有l(wèi)ncRNA-mRNA調(diào)控網(wǎng)絡(luò)的46.9%，lncRNA與編碼基因間通過正向和反向調(diào)控構(gòu)成了遠超本身數(shù)量的復(fù)雜機體微調(diào)網(wǎng)絡(luò)，再次說明同一共有l(wèi)ncRNA廣泛參與到捻轉(zhuǎn)血矛線蟲伊維菌素耐藥性的產(chǎn)生過程。

KEGG富集結(jié)果顯示，被富集到較多的通路是代謝和轉(zhuǎn)運相關(guān)通路，其中顯著差異基因約有三分之二富集到代謝過程。Siddiqui等[28]發(fā)現(xiàn)在大多數(shù)給藥后或耐藥生物的通路富集結(jié)果中其他藥物代謝酶、藥物代謝細胞色素 P450、谷胱甘肽代謝、甘油磷脂代謝、糖酵解/糖異生和異生物質(zhì)代謝等通路均存在富集情況，這些代謝途徑與信號通路存在緊密的聯(lián)系，不同的信號通路都發(fā)揮著不同的調(diào)控作用。在腫瘤細胞的研究中耐藥細胞與藥物敏感細胞相比，許多脂質(zhì)都發(fā)生了顯著變化[29]，本研究中甘油磷脂代謝等相關(guān)脂質(zhì)通路被顯著富集，甘油磷脂代謝作為重要的藥物代謝通路已被證實在蠕蟲中藥物會通過氨基酸和甘油磷脂產(chǎn)生的小分子進行代謝[30]。根據(jù)功能預(yù)測分析，甘油磷脂代謝通路中的 F09G2.8作為一類磷脂水解酶酶和磷脂酶D病毒包膜結(jié)構(gòu)域蛋白在MSTRG.14103.3的下調(diào)影響下低豐度表達，由此推測MSTRG.14103.3可靶向調(diào)節(jié)脂類小分子蛋白參與捻轉(zhuǎn)血矛線蟲伊維菌素的代謝過程從而影響耐藥性。通過對lncRNA富集情況進行統(tǒng)計發(fā)現(xiàn)cAMP信號通路和MAPK信號通路同時富集到MSTRG.3141.1及其靶基因PRKACA。Moody等[31]發(fā)現(xiàn)除了MAPK和PI3K-AKt等信號通路外，PRKACA在曲妥珠單抗耐藥的乳腺癌樣本中上調(diào)表達，并在人表皮生長因子受體2的治療中拯救了細胞，從而使細胞對拉帕替尼治療敏感，其中PRKACA的表達情況與本研究中的結(jié)果一致。在人類肺癌對鉑類藥物的抗藥性研究中表明，PRKACA作為一類cAMP依賴性蛋白，ATP轉(zhuǎn)化為的環(huán)AMP可以與cAMP依賴蛋白結(jié)合調(diào)節(jié)細胞轉(zhuǎn)運激酶、鳥嘌呤交換因子和感覺信號向環(huán)苷酸門控通道轉(zhuǎn)導(dǎo)信號，進而提高細胞對藥物的敏感性[32]。同樣被注釋到的MAPK信號通路可以將細胞外信號傳導(dǎo)至細胞核內(nèi)，在脅迫應(yīng)答和細胞增殖與凋亡等方面起到關(guān)鍵作用[33]，由此可以推測在伊維菌素耐藥性產(chǎn)生過程中MSTRG.3141.1的下調(diào)可能激活下游轉(zhuǎn)錄因子與cAMP依賴性蛋白啟動子的近端區(qū)域相互作用，并將信號傳遞到細胞內(nèi)以促進靶基因PRKACA表達。

ABC轉(zhuǎn)運蛋白是一類專一性外排轉(zhuǎn)運蛋白，通過跨細胞膜轉(zhuǎn)運分子逆濃度梯度調(diào)節(jié)內(nèi)源性和外源性物質(zhì)的細胞水平，并參與多種生物過程。ABC轉(zhuǎn)運蛋白的過度表達與大環(huán)內(nèi)酯類藥物的抗性存在密切聯(lián)系[34]。P-糖蛋白(Pgp)是ABC轉(zhuǎn)運蛋白超家族中重要的功能蛋白，在人類的多種疾病和寄生蟲的耐藥性研究中均出現(xiàn)明顯異常表達，是重要的耐藥相關(guān)蛋白[35-36]。在捻轉(zhuǎn)血矛線蟲耐伊維菌素的研究中發(fā)現(xiàn)，Pgp基因均出現(xiàn)表達情況的異常改變及SNP位點的突變[37]，本研究中ABC轉(zhuǎn)運途徑中富集到Pgp-1和ced-7等多個基因，但通過ceRNA結(jié)合能力預(yù)測Pgp-1并未靶向到lncRNA和miRNA，而ced-7卻靶向到大量miRNA，由此推測機體可能通過ABC轉(zhuǎn)運蛋白家族中的其他轉(zhuǎn)運分子進行ceRNA網(wǎng)絡(luò)調(diào)節(jié)進而達到藥物逃逸的作用。lncRNA作為研究起步較晚的一類非編碼RNA，在捻轉(zhuǎn)血矛線蟲伊維菌素耐藥相關(guān)的研究中尚未見相關(guān)報道，由于各物種間lncRNA序列的保守性較差，無法確定lncRNA的轉(zhuǎn)錄確切位置。本研究通過編碼潛能預(yù)測獲取的lncRNA表達譜并與靶基因和miRNA進行ceRNA網(wǎng)絡(luò)構(gòu)建，盡可能去除無效的非編碼RNA分子。對于篩選出的耐藥相關(guān)且特異性表達的lncRNA及其調(diào)控組合，后續(xù)還需對蟲體進行干擾和過表達等方式進行生物學(xué)功能驗證，以期進一步闡述非編碼RNA在捻轉(zhuǎn)血矛線蟲耐藥性產(chǎn)生過程中發(fā)揮的重要作用。

4 結(jié) 論

通過RNA-seq技術(shù)，獲取到捻轉(zhuǎn)血矛線蟲伊維菌素敏感蟲株和耐藥蟲株lncRNA表達譜并進行全面分析，解析在耐藥性產(chǎn)生過程中l(wèi)ncRNA的數(shù)量、種類以及表達譜間存在的明顯差異，揭示部分lncRNA可通過順式調(diào)控作用參與ceRNA競爭調(diào)節(jié)機制，并廣泛參與到藥物代謝和轉(zhuǎn)運等生物學(xué)過程。研究結(jié)果可為耐藥性關(guān)鍵調(diào)控分子的篩選和功能性研究提供重要參考，也為深度開發(fā)捻轉(zhuǎn)血矛線蟲耐藥性關(guān)鍵lncRNA及其潛在ceRNA調(diào)控元件提供數(shù)據(jù)支持。