王曉陽(yáng)，彭振,2，邢愛(ài)雙，趙盈睿，馬欣麗，劉方,2，杜雄明,2，何守樸,2

亞洲棉短纖維發(fā)育相關(guān)長(zhǎng)鏈非編碼RNA的鑒定及表達(dá)

王曉陽(yáng)1，彭振1,2，邢愛(ài)雙1，趙盈睿1，馬欣麗1，劉方1,2，杜雄明1,2，何守樸1,2

1中國(guó)農(nóng)業(yè)科學(xué)院棉花研究所/棉花生物育種與綜合利用全國(guó)重點(diǎn)實(shí)驗(yàn)室，河南安陽(yáng) 455000；2鄭州大學(xué)農(nóng)學(xué)院/棉花生物育種與綜合利用全國(guó)重點(diǎn)實(shí)驗(yàn)室鄭州基地，鄭州 450001

【目的】長(zhǎng)鏈非編碼RNA（long non-coding RNAs，lncRNAs）是一類無(wú)蛋白質(zhì)編碼能力，但參與許多重要生命活動(dòng)調(diào)控過(guò)程的長(zhǎng)度大于200 nt的RNA。通過(guò)對(duì)亞洲棉無(wú)短纖維突變體（GA0149）和野生型（GA0146）纖維發(fā)育早期的轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行分析，挖掘調(diào)控短纖維發(fā)育的lncRNA，并明確其調(diào)控網(wǎng)絡(luò)，為進(jìn)一步解析棉花纖維發(fā)育機(jī)制奠定基礎(chǔ)?！痉椒ā窟x擇GA0146和GA0149 2個(gè)材料在開(kāi)花后當(dāng)天（0 DPA）及花后3 d（3 DPA）、5 d（5 DPA）和8 d（8 DPA）的胚珠和纖維為材料進(jìn)行轉(zhuǎn)錄組測(cè)序。鑒定lncRNA并預(yù)測(cè)其調(diào)控的靶基因；通過(guò)mRNA和lncRNA的差異表達(dá)分析，比較2個(gè)材料在不同纖維發(fā)育時(shí)期的差異。進(jìn)一步利用KOBAS軟件預(yù)測(cè)對(duì)差異lncRNA的靶基因進(jìn)行富集分析并預(yù)測(cè)其參與的生物過(guò)程；最后通過(guò)實(shí)時(shí)熒光定量（RT-qPCR）技術(shù)對(duì)25個(gè)差異表達(dá)的lncRNA轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行驗(yàn)證。【結(jié)果】共鑒定獲得15 339個(gè)lncRNA，其中11 595個(gè)lncRNA位于基因間區(qū)，包括2 428個(gè)反義lncRNA、350個(gè)內(nèi)含子lncRNA及966個(gè)正義lncRNA。共有1 932個(gè)差異表達(dá)lncRNA（DE-lncRNA），它們所對(duì)應(yīng)的8 134個(gè)靶基因中，有788個(gè)為差異表達(dá)基因（DE-mRNA）。KEGG代謝通路富集分析表明，DE-mRNA主要參與植物激素信號(hào)轉(zhuǎn)導(dǎo)（plant hormone signal transduction）和內(nèi)質(zhì)網(wǎng)中蛋白質(zhì)加工過(guò)程（protein processing in endoplasmic reticulum）。共表達(dá)調(diào)控網(wǎng)絡(luò)分析顯示，表達(dá)量差異比較顯著的lncRNA（）和其所調(diào)控的靶基因表達(dá)趨勢(shì)一致，僅在野生型（GA0146）短纖維發(fā)育早期胚珠中特異表達(dá)；而lncRNA（）與其調(diào)控的靶基因表達(dá)趨勢(shì)相反，在突變體（GA0149）中的表達(dá)量顯著高于野生型。RT-qPCR結(jié)果證實(shí)了轉(zhuǎn)錄組數(shù)據(jù)的真實(shí)性?！窘Y(jié)論】鑒定了26個(gè)與亞洲棉短纖維發(fā)育相關(guān)的lncRNA，其通過(guò)調(diào)控植物激素信號(hào)轉(zhuǎn)導(dǎo)途徑相關(guān)的吲哚乙酸合成酶基因（）和生長(zhǎng)素響應(yīng)蛋白基因（）的表達(dá)而影響短纖維的發(fā)育。

亞洲棉；短纖維突變體；長(zhǎng)鏈非編碼RNA；調(diào)控網(wǎng)絡(luò)；熒光定量PCR

0 引言

【研究意義】棉花纖維起始于種子表皮的一個(gè)單細(xì)胞，是紡織工業(yè)的主要原料。棉花纖維的產(chǎn)量主要取決于起始階段，該時(shí)期決定了胚珠纖維的數(shù)量。盡管所有棉花胚珠表皮細(xì)胞都具有分化成纖維的潛在能力，但最終只有30%的細(xì)胞能夠分化成纖維[1]。因此，可通過(guò)調(diào)節(jié)纖維起始過(guò)程表皮細(xì)胞數(shù)量來(lái)提高纖維產(chǎn)量。根據(jù)棉花成熟纖維長(zhǎng)度可以分為長(zhǎng)纖維（Lint）和短纖維（Fuzz），具有重要經(jīng)濟(jì)作用的長(zhǎng)纖維起始于開(kāi)花當(dāng)天（day post anthesis，0 DPA），最終長(zhǎng)度能達(dá)到2.5—3.5 cm，短纖維起始于開(kāi)花后3—5 DPA，長(zhǎng)度為5—10 mm[2]。雖然短纖維的起始要晚于長(zhǎng)纖維的起始，但是長(zhǎng)纖維和短纖維在分化過(guò)程可能經(jīng)歷相似的分子機(jī)制。因此，研究亞洲棉短纖維起始及其相關(guān)lncRNA的表達(dá)調(diào)控對(duì)解析棉花纖維發(fā)育具有重要的理論意義和實(shí)踐指導(dǎo)作用。【前人研究進(jìn)展】天然纖維突變體由于其特殊的性狀，是解析長(zhǎng)纖維和短絨發(fā)育機(jī)制的重要種質(zhì)資源[3-7]。目前，科學(xué)家對(duì)四倍體陸地棉纖維突變體進(jìn)行了大量研究，而對(duì)二倍體亞洲棉突變體的研究較少。與四倍體棉花相比，二倍體亞洲棉基因組結(jié)構(gòu)較為簡(jiǎn)單，因此，適合用于挖掘光籽性狀關(guān)鍵基因及闡明短絨發(fā)育的分子機(jī)制[8-10]。前期杜雄明課題組利用215份亞洲棉重測(cè)序數(shù)據(jù)進(jìn)行短纖維性狀的全基因組關(guān)聯(lián)分析以及轉(zhuǎn)基因試驗(yàn)，確定（）可能是控制亞洲棉短纖維發(fā)育的候選基因[11]。隨后，F(xiàn)eng等[12-13]利用QTL定位和轉(zhuǎn)錄組分析方法也鑒定出（）及編碼Fasciclin-like阿拉伯半乳聚糖蛋白18（FLA18）和轉(zhuǎn)運(yùn)蛋白的基因?yàn)榭刂贫汤w維性狀的候選基因。盡管已經(jīng)鑒定出許多參與纖維發(fā)育過(guò)程的基因，但是，控制纖維起始的整個(gè)調(diào)控網(wǎng)絡(luò)目前還未解析。長(zhǎng)鏈非編碼RNA（long noncoding RNAs，lncRNAs）是一類長(zhǎng)度大于200 bp，不具備蛋白編碼能力的RNA分子，與mRNA類似，都是由RNA聚合酶Ⅱ（RNA pol Ⅱ）進(jìn)行反轉(zhuǎn)錄，在結(jié)構(gòu)上具有5′加帽、3′聚腺苷酸化等特點(diǎn)[14-15]。與mRNA不同的是，lncRNA通常是從基因組的基因間區(qū)域轉(zhuǎn)錄的，而一些lncRNA來(lái)源于編碼基因的反義鏈[16]。lncRNA在不同植物的基因表達(dá)和維持染色體穩(wěn)定性方面起重要調(diào)控作用，參與廣泛的生物調(diào)節(jié)過(guò)程，例如擬南芥開(kāi)花過(guò)程、水稻光敏雄性不育、沙棘果維生素C合成、光誘導(dǎo)的蘋(píng)果種皮花青素合成，以及植物相應(yīng)逆境和非逆境生物脅迫等[17-26]。隨著多個(gè)棉花品種的高質(zhì)量基因組數(shù)據(jù)的發(fā)布[27-28]，前人已運(yùn)用RNA-seq和生物信息分析方法鑒定出一批棉花lncRNA。Wang等[29]和閆飛林[30]通過(guò)對(duì)海島棉中l(wèi)ncRNA進(jìn)行系統(tǒng)分析，鑒定出30 550個(gè)基因間lncRNA、4 718個(gè)反義lncRNA和35 268個(gè)反義lncRNA（lncNATs），經(jīng)組織特異性分析，發(fā)現(xiàn)約76%的同源1ncRNA表現(xiàn)出At或Dt亞組的偏好性表達(dá)模式；通過(guò)調(diào)查轉(zhuǎn)基因RNAi株系表型，發(fā)現(xiàn)lncRNA棉花轉(zhuǎn)基因株系的纖維均有變短趨勢(shì)，說(shuō)明lncRNA可能負(fù)調(diào)控棉纖維的起始和伸長(zhǎng)。通過(guò)對(duì)亞洲棉纖維起始期和快速伸長(zhǎng)期的轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行分析，發(fā)現(xiàn)51.9%的基因間lncRNA（lincRNAs）和54.5%的lncNAT傾向于在纖維發(fā)育的一個(gè)特定階段優(yōu)先表達(dá)[31]。通過(guò)分析Xuzhou142及其纖維突變體，以及由這兩個(gè)材料雜交產(chǎn)生的3個(gè)不同衣分產(chǎn)量株系0 DPA和5 DPA的轉(zhuǎn)錄組數(shù)據(jù)，鑒定出35 802個(gè)lncRNA，其中645個(gè)lncRNAs在中優(yōu)先表達(dá)，病毒誘導(dǎo)的基因沉默試驗(yàn)證實(shí)XLOC_545639和XLOC_039050能增加胚珠細(xì)胞起始數(shù)目[32]。通過(guò)對(duì)長(zhǎng)纖維突變體（）和野生型開(kāi)花0 DPA和8 DPA的轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行分析，鑒定出18 333個(gè)lncRNA，差異分析和GO富集通路分析發(fā)現(xiàn)，差異lncRNA的靶基因主要參與脂肪酸生物合成和延伸途徑[33]。同時(shí)，一些參與陸地棉鹽脅迫過(guò)程的lncRNA也被鑒定，擬南芥轉(zhuǎn)基因和病毒誘導(dǎo)的基因沉默試驗(yàn)證明，lncRNA973能提高陸地棉的耐鹽性[34]。目前，利用多策略RNA-seq技術(shù)對(duì)亞洲棉中的lncRNA進(jìn)行全面地功能注釋和鑒定，并對(duì)其轉(zhuǎn)錄調(diào)控機(jī)制進(jìn)行了解析，為理解棉花lncRNA的功能提供了寶貴研究資源[35]?！颈狙芯壳腥朦c(diǎn)】目前，盡管關(guān)于棉花短纖維候選基因的挖掘已經(jīng)取得了很大進(jìn)展，但其發(fā)育機(jī)制和調(diào)控網(wǎng)絡(luò)尚不清楚。lncRNA在調(diào)控植物生長(zhǎng)發(fā)育中起不可或缺的作用，但對(duì)其調(diào)控棉花纖維起始和伸長(zhǎng)的機(jī)制還知之甚少?！緮M解決的關(guān)鍵問(wèn)題】本研究通過(guò)分析亞洲棉無(wú)短纖維突變體（GA0149）和其野生型（GA0146）纖維發(fā)育早期的轉(zhuǎn)錄組數(shù)據(jù)，鑒定差異表達(dá)lncRNA，明確其功能性特征及表達(dá)水平，為進(jìn)一步研究lncRNA在棉花纖維發(fā)育中的作用提供依據(jù)。

