謝留威，李春義，巴恒星※

（1.中國(guó)農(nóng)業(yè)科學(xué)院特產(chǎn)研究所吉林省特種經(jīng)濟(jì)動(dòng)物分子生物學(xué)省部共建國(guó)家重點(diǎn)實(shí)驗(yàn)室，吉林 長(zhǎng)春130112；2.長(zhǎng)春科技學(xué)院，吉林 長(zhǎng)春130600）

通常情況下，每個(gè)哺乳動(dòng)物細(xì)胞含有大約100個(gè)線粒體，每個(gè)線粒體有2~10個(gè)線粒體DNA（mtDNA）拷貝[1-3]。在沒有組蛋白保護(hù)的情況下，mtDNA是高度易變的，其突變率是核DNA的6~17倍[4]。突變?cè)诩?xì)胞線粒體內(nèi)累積，與野生型共存，稱為mtDNA異質(zhì)性。目前，mtDNA異質(zhì)性已成為線粒體遺傳疾病的研究熱點(diǎn)[5]。許多研究認(rèn)為，發(fā)育、衰老和世代演進(jìn)的高頻率mtDNA異質(zhì)性是人類遺傳疾病的核心[6-7]。

對(duì)于已知的mtDNA異質(zhì)性突變，基于幾種PCR方法可以確定樣本中的突變位點(diǎn)及其相對(duì)頻率[8-11]。然而，對(duì)于發(fā)現(xiàn)未知突變，大規(guī)模的高通量測(cè)序技術(shù)更適合。但最近研究報(bào)道顯示，mtDNA異質(zhì)性檢測(cè)正從掃描已知的有限數(shù)量突變轉(zhuǎn)變?yōu)槿€粒體基因組篩查[12-13]。

ABI SOLiD測(cè)序平臺(tái)使用短讀長(zhǎng)（35~75 bp），通過一種顏色編碼框架（color space），測(cè)序不同堿基之間顏色信號(hào)。每一個(gè)測(cè)序堿基都依賴于2個(gè)連續(xù)的顏色編碼值，當(dāng)短讀長(zhǎng)與參考基因組比對(duì)后，如出現(xiàn)測(cè)序錯(cuò)誤（即單色差異）就可以與正確的堿基變化區(qū)別分開，這為重測(cè)序研究提供了新的非常具有優(yōu)勢(shì)的判定方法。這種獨(dú)特的編碼方式具有前后堿基校對(duì)功能，使堿基測(cè)序準(zhǔn)確率高達(dá)99.94%[14]。本研究通過對(duì)一只雄性梅花鹿血細(xì)胞全基因組SOLiD測(cè)序數(shù)據(jù)進(jìn)行分析，檢測(cè)mtDNA基因組中異質(zhì)性突變，以期為充分利用SOLiD全基因組測(cè)序數(shù)據(jù)鑒定未知的mtDNA異質(zhì)性點(diǎn)變異提供基本生物信息分析框架。

1 材料與方法

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

利用ABI SOLiD測(cè)序平臺(tái)對(duì)一只雄性梅花鹿血細(xì)胞進(jìn)行全基因組雙末端測(cè)序，讀段長(zhǎng)度為50 bp。另外，從GenBank中下載與全基因組測(cè)序梅花鹿屬于同一亞種的mtDNA全長(zhǎng)序列（登錄號(hào)：KX689229）作為SOLiD測(cè)序數(shù)據(jù)比對(duì)的線粒體參考序列。

1.2 分析方法

低頻變異檢測(cè)需要高質(zhì)量的測(cè)序數(shù)據(jù)，過濾掉包含未知堿基或平均質(zhì)量分?jǐn)?shù)低于20的讀段。使用PerM v0.4軟件[15]，選項(xiàng)設(shè)置為“-A-E-seed F4-V5-L400-e”，將高質(zhì)量的讀段與mtDNA參考序列進(jìn)行比對(duì)。在許多真核生物中，mtDNA經(jīng)常被轉(zhuǎn)移到核基因組中，產(chǎn)生了核線粒體假基因序列（nuclear mitochondrial DNA segments，NUMTs）。不同大小NUMTs片段與mtDNA序列存在高度差異[16-17]。根據(jù)NUMTs高變異率的特點(diǎn)，排除包含大于5個(gè)不匹配顏色值，相當(dāng)于2~3個(gè)不匹配堿基（相似度≥94%）的讀段。由于NUMTs平均長(zhǎng)度為240 bp[18]，SOLiD成對(duì)末端測(cè)序文庫(kù)的平均雙末端長(zhǎng)度為831 bp（圖1），這可以在一定程度上減少NUMTs的影響，同時(shí)去掉單端映射讀段。最后，利用Samtools v1.2[19]中的mpileup工具和Bcftools v1.3.1軟件[20]聯(lián)合檢測(cè)mtDNA異質(zhì)性點(diǎn)變異，在Tablet v1.18軟件中對(duì)點(diǎn)變異可視化顯示。

