郭思武，劉玉芳※，儲明星

（1.中國農(nóng)業(yè)科學院北京畜牧獸醫(yī)研究所，農(nóng)業(yè)農(nóng)村部動物遺傳育種與繁殖重點實驗室，北京 100193；2.河北工程大學生命科學與食品工程學院，河北 邯鄲 056001）

胴體屠宰率是綿羊生產(chǎn)中的重要經(jīng)濟性狀。胸椎數(shù)和胴體長是提高胴體屠宰率的重要指標，但其遺傳機制尚不明晰。因此，挖掘相關候選基因多態(tài)性位點對提高綿羊胴體產(chǎn)量的分子標記輔助選育具有重要意義。HERC 蛋白家族由4 個蛋白組成，分為2 個亞家族（大型HERC 和小型HERC），其中大型HERCs（包括HERC1 和HERC2 兩種蛋白）是長度近5 000 個氨基酸殘基的巨型蛋白質，人中位于15 號染色體，綿羊中位于7 號染色體[1]。而小型HERC（包括HERC3 和HERC4 兩種蛋白）的大小不到大型HERCs 的1/4[2]。HERC1（HECT and RLD domain containing E3 ubiquitin protein ligase family member 1）是該家族中最大的成員（532 kDa），包含許多保守的序列特征，這些特征被認為在蛋白質的整體功能中起著多種不同的作用[3]。

醛脫氫酶（ALDHs）是一個具有類似主要結構的NAD(P)+-依賴性酶的超家族，催化氧化廣泛的內(nèi)源性和外源性脂肪族和芳香族醛類[4]。到目前為止，在人類基因組中已經(jīng)發(fā)現(xiàn)了16 個ALDH 家族基因，分布在染色體的不同位置，其中，ALDH3A2、ALDH4A1、ALDH5A1 和ALDH6A1 的多態(tài)性與通常以神經(jīng)系統(tǒng)并發(fā)癥為特征的代謝性疾病相關[4]。ALDH 基因超家族包括原核生物和真核生物物種中的250 多個ALDH基因[5]。迄今為止，真核生物ALDH 基因超家族至少由18 個家族組成，包括37 個亞家族[6]。ALDH6A1 參與纈氨酸和嘧啶的降解途徑，分別催化丙二酸和丙二酸甲酯不可逆氧化脫羧為乙酰和丙酰輔酶A[7]。在大鼠的研究中發(fā)現(xiàn)，ALDH6A1 在肝臟、腎臟、心臟、肌肉和大腦中均有表達[7]。ALDH6A1 缺陷與發(fā)育遲緩有關，但目前尚不清楚這種缺陷的分子基礎[8]。

迄今為止，未有HERC1 和ALDH6A1 突變位點與家畜胸椎數(shù)性狀相關的文獻報道，相關功能多集中在骨髓增生腫瘤[9]、乳腺癌[10]和腦部生長發(fā)育[11,12]等疾病研究領域。因此，本研究針對綿羊胸椎數(shù)變異相關候選基因HERC1 和ALDH6A1 SNPs 位點與胸椎數(shù)等性狀進行關聯(lián)分析，并通過對比突變前后蛋白理化性質和空間結構的差異，探討HERC1 和ALDH6A1基因多態(tài)位點與胸椎數(shù)性狀的關聯(lián)性，為進一步拓展綿羊胸椎數(shù)變異的研究和分子標記輔助育種提供理論基礎。

1 材料與方法

1.1 樣品采集

實驗選取天津市畜牧獸醫(yī)研究所畜禽繁育基地飼養(yǎng)的383 只綿羊。其中蘇尼特羊（SNT）來自內(nèi)蒙古巴彥淖爾市烏拉特中旗，小尾寒羊（STH）來自山東鄆城。胸椎數(shù)13 根（T13）樣本259 只（SNT 122 只，STH 137只），胸椎數(shù)14 根（T14）樣本124 只（SNT 66 只，STH 58 只）。綿羊健康狀況良好，飼養(yǎng)環(huán)境一致。綿羊屠宰后，用2 mL RNase-Free 離心管快速采集新鮮背最長肌肌肉組織，用液氮冷凍運輸帶回實驗室轉移至-80℃超低溫冰箱保存?zhèn)溆?，記錄小尾寒羊胴體長度數(shù)據(jù)。

