軍曹魚響應(yīng)低氧脅迫轉(zhuǎn)錄組SNP位點(diǎn)鑒定及其功能注釋分析

楊二軍，楊林桐，王維政，黃建盛,2,3*，張健東,2,3，王忠良,2,3，陳剛,2,3*

( 1. 廣東海洋大學(xué) 水產(chǎn)學(xué)院，廣東 湛江 524088；2. 南方海洋科學(xué)與工程廣東省實(shí)驗(yàn)室（湛江），廣東 湛江 524025；3. 廣東省水產(chǎn)經(jīng)濟(jì)動物病原生物學(xué)及流行病學(xué)重點(diǎn)實(shí)驗(yàn)室，廣東 湛江 524088)

1 引言

溶解氧（Dissolved Oxygen, DO）對水生動物的發(fā)育和健康至關(guān)重要[1]，同時也是水生動物生存不可缺少的環(huán)境因子之一。溶解氧不足時，機(jī)體的生長發(fā)育、呼吸、新陳代謝、繁殖及免疫功能等生化反應(yīng)和生理功能會受到影響，嚴(yán)重時會發(fā)生死亡[2]。隨著全球氣候變暖和水體富營養(yǎng)化引發(fā)的一系列生態(tài)問題，水體低氧化現(xiàn)象越來越嚴(yán)重。調(diào)查發(fā)現(xiàn)，我國近岸海域出現(xiàn)的低氧層（DO濃度約為3 mg/L）水域范圍越來越大，對生活在其中的海洋生物造成了極大的威脅，低氧現(xiàn)象已成為一個被高度重視和廣泛關(guān)注的問題[3-4]。

軍曹魚(Rachycentron canadum)，隸屬于鱸形目(Periformes)，軍曹魚科(Rachycentridae)，軍曹魚屬(Rachycentron)，又名海鱺、海龍魚等，是熱帶和亞熱帶海域的廣鹽性、肉食性洄游魚類。因其具有生長快、抗病能力強(qiáng)、營養(yǎng)豐富、肉質(zhì)鮮美等特點(diǎn)，使其成為我國南方沿海地區(qū)深水網(wǎng)箱養(yǎng)殖的重要品種之一[5]。近年來，隨著高密度集約化養(yǎng)殖模式的推廣，軍曹魚近海浮動式網(wǎng)箱養(yǎng)殖和深海網(wǎng)箱養(yǎng)殖發(fā)展迅速，然而受養(yǎng)殖密度大、管理不善、養(yǎng)殖生態(tài)環(huán)境惡化、天氣、溫度等因素的影響，軍曹魚在養(yǎng)殖過程中容易出現(xiàn)缺氧、抗病能力下降等情況，導(dǎo)致病害日趨頻繁，嚴(yán)重制約了軍曹魚的健康養(yǎng)殖發(fā)展[6-7]。針對軍曹魚養(yǎng)殖過程中出現(xiàn)的缺氧問題，本實(shí)驗(yàn)室前期開展了低氧脅迫對軍曹魚幼魚生長和血清生化指標(biāo)、氧化應(yīng)激、能量代謝以及腸道菌群影響的研究，結(jié)果表明，低氧脅迫對軍曹魚幼魚的生理狀況有顯著影響[7-12]。因此，開展軍曹魚的遺傳改良和優(yōu)良性狀品種培育的工作，對促進(jìn)我國軍曹魚養(yǎng)殖業(yè)健康快速的發(fā)展具有現(xiàn)實(shí)意義。

單核苷酸多態(tài)性（Single Nucleotide Polymorphism，SNP）是指因單個核苷酸變異引起的DNA序列多態(tài)性，其形式包括單個堿基的轉(zhuǎn)換、顛換、插入和缺失，是一種常見的基因突變[13]。SNP作為第三代分子標(biāo)記方法，是基因型和表型之間最好的研究載體，相比于前幾代的標(biāo)記方法，它具有位點(diǎn)豐富、分布廣泛、遺傳穩(wěn)定性強(qiáng)、檢測快速、自動化分析便利和成本低等特點(diǎn)，已成為最具前景的遺傳標(biāo)記方法之一，被廣泛應(yīng)用于生物、農(nóng)學(xué)、醫(yī)學(xué)等眾多領(lǐng)域[14-15]。近年來，SNP分子標(biāo)記技術(shù)在魚類[16-17]、甲殼類[18-19]和貝類[20-21]等水產(chǎn)動物的研究中也取得了重要進(jìn)展。

