鰻弧菌(Vibrio anguillarum)核酸適配體篩選中高頻序列的進化多樣性研究*

許 凈 劉慧敏 林筱鈞 鄭 江,,3① 江興龍,鄢慶枇 范云庭 湯學敏

(1. 集美大學水產學院 廈門 361021; 2. 鰻鱺現代產業技術教育部工程研究中心 廈門 361021; 3. 福建省水產生物育種與健康養殖工程研究中心 廈門 361021; 4. 福建省特種水產配合飼料重點實驗室 福清 350308)

核酸適配體是通過指數富集配體的系統進化技術(Systematic Evolution of Ligands by Exponential Enrichment), 即SELEX 篩選技術, 獲得的對靶目標有較好親和特異性的寡核苷酸分子(Ellingtonet al,1990)。核酸適配體具有靶目標范圍廣、體積小、易合成、親和力高、特異性強等優點, 在癌癥治療、微生物檢測、食品安全等領域中呈現出廣闊的應用前景(Hamulaet al, 2011; 于寒松等, 2015; Yazdian-Robatiet al, 2017; Liuet al, 2020)。SELEX 篩選中, 需要對獲得的篩選產物進行測序, 然后從測序結果中選擇部分序列進行親和特異性的驗證, 親和特異性較好的序列則被確認為核酸適配體序列。然而, 通過傳統克隆測序獲得的序列有成百上千條, 通過高通量測序獲得的序列更是高達數萬條, 因此, 如何從如此眾多的序列中精準、高效的選擇出相應的候選序列進行核酸適配體的驗證, 是一個亟待解決的問題。

高頻序列是測序結果中出現的頻率大于等于2次的序列, 研究表明, 相當多的高頻序列都是在篩選進化中具有競爭優勢的序列, 且大概率是對靶目標有較好親和力的核酸適配體序列(鄭江等, 2014), 而高頻序列的數量在總序列中的占比常常很低, 因此,對高頻序列進行分析研究, 研究其在篩選進化過程中的優勢序列及其多樣性變化, 將能有針對性地尋找到相應的核酸適配體序列, 從而能夠降低篩選的盲目性, 提高篩選的效率。

鰻弧菌(Vibrio anguillarum)是多種淡海水養殖品種的致病菌(饒穎竹等, 2016; 王鳳青等, 2018;Mohamadet al, 2019), 由其導致的疾病給水產養殖業造成巨大的經濟損失, 而對鰻弧菌進行快速的分析檢測則是其病害防治的前提和基礎。核酸適配體因其具有親和力高、特異性強等優點, 已廣泛應用于微生物檢測(Tanget al, 2013; Zhenget al, 2015; Yuet al,2019; Sadsriet al, 2020), 因此篩選鰻弧菌的核酸適配體, 對于鰻弧菌的快速檢測及其病害防治都具有積極意義。而研究其篩選過程中高頻序列的進化特點及其多樣性變化, 不僅能快速準確地獲得鰻弧菌核酸適配體的候選序列, 對于理解掌握篩選文庫的進化特點, 提高相關水產病原菌核酸適配體的篩選效率都具有積極意義。

本文選擇了13 輪鰻弧菌核酸適配體的篩選產物進行高通量測序, 采用生態學中常用的豐富度指數、均勻度指數、多樣性指數和優勢度指數對其中的高頻序列進行了多樣性分析, 并采用相對重要性指數(IRI)對高頻序列進行了歸納分類, 最后利用在線網站對IRI 較高的7 種核酸適配體的二級結構進行了預測模擬, 相關研究對于核酸適配體的篩選以及鰻弧菌的檢測防治都具有重要意義。

1 材料與方法

1.1 序列的來源及其種類、個體數的定義

