海南黑山羊GH1基因nsSNPs功能性預測及生物信息學分析

管凇 周漢林 徐鐵山 孫衛(wèi)平 胡海超

（中國熱帶農業(yè)科學院熱帶作物品種資源研究所 571101）

生物信息學技術是現階段篩選候選功能位點最有效的技術之一，因此實驗研究中可利用該技術預測未知表型的nsSNPs的功能[1]。已知，現階段研究人員利用各類生物信息分析工具對人類或動物表型和性狀基因功能性突變進行了大量的分析和預測，且取得了突破性進展[2-6]。本試驗主要通過利用生物信息學技術預測海南黑山羊GH1 基因編碼區(qū)nsSNPs 的功能。旨在篩選可顯著影響海南黑山羊生長的功能性突變位點，并分析該突變位點的作用機理，以期更好的理解黑山羊GH1 基因的作用機理，同時也對應用標記輔助選擇改善山羊生長發(fā)育提供一定的參考。

1 材料和方法

1.1 樣品采集和DNA 處理

海南黑山羊樣品采自品資所山羊保種場山羊耳組織。對104 只成年山羊的體長、體高、體重進行測量。剪取耳組織，置于裝有75%乙醇的離心管中，4℃保存。采用DNA 提取試劑盒提取耳組織總DNA，用紫外分光光度計檢驗每個樣品DNA質量濃度。

1.2 引物設計

從NCBI 數據庫中查詢山羊GH1 和IGF2 基因DNA 序列（GenBank 登錄號分別為：NC_030826 和NC_030836），根據單核苷酸多態(tài)性位點的序列信息，利用引物設計軟件Primer Plex 2 設計引物，使用primer3plus 在線設計單堿基延伸引物，引物信息見表1。

表1 引物信息

1.3 SNP 位點的篩選和分型

對提取的DNA 基因組序列進行擴增，PCR 擴增反應總體積為10μL：5μL5×TaqPCRmastermix，上下游引物各0.5μL（20μmol/L），模板DNA1μL (20ng/ul），加ddH2O 至10μL。PCR 擴增程序：94℃預變性5min；45 個循環(huán)（94℃，20s；60℃，30s；72℃，30min）；72℃再延伸3min。將PCR 產物去除體系中游離的dNTPs，再進行單堿基延伸反應，反應體系總體積5μl。PCR 產物委托武漢天一輝遠生物科技有限公司采用SnaPshot 方法，用Gene marker 分析得到的.fsa 結果，導出峰圖及表格文件，并統(tǒng)計出每個樣品的SNP 位點的基因型。

1.4 數據分析

用 ROVEN（Protein Variation Effect Analyzer）http://provean.jcvi.org/index.php 和SITF（Sorting Intolerant From Tolerant）https://sift.bii.a-star.edu.sg/和分析nsSNPs 對蛋白功能的影響，利用http://rna.tbi.univie.ac.at/cgi-bin/RNAfold.cgi 軟件進行mRNA 二級結構預測[7]。利用Sopma（https://npsa-prabi.ib -cp.fr/cgi -bin/npsa_automat.pl?page=npsa_sop -ma.htm 和CPHmodels（http://www.cbs.dtu.dk/services/CPHmodels/）軟件分別對突變體的二級結構和三級結構分別進行預測和分析，進而評估nsSNPs 可能造成的影響，三級結構3D 結果展示使用軟件BIOVIA（https://www.3dsbiovia.com）。

2 結果與分析

2.1 單核苷酸多態(tài)性位點分析

對海南黑山羊GH1 基因的測序比對結果顯示，在GH1 基因外顯子1、2、3、4 區(qū)域檢測到10 個SNP 位點，有4 個錯義突變位點，4 個錯義突變突變位點測序峰圖見圖1 和表2。在第一外顯子350bp 處發(fā)生G→A 錯義突變，突變位點導致氨基酸由原來的甘氨酸（Gly）變?yōu)榻z氨酸（Ser），第二外顯子739bp 處發(fā)生G→A 錯義突變，突變位點導致氨基酸由原來甘氨酸（Gly）變?yōu)榻z氨酸（Ser），第三外顯子1101bp 處發(fā)生A→G 錯義突變，突變位點導致氨基酸由原來天冬氨酸（Asp）變?yōu)楦拾彼幔℅ly），第三外顯子1154bp 處發(fā)生C→T 錯義突變，突變位點導致氨基酸由原來精氨酸（Arg）變?yōu)樯彼幔═rp），其余6 個同義突變位點沒有導致氨基酸改變。GH1 基因Exon2-729G→A、Exon2-739G→A 兩個位點在群體中分布的等位基因型頻率突變前后相同，兩個位點連鎖，Exon3-1101A→G、xon3-1120G→A、Exon3-1154C→T3 個位個等位基因型頻率突變前后相同，完全連鎖。

