王憲宗，宋 晶，楊春娟

（1. 山西農(nóng)業(yè)大學(xué)動(dòng)物科學(xué)學(xué)院，山西太谷 030801；2. 山西省水產(chǎn)技術(shù)推廣站，山西太原 030002）

對(duì)生物有機(jī)體來說，穩(wěn)定的氧供應(yīng)是維持細(xì)胞內(nèi)穩(wěn)態(tài)、保證各種組織正常發(fā)揮功能必不可少的因素。但在胚胎發(fā)育過程中細(xì)胞快速分裂的時(shí)期，腫瘤細(xì)胞急劇生長時(shí)，以及循環(huán)系統(tǒng)功能異常等情況下，細(xì)胞內(nèi)的氧平衡會(huì)被打破，進(jìn)而威脅有機(jī)體的生存[1‐4]。為了適應(yīng)低氧壓力，細(xì)胞和組織需要激活一系列參與血管生成、細(xì)胞分化或存活，以及葡萄糖和鐵代謝的基因的表達(dá)[4‐5]。在真核生物中，這一過程主要由低氧誘導(dǎo)因子-1a（Hypoxia‐inducible factor‐1a，Hif‐1a）介導(dǎo)。Hif-1a 的穩(wěn)定性受到氧濃度的嚴(yán)格調(diào)控，其在常氧條件下極易降解，但在低氧條件下由于FIH-1（Factor inhibiting HIF‐1）和pVHL（Von hippel-lindau protein）抑制作用的消失則可以大量積累。穩(wěn)定性得到增強(qiáng)的Hif-1a 會(huì)通過其C 末端的C-TAD 結(jié)構(gòu)域與EP300/CREBBP 蛋白結(jié)合，并通過其N 末端的bHLH 結(jié)構(gòu)域與Hif-1b 發(fā)生二聚化，形成有活性的轉(zhuǎn)錄激活復(fù)合物，最終與靶基因啟動(dòng)子中的低氧響應(yīng)元件（Hypoxia response elements，HRE）結(jié)合，并激活其表達(dá)[4，6]。

魚類生存所依賴的水環(huán)境中溶解氧的含量遠(yuǎn)低于陸地且波動(dòng)很大，因此與陸生動(dòng)物相比，魚類更易受到低氧脅迫[7‐9]。 冷水性的虹鱒（Oncorhynchus mykiss）對(duì)溶氧的要求高于多數(shù)魚類[10]，這就導(dǎo)致其養(yǎng)殖受自然條件的制約非常大，不利于這一優(yōu)良經(jīng)濟(jì)魚類[11‐13]的進(jìn)一步推廣。因此，對(duì)虹鱒hif-1a基因家族的拷貝數(shù)及結(jié)構(gòu)、功能進(jìn)行深入研究，可以對(duì)其不耐低氧的分子機(jī)制有更深入認(rèn)識(shí)，以期為虹鱒的遺傳改良提供新思路。

1 材料和方法

1.1 相關(guān)物種的Hif家族成員序列

本研究選擇斑馬魚（Danio rerio）的2 個(gè)Hif-1a成員作為查詢序列，選擇斑馬魚、小鼠（Mus musculus）和原雞（Gallus gallus）的Hif 家族成員作為重構(gòu)進(jìn)化樹的外群序列。通過檢索NCBI 的基因數(shù)據(jù)庫獲得已被注釋的3個(gè)物種的hif-1、hif-2等基因家族成員，各基因的詳細(xì)信息見表1。

1.2 BLAST搜索及結(jié)果提取

使用表1 中所列舉的2 條斑馬魚Hif-1a 蛋白序列進(jìn)行在線TBLASTN 搜索和BLASTP 搜索，參考數(shù)據(jù)庫分別是非冗余基因組數(shù)據(jù)庫（refseq_genomes）和非冗余蛋白質(zhì)數(shù)據(jù)庫（refseq_protein），物種均限定為虹鱒。搜索參數(shù)中max target sequences 均設(shè)置為5 000，e-value 則設(shè)置為1e-5。編寫python 腳本從XML格式的結(jié)果文件中提取搜索結(jié)果，對(duì)覆蓋度低于30%的hits進(jìn)行過濾。BLASTP搜索結(jié)果中，將具有顯著性的hits 所對(duì)應(yīng)基因的詳細(xì)信息，從gene2accession 文件（下載自https://ftp.ncbi.nlm.nih.gov/gene/DATA/）中提取。

