MutPrimerDesign：用于人類基因編碼區(qū)域突變位點(diǎn)的引物設(shè)計(jì)程序

曹英豪，彭公信

(中國(guó)醫(yī)學(xué)科學(xué)院 基礎(chǔ)醫(yī)學(xué)研究所 & 北京協(xié)和醫(yī)學(xué)院 基礎(chǔ)醫(yī)學(xué)院，北京 100730)

DNA突變不僅是導(dǎo)致生物多樣性和生物進(jìn)化的重要原因之一，同時(shí)也與人類的健康和疾病發(fā)生有著密切的關(guān)系。比如引起腫瘤發(fā)生最重要的原癌基因和抑癌基因，其大部分功能與細(xì)胞的生長(zhǎng)、分裂、錯(cuò)配修復(fù)和細(xì)胞凋亡等相關(guān)。如果這兩類基因在蛋白編碼區(qū)發(fā)生突變，則可能導(dǎo)致蛋白質(zhì)功能的改變或者喪失，造成細(xì)胞的非正常生長(zhǎng)，從而導(dǎo)致腫瘤的產(chǎn)生。因此，檢測(cè)基因突變并尋找可用于診斷和治療的突變位點(diǎn)已經(jīng)成為腫瘤研究中的常規(guī)內(nèi)容，具有重要的臨床和生物學(xué)意義[1]。

發(fā)生在基因編碼區(qū)的DNA突變，根據(jù)其對(duì)蛋白產(chǎn)物的影響可以分為移碼突變、非同義突變、同義突變等。其中前兩種突變由于引起了蛋白質(zhì)的序列變化，并可能影響其正常功能的發(fā)揮，因此，對(duì)這兩種突變的檢測(cè)和研究尤為重要。近年來，以高通量測(cè)序技術(shù)為主的檢測(cè)技術(shù)得到了迅猛的發(fā)展，比如目標(biāo)區(qū)域捕獲技術(shù)(Panel)、全外顯子組測(cè)序技術(shù)(WES)、全基因組測(cè)序技術(shù)(WGS)等[2]。尤其是Panel測(cè)序技術(shù)，通過特異性捕獲探針對(duì)多個(gè)基因目標(biāo)區(qū)域DNA片段進(jìn)行捕獲測(cè)序，具有成本低、測(cè)序深度高的優(yōu)勢(shì)，能夠準(zhǔn)確給出重要基因的突變信息，在臨床上得到了廣泛的應(yīng)用[3-5]。