隨著高通量測序技術(shù)的快速發(fā)展以及測序成本的下降，利用轉(zhuǎn)錄組測序篩選目標(biāo)性狀SNP位點(diǎn)已成為研究遺傳多樣性、分子育種等分子標(biāo)記的常用手段[22]。研究發(fā)現(xiàn)，利用轉(zhuǎn)錄組測序技術(shù)和物種基因組信息，更易獲得與目標(biāo)性狀相關(guān)的SNP[23]。目前，通過轉(zhuǎn)錄組測序來篩選大量SNP位點(diǎn)的研究已在多種水產(chǎn)動物中被報道，如馬氏珠母貝(Pinctada fucata)[20]、方斑東風(fēng)螺(Babylonia areolata)[23]、羅氏沼蝦(Macrobrachium rosenbergii)[24]、團(tuán)頭魴(Megalobrama amblyocephala)[25]、鳙魚(Hypophthalmichthys nobilis)[26]等。關(guān)于軍曹魚群體遺傳的研究資料較為匱乏，曹丹煜[27]曾利用轉(zhuǎn)錄組測序結(jié)果篩選與軍曹魚鹽度適應(yīng)相關(guān)SNP位點(diǎn)，但未對其響應(yīng)低氧脅迫SNP位點(diǎn)挖掘和功能注釋進(jìn)行研究。本研究基于低氧脅迫后的軍曹魚幼魚腸道轉(zhuǎn)錄組測序結(jié)果，篩選出大量SNP位點(diǎn)，對SNP所在基因的SNP-Unigene進(jìn)行功能注釋，側(cè)重對免疫防御及低氧通路中的相關(guān)基因進(jìn)行SNP的挖掘與篩選，以期為今后軍曹魚耐低氧和抗病品系的遺傳選育、分子輔助育種研究提供基礎(chǔ)資料。

2 材料和方法

2.1 實(shí)驗(yàn)用魚及飼養(yǎng)

軍曹魚幼魚由廣東海洋大學(xué)湛江東海島生物研究基地提供。實(shí)驗(yàn)在廣東恒興飼料實(shí)業(yè)股份有限公司863基地進(jìn)行，養(yǎng)殖設(shè)施為室內(nèi)24 h持續(xù)充氣的流水養(yǎng)殖系統(tǒng)，水體為500 L。選取大小規(guī)格整齊、活力好、體格健壯的軍曹魚210尾（體重為（50.44±2.78）g），先在養(yǎng)殖水槽中暫養(yǎng)1周使其適應(yīng)新環(huán)境后再進(jìn)行低氧脅迫實(shí)驗(yàn)，每天投喂2次，及時用虹吸管清理糞便和剩余餌料，暫養(yǎng)期間水體溶解氧在6 mg/L以上，水溫為（29±1）℃，pH為7.8～8.0，鹽度為28～30，氨氮濃度小于0.02 mg/L，自然光照。

2.2 實(shí)驗(yàn)設(shè)計及采樣