1.3 系統(tǒng)發(fā)育分析

基于對(duì)BLASTP 搜索結(jié)果的分析，下載符合篩選標(biāo)準(zhǔn)的虹鱒候選Hif-1a 蛋白序列，使用MAFFT軟件[14]對(duì)這些序列及表1 中所列舉的斑馬魚、小鼠和原雞的Hif 蛋白序列進(jìn)行多重比對(duì)，參數(shù)選擇最精確的L-INS-i 模式。由于Hif-1a 序列較短，為避免信息過度丟失，不對(duì)多重比對(duì)結(jié)果進(jìn)行修剪，而是直接使用RAxML 8.2.8 軟件[15]重構(gòu)最大似然樹，采用GAMMA 速率異質(zhì)性模型，自動(dòng)選擇氨基酸替代模型，自展抽樣次數(shù)設(shè)為500。

1.4 基因表達(dá)及功能富集分析

NCBI 登記的PRJEB37848 項(xiàng)目包含了對(duì)虹鱒和斑馬魚10余個(gè)組織進(jìn)行轉(zhuǎn)錄組測序的原始數(shù)據(jù)，下載這2 種魚的全部測序數(shù)據(jù)，使用sratoolkit 工具包中的fasterq-dump 進(jìn)行數(shù)據(jù)格式轉(zhuǎn)換?；谏鲜鲈紨?shù)據(jù)使用Salmon[16]以SAF genomic 索引模式對(duì)2 種魚非冗余轉(zhuǎn)錄本（下載自https://ftp.ncbi.nlm.nih.gov/genomes/refseq/vertebrate_other/）的組織表達(dá)水平進(jìn)行計(jì)數(shù)。Salmon 衡量一個(gè)轉(zhuǎn)錄本相對(duì)表達(dá)水平的單位是TPM（Transcripts per million），即每100萬個(gè)測序生成的原始reads里有多少個(gè)來自于該轉(zhuǎn)錄本。基因的TPM 則是將該基因所有轉(zhuǎn)錄本的TPM 相加的結(jié)果，基因與轉(zhuǎn)錄本的對(duì)應(yīng)關(guān)系同樣從gene2accession文件獲取。

使用皮爾遜相關(guān)系數(shù)（Pearson correlation coefficient）衡量不同基因之間組織表達(dá)譜的相關(guān)性，由基于python 語言的scipy.stats 模塊[17]的pearsonr 函數(shù)計(jì)算而來，該函數(shù)同時(shí)會(huì)返回衡量非相關(guān)性的P值。設(shè)置合理閾值（r>0.9，同時(shí)要求P<0.01）篩選各hif-1a成員的共表達(dá)基因。將虹鱒hif-1a共表達(dá)基因編碼的蛋白質(zhì)序列參考斑馬魚非冗余的蛋白質(zhì)數(shù)據(jù)庫進(jìn)行本地BLASTP 搜索，基于gene2accession 文件中基因和蛋白質(zhì)序列的對(duì)應(yīng)關(guān)系將虹鱒的基因轉(zhuǎn)換為直系同源的斑馬魚基因，然后使用GOATOOLs[18]對(duì)轉(zhuǎn)換后的共表達(dá)基因進(jìn)行功能富集分析（斑馬魚的共表達(dá)基因則無需轉(zhuǎn)換，可以直接進(jìn)行富集分析）。

1.5 保守結(jié)構(gòu)域搜索及三維結(jié)構(gòu)建模