圖1雙末端長(zhǎng)度分布Fig.1 Insert size distribution

2 結(jié)果與分析

2.1 比對(duì)結(jié)果

過濾低質(zhì)量數(shù)據(jù)后，共獲得410 662 782（38 Gb）成對(duì)短讀段，相當(dāng)于全基因組測(cè)序深度約13（38 Gb/3 Gb）。大約0.035%的成對(duì)讀段比對(duì)到mtDNA參考序列，覆蓋率為99.62%，在2 265~2 329（65 bp）區(qū)域未被短讀段覆蓋（圖2）。線粒體基因組的平均測(cè)序深度約為892，一方面表明血細(xì)胞中mtDNA的拷貝數(shù)較高，另一方面表明部分源于NUMTs的讀段被比對(duì)。然而，動(dòng)物中NUMTs的總量約為核基因組的0.1%。在梅花鹿基因組中，NUMTs的比例（0.035%）相對(duì)較低，表明高度可變的NUMTs讀段被排除掉，這進(jìn)一步提高了鑒定mtDNA異質(zhì)性點(diǎn)突變的準(zhǔn)確性。

圖2 mtDNA上2 265~2 329（65 bp）區(qū)域未被短讀段覆蓋Fig.2 mtDNA region 2 265-2 329（65 bp）uncovered by short reads

2.2 mtDNA異質(zhì)性點(diǎn)變異鑒定

本研究共檢測(cè)到8個(gè)點(diǎn)突變（圖3），包括4種轉(zhuǎn)換和4種顛換（表1）。其中，5個(gè)點(diǎn)突變位于蛋白質(zhì)編碼基因（COX1、COX2、ND4和ND5），2個(gè)位于tRNASer，1個(gè)存在于12SrRNA中。為了驗(yàn)證NUMTs讀段是否影響點(diǎn)突變的鑒定，通過檢索GenBank中18個(gè)梅花鹿線粒體全基因組中相應(yīng)位點(diǎn)的變化，證實(shí)這8個(gè)點(diǎn)突變都存在于18個(gè)mtDNA基因組中（表1），表明它們是mtDNA特有的單核苷酸多態(tài)位點(diǎn)，而不是由NUMTs引入的。重要的是，C6180T、T7481A和A10909T這3個(gè)位點(diǎn)分別導(dǎo)致氨基酸的變化，即：絲氨酸（Ser）脯氨酸（Phe）→亮氨酸（Leu）→脯氨酸（Phe）、亮氨酸（Leu）→谷氨酰胺（Gln）。mtDNA異質(zhì)性點(diǎn)突變的選擇性壓力與組織特異性代謝率、細(xì)胞周期和生物能量需求有關(guān)[18]，這暗示梅花鹿不同類型的血細(xì)胞可能承受不同的選擇壓力。

圖3 8個(gè)mtDNA異質(zhì)性點(diǎn)突變可視化比對(duì)Fig.3 Visual comparison of 8 mt DNA heteroplasmic point variations

3 結(jié)論

目前，盡管SOLiD測(cè)序平臺(tái)已不常用，但其已產(chǎn)生大量原始基因組測(cè)序數(shù)據(jù)，這些數(shù)據(jù)仍然具有較高的科學(xué)價(jià)值，特別是SOLiD平臺(tái)對(duì)測(cè)序堿基進(jìn)行顏色編碼校對(duì)特性，在利用重測(cè)序進(jìn)行點(diǎn)變異檢測(cè)方面具有固有優(yōu)勢(shì)。本研究通過生物信息學(xué)分析方法充分挖掘一只雄性梅花鹿血細(xì)胞全基因組SOLiD測(cè)序中短讀長(zhǎng)數(shù)據(jù)，在mtDNA基因組中鑒定了8個(gè)異質(zhì)性點(diǎn)突變，將為利用全基因組測(cè)序數(shù)據(jù)鑒定未知的mtDNA異質(zhì)性點(diǎn)變異提供一個(gè)基本分析框架。