真核模式生物核小體分布參數(shù)模型分析

豐繼華, 郭亞茹, 牟 錦, 黃月月， 劉 珂， 范力棟

(云南民族大學(xué) 電氣信息工程學(xué)院, 昆明 650504)

在生物遺傳和進(jìn)化過程中，除了以DNA編碼為代表的“硬”遺傳物質(zhì)外，還存在著另一種重要的“軟”遺傳物質(zhì)，即基因組中的功能性蛋白質(zhì)及其化學(xué)修飾。以染色體為例，其在微觀上組織得井然有序，使得DNA序列在空間結(jié)構(gòu)上非常緊湊，但又互不干擾，尤其是一些重要的編碼片段會(huì)受到特殊“保護(hù)”?；蚪M實(shí)現(xiàn)這些功能主要依賴于一種染色質(zhì)基本亞基——核小體(Nucleosome)[1]。

核小體在基因組中擔(dān)負(fù)著兩個(gè)重要作用[2]：1)壓縮和折疊DNA，以適應(yīng)細(xì)胞核空間的大??；2)限制DNA的易接近性。細(xì)胞廣泛利用后一種功能來調(diào)控代謝過程，其中一個(gè)關(guān)鍵步驟是通過核小體定位對(duì)基因表達(dá)進(jìn)行精細(xì)的調(diào)節(jié)[3-4]。核小體沿DNA的定位控制了遺傳信息傳遞和獲取，僅相隔幾個(gè)核苷酸的替代位置就可能對(duì)基因表達(dá)產(chǎn)生深遠(yuǎn)的影響。因此，核小體作為染色體的基本單元，被認(rèn)為是影響表觀遺傳狀態(tài)的主要驅(qū)動(dòng)因素。然而，單細(xì)胞生物和多細(xì)胞生物核小體在結(jié)構(gòu)上的差異性仍然不明確[5-6]，影響了對(duì)基因組的認(rèn)識(shí)進(jìn)程。因此，如何定量獲取核小體空間尺度特征尤為迫切。

研究發(fā)現(xiàn)，真核生物基因轉(zhuǎn)錄起始區(qū)域的核小體分布具有高度保守性[7-8]。無論是單細(xì)胞的酵母，還是屬于高級(jí)哺乳動(dòng)物的人類，其核小體在基因轉(zhuǎn)錄區(qū)和編碼區(qū)的定位圖譜總體上都呈現(xiàn)出一種周期性振蕩衰減趨勢(shì)[9-10]。但深入觀察后會(huì)發(fā)現(xiàn)：不同物種的核小體占位圖譜在細(xì)節(jié)上是有明顯差別的，這種差異性可能代表了物種在染色質(zhì)結(jié)構(gòu)和功能上的進(jìn)化印跡[11]。在此提出構(gòu)建一個(gè)高精度數(shù)學(xué)模型對(duì)酵母和果蠅核小體分布進(jìn)行空間尺度和頻域分析。

1 核小體定位模型構(gòu)建

1.1 酵母與果蠅核小體定位圖譜

1.1.1 數(shù)據(jù)來源

本文使用的生物實(shí)驗(yàn)數(shù)據(jù)主要來源于兩個(gè)研究團(tuán)隊(duì)。酵母數(shù)據(jù)來源于Lee等[8]于2007發(fā)布的高分辨率酵母核小體定位率實(shí)驗(yàn)數(shù)據(jù)；果蠅數(shù)據(jù)則來自于Mavrich等[12]獲得的果蠅胚胎期核小體定位實(shí)驗(yàn)數(shù)據(jù)。

由于真核生物在基因轉(zhuǎn)錄起始區(qū)域周圍有著共同的性質(zhì)，本文選取的兩個(gè)實(shí)驗(yàn)數(shù)據(jù)其實(shí)驗(yàn)精度已達(dá)到我們的實(shí)驗(yàn)要求。

1.1.2 轉(zhuǎn)錄起始位點(diǎn)周圍核小體定位圖譜

由于基因轉(zhuǎn)錄起始位點(diǎn)(TSS)周圍核小體占位具有高度的保守性和代表性，我們選取了酵母和果蠅基因轉(zhuǎn)錄起始位點(diǎn)上、下游1000 bp范圍，對(duì)齊平均后的核小體定位圖譜(如圖1)，并基于定位圖譜分別建立擬合模型。

圖1 中，縱坐標(biāo)是以兩個(gè)物種基因轉(zhuǎn)錄起始位點(diǎn)(TSS)為中心，對(duì)齊平均并歸一化的核小體占位率曲線；橫坐標(biāo)是DNA的長度計(jì)量單位bp(核苷酸堿基對(duì)，Base Pair)。

1.2 擬合函數(shù)選取