基于對(duì)hif-1a基因不同拷貝表達(dá)分析的結(jié)果，確定各拷貝表達(dá)量最高的轉(zhuǎn)錄本，提取其編碼的蛋白質(zhì)序列，使用CDD search 在線程序（https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi）[19‐20]搜索其中所包含的保守結(jié)構(gòu)域。根據(jù)CDD search 搜索結(jié)果提取各序列中bHLH 結(jié)構(gòu)域和C-TAD 結(jié)構(gòu)域所對(duì)應(yīng)的子序列，使用部署在本地的I-TASSER 5.1 軟件包[21]對(duì)這2 個(gè)結(jié)構(gòu)域進(jìn)行三維結(jié)構(gòu)建模，建模得到的最佳模型使用PyMOL[22]進(jìn)行可視化。

2 結(jié)果與分析

2.1 虹鱒hif-1a家族的拷貝數(shù)

TBLASTN 搜索結(jié)果顯示，2 條斑馬魚Hif-1a 查詢序列在虹鱒染色體上有6 個(gè)較長的匹配區(qū)間，它們的長度在6 800～36 464 nt（表2）。由于這些匹配區(qū)間對(duì)查詢序列的覆蓋度較低，最高僅54.18%，難以確定基因座。通過BLASTP 搜索得到8 條覆蓋度較高的蛋白質(zhì)序列（hits），其中6條序列所屬基因座恰好與TBLASTN 搜索得到的6 個(gè)匹配區(qū)間一致（表3）。將上述的8 條序列與斑馬魚、小鼠和原雞的Hif-1a 蛋白序列進(jìn)行比對(duì)，并重構(gòu)最大似然樹，結(jié)果顯示，有2 條序列與斑馬魚的Hif-1aa 聚在一起，但沒有序列與斑馬魚的Hif-1ab 聚在一起（圖1）。最大似然樹的拓?fù)浣Y(jié)構(gòu)表明，虹鱒只有2個(gè)與斑馬魚hif-1aa直系同源的基因，但沒有與斑馬魚hif-1ab直系同源的基因。根據(jù)通行的命名規(guī)則，將這2個(gè)成員分別命名為hif-1aaa（Hit gene：100135944）和hif-1aab（Hit gene：110505629）。

表2 斑馬魚Hif-1a參考虹鱒基因組的TBLASTN搜索結(jié)果Tab.2 TBLASTN search results of D.rerio Hif-1a against O.mykiss genome

表3 斑馬魚Hif-1a參考虹鱒蛋白質(zhì)數(shù)據(jù)庫的BLASTP搜索結(jié)果Tab.3 BLASTP search results of D.rerio Hif-1a against O.mykiss protein database

續(xù)表3 斑馬魚Hif-1a參考虹鱒蛋白質(zhì)數(shù)據(jù)庫的BLASTP搜索結(jié)果Tab.3（Continued） BLASTP search results of D.rerio Hif-1a against O.mykiss protein database

2.2 虹鱒hif-1a家族不同成員的組織表達(dá)特征

基于對(duì)PRJEB37848 項(xiàng)目轉(zhuǎn)錄組測序原始數(shù)據(jù)的分析，發(fā)現(xiàn)無論是基于12個(gè)組織還是基于和斑馬魚相同的8 個(gè)組織，虹鱒2 個(gè)拷貝的表達(dá)譜均不存在顯著相關(guān)性（r<0.5，P>0.1）。Hif-1aaa的最高表達(dá)量出現(xiàn)在心臟，而hif-1aab的最高表達(dá)量出現(xiàn)在鰓，二者在其他組織中也有較高的表達(dá)水平（圖2a和圖2b）。斑馬魚2 個(gè)拷貝的表達(dá)譜也不存在顯著相關(guān)性，其中hif-1aa僅在心臟中有較高水平的表達(dá)，而hif-1ab則在眼、腦、鰓、皮膚和心臟中都有較高水平的表達(dá)（圖2c和圖2d）。顯然，無論是虹鱒還是斑馬魚，它們所擁有的2 個(gè)拷貝之間均發(fā)生了較為明顯的功能分化。虹鱒的2 個(gè)拷貝與斑馬魚的hif-1aa均是直系同源的關(guān)系，但僅有hif-1aaa與斑馬魚hif-1aa的表達(dá)譜存在顯著相關(guān)性，而且這種顯著相關(guān)性是基于二者均在心臟中有最高水平的表達(dá)。當(dāng)排除掉心臟的表達(dá)數(shù)據(jù)后，二者的表達(dá)譜也不再顯著相關(guān)，表明虹鱒和斑馬魚的hif-1a成員之間也發(fā)生了功能分化。

