虹鳟hif-1a 基因家族的拷贝数变异及分歧进化

王宪宗，宋 晶，杨春娟

（1. 山西农业大学动物科学学院，山西太谷 030801；2. 山西省水产技术推广站，山西太原 030002）

对生物有机体来说，稳定的氧供应是维持细胞内稳态、保证各种组织正常发挥功能必不可少的因素。但在胚胎发育过程中细胞快速分裂的时期，肿瘤细胞急剧生长时，以及循环系统功能异常等情况下，细胞内的氧平衡会被打破，进而威胁有机体的生存[1‐4]。为了适应低氧压力，细胞和组织需要激活一系列参与血管生成、细胞分化或存活，以及葡萄糖和铁代谢的基因的表达[4‐5]。在真核生物中，这一过程主要由低氧诱导因子-1a（Hypoxia‐inducible factor‐1a，Hif‐1a）介导。Hif-1a 的稳定性受到氧浓度的严格调控，其在常氧条件下极易降解，但在低氧条件下由于FIH-1（Factor inhibiting HIF‐1）和pVHL（Von hippel-lindau protein）抑制作用的消失则可以大量积累。稳定性得到增强的Hif-1a 会通过其C 末端的C-TAD 结构域与EP300/CREBBP 蛋白结合，并通过其N 末端的bHLH 结构域与Hif-1b 发生二聚化，形成有活性的转录激活复合物，最终与靶基因启动子中的低氧响应元件（Hypoxia response elements，HRE）结合，并激活其表达[4，6]。

鱼类生存所依赖的水环境中溶解氧的含量远低于陆地且波动很大，因此与陆生动物相比，鱼类更易受到低氧胁迫[7‐9]。 冷水性的虹鳟（Oncorhynchus mykiss）对溶氧的要求高于多数鱼类[10]，这就导致其养殖受自然条件的制约非常大，不利于这一优良经济鱼类[11‐13]的进一步推广。因此，对虹鳟hif-1a基因家族的拷贝数及结构、功能进行深入研究，可以对其不耐低氧的分子机制有更深入认识，以期为虹鳟的遗传改良提供新思路。

1 材料和方法

1.1 相关物种的Hif家族成员序列

本研究选择斑马鱼（Danio rerio）的2 个Hif-1a成员作为查询序列，选择斑马鱼、小鼠（Mus musculus）和原鸡（Gallus gallus）的Hif 家族成员作为重构进化树的外群序列。通过检索NCBI 的基因数据库获得已被注释的3个物种的hif-1、hif-2等基因家族成员，各基因的详细信息见表1。

1.2 BLAST搜索及结果提取

使用表1 中所列举的2 条斑马鱼Hif-1a 蛋白序列进行在线TBLASTN 搜索和BLASTP 搜索，参考数据库分别是非冗余基因组数据库（refseq_genomes）和非冗余蛋白质数据库（refseq_protein），物种均限定为虹鳟。搜索参数中max target sequences 均设置为5 000，e-value 则设置为1e-5。编写python 脚本从XML格式的结果文件中提取搜索结果，对覆盖度低于30%的hits进行过滤。BLASTP搜索结果中，将具有显著性的hits 所对应基因的详细信息，从gene2accession 文件（下载自https://ftp.ncbi.nlm.nih.gov/gene/DATA/）中提取。

1.3 系统发育分析

基于对BLASTP 搜索结果的分析，下载符合筛选标准的虹鳟候选Hif-1a 蛋白序列，使用MAFFT软件[14]对这些序列及表1 中所列举的斑马鱼、小鼠和原鸡的Hif 蛋白序列进行多重比对，参数选择最精确的L-INS-i 模式。由于Hif-1a 序列较短，为避免信息过度丢失，不对多重比对结果进行修剪，而是直接使用RAxML 8.2.8 软件[15]重构最大似然树，采用GAMMA 速率异质性模型，自动选择氨基酸替代模型，自展抽样次数设为500。

1.4 基因表达及功能富集分析

NCBI 登记的PRJEB37848 项目包含了对虹鳟和斑马鱼10余个组织进行转录组测序的原始数据，下载这2 种鱼的全部测序数据，使用sratoolkit 工具包中的fasterq-dump 进行数据格式转换。基于上述原始数据使用Salmon[16]以SAF genomic 索引模式对2 种鱼非冗余转录本（下载自https://ftp.ncbi.nlm.nih.gov/genomes/refseq/vertebrate_other/）的组织表达水平进行计数。Salmon 衡量一个转录本相对表达水平的单位是TPM（Transcripts per million），即每100万个测序生成的原始reads里有多少个来自于该转录本。基因的TPM 则是将该基因所有转录本的TPM 相加的结果，基因与转录本的对应关系同样从gene2accession文件获取。