1 材料與方法

1.1 試驗(yàn)材料

亞洲棉野生型中美棉GA0146既有長(zhǎng)纖維又有短纖維，其近等基因系短纖維突變體GA0149有長(zhǎng)纖維沒(méi)有短纖維，二者來(lái)自于中國(guó)農(nóng)業(yè)科學(xué)院棉花研究所國(guó)家棉花種質(zhì)資源中期庫(kù)，經(jīng)過(guò)大約連續(xù)10代自交，均沒(méi)有出現(xiàn)短纖維性狀的分離，其他農(nóng)藝性狀接近。所有材料均在中國(guó)農(nóng)業(yè)科學(xué)院棉花研究所安陽(yáng)試驗(yàn)田種植，并按照標(biāo)準(zhǔn)的田間生產(chǎn)管理方法。開(kāi)花前一天自交，同時(shí)掛牌記錄開(kāi)花時(shí)間，即為-1 DPA（開(kāi)花前一天）。在上午分別取每個(gè)材料0 DPA（開(kāi)花當(dāng)天）、3 DPA、5 DPA和8 DPA的胚珠或纖維，液氮速凍，-80 ℃冷凍保存，用于RNA的提取和轉(zhuǎn)錄組測(cè)序。每個(gè)時(shí)期設(shè)置3個(gè)生物學(xué)重復(fù)。

1.2 RNA提取和轉(zhuǎn)錄組測(cè)序

使用多糖多酚植物總RNA提取試劑盒（No. DP441）提取GA0146和GA0149短纖維發(fā)育早期0、3和5 DPA的胚珠及8 DPA的胚珠和纖維混合物，經(jīng)1.2%瓊脂糖凝膠電泳檢測(cè)完整性，用NanoDrop-2000（Thermo Scientific）分光度計(jì)測(cè)定濃度，然后用于RT-qPCR分析和cDNA文庫(kù)構(gòu)建。利用Illumina HiSeq2000平臺(tái)進(jìn)行雙末端測(cè)序，測(cè)序讀長(zhǎng)為PE150，獲得原始數(shù)據(jù)。

1.3 亞洲棉lncRNA的生物信息學(xué)預(yù)測(cè)和表達(dá)分析

將原始數(shù)據(jù)通過(guò)FastQC軟件去掉接頭和低質(zhì)量的數(shù)據(jù)，獲得1 064.24 Gb高質(zhì)量clean reads（Q＞20），各樣品Clean Data均達(dá)到16.45 Gb，Q30堿基百分比在86.17%及以上。利用TopHat2軟件[36]將clean reads與亞洲棉石系亞1號(hào)三代基因組[28]比對(duì)。利用HTSeq v0.6.1軟件與外顯子上的reads數(shù)進(jìn)行比對(duì)、計(jì)數(shù)[37]。通過(guò)Cufflinks 2.0程序組裝每個(gè)轉(zhuǎn)錄本。然后，利用Cuffdiff組裝來(lái)自BAM輸出文件的所有轉(zhuǎn)錄本，使用每千堿基轉(zhuǎn)錄物每百萬(wàn)映射讀取的片段（fragments per kilobase of transcript per million fragments mapped，F(xiàn)PKM）值來(lái)確定轉(zhuǎn)錄本豐度[38]。采用嵌入到CUFFLINKS軟件中的多次讀數(shù)功能和片段偏差校正方法來(lái)提高基因表達(dá)水平估計(jì)的準(zhǔn)確性。

根據(jù)起始位置，lncRNA可以分為基因內(nèi)lncRNA（intronic long non-coding RNAs，intronic-lncRNAs）、反義lncRNA（antisense long noncoding RNAs，antisense- lncRNAs）、正義lncRNA（sense lncRNAs）和基因間lncRNA（intergenic long non-coding RNAs，lincRNAs）等4種類型[39]。通過(guò)以下6個(gè)步驟鑒定lncRNA的類型：（1）丟棄在少于2個(gè)樣本中被檢測(cè)到的同時(shí)用2種方法組裝的轉(zhuǎn)錄本；（2）刪除覆蓋率小于轉(zhuǎn)錄本長(zhǎng)度一半的轉(zhuǎn)錄本；（3）去除由rRNA和tRNA延伸而來(lái)的轉(zhuǎn)錄本（cutoff E-value＜0.001）；（4）去除長(zhǎng)度小于200 bp的轉(zhuǎn)錄組；（5）利用CPC（Coding Potential Calculator）和CNCI（Coding-Non- Coding Index）軟件排除具有編碼能力的轉(zhuǎn)錄本[40]；（6）將剩余的轉(zhuǎn)錄本在Swiss-Prot和Pfam數(shù)據(jù)庫(kù)中進(jìn)行比對(duì)，以消除具有蛋白編碼能力和蛋白結(jié)構(gòu)域的轉(zhuǎn)錄本。

1.4 靶標(biāo)基因預(yù)測(cè)

為了研究lncRNA的功能，對(duì)其作用于相鄰基因的靶基因進(jìn)行預(yù)測(cè)。在lncRNA上下游10 kb內(nèi)搜索每個(gè)lncRNA周圍相鄰的編碼基因，研究它們的功能。通過(guò)基因表達(dá)水平來(lái)預(yù)測(cè)lncRNA和其靶基因的相互作用。使用R工具，應(yīng)用皮爾遜相關(guān)系數(shù)來(lái)闡明lncRNA和mRNA之間的相對(duì)表達(dá)關(guān)系，選擇皮爾遜相關(guān)系數(shù)（2）≥0.95。然后，利用lncRNA和mRNA的相關(guān)參數(shù)構(gòu)建共表達(dá)網(wǎng)絡(luò)，揭示lncRNA和mRNA之間的相互作用，并用Cytoscape軟件（3.8.2）將共表達(dá)網(wǎng)絡(luò)可視化；利用goseq和KOBAS軟件對(duì)lncRNA的靶基因功能進(jìn)行富集[41-42]。

1.5 差異表達(dá)基因（mRNA）和lncRNA分析

運(yùn)用“Deseq2”R軟件包對(duì)野生型和突變體同一發(fā)育時(shí)期的mRNA和lncRNA的FPKM和count值進(jìn)行差異表達(dá)分析，同時(shí)，過(guò)濾掉FPKM＜0.5的lncRNA[43]，篩選差異表達(dá)基因（differentially expressed mRNA，DE-mRNA）和lncRNA（DE-lncRNA）。將|log2(fold change)|＞1和-value＜0.05作為篩選DE-lncRNA的標(biāo)準(zhǔn)。

1.6 KEGG功能富集分析

為了解析差異基因的功能，將獲得的DE-lncRNA的靶基因在數(shù)據(jù)庫(kù)中進(jìn)行KEGG（kyoto encyclopedia of genes and genomes）富集通路分析。使用KOBAS軟件檢驗(yàn)KEGG途徑中DE基因或lncRNA靶基因的富集因子，選擇-value＜0.05的富集通路。

1.7 轉(zhuǎn)錄組數(shù)據(jù)的RT-qPCR分析