2.3 虹鱒hif-1a家族不同成員共表達(dá)基因的功能富集分析

以表達(dá)譜相關(guān)系數(shù)大于0.9為閾值（同時(shí)要求P值小于0.01），發(fā)現(xiàn)虹鱒的2 個(gè)hif-1a拷貝都有與其表達(dá)譜顯著相關(guān)的基因，而且這2 個(gè)基因集合之間并不存在重復(fù)，這與2 個(gè)拷貝自身表達(dá)譜不存在顯著相關(guān)是一致的。對(duì)2個(gè)拷貝的共表達(dá)基因分別進(jìn)行功能富集分析，發(fā)現(xiàn)它們都能富集到一定數(shù)量的與血管生成、細(xì)胞分化、物質(zhì)運(yùn)輸及代謝相關(guān)的GO條目（表4），這四大類GO 條目代表了低氧條件下Hif-1a 蛋白的核心功能?？梢?，每個(gè)拷貝都仍具有其原始功能，但共表達(dá)基因集合的分離表明它們所作用的下游基因發(fā)生了分化。hif-1aab的共表達(dá)基因數(shù)量是hif-1aaa的2 倍以上，但它所富集到的四大類GO 條目并不總是多于后者，表明它們二者所繼承的原始功能各有側(cè)重。斑馬魚上這一特征更為明顯，hif-1aa的總體表達(dá)水平比hif-1ab低（圖2），而且共表達(dá)基因數(shù)量也更少，但它所富集的四大類GO 條目中有三大類的數(shù)量都多于后者。hif-1aa的共表達(dá)基因所富集到的與血管生成相關(guān)的GO 條目遠(yuǎn)多于hif-1ab，這與它僅在心臟中有較高水平的表達(dá)相一致。

表4 虹鱒和斑馬魚hif-1a家族成員共表達(dá)基因的功能富集結(jié)果Tab.4 GO enrichment analysis of co-expressed genes of hif-1a members of O.mykiss and D.rerio個(gè)

2.4 虹鱒Hif-1a家族不同拷貝的結(jié)構(gòu)特征

使用CDD search 在線程序?qū)琪V、斑馬魚及小鼠的Hif-1a 成員的保守結(jié)構(gòu)域進(jìn)行搜索，發(fā)現(xiàn)5 條序列中均有5 個(gè)保守結(jié)構(gòu)域（圖3），對(duì)應(yīng)Hif-1a 蛋白中的1 個(gè)堿性螺旋-環(huán)-螺旋（Basic helix‐loop‐helix，bHLH）結(jié)構(gòu)域、2 個(gè)PAS（Per/ARNT/Sim）結(jié)構(gòu)域，以及2 個(gè)轉(zhuǎn)錄激活結(jié)構(gòu)域N-TAD 和C-TAD（N/C‐terminal activation domain）。CDD search 搜索的結(jié)果還顯示，虹鱒的2 條序列各有1 個(gè)結(jié)構(gòu)域的特征不夠明顯，分別是N 末端的bHLH 結(jié)構(gòu)域和C 末端的C-TAD 結(jié)構(gòu)域；斑馬魚的2 條序列中這2 個(gè)結(jié)構(gòu)域也存在特征不夠明顯的情況（圖3）。這2 個(gè)結(jié)構(gòu)域所發(fā)生的結(jié)構(gòu)變異很可能是虹鱒和斑馬魚中2個(gè)拷貝功能分化的基礎(chǔ)，因此，使用I-TASSER 5.1 軟件對(duì)5條序列的bHLH 和C-TAD 結(jié)構(gòu)域分別進(jìn)行結(jié)構(gòu)建模。對(duì)于bHLH 結(jié)構(gòu)域，從表5 可以看出，虹鱒Hif-1aaa 最佳模型的C-score 和TM-score 都高于其他序列，表明它的結(jié)構(gòu)特征很可能更為突出，從而在結(jié)構(gòu)搜索時(shí)更易收斂。從圖4 可以看出，虹鱒Hif-1aaa 的2條長螺旋都非常完整，而其他4 條序列，即使是小鼠的，也有1條長螺旋存在中間斷開的情況。對(duì)于C-TAD 結(jié)構(gòu)域，綜合5 條序列的建模結(jié)果，可以確定該結(jié)構(gòu)域原本應(yīng)該只在N 端存在1 個(gè)由十幾個(gè)氨基酸所形成的螺旋，但虹鱒的2 條序列均在C 端存在1個(gè)短螺旋，表明虹鱒在此結(jié)構(gòu)域很可能發(fā)生了結(jié)構(gòu)創(chuàng)新（圖5）。斑馬魚Hif-1aa 的CTAD 結(jié)構(gòu)域中螺旋明顯短于其他序列（圖5c），很可能發(fā)生了結(jié)構(gòu)退化。