使用皮尔逊相关系数（Pearson correlation coefficient）衡量不同基因之间组织表达谱的相关性，由基于python 语言的scipy.stats 模块[17]的pearsonr 函数计算而来，该函数同时会返回衡量非相关性的P值。设置合理阈值（r>0.9，同时要求P<0.01）筛选各hif-1a成员的共表达基因。将虹鳟hif-1a共表达基因编码的蛋白质序列参考斑马鱼非冗余的蛋白质数据库进行本地BLASTP 搜索，基于gene2accession 文件中基因和蛋白质序列的对应关系将虹鳟的基因转换为直系同源的斑马鱼基因，然后使用GOATOOLs[18]对转换后的共表达基因进行功能富集分析（斑马鱼的共表达基因则无需转换，可以直接进行富集分析）。

1.5 保守结构域搜索及三维结构建模

基于对hif-1a基因不同拷贝表达分析的结果，确定各拷贝表达量最高的转录本，提取其编码的蛋白质序列，使用CDD search 在线程序（https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi）[19‐20]搜索其中所包含的保守结构域。根据CDD search 搜索结果提取各序列中bHLH 结构域和C-TAD 结构域所对应的子序列，使用部署在本地的I-TASSER 5.1 软件包[21]对这2 个结构域进行三维结构建模，建模得到的最佳模型使用PyMOL[22]进行可视化。

2 结果与分析

2.1 虹鳟hif-1a家族的拷贝数

TBLASTN 搜索结果显示，2 条斑马鱼Hif-1a 查询序列在虹鳟染色体上有6 个较长的匹配区间，它们的长度在6 800～36 464 nt（表2）。由于这些匹配区间对查询序列的覆盖度较低，最高仅54.18%，难以确定基因座。通过BLASTP 搜索得到8 条覆盖度较高的蛋白质序列（hits），其中6条序列所属基因座恰好与TBLASTN 搜索得到的6 个匹配区间一致（表3）。将上述的8 条序列与斑马鱼、小鼠和原鸡的Hif-1a 蛋白序列进行比对，并重构最大似然树，结果显示，有2 条序列与斑马鱼的Hif-1aa 聚在一起，但没有序列与斑马鱼的Hif-1ab 聚在一起（图1）。最大似然树的拓扑结构表明，虹鳟只有2个与斑马鱼hif-1aa直系同源的基因，但没有与斑马鱼hif-1ab直系同源的基因。根据通行的命名规则，将这2个成员分别命名为hif-1aaa（Hit gene：100135944）和hif-1aab（Hit gene：110505629）。

表2 斑马鱼Hif-1a参考虹鳟基因组的TBLASTN搜索结果Tab.2 TBLASTN search results of D.rerio Hif-1a against O.mykiss genome

表3 斑马鱼Hif-1a参考虹鳟蛋白质数据库的BLASTP搜索结果Tab.3 BLASTP search results of D.rerio Hif-1a against O.mykiss protein database

续表3 斑马鱼Hif-1a参考虹鳟蛋白质数据库的BLASTP搜索结果Tab.3（Continued） BLASTP search results of D.rerio Hif-1a against O.mykiss protein database

2.2 虹鳟hif-1a家族不同成员的组织表达特征

基于对PRJEB37848 项目转录组测序原始数据的分析，发现无论是基于12个组织还是基于和斑马鱼相同的8 个组织，虹鳟2 个拷贝的表达谱均不存在显著相关性（r<0.5，P>0.1）。Hif-1aaa的最高表达量出现在心脏，而hif-1aab的最高表达量出现在鳃，二者在其他组织中也有较高的表达水平（图2a和图2b）。斑马鱼2 个拷贝的表达谱也不存在显著相关性，其中hif-1aa仅在心脏中有较高水平的表达，而hif-1ab则在眼、脑、鳃、皮肤和心脏中都有较高水平的表达（图2c和图2d）。显然，无论是虹鳟还是斑马鱼，它们所拥有的2 个拷贝之间均发生了较为明显的功能分化。虹鳟的2 个拷贝与斑马鱼的hif-1aa均是直系同源的关系，但仅有hif-1aaa与斑马鱼hif-1aa的表达谱存在显著相关性，而且这种显著相关性是基于二者均在心脏中有最高水平的表达。当排除掉心脏的表达数据后，二者的表达谱也不再显著相关，表明虹鳟和斑马鱼的hif-1a成员之间也发生了功能分化。