根據(jù)實(shí)驗(yàn)核小體的定位特征，我們分別使用了多項(xiàng)式、傅里葉級(jí)數(shù)、高斯函數(shù)和正弦函數(shù)構(gòu)建定位模型。表1為以上4種函數(shù)對(duì)酵母和果蠅核小體定位曲線的擬合性能指標(biāo)，擬合結(jié)果如圖2所示。

圖1 兩種模式生物核小體定位圖譜

物種擬合函數(shù)和方差擬合優(yōu)度標(biāo)準(zhǔn)差自由度校正決定系數(shù)酵母多項(xiàng)式34.80370.58660.13221991 0.5847傅里葉6.37800.92420.05671983 0.9236高斯7.24460.91390.0605 1977 0.9129正弦函數(shù)1.81380.99170.01881974 0.9916 果蠅多項(xiàng)式27.19420.70060.11691991 0.6992傅里葉6.61830.92710.05781983 0.9265 高斯2.10640.97680.03261977 0.9765正弦函數(shù)1.08940.98800.02351974 0.9878

綜合比較4種函數(shù)的擬合性能和擬合效果后，發(fā)現(xiàn)由正弦復(fù)合函數(shù)構(gòu)成的模型擬合精度最高，如圖2-g、h所示，模型幾乎能在整個(gè)區(qū)域跟蹤核小體的定位特征；高斯函數(shù)次之；如圖2-f所示只是在遠(yuǎn)離TSS區(qū)域出現(xiàn)了誤差；而圖2-e中除在遠(yuǎn)離TSS區(qū)域出現(xiàn)了誤差外，TSS附近區(qū)域也存在較大誤差。傅里葉函數(shù)能基本擬合占位圖譜趨勢(shì)，但誤差較大(圖2-c、d)，其精度不能滿足后續(xù)分析需要；由多項(xiàng)式構(gòu)建的模型則擬合度最低，無法捕獲定位特征的細(xì)節(jié)特征(圖2-a、b)。根據(jù)以上結(jié)論，最終確定以正弦函數(shù)構(gòu)建核小體定位模型，并基于該模型進(jìn)行后續(xù)的尺度和頻域分析。

圖2 4種函數(shù)的擬合曲線

1.3 定位模型及參數(shù)

本文構(gòu)建的正弦復(fù)合函數(shù)如下：

(1)

(1)式由9個(gè)正弦函數(shù)線性組合而成，其中An為正弦函數(shù)的振幅，ωn為正弦函數(shù)的角速度，φn為正弦函數(shù)的初始相位。經(jīng)過對(duì)酵母和果蠅核小體占位率曲線擬合后，得到模型參數(shù)如表2所示。

表2 核小體定位模型參數(shù)

2 定位模型尺度計(jì)算

利用經(jīng)典信號(hào)處理方法挖掘基因組數(shù)據(jù)的內(nèi)在機(jī)理，在生物信息學(xué)領(lǐng)域得到廣泛應(yīng)用[13]。為了研究核小體定位的動(dòng)態(tài)與靜態(tài)特性，在此分別對(duì)兩種模式生物的定位模型(式1)進(jìn)行微積分計(jì)算。

2.1 模型微分特性

對(duì)式(1)求導(dǎo)后，得：

(2)

f′(x)表示核小體分布的動(dòng)態(tài)特性(占位變化)。令定位模型的導(dǎo)函數(shù)f′(x)=0，可得到核小體定位模型f(x)在TSS周圍的各極值點(diǎn)(表3)。其中，極大值表示單個(gè)核小體的分布中心，極小值表示連接DNA(核小體之間的DNA)分布中心(圖3)。

為便于分析，我們根據(jù)核小體分布中心位置與TSS距離的遠(yuǎn)近，對(duì)其依次進(jìn)行了編號(hào)(圖3和表3)。

圖3 兩種模式生物定位模型求導(dǎo)函數(shù)曲線

表3 核小體與連接DNA的分布中心位置(以TSS為坐標(biāo)原點(diǎn)，單位為bp)

從表3可以觀察到，盡管酵母和果蠅在TSS附近均形成一個(gè)最明顯的核小體缺失區(qū)域，但區(qū)別在于，果蠅核小體缺失區(qū)域范圍更大，寬度近400 bp，其特征更復(fù)雜，且存在兩個(gè)較弱的分布峰。與之對(duì)應(yīng)，酵母核小體缺失區(qū)域?qū)挾葍H300 bp左右，特征較為單一。

2.2 模型積分特性

為了確定兩個(gè)物種在TSS附近單核小體的定位穩(wěn)定性，我們以連接DNA分布中心作為積分區(qū)間，對(duì)式(1)計(jì)算定積分。

(3)

式(3)中，積分上、下限a、b取值為相鄰連接DNA分布中心坐標(biāo)，每個(gè)積分區(qū)間都包含有一個(gè)核小體分布中心。在此，我們計(jì)算了TSS上、下游共10個(gè)單核小體分布區(qū)間積分值(圖4)，An、ωn、φn取值為模型參數(shù)(表2)。