圖3 虹鱒、斑馬魚和小鼠Hif-1a成員的CDD search搜索結(jié)果Fig.3 CDD search results of Hif-1a members of O.mykiss，D.rerio and M.musculus

圖4 虹鱒、斑馬魚和小鼠Hif-1a成員bHLH結(jié)構(gòu)域的最佳模型Fig.4 Best models of bHLH domain of Hif-1a members of O.mykiss，D.rerio and M.musculus

圖5 虹鱒、斑馬魚和小鼠Hif-1a成員C-TAD結(jié)構(gòu)域的最佳模型Fig.5 Best models of C-TAD domain of Hif-1a members of O.mykiss，D.rerio and M.musculus

表5 虹鱒、斑馬魚和小鼠Hif-1a成員各結(jié)構(gòu)域最佳模型的準(zhǔn)確度指標(biāo)Tab.5 Estimated accuracy of best model of each domain of Hif-1a members of O.mykiss，D.rerio and M.musculus

3 結(jié)論與討論

真骨魚類共同祖先從脊椎動(dòng)物的祖先分支中分化出來后在3.2 億～3.5 億年前經(jīng)歷過1 次全基因組重復(fù)（Whole genome duplication，WGD）事件[23]，多數(shù)基因都在較短時(shí)間內(nèi)丟失了其中1 個(gè)拷貝，重新變?yōu)閱慰截悹顟B(tài)[24]，但一些較為重要的調(diào)節(jié)基因卻能夠在較長的進(jìn)化時(shí)期保持多拷貝的狀態(tài)[25]，如斑馬魚的hif-1a仍保留有2 個(gè)拷貝。虹鱒較近的祖先物種還經(jīng)歷過1 次額外的全基因組重復(fù)事件（約1億年前）[26‐27]，其hif-1a基因家族理論上應(yīng)該有4 個(gè)拷貝。本研究綜合TBLASTN、BLASTP 搜索及重構(gòu)的最大似然樹結(jié)果得出，虹鱒hif-1a家族只有2 個(gè)成員，這有2種解釋：一是虹鱒現(xiàn)有的基因組測序及拼接精度欠佳，導(dǎo)致大段區(qū)間缺失，而未被發(fā)現(xiàn)的hif-1a成員恰好在這些缺失的區(qū)間；二是虹鱒hif-1a家族發(fā)生過基因丟失。本研究中，虹鱒Arlee品系基因組序列由PacBio 平臺(tái)生成的長測序片段拼接而來，原始數(shù)據(jù)及拼接結(jié)果的各項(xiàng)指標(biāo)均遠(yuǎn)高于之前的Swanson 品系，大段區(qū)間缺失的可能性非常小[28]。而且筆者同期開展的另一項(xiàng)研究表明，虹鱒的ep300/crebbp家族有8 個(gè)成員，恰好是斑馬魚上該家族成員數(shù)量的2倍（數(shù)據(jù)尚未發(fā)表）。綜合來看，是基因丟失導(dǎo)致了虹鱒的hif-1a家族現(xiàn)在只有2個(gè)成員。虹鱒現(xiàn)有的2 個(gè)hif-1a成員均是斑馬魚hif-1aa的直系同源基因，表明其hif-1ab基因很可能在額外的全基因組重復(fù)事件之前已經(jīng)丟失。