2.3 虹鳟hif-1a家族不同成员共表达基因的功能富集分析

以表达谱相关系数大于0.9为阈值（同时要求P值小于0.01），发现虹鳟的2 个hif-1a拷贝都有与其表达谱显著相关的基因，而且这2 个基因集合之间并不存在重复，这与2 个拷贝自身表达谱不存在显著相关是一致的。对2个拷贝的共表达基因分别进行功能富集分析，发现它们都能富集到一定数量的与血管生成、细胞分化、物质运输及代谢相关的GO条目（表4），这四大类GO 条目代表了低氧条件下Hif-1a 蛋白的核心功能。可见，每个拷贝都仍具有其原始功能，但共表达基因集合的分离表明它们所作用的下游基因发生了分化。hif-1aab的共表达基因数量是hif-1aaa的2 倍以上，但它所富集到的四大类GO 条目并不总是多于后者，表明它们二者所继承的原始功能各有侧重。斑马鱼上这一特征更为明显，hif-1aa的总体表达水平比hif-1ab低（图2），而且共表达基因数量也更少，但它所富集的四大类GO 条目中有三大类的数量都多于后者。hif-1aa的共表达基因所富集到的与血管生成相关的GO 条目远多于hif-1ab，这与它仅在心脏中有较高水平的表达相一致。

表4 虹鳟和斑马鱼hif-1a家族成员共表达基因的功能富集结果Tab.4 GO enrichment analysis of co-expressed genes of hif-1a members of O.mykiss and D.rerio个

2.4 虹鳟Hif-1a家族不同拷贝的结构特征

使用CDD search 在线程序对虹鳟、斑马鱼及小鼠的Hif-1a 成员的保守结构域进行搜索，发现5 条序列中均有5 个保守结构域（图3），对应Hif-1a 蛋白中的1 个碱性螺旋-环-螺旋（Basic helix‐loop‐helix，bHLH）结构域、2 个PAS（Per/ARNT/Sim）结构域，以及2 个转录激活结构域N-TAD 和C-TAD（N/C‐terminal activation domain）。CDD search 搜索的结果还显示，虹鳟的2 条序列各有1 个结构域的特征不够明显，分别是N 末端的bHLH 结构域和C 末端的C-TAD 结构域；斑马鱼的2 条序列中这2 个结构域也存在特征不够明显的情况（图3）。这2 个结构域所发生的结构变异很可能是虹鳟和斑马鱼中2个拷贝功能分化的基础，因此，使用I-TASSER 5.1 软件对5条序列的bHLH 和C-TAD 结构域分别进行结构建模。对于bHLH 结构域，从表5 可以看出，虹鳟Hif-1aaa 最佳模型的C-score 和TM-score 都高于其他序列，表明它的结构特征很可能更为突出，从而在结构搜索时更易收敛。从图4 可以看出，虹鳟Hif-1aaa 的2条长螺旋都非常完整，而其他4 条序列，即使是小鼠的，也有1条长螺旋存在中间断开的情况。对于C-TAD 结构域，综合5 条序列的建模结果，可以确定该结构域原本应该只在N 端存在1 个由十几个氨基酸所形成的螺旋，但虹鳟的2 条序列均在C 端存在1个短螺旋，表明虹鳟在此结构域很可能发生了结构创新（图5）。斑马鱼Hif-1aa 的CTAD 结构域中螺旋明显短于其他序列（图5c），很可能发生了结构退化。

图3 虹鳟、斑马鱼和小鼠Hif-1a成员的CDD search搜索结果Fig.3 CDD search results of Hif-1a members of O.mykiss，D.rerio and M.musculus

图4 虹鳟、斑马鱼和小鼠Hif-1a成员bHLH结构域的最佳模型Fig.4 Best models of bHLH domain of Hif-1a members of O.mykiss，D.rerio and M.musculus

图5 虹鳟、斑马鱼和小鼠Hif-1a成员C-TAD结构域的最佳模型Fig.5 Best models of C-TAD domain of Hif-1a members of O.mykiss，D.rerio and M.musculus

表5 虹鳟、斑马鱼和小鼠Hif-1a成员各结构域最佳模型的准确度指标Tab.5 Estimated accuracy of best model of each domain of Hif-1a members of O.mykiss，D.rerio and M.musculus