1.2 DNA提取

使用DNA 提取試劑盒（天根生物科技有限公司，北京）對組織樣品進行DNA 提取。提取后樣品DNA用NanoDrop2000 和1.5%瓊脂糖凝膠電泳進行濃度和質量評估，合格的DNA 樣品供后續(xù)使用。

1.3 基因分型

根據(jù)綿羊HERC1 和ALDH6A1 SNPs 位點序列信息，使用Assay design3.1（Sequenom 公司）引物設計軟件設計PCR 反應和單堿基擴展引物合成，引物信息見表1。采用Sequenom MassARRAY SNP[13]分型技術對HERC1 c.8950G ＞A（rs428482815）、g.43407058T ＞G、c.11902A ＞G（rs590493234）和ALDH6A1 c.1126A＞G（rs429743276）進行檢測分型，分型樣品為DNA，每個反應需要20L，DNA 濃度為40~80 ng/L。

表1 HERC1 與ALDH6A1 分型引物序列Table 1 Primer sequences for HERC1 and ALDH6A1 typing

1.4 數(shù)據(jù)統(tǒng)計分析

利用Excel 計算位點基因頻率、基因型頻率、有效等位基因數(shù)（Ne）、多態(tài)信息含量（PIC）、雜合度（He）及Hardy-Weinberg 平衡檢驗。計算PIC、He、有效等位基因數(shù)和Hardy-Weinberg 平衡檢驗等群體遺傳學參數(shù)方法詳見參考文獻[14]。利用SPSS19.0 運用一般線性模型單變量分析與單因素方差分析LSD 對HERC1和ALDH6A1 不同基因型綿羊群體進行胸椎數(shù)關聯(lián)分析。模型采用以下公式：

1.5 生物信息學分析

采用ProtParam（https://web.expasy.org/protparam/）分析目的蛋白理化性質；ProtScale（https://web.expasy.org/protscale/）分析目的蛋白疏水性；TMHMM Server v.2.0（http://www.cbs.dtu.dk/services/TMHMM/）分析目的蛋白跨膜區(qū)；運用SOPMA（https://npsa-prabi.ibcp.fr/cgi-bin/npsa_automat.pl? page=npsa_sopma.html）分析目的蛋白二級結構，Phyre2（https://npsa-prabi.ibcp.fr/cgi-bin/npsa_automat.pl? page=npsa_sopma.html）分析目的蛋白三級結構；STRING 11.0（https://string-db.org/cgi/input.pl）對目的蛋白進行蛋白互作分析。

2 結果與分析

2.1 HERC1 和ALDH6A1 多態(tài)性分析

383 只綿羊群體基因型分型（圖1）分別是：HERC1 c.8950G ＞A 位點GG 型365 只、GA 型15 只，g.43407058T ＞G 位點TT 型376 只、GG 型1 只、GT型4 只，c.11902A ＞G 位點AA 型364 只、GG 型3只、GT 型13 只；ALDH6A1 c.1126A ＞G 位點TT 型364 只、CT 型18 只、CC 型1 只。群體遺傳學分析表明（表2 和表3），STH 群體中HERC1 c.8950G ＞A、g.43407058T ＞G、c.11902A ＞G 位點的多態(tài)信息含量分別是0.035、0.030 和0.021，屬低度多態(tài)（PIC＜0.25）；SNT 群體中HERC1 基因c.8950G ＞A、c.11902A ＞G位點的多態(tài)信息含量分別是0.041、0.075，屬低度多態(tài)（PIC ＜0.25），SNT 群體中無g.43407058T ＞G 突變位點。STH 與SNT 群體中ALDH6A1 基因c.1126A ＞G位點的多態(tài)信息含量分別是0.084 及0.011，均屬低度多態(tài)（PIC ＜0.25）。群體中SNPs 位點卡方檢驗值表明：STH 群體中HERC1 c.8950G ＞A、g.43407058T ＞G、c.11902A ＞G 與ALDH6A1 c.1126A ＞G 分別是0.798、0、0.885 和0.342；SNT群體中HERC1c.8950G＞A、c.11902A ＞G 與ALDH6A1 c.1126A ＞G 分別是0.764、0 和0.941，SNT 群體中無g.43407058T ＞G 突變位點，不計算其群體遺傳學信息。結果表明，HERC1 基因g.43407058T ＞G 在STH 群體中不符合Hardy-Weinberg 平衡狀態(tài)，且c.11902A ＞G 位點在SNT 群體中也不符合Hardy-Weinberg 平衡狀態(tài)。其余SNPs位點在兩種綿羊群體中均處于Hardy-Weinberg平衡狀態(tài)（P ＞0.05）。

