摘要：為了明晰沙鞭（Psammochloa villosa）的轉(zhuǎn)錄組特征，本研究利用PacBio Sequel測(cè)序平臺(tái)首次對(duì)其進(jìn)行全長(zhǎng)轉(zhuǎn)錄組測(cè)序和數(shù)據(jù)分析，結(jié)果共獲得323 309個(gè)clean reads，自我矯正產(chǎn)生環(huán)形一致性序列（CCS）673 540個(gè)，預(yù)測(cè)得到蛋白編碼區(qū)（CDS）序列28 447個(gè)；MISA軟件共搜索得到沙鞭93 563個(gè)簡(jiǎn)單重復(fù)序列（SSR），分布于56 824條unigene上；轉(zhuǎn)錄本基因功能注釋，共有166 541條序列得到NR注釋，結(jié)果顯示沙鞭與二穗短柄草（Brachypodium distachyon）親緣最近；KOG數(shù)據(jù)庫(kù)比對(duì)將97 892個(gè)unigene分為25個(gè)功能類別，其中注釋較多的功能為翻譯后修飾、蛋白質(zhì)周轉(zhuǎn)和分子伴侶；GO數(shù)據(jù)庫(kù)比對(duì)共注釋到87 930個(gè)unigene，分為生物過(guò)程、細(xì)胞組成和分子功能3大類及62個(gè)亞類分支；KEGG結(jié)果表明碳水化合物代謝、信號(hào)轉(zhuǎn)導(dǎo)、能量代謝通路中注釋基因較多。本研究結(jié)果豐富了沙鞭的遺傳信息，為今后沙鞭關(guān)鍵耐旱基因的挖掘提供了理論依據(jù)。

關(guān)鍵詞：禾本科；沙鞭；全長(zhǎng)轉(zhuǎn)錄組；功能注釋；SSR

中圖分類號(hào)：Q346+.5；Q949.95

文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1007-0435（2023）06-1673-09

Full-length Transcriptome Sequencing and Bioinformatic Analysis of

Psammochloa villosa （Poaceae）

MAO Xuan-rui¹， SU Xu^{1，2，3}， LIU Yu-ping^{1，2，4*}， LIU Tao¹， CHEN Jin-yuan^1，3，

HU Xia-yu¹， YANG Ping¹， LI Xiao-li¹， YUAN Yuan¹

（1. School of Life Sciences， Qinghai Normal University， Xining， Qinghai Province 810008， China; 2. Academy of Plateau Science and

Sustainability， Qinghai Normal University， Xining， Qinghai Province 810016， China; 3. Key Laboratory of Biodiversity Formation Mechanism

and Comprehensive Utilization of the Qinghai-Tibet Plateau in Qinghai Province， Qinghai Normal University， Xining， Qinghai Province

810008， China; 4. Key Laboratory of Medicinal Animal and Plant Resources of the Qinghai-Tibet Plateau in Qinghai Province， Qinghai

Normal University， Xining， Qinghai Province 810008， China）

Abstract：In order to explore the transcriptome characteristics of Psammochloa villosa （Poaceae），the full-length transcriptome was firstly sequenced，then the genomic library of P. villosa were established by using the PacBio Sequel platform in this study. The result showed that 323 309 clean reads were obtained from that transcriptome sequence，which could form 673 540 circular consensus reads （CCS） through self-correction and a total of 28 447 coding sequences （CDS） in prediction. Moreover，93 563 SSRs in the transcriptome of P. villosa were searched and found out by MISA software，which were distributed in 56 824 unigenes. The transcripts of P.villosa was matched to the counterparts in NR database for the gene function，and it was identified out that 166 541 unigenes of P. villosa were annotated and that P. villosa had the closest genetic relationship with Brachypodium distachyon. Based on KOG database，97 892 unigenes within the transcriptome of P. villosa were classified into 25 functional groups，most of which were annotated in the participation of the posttranslational modification，protein turnover and chaperones. Similarly，87 930 unigenes within that transcriptome were annotated to GO database，whose functions were divided into 62 subclades of three clades of the biological processes，cellular components，and molecular functions. Besides，the functional analysis of unigenes referred to the KEGG database indicated that the functions of annotated genes in that transcriptome were related with carbohydrate metabolism，signal transduction and energy metabolism. These results above enriched the genetic information of P. villosa，and could provide a theoretical basis for excavating the key drought-resistant genes of P. villosa in coming days.

Key words：Poaceae;Psammochloa villosa;Full-length transcriptome;Function annotation;SSR