圖1 錯義突變位點測序峰圖

表2 突變位點信息

2.2 不同SNP 位點與性狀關聯分析

利用SPSS 軟件統(tǒng)計各基因的不同SNP 對應的表型（體長，體高和體重）的平均值及標準誤差，詳細結果見下表4。利用SPSS 中的一般線性模型對各基因的SNPs 位點與表型進行關聯分析，GH1 基因外顯子1、2、3、4 檢測到10 個SNP 位點，每個位點均只有2 種基因型。一般線性模型的主體間效應檢測結果顯示，所有SNP 位點基因型與表型（體長，體高和體重）均未有顯著差異。

表4 突變位點與生長性狀的關聯性分析

2.3 SNP 位點的生物信息學分析

2.3.1 nsSNPs 對蛋白功能的影響

2.3.1.1 ROVEN 分析

用ROVEN 分析nsSNPs 對蛋白功能的影響，當SNP 的PROVEAN score 分數等于或低于-2.5 的變體被認為是“有害的”。當SNP 的PROVEAN score 分數高于-2.5 的變體被視為“中性”突變。分析結果顯示，D129G（天冬氨酸變?yōu)楦拾彼幔┖蚏147W（精氨酸變?yōu)樯彼幔㏒NP 的PROVEAN score 分數低于-2.5，該突變位點屬于有害突變，對蛋白功能有較大影響見表5。

表5 nsSNPs 的score 評判分值

2.3.1.2 SITF 分析

用SITF 分析nsSNPs 對蛋白功能的影響，評級為4 個標準：tolerated、low confidence tolerated、low confidence deleterious、deleterious。當SNP 的score 大于0.5 時，說明該SNP 位點對蛋白質功能沒有影響或者影響較小；當SNP 的score 小于0.5 時，說明這個突變的SNP 位點對蛋白質功能有較大影響。一個SNP 的得分越低，危害性越大，分析結果顯示，G85S（甘氨酸變?yōu)榻z氨酸）、D129G（天冬氨酸變?yōu)楦拾彼幔?、R147W（精氨酸變?yōu)樯彼幔┩蛔兾稽cscore＜0.5，該突變位點屬于有害突變，對蛋白功能有較大影響，見表6。

表6 nsSNPs 的score 評判分值

2.3.2 mRNA 二級結構分析

對引起錯義突變的GH1 基因Exon1-350G→A、Exon2-739G →A；Exon3-1101A →G、Exon3-1154C →T位點進行mRNA 二級結構預測，結果顯示，4 個突變位點突變前后均引起mRNA 二級結構和自由能發(fā)生變化，Exon1-350G→A 位點的自由能和結構變化較明顯。說明基因序列SNPs 突變前和突變后mRNA 自由能的改變、二級結構的改變均導致其結構的穩(wěn)定性發(fā)生了改變，見圖2 和表7。

表7 mRNA 二級結構突變前后自由能變化

圖2 mRNA 二級結構突變前后的變化結果

2.3.3 蛋白質二級結構分析