表2 HERC1 群體遺傳學分析Table 2 Population genetic analysis of HERC1

表3 ALDH6A1 群體遺傳學分析Table 3 Population genetic analysis of ALDH6A1

2.2 HERC1 和ALDH6A1 SNPs與胸椎數(shù)的關聯(lián)分析

為探究HERC1 和ALDH6A1 SNPs 對SNT 和STH 胸椎數(shù)的影響，利用SPSS19.0 單因素方差分析分別在不同群體中分析SNPs 對胸椎數(shù)的影響。結果發(fā)現(xiàn)，HERC1 c.8950G ＞A、g.43407058T ＞G、c.11902A＞G 在SNT 與STH 中均與胸椎數(shù)呈不顯著相關（P ＞0.05）；ALDH6A1 c.1126A > G 在SNT 與STH 中與胸椎數(shù)呈不顯著相關（P ＞0.05）（表4 和表5）。采用一般線性模型單變量分析品種和基因型對綿羊不同胸椎數(shù)的影響。結果表明，在考慮品種效應時，基因位點基因型對于綿羊胸椎數(shù)的影響均不顯著（P ＞0.05）（表6）。

表4 SNPs 位點基因型在不同品種中與胸椎數(shù)的關聯(lián)分析結果Table 4 Results of association analysis between genotypes of SNPs loci in different varieties and the number of thoracic vertebrae

表5 SNPs 位點基因型小尾寒羊和蘇尼特羊胸椎數(shù)最小二乘均值及標準誤Table 5 Least squares means and standard errors for the number of thoracic vertebrae of genotypes of SNPs loci in Small Tail Han sheep and Sunit sheep

表6 SNPs 位點基因型、品種與胸椎數(shù)關聯(lián)分析結果Table 6 Results of association analysis between genotype,variety and number of thoracic vertebrae at SNPs loci

2.3 HERC1 和ALDH6A1 SNPs位點與STH 胴體長的關聯(lián)分析

在對SNT 和STH 群體HERC1 和ALDH6A1 突變位點多態(tài)性檢測的基礎上，通過單因素方差分析對記錄胴體長數(shù)據(jù)的90只STH 進行SNPs 位點與胴體長的關聯(lián)分析（表7）。結果顯示，HERC1c.8950G ＞A 對STH 胴體 長 有 顯 著 影 響（P ＜ 0.05），HERC1 g.43407058T ＞G、c.11902A ＞G 與ALDH6A1 c.1126A ＞G 在STH 胴體長性狀中均不存在顯著相關（P ＞0.05）。

表7 SNPs 位點基因型與STH 胴體體長的關聯(lián)分析（平均值±標準誤）Table 7 Association analysis of genotypes of SNPs loci with STH carcass length (mean±standard error )

2.4 HERC1生物信息學分析

2.4.1 HERC1 蛋白理化性質分析 針對HERC1 c.8950G ＞A，經(jīng)ProtParam 在線軟件預測得知（表8），HERC1 蛋白共4 882 個氨基酸，分子量由突變前的534 217.22 增加到534 231.25，等電點5.7，其中亮氨酸（Leu）所占比例最高（11.50%），色氨酸（Trp）所占比例最低（1.60%），不穩(wěn)定性指數(shù)為47.42，屬于不穩(wěn)定蛋白。

表8 氨基酸數(shù)量與比例Table 8 Number and proportion of amino acids