將實(shí)驗(yàn)用魚分為對照組（CT組）和低氧組（HT組），每組設(shè)置3個重復(fù)實(shí)驗(yàn)，每個重復(fù)實(shí)驗(yàn)有35尾魚。HT組溶解氧濃度為（3.15±0.21）mg/L，通過覆蓋塑料薄膜、調(diào)整充氣量和水流速度的大小來降低水體溶解氧濃度。水體溶解氧采用碘量法(GB 7489-87)進(jìn)行測定，CT組無任何處理，溶解氧濃度為 （6.18±0.23）mg/L。實(shí)驗(yàn)處理周期為28 d。軍曹魚采樣前24 h禁止飼喂，采樣時先經(jīng)丁香酚麻醉，測定魚體的生物學(xué)指標(biāo)，然后置于冰上解剖，取出整個腸道并剔除脂肪、結(jié)締組織和內(nèi)容物，用無菌生理鹽水沖洗干凈后，放置于2.0 mL的凍存管中，投入液氮速凍，運(yùn)回實(shí)驗(yàn)室后轉(zhuǎn)入-80℃冰箱保存進(jìn)行后續(xù)實(shí)驗(yàn)。

2.3 轉(zhuǎn)錄組測序數(shù)據(jù)

實(shí)驗(yàn)樣品送至廣州基迪奧生物科技有限公司進(jìn)行轉(zhuǎn)錄組測序工作，使用TRIzol法提取總RNA，純度（NanoDrop 2000）和濃度（Agilent 2100）檢測合格后，每組建3個轉(zhuǎn)錄組文庫，每個文庫以2尾魚的腸道RNA等效混合進(jìn)行建庫，利用Illumina Hiseq-2000進(jìn)行測序[23]。轉(zhuǎn)錄組測序共獲得302 588 246條原始讀數(shù)（Raw Reads），經(jīng)去除含有接頭、重復(fù)及測序質(zhì)量較低的原始讀數(shù)后共獲得293 736 220條干凈讀數(shù)（Clean Reads）（數(shù)據(jù)未發(fā)表）。

2.4 SNP位點(diǎn)檢測及所在Unigene的注釋與功能分析

使用hisat2（http://ccb.jhu.edu/software/hisat2/index.shtml）全局比對方式，將不同長度的干凈讀數(shù)比對到本實(shí)驗(yàn)室軍曹魚參考基因組（NCBI Bioproject ID:PRJNA634421）獲得比對文件，使用samtools sort對其進(jìn)行排序[28]；將CT組和HT組中所有樣本排序后的比對文件和相應(yīng)的參考基因組序列輸入到bcftools（https://github.com/samtools/bcftools）軟件中，使用bcftools mpileup命令將兩組軍曹魚比對后的bam文件轉(zhuǎn)換成vcf文件，進(jìn)行SNP分型；使用bcftools filter對vcf文件進(jìn)行過濾，將質(zhì)量值小于20的低質(zhì)量突變位點(diǎn)過濾掉，保留質(zhì)量值不小于20的結(jié)果；使用annovar（http://www.openbioinformatics.org/annovar/）軟件對過濾后的vcf文件進(jìn)行注釋；使用Excel統(tǒng)計SNP位點(diǎn)的數(shù)量及分布情況[29-31]。再將包含SNP的基因序列比對到GO、KOG和KEGG數(shù)據(jù)庫進(jìn)行功能注釋。同時，從轉(zhuǎn)錄組數(shù)據(jù)中篩選免疫和低氧相關(guān)的KEGG通路，再由其通路篩選出差異表達(dá)的免疫防御和低氧相關(guān)基因，進(jìn)行SNP位點(diǎn)分析[23]。

3 結(jié)果與分析

3.1 SNP位點(diǎn)分析

3.1.1 SNP位點(diǎn)數(shù)量

利用軟件SOAPsnp將CT組和HT組的轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行分析，26 120條Unigene中共存在431 845個SNP位點(diǎn)，其中CT組有215 953個SNP位點(diǎn)，HT組有215 892個SNP位點(diǎn)。所有SNP位點(diǎn)中，純合SNP位點(diǎn)有152 642個（CT組76 550個，HT組76 092個）（表1）。SNP發(fā)生頻率約為1/171 bp（約171 bp有1個SNP位點(diǎn)）。SNP密度頻數(shù)分析顯示，每1 000 bp含有0～5個SNP位點(diǎn)的Unigene出現(xiàn)的頻率最高（圖1）。

表1 SNP位點(diǎn)數(shù)量概況Table 1 The number of the SNP sites

圖1 SNP密度頻數(shù)分布Fig. 1 The density frequency distribution of SNP

3.1.2 SNP的分布及轉(zhuǎn)換、顛換信息

對于Unigene中SNP位點(diǎn)數(shù)目統(tǒng)計得出，CT組和HT組中分別含有1個SNP位點(diǎn)的Unigene有4 337和4 384條，含有2～10個SNP位點(diǎn)的Unigene有8 299和8 258條；含有10個以上SNP位點(diǎn)的Unigene有400條和402條（圖2）。

圖2 SNP的分布統(tǒng)計Fig. 2 The statistics distribution of SNP

對軍曹魚SNP突變類型的數(shù)量分析表明，在純合SNP位點(diǎn)中，轉(zhuǎn)換數(shù)量明顯高于顛換數(shù)量，即CT組SNP位點(diǎn)中，轉(zhuǎn)換類型131 054個（60.69%），顛換類型84 899個（39.31%）；HT組SNP位點(diǎn)中，轉(zhuǎn)換類型131 587個（60.95%），顛換類型84 305個（39.05%）。在核苷酸的變異類型中，C-T發(fā)生轉(zhuǎn)換的頻率最高，分別占純合SNP數(shù)的44.60%（CT組34 145個）和45.21%（HT組34 400個）；T-A發(fā)生顛換的頻率最高，分別占純合SNP數(shù)的16.45%（CT組12 591個）和16.37%（HT組12 458個）（圖3）。

圖3 SNP突變類型數(shù)量統(tǒng)計Fig. 3 The numbers of different mutation types statistics of SNP

3.1.3 SNP位置分析及突變類型分析

對軍曹魚腸道轉(zhuǎn)錄組SNP序列位置變異類型分析結(jié)果顯示（圖4），CT組和HT組中SNP位點(diǎn)總數(shù)分別為215 953和215 892個，位于內(nèi)含子和基因間區(qū)的SNP位點(diǎn)較多。在CT組和HT組中，位于內(nèi)含子的SNP位點(diǎn)分別占45.13%和44.62%；位于基因間區(qū)的SNP位點(diǎn)分別占24.66%和25.17%；位于外顯子的SNP位點(diǎn)分別占20.53%和20.43%；位于基因下游的SNP位點(diǎn)分別占13.54%和13.71%；位于3′非編碼區(qū)的SNP位點(diǎn)分別占8.74%和8.76%；位于基因上游的SNP位點(diǎn)分別占4.45%和4.44%；位于5′非編碼區(qū)的SNP位點(diǎn)分別占2.20%和2.15%；位于基因上游和下游的SNP位點(diǎn)分別為2 031個和2 028個; 位于剪切位點(diǎn)的SNP位點(diǎn)分別為523個和518個。

圖4 SNP位置鑒定Fig. 4 Identification of SNP location

對SNP的突變類型分析結(jié)果顯示（圖5），突變的SNP位點(diǎn)數(shù)量較少，其中同義突變SNP位點(diǎn)數(shù)量占比最高，CT組為10.10%，HT組為10.04%；非同義突變SNP位點(diǎn)數(shù)量在CT組中占比為9.15%，HT組占比為9.11%。其他少量突變類型包括移碼缺失、無義突變和非移碼缺失等。

圖5 SNP突變類型鑒定Fig. 5 Identification of SNP variants type

3.2 SNP-Unigene的注釋與功能

將26 120條SNP-Unigene與GO數(shù)據(jù)庫比對，結(jié)果顯示，分布于生物學(xué)過程（Biological Process，BP）、分子功能（Molecular Function，MF）和細(xì)胞組分（Cellular Component，CC）3大功能分類的GO條目分別為8 043條、7 537條和7 854條，其中參與生物學(xué)過程的SNP-Unigene數(shù)量最多（圖6）。將SNP-Unigene序列比對到KOG數(shù)據(jù)庫進(jìn)行注釋和分類，結(jié)果顯示，共有14 334條SNP-Unigene匹配到相應(yīng)的注釋信息，根據(jù)功能分為25類，其中以“信號轉(zhuǎn)導(dǎo)機(jī)制”和“一般功能預(yù)測”類較多，分別包含2 822條和2 548條SNPUnigene（圖7）。將26 120條SNP-Unigene序列進(jìn)行KEGG富集分析，發(fā)現(xiàn)有5 160條SNP-Unigene富集到349條KEGG代謝通路中，其中以“信號轉(zhuǎn)導(dǎo)” “全球概覽地圖”和“傳染病”通路中富集的SNP-Unigene較多，分別為1 541條、1 123條和1 153條（圖8）。

圖6 SNP-Unigene GO功能注釋分析Fig. 6 GO annotation analysis of SNP-Unigene

圖7 SNP-Unigene KOG功能注釋分析Fig. 7 KOG function annotation analysis of SNP-Unigene

圖8 SNP-Unigene的KEGG功能注釋分析Fig. 8 KEGG classification of SNP-Unigene

3.3 免疫防御相關(guān)SNP-Unigene的富集及特異SNP位點(diǎn)

根據(jù)KEGG信號通路的富集分析，發(fā)現(xiàn)共有3 417條免疫防御相關(guān)的SNP-Unigene被注釋到“MAPK信號通路”等35條與免疫相關(guān)的通路中，其中以“MAPK信號通路”中注釋的SNP-Unigene最多（262條），其次分別為“神經(jīng)活性配體-受體相互作用”（228條）、“Rap1信號通路”（190條）、“Ras信號通路”（164條）等信號通路（表2）。

表2 免疫防御相關(guān)SNP-Unigene的KEGG富集分析Table 2 KEGG enrichment analysis of immune-related SNP-Unigene

3.4 差異表達(dá)免疫防御及低氧相關(guān)基因SNP位點(diǎn)分析

基于軍曹魚低氧脅迫后腸道轉(zhuǎn)錄組分析，得到大量的差異基因（篩選條件為錯誤發(fā)生率（False Discovery Rate，F(xiàn)DR）小于0.05且|log2FC|大于1，其中FC為差異倍數(shù)）。對“MAPK信號通路”等7個免疫通路及“HIF-1信號通路”中差異基因進(jìn)行了SNP位點(diǎn)分析（表3）。

表3 CT和HT轉(zhuǎn)錄組中差異表達(dá)免疫防御及低氧相關(guān)基因SNP位點(diǎn)分析Table 3 SNP identified in differentially expressed immune and hypoxia-related genes in CT and HT transcriptomes

4 討論

4.1 低氧脅迫下軍曹魚轉(zhuǎn)錄組數(shù)據(jù)中SNP位點(diǎn)鑒定及功能注釋

SNP作為一種分子遺傳標(biāo)記，在水生動物的遺傳選育中發(fā)揮著重要作用[25]。本研究從CT轉(zhuǎn)錄組、HT轉(zhuǎn)錄組數(shù)據(jù)中分別獲得大量的SNP位點(diǎn)，SNP發(fā)生頻率依次為0.00 586和0.005 858（1/171 bp），高于鹽度脅迫對軍曹魚（1/730 bp）[27]、半滑舌鰨（Cynoglossus semilaevis）（1/491 bp）[32]、大 黃 魚（Larimichthy crocea）（1/506 bp）[33]、羅氏沼蝦（Macrobrachium rosenbergii）（1/1 152 bp）[24]、大西洋鮭（Salmo salar）（1/641 bp）[34]SNP發(fā)生頻率，但低于馬氏珠母貝（1/100 bp）[20]。有研究表明，不同物種中SNP分布頻率差異很大，可能與研究物種、遺傳背景和生存環(huán)境有關(guān)[19]，也可能是因?yàn)樵诒磉_(dá)序列標(biāo)簽中篩選SNP的要求和條件不同[21]。SNP堿基替換是產(chǎn)生基因突變和推動進(jìn)化的原始動力[35]。堿基替換類型分為轉(zhuǎn)換(A-G和T-C的相互置換)和顛換(A-T、C-G、A-C、T-G的相互置換)。理論上，生物中發(fā)生顛換與轉(zhuǎn)換比例為1∶2[23]。本研究中，SNP轉(zhuǎn)換類型的比例約為60%，顛換類型的比例約為40%，轉(zhuǎn)換大于顛換，符合“轉(zhuǎn)換偏差”原理。在核苷酸變異類型中，C-T發(fā)生轉(zhuǎn)換的比例最高，與大多數(shù)水產(chǎn)動物的研究結(jié)果相一致[19,36]，原因可能是CpG二核苷酸上甲基化的胞嘧啶殘基容易脫去氨基形成胸腺嘧啶[37]，或受物種進(jìn)化過程中承受的壓力影響[38]。

根據(jù)SNP在基因組的分布位置可分為3類：基因編碼區(qū)SNPs（cSNPs）、基因間SNPs（iSNPs）和基因周邊SNPs（pSNPs）[39]。其中cSNP和pSNP會對基因的表達(dá)產(chǎn)生影響[22]。根據(jù)對生物遺傳性狀的影響，cSNPs又分為同義突變和非同義突變，前者編碼序列的變化不會影響翻譯后的氨基酸序列，故不會導(dǎo)致蛋白質(zhì)的功能變化；后者編碼序列的變化會引起翻譯后氨基酸序列的變化，改變蛋白質(zhì)的生物學(xué)功能[40]。由于受到較強(qiáng)的選擇壓力，非同義SNP位點(diǎn)的突變經(jīng)常會影響生物的表型，而關(guān)于非同義SNP位點(diǎn)的研究已在團(tuán)頭魴[25]、大黃魚[41]有相關(guān)報道。對軍曹魚腸道轉(zhuǎn)錄組序列SNP位置和突變類型分析發(fā)現(xiàn)，位于基因編碼區(qū)的SNP位點(diǎn)以及非同義突變SNP位點(diǎn)占比較多。因此，研究基因編碼區(qū)非同義SNP可能對軍曹魚具有重要的生物學(xué)意義。此外，GO、KOG及KEGG功能注釋結(jié)果表明，大部分SNP-Unigene涉及信號轉(zhuǎn)導(dǎo)、傳染病、代謝、癌癥、內(nèi)分泌系統(tǒng)和免疫系統(tǒng)等多種途徑，表明軍曹魚幼魚遭受低氧脅迫后，這些SNP-Unigene參與各種生命活動，并影響機(jī)體的各種生物學(xué)性狀。

4.2 低氧脅迫下軍曹魚免疫相關(guān)基因SNP位點(diǎn)分析

免疫相關(guān)基因的SNP在水產(chǎn)遺傳育種中的應(yīng)用對魚類抗逆性性狀的選育具有重要意義[41]。在本研究中，為了解軍曹魚免疫防御機(jī)制，對轉(zhuǎn)錄組測序數(shù)據(jù)進(jìn)一步分析，結(jié)果發(fā)現(xiàn)，3 417條SNP-Unigene被注釋到35條與免疫防御相關(guān)的通路中，大部分SNP-

Unigene參與“MAPK信號通路” “神經(jīng)活性配體-受體相互作用” “Rap1信號通路” “Ras信號通路”等免疫通路。我們從相關(guān)免疫通路中篩選出了大量差異表達(dá)的免疫相關(guān)基因及其SNP位點(diǎn)。因此，這里對部分關(guān)鍵基因進(jìn)行后續(xù)的研究分析。

IL1R1是細(xì)胞表面重要的免疫基因，是白細(xì)胞介素-1β（IL-1β）的受體，在機(jī)體炎癥反應(yīng)和免疫信號轉(zhuǎn)導(dǎo)過程中非常重要[42]。Lam等[43]研究報道，慢性缺氧時，大鼠頸動脈體細(xì)胞中促炎細(xì)胞因子及其受體（IL-1β和IL1R1）的表達(dá)量顯著增加，表明其在頸動脈體局部炎癥中起著重要功能。在本研究中，軍曹魚受低氧脅迫后，腸道IL1R1基因的表達(dá)量顯著上調(diào)，提示IL1R1基因可能在軍曹魚免疫調(diào)節(jié)和抵御低氧應(yīng)激中發(fā)揮著重要作用。此外，在基因多態(tài)性方面，已發(fā)現(xiàn)IL1R1基因的SNP與多種人類疾病相關(guān)[44]。胡亮[45]對崗巴綿羊的候選特異性SNP位點(diǎn)進(jìn)行基因注釋和富集分析，發(fā)現(xiàn)IL1R1和IL1R2基因參與炎癥反應(yīng)及自身免疫性調(diào)節(jié)。因此，暗示IL1R1基因可作為軍曹魚免疫防御的候選基因，這為今后進(jìn)行IL1R1基因的多態(tài)性與疾病關(guān)聯(lián)研究提供了參考。

Wnt信號通路是一個高度保守的系統(tǒng)，通過Wnt蛋白與細(xì)胞表面受體相互作用，調(diào)控著細(xì)胞增殖、分化和黏附等多種復(fù)雜的生物學(xué)過程[23,46]。Wnt配體蛋白與跨膜受體卷曲蛋白(FZD)和低密度脂蛋白受體相關(guān)蛋白(LRP5/6)形成的細(xì)胞表面受體復(fù)合物結(jié)合，從而激活細(xì)胞內(nèi)的Wnt信號通路[47]。在本研究中，軍曹魚受低氧脅迫后，Wnt信號通路中FZD3和LRP5基因表達(dá)量顯著下調(diào)，且存在多個SNP位點(diǎn)。目前有關(guān)魚類中FZDs家族基因的研究主要涉及基因克隆、組織表達(dá)及在魚類發(fā)育過程中的生理作用[48]。張麗晗等[48]研究了銅脅迫下黃顙魚（Pelteobagrus fulvidraco）FZD家 族 基 因（FZD2、FZD3A、FZD4和FZD10）的mRNA水平。FZD3在小鼠和人的毛囊發(fā)育、細(xì)胞遷移等方面發(fā)揮著關(guān)鍵作用[49]，已有研究者對該基因多態(tài)性位點(diǎn)進(jìn)行分析，篩選影響性狀的SNP位點(diǎn)。Zhao等[49]對中國美利奴羊FZD3基因的多態(tài)性進(jìn)行檢測，發(fā)現(xiàn)SNP1和SNP2與羊毛性狀有顯著的相關(guān)性，可作為綿羊育種的潛在遺傳標(biāo)記。LRP5基因在Wnt信號通路中扮演著重要角色，其突變與骨質(zhì)疏松、癌癥等疾病的發(fā)生密切相關(guān)[50]。已有研究表明[50]，LRP5基因多態(tài)性與2型糖尿病患者的β細(xì)胞功能和脂質(zhì)代謝有關(guān)。由此推測，F(xiàn)ZD3和LRP5基因可能參與調(diào)控軍曹魚低氧應(yīng)激過程中的免疫反應(yīng)，具體機(jī)制有待進(jìn)一步探討，而研究這些免疫基因的SNP對低氧脅迫下軍曹魚免疫防御相關(guān)SNP位點(diǎn)的開發(fā)也具有潛在意義。

本研究結(jié)果表明，參與mTOR信號通路的差異基因LIPIN1表達(dá)水平上調(diào)，也存在多個SNP位點(diǎn)。Lipin1是由LIPIN1基因編碼的磷脂酸磷脂酶，具有雙向調(diào)控脂肪代謝的功能，既能催化磷脂酸去磷酸化得到二脂酰甘油，從而促進(jìn)甘油三酯和磷脂的合成；還具有轉(zhuǎn)錄激活功能，參與調(diào)控脂肪代謝相關(guān)基因表達(dá)[51-52]。關(guān)于LIPIN蛋白家族的研究主要集中在蛋白功能和基因表達(dá)方面，多態(tài)性方面以小鼠、豬、牛等哺乳動物為主，水產(chǎn)經(jīng)濟(jì)動物較為匱乏。有報道稱[53]，LIPIN1基因受HIF-1調(diào)控，其表達(dá)水平上調(diào)與非脂肪細(xì)胞中的甘油三酯和脂滴的積累有關(guān)，提示LIPIN1是細(xì)胞適應(yīng)低氧應(yīng)激的重要中介物。He等[54]在豬的LIPIN1基因的編碼區(qū)和UTR檢測到5個SNP位點(diǎn)與豬的脂肪沉積性狀存在顯著相關(guān)。這些研究表明，LIPIN1基因可能對軍曹魚低氧脅迫過程中脂質(zhì)代謝有一定程度的影響，可為LIPIN1基因多態(tài)性與目標(biāo)性狀關(guān)聯(lián)研究提供了新思路。

4.3 低氧脅迫下軍曹魚HIF-1信號通路中差異基因SNP位點(diǎn)分析

除了免疫防御相關(guān)信號通路所涉及的差異基因外，還有多種基因是通過其他信號通路來共同應(yīng)答魚類的低氧脅迫。HIF-1信號通路是魚類遭受低氧應(yīng)激時的關(guān)鍵通路之一，通路中的相關(guān)基因不僅參與機(jī)體代謝的調(diào)控，而且還會通過抑制血紅細(xì)胞的增殖和凋亡，提高血氧親和力，以此來適應(yīng)低氧應(yīng)激[55]。在本研究結(jié)果中發(fā)現(xiàn)，HIF-1信號通路中存在PIK3CA、HK1、ANGPT1、EPAS1、GAPDH、ANGPT2、EDN1和ENO3等8個差異基因也同樣含有多個SNP位點(diǎn)。

PIK3CA基因是編碼I類磷脂酰肌醇3激酶（PI3Ks）的催化亞單位，其突變能激活PI3K/AKT信號通路，進(jìn)而參與細(xì)胞增殖、凋亡和蛋白質(zhì)合成等多種生理過程[56-57]。Zhang等[56]研究發(fā)現(xiàn)，SmPIK3CA基因參與大菱鲆的免疫反應(yīng)，并在該基因中篩選出4個抗病相關(guān)的SNPs。EPAS1作為激活低氧調(diào)控基因表達(dá)的關(guān)鍵轉(zhuǎn)錄因子[58]，在細(xì)胞代謝和氧化還原反應(yīng)中起著重要作用。關(guān)于EPAS1基因SNP位點(diǎn)的研究多集中在藏族人和高原哺乳動物，在魚類中的研究較少。在藏族人和高原動物低氧適應(yīng)研究中發(fā)現(xiàn)，EPAS1基因都受到選擇且信號顯著。在牦牛中發(fā)現(xiàn)EPAS1基因多態(tài)性與血紅蛋白(HGB)關(guān)系緊密[58]。在藏綿羊EPAS1基因與低氧適應(yīng)相關(guān)性研究中發(fā)現(xiàn)，該基因第8內(nèi)含子C871T 變異位點(diǎn)的等位基因(T)突變頻率與海拔高度呈現(xiàn)正相關(guān)，且TT基因型能顯著增加血紅蛋白濃度[59]。由此推測，HIF-1信號通路中篩選出來的PIK3CA等8個顯著差異基因可能對軍曹魚低氧相關(guān)SNP位點(diǎn)的挖掘具有重要的參考價值。

綜上所述，本研究基于低氧脅迫下軍曹魚腸道轉(zhuǎn)錄組測序數(shù)據(jù)，獲得了大量與免疫、低氧相關(guān)差異基因的SNP位點(diǎn)，豐富了軍曹魚分子標(biāo)記的數(shù)據(jù)，這些候選標(biāo)記將為今后深入了解軍曹魚響應(yīng)低氧脅迫的分子機(jī)制、SNP分子標(biāo)記開發(fā)和遺傳連鎖圖譜構(gòu)建等研究提供重要的理論依據(jù)。