對引起錯義突變的GH1 基因Exon1-350G→A、Exon2-739G→A；Exon3-1120G→A、Exon3-1154C→T 位點進行進行蛋白質二級結構分析，結果表明：Exon1-350G→A（G31S）位點突變前后的無規(guī)則卷曲、α-螺旋、β-轉角、擴展鏈均無變化，該突變位點對蛋白質二級結構沒有發(fā)生影響；Exon2-739G→A（G85S）位點突變前后的無規(guī)則卷曲和擴展鏈發(fā)生改變，Exon3 -1120G →A（D129G）、Exon3 -1154C →T（R147W）位點突變前后的α-螺旋和無規(guī)則卷曲發(fā)生改變，這3 個突變位點對蛋白質二級結構發(fā)生影響。統(tǒng)計結果見表8。

表8 蛋白質突變前后二級結構預測分析

2.3.4 位點突變前后蛋白質三級結構分析

突變位點Exon1-350G →A（G31S）、Exon2-739G →A（G85S）、Exon3-1120G →A（D129G）、Exon3-1154C →T（R147W）突變前后蛋白質三級結構進行預測分析，結果三級結構與二級結構變化一致，Exon1-350G→A（G31S）位點突變前后三級結構無變化，該突變位點對蛋白質三級結構沒有發(fā)生影響；Exon2-739G→A（G85S），Exon3-1120G→A（D129G）、Exon3-1154C→T（R147W）位點突變前后三級結構發(fā)生變化（見圖3），預測這3 個突變位點對蛋白質的結構功能可能會造成影響。

圖3 突變前后蛋白質三級結構

3 討論

隨著生物信息學的不斷發(fā)展，生物信息學工具已廣泛應用于人類醫(yī)學和農業(yè)領域中對功能性SNP 研究中。該工具可以顯著降低SNPs 的數目，且能篩選出最有可能的功能性nsSNPs，為后續(xù)功能試驗驗證提供一定的參考依據。因此，使用該工具可顯著降低研究成本和工作量[7]。目前較常用的工具為SIFT和ROVEN，建議聯合使用，以擴大數據庫，提高檢測的靈敏度和準確度。單獨使用某一工具時預測結果不同，SIFT 預測結果假陽性率更低[8-9]。ROVEN 對高風險位點的預測結果更準確，但對中風險位點的預測準確度較低[10]，本研究同時采用SIFT 和ROVEN 兩種方法對海南黑山羊GH1 基因的外顯子區(qū)的功能性nsSNPs 進行了預測（表5-6）。預測G85S、D129G、R147W3 個位點可能是有害突變位點。D129G、R147W 兩個位點突變前后二級結構的α-螺旋比列下降，無規(guī)則卷曲比例有上升，G85S 位點的擴展鏈發(fā)生變化。通過nsSNPs、mRNA二級結構、蛋白質二級結構和蛋白質三級結構分析結果，預測G85S、D129G、R147W3 個突變位點可能是導致海南黑山羊個體矮小，生長緩慢的功能性SNP。G85S、D129G、R147W3 個突變位點分別使氨基酸由甘氨酸（Gly）變?yōu)榻z氨酸（Ser），天冬氨酸（Asp）變甘氨酸（Gly），精氨酸（Arg）變?yōu)樯彼幔═rp）。甘氨酸是在氨基酸系列中結構最為簡單的非極性氨基酸，在中樞神經系統(tǒng)，尤其是在脊椎里，甘氨酸是抑制性神經遞質；絲氨酸是一種非必需氨基酸，在脂肪和脂肪酸的新陳代謝及肌肉的生長中發(fā)揮著作用；天冬氨酸是一種α-氨基酸，對肝臟功能有增強作用，消除疲勞；精氨酸為堿性氨基酸，是維持嬰幼兒生長發(fā)育必不可少的氨基酸；色氨酸對泌乳期的乳牛和母豬有促進泌乳作用，當畜禽缺乏色氨酸時，生長停滯，體重下降，脂肪積累降低。本文預測到G85S、D129G、R147W3 個突變位點屬于有害突變，對蛋白功能有較大影響，可能是導致海南黑山羊個體矮小，生長緩慢的功能性SNP。