摘 要 為探究低、中、高3個不同海拔梯度下喜馬拉雅旱獺（Marmota himalayana）基因表達是否存在差異，篩選該物種與高原低氧及強紫外輻射環(huán)境適應(yīng)的相關(guān)候選基因，采用Illumina HiSeq 2500高通量測序平臺對喜馬拉雅旱獺肝臟組織進行無參轉(zhuǎn)錄組學(xué)分析，通過對質(zhì)控后測序數(shù)據(jù)與參考基因組的比較分析，篩選差異表達基因后進行GO功能注釋和KEGG富集分析。結(jié)果顯示：樂都樣品（低海拔組，LD）、玉樹樣品（高海拔組，YS）喜馬拉雅旱獺肝臟組織中差異表達基因差異較大，剛察樣品（中海拔組，GC）與其他兩組間均具有較小的基因表達差異模式。上調(diào)表達基因在19個GO Terms和6個KEGG通路上有顯著富集，氧化磷酸化相關(guān)基因（HiF-1α、MPKA3）、礦物質(zhì)吸收及脂質(zhì)代謝相關(guān)基因（PPARs）在高海拔組中呈現(xiàn)上調(diào)表達，推測這些基因可能直接或間接在喜馬拉雅旱獺適應(yīng)高海拔低氧等極端環(huán)境中，通過調(diào)整細胞代謝、提高氧化還原能力以及增強免疫反應(yīng)來維持自身的生理功能和生存能力，為后續(xù)深入研究高原野生動物基因多樣性以及喜馬拉雅旱獺低氧適應(yīng)性提供基因組學(xué)數(shù)據(jù)。

關(guān)鍵詞：喜馬拉雅旱獺；低氧適應(yīng)性；轉(zhuǎn)錄組學(xué)；肝臟

中圖分類號：Q953 文獻標識碼：A 文章編號：2310 - 1490（2024）- 02 - 0305 - 09

DOI：10.12375/ysdwxb.20240209

喜馬拉雅旱獺（Marmota himalayana）隸屬于哺乳綱（Mammalia）嚙齒目（Rodentia）松鼠科（Sciuri?dae）旱獺屬（Marmota）。目前，世界范圍內(nèi)有2個亞種，分別是指名亞種（M. h. himalayana）和川西亞種（M. h. rubustus），指名亞種主要分布在我國甘肅西部、青海（柴達木盆地）及新疆等地，川西亞種主要分布在我國四川西部和西北部等地［1］。作為青藏高原特有野生動物資源，其在生物醫(yī)學(xué)領(lǐng)域中具有廣泛應(yīng)用［2］。

2012年，王忠東［3］首次通過旱獺活體肝臟穿刺術(shù)，成功建立實驗技術(shù)與室內(nèi)飼養(yǎng)管理技術(shù)，隨后，劉海青等［4］成功采用人工批量化繁育技術(shù)并建立人工繁育種群，為后續(xù)喜馬拉雅旱獺實驗動物化研究的開展提供基礎(chǔ)技術(shù)支持。研究人員通過觀測并分析該物種骨骼?。?］與心?。?］等組織學(xué)相關(guān)指標，發(fā)現(xiàn)經(jīng)過長期的自然選擇和進化，喜馬拉雅旱獺在高原環(huán)境中表現(xiàn)出良好的適應(yīng)性，獲得適應(yīng)低氧環(huán)境且穩(wěn)定遺傳的生物學(xué)特性。隨著分子生物學(xué)技術(shù)的發(fā)展，李優(yōu)等［7］利用dd-RAD基因組學(xué)測序技術(shù)對青海省44例喜馬拉雅旱獺進行遺傳進化分析，揭示不同地區(qū)喜馬拉雅旱獺種群之間存在遺傳差異。盡管線粒體基因組在研究物種進化和親緣關(guān)系等方面發(fā)揮了重要作用，但無法提供整體基因組完整信息［8］。

現(xiàn)有喜馬拉雅旱獺相關(guān)研究大多集中于傳染病學(xué)及生物學(xué)等方面，包括病原體鑒定［9］、疾病傳播途徑［10］及遺傳進化差異［7］等。隨著深度測序的市場化，RNA測序（RNA-Seq）利用高通量測序平臺更精準全面地探索特定生物在特定條件下的基因表達和生物代謝，具有低成本、全面及高準確度的優(yōu)勢，現(xiàn)已成為探究物種進化調(diào)控機制的重要技術(shù)手段［11］，然而目前對不同海拔下喜馬拉雅旱獺適應(yīng)低氧的調(diào)控機制及整體性分析研究尚不深入。本研究基于轉(zhuǎn)錄組學(xué)測序技術(shù)，通過對3個不同海拔梯度下的喜馬拉雅旱獺肝臟組織進行無參轉(zhuǎn)錄組學(xué)分析，初步篩選出喜馬拉雅旱獺與低氧適應(yīng)有關(guān)的候選基因，為進一步挖掘高原土著動物低氧適應(yīng)機制及相關(guān)功能基因奠定基礎(chǔ)。

1 材料與方法

1. 1 樣本采集

在海東市樂都區(qū)引勝鄉(xiāng)［海拔1 992 m；36°35' N，102°24' E；3例；記為低海拔組（LD組）］、海北藏族自治州剛察縣［海拔3 127 m；35°29' N，102°16' E；3例；記為中海拔組（GC組）］和玉樹州曲麻萊縣麻多鄉(xiāng)［海拔4 101 m；35°1' N，96°22' E；3例；記為高海拔組（YS組）］3個地區(qū)采集喜馬拉雅旱獺新鮮肝臟組織，共9例（圖1），分別裝入2 mL凍存管內(nèi)（n=3），液氮速凍后按照分組標記，置于冰箱中-80 ℃保存?zhèn)溆?。取樣過程與青海省地方病預(yù)防控制所合作，遵循動物倫理相關(guān)法規(guī)（PJ20230306）。

1. 2 RNA 提取與質(zhì)控

使用TRIzol法提取喜馬拉雅旱獺肝臟組織中的總RNA，1%瓊脂糖凝膠電泳初步分析RNA降解及污染程度，采用Agilent 2100生物分析儀與NanoDrop設(shè)備檢測RNA 的完整性（RNA integrity number，RIN）及純度。建庫要求為動物組織樣本不得少于250 mg，OD260/280為1. 8～2. 2，RIN值≥7，28S∶18S≥7。

1. 3 文庫構(gòu)建與基礎(chǔ)分析

利用Oligo（dT）磁珠富集法富集mRNA，使用Ragmentation Buffer 隨機打斷mRNA，通過AMPureXP beads 純化PCR 產(chǎn)物后，構(gòu)建cDNA 文庫。采用Illumina HiSeq 2500高通量測序平臺同時進行合成與測序。測序原始數(shù)據(jù)（raw reads）經(jīng)過濾后得到優(yōu)質(zhì)測序數(shù)據(jù)（clean reads）。使用Trinity軟件對cleanreads 進行拼接，得到Unigenes 轉(zhuǎn)錄本，獲取后續(xù)分析參考序列。

1. 4 差異表達基因分析

采用DESeq2（1. 6. 3）軟件篩選差異表達基因（differentially expressed genes，DEGs），顯著差異表達閾值為plt;0. 05且 | log2FoldChange | 值大于1。為分析差異表達基因的功能特性和參與的生理生化代謝途徑，對DEGs進行GO功能注釋和KEGG富集分析，展現(xiàn)3組喜馬拉雅旱獺種群間差異表達基因序列在兩數(shù)據(jù)庫中的分布狀況，在DEGs中顯著富集的閾值（q）不超過0. 05。所有統(tǒng)計分析均用SPSS 18. 0軟件進行，使用兩兩獨立樣本t 檢驗進行比較，同時用GraphPad Prism 8. 0 軟件繪圖，plt;0. 05，差異具有統(tǒng)計學(xué)意義。

2 結(jié)果與分析

2. 1 總RNA 質(zhì)檢結(jié)果

從LD組、GC組及YS組喜馬拉雅旱獺肝臟組織中提取的總RNA無DNA等雜質(zhì)污染，純度（OD260/280為2. 02～2. 08，OD260/230 為2. 33～2. 35）、濃度（545. 9～584. 1 ng/μL）和檢測的RIN 值（8. 8～9. 6）滿足建庫及測序要求。

2. 2 功能注釋分析

質(zhì)控后轉(zhuǎn)錄組分析數(shù)據(jù)總量為20. 38 Gb，Q20≥96. 45%，GC 比例介于47. 63%～48. 20%，堿基錯誤率低，數(shù)據(jù)質(zhì)量合格，表明測序數(shù)據(jù)整體質(zhì)量較好，可用于后續(xù)分析（表1）。將質(zhì)控后數(shù)據(jù)與NT、NR、KO、PFAM、KOG和Swiss Prot 六個數(shù)據(jù)庫比對，統(tǒng)計各數(shù)據(jù)庫占比情況并可視化處理（圖2），注釋結(jié)果表明：NT數(shù)據(jù)庫注釋信息較為豐富，轉(zhuǎn)錄本平均占比最高，KOG數(shù)據(jù)庫信息匹配率較低，轉(zhuǎn)錄本平均占比最低。

2. 3 差異表達基因的篩選及聚類分析

對LD、GC、YS兩兩組間交集差異表達基因關(guān)聯(lián)后，將信息整合繪制Venn圖（圖3），其中3組共表達基因有991個。對DEGs進行表達模式聚類分析（圖4），發(fā)現(xiàn)LD、YS組喜馬拉雅旱獺的表達模式區(qū)別顯著，推測是相鄰海拔梯度間環(huán)境氧含量變化程度較為溫和所致。

2. 4 差異表達基因的GO 富集分析

將喜馬拉雅旱獺差異表達基因與GO數(shù)據(jù)庫進行比對，差異基因注釋結(jié)果被分為生物過程（BP）、細胞組分（CC）及分子功能（MF）三大功能類別，主要顯著富集在外切脫氧核酸酶Ⅶ復(fù)合物部分DNA降解、DNA修飾、芳香氨基酸家族代謝和賴氨酸代謝、脫氧核糖核酸酶以及結(jié)合血紅素部分活性等生物學(xué)過程（圖5）。

上調(diào)差異表達基因主要富集在線粒體電子傳遞（GO：0006120）、礦物質(zhì)吸收（GO：0055072）以及DNA降解（GO：0006308）等生物學(xué)過程，表明機體通過上調(diào)線粒體電子傳遞鏈的活性，增加線粒體數(shù)量來維持低氧環(huán)境下HIF-1α 基因的穩(wěn)定，以優(yōu)化氧氣的吸收和運輸能力，從而更有效地利用有限的氧氣資源，保持正常的能量代謝。下調(diào)基因顯著富集在分子功能（MF）分類中，主要是DNA依賴的轉(zhuǎn)錄調(diào)控（GO：0006355）、小GTP酶介導(dǎo)的信號轉(zhuǎn)導(dǎo)（GO：0007264）、DNA定向DNA聚合酶活性（GO：0003887）以及DNA結(jié)合（GO：0003677）等方面，說明可能在低氧低壓條件下機體通過抑制DNA損傷過程，以及下調(diào)小GTP活性狀態(tài)降低高原環(huán)境中喜馬拉雅旱獺的一些基因表達以及生理代謝過程（表2）。

2. 5 差異表達基因的KEGG 富集分析

為探究差異基因參與影響的信號轉(zhuǎn)導(dǎo)通路或生物學(xué)功能，進行KEGG富集分析后發(fā)現(xiàn)這些基因顯著富集在乙酰乙酸和二羧酸代謝通路，纈氨酸、亮氨酸、異亮氨酸的降解代謝，胰島素信號通路，真核生物中核糖體生物合成等通路（圖6）。

3 討論

青藏高原的低氧環(huán)境對本土動物施加選擇壓力，使其在分子水平上表現(xiàn)出趨同進化等遺傳特性，在長期低氧及強紫外輻射等極端氣候環(huán)境中，高原土著動物形成穩(wěn)定可遺傳的適應(yīng)性機制［12］。KEGG通路富集分析表明，不同海拔梯度下喜馬拉雅旱獺差異表達基因顯著富集在氧化磷酸化（ko00190）、細胞周期（ko04110）、鐵調(diào)控信號通路（ko04927）、PI3K-Akt 信號通路（ko04151）、脂質(zhì)代謝（ko00061）以及細胞凋亡（ko04210）等通路（表3）。

線粒體作為“能量工廠”為機體提供ATP，支撐高原哺乳動物的正常生命活動，能量產(chǎn)生方式主要依賴于氧化磷酸化及三羧酸循環(huán)［13］。文禹粱［14］在藏綿羊?qū)嶒炛邪l(fā)現(xiàn)，氧化磷酸化的上調(diào)表達可維持藏綿羊日常生理需求使機體表現(xiàn)出低氧適應(yīng)性。HIF-1α 作為缺氧時下游基因的主要激活因子，通過調(diào)控多個靶基因的轉(zhuǎn)錄活性并上調(diào)表達增加機體對ATP的依賴性，形成低氧適應(yīng)性機制［15］，在藏綿羊［16］、藏雞［17］等高原土著動物的研究中該基因均受到正選擇影響。此外，研究發(fā)現(xiàn)MAPK信號通路和PI3K-Akt信號轉(zhuǎn)導(dǎo)通路等與低氧適應(yīng)有關(guān)聯(lián)的通路顯著富集。細胞通過MAPK信號轉(zhuǎn)導(dǎo)通路響應(yīng)外部刺激，參與細胞增殖、分化及凋亡等多種生物學(xué)反應(yīng)［18］。PI3K-Akt信號通路是調(diào)節(jié)細胞增殖、分化和凋亡的重要細胞內(nèi)信號通路，具有在缺氧條件下保護機體免受細胞凋亡影響的作用［19］。張倩等［20］在牦牛大腦皮質(zhì)組織轉(zhuǎn)錄學(xué)研究中證明，MAPK 信號通路和PI3K-Akt信號轉(zhuǎn)導(dǎo)通路促進低氧適應(yīng)機制的形成。青藏地區(qū)氧氣稀薄且海拔較高，大氣層厚度變薄，對入射輻射的吸收和散射能力減弱［21］，Song等［22］研究表明，高原土著生物通過修復(fù)基因的高表達，消除紫外誘導(dǎo)下DNA損傷等不良反應(yīng)，獲得強紫外輻射適應(yīng)性優(yōu)勢。

綜上，本研究利用RNA-Seq技術(shù)對LD、GC和YS三組不同海拔梯度下喜馬拉雅旱獺肝臟組織進行轉(zhuǎn)錄組學(xué)分析，檢測出大量與能量代謝及細胞維穩(wěn)相關(guān)的DEGs，KEGG富集分析結(jié)果表明這些基因主要富集在氧化磷酸化、DNA損傷感應(yīng)與信號轉(zhuǎn)導(dǎo)通路及PI3K-Akt等途徑，GO注釋分析中也發(fā)現(xiàn)了相似的結(jié)果，這可能提示成年喜馬拉雅旱獺通過增強氧化應(yīng)激、維護細胞基因組穩(wěn)定性等生理機制維持自身穩(wěn)定，獲得高原環(huán)境適應(yīng)優(yōu)勢。以上結(jié)果有助于更好地理解喜馬拉雅旱獺在高海拔環(huán)境中的適應(yīng)生存策略，進一步對研究高原生物適應(yīng)性機制提供更多的參考，為高原野生動物基因多樣性的保護和管理奠定基礎(chǔ)。

參考文獻：

［1］ 魏輔文， 楊奇森， 吳毅， 等. 中國獸類分類與分布［M］. 北京：科學(xué)出版社， 2022： 236.

WEI F W， YANG Q S， WU Y， et al. Taxonomy and distributionof mammals in China［M］. Beijing： Science Press， 2022： 236.

［2］ 劉海青， 范微， 張靜宵， 等. 喜馬拉雅旱獺實驗動物化的研究進展［J］. 中國比較醫(yī)學(xué)雜志， 2015， 25（11）： 64-68.

LIU H Q， FAN W， ZHANG J X， et al. Research progress on do?mestication of Marmota himalayana［J］. Chinese Journal of Com?parative Medicine， 2015， 25（11）： 64-68.

［3］ 王忠東. 喜馬拉雅旱獺動物實驗技術(shù)與飼養(yǎng)管理［J］. 醫(yī)學(xué)動物防制， 2012， 28（6）： 653-656.

WANG Z D. Feeding management and animal experiments of Mar?mota himalayana［J］. Journal of Medical Pest Control， 2012， 28（6）： 653-656.

［4］ 劉海青， 范微， 張靜宵， 等. 喜馬拉雅旱獺的人工繁育技術(shù)［J］. 中國比較醫(yī)學(xué)雜志， 2015， 25（11）： 59-63；9.

LIU H Q， FAN W， ZHANG J X， et al. Scale breeding and the es?tablishment of breeding population of Marmota himalayana［J］.Chinese Journal of Comparative Medicine， 2015， 25（11）：59-63；9.

［5］ 秦鴻楠， 周娟， 李琳， 等. 喜馬拉雅旱獺骨骼肌低氧適應(yīng)的組織學(xué)特點［J］. 中國獸醫(yī)雜志， 2022， 58（9）： 23-28.

QIN H N， ZHOU J， LI L， et al. Histological characteristics of hy?poxic adaptation of skeletal muscle of Himalayan marmot［J］. Chi?nese Journal of Veterinary Medicine， 2022， 58（9）： 23-28.

［6］ 李琳， 秦鴻楠， 周娟， 等. 喜馬拉雅旱獺心肌低氧適應(yīng)的組織學(xué)特點［J］. 黑龍江畜牧獸醫(yī)， 2023（5）： 127-131；140.

LI L， QIN H N， ZHOU J， et al. Histological characteristics ofmyocardial hypoxic adaptation in Marmota himalayana［J］. Hei?longjiang Animal Science and Veterinary Medicine， 2023（5）：127-131；140.

［7］ 李優(yōu)， 南新營， 趙坤鈺， 等. 基于簡化基因組測序技術(shù)的喜馬拉雅旱獺遺傳多樣性分析［J］. 野生動物學(xué)報， 2023， 44（1）：68-75.

LI Y， NAN X Y， ZHAO K Y， et al. Genetic diversity of Marmotahimalayana based on simplified genome sequencing technique［J］. Chinese Journal of Wildlife， 2023， 44（1）： 68-75.

［8］ SMITH D R. The past， present and future of mitochondrial genom?ics： Have we sequenced enough mtDNAs？［J］. Briefings in Func?tional Genomics， 2016， 15（1）： 47-54.

［9］ 李翔， 趙海紅， 徐小青， 等. 青海省烏蘭縣與同德縣鼠疫噬菌體分離及流行病學(xué)比較分析［J］. 醫(yī)學(xué)動物防制， 2024， 40（2）： 113-116；121.

LI X， ZHAO H H， XU X Q， et al. Comparative analysis of isola?tion and epidemiology of Yersinia pestis phage in Wulan Countyand Tongde County of Qinghai province［J］. Journal of MedicalPest Control， 2024， 40（2）： 113-116；121.

［10］ 李勝， 靳娟， 何建， 等. 我國喜馬拉雅旱獺鼠與南方家鼠鼠疫疫源地鼠疫菌遺傳特征研究［J］. 中國熱帶醫(yī)學(xué)， 2023， 23（9）： 916-921.

LI S， JIN J， HE J， et al. Genetic study of Yersinia pestis strainsisolated from the Himalayan marmot natural focus area and do?mestic rat plague focus area in southern China［J］. China Tropi?cal Medicine， 2023， 23（9）： 916-921.

［11］ HRDLICKOVA R， TOLOUE M， TIAN B. RNA-Seq methods fortranscriptome analysis［J］. Wiley Interdisciplinary ReviewsRNA， 2017， 8（1）. DOI： 10. 1002/wrna. 1364.

［12］ LYU T S， ZHOU S Y， FANG J H， et al. Convergent genomicsignatures of high-altitude adaptation among six independentlyevolved mammals［J］. Animals， 2022， 12（24）： 3572.

［13］ 李勃深， 張宇軒， 范容暉， 等. 高原對線粒體功能及能量代謝影響研究進展［J］. 中國比較醫(yī)學(xué)雜志， 2023， 33（10）：106-113.

LI B S， ZHANG Y X， FAN R H， et al. Research progresses onthe effects of high altitude on mitochondrial function and energymetabolism［J］. Chinese Journal of Comparative Medicine，2023， 33（10）： 106-113.

［14］ 文禹粱. 基于mRNA/miRNA調(diào)控網(wǎng)絡(luò)、線粒體動力學(xué)和功能解析藏綿羊低氧適應(yīng)機制［D］ . 蘭州： 甘肅農(nóng)業(yè)大學(xué)， 2022.

WEN Y L. Adaptation mechanism of hypoxia based on mRNA/miRNA regulation network， mitochondrial dynamics and func?tions on Qinghai-Tibet Plateau in Tibetan sheep［D］. Lanzhou：Gansu Agricultural University， 2022.

［15］ XIAO B， WANG S Y， YANG G S， et al. HIF-1α contributes tohypoxia adaptation of the naked mole rat［J］. Oncotarget， 2017，8（66）： 109941-109951.

［16］ ZHAO P F， HE Z H， XI Q M， et al. Variations in HIF-1α con?tributed to high altitude hypoxia adaptation via affected oxygenmetabolism in Tibetan sheep［J］. Animals， 2021， 12（1）： 58.

［17］ ZHANG Y， ZHANG H L， ZHANG B， et al. Identification ofkey HIF-1α target genes that regulate adaptation to hypoxic con?ditions in Tibetan chicken embryos［J］. Gene， 2020， 729：144321.

［18］ 高黎， 趙志英. MAPK信號轉(zhuǎn)導(dǎo)通路及其意義［J］. 世界最新醫(yī)學(xué)信息文摘， 2015， 15（53）： 25-27.

GAO L， ZHAO Z Y. MAPK signal transduction pathway and itssignificance［J］. World Latest Medicine Information， 2015， 15（53）： 25-27.

［19］ 鄧軍， 陳中葉， 拉門措， 等. 青藏高原強紫外輻射適應(yīng)機制研究進展［J］. 青海草業(yè)， 2023， 32（2）： 41-43；52.

DENG J， CHEN Z Y， LA M C， et al. Advances in research onthe adaptation mechanism of strong ultraviolet radiation on theQinghai-Tibet Plateau［J］. Qinghai Prataculture， 2023， 32（2）：41-43；52.

［20］ 張倩， 崔燕， 余四九， 等. 初生及成年牦牛大腦皮質(zhì)的轉(zhuǎn)錄組分析［J］. 畜牧獸醫(yī)學(xué)報， 2023， 54（11）： 4605-4614.

ZHANG Q， CUI Y， YU S J， et al. Transcriptome analysis of thecerebral cortex in newborn and adult yaks （Bos grunniens）［J］.Acta Veterinaria et Zootechnica Sinica， 2023， 54（11）： 4605-4614.

［21］ 吳娟. 紫外輻射對不同海拔高原林蛙的影響及其生理機制研究［D］. 蘭州： 蘭州大學(xué)， 2023.

WU J. A study on the effects of ultraviolet radiation on Rana ku?kunoris at different altitudes and its physiological mecanisms［D］. Lanzhou： Lanzhou University， 2023.

［22］ SONG H P， CHU Z G， ZHANG D X， et al. PI3K-AKT pathwayprotects cardiomyocytes against hypoxia-induced apoptosis byMitoKATP-mediated mitochondrial translocation of pAKT［J］.Cellular Physiology and Biochemistry， 2018， 49（2）： 717-727.

基金項目：青海省科技計劃項目（2018-ZJ-775）