為了驗(yàn)證轉(zhuǎn)錄組數(shù)據(jù)的準(zhǔn)確性，挑選25個(gè)在4個(gè)時(shí)期共有的差異lncRNA進(jìn)行RT-qPCR驗(yàn)證。利用在線軟件NCBI-Primer（http://www.ncbi.nlm.nih.gov/ tools/primer-blast/）進(jìn)行l(wèi)ncRNA qRT-PCR引物設(shè)計(jì)（電子附表1）。使用EvoM-MLV預(yù)混型（AG11728）反轉(zhuǎn)錄試劑盒（含去除gDNA試劑，用于qPCR）將1 μg總RNA合成cDNA，以為內(nèi)參，進(jìn)行實(shí)時(shí)熒光定量PCR驗(yàn)證。每個(gè)基因至少設(shè)置3個(gè)技術(shù)重復(fù)，3個(gè)生物學(xué)重復(fù)。以GA0146材料0 DPA樣品的基因表達(dá)量為1，用2-ΔΔCT方法[44]分析突變體和野生型各個(gè)時(shí)期的lncRNA相對(duì)表達(dá)量。

2 結(jié)果

2.1 亞洲棉短纖維突變體（GA0149）和其野生型（GA0146）中的lncRNA鑒定

前期通過(guò)觀察亞洲棉短纖維起始早期的表型觀察，在3 DPA時(shí)突變體胚珠表皮未發(fā)現(xiàn)細(xì)胞突起，而野生型有大量細(xì)胞突起，說(shuō)明短纖維可能起始于3 DPA[45]。當(dāng)胚珠長(zhǎng)至成熟期，對(duì)其種子表型進(jìn)行調(diào)查，發(fā)現(xiàn)野生型GA0146成熟種子表皮同時(shí)附著短纖維和長(zhǎng)纖維，而其短纖維突變體GA0149成熟種子表皮則僅附著長(zhǎng)纖維（圖1-A）。為了進(jìn)一步了解造成種子短纖維表型差異的原因，選取亞洲棉野生型和其突變體短纖維發(fā)育關(guān)鍵時(shí)期（0 DPA、3 DPA、5 DPA和8 DPA）的胚珠及纖維通過(guò)Illumina平臺(tái)進(jìn)行全轉(zhuǎn)錄組測(cè)序，鑒定調(diào)控短纖維發(fā)育的關(guān)鍵lncRNA。將獲得數(shù)據(jù)與亞洲棉參考基因組比對(duì)，共獲得15 339個(gè)lncRNA。根據(jù)與編碼基因的相對(duì)位置，鑒定出11 595個(gè)基因間lncRNA、2 428個(gè)反義lncRNA、350個(gè)基因內(nèi)lncRNA和966個(gè)正義lncRNA（圖2-B）。

A：種子成熟期的亞洲棉野生型和短纖維突變體表型；B：4種類型lncRNA的分布

GC含量與四倍體棉花基因組中亞基因組間的基因偏向性表達(dá)有顯著相關(guān)性[29]。通過(guò)對(duì)lncRNA和mRNA的GC含量、轉(zhuǎn)錄本長(zhǎng)度以及外顯子數(shù)目進(jìn)行分析，與mRNA相比，lncRNA的GC含量更低（圖2-A）；同時(shí)，發(fā)現(xiàn)mRNA的平均外顯子長(zhǎng)度要大于lncRNA的外顯子長(zhǎng)度，lncRNA外顯子變化范圍為200—11 619 bp，mRNA的變化范圍為200—25 929 bp（圖2-B）；大部分lncRNAs含有2個(gè)外顯子，而mRNA含有多個(gè)外顯子（圖2-C）。lncRNA對(duì)mRNA的調(diào)控主要分為順式和反式調(diào)節(jié)，通過(guò)分析調(diào)控方式發(fā)現(xiàn)，大部分lncRNA對(duì)mRNA起順式調(diào)控作用，少部分為反式調(diào)節(jié)或者為2種方式的調(diào)節(jié)（圖2-D）。

A：lncRNA和mRNA中的GC含量；B：lncRNA和mRNA中的外顯子長(zhǎng)度分布情況；C：每個(gè)lncRNA和mRNA轉(zhuǎn)錄本的外顯子數(shù)目；D：lncRNA和mRNA之間的相互作用方式

2.2 亞洲棉短纖維起始期差異lncRNA的表達(dá)分析

為了鑒定出參與短纖維發(fā)育的候選lncRNA，對(duì)亞洲棉野生型和其突變體4個(gè)纖維發(fā)育時(shí)期（GA0149 vs GA0146 0 DPA、3 DPA、5 DPA和8 DPA）的lncRNA轉(zhuǎn)錄本表達(dá)量進(jìn)行分析，共鑒定出1 932個(gè)差異表達(dá)lncRNA，其中，612個(gè)來(lái)自于短纖維發(fā)育5 DPA，520個(gè)來(lái)自于長(zhǎng)纖維發(fā)育0 DPA，570個(gè)來(lái)自于短纖維發(fā)育3 DPA，440個(gè)來(lái)自于短纖維發(fā)育8 DPA，短纖維發(fā)育3和5 DPA共有的DE-lncRNA數(shù)目為167個(gè)（圖3和圖4-A）?；鹕綀D（圖3）顯示，4個(gè)比較組中，DE-lncRNA的基因表達(dá)譜，上調(diào)表達(dá)和下調(diào)表達(dá)的lncRNA具有明確的分布模式，表明不同比較組中存在不同的表達(dá)模式。4個(gè)時(shí)期共有的DE-lncRNA數(shù)目為65個(gè)（圖4-B、表1和電子附表2）。其中，有28個(gè)lncRNA在突變體中發(fā)生上調(diào)表達(dá)，37個(gè)lncRNA在突變體中發(fā)生了下調(diào)表達(dá)（表1）。說(shuō)明3和5 DPA共有的DE-lncRNA可能調(diào)控短纖維的起始和發(fā)育；在短纖維起始和發(fā)育階段，4個(gè)時(shí)期共有的65個(gè)DE-lncRNA的表達(dá)在2個(gè)材料中保持趨勢(shì)一致，同時(shí)，可以證明這些共有的DE-lncRNA可能調(diào)控長(zhǎng)纖維起始和發(fā)育。

2.3 與短纖維起始相關(guān)的差異lncRNA的功能分析

為了解析lncRNA的表達(dá)模式與短纖維發(fā)育的關(guān)系，分析野生型和突變體同一纖維發(fā)育時(shí)期的lncRNA表達(dá)情況（圖5-A—B）。在長(zhǎng)纖維起始期（0 DPA），突變體中有261個(gè)lncRNA發(fā)生了上調(diào)表達(dá)，259個(gè)lncRNA發(fā)生了下調(diào)表達(dá)，說(shuō)明有大量lncRNA參與了長(zhǎng)纖維的起始；在短纖維起始期（3 DPA），有316個(gè)lncRNA在突變體中發(fā)生了上調(diào)表達(dá)，254個(gè)lncRNA在突變體中發(fā)生了下調(diào)表達(dá)；在短纖維發(fā)育期（5 DPA），有300個(gè)lncRNA在突變體中發(fā)生了上調(diào)表達(dá)，312個(gè)lncRNA在突變體中發(fā)生了下調(diào)表達(dá)；在纖維伸長(zhǎng)期（8 DPA），有179個(gè)lncRNA在突變體中發(fā)生了上調(diào)表達(dá)，261個(gè)lncRNA在突變體中發(fā)生了下調(diào)表達(dá)。說(shuō)明在短纖維起始期發(fā)生上調(diào)表達(dá)或下調(diào)表達(dá)的差異lncRNA多于長(zhǎng)纖維起始期和伸長(zhǎng)期的差異lncRNA。

圖3 火山圖展示lncRNA表達(dá)量變化和FDR值之間的分布

表1 亞洲棉短纖維突變體（GA0149）和其野生型（GA0146）短纖維發(fā)育早期共有差異lncRNA及其表達(dá)量變化倍數(shù)

續(xù)表1 Continued table 1

A：2個(gè)材料4個(gè)纖維發(fā)育時(shí)期共有差異lncRNA的韋恩圖；B：65個(gè)共有的差異lncRNA表達(dá)量熱圖

為了解短纖維發(fā)育早期各時(shí)期lncRNA表達(dá)情況，比較2個(gè)材料相鄰2個(gè)時(shí)期lncRNA的表達(dá)變化（圖5-A—B）。在GA0146中，在長(zhǎng)纖維起始期至短纖維起始期（0—3 DPA），291個(gè)lncRNA發(fā)生上調(diào)表達(dá)，201個(gè)下調(diào)表達(dá)；在短纖維發(fā)育初期（3—5 DPA），254個(gè)lncRNA上調(diào)表達(dá)，189個(gè)下調(diào)表達(dá)；在短纖維發(fā)育初期至長(zhǎng)纖維伸長(zhǎng)期（5—8 DPA），334個(gè)lncRNA上調(diào)表達(dá)，426個(gè)下調(diào)表達(dá)。而在短纖維突變體GA0149中，從0—3 DPA時(shí)期，344個(gè)lncRNA上調(diào)表達(dá)，199個(gè)下調(diào)表達(dá)；在3—5 DPA時(shí)期，271個(gè)lncRNA上調(diào)表達(dá)，284個(gè)下調(diào)表達(dá)；在5—8 DPA時(shí)期，435個(gè)lncRNA上調(diào)表達(dá)，805個(gè)下調(diào)表達(dá)?？傊?，在長(zhǎng)纖維起始期至短纖維起始期，上調(diào)表達(dá)的lncRNA多于下調(diào)表達(dá)的lncRNA，說(shuō)明大部分lncRNA可能正調(diào)控長(zhǎng)纖維起始；在短纖維發(fā)育初期至長(zhǎng)纖維伸長(zhǎng)期，差異lncRNA較少；在纖維伸長(zhǎng)期，差異lncRNA最多，并且下調(diào)表達(dá)的lncRNA多于上調(diào)表達(dá)的lncRNA，表明大部分lncRNA可能負(fù)調(diào)控纖維伸長(zhǎng)，同時(shí)，由于伸長(zhǎng)期的纖維細(xì)胞生長(zhǎng)最快，并且該時(shí)期主要影響纖維長(zhǎng)度、強(qiáng)度和細(xì)度等品質(zhì)性狀，因此，在纖維伸長(zhǎng)期，lncRNA的表達(dá)水平可能發(fā)生了較大變化。