沙鞭（Psammochloa villosa）是禾本科（Poaceae）、沙鞭屬（Psammochloa）的一種多年生草本植物，主要生長(zhǎng)于海拔900～2 900 m的沙地和荒漠中^［1］，主產(chǎn)青海、甘肅、新疆、陜西北部、內(nèi)蒙古等省區(qū)。沙鞭具有極強(qiáng)的抗風(fēng)沙、抗旱和抗寒能力，根莖十分發(fā)達(dá)，是良好的固沙植物^［2-3］。目前，國(guó)內(nèi)外關(guān)于沙鞭的研究主要集中于形態(tài)性狀^［4-8］、生理特征^［9-10］、生態(tài)適應(yīng)性^［11-13］、遺傳多樣性^［14-16］及種質(zhì)資源^［17］等領(lǐng)域。譬如，呂婷等^［7］通過(guò)觀察沙鞭50個(gè)代表種群的葉表皮微形態(tài)特征，發(fā)現(xiàn)沙鞭種群葉表皮氣孔密度與環(huán)境干旱程度呈正相關(guān)，推測(cè)這種特征有利于沙鞭減少體內(nèi)的水分散失，旨在更好地適應(yīng)荒漠環(huán)境；黃振英等^［9］比較了不同深度沙層對(duì)沙鞭種子萌發(fā)和幼苗出土特性的影響，發(fā)現(xiàn)隨沙層深度的增加，沙鞭種子的萌發(fā)率和幼苗出土率均呈下降趨勢(shì)；賀學(xué)禮等^［11］和徐翥骍等^［12］通過(guò)探討叢枝菌根（Arbuscular mycorrhiza，AM）真菌的時(shí)空分布特征以及侵入群落間隔空地前后對(duì)沙鞭AM真菌的影響，發(fā)現(xiàn)沙鞭與AM真菌可形成良好的共生關(guān)系。然而，目前國(guó)內(nèi)外對(duì)沙鞭進(jìn)行全長(zhǎng)轉(zhuǎn)錄組測(cè)序及生物學(xué)信息學(xué)分析的研究尚屬空白。

轉(zhuǎn)錄組是指某一功能狀態(tài)下特定細(xì)胞轉(zhuǎn)錄出的總RNA^［18］。轉(zhuǎn)錄組測(cè)序是挖掘基因結(jié)構(gòu)和功能的關(guān)鍵技術(shù)手段^［19］。其中，全長(zhǎng)轉(zhuǎn)錄組測(cè)序（Full Length Transcriptome Sequencing）是基于第三代測(cè)序平臺(tái)（PacBio Sequel）對(duì)特定細(xì)胞或組織在某個(gè)時(shí)期或某一特定狀態(tài)下轉(zhuǎn)錄出的總mRNA進(jìn)行研究。目前，該全長(zhǎng)轉(zhuǎn)錄組測(cè)序技術(shù)已被廣泛應(yīng)用于醫(yī)學(xué)、生物學(xué)、藥物研發(fā)等諸多領(lǐng)域，具有讀長(zhǎng)長(zhǎng)、可對(duì)反轉(zhuǎn)錄全長(zhǎng)cDNA直接讀取等優(yōu)點(diǎn)，能夠獲得高質(zhì)量的單個(gè)RNA分子的全部序列，實(shí)現(xiàn)挖掘新基因的目的；同時(shí)，該測(cè)序技術(shù)也是分子標(biāo)記和基因資源開(kāi)發(fā)的一個(gè)重要途徑，現(xiàn)已為多種植物的轉(zhuǎn)錄組建立了豐富的遺傳數(shù)據(jù)庫(kù)^［20］。據(jù)此，本研究以沙鞭為材料，利用高通量測(cè)序技術(shù)對(duì)其三代全長(zhǎng)轉(zhuǎn)錄組進(jìn)行測(cè)序和數(shù)據(jù)分析，旨在豐富該物種的遺傳信息，并為后續(xù)挖掘沙鞭的關(guān)鍵耐旱基因提供理論參考。

1 材料與方法

1.1 試驗(yàn)材料

試驗(yàn)樣品采自甘肅省武威市民勤縣（103.65° N，39.17° E，海拔1 250 m），即從沙鞭生長(zhǎng)良好且健康的植株上選取新鮮的葉片、嫩莖和根進(jìn)行混合，然后用錫箔紙包裹后迅速保存于液氮中，用于后續(xù)RNA提取。憑證標(biāo)本存放于中國(guó)科學(xué)院西北高原生物研究所青藏高原生物標(biāo)本館（HNWP）。

1.2 試驗(yàn)方法

1.2.1 總RNA提取和檢測(cè) 試驗(yàn)樣品于液氮預(yù)冷的研缽中研磨并不斷補(bǔ)充液氮，直至成粉末狀且無(wú)明顯可見(jiàn)顆粒為止，隨后轉(zhuǎn)移至裝有裂解液的EP管中，采用Trizol法提取總RNA^［21］；利用瓊脂糖凝膠電泳檢測(cè)RNA純度，采用Nanodrop測(cè)定RNA的D260/280比值，運(yùn)用Qubit精準(zhǔn)定量RNA濃度，使用Agilent 2100分析RNA完整性，純度和濃度檢測(cè)合格的RNA樣品用于全長(zhǎng)轉(zhuǎn)錄組測(cè)序和文庫(kù)構(gòu)建。

1.2.2 文庫(kù)構(gòu)建和數(shù)據(jù)評(píng)估 采用SMARTer PCR cDNA Synthesis Kit（Clontech，634925）將質(zhì)檢合格的總RNA反轉(zhuǎn)錄成cDNA，然后選取部分cDNA利用BluePippin篩選和富集大于4 kb的片段；篩選的片段進(jìn)行PCR擴(kuò)增、全長(zhǎng)cDNA純化、末端修復(fù)和連接接頭后構(gòu)建文庫(kù)；最后，利用PacBio Sequel測(cè)序平臺(tái)對(duì)質(zhì)檢合格的文庫(kù)進(jìn)行全長(zhǎng)轉(zhuǎn)錄組測(cè)序。下機(jī)得到的原始數(shù)據(jù)采用SMRTlink軟件的Iso-Seq3程序評(píng)估質(zhì)量，即獲得的subreads序列首先提取環(huán)形一致性序列（Circular Consensus Sequences，CCS），然后根據(jù)CCS是否含有3′引物、5′引物和PloyA尾將其劃分為全長(zhǎng)非嵌合序列（Full-Length Non-Chimeric Read，F(xiàn)LNC）和非全長(zhǎng)序列（Non-full-length Read，NFL），接著對(duì)FLNC相似的序列聚類得到一致性序列，對(duì)其進(jìn)行矯正（Polishing）獲得高質(zhì)量的全長(zhǎng)轉(zhuǎn)錄本序列（Polished High-quality Isoforms），最后使用CD-HIT軟件通過(guò)序列比對(duì)選擇相似性gt;95%的高質(zhì)量全長(zhǎng)轉(zhuǎn)錄本序列去冗余，獲得沙鞭的三代全長(zhǎng)轉(zhuǎn)錄組unigenes序列。