2.4.2 HERC1 蛋白親/疏水性分析 運用ExPASy 服務器上的ProtScale 程序預測HERC1 蛋白突變前后的疏水性變化，結果顯示，HERC1 蛋白在1 766 位蘇氨酸（T）親水分數(shù)最大，達3.056，在2422位賴氨酸（K）親水分數(shù)最小，僅為-3.756，總平均親水性為-0.250，屬于親水蛋白（圖2），突變前后并未引起明顯的疏水性變化。

圖2 親水性分析Fig.2 Hydrophilic analysis

2.4.3 HERC1 蛋白跨膜結構預測與分析 運用TMHMM Serverv.2.0 對HERC1 蛋白跨膜結構進行分析（圖3），結果顯示，HERC1 蛋白突變前后均不存在跨膜螺旋區(qū)。

圖3 蛋白質跨膜結構預測（突變前后一致）Fig.3 Protein transmembrane structure prediction(consistent before and after mutation)

2.4.4 HERC1 蛋白二級結構預測 通過ExPASy 的SOPMA 程序，針對HERC1 c.8950G ＞A 突變位置，預測其蛋白二級結構（圖4），結果表明，突變前至突變后HERC1 蛋白-螺旋區(qū)域由1 720 個增加至1 721個，占比由37.22%增加至37.24%；-轉角區(qū)域由327個增加至328 個，占比由7.08% 增加至7.10%；無規(guī)卷曲區(qū)域由1 780 個減少至1 778 個，占比由38.52%減少至38.48%；延伸鏈區(qū)域并無變化，共794 個，占17.18%。

圖4 HERC1 蛋白質二級結構預測（截取部分，其他區(qū)域完全一致）Fig.4 Predicted secondary structure of HERC1 protein (intercepted part,other regions are identical)

2.4.5 HERC1 蛋白三級結構預測 采用Phyre2 軟件對HERC1 蛋白突變前后進行同源建模獲得三級結構模型。結果顯示，HERC1 蛋白三級結構與二級結構一致（圖5~6）。

圖5 HERC1 蛋白質三級結構預測Fig.5 Predicted tertiary structure of HERC1 protein

圖6 HERC1 蛋白質三級結構預測（氨基酸變化）Fig.6 Predicted tertiary structure of HERC1 protein (amino acid changes)

2.4.6 HERC1蛋白質的互作分析 采用STRING 11.0軟件對HERC1 蛋白進行蛋白質互作分析，結果如圖7 所示。

圖7 蛋白質互作分析Fig.7 Protein interactions analysis

3 討論

對綿羊等脊椎動物來說，脊椎的發(fā)育是從原腸胚中胚層開始的，其控制脊椎生長發(fā)育的初級神經(jīng)胚主要受到骨形成蛋白調(diào)控[15]。其中，各階段發(fā)育極為復雜，出現(xiàn)任何錯誤都可能造成脊椎分化異常[16]。由于豬繁殖快、產(chǎn)仔多的特點，使其在多脊椎性狀相關候選基因研究中較其他家畜成果更多、進展更快。因此，可以將家豬多脊椎性狀的研究作為參照，從遺傳育種的角度發(fā)掘更多對綿羊胸椎數(shù)有顯著影響的相關候選基因。目前，在綿羊多脊椎變異的研究中，已發(fā)現(xiàn)多個能夠引起脊椎數(shù)變異的有效基因。劉建明等[17]在STH脊椎數(shù)變異的關聯(lián)研究中指出，ActRⅡB 基因第85 位點發(fā)現(xiàn)的C ＞T 單核苷酸突變引發(fā)酶切位點消失，導致STH 第1 腰椎變異為第14 胸椎起到一定的作用，而該基因內(nèi)含子上的1 個SNP 也對小尾寒羊脊椎發(fā)育存在影響。張立嶺等[18]研究發(fā)現(xiàn)多脊椎畜禽在脊椎數(shù)目上比正常個體多1~2 個甚至更多，其中以蒙古羊和烏珠穆沁羊為首，多胸、腰椎個體的比例明顯高于其他品種。Krumlauf等[19]也指出動物脊椎數(shù)的變化是由調(diào)控脊椎發(fā)育的同源框基因（簡寫Hox）發(fā)生突變導致的。雖然已存在許多綿羊胸椎數(shù)變異等性狀的報道，但胸椎數(shù)變異受到多個QTLs的共同作用影響，其控制多脊椎性狀的分子機制仍不明晰，從遺傳育種角度有必要進一步發(fā)掘更多的對綿羊胸椎數(shù)有顯著影響的候選基因。