由全基因組重復(fù)所產(chǎn)生的多拷貝基因，其各個(gè)拷貝的進(jìn)化速率通常都存在差異，隨著時(shí)間的推移往往會(huì)以新功能化或亞功能化的形式發(fā)生分歧，最終的結(jié)局是形成新的基因或者部分拷貝假基因化[29]。本研究發(fā)現(xiàn)，無論是虹鱒還是斑馬魚，它們各自的2個(gè)hif-1a拷貝在組織表達(dá)譜上都發(fā)生了顯著的分化，表明它們在不同組織中所發(fā)揮的生物學(xué)功能也已經(jīng)發(fā)生了一定程度的分化。Hif-1a 通路所介導(dǎo)的低氧適應(yīng)需要激活上百個(gè)基因的表達(dá)，這些基因所參與的功能包括血管生成、細(xì)胞分化/存活和糖代謝等，即通過增強(qiáng)循環(huán)系統(tǒng)來增加氧的供應(yīng)或者增強(qiáng)無氧代謝的強(qiáng)度以獲得足夠的能量[4]。虹鱒或斑馬魚2個(gè)拷貝的共表達(dá)基因都能夠富集到一定數(shù)量的與血管生成、細(xì)胞分化、物質(zhì)運(yùn)輸及物質(zhì)代謝相關(guān)的GO 條目，可見它們大體上仍保留了原始基因的核心功能。但2個(gè)拷貝的共表達(dá)基因并無交集，表明信號(hào)通路本身也已經(jīng)發(fā)生了分化，2個(gè)拷貝分歧的形式主要是亞功能化。

虹鱒Hif-1aaa 和Hif-1aab 在bHLH 結(jié)構(gòu)域上存在較大結(jié)構(gòu)差異，同樣支持它們發(fā)生功能分歧的結(jié)論。本研究還發(fā)現(xiàn)，虹鱒Hif-1aaa 和Hif-1aab 同樣存在結(jié)構(gòu)創(chuàng)新。就C-TAD結(jié)構(gòu)域而言，小鼠和斑馬魚的序列建模都只能得到1 條螺旋，而虹鱒的2 條序列除了能得到比前2 個(gè)物種更長的螺旋外，還能得到1 條較短的螺旋。低氧狀態(tài)下，C-TAD 結(jié)構(gòu)域的主要功能是募集轉(zhuǎn)錄激活因子Ep300/Crebbp，而后者自身也擁有非常多的結(jié)合結(jié)構(gòu)域[30]，C-TAD 結(jié)構(gòu)域中額外的較短螺旋很可能有助于它們與Ep300/Crebbp 結(jié)合。虹鱒的Ep300/Crebbp 家族未發(fā)生過拷貝丟失，但這些拷貝之間必然也會(huì)發(fā)生功能分歧[31]，C-TAD 結(jié)構(gòu)域的改變很可能就是為了適應(yīng)不同Ep300/Crebbp 拷貝發(fā)生的結(jié)構(gòu)變異[32]，從而實(shí)現(xiàn)更精細(xì)調(diào)控，其主要作用仍然是促進(jìn)2 個(gè)拷貝的亞功能化。

綜合來看，虹鱒的hif-1a基因家族在最近的全基因組重復(fù)事件之前丟失過拷貝，其祖先物種在某段時(shí)期內(nèi)僅有1 個(gè)拷貝，這可能是它不耐低氧的原因之一。虹鱒和斑馬魚的2 個(gè)hif-1a拷貝之間較大的功能分歧表明，該基因家族的進(jìn)化速度較快，未來的研究可以關(guān)注虹鱒不同地理種群中2 個(gè)hif-1a拷貝的遺傳變異，特別是這些遺傳變異所包含的結(jié)構(gòu)與功能改變，從而找到有助于提高其耐低氧性能的基因型。