為了應(yīng)用于疾病的診斷和治療，通過測(cè)序平臺(tái)獲取的或者文獻(xiàn)中已發(fā)表的突變信息，需要進(jìn)一步擴(kuò)大人群規(guī)模來驗(yàn)證該突變的可靠性。利用qPCR及衍生技術(shù)檢測(cè)突變位點(diǎn)，在降低成本的同時(shí)又保證其準(zhǔn)確性，是一種簡(jiǎn)單而可行的方法[6-8]。設(shè)計(jì)引物是qPCR的首要問題，目前已有許多成熟的引物設(shè)計(jì)軟件可以使用。這其中除了Primer3[9]、Oligo[10]等綜合性的引物設(shè)計(jì)軟件包之外，還有一些專用的引物設(shè)計(jì)軟件，比如整合了Primer3和BLAST的模板特異性引物設(shè)計(jì)工具Primer-BLAST[11]，突變引物設(shè)計(jì)工具WASP[12]，PrimerX(http://www.bioinformatics.org/primerx/)，MutScreener[13]，Edesign[14]以及適用于克隆表達(dá)載體的PrimerCE[15]等。但由于測(cè)序平臺(tái)的突變檢測(cè)報(bào)告以及文獻(xiàn)中給定的突變信息，大部分均是以HGNC(HUGO Gene Nomenclature Committee)[16]的基因名稱(Gene Symbol)命名的基因編碼區(qū)域的相關(guān)坐標(biāo)，因此在設(shè)計(jì)引物時(shí)，這些軟件或多或少有些不適用。主要有兩個(gè)問題：(1)人類基因普遍有多個(gè)轉(zhuǎn)錄本，需要通過不同的數(shù)據(jù)庫尋找對(duì)應(yīng)的標(biāo)準(zhǔn)轉(zhuǎn)錄本序列。(2)由于內(nèi)含子的影響，對(duì)于外顯子編碼區(qū)的坐標(biāo)，需要根據(jù)基因結(jié)構(gòu)進(jìn)行基因組坐標(biāo)的手工轉(zhuǎn)換。另外，為了能夠進(jìn)一步直接從擴(kuò)增子中檢測(cè)序列是否發(fā)生突變，針對(duì)突變位點(diǎn)可增加探針序列的設(shè)計(jì)。

利用Python語言編寫了MutPrimerDesign，通過對(duì)人類基因組序列數(shù)據(jù)庫以及基因注釋信息的解析，可識(shí)別多種數(shù)據(jù)庫的基因名稱編號(hào)，轉(zhuǎn)換基因編碼區(qū)坐標(biāo)為基因組坐標(biāo)。同時(shí)，通過對(duì)primer3-py(https://github.com/libnano/primer3-py)程序包的整合，可批量設(shè)計(jì)突變位點(diǎn)附近的引物及探針序列。軟件的下載地址為:https://github.com/bioinfo-ibms-pumc/MutPrimerDesign。

1 方 法

MutPrimerDesign的流程軟件主要分為三步，第一步接收輸入信息，并將基因名稱轉(zhuǎn)換為HGNC中標(biāo)準(zhǔn)基因名稱，同時(shí)獲取標(biāo)準(zhǔn)轉(zhuǎn)錄本以及相對(duì)應(yīng)的RefSeq[17]編號(hào)。在檢測(cè)輸入數(shù)據(jù)時(shí)，MutPrimerDesign可以接收三種不同的基因命名，包括基因名稱、Ensembl編號(hào)或者HGNC編號(hào)。第二步轉(zhuǎn)換編碼區(qū)坐標(biāo)為基因組坐標(biāo)，并提取突變位點(diǎn)附近基因組序列。在轉(zhuǎn)換坐標(biāo)前，通過對(duì)人類基因組的注釋文件進(jìn)行解析，構(gòu)建所有基因?qū)?yīng)的基因結(jié)構(gòu)特征，并根據(jù)該結(jié)構(gòu)，將輸入的標(biāo)準(zhǔn)轉(zhuǎn)錄本的編碼區(qū)坐標(biāo)，轉(zhuǎn)換為對(duì)應(yīng)基因組上的坐標(biāo)，并進(jìn)一步利用python的BioSeq程序包對(duì)人類基因組序列進(jìn)行解析，提取相應(yīng)的擴(kuò)增子序列(Amplicon)。第三步調(diào)用python的primer3-py程序包，進(jìn)行引物設(shè)計(jì)，輸出報(bào)告。為了方便在引物擴(kuò)增序列的同時(shí)能夠區(qū)分野生型和突變型位點(diǎn)，軟件中加入了探針序列的設(shè)計(jì)，并默認(rèn)設(shè)定探針序列的Tm值比引物序列的Tm值低5度。具體流程(見圖1)。

1.1 數(shù)據(jù)源下載及安裝

通過git命令可以克隆下載軟件全部代碼，命令如下：

git clone https://github.com/bioinfo-ibms-pumc/MutPrimerDesign.git

下載后的目錄中主要包括4個(gè)文件，分別為MutPrimerDesign.py，gffLoader.py，fastaLoc.py和namedb文件。主程序MutPrimerDesign.py用于接收基因突變信息并設(shè)計(jì)引物，gffLoader.py用于解析基因注釋文件，fastaLoc.py用于提取基因組序列，namedb文件為基因名稱數(shù)據(jù)庫。MutPrimerDesign依賴BioPython、Pandas、Primer3-Py數(shù)據(jù)包，這些python依賴包均可使用pip命令進(jìn)行標(biāo)準(zhǔn)安裝。依賴包的安裝命令如下：