本文中的高頻序列來自于以鰻弧菌為靶目標的SELEX 篩選, 具體流程如下: (1) 構建并合成初始隨機寡核苷酸文庫: 5′-TCAGTCGCTTCGCCGTCTCC TTC(N35)GCACAAGAG GGAGACCCCAGAGGG-3′,N35 為含35 個隨機堿基的寡核苷酸序列。(2) 孵育:將上述寡核苷酸文庫與鰻弧菌共同孵育30 min—2 h,使寡核苷酸分子與鰻弧菌結合; (3) 分離: 利用離心、洗滌、熱變性等分離技術, 將未結合的或結合較弱的寡核苷酸分子與鰻弧菌分離, 獲得能與鰻弧菌較好結合的寡核苷酸分子; (4) 擴增: 將能與鰻弧菌結合的寡核苷酸分子進行不對稱PCR 擴增, 此為完成一輪篩選, 該PCR 產物即為該輪的篩選產物, 也是下一輪篩選的次級文庫; (5) 重復: 重復步驟(1)到(4), 即次級文庫再次與鰻弧菌進行孵育、分離、擴增, 獲得相應的篩選產物; (6) 測序: 對篩選產物進行高通量測序, 從中選擇核酸適配體的候選序列進行后續的驗證和鑒定。我們共進行了15 輪的篩選, 選擇了其中的13 輪篩選產物, 即第1、2、3、4、5、6、7、9、11、12、13、14、15 輪的PCR 產物, 進行高通量測序, 選擇測序結果中的高頻序列進行分析研究。高通量測序由生工生物工程有限公司(上海)完成。

為了方便理解, 我們對本文中核酸序列的種類和個體數作如下定義或表述: 不同種類的序列是指堿基的組成、數量或排列順序不同的序列, 通常被簡稱為不同序列; 同一種序列在測序結果中如果出現多次, 它的出現頻率可以看作是這種序列的個體數,而每輪測序結果中的總序列數則可看作該輪測序獲得的總個體數。

1.2 高頻序列在總序列中的占比

統計每輪測序結果中高頻序列的個體數和每種高頻序列出現的頻率(即該高頻序列的個體數), 則該輪高頻序列個體數的占比(%)=(該輪高頻序列個體數/該輪測序序列總個體數)×100=(該輪所有高頻序列的頻率之和/該輪測序序列總個體數)×100, 分別計算出13 個測序輪中高頻序列的占比, 比較分析其變化規律。

1.3 高頻序列的多樣性分析

參考王旭娜等(2018)和黃雅琴等(2020)的方法,分別選用豐富度指數、均勻度指數、多樣性指數和優勢度指數, 來研究高頻序列的多樣性, 這些指數的原計算公式如下:

原公式中,S為所有物種的種類數,N為每輪物種的個體數量的總和,Ni為第i個物種的個體數,Pi=Ni/N,Nmax為優勢種的總個體數,Nt為各物種的總個體數。將上述公式應用于本文的高頻序列, 則公式中的S就是某測序輪中高頻序列的種類數,N為該輪中高頻序列的總個體數,Ni是該輪第i個高頻序列的個體數, 即第i個高頻序列的出現次數或頻率,Pi=Ni/N,Nmax為該輪中優勢高頻序列的個體數,Nt是該輪中高頻序列的總出現次數。因此, 根據上述公式, 采用Primer 6.0 軟件可計算出每一個測序輪中的上述四種指數。

1.4 高頻序列的相對重要性分析

優勢種的研究通常是先選擇一系列的站點, 然后計算每個站點中的各個物種的相對重要性指數(Index of relative importance, IRI), 再根據該指數對相應的物種進行判定分類(王雪輝等, 2010; 劉曉霞等,2016; 張鵬弛等, 2017; 粟麗等, 2018)。IRI 的原始計算公式為: IRI = (N+W)×F×10 000 (沃佳等, 2017; 丁朋朋等, 2019), 式中,N為某一種類的個體數占總個體數的百分率;W為某一種類的重量占總重量的百分率;F為某一種類出現的站數占調查總站數的百分率。在本文高頻序列的研究中, 我們將各個篩選輪或測序輪看作是原公式中的各個站點, 則總共13 個測序輪就對應著13 個站點, 某種高頻序列出現的次數或頻率就相當于該高頻序列的個體數, 則公式中N就是某種高頻序列在13 個測序輪中出現的總次數占13輪中全部高頻序列的出現總次數的百分比,F則為該高頻序列出現的輪數占全部輪數(即13 輪)的百分比。另外, 由于核酸序列重量與其序列長度成正比, 本文中的高頻序列的長度大多在82 nt 左右, 差別不大,因此公式中的重量百分比W就近似等于其數量百分比N, 因此IRI 公式就無需加和, 可只取數量百分比N, 由此相對重要性指數IRI 的公式可修正為IRI =N×F×10 000。然后參考相關文獻的判定方法(陳國寶等, 2007; 李顯森等, 2013; 粟麗等, 2018), 依據IRI的值對不同序列進行判斷分類, 具體如下: 當該序列的IRI≥500時, 該序列為優勢序列; 100≤IRI<500時,該序列為重要序列; 10≤IRI<100 時, 該序列為常見序列, 1≤IRI<10 時, 該序列為一般序列; IRI<1 時, 該序列為少見序列。

1.5 序列的二級結構模擬

利用RNA structure 網站(http://rna.urmc.rochester.edu/RNAstructureWeb/)進行在線預測, 核酸類型選擇為DNA, 其他參數為默認值, 最后選擇最大可能性的二級結構, 即MaxExpect 二級結構作為各序列的二級結構。

2 結果

2.1 篩選過程中高頻序列的種類數及其在總序列中的占比變化

高頻序列種類數的變化(圖1)與其占比的變化(圖2)相類似, 都是在第4、5 輪出現突然性的大幅上升,在第4 輪達到最大值, 之后快速回落后出現穩步增加,在最后一輪達到次高點。從總趨勢上看, 高頻序列的種類數在前6 輪出現較劇烈的波動, 之后隨著篩選的進行呈穩步增加的趨勢, 高頻序列的占比也呈現類似變化, 不過高頻序列占比的波動要更劇烈些, 后續的增加趨勢也要更明顯。

圖1 篩選過程中高頻序列的種類數變化Fig.1 Changes of high frequency sequences in the species during the selection

圖2 篩選過程中高頻序列在總序列中的占比變化Fig.2 The proportion changes of high frequency sequence in the total sequence during the selection

2.2 高頻序列的多樣性分析

反映多樣性的Shannon-Wiener 指數和豐富度指數與高頻序列的種類數和占比變化(如圖1 和圖2)也類似, 也是在篩選的第4 輪達到最高, 之后快速回落,然后再逐步上升, 說明隨著篩選進行, 高頻序列也在進化, 其豐富度和多樣性在經歷了前6 輪的劇烈波動后, 都隨著篩選進化過程而逐步提高。

Pielou 均勻度指數反映的是各物種個體數目分配的均勻程度, 指數越高, 說明個體間數目越均勻。圖3 顯示, 該指數在前6 輪在一定區間內呈現較大波動, 之后則趨于穩定, 并呈現逐漸降低趨勢, 說明隨著篩選進行, 均勻度在經歷了前6 輪的劇烈波動后,呈現越來越低的趨勢, 這也進一步說明篩選進化是向著某一方向進行的, 適應這一進化方向的高頻序列, 其個體數就會越來越多, 不適應這一方向的高頻序列, 其個體數就會越來越少, 從而呈現出越來越大的差異和不均勻性。

優勢度指數反映的是優勢種的優勢程度, 優勢度指數越高, 說明該輪優勢種的優勢越突出。圖3 顯示, 該指數在前6 輪波動較為劇烈, 說明優勢種的優勢度并不穩定, 其他種的高頻序列可能對優勢種呈現較大的競爭壓力, 之后隨著篩選進行, 優勢度指數也逐漸走高, 說明優勢種在篩選進化中的競爭力越來越強, 其在進化競爭中的優勢越來越明顯, 這也反映了篩選效果越來越好。

圖3 篩選過程中Margalef 豐富度指數、Shannon-Wiener 多樣性指數、Pielou 均勻度指數、Berger-Parker 優勢度指數的變化Fig.3 Changes of the Margalef richness index, Shannon-Wiener diversity index, Pielou evenness index, and Berger-Parker dominance index during the selection

2.3 高頻序列的相對重要性分析及其分類

對13 輪中的79 個高頻序列進行了相對重要性指數(IRI)的計算, 結果如表1。依據IRI 的大小可將這些高頻序列分為以下五種類型, 即優勢序列(IRI≥500)、重要序列(100 ≤IRI<500)、常見序列(10 ≤IRI<100)、一般序列(1≤IRI<10)和少見序列(IRI<1)。從表中可以看出, 優勢序列和重要序列分別只有1 種,約占全部高頻序列的1.27%, 常見序列、一般序列和少見序列, 則分別約占全部高頻序列的 8.86%、18.99%和69.62%。另外, 優勢序列的IRI 值高達8 707.25, 是重要序列的19 倍多, 這表明優勢序列在篩選進化中的競爭優勢是較為突出的。

通過表1 的分類可看出, 優勢序列、重要序列和常見序列應該是篩選進化中的相對優勢種群, 它們對靶目標鰻弧菌的競爭力大概率要高于一般序列和少見序列, 因此, 按照IRI 指數從高到低, 優先從優勢序列、重要序列中挑選相應的序列進行核酸適配體的驗證, 將能大大降低驗證的盲目性和隨機性, 提高篩選效率。后續選擇了IRI 最高的7 個高頻序列(H5、H1、H6、H12、H17、H28、H21)進行了親和特異性驗證, 發現這7 個高頻序列均是親和特異性較好的核酸適配體(相關數據另文發表)。

表1 高頻序列的相對重要性指數及其分類Tab.1 Index of relative importance of high frequency sequences and their classification

2.4 核酸適配體序列的二級結構

選擇已驗證確認的7 個核酸適配體(H5、H1、H6、H12、H17、H28、H21)進行了二級結構的模擬(圖4), 從圖中可看出, 這些序列均形成了多個大小不一的閉合環狀結構, 表2 對這些核酸適配體序列中的閉合環進行總結, 可以發現, 每個序列都含有2—6 個閉合環, 每個環由5—25 核苷酸組成, 根據每個環中的核苷酸數量可以把這些環大致分為三種:小型環(由5—8 個核苷酸構成)、中型環(由9—13 個核苷酸構成)、大型環(由17 個以上的核苷酸構成)。這些環的大小和數量很可能與核酸適配體的結合位點相匹配。

表2 核酸適配體序列的IRI 值及其二級結構中環的數量和大小Tab.2 Values of index of relative importance of aptamers and the numbers and sizes of their rings in their secondary structures

圖4 核酸適配體二級結構Fig.4 Secondary structures of the aptamers

3 討論

研究鰻弧菌核酸適配體篩選中高頻序列的進化多樣性對于揭示核酸適配體的篩選進化規律、提高核酸適配體的篩選效率都具有重要意義。本文分別從高頻序列的種類數、多樣性、相對重要性指數(IRI)和二級結構等四個方面對高頻序列的篩選進化特點進行了研究, 下面也從這四個方面對其進行分析討論。

3.1 篩選過程中高頻序列的種類數和占比變化

理論上, 高頻序列是在篩選競爭中占優勢的序列, 其種類數量和在總序列中的占比應該隨著篩選輪數的增加而逐步增加。但我們的研究卻發現, 在篩選的中間階段出現了高頻序列大量爆增、劇烈波動的現象, 其種類和占比都超過了最后一輪。出現這種劇烈波動現象的原因很可能和篩選條件有關。前6 輪的篩選是在較為寬松的條件下進行的(結合1—2 h, 洗滌1—2 次), 所以, 雖然第4、5 輪出現了高頻序列的大爆發, 但大部分高頻序列的親和力可能并不高, 第7 輪后為了去除親和力較弱的序列, 篩選條件逐漸嚴苛, 結合時間縮短為0.5 h, 洗滌次數增加為3 次, 從而導致相當多親和力不高的高頻序列被淘汰, 最后雖然高頻序列的數量和占比沒有達到最高, 但其中高親和力序列的數量大大增加, 文庫中高頻序列的質量得到了較大的改善。

3.2 篩選過程中高頻序列的進化多樣性

Margalef 豐富度指數、Shannon-Wiener 多樣性指數、Pielou 均勻度指數和Berger-Parker 優勢度指數,經常被用于動植物的多樣性研究。作為一種有效的分析和描述方法, 這些指數還被用于手足口病的流行性特征分析(Yanget al, 2020)、口腔微生物菌群的組成分析(Wolffet al, 2019)、以及土壤微生物多樣性的研究(夏圍圍等, 2014)。本文則采用這四個指數研究了核酸適配體篩選過程中高頻寡核苷酸序列的進化多樣性, 結果表明, 隨著篩選進行, 高頻序列的豐富度、多樣性和優勢度都呈現逐步增加趨勢, 而均勻度出現降低。這一結果很可能是因為在篩選競爭中, 在總序列中占比90%以上的、大量親和力較弱的非高頻序列逐漸被淘汰, 而親和力較強的非高頻序列則逐漸進化為高頻序列, 從而導致含有大量核酸適配體的、有較高親和力的高頻序列在篩選競爭中逐漸占據優勢, 其種類數量不斷增加, 競爭優勢也越來越強,并呈現出強者越強的定向進化趨勢。這些變化趨勢與核酸適配體的篩選進化特點也是一致的。

3.3 相對重要性指數的意義

核酸適配體的篩選中, 需要對篩選產物進行測序, 然后再從眾多的測序序列中挑選出可能的候選序列, 進行親和特異性的驗證, 只有親和特異性較好的序列才被認定為核酸適配體。如何挑選候選序列,目前尚沒有統一的方法。早期是采用隨機的方法從測序結果中隨機挑選數個序列去進行親和特異性驗證(李衛濱等, 2007), 或者利用同源性將測序的序列分為數個家族, 然后從數個家族中隨機挑選序列進行驗證(Liuet al, 2014), 這種隨機方法要進行大量的驗證實驗, 工作量較大。另外, 也有學者是按照序列的頻率和占比來挑選候選序列的(Zhuet al, 2021)。而我們則綜合考慮了頻率和出現輪數, 采用相對重要性指數(IRI)這個指標, 按照IRI 從大到小來挑選相應的候選序列, 挑選出的7 個高頻序列后續均被確認是親和特異性較好的核酸適配體, 說明采用IRI 指數來挑選候選核酸適配體序列也是可行的。

3.4 高頻序列的二級結構

核酸適配體通常借助其結構中的堆積作用、形狀互補、靜電作用和氫鍵等相互作用, 通過頸環等結構,與靶目標或嵌合或包被, 從而形成穩定的復合物(Ellingtonet al, 1992; Hermannet al, 2000)。有研究發現, 適配體二級結構中的不同環區可以聯合參與靶標的識別與結合(孫羽菡, 2019), 還有研究表明, 兩個頸環的適配體要比一個頸環的適配體有更好的穩定性和親和特異性(Bhardwajet al, 2019), 我們的研究也表明二級結構中環的數量和大小很可能是與其結合位點相匹配的, 進而影響其與靶目標的結合和穩定性。

4 結論

利用Margalef 豐富度指數、Shannon-Wiener 多樣性指數、Pielou 均勻度指數和Berger-Parker 優勢度指數, 研究了鰻弧菌核酸適配體的篩選中高頻序列的進化多樣性, 結果發現, 在篩選進化過程中, 高頻序列的豐富度、多樣性、優勢度和均勻度都在前6 輪出現劇烈的波動, 之后隨著篩選的進行, 豐富度、多樣性和優勢度都呈逐步增加趨勢, 而均勻度則呈現逐步降低的趨勢, 說明含有大量核酸適配體的高頻序列在篩選競爭中逐漸占據優勢, 其種類數量不斷增加, 競爭優勢越來越強; 利用相對重要性指數(IRI)可將高頻序列分成優勢序列、重要序列、常見序列、一般序列和少見序列, IRI 越高的序列在篩選進化中的競爭優勢也越大, 也越有可能成為高親和力適配體;最后, 對IRI 較高的7 個核酸適配體的二級結構進行了模擬分析, 發現它們的二級結構中都存在有大、中、小三種不同尺寸的環狀結構, 推測這些環的大小和數量是與靶目標上的結合位點相對應的。