HERC1 c.8950G ＞A 和c.11902A ＞G 在兩種綿羊群體中表現(xiàn)為低度多態(tài)（PIC ＜0.25），表明這兩個位點存在較強的選擇潛力，這與Ensembl 數(shù)據(jù)庫發(fā)現(xiàn)世界山羊群體中HERC1 c.8950G ＞A 位點G 等位基因與c.11902A ＞G 位點A 等位基因為優(yōu)勢等位基因（G 等位基因頻率90%，A等位基因頻率97%）相一致。有趣的是，HERC1 g.43407058T ＞G 是一個新位點，在SNT 群體中并未發(fā)現(xiàn)該突變位點，且本研究發(fā)現(xiàn)，該位點的突變雜合型（GT）與突變純和型（GG）均高于野生型（GG），說明該基因突變位點在不同綿羊品種中可能存在選擇性，具有成為影響胸椎數(shù)性狀的候選基因潛力，這仍需進一步驗證研究。另外，STH 群體HERC1 g.43407058T ＞G 位點以及SNT群體HERC1c.11902A ＞G 位點處于Hardy-Weinberg不平衡狀態(tài)，可能是因為自然或人工選擇干預，也可能與樣品數(shù)量不足有關。關聯(lián)分析結果顯示，ALDH6A1 c.1126A>G 位點不同基因型與STH 及SNT 胸椎數(shù)和胴體長均無顯著關聯(lián)（P ＞0.05），但在胴體長關聯(lián)分析中，該位點突變雜合型（CT）胴體長長度均值比野生型（TT）多2 cm 左右，說明該位點仍有可能影響胴體長性狀，其中影響機制尚待研究發(fā)掘。而STH 群體中HERC1 c.8950G ＞A 位點在胴體長性狀上GG 型顯著高于GA 型，推測該位點突變在一定程度上降低了小尾寒羊的胴體產(chǎn)量。

通過HERC1 蛋白二級結構預測發(fā)現(xiàn)該蛋白沒有跨膜結構域，這與HERC1 定位在高爾基體和胞漿等細胞內(nèi)膜中[20]的結果相一致，表明該蛋白在細胞內(nèi)膜轉運中具有重要功能。而HERC1 c.8950G ＞A 位點的改變影響了HERC1 蛋白的二級與三級結構，使得突變前蛋白的兩個無規(guī)則卷曲區(qū)域轉變成一個-螺旋區(qū)域和一個-轉角區(qū)域，且使蛋白氨基酸序列中突變位點的纈氨酸轉變?yōu)橥蛔兒蟮漠惲涟彼?，推測其蛋白功能發(fā)生變化。同時，在群體遺傳學分析中也發(fā)現(xiàn)該突變與胴體長顯著相關。通過進一步對HERC1 蛋白進行蛋白質的互作分析發(fā)現(xiàn)，與其關聯(lián)的SOX10 在神經(jīng)發(fā)育及幼兒生長發(fā)育方面起著重要作用[21,22]。推測HERC1 可能與脊椎發(fā)育的起始原腸胚時期中胚層的初級神經(jīng)胚的發(fā)育有關，這需要進一步探索研究。

4 結論

綜上所述，HERC1 c.8950G ＞A、g.43407058T ＞G、c.11902A ＞G 和ALDH6A1 c.1126A>G 在小尾寒羊和蘇尼特羊中均與胸椎數(shù)不存在顯著相關（P ＞0.05）。小尾寒羊HERC1 c.8950G ＞A 位點與胴體長顯著相關（P ＜0.05）。綿羊HERC1 為親水蛋白，蛋白結構不穩(wěn)定，無跨膜結構域，可能不參與跨膜間生化反應。HERC1 c.8950G ＞A 位點可作為潛在影響綿羊胴體長度的分子輔助選擇標記。