pip3 install biopython pandas primer3-py

另外，由于HGNC中基因名稱對(duì)應(yīng)于RefSeq的標(biāo)準(zhǔn)序列，因此，為了確保序列坐標(biāo)的準(zhǔn)確性，MutPrimerDesign目前只接收由NCBI網(wǎng)址提供的人類基因組序列文件以及相對(duì)應(yīng)版本的基因注釋文件。人類GRCh38版本基因組及基因注釋文件下載及解壓縮命令如下：

wget-c ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.39_GRCh38.p13/GCF_000001405.39_GRCh38.p13_genomic.fna.gz

wget-c ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.39_GRCh38.p13/GCF_000001405.39_GRCh38.p13_genomic.gff.gz

gunzip-c GCF_000001405.39_GRCh38.p13_genomic.fna.gz > genome.fa

gunzip-c GCF_000001405.39_GRCh38.p13_genomic.gff.gz > gene.gff

圖1 MutPrimerDesign引物設(shè)計(jì)流程圖及原理Fig.1 Flow chart of MutPrimerDesign and its principle

1.2 命令行說明

MutPrimerDesign主要由兩部分命令組成，一部分為程序主要參數(shù)，包括輸入輸出以及各種依賴數(shù)據(jù)庫的定義；另一部分為引物參數(shù)，主要包括引物的Tm值，長(zhǎng)度，GC含量等(見圖2)。

圖2 MutPrimerDesign的命令行說明截圖Fig.2 Screen shot for commands of MutPrimerDesign

MutPrimerDesign接收的輸入?yún)?shù)有兩種類型。第一種，是“-i”命令，直接接收基因突變信息。每個(gè)基因突變的信息分為三個(gè)部分，以冒號(hào)分開，內(nèi)容依次分別為基因名稱、突變?cè)谠摶虻腃DS區(qū)域坐標(biāo)、突變長(zhǎng)度。不同突變信息以逗號(hào)隔開，比如“KRAS:24:2,MET:1124:1,EGFR:2573:1”。另一種是“-f”命令，直接接收文本文件。文件中每一行為一個(gè)基因突變信息，內(nèi)容與-i中突變信息格式一致。

MutPrimerDesign除了輸入的基因信息外，還需要三個(gè)輸入文件，分別為基因名稱轉(zhuǎn)換文件namedb，基因組序列文件genome.fa以及相對(duì)應(yīng)版本的基因注釋文件gene.gff。其中namedb為MutPrimerDesign預(yù)處理后的文件，另外兩個(gè)文件可從NCBI網(wǎng)站上下載解壓后得到。

MutPrimerDesign的輸出分為三個(gè)文本文件，通過“-o”參數(shù)來控制。主要包括擴(kuò)增子序列文件，另外還有一個(gè)以“.primer”后綴結(jié)尾的引物序列文件，以及一個(gè)以“.detail”后綴結(jié)尾的詳細(xì)報(bào)告文件。

2 結(jié) 果

2.1 不同數(shù)據(jù)庫來源的引物設(shè)計(jì)舉例

經(jīng)統(tǒng)計(jì)，50%的黑色素瘤中原癌基因BRAF均發(fā)生了突變，其中90%的突變發(fā)生在第600位纈氨酸的位置上[18]。該位點(diǎn)是一個(gè)非常重要的免疫檢查點(diǎn)，位點(diǎn)的突變會(huì)影響腫瘤用藥效果。對(duì)此，分別使用BRAF基因的基因名稱、Ensembl ID號(hào)、HGNC ID號(hào)作為輸入，對(duì)其V600突變位點(diǎn)進(jìn)行引物設(shè)計(jì)，運(yùn)行命令如下所示：

(1)python3 MutPrimerDesign.py -g genome.fa -G gene.gff -k namedb -o result.fa -i “BRAF:1799:1”

(2)python3 MutPrimerDesign.py -g genome.fa -G gene.gff -k namedb -o result.fa -i “ENSG00000157764:1799:1” -s ensembl_id