為了研究lncRNA的功能，根據(jù)在染色體上的位置，預(yù)測(cè)差異lncRNA所調(diào)控的靶基因。共篩選到8 134個(gè)mRNA，同時(shí)按照差異基因分析標(biāo)準(zhǔn)，對(duì)mRNA的表達(dá)量進(jìn)行計(jì)算以及差異性分析，最終鑒定到788個(gè)差異表達(dá)的靶基因。經(jīng)KEGG代謝通路富集顯示（圖5-C），差異靶基因主要參與的通路如：植物激素信號(hào)轉(zhuǎn)導(dǎo)（plant hormone signal transduction）和內(nèi)質(zhì)網(wǎng)中蛋白質(zhì)加工過(guò)程（protein processing in endoplasmic reticulum）。對(duì)富集到植物激素信號(hào)轉(zhuǎn)導(dǎo)通路中的14個(gè)差異基因表達(dá)量進(jìn)行分析（表2），發(fā)現(xiàn)8個(gè)差異基因參與生長(zhǎng)素響應(yīng)和生物合成過(guò)程，與野生型相比，生長(zhǎng)素響應(yīng)蛋白基因在突變體的開(kāi)花當(dāng)天上調(diào)表達(dá)，而在3、5和8 DPA表達(dá)水平?jīng)]有發(fā)生變化；吲哚乙酸合成酶基因在突變體的5 DPA表達(dá)量上升，其他時(shí)期的表達(dá)無(wú)顯著差異；生長(zhǎng)素響應(yīng)蛋白編碼基因和，在突變體的5 DPA表達(dá)量上升，8 DPA表達(dá)量下降，其他時(shí)期表達(dá)量沒(méi)有變化，表明生長(zhǎng)素對(duì)短纖維的起始可能沒(méi)有影響，對(duì)纖維的伸長(zhǎng)起正調(diào)控作用。參與茉莉酸甲酯合成途徑的轉(zhuǎn)錄因子MYC基因（）在突變體5 DPA的表達(dá)量發(fā)生了顯著上調(diào)表達(dá)，說(shuō)明該轉(zhuǎn)錄因子可能負(fù)調(diào)控短纖維的起始。參與乙烯合成途徑的基因以及參與脫落酸（ABA）合成途徑的基因在突變體的5 DPA表達(dá)量上升，8 DPA表達(dá)量顯著下降，說(shuō)明過(guò)量的乙烯和ABA能抑制短纖維起始，少量的乙烯能促進(jìn)纖維伸長(zhǎng)。同時(shí)，也發(fā)現(xiàn)水楊酸（SA）和赤霉素（GA）合成途徑的基因也參與了纖維發(fā)育過(guò)程。說(shuō)明植物激素在纖維起始和伸長(zhǎng)過(guò)程中起重要作用，這些可能的靶基因?yàn)閘ncRNA在棉花纖維發(fā)育中的作用提供了新見(jiàn)解。

A：短纖維發(fā)育早期時(shí)期上調(diào)表達(dá)和下調(diào)表達(dá)差異lncRNA數(shù)目，紅色代表上調(diào)的lncRNA，藍(lán)色代表下調(diào)的lncRNA；B：4個(gè)纖維發(fā)育時(shí)期差異lncRNA表達(dá)量熱圖；C：差異lncRNA的靶基因KEGG富集通路

表2 植物激素信號(hào)通路的差異靶基因及其表達(dá)情況

2.4 lncRNA和其靶基因的調(diào)控網(wǎng)絡(luò)分析

lncRNA的一個(gè)主要功能是順式調(diào)控同一等位基因上的相鄰基因，從而調(diào)節(jié)基因的轉(zhuǎn)錄或轉(zhuǎn)錄后表達(dá)。因?yàn)轫樖阶饔玫膌ncRNA靶向鄰近基因，在所有已鑒定的lncRNA上、下游10 kb區(qū)域搜索了編碼基因，同時(shí)預(yù)測(cè)它們的功能。通過(guò)分析具有順式調(diào)控作用的lncRNA，構(gòu)建lncRNA和mRNA的共表達(dá)網(wǎng)絡(luò)。大多數(shù)mRNA和lncRNA是一對(duì)一的匹配。然而，大部分lncRNA也可以調(diào)控多個(gè)mRNA。考慮到表格無(wú)法顯示lncRNA和mRNA之間的大量網(wǎng)絡(luò)信息，選擇了一些mRNA與lncRNA以繪制它們之間的網(wǎng)絡(luò)圖（圖6）。結(jié)果顯示，lncRNA（）可以同時(shí)調(diào)控12個(gè)mRNA的表達(dá)，而則調(diào)控18個(gè)mRNA的表達(dá)；同時(shí)發(fā)現(xiàn)一個(gè)mRNA也可以受多個(gè)lncRNA的調(diào)控，如6個(gè)lncRNA調(diào)控編碼基因的表達(dá)，則受22個(gè)lncRNA的調(diào)控。

位于蛋白質(zhì)編碼基因上游的lncRNA可能與啟動(dòng)子區(qū)域或者順式調(diào)節(jié)元件有重疊，并可能在轉(zhuǎn)錄或轉(zhuǎn)錄后水平調(diào)節(jié)其附近基因的表達(dá)。位于蛋白質(zhì)編碼基因下游的lncRNA可以從基因的3′UTR或下游區(qū)域啟動(dòng)轉(zhuǎn)錄，并可能參與基因間相互調(diào)節(jié)作用。為了進(jìn)一步解析差異lncRNA對(duì)其下游靶基因表達(dá)的影響，隨機(jī)選擇4個(gè)差異表達(dá)明顯的lncRNA及其所對(duì)應(yīng)的靶基因進(jìn)行網(wǎng)絡(luò)調(diào)控分析（圖7）。lncRNA（）僅在野生型短纖維發(fā)育早期表達(dá)，在短纖維突變體的各個(gè)時(shí)期都不表達(dá)，而其所調(diào)控的靶基因（功能未知蛋白）、（功能未知蛋白）和（含五肽重復(fù)序列的蛋白）也是在野生型中顯示高表達(dá)，在突變體中幾乎不表達(dá)；（功能未知蛋白）和（功能未知蛋白）則與lncRNA的表達(dá)趨勢(shì)相反，在突變體中的表達(dá)水平高于其在野生型中的表達(dá)水平，說(shuō)明MSTRG.454250.3通過(guò)正調(diào)控、和的表達(dá)，負(fù)調(diào)控和的表達(dá)，進(jìn)而抑制纖維的發(fā)育。僅在野生型中表達(dá)的lncRNAs（）與其下游靶基因的表達(dá)趨勢(shì)一致，正調(diào)控下游靶基因的表達(dá)。正調(diào)控下游靶基因的表達(dá)，在野生型中的表達(dá)量高于突變體中的表達(dá)量。與上述3個(gè)lncRNA的表達(dá)趨勢(shì)相反，在突變體中的表達(dá)量高于野生型中的表達(dá)量，同時(shí)正調(diào)控下游靶基因的表達(dá)，說(shuō)明該lncRNA在短纖維起始和發(fā)育過(guò)程起負(fù)調(diào)控作用。

圖6 lncRNA與mRNA之間的相互作用關(guān)系

2.5 參與短纖維起始的lncRNA表達(dá)量驗(yàn)證

長(zhǎng)鏈非編碼RNA在許多植物的生長(zhǎng)、發(fā)育和適應(yīng)機(jī)制中發(fā)揮著重要作用，如花發(fā)育、有性繁殖、果實(shí)成熟、纖維發(fā)育、生物脅迫和非生物脅迫反應(yīng)。為了揭示lncRNA在亞洲棉短纖維發(fā)育中的作用，結(jié)合RNA-seq中的lncRNA表達(dá)水平，隨機(jī)選擇25個(gè)lncRNA在亞洲棉野生型（GA0146）和其短纖維突變體（GA0149）的0、3、5和8 DPA胚珠及纖維中進(jìn)行熒光定量PCR驗(yàn)證。轉(zhuǎn)錄組數(shù)據(jù)與RT-qPCR結(jié)果一致（圖8），一些lncRNA（如、、、和）傾向于在突變體4個(gè)纖維發(fā)育時(shí)期（GA0149）中高表達(dá)，在野生型（GA0146）中不表達(dá)，說(shuō)明這些lncRNA可能在棉花短纖維發(fā)育中起負(fù)調(diào)控作用。另外一些lncRNA（如、、、和）在亞洲棉野生型（GA0146）的0、3、5和8 DPA時(shí)期高于其在突變體中的表達(dá)量，說(shuō)明這些lncRNA正調(diào)控棉花短纖維的發(fā)育。

圖7 亞洲棉野生型和短纖維突變體4個(gè)纖維發(fā)育時(shí)期的差異lncRNA和其差異靶基因的共調(diào)控網(wǎng)絡(luò)