3 结论与讨论

真骨鱼类共同祖先从脊椎动物的祖先分支中分化出来后在3.2 亿～3.5 亿年前经历过1 次全基因组重复（Whole genome duplication，WGD）事件[23]，多数基因都在较短时间内丢失了其中1 个拷贝，重新变为单拷贝状态[24]，但一些较为重要的调节基因却能够在较长的进化时期保持多拷贝的状态[25]，如斑马鱼的hif-1a仍保留有2 个拷贝。虹鳟较近的祖先物种还经历过1 次额外的全基因组重复事件（约1亿年前）[26‐27]，其hif-1a基因家族理论上应该有4 个拷贝。本研究综合TBLASTN、BLASTP 搜索及重构的最大似然树结果得出，虹鳟hif-1a家族只有2 个成员，这有2种解释：一是虹鳟现有的基因组测序及拼接精度欠佳，导致大段区间缺失，而未被发现的hif-1a成员恰好在这些缺失的区间；二是虹鳟hif-1a家族发生过基因丢失。本研究中，虹鳟Arlee品系基因组序列由PacBio 平台生成的长测序片段拼接而来，原始数据及拼接结果的各项指标均远高于之前的Swanson 品系，大段区间缺失的可能性非常小[28]。而且笔者同期开展的另一项研究表明，虹鳟的ep300/crebbp家族有8 个成员，恰好是斑马鱼上该家族成员数量的2倍（数据尚未发表）。综合来看，是基因丢失导致了虹鳟的hif-1a家族现在只有2个成员。虹鳟现有的2 个hif-1a成员均是斑马鱼hif-1aa的直系同源基因，表明其hif-1ab基因很可能在额外的全基因组重复事件之前已经丢失。

由全基因组重复所产生的多拷贝基因，其各个拷贝的进化速率通常都存在差异，随着时间的推移往往会以新功能化或亚功能化的形式发生分歧，最终的结局是形成新的基因或者部分拷贝假基因化[29]。本研究发现，无论是虹鳟还是斑马鱼，它们各自的2个hif-1a拷贝在组织表达谱上都发生了显著的分化，表明它们在不同组织中所发挥的生物学功能也已经发生了一定程度的分化。Hif-1a 通路所介导的低氧适应需要激活上百个基因的表达，这些基因所参与的功能包括血管生成、细胞分化/存活和糖代谢等，即通过增强循环系统来增加氧的供应或者增强无氧代谢的强度以获得足够的能量[4]。虹鳟或斑马鱼2个拷贝的共表达基因都能够富集到一定数量的与血管生成、细胞分化、物质运输及物质代谢相关的GO 条目，可见它们大体上仍保留了原始基因的核心功能。但2个拷贝的共表达基因并无交集，表明信号通路本身也已经发生了分化，2个拷贝分歧的形式主要是亚功能化。

虹鳟Hif-1aaa 和Hif-1aab 在bHLH 结构域上存在较大结构差异，同样支持它们发生功能分歧的结论。本研究还发现，虹鳟Hif-1aaa 和Hif-1aab 同样存在结构创新。就C-TAD结构域而言，小鼠和斑马鱼的序列建模都只能得到1 条螺旋，而虹鳟的2 条序列除了能得到比前2 个物种更长的螺旋外，还能得到1 条较短的螺旋。低氧状态下，C-TAD 结构域的主要功能是募集转录激活因子Ep300/Crebbp，而后者自身也拥有非常多的结合结构域[30]，C-TAD 结构域中额外的较短螺旋很可能有助于它们与Ep300/Crebbp 结合。虹鳟的Ep300/Crebbp 家族未发生过拷贝丢失，但这些拷贝之间必然也会发生功能分歧[31]，C-TAD 结构域的改变很可能就是为了适应不同Ep300/Crebbp 拷贝发生的结构变异[32]，从而实现更精细调控，其主要作用仍然是促进2 个拷贝的亚功能化。

综合来看，虹鳟的hif-1a基因家族在最近的全基因组重复事件之前丢失过拷贝，其祖先物种在某段时期内仅有1 个拷贝，这可能是它不耐低氧的原因之一。虹鳟和斑马鱼的2 个hif-1a拷贝之间较大的功能分歧表明，该基因家族的进化速度较快，未来的研究可以关注虹鳟不同地理种群中2 个hif-1a拷贝的遗传变异，特别是这些遗传变异所包含的结构与功能改变，从而找到有助于提高其耐低氧性能的基因型。