(3)python3 MutPrimerDesign.py -g genome.fa -G gene.gff -k namedb -o result.fa -i “1097:1799:1” -s hgnc_id

BRAF基因編碼區(qū)第1 799位引物設(shè)計(jì)的部分結(jié)果輸出截圖(見圖3)。輸出內(nèi)容中，每一個(gè)擴(kuò)增子信息以Amplicon開始，后接設(shè)計(jì)順序及引物編號(hào)，第一行為擴(kuò)增子編號(hào)，第二行包括基因突變位點(diǎn)信息，依次為輸入突變信息、基因、Refseq編號(hào)、染色體編號(hào)、擴(kuò)增子在基因組上坐標(biāo)區(qū)間、突變位點(diǎn)在基因組上坐標(biāo)、突變位點(diǎn)在編碼區(qū)坐標(biāo)、參考?jí)A基及長(zhǎng)度和引物編號(hào)。第三行為擴(kuò)增子長(zhǎng)度。第四行到第七行為primer3-py設(shè)計(jì)的引物及探針序列信息。其中，第四行為根據(jù)熱力學(xué)方法計(jì)算引物對(duì)整體形成二級(jí)結(jié)構(gòu)的可能性以及3’末端形成二級(jí)結(jié)構(gòu)的可能性。第五行為5’端引物，主要信息包括：引物序列，引物位置，Tm值，GC含量，引物自身形成二級(jí)結(jié)構(gòu)的可能性，3’末端形成二級(jí)結(jié)構(gòu)的可能性，形成發(fā)夾結(jié)構(gòu)的可能性，3’末端的穩(wěn)定性(delta G, kcal/mol)。同樣，第六行和第七行分別為3’端引物和探針序列信息。第八行為擴(kuò)增子序列，第九行為引物、探針、突變位點(diǎn)在擴(kuò)增子中的位置。由圖3中可以看出，無論是采用基因名稱，Ensembl編號(hào)或者HGNC編號(hào)，軟件均可以完成相同的引物設(shè)計(jì)。

2.2 批量引物設(shè)計(jì)舉例

Maki-Nevala S等通過對(duì)425例非小細(xì)胞型肺癌病人樣本進(jìn)行測(cè)序，找到了9個(gè)復(fù)發(fā)突變熱點(diǎn)[19]。利用MutPrimerDesign的命令行輸入“-i”和文件輸入“-f”兩種方式進(jìn)行引物設(shè)計(jì)，引物及探針設(shè)計(jì)結(jié)果(見表1)。運(yùn)行命令如下所示：

(1)python3 MutPrimerDesign.py -g genome.fa -G gene.gff -k namedb -o lung.fa -i "KRAS:34:2,MET:1124:1,MET:3029:1,EGFR:2573:1,TP53:820:1,PIK3CA:1624:1,EGFR:2235:1,TP53:469:1"

(2)python3 MutPrimerDesign.py -g genome.fa -G gene.gff -k namedb -o lung.fa -f samples.file

2.3 引物重新設(shè)計(jì)舉例

對(duì)于已經(jīng)完成引物設(shè)計(jì)的基因位點(diǎn)，并且擴(kuò)增子序列已經(jīng)保存在result.fa文件中時(shí)，如果需要修改引物參數(shù)進(jìn)行重新設(shè)計(jì)，則可以使用-p參數(shù)忽略序列提取過程，僅進(jìn)行引物重新設(shè)計(jì)，運(yùn)行命令如下所示：

(1)python3 MutPrimerDesign.py -g genome.fa -G gene.gff -k namedb -o result.fa -p --primerNum 3 --primer_opt_tm 61 --primer_max_tm 64 --primer_opt_size 23

2.4 軟件運(yùn)行效率測(cè)試