2.6 二倍體和四倍體棉花中l(wèi)ncRNA的鑒定

通過(guò)比較二倍體和四倍體棉種間lncRNA序列以及使用共線分析方法，在高置信水平上鑒定具有保守性低的同源lncRNA位點(diǎn)，找出亞洲棉和陸地棉中共有及特有的lncRNA，結(jié)果顯示，在2個(gè)棉種中共有的lncRNA有9 260個(gè)，陸地棉特有的lncRNA為24 856個(gè)，亞洲棉特有的lncRNA為6 076個(gè)。為進(jìn)一步探討lncRNA在2個(gè)棉種中對(duì)纖維發(fā)育機(jī)制的調(diào)控異同，對(duì)這些共有和特有的lncRNA進(jìn)行KEGG富集通路分析。發(fā)現(xiàn)亞洲棉特有的lncRNA主要參與植物激素信號(hào)轉(zhuǎn)導(dǎo)途徑、糖類代謝及代謝、苯丙烷生物合成及代謝途徑；2個(gè)棉種共有的lncRNA主要參與花生四烯酸代謝、植物激素信號(hào)轉(zhuǎn)導(dǎo)和內(nèi)黃酮生物合成途徑等；陸地棉共有的lncRNA主要參與肌醇磷酸鹽代謝、脂質(zhì)代謝、氨基酸降解、蛋白質(zhì)加工和角質(zhì)、軟木質(zhì)和蠟質(zhì)生物合成途徑等（圖9）。

3 討論

纖維起始和快速伸長(zhǎng)階段是影響棉花纖維產(chǎn)量和質(zhì)量的關(guān)鍵步驟。因此，了解棉纖維發(fā)育過(guò)程中的各種分子途徑至關(guān)重要。在過(guò)去的10年里，通過(guò)遺傳學(xué)篩選，挖掘出許多調(diào)控棉纖維發(fā)育的基因[46]。然而，調(diào)控胚珠表皮分化和纖維伸長(zhǎng)過(guò)程的機(jī)制目前還尚不清楚。本研究系統(tǒng)地比較分析了2個(gè)短纖維性狀差異顯著的亞洲棉lncRNA。通過(guò)分析胚珠和纖維起始、分化和伸長(zhǎng)階段的轉(zhuǎn)錄組數(shù)據(jù)，鑒定到一批與短纖維發(fā)育相關(guān)的新lncRNA。發(fā)現(xiàn)差異lncRNA的靶基因主要富集在植物激素信號(hào)轉(zhuǎn)導(dǎo)和蛋白質(zhì)加工過(guò)程。建立了lncRNA-mRNA調(diào)控短纖維起始和形成的網(wǎng)絡(luò)圖，首次發(fā)現(xiàn)大量lncRNA在短纖維發(fā)育階段特異性表達(dá)。這些數(shù)據(jù)為進(jìn)一步研究lncRNA調(diào)控棉花纖維發(fā)育提供了良好的基礎(chǔ)。

A：熱圖顯示亞洲棉野生型（GA0146）和突變體（GA0149）lncRNA的轉(zhuǎn)錄差異；B：25個(gè)lncRNA在亞洲棉野生型和短纖維突變體的熒光定量結(jié)果，平均值±SD（n=3），*P＜0.05，**P＜0.01，****P＜0.0001

A：lncRNA同源基因在亞洲棉（Ga）和陸地棉（Gh）中的分布；B：lncRNA同源基因在亞洲棉和陸地棉中的富集通路比較

3.1 亞洲棉短纖維發(fā)育過(guò)程中的lncRNA鑒定

目前，lncRNA已經(jīng)在許多物種中鑒定出來(lái)，在植物生長(zhǎng)發(fā)育和逆境脅迫過(guò)程起重要作用。lncRNA通過(guò)調(diào)控基因的甲基化水平進(jìn)而影響水稻的育性[47]，參與擬南芥的非生物脅迫和低溫春化過(guò)程[21]，參與玉米的干旱脅迫和與赤霉素相關(guān)的植物激素途徑[48]。此外lncRNA也參與棉花干旱和鹽脅迫等非生物逆境過(guò)程以及棉纖維的發(fā)育過(guò)程[31, 49-50]。上述結(jié)果說(shuō)明，lncRNA是植物中參與轉(zhuǎn)錄調(diào)控過(guò)程的一類重要媒介。另外，有研究報(bào)道，在轉(zhuǎn)錄水平上，lncRNA通過(guò)干擾編碼mRNA和非編碼RNA的轉(zhuǎn)錄、干擾蛋白質(zhì)結(jié)合復(fù)合物的形成，通過(guò)順式作用元件調(diào)節(jié)鄰近基因的表達(dá)；在轉(zhuǎn)錄后水平上lncRNA通過(guò)降解mRNA、調(diào)控翻譯后修飾和可變剪接過(guò)程進(jìn)而影響植物的生長(zhǎng)發(fā)育。盡管越來(lái)越多的報(bào)告表明lncRNA在哺乳動(dòng)物和植物發(fā)育過(guò)程中發(fā)揮作用，但在植物中鑒定的lncRNA尚處于起步階段，僅有少部分lncRNA被明確鑒定在調(diào)節(jié)植物發(fā)育過(guò)程中起重要作用[39]。本研究通過(guò)對(duì)亞洲棉4個(gè)纖維發(fā)育時(shí)期的RNA-Seq數(shù)據(jù)進(jìn)行分析，鑒定出15 339個(gè)lncRNA，大部分lncRNA為基因間lncRNA（lincRNAs）。為了研究亞洲棉中l(wèi)ncRNAs的結(jié)構(gòu)特點(diǎn)，對(duì)lncRNA和蛋白編碼基因（protein coding genes，PCGs）外顯子長(zhǎng)度、數(shù)目、GC含量和表達(dá)水平進(jìn)行分析，結(jié)果發(fā)現(xiàn)，與編碼基因相比，lncRNA含有較低的GC含量、較長(zhǎng)的外顯子長(zhǎng)度和較低的表達(dá)水平[35]。本研究序列結(jié)構(gòu)分析發(fā)現(xiàn)了同樣的結(jié)果，lncRNA序列長(zhǎng)度短于編碼基因的序列長(zhǎng)度，lncRNA中的GC含量低于mRNA中的GC含量。lncRNA主要通過(guò)順式和反式作用機(jī)制調(diào)控相關(guān)基因表達(dá)，反式作用的lncRNA可以作為染色質(zhì)的信號(hào)，遠(yuǎn)距離地調(diào)節(jié)相同染色體域甚至不同染色體中的靶基因表達(dá)；順式作用lncRNA位于編碼蛋白的上游和下游，與啟動(dòng)子或共表達(dá)基因的順式核心元件相互作用，從而在轉(zhuǎn)錄或轉(zhuǎn)錄后水平調(diào)節(jié)基因表達(dá)。本研究根據(jù)lncRNA與mRNA的物理位置和表達(dá)關(guān)系預(yù)測(cè)了其潛在的順式和反式靶基因，發(fā)現(xiàn)大部分lncRNA對(duì)其靶基因的調(diào)控方式為順式調(diào)節(jié)。Zheng等[35]同樣在亞洲棉中發(fā)現(xiàn)許多l(xiāng)ncRNA以順式作用方式調(diào)控鄰近編碼蛋白基因的表達(dá)。上述結(jié)果說(shuō)明亞洲棉lncRNA具有獨(dú)特的結(jié)構(gòu)特點(diǎn)。

3.2 差異表達(dá)lncRNA在亞洲棉纖維發(fā)育過(guò)程中的作用

通常認(rèn)為lncRNA在各種物種中具有高度的組織特異性，lncRNA的這一特性可能意味著它們可以在維持組織特性及組織發(fā)育和分化中發(fā)揮重要作用。為了進(jìn)一步了解lncRNA對(duì)亞洲棉野生型和短纖維突變體材料的短纖維發(fā)育影響，分析了這兩個(gè)材料4個(gè)短纖維發(fā)育早期的基因表達(dá)譜。結(jié)果顯示，大部分差異表達(dá)的lncRNA表現(xiàn)出組織特異性表達(dá)趨勢(shì)，暗示在纖維發(fā)育過(guò)程中l(wèi)ncRNA出現(xiàn)了功能分化。研究報(bào)道，lncRNA的表達(dá)具有時(shí)空特異性、高度保守性，在棉花胚珠發(fā)育過(guò)程中起重要作用[51]。本研究發(fā)現(xiàn)、和等僅在野生型的胚珠中特異性表達(dá)，在突變體中不表達(dá)，說(shuō)明這些lncRNA對(duì)短纖維的起始可能起正調(diào)控作用。和在突變體胚珠中的表達(dá)水平顯著高于其在野生型中的表達(dá)水平，說(shuō)明這兩個(gè)lncRNA可能抑制亞洲棉短纖維的起始。迄今為止，大多數(shù)lncRNA的功能還未完全解析。通過(guò)分析lncRNA對(duì)編碼基因的調(diào)控方式，構(gòu)建lncRNA-mRNA共表達(dá)網(wǎng)絡(luò)，以進(jìn)一步鑒定lncRNA與mRNA之間的相互調(diào)控關(guān)系，有助于lncRNA的功能預(yù)測(cè)[52]。大部分lncRNA與mRNA屬于一對(duì)多的調(diào)控方式，同樣的一個(gè)編碼基因也可以受多個(gè)lncRNAs的調(diào)控。例如，、和正調(diào)控其下游靶基因的表達(dá)，進(jìn)而影響短纖維的發(fā)育。上述結(jié)果說(shuō)明lncRNA可能通過(guò)抑制或激活其下游靶基因的表達(dá)，進(jìn)而影響亞洲棉短纖維發(fā)育過(guò)程。但是這些lncRNA的功能研究需要進(jìn)一步通過(guò)轉(zhuǎn)基因和病毒誘導(dǎo)的基因沉默等試驗(yàn)驗(yàn)證。