圖4是歸一化后的單核小體分布區(qū)間定積分值，一定程度上代表了相應(yīng)位置核小體的穩(wěn)定性(定位強(qiáng)度)。從計(jì)算結(jié)果可看到，酵母和果蠅核小體在TSS下游的定位強(qiáng)度總體要高于上游區(qū)域。但有趣的是，在+1和+2核小體分布區(qū)域，酵母的+2核小體定位強(qiáng)度要略高于+1核小體，而果蠅在總體上呈現(xiàn)出從+1到+5核小體定位強(qiáng)度依次遞減的規(guī)律。這一進(jìn)化上的差異，再次證明+1核小體在多細(xì)胞生物中承擔(dān)了更為關(guān)鍵的作用[14]。

圖4 轉(zhuǎn)錄起始位點(diǎn)周圍核小體的定位穩(wěn)定性

3 定位模型頻譜分析

3.1 模型傅里葉變換

根據(jù)定位模型式(1)，我們對(duì)其進(jìn)行了周期信號(hào)傅里葉變換，用于分析核小體定位模型的頻域特征。

(4)

3.2 模型的功率譜與相位譜

在式(4)基礎(chǔ)上，可得到模型的功率譜表達(dá)式。

(5)

圖5 兩種模式生物功率譜和相位譜

圖5橫坐標(biāo)為歸一化頻率，其中圖5-a觀察到兩種模式生物的核小體定位頻率(體現(xiàn)了占位變化速率)主要集中于0～0.007 Hz低頻段，表明核小體在基因轉(zhuǎn)錄區(qū)的位置較穩(wěn)定，進(jìn)一步證實(shí)了這一區(qū)域核小體組織的保守性和穩(wěn)健性，這一性質(zhì)對(duì)于細(xì)胞的生存具有重要意義[15]。觀察兩個(gè)物種頻譜后，發(fā)現(xiàn)單細(xì)胞酵母核小體占位變化只集中在頻點(diǎn)0～0.002 Hz附近，而多細(xì)胞生物果蠅的核小體占位變化頻點(diǎn)則在接近0.005～0.006 Hz附近出現(xiàn)一個(gè)明顯峰值，這一特征可能代表了二者在進(jìn)化上的差異，多細(xì)胞生物的核小體定位機(jī)制更為復(fù)雜。

圖5-b、c為兩種模式生物核小體定位模型的相位譜?？傮w上，相位譜隨著頻率升高呈現(xiàn)出周期性振蕩趨勢(shì)。在0～1.25×10-11Hz區(qū)間，酵母核小體的相位呈現(xiàn)出近似線性增長趨勢(shì)，而果蠅則呈現(xiàn)出復(fù)雜的非線性變化。兩個(gè)物種的相位譜在1.25×10-11Hz附近同時(shí)出現(xiàn)了一個(gè)拐點(diǎn)，但隨著頻率的增加，酵母核小體定位的相位向負(fù)方向急劇變化后，呈現(xiàn)出緩慢振蕩上行的趨勢(shì)；與之相反，果蠅核小體定位的相位卻在向正方向急劇變化后，呈現(xiàn)出了緩慢振蕩下行的趨勢(shì)。這一現(xiàn)象與圖3中的定位模型曲線是一致的，即以+1核小體位置為參照，在轉(zhuǎn)錄起始位點(diǎn)下游區(qū)域，果蠅核小體出現(xiàn)的位置要明顯滯后于酵母核小體出現(xiàn)位置，其中的生物意義還有待進(jìn)一步研究。

4 結(jié)論

本文研究了兩種真核模式生物核小體定位模型的構(gòu)建問題?；诙ㄎ荒Ｐ偷目臻g尺度和頻譜特征，比較了兩種生物在染色體結(jié)構(gòu)上的差異性，探索了利用數(shù)學(xué)模型對(duì)核小體定位過程進(jìn)行定量分析的方法，在此基礎(chǔ)上得出以下結(jié)論：1)利用數(shù)學(xué)模型對(duì)酵母和果蠅核小體定位曲線進(jìn)行擬合，為應(yīng)用經(jīng)典信號(hào)處理方法分析核小體定位特征提供了一種新思路；2)正弦作為一種單一頻率的基本信號(hào)與作為染色質(zhì)基本單元的核小體定位特征高度擬合，其背后可能隱藏著未知的自然規(guī)律;3)建立的核小體定位模型在一定程度上解釋了兩個(gè)物種在進(jìn)化過程中存在的差異性和共性(物種間的保守性)，為后續(xù)研究其他生物核小體定位及染色質(zhì)進(jìn)化問題做了鋪墊。

當(dāng)然，本文使用核小體定位模型和時(shí)頻分析方法所獲得的結(jié)論，只是針對(duì)核小體定位的理論分析，其實(shí)際機(jī)理非常復(fù)雜，其間還有諸多未解決的問題需要進(jìn)一步探索。