MutPrimerDesign在設(shè)計(jì)引物時(shí)消耗內(nèi)存及計(jì)算資源較小，普通的臺(tái)式機(jī)或者筆記本環(huán)境下均可以運(yùn)行。為了測(cè)試軟件的運(yùn)行性能和穩(wěn)定性，以Linux操作系統(tǒng)為例(系統(tǒng)版本：Linux Mint 19.3 Cinnamon；CPU：Intel Core i7-9700；內(nèi)存：24 GB；存儲(chǔ)：1T機(jī)械硬盤)，采用MutPrimerDesign軟件默認(rèn)參數(shù)進(jìn)行了性能測(cè)試：系統(tǒng)開機(jī)后，首次完成8個(gè)突變位點(diǎn)的引物設(shè)計(jì)耗時(shí)為1 m 33 s。而后當(dāng)系統(tǒng)自動(dòng)加載了基因組信息緩存后，無論是多基因還是單基因的100個(gè)突變位點(diǎn)，設(shè)計(jì)引物只需要30 s即可完成。

3 結(jié) 論

通過Python，開發(fā)了MutPrimerDesign引物設(shè)計(jì)軟件，可專門用于人類基因編碼區(qū)的突變位點(diǎn)引物設(shè)計(jì)。通過對(duì)基因名稱進(jìn)行解析，MutPrimerDesign將自動(dòng)選擇HGNC給定的標(biāo)準(zhǔn)轉(zhuǎn)錄本信息，完成基因編碼區(qū)坐標(biāo)到基因組坐標(biāo)的轉(zhuǎn)換以及基因組序列的提取，從而能夠批量、自動(dòng)化地完成基因編碼區(qū)突變位點(diǎn)的引物及探針序列的設(shè)計(jì)。目前軟件的不足之處是輸出格式較為簡(jiǎn)單，輸出結(jié)果可以更加豐富全面。

4 討 論

以Primer3, Oligo為首的引物設(shè)計(jì)軟件的出現(xiàn)，使得引物設(shè)計(jì)相關(guān)工作變得越來越簡(jiǎn)單。但是，伴隨著生物技術(shù)的不斷發(fā)展，生物學(xué)問題越來越多，引物設(shè)計(jì)的需求也越來越細(xì)化。比如高通量測(cè)序平臺(tái)的不斷發(fā)展，突變位點(diǎn)檢測(cè)等相關(guān)的下游應(yīng)用將會(huì)變得越來越多，對(duì)此，WASP，MutScreener等特定的引物設(shè)計(jì)軟件一定程度下緩解了這些需求。但針對(duì)突變位點(diǎn)設(shè)計(jì)引物時(shí)，用戶輸入基因的名稱可能來自不同數(shù)據(jù)庫，如何自動(dòng)地轉(zhuǎn)換基因名稱，對(duì)含有多個(gè)轉(zhuǎn)錄本的基因進(jìn)行標(biāo)準(zhǔn)轉(zhuǎn)錄本的選擇，以及如何根據(jù)基因的編碼區(qū)位置自動(dòng)定位到基因組上相應(yīng)的位置，這些看似簡(jiǎn)單的問題卻沒有完善及考慮。MutPrimerDesign通過對(duì)不同數(shù)據(jù)庫的整合以及基因結(jié)構(gòu)的解析解決了這些特定的問題。另外，引物設(shè)計(jì)開發(fā)的一個(gè)重要指標(biāo)是如何能夠更合理地減少用戶的工作量，輔助用戶完成批量、自動(dòng)化的引物設(shè)計(jì)。MutPrimerDesign試圖通過最簡(jiǎn)潔的輸入來完成多種情況下突變位點(diǎn)的引物設(shè)計(jì)，也是對(duì)該指標(biāo)的一次嘗試。在將來，MutPrimerDesign將會(huì)考慮整合HGVS(Human Genome Variation Society)[20]的突變標(biāo)準(zhǔn)規(guī)則，進(jìn)一步直接對(duì)接測(cè)序平臺(tái)的突變檢測(cè)報(bào)告，豐富并完善突變序列與參考序列的對(duì)比信息，更加方便、合理地完成突變位點(diǎn)檢測(cè)的引物設(shè)計(jì)工作。另外，MutPrimerDesign將會(huì)考慮應(yīng)用到更多的物種，實(shí)現(xiàn)多物種的突變引物設(shè)計(jì)。