3.3 參與植物激素信號(hào)轉(zhuǎn)導(dǎo)過(guò)程的lncRNA調(diào)控亞洲棉短纖維的發(fā)育

植物激素生長(zhǎng)素、乙烯、赤霉素和茉莉酸甲酯在纖維發(fā)育過(guò)程起重要作用。近年來(lái)，許多研究表明，植物主要生長(zhǎng)素吲哚-3-乙酸（IAA）通過(guò)在纖維細(xì)胞中建立濃度梯度來(lái)決定棉花纖維的起始[53]。在纖維起始階段，利用花結(jié)合蛋白基因7（）啟動(dòng)子驅(qū)動(dòng)IAA生物合成基因在棉花中特異性表達(dá)，顯著增加了轉(zhuǎn)基因株系胚珠表皮中IAA的積累和起始纖維細(xì)胞數(shù)量，進(jìn)而使棉纖維的產(chǎn)量得到了提高[54]。吲哚乙酸誘導(dǎo)蛋白16基因（）在野生型的胚珠中表達(dá)量非常低，而在棉纖維突變體（）開(kāi)花后的胚珠中表達(dá)量迅速升高，說(shuō)明該基因可能負(fù)調(diào)控纖維的起始和伸長(zhǎng)[55]。本研究通過(guò)對(duì)差異lncRNA的靶基因進(jìn)行功能富集發(fā)現(xiàn)，大部分靶基因都匯集到植物激素信號(hào)轉(zhuǎn)導(dǎo)通路，吲哚乙酸合成酶基因（）和生長(zhǎng)素響應(yīng)蛋白基因（）在短纖維突變體5 DPA的表達(dá)量高于野生型中的表達(dá)量，說(shuō)明這兩個(gè)基因可能負(fù)調(diào)控短纖維的起始。茉莉酸甲酯（JA）在纖維發(fā)育過(guò)程也發(fā)揮著一定作用，例如當(dāng)用JA的抑制劑GhJAZ2在棉花中超表達(dá)后，顯著抑制了長(zhǎng)纖維和短纖維的起始和伸長(zhǎng)[56]。向體外離體培養(yǎng)的胚珠中加入2.5 mmol·L-1的JA，能顯著抑制纖維的伸長(zhǎng)，而當(dāng)加入0.05 mmol·L-1的JA時(shí)，纖維長(zhǎng)度顯著增加。說(shuō)明體外離體培養(yǎng)時(shí)，高濃度的JA能抑制纖維起始，而適當(dāng)?shù)蜐舛鹊腏A能促進(jìn)纖維伸長(zhǎng)。這與本研究結(jié)果相符，如參與JA途徑的轉(zhuǎn)錄因子MYC4（）和茉莉酸甲酯合成酶基因JAR1（）在突變體5DPA表達(dá)量高，8 DPA表達(dá)量低。通過(guò)對(duì)棉花胚珠中ABA含量進(jìn)行測(cè)定發(fā)現(xiàn)，長(zhǎng)纖維突變體開(kāi)花當(dāng)天（0 DPA）胚珠中ABA的含量高于野生型胚珠中的ABA含量，說(shuō)明ABA可能作為一個(gè)負(fù)調(diào)節(jié)因子參與棉纖維的發(fā)育[57]。本研究發(fā)現(xiàn)ABA信號(hào)因子受體ABI2（）在突變體5 DPA中的表達(dá)量顯著高于在野生型中的表達(dá)量，這與前人結(jié)果相符合。但是上述結(jié)果需要進(jìn)一步的試驗(yàn)證明。

3.4 lncRNA調(diào)控亞洲棉和陸地棉纖維機(jī)制的異同

異源四倍體栽培種棉陸地棉（，A1A1D1D1）和海島棉（，A2A2D2D2）是由二倍體A基因組棉（亞洲棉或者草棉）與二倍體D基因組（雷蒙德氏棉）通過(guò)天然雜交加倍而形成。后來(lái)又經(jīng)歷大約五千年的獨(dú)立馴化和人工選擇，以及環(huán)境適應(yīng)，形成現(xiàn)代棉花栽培品種[58]。二倍體棉種和四倍體棉種的棉纖維性狀不相同，亞洲棉纖維粗短且產(chǎn)量低，陸地棉具有纖維品質(zhì)好、產(chǎn)量高和適應(yīng)性強(qiáng)等特點(diǎn)，是主要的栽培種。為解析異源四倍體陸地棉和二倍體亞洲棉纖維發(fā)育機(jī)制的異同。Zhao等[59]從表觀遺傳學(xué)方面對(duì)亞洲棉、雷蒙德氏棉、陸地棉及由亞洲棉和雷蒙德氏棉雜交F1的胚珠和葉片中的lncRNA鑒定，發(fā)現(xiàn)80% lncRNA在異源四倍體棉花基因組中等位表達(dá)，基因組沖擊會(huì)產(chǎn)生新的非編碼轉(zhuǎn)錄本。通過(guò)對(duì)二倍體亞洲棉和四倍體陸地棉lncRNA進(jìn)行鑒定，顯示由二倍體進(jìn)化到四倍體的過(guò)程可能由于基因組重排產(chǎn)生了新的lncRNA[60]。本研究通過(guò)對(duì)2個(gè)棉種lncRNA進(jìn)行富集通路分析發(fā)現(xiàn)，亞洲棉及2個(gè)棉種共有的lncRNA主要富集在植物信號(hào)轉(zhuǎn)導(dǎo)途徑，陸地棉lncRNA主要參與脂肪酸合成代謝途徑。說(shuō)明由于地理位置的不同及基因組重組和自然選擇的關(guān)系，在多倍化形成過(guò)程中l(wèi)ncRNA產(chǎn)生了新的功能分化。

總之，利用生物信息學(xué)方法，鑒定出許多與亞洲棉花纖維發(fā)育相關(guān)的lncRNA，為未來(lái)研究lncRNA在棉花纖維中的調(diào)控作用和功能分析提供了數(shù)據(jù)基礎(chǔ)。未來(lái)的工作將旨在解析lncRNA與棉花纖維發(fā)育相關(guān)的生物學(xué)功能及其提高纖維產(chǎn)量和品質(zhì)的遺傳學(xué)機(jī)制。

4 結(jié)論

從亞洲棉短纖維突變體（GA0149）和其野生型（GA0146）4個(gè)短纖維發(fā)育早期（0、3、5和8 DPA）的轉(zhuǎn)錄組數(shù)據(jù)中共鑒定出15 329個(gè)lncRNA和1 932個(gè)具有差異表達(dá)的lncRNA。同時(shí)，發(fā)現(xiàn)lncRNA的GC含量低于mRNA的GC含量、序列長(zhǎng)度小于mRNA、表達(dá)方式具有組織特異性。單個(gè)lncRNA可以同時(shí)調(diào)控多個(gè)靶基因的表達(dá)，而單個(gè)mRNA也可以受多個(gè)lncRNA的調(diào)控。差異表達(dá)的lncRNA通過(guò)調(diào)控植物激素信號(hào)轉(zhuǎn)導(dǎo)通路基因的表達(dá)，影響亞洲棉短纖維的發(fā)育。大部分差異表達(dá)的lncRNA的靶基因都匯集到植物激素信號(hào)轉(zhuǎn)導(dǎo)通路，其中，吲哚乙酸合成酶基因（）和生長(zhǎng)素響應(yīng)蛋白基因（）可能抑制纖維的起始。

[1] KIM H J, TRIPLETT B A. Cotton fiber growth in planta and in vitro. Models for plant cell elongation and cell wall biogenesis. Plant Physiology, 2001, 127(4): 1361-1366.

[2] STEWART J M. Fiber initiation on the cotton ovule (). American Journal of Botany, 1975, 62(7): 723-730.

[3] WU H T, TIAN Y, WAN Q, FANG L, GUAN X Y, CHEN J D, HU Y, YE W X, ZHANG H, GUO W Z, CHEN X Y, Zhang T Z. Genetics and evolution of MIXTA genes regulating cotton lint fiber development. The New Phytologist, 2018, 217(2): 883-895.

[4] WAN Q, GUAN X Y, YANG N N, WU H T, PAN M Q, LIU B L, FANG L, YANG S P, HU Y, YE W X, ZHANG H, MA P Y, CHEN J D, WANG Q, MEI G F, CAI C P, YANG D L, WANG J W, GUO W Z, ZHANG W H, CHEN X Y, ZHANG T Z. Small interfering RNAs from bidirectional transcripts ofregulate cotton fiber development. The New Phytologist, 2016, 210(4): 1298-1310.

[5] CHEN W, LI Y, ZHU S H, FANG S T, ZHAO L J, GUO Y, WANG J Y, YUAN L, LU Y J, LIU F, YAO J B, ZHANG Y S. A retrotransposon insertion inis likely responsible for the lintless locus liof tetraploid cotton. Frontiers in Plant Science, 2020, 11: 593679.

[6] NAOUMKINA M, THYSSEN G N, FANG D D, LI P, FLORANE C B. Elucidation of sequence polymorphism in fuzzless-seed cotton lines. Molecular Genetics and Genomics, 2021, 296(1): 193-206.

[7] ZHU Q H, STILLER W, MONCUQUET P, GORDON S, YUAN Y M, BARNES S, WILSON I. Genetic mapping and transcriptomic characterization of a new fuzzless-tufted cottonseed mutant. G3 Genes|Genomes|Genetics, 2021, 11(1): jkaa042.

[8] DING M Q, CAO Y F, HE S E, SUN J, DAI H Q, ZHANG H, SUN C D, JIANG Y R, PATERSON A H, RONG J K., a candidate gene for theSMA-4 mutant, promotes trichome and fiber initiation by cellular H2O2and Ca2+signals. Plant Molecular Biology,2020, 103(4/5): 409-423.

[9] LIU X Y, MONCUQUET P, ZHU Q H, STILLER W, ZHANG Z S, WILSON I. Genetic identification and transcriptome analysis of lintless and fuzzless traits inL.. International Journal of Molecular Sciences, 2020, 21(5): 1675-1696.

[10] DU S J, DONG C J, ZHANG B, LAI T F, DU X M, LIU J Y. Comparative proteomic analysis reveals differentially expressed proteins correlated with fuzz fiber initiation in diploid cotton (L.). Journal of Proteomics, 2013, 82: 113-129.

[11] WANG X Y, MIAO Y C, CAI Y F, SUN G F, JIA Y H, SONG S, PAN Z E, ZHANG Y M, WANG L Y, FU G Y, GAO Q, JI G X, WANG P P, CHEN B J, PENG Z, ZHANG X M, WANG X, DING Y, HU D W, GENG X L, WANG L R, PANG B Y, GONG W F, HE S P, DU X M. Large-fragment insertion activates gene() and is associated with the fuzz and trichome reduction in cotton (). Plant Biotechnology Journal, 2021, 19(6): 1110-1124.

[12] FENG X X, CHENG H L, ZUO D Y, ZHANG Y P, WANG Q L, LIU K, ASHRAF J, YANG Q H, LI S M, CHEN X Q, SONG G L. Fine mapping and identification of the fuzzless genein DPL972 (). Theoretical and Applied Genetics, 2019, 132(8): 2169-2179.

[13] FENG X X, LIU S, CHENG H L, ZUO D Y, ZHANG Y P, WANG Q L, LV L M, SONG G L. Weighted gene co-expression network analysis reveals hub genes contributing to fuzz development in. Genes, 2021, 12(5): 753-768.

[14] LIU J, JUNG C, XU J, WANG H, DENG S L, BERNAD L, ARENAS-HUERTERO C, CHUA N H. Genome-wide analysis uncovers regulation of long intergenic noncoding RNAs in. The Plant Cell, 2012, 24(11): 4333-4345.

[15] GUTTMAN M, RUSSELL P, INGOLIA N T, WEISSMAN J S, LANDER E S. Ribosome profiling provides evidence that large noncoding RNAs do not encode proteins.Cell, 2013, 154(1): 240-251.

[16] PONTING C P, OLIVER P L, REIK W. Evolution and functions of long noncoding RNAs. Cell, 2009, 136(4): 629-641.

[17] WANG H L V, CHEKANOVA J A. Long noncoding RNAs in plants. Advances in Experimental Medicine and Biology. Singapore: Springer Singapore, 2017, 1008: 133-154.

[18] SUN X, ZHENG H X, SUI N. Regulation mechanism of long noncoding RNA in plant response to stress. Biochemical and Biophysical Research Communications, 2018, 503(2): 402-407.

[19] ZHANG X P, SHEN J, XU Q J, DONG J, SONG L R, WANG W, SHEN F F. Long noncoding RNA lncRNA354 functions as a competing endogenous RNA of miR160b to regulategenes in response to salt stress in upland cotton. Plant, Cell & Environment, 2021, 44(10): 3302-3321.

[20] ZHANG L, LIU J L, CHENG J R, SUN Q, ZHANG Y, LIU J G, LI H M, ZHANG Z, WANG P, CAI C W, CHU Z Y, ZHANG X, YUAN Y L, SHI Y Z, CAI Y F.andmodulate cell wall defense genes to regulate cotton resistance towilt. Plant Physiology, 2022, 189(1): 264-284.

[21] SUN Q W, CSORBA T, SKOURTI-STATHAKI K, PROUDFOOT N J, DEAN C. R-loop stabilization represses antisense transcription at thelocus. Science, 2013, 340(6132): 619-621.

[22] JABNOUNE M, SECCO D, LECAMPION C, ROBAGLIA C, SHU Q Y, POIRIER Y. A rice-natural antisense RNA acts as a translational enhancer for its cognate mRNA and contributes to phosphate homeostasis and plant fitness. The Plant Cell, 2013, 25(10): 4166-4182.

[23] DING J H, SHEN J Q, MAO H L, XIE W B, LI X H, ZHANG Q F. RNA-directed DNA methylation is involved in regulating photoperiod- sensitive male sterility in rice. Molecular Plant, 2012, 5(6): 1210-1216.

[24] ZHANG Y C, LIAO J Y, LI Z Y, YU Y, ZHANG J P, LI Q F, QU L H, SHU W S, CHEN Y Q. Genome-wide screening and functional analysis identify a large number of long noncoding RNAs involved in the sexual reproduction of rice. Genome Biology, 2014, 15(12): 512.

[25] ZHANG G Y, CHEN D G, ZHANG T, DUAN A G, ZHANG J G, HE C Y. Transcriptomic and functional analyses unveil the role of long non-coding RNAs in anthocyanin biosynthesis during sea buckthorn fruit ripening. DNA Research, 2018, 25(5): 465-476.

[26] YANG T, MA H Y, ZHANG J, WU T, SONG T T, TIAN J, YAO Y C. Systematic identification of long noncoding RNAs expressed during light-induced anthocyanin accumulation in apple fruit. The Plant Journal, 2019, 100(3): 572-590.

[27] YANG Z E, GE X Y, YANG Z R, QIN W Q, SUN G F, WANG Z, LI Z, LIU J, WU J, WANG Y, LU L L, WANG P, MO H J, ZHANG X Y, LI F G. Extensive intraspecific gene order and gene structural variations in upland cotton cultivars. Nature Communications, 2019, 10(1): 1-13.

[28] DU X M, HUANG G, HE S P, YANG Z E, SUN G F, MA X F, LI N, ZHANG X Y, SUN J L, LIU M, JIA Y H, PAN Z E, GONG W F, LIU Z H, ZHU H Q, MA L, LIU F Y, YANG D G, WANG F, FAN W, GONG Q, PENG Z, WANG L R, WANG X Y, XU S J, SHANG H H, LU C R, ZHENG H K, HUANG S W, LIN T, ZHU Y X, LI F G. Resequencing of 243 diploid cotton accessions based on an updated a genome identifies the genetic basis of key agronomic traits. Nature Genetics, 2018, 50(6): 796-802.

[29] WANG M J, YUAN D J, TU L L, GAO W H, HE Y H, HU H Y, WANG P C, LIU N, LINDSEY K, ZHANG X L. Long noncoding RNAs and their proposed functions in fibre development of cotton (). The New Phytologist, 2015, 207(4): 1181-1197.

[30] 閆飛林. 棉花纖維起始相關(guān)長(zhǎng)鏈非編碼RNA的功能鑒定[D]. 武漢: 華中農(nóng)業(yè)大學(xué), 2019.

YAN F L. The functional identification of long non-coding RNA related to cotton fiber initial[D]. Wuhan: huazhong agricultural university, 2019. (in Chinese)

[31] ZOU C S, WANG Q L, LU C R, YANG W C, ZHANG Y P, CHENG H L, FENG X X, PROSPER M A, SONG G L. Transcriptome analysis reveals long noncoding RNAs involved in fiber development in cotton (). Science China Life Sciences, 2016, 59(2): 164-171.

[32] HU H Y, WANG M J, DING Y H, ZHU S T, ZHAO G N, TU L L, ZHANG X L. Transcriptomic repertoires depict the initiation of lint and fuzz fibres in cotton (L.). Plant Biotechnology Journal, 2018, 16(5): 1002-1012.

[33] SALIH H, GONG W F, HE S P, XIA W, ODONGO M R, DU X M. Long non-coding RNAs and their potential functions in Ligon- lintless-1 mutant cotton during fiber development. BMC Genomics, 2019, 20(1): 661-678.

[34] ZHANG X P, DONG J, DENG F N, WANG W, CHENG Y Y, SONG L R, HU M J, SHEN J, XU Q J, SHEN F F. The long non-coding RNA lncRNA973 is involved in cotton response to salt stress. BMC Plant Biology, 2019, 19(1): 459-475.

[35] ZHENG X M, CHEN Y J, ZHOU Y F, SHI K K, HU X, LI D Y, YE H Z, ZHOU Y, WANG K. Full-length annotation with multistrategy RNA-seq uncovers transcriptional regulation of lncRNAs in cotton. Plant Physiology, 2021, 185(1): 179-195.

[36] KIM D, PERTEA G, TRAPNELL C, PIMENTEL H, KELLEY R, SALZBERG S L. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biology, 2013, 14(4): R36.

[37] ANDERS S, PYL P T, HUBER W. HTSeq–a python framework to work with high-throughput sequencing data. Bioinformatics, 2015, 31(2): 166-169.

[38] TRAPNELL C, ROBERTS A, GOFF L, PERTEA G, KIM D, KELLEY D R, PIMENTEL H, SALZBERG S L, RINN J L, PACHTER L. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nature Protocols, 2012, 7(3): 562-578.

[39] DATTA R, PAUL S. Long non-coding RNAs: Fine-tuning the developmental responses in plants. Journal of Biosciences, 2019, 44(4): 77-88.

[40] SUN L, LUO H T, BU D C, ZHAO G G, YU K T, ZHANG C H, LIU Y N, CHEN R S, ZHAO Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Research, 2013, 41(17): e166.

[41] YOUNG M D, WAKEFIELD M J, SMYTH G K, OSHLACK A. Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biology, 2010, 11(2): R14.

[42] WU J M, MAO X Z, CAI T, LUO J C, WEI L P. KOBAS server: a web-based platform for automated annotation and pathway identification. Nucleic Acids Research, 2006, 34(suppl_2): W720-W724.

[43] 王驍. 陸地棉纖維品質(zhì)和產(chǎn)量雜種優(yōu)勢(shì)轉(zhuǎn)錄本的動(dòng)態(tài)表達(dá)與調(diào)控模式分析[D]. 北京: 中國(guó)農(nóng)業(yè)科學(xué)院, 2021.

WANG X. Study on dynamic expression and regulation pattern of transcripts in the heterosis of fiber quality and yield ofL[D]. Beijing:Chinese Academy of Agricultural Sciences, 2021. (in Chinese)

[44] LIVAK K J, SCHMITTGEN T D. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T))method. Methods, 2001, 25(4): 402-408.

[45] 王曉陽(yáng), 王麗媛, 潘兆娥, 何守樸, 王驍, 龔文芳, 杜雄明. 亞洲棉短絨突變體纖維發(fā)育及其差異基因表達(dá)分析. 作物學(xué)報(bào), 2020, 46(5): 645-660.

WANG X Y, WANG L Y, PAN Z E, HE S P, WANG X, GONG W F, DU X M. Analysis of differentially expressed genes and fiber development infuzzless mutant. Acta Agronomica Sinica, 2020, 46(5): 645-660. (in Chinese)

[46] PEI Y. The homeodomain-containing transcription factor, GhHOX3, is a key regulator of cotton fiber elongation. Science China Life Sciences, 2015, 58(3): 309-310.

[47] DING J H, LU Q, OUYANG Y D, MAO H L, ZHANG P B, YAO J L, XU C G, LI X H, XIAO J H, ZHANG Q F. A long noncoding RNA regulates photoperiod-sensitive male sterility, an essential component of hybrid rice. Proceedings of the National Academy of Sciences of the United States of America, 2012, 109(7): 2654-2659.

[48] ZHANG W, HAN Z X, GUO Q L, LIU Y, ZHENG Y X, WU F L, JIN W B. Identification of maize long non-coding RNAs responsive to drought stress. Plos One, 2014, 9(6): e98958.

[49] DENG F N, ZHANG X P, WANG W, YUAN R, SHEN F F. Identification of Gossypium hirsutum long non-coding RNAs (lncRNAs) under salt stress. BMC Plant Biology, 2018, 18(1): 23-37.

[50] ZHANG X P, DONG J, DENG F N, WANG W, CHENG Y Y, SONG L R, HU M J, SHEN J, XU Q J, SHEN F F. The long non-coding RNA lncRNA973 is involved in cotton response to salt stress. BMC Plant Biology, 2019, 19(1): 459-475.

[51] NECSULEA A, SOUMILLON M, WARNEFORS M, LIECHTI A, DAISH T, ZELLER U, BAKER J C, GRüTZNER F, KAESSMANN H. The evolution of lncRNA repertoires and expression patterns in tetrapods. Nature, 2014, 505(7485): 635-640.

[52] CUI J, LUAN Y S, JIANG N, BAO H, MENG J. Comparative transcriptome analysis between resistant and susceptible tomato allows the identification of lncRNA16397 conferring resistance toby co-expressing glutaredoxin. The Plant Journal, 2017, 89(3): 577-589.

[53] ZHANG M, ZENG J Y, LONG H, XIAO Y H, YAN X Y, PEI Y. Auxin regulates cotton fiber initiation via GhPIN-mediated auxin transport. Plant and Cell Physiology, 2017, 58(2): 385-397.

[54] ZHANG M, ZHENG X L, SONG S Q, ZENG Q W, HOU L, LI D M, ZHAO J, WEI Y, LI X B, LUO M, XIAO Y H, LUO X Y, ZHANG J F, XIANG C B, PEI Y. Spatiotemporal manipulation of auxin biosynthesis in cotton ovule epidermal cells enhances fiber yield and quality. Nature Biotechnology, 2011, 29(5): 453-458.

[55] HAN X Y, XU X Y, FANG D D, ZHANG T Z, GUO W Z. Cloning and expression analysis of novelfamily genes in. Gene, 2012, 503(1): 83-91.

[56] HU H Y, HE X, TU L L, ZHU L F, ZHU S T, GE Z H, ZHANG X L. GhJAZ2 negatively regulates cotton fiber initiation by interacting with the R2R3-MYB transcription factor GhMYB25-like. The Plant Journal, 2016, 88(6): 921-935.

[57] GILBERT M K, BLAND J M, SHOCKEY J M, CAO H P, HINCHLIFFE D J, FANG D D, NAOUMKINA M. A transcript profiling approach reveals an abscisic acid-specific glycosyltransferase (UGT73C14) induced in developing fiber ofmutant of cotton (L.). PLoS One, 2013, 8(9): e75268.

[58] FLAGEL L E, WENDEL J F, UDALL J A. Duplicate gene evolution, homoeologous recombination, and transcriptome characterization in allopolyploid cotton. BMC genomics, 2012, 13: 302-315.

[59] ZHAO T, TAO X Y, FENG S L, WANG L Y, HONG H, MA W, SHANG G D, GUO S S, HE Y X, ZHOU B L, GUAN X Y. LncRNAs in polyploid cotton interspecific hybrids are derived from transposon neofunctionalization. Genome Biology, 2018, 19(1): 195-212.

[60] ZOU X, ALI F, JIN S, LI F, WANG Z. RNA-Seq with a novel glabrous-ZM24reveals some key lncRNAs and the associated targets in fiber initiation of cotton. BMC Plant Biology,2022, 22: 61-77.

Identification and Expression Analysis of Fuzz Fiber Development Related Long Noncoding RNAs in

WANG Xiaoyang1, PENG Zhen1,2, XING Aishuang1, ZHAO Yingrui1, MA Xinli1, LIU Fang1,2, DU Xiongming1, HE Shoupu1,2

1Institute of Cotton Research of Chinese Academy of Agricultural Sciences/National Key Laboratory of Cotton Biological Breeding and Comprehensive Utilization, Anyang 455000, Henan;2School of Agricultural Sciences, Zhengzhou University/Zhengzhou Research Base, National Key Laboratory of Cotton Biological Breeding and Comprehensive Utilization, Zhengzhou 450001

【Objective】Long non-coding RNAs(lncRNAs) are a group of RNA molecules longer than 200 bp with no protein coding capacity, which are involved in various biological regulatory processes. In this study, we aim to analyze the RNA-sequencing data of twoisogenic lines, a fuzzless mutant (GA0149) and its wildtype (GA0146), to identify the lncRNA involved in early fuzz fiber development, providing a foundation for investigation the mechanism of fiber development. 【Method】We collected 0 DPA, 3 DPA and 5 DPA ovule and 8 DPA ovule and fiber from thefuzzless mutant GA0149 and its isogenic line GA0146 with normal fuzz and lint fibers, were used for RNA-seq to identify lncRNA and predict their target genes. Differentially expressed mRNA (DE-mRNA) and lncRNA(DE-lncRNAs) between the samples were identified. The KOBAS software was used to predict the KEGG enrichment pathways which DE-lncRNAs targets were involved in. To ensure the quality of high-through sequencing, 25 DE-lncRNAs were selected for RT-qPCR detection. 【Result】We identified 15 339 lncRNA-encoding transcripts that 11 595 lncRNAs were located to intergenic regions, 2 428 lncRNAs were classified as antisense lncRNAs, 350 were categorized as intronic lncRNAs and 966 belonged to sense lncRNAs. Compared to mRNAs, lncRNAs in Asian cotton showed shorter exons and lower GC content. Most of lncRNAs had cis-regulatory effects on their neighboring mRNAs. We identified 1 932 differentially expressed (DE) lncRNAs, with 8 134 predicted DE-lncRNA target genes. Further analysis showed that 788 genes (mRNA) were differentially expressed (DE-genes) during four fiber development stages. KEGG enrichment pathways analysis showed that DE-target-mRNAs were mainly enriched in plant hormone signal transduction and protein processing in endoplasmic reticulum. Co-expression network analysis revealed that lncRNA () and its associated target genes showed identical expression trends during four fuzz fiber development stages, while lncRNAs () and its associated target genes showed contrary expression tendency, exhibiting dramatic higher expression in fuzzless GA0149 compared to wildtype GA0146. The results of RT-qPCR analysis confirmed the authenticity of our RNA-seq data.【Conclusion】A total of 26 specifically expressed lncRNAs were identified which related to cotton fuzz fiber development process. We further confirmed that these lncRNAs affected the fuzz fiber development by regulating the expression of indole-3-acetic acid-amido synthetase () and Auxin-responsive protein () in the plant hormone signal transduction pathway.

; fuzzless mutant; long non-coding RNAs; regulation network; RT-qPCR

10.3864/j.issn.0578-1752.2023.23.002

2022-12-30；

2023-03-16

國(guó)家自然科學(xué)基金（32201875）

王曉陽(yáng)，E-mail：wangxiaoyang198806@126.comwangxiaoyang198806@126.com。通信作者何守樸，E-mail：heshoupu@caas.cnheshoupu@caas.cn。通信作者杜雄明，E-mail：dxm630723@163.comdxm630723@163.com。通信作者劉方，E-mail：liufcri@163.com

（責(zé)任編輯 李莉）