1.2.3 CDS預(yù)測(cè)及SSR分析 使用TransDecoder軟件中的TransDecoder. LongOrfs對(duì)開(kāi)放閱讀框（Open Reading Frame，ORF）進(jìn)行搜索^［24］，鑒定出最長(zhǎng)的ORF序列，然后選擇500個(gè)最長(zhǎng)的ORF序列作為訓(xùn)練集參數(shù)，利用Markov模型預(yù)測(cè)并打分。若假定的ORF編碼結(jié)構(gòu)分值為正且高于其他假定錯(cuò)誤的ORF編碼結(jié)構(gòu)分值，則該編碼結(jié)構(gòu)確定為ORF；若高分值的ORF被長(zhǎng)ORF的其他不同閱讀框全部包含，則該ORF排除；若預(yù)測(cè)的氨基酸序列有Pfam domain或SwissProt數(shù)據(jù)注釋，則該ORF優(yōu)先篩選。根據(jù)上述標(biāo)準(zhǔn)，對(duì)沙鞭的每個(gè)unigene序列篩選出一個(gè)最佳的ORF預(yù)測(cè)結(jié)果，然后利用TransDecoder. Predict在最佳ORF結(jié)果中預(yù)測(cè)CDS序列。使用MISA（Microsatellite Identification Tool）軟件按照單核苷酸至六核苷酸的最少重復(fù)次數(shù)分別不少于10，6，5，5，5和5搜索轉(zhuǎn)錄本unigene序列中的SSR位點(diǎn)^［25-27］。

1.2.4 全長(zhǎng)轉(zhuǎn)錄組序列功能注釋 利用BLAST軟件將沙鞭的unigene序列與非冗余蛋白數(shù)據(jù)庫(kù)（Non-Redundant protein database，NR）、京都基因與基因組百科全書(shū)（Kyoto encyclopedia of genes and genomes，KEGG）^［31］、基因本體論數(shù)據(jù)庫(kù)（Gene ontology，GO）^［32］、蛋白質(zhì)真核同源數(shù)據(jù)庫(kù)（Eukaryotic ortholog groups，KOG）^［33］、蛋白質(zhì)家族數(shù)據(jù)庫(kù)（Protein families database，Pfam）^［34］、蛋白質(zhì)序列數(shù)據(jù)庫(kù)（SwissProt protein sequence database，SwissProt）^［35］等公共數(shù)據(jù)庫(kù)分別進(jìn)行比對(duì)和基因功能注釋^［28-30］。

2 結(jié)果與分析

2.1 轉(zhuǎn)錄組測(cè)序和組裝

2.1.1 測(cè)序結(jié)果 利用第三代測(cè)序平臺(tái)PacBio Sequel對(duì)沙鞭根、莖和葉的混合樣品進(jìn)行全長(zhǎng)轉(zhuǎn)錄組測(cè)序，共獲得673 786個(gè)clean reads，平均序列長(zhǎng)度為30 433 bp，N50為53 500 bp。使用SMRTLink軟件Iso-Seq3流程從subreads中獲得環(huán)形一致性序列（CCS）673 540個(gè)，平均序列長(zhǎng)度為2 197 bp，其中全長(zhǎng)序列（Full-length Read）共478 972個(gè)，占CCS序列總數(shù)的71.11%；全長(zhǎng)非嵌合體序列（Full-length non-chimeric Read）共463 338個(gè)（68.79%）；非全長(zhǎng)序列（Non-full-length Read）共194 314個(gè)（28.85%）；對(duì)FLNC序列進(jìn)行校正，獲得37 890個(gè)（Polished High-quality Isoforms）；最終利用CD-HIT軟件對(duì)高質(zhì)量的全長(zhǎng)轉(zhuǎn)錄本序列去冗余后共獲得unigene數(shù)目為184 076個(gè)，平均長(zhǎng)度2 461.22 bp，N50長(zhǎng)度3 292 bp。轉(zhuǎn)錄組測(cè)序數(shù)據(jù)已上傳至NCBI的SRA數(shù)據(jù)庫(kù)，登錄號(hào)PRJNA730046。

2.1.2 CDS預(yù)測(cè) CDS（Coding Sequence）是編碼一段蛋白產(chǎn)物的序列。本研究利用編碼預(yù)測(cè)軟件TransDecoder對(duì)獲得的基因片段進(jìn)行編碼區(qū)預(yù)測(cè)，獲得蛋白編碼區(qū)序列28 447個(gè)，長(zhǎng)度范圍為200～7 896 bp，主要分布于200～1 600 bp，其中長(zhǎng)度≥500 bp的序列占總序列數(shù)的69.08%，≥1 000 bp的序列占比31.18%（圖1）。

2.1.3 SSR和LncRNA分析 簡(jiǎn)單重復(fù)序列（Simple Sequence Repeat，SSR）又稱微衛(wèi)星序列（Microsatellite），是研究動(dòng)植物遺傳多樣性的一種重要分子標(biāo)記。本研究通過(guò)對(duì)沙鞭184 076條轉(zhuǎn)錄本的SSR檢索，共獲得93 563個(gè)SSR重復(fù)序列，其分布于56 824條unigene上，且重復(fù)類型以單核苷酸、二核苷酸和三核苷酸為主，其中單核苷酸重復(fù)最多（50 027）、占比53.47%，三核苷酸次之（26 112）、占比27.91%，二核苷酸最少（15 609）、占比16.68%，四、五、六核苷酸重復(fù)序列類型極少，分別為1.10%，0.34%和0.50%（表1）。SSR搜索結(jié)果統(tǒng)計(jì)分析，發(fā)現(xiàn)二核苷酸重復(fù)有4種基元類型，其中（AG/CT）n重復(fù)次數(shù)最多，形成8 691個(gè)SSR位點(diǎn)，（CG/CG）n重復(fù)次數(shù)最少（862）；三核苷酸重復(fù)有10種基元類型，（CCG/CGG）n最多（9 175），（ACT/AGT）n最少（245）；四、五、六核苷酸重復(fù)各有31，42，79種，占比1.94%（圖2）。此外，長(zhǎng)鏈非編碼RNA（LncRNA）是一類不編碼蛋白質(zhì)的RNA分子。采用轉(zhuǎn)錄本蛋白編碼潛能預(yù)測(cè)工具PLEK通過(guò)預(yù)測(cè)沙鞭轉(zhuǎn)錄本的LncRNA，并根據(jù)Transdecoder預(yù)測(cè)結(jié)果對(duì)候選的LncRNA進(jìn)行過(guò)濾，本研究共獲得LncRNA 14 569個(gè)，最短長(zhǎng)度為201 bp，最長(zhǎng)為8 742 bp，平均長(zhǎng)度748.15 bp。

2.2 Unigene功能注釋

2.2.1 NR功能注釋和分類 沙鞭三代全長(zhǎng)轉(zhuǎn)錄組數(shù)據(jù)與NR數(shù)據(jù)庫(kù)比對(duì)，結(jié)果發(fā)現(xiàn)184 076個(gè)unigene中有166 541條得到注釋，與其相似的前10個(gè)物種分別是二穗短柄草（Brachypodium dis-tachyon，72 093個(gè)，47.7%）、節(jié)節(jié)麥（Aegilops tauschii，12 811個(gè)，8.5%）、粳稻（Oryza sativa subsp. japonica，12 803個(gè)，8.5%）、烏拉圖小麥（Triticum urartu，8 359個(gè)，5.5%）、谷子（Setaria italic，6 726個(gè)，4.5%）、短花藥野生稻（Oryza brachyantha，6 191個(gè)，4.1%）、小麥（Triticum aestivum，5 402個(gè)，3.6%）、玉米（Zea mays，3 546個(gè)，2.3%）、高粱（Sorghum bicolor，3 159個(gè)，2.1%）、秈稻（Oryza sativa subsp. indica，2 941個(gè)，1.9%），其他16 973個(gè)，11.2%（圖3A）；從注釋物種看，二穗短柄草基因數(shù)最多，說(shuō)明沙鞭與二穗短柄草同源性較高，親緣關(guān)系較近。同時(shí)，與NR，COG，SwissProt，KEGG和Pfam數(shù)據(jù)庫(kù)比對(duì)后，注釋到的unigene分別為166 541，97 892，130 885，89 086和111 373個(gè)，特有unigene依次為17 640，126，141，1 217和353個(gè)，其中NR與SwissProt的相同unigene最多，有130 619個(gè)，關(guān)系密切；與KEGG的相同unigene最少，為87 633個(gè)；五個(gè)數(shù)據(jù)庫(kù)均注釋到的轉(zhuǎn)錄本為53 806個(gè)（圖3B）。

2.2.2 KOG功能注釋和分布 KOG/COG數(shù)據(jù)庫(kù)（Eukaryotic Ortholog Groups/Cluster of Orthologous Groups of Proteins）是基于真核生物、細(xì)菌及藻類的完整基因組編碼蛋白構(gòu)建的系統(tǒng)進(jìn)化關(guān)系。本研究將沙鞭轉(zhuǎn)錄本與KOG數(shù)據(jù)庫(kù)比對(duì)，共發(fā)現(xiàn)97 892個(gè)unigene得到KOG功能注釋，可以劃分為25個(gè)功能類別，其中注釋最多的功能類別為翻譯后修飾、蛋白質(zhì)周轉(zhuǎn)和分子伴侶，有12 246個(gè)，占KOG注釋unigene總數(shù)的12.5%；其次是一般功能預(yù)測(cè)和信號(hào)轉(zhuǎn)導(dǎo)機(jī)制，有11 319個(gè)，占比11.6%；細(xì)胞運(yùn)動(dòng)（63個(gè)）和細(xì)胞外結(jié)構(gòu)（33個(gè)）數(shù)量最少（圖4）。

2.2.3 GO功能注釋和分類 將沙鞭的unigene序列與GO數(shù)據(jù)庫(kù)比對(duì)，結(jié)果表明有87 930個(gè)unigene得到注釋，功能分為生物過(guò)程（Biological Processes，BP）（60 301個(gè)）、細(xì)胞組成（Cellular Component，CC）（50 618個(gè)）和分子功能（Molecular Function，MF）（70 552個(gè)）3大類和62個(gè)亞類分支。其中，細(xì)胞組分有18個(gè)亞類分支，以細(xì)胞和細(xì)胞部分最多，注釋到42 929個(gè)unigene，其次為細(xì)胞器（34 520個(gè)）、膜結(jié)構(gòu)（22 942個(gè)）和細(xì)胞器部分（11 981個(gè)），病毒體（16個(gè)）和其他細(xì)胞器（12個(gè)）最少；分子功能分為13個(gè)亞類分支，催化活性基因最多，為47 510個(gè)，與分子載體調(diào)節(jié)（5個(gè)）和蛋白質(zhì)標(biāo)記物（1個(gè)）有關(guān)基因最少；生物過(guò)程涉及31個(gè)亞類分支，代謝過(guò)程（46 512個(gè)）最多，細(xì)胞過(guò)程（40 541個(gè)）次之，細(xì)胞殺傷（1個(gè)）最少（表2）。

2.2.4 KEGG代謝通路分析 沙鞭KEGG代謝通路有89 086個(gè)unigene得到注釋，分為細(xì)胞加工（Cellular Processes）、環(huán)境信息處理（Environmental Information Processing）、遺傳信息處理（Genetic Information Processing）、代謝（Metabolism）等4大類、23小類，其中碳水化合物代謝（Carbohydrate Metabolism）（23 460個(gè)）、信號(hào)轉(zhuǎn)導(dǎo)（Signal Transduction）（22 642個(gè)）、能量代謝（Energy Metabolism）（18 617個(gè)）及全局和概覽圖（Global and Overview Maps）（18 267個(gè)）為基因注釋數(shù)量較多的代謝通路（表3）。

2.2.5 轉(zhuǎn)錄因子分析 基于植物轉(zhuǎn)錄因子數(shù)據(jù)庫(kù)PlantTFDB，本研究從沙鞭中共鑒定到1 110個(gè)轉(zhuǎn)錄因子，隸屬于35個(gè)轉(zhuǎn)錄因子家族，其中數(shù)量最多的前10位轉(zhuǎn)錄因子家族依次為FAR1（95，8.56%），MYB superfamily（92，8.29%），bHLH（87，7.84%），bZIP（79，7.12%），GARP（70，6.31%），C3H（66，5.95%），B3 superfamily（64，5.77%），C2C2（62，5.59%），WRKY（62，5.59%）和AP2/ERF（53，4.77%）（圖5）。

3 討論

全長(zhǎng)轉(zhuǎn)錄組測(cè)序是基于PacBio Sequel平臺(tái)的一種單分子實(shí)時(shí)測(cè)序技術(shù)，具有無(wú)需打斷和PCR擴(kuò)增、讀長(zhǎng)長(zhǎng)、可對(duì)反轉(zhuǎn)錄全長(zhǎng)cDNA直接讀取等特點(diǎn)，能夠獲得高質(zhì)量單個(gè)RNA分子的全部序列，精準(zhǔn)得到mRNA的同源異構(gòu)體、單堿基變異、可變剪接、同源基因和等位基因等信息^［36］，是轉(zhuǎn)錄組從頭測(cè)序的首選^［20］。本研究利用PacBio Sequel測(cè)序平臺(tái)對(duì)沙鞭的全長(zhǎng)轉(zhuǎn)錄組測(cè)序和文庫(kù)構(gòu)建，結(jié)果表明沙鞭去冗余后的轉(zhuǎn)錄本為184 076個(gè)、平均長(zhǎng)度2 461.22 bp、N50長(zhǎng)度3 292 bp。N50長(zhǎng)度是評(píng)價(jià)轉(zhuǎn)錄組數(shù)據(jù)組裝質(zhì)量的一個(gè)重要指標(biāo)，沙鞭轉(zhuǎn)錄組的N50長(zhǎng)度大于1 000 bp，比先前報(bào)道的薏苡（Coix lacryma-jobi）（1 700 bp）^［37］、大青木（Clerodendron cyrtophyllum）（1 055 bp）^［38］、披堿草（Elymus dahuricus）（481 bp）^［39］和茶樹(shù)（Camellia sinensis）（1 081 bp）^［40］的N50長(zhǎng)度長(zhǎng)，且GC含量穩(wěn)定，說(shuō)明沙鞭轉(zhuǎn)錄組的組裝質(zhì)量完整性高，完全能夠確保后續(xù)SSR搜索和基因功能注釋分析。同時(shí)，本研究還檢測(cè)到分布于56 824條unigene序列上的93 563個(gè)SSR重復(fù)序列，SSR發(fā)生頻率（含SSR位點(diǎn)unigene數(shù)占unigene總數(shù)比）為30.87%，明顯高于先前學(xué)者認(rèn)為的青稞（Hordeum vulgare）（16.49%）^［41］、甘蔗（Sacchrum spp.）（14.02%）^［42］和蘇丹草（Sorghum sudanense）（16.82%）^［43］，表明沙鞭轉(zhuǎn)錄組中SSR位點(diǎn)數(shù)量豐富，具有較高的發(fā)生頻率和開(kāi)發(fā)潛能。尤為重要的是，沙鞭轉(zhuǎn)錄組SSR重復(fù)單元中單核苷酸重復(fù)數(shù)量最多（50 027個(gè)），占比53.47%，超過(guò)總SSR的一半，二核苷酸（15 609個(gè)，16.68%）和三核苷酸（26 112，27.91%）重復(fù)數(shù)均在1.5萬(wàn)以上，說(shuō)明這些重復(fù)單元均可為沙鞭SSR分子標(biāo)記開(kāi)發(fā)提供足夠豐富的信息位點(diǎn)。

沙鞭轉(zhuǎn)錄組的功能注釋表明，168 650個(gè)unigene至少在1個(gè)數(shù)據(jù)庫(kù)中得到注釋，15 426個(gè)unigene在任意一個(gè)數(shù)據(jù)庫(kù)中得以注釋，我們認(rèn)為這可能是由于先前從未對(duì)沙鞭進(jìn)行過(guò)任何的測(cè)序和分析，數(shù)據(jù)庫(kù)中缺乏相關(guān)的注釋信息，或者是由沙鞭的轉(zhuǎn)錄組序列為非編碼或不完整序列所導(dǎo)致。就NR數(shù)據(jù)庫(kù)注釋而言，注釋最多的前10個(gè)物種均屬于禾本科，與沙鞭具有相對(duì)較近的親緣關(guān)系，其中比對(duì)到二穗短柄草的基因數(shù)最多，說(shuō)明沙鞭與二穗短柄草的同源性較高，兩者親緣關(guān)系最近，這為禾本科植物基因表達(dá)差異分析提供了基礎(chǔ)數(shù)據(jù)。從KOG數(shù)據(jù)庫(kù)來(lái)看，沙鞭的unigene序列依據(jù)功能不同劃分為25類，其中注釋較多的功能類別為翻譯后修飾、蛋白質(zhì)周轉(zhuǎn)和分子伴侶；同樣，根據(jù)GO數(shù)據(jù)庫(kù)，沙鞭的unigene分為3大類、62亞類分支，其中以催化活性和代謝過(guò)程最多；KEGG數(shù)據(jù)庫(kù)則注釋到4大類、23小類，其中涉及碳水化合物代謝途徑的unigene最多，可能是由于碳水化合物在沙鞭生長(zhǎng)發(fā)育和耐旱性方面有重要作用，為今后探討沙鞭耐旱性代謝通路提供了理論依據(jù)。

先前研究表明，F(xiàn)AR1是由轉(zhuǎn)座酶衍生而來(lái)的一類轉(zhuǎn)錄因子家族，通過(guò)調(diào)控植物下游基因表達(dá)，在光響應(yīng)、淀粉合成代謝、生長(zhǎng)發(fā)育及逆境脅迫應(yīng)答等方面具有重要的調(diào)控作用^［44］，其大多數(shù)成員能被NaCl和PEG6000脅迫誘導(dǎo)而被ABA抑制，從而使植物達(dá)到一種動(dòng)態(tài)平衡以抵御外界環(huán)境脅迫^［45］。譬如，羅文巧^［46］通過(guò)對(duì)二穗短柄草的研究，發(fā)現(xiàn)FAR轉(zhuǎn)錄因子家族參與合成木栓質(zhì)中的初級(jí)脂肪醇，木栓質(zhì)的大量積累能顯著提高其保水能力，提高了二穗短柄草的抗旱能力；同樣，MYB superfamily和bHLH也是植物中數(shù)量龐大的轉(zhuǎn)錄因子家族，其已被充分證實(shí)對(duì)植物生長(zhǎng)發(fā)育和各種非生物脅迫發(fā)揮著重要功能，可以提高植物的抗逆性^［47-49］。本研究對(duì)沙鞭轉(zhuǎn)錄因子的分析鑒定出35個(gè)家族、1 110個(gè)轉(zhuǎn)錄因子，其中FAR1、MYB superfamily和bHLH位居前三，這些轉(zhuǎn)錄因子家族的存在可能與沙鞭有較強(qiáng)的抗旱和抗寒等生理特性密切相關(guān)，從而進(jìn)一步驗(yàn)證了先前研究的正確性和合理性^［46-49］，為將來(lái)沙鞭及近緣種耐旱分子機(jī)制的研究提供了理論依據(jù)?？傊狙芯客ㄟ^(guò)對(duì)沙鞭轉(zhuǎn)錄組測(cè)序和數(shù)據(jù)分析，彌補(bǔ)了先前沙鞭基因組遺傳信息的不足，尤其轉(zhuǎn)錄組中富集到的SSR位點(diǎn)為將來(lái)沙鞭SSR引物開(kāi)發(fā)、遺傳多樣性分析和親緣關(guān)系探討乃至新基因發(fā)現(xiàn)等均提供了重要參考依據(jù)。

4 結(jié)論

基于沙鞭全長(zhǎng)轉(zhuǎn)錄組數(shù)據(jù)分析，本研究共獲得184 076個(gè)去冗余后的轉(zhuǎn)錄本、93 563個(gè)SSR重復(fù)序列，表現(xiàn)出較高的多態(tài)性潛能；NR注釋結(jié)果前10位均屬于禾本科，比對(duì)到二穗短柄草的基因數(shù)最多，說(shuō)明沙鞭與二穗短柄草的同源性最高，親緣關(guān)系最近；基因功能主要以催化活性和代謝過(guò)程、翻譯后修飾、蛋白質(zhì)周轉(zhuǎn)和分子伴侶、碳水化合物代謝途徑為主，推測(cè)可能是碳水化合物在沙鞭生長(zhǎng)發(fā)育和耐旱性方面具有重要調(diào)節(jié)作用；轉(zhuǎn)錄因子家族分析發(fā)現(xiàn)位居前三的FAR1、MYB superfamily和bHLH均涉及沙鞭生長(zhǎng)發(fā)育及逆境脅迫應(yīng)答調(diào)控。

參考文獻(xiàn)

［1］ 郭本兆. 中國(guó)植物志［M］. 北京：科學(xué)出版社，1997：309

［2］ YU F H，DONG M，BERTIL K. Clonal integration helps Psammochloa villosa survive sand burial in an inland dune ［J］. New Phytologist，2004，162（3）：697-904

［3］ 劉峰，馬子蘭，劉濤，等. 不同生態(tài)環(huán)境條件沙生植物沙鞭的結(jié)實(shí)研究［J］. 廣西植物，2021，41（9）：1457-1464

［4］ 馬毓泉. 內(nèi)蒙古植物志［M］. 呼和浩特：內(nèi)蒙古人民出版社，1997：115-152

［5］ 董鳴. 切斷根莖對(duì)根莖禾草沙鞭和賴草克隆生長(zhǎng)的影響［J］. 植物學(xué)報(bào)，1999（2）：194-198

［6］ 董鳴，阿拉騰寶，邢雪榮，等. 根莖禾草沙鞭的克隆基株及分株種群特征［J］. 植物生態(tài)學(xué)報(bào)，1999（4）：15-23

［7］ 呂婷，劉玉萍，周勇輝，等. 沙生植物沙鞭葉表皮微形態(tài)特征及其生態(tài)適應(yīng)性研究［J］. 植物研究，2018，38（3）：330-337

［8］ 呂婷，劉濤，梁瑞芳，等. 沙生植物沙鞭不同居群形態(tài)變異研究［J］. 植物研究，2021，41（1）：60-66

［9］ 黃振英. 鄂爾多斯高原固沙禾草沙鞭種子休眠和萌發(fā)與環(huán)境的關(guān)系［J］. 西北植物學(xué)報(bào)，2003（7）：72-77

［10］黃振英，董鳴，張淑敏. 沙鞭（禾本科）種子在沙丘上的萌發(fā)策略及幼苗的耐干燥特性［J］. 生態(tài)學(xué)報(bào)，2005（2）：298-303

［11］賀學(xué)禮，李英鵬，趙麗莉，等. 毛烏素沙地克隆植物沙鞭生長(zhǎng)對(duì)AM真菌生態(tài)分布的影響［J］. 生態(tài)學(xué)報(bào)，2010，30（3）：751-758

［12］徐翥骍，賀學(xué)禮，郭輝娟，等. 內(nèi)蒙古農(nóng)牧交錯(cuò)區(qū)沙鞭和羊柴AM真菌侵染及其土壤因子［J］. 河北農(nóng)業(yè)大學(xué)學(xué)報(bào)，2011，34（1）：56-61

［13］CHEN M Y，LIN L，ZHANG Y M. Genome sequence of Klebsiella oxytoca SA2，an endophytic nitrogen-fixing bacterium isolated from the pioneer grass Psammochloa villosa ［J］. Genome Announcements，2013，1（4）：e00601

［14］王可青，葛頌，董鳴. 根莖禾草沙鞭的等位酶變異及克隆多樣性［J］. 植物學(xué)報(bào)，1999，41（5）：537-540

［15］LI A，GE S. Genetic variation and clonal diversity of Psammochloa villosa （Poaceae） detected by ISSR markers ［J］. Annals of Botany，2001（87）：585-590

［16］呂婷. 沙生牧草沙鞭的遺傳多樣性研究［D］. 西寧：青海師范大學(xué)，2019：2-62

［17］呂婷，劉玉萍，周勇輝，等. 荒漠植物沙鞭（Psammochloa villosa）種質(zhì)資源收集及譜系遺傳分化初探［J］. 草地學(xué)報(bào)，2018，26（3）：733-740

［18］WANG Z，GERSTEIN M，SNYDER M. RNA-Seq：a revolutionary tool for transcriptomics ［J］. Nature Reviews Genetics，2009，10（1）：57-63

［19］LIAO X，LIU Y Q，HAN T Y，et al. Full-Length transcriptome sequencing reveals tissue-specific gene expression profile of mangrove clam Geloina erosa ［J］. Frontiers in Physiology，2022（13）：851957

［20］FRANSSEN S U，SHRESTHA R P，BRAUTIGAM A，et al. Comprehensive transcriptome analysis of the highly complex Pisum sativum genome using next generation sequencing［J］. BMC Genomics，2011，12（1）：227

［21］劉春曉，黃小慶，劉自廣，等. 十字花科植物種子低分子RNA提取方法比較［J］. 基因組學(xué)與應(yīng)用生物學(xué)，2019，38（3）：1236-1241

［22］崔凱，吳偉偉，刁其玉. 轉(zhuǎn)錄組測(cè)序技術(shù)的研究和應(yīng)用進(jìn)展［J］. 生物技術(shù)通報(bào)，2019，35（7）：1-9

［23］孫銘陽(yáng)，徐世強(qiáng)，顧艷，等. 穿心蓮全長(zhǎng)轉(zhuǎn)錄組測(cè)序及特性分析［J］. 中國(guó)農(nóng)學(xué)通報(bào)，2021，37（27）：82-89

［24］潘敏，于旭東，蔡澤坪，等. 菠蘿蜜莖葉全長(zhǎng)轉(zhuǎn)錄組分析［J］. 熱帶作物學(xué)報(bào)，2020，41（7）：1288-1297

［25］THIEL T，MICHALEK W，VARSHNEY R K，et al. Exploiting EST databases for the development and characterization of gene derived SSR-markers in barley （Hordeum vulgare L.） ［J］. Theoretical and Applied Genetics，2003，106（3）：411-422

［26］毛軒睿，劉玉萍，蘇旭，等. 沙鞭轉(zhuǎn)錄組簡(jiǎn)單重復(fù)序列（SSR）位點(diǎn)特征分析［J］. 草地學(xué)報(bào)，2022，30（8）：1990-2001

［27］富貴，劉玉萍，蘇旭. 基于轉(zhuǎn)錄組數(shù)據(jù)的密花香薷SSR位點(diǎn)特征分析［J］. 西北植物學(xué)報(bào)，2021，41（4）：654-663

［28］ALTSCHUL S F，GISH W，MILLER W，et al. Basic local alignment search tool［J］. Journal of Molecular Biology，1990，215（3）：403-410

［29］陳春光，米華玲. 天麻素生物合成途徑相關(guān)基因的分析［J］. 植物生理學(xué)報(bào)，2021，57（9）：1819-1828

［30］周延清，邵露營(yíng)，李慧敏，等. 地黃全長(zhǎng)轉(zhuǎn)錄組分析及毛蕊花糖苷生物合成下游途徑關(guān)鍵酶基因的挖掘與解析［J］. 植物生理學(xué)報(bào)，2022，58（1）：197-206

［31］WIXON J，KELL D. The kyoto encyclopedia of genes and genomes-KEGG ［J］. Yeast，2000，17（1）：48-55

［32］MICHAEL A，CATHERINE A B，JUDITH A B，et al. Gene ontology：tool for the unification of biology ［J］. Nature Genetics，2000，25（1）：25-29

［33］TATUSOV R，F(xiàn)EDOROVA N，JACKSON J，et al. The COG database：an updated version includes eukaryotes ［J］. BMC Bioinformatics，2003，4（1）：41

［34］BATEMAN A，COIN L，DURBIN R，et al. The Pfam protein family’s database ［J］. Nucleic Acids Research，2004，32（1）：138-141

［35］BAIROCH A，APWEILER R. The SWISS-PROT protein sequence data bank and its supplement TrEMBL in 2000 ［J］. Nucleic Acids Research，2000，28（1）：45-48

［36］夏麗飛，孫云南，宋維希，等. 基于PacBio平臺(tái)的紫娟茶樹(shù)全長(zhǎng)轉(zhuǎn)錄組分析［J］. 基因組學(xué)與應(yīng)用生物學(xué)，2020，39（6）：2646-2658

［37］歐陽(yáng)雨晴，李玲玲，石好宇，等. 薏苡轉(zhuǎn)錄組測(cè)序及基因功能注釋［J］. 中南藥學(xué)，2021，19（7）：1286-1293

［38］劉敏，黃煒忠，何孟璐，等. 大青轉(zhuǎn)錄組測(cè)序及生物信息學(xué)分析［J］. 廣州中醫(yī)藥大學(xué)學(xué)報(bào)，2022，39（1）：177-183

［39］馬玉金. 七個(gè)禾本科牧草轉(zhuǎn)錄組測(cè)序（RNA-Seq）及轉(zhuǎn)錄組數(shù)據(jù)分析與基因發(fā)掘［D］. 昆明：昆明理工大學(xué)，2016：33-38

［40］龐丹丹，張亞真，徐禮羿，等. 基于轉(zhuǎn)錄組測(cè)序的茶樹(shù)兒茶素合成調(diào)控相關(guān)基因的挖掘［J］. 分子植物育種，2019，17（17）：5604-5612

［41］徐金青，夏騰飛，王蕾，等. 青稞轉(zhuǎn)錄組SSR位點(diǎn)及其基因功能分析［J］. 麥類作物學(xué)報(bào)，2017，37（2）：175-184

［42］王恒波，祁舒婷，陳姝琦，等. 甘蔗栽培種單倍體基因組SSR位點(diǎn)的發(fā)掘與應(yīng)用［J］. 作物學(xué)報(bào)，2020，46（4）：631-642

［43］朱永群，彭丹丹，林超文，等. 蘇丹草轉(zhuǎn)錄組SSR分子標(biāo)記開(kāi)發(fā)及遺傳多樣性評(píng)價(jià)［J］. 草業(yè)學(xué)報(bào)，2018，27（5）：178-189

［44］陳渝，鄧潔，陳君愉，等. 番茄FAR1/FHY3轉(zhuǎn)錄因子家族的全基因組鑒定及表達(dá)分析［J］. 植物生理學(xué)報(bào)，2021，57（10）：1983-1995

［45］袁娜，王彤，劉廷利，等. 棉花FAR1/FHY3基因家族的全基因組分析［J］. 棉花學(xué)報(bào)，2018，30（1）：1-11

［46］羅文巧. 二穗短柄草根中木栓質(zhì)脂肪醇合成基因的克隆與功能驗(yàn)證［D］. 楊凌：西北農(nóng)林科技大學(xué)，2018：46-49

［47］陳柳君. bHLH轉(zhuǎn)錄因子在植物低溫脅迫中的研究進(jìn)展［J］. 江西農(nóng)業(yè)學(xué)報(bào)，2022，34（1）：79-87

［48］何江峰，王力偉，房永雨，等. 干旱脅迫和復(fù)水處理后梭梭轉(zhuǎn)錄因子的轉(zhuǎn)錄組分析［J］. 華北農(nóng)學(xué)報(bào)，2020，35（1）：36-43

［49］關(guān)淑艷，焦鵬，蔣振忠，等. MYB轉(zhuǎn)錄因子在植物非生物脅迫中的研究進(jìn)展［J］. 吉林農(nóng)業(yè)大學(xué)學(xué)報(bào)，2019，41（3）：253-260

（責(zé)任編輯 劉婷婷）