基于RNA測序的尼古丁影響下大鼠前額葉皮質的基因表達分析*

張莉華 李明定 王 舉,△

實驗研究

基于RNA測序的尼古丁影響下大鼠前額葉皮質的基因表達分析*

張莉華1李明定2王 舉1,2△

目的 研究尼古丁處理對大鼠前額葉皮質（PFC）中基因表達的影響。方法24只F344大鼠隨機分為實驗組和對照組，每組12只，分別注射尼古丁水溶液或生理鹽水17 d。藥物處理結束后，從每只大鼠大腦采集PFC組織。然后采用RNA測序（RNA-Seq）技術對24只大鼠大腦PFC區域的mRNA進行測序，測出的讀段利用Bowtie、Tophat及Cufflinks等生物信息學工具比對到參考基因組上并測量轉錄子表達水平。隨后采用統計學方法確定表達水平在實驗組和對照組發生顯著變化的基因，并進一步分析了這些基因的功能類別及與其相關的生化通路。結果通過分析每個轉錄子在處理組和對照組的表達水平，確定了在尼古丁影響下的表達水平發生顯著變化的基因885個。在這些基因中，與多個Gene Ontology生物學過程如G蛋白偶聯受體信號通路、蛋白質泛素化和神經遞質吸收等，以及生化通路如胰島素信號、阿爾茨海默病和長時程增強等相關的基因均顯著富集。結論尼古丁處理可能調節大腦PFC區域神經信號傳導和能量代謝等生物學過程。

尼古丁；神經系統；前額葉皮質；序列分析,RNA；基因表達；信號傳導；能量代謝；大鼠,近交F344

吸煙引起的疾病是當今最為普遍而又危害巨大的公共衛生問題之一[1]。吸煙能加速腦萎縮和退化改變，導致中樞神經系統、大腦皮質和海馬區發育異常，影響與突觸、神經遞質和神經營養因子相關的基因的表達[2]。尼古丁是煙草中的主要致癮性物質，研究尼古丁對神經系統的作用對于認識尼古丁成癮的分子機制，開發更有效的戒煙療法具有重要意義。RNA深度測序（RNA-Seq）是近年來快速發展的高通量測序技術，具有信噪比高、分辨率高、應用范圍廣等優勢[3]。本研究利用RNA-Seq技術，分析尼古丁處理組和鹽水對照組的F344大鼠前額葉皮質（prefrontal cortex,PFC）大腦區域mRNA表達水平，通過對基因的差異表達分析和富集分析，探討尼古丁處理對大鼠PFC區域的影響及其分子機制。

1 材料與方法

1.1 實驗動物與樣品采集 24只F344大鼠購于美國Harlan公司，雄性，8周齡，SPF級，采用標準條件飼養，并對其生理狀況和行為進行監控[4]。大鼠按隨機數字表法分為2組，每組12只，實驗組按照0.25 mg/kg體質量劑量注射尼古丁溶液（尼古丁預先溶解在生理鹽水中），對照組則按照1.0 mL/kg體質量劑量注射生理鹽水，采用皮下注射方式連續17 d。實驗結束后，對大鼠大腦進行組織切片，得到大約1 mm厚的PFC區域大腦切片。利用3.00 mm的Harris微型打孔器（美國GE公司）對PFC區域進行雙邊組織提取，存放于-80℃環境，用于提取組織總RNA[5]。

1.2 RNA提取和樣本制備 TRIzol法提取PFC區域的總RNA，Qubit RNA BR分析工具（美國Life Technologies公司）檢測RNA濃度，最后由生物分析儀（美國Agilent）進行質量評估。

1.3 構建RNA測序文庫及深度測序 樣品提取總RNA后，利用TruSeq RNA樣本制備工具構建RNA樣本測序文庫。使用鏈霉親和素磁珠從1 μg RNA中提取全部帶poly（A）尾的RNA，隨機打斷成片段，以mRNA片段為模板進行反轉錄并轉換為雙鏈cDNA。通過末端修復、加堿基A，加測序接頭成為DNA片段。再經過2%的瓊脂糖收集200～250 bp的片段，進行15個循環的PCR擴增并利用QIAquick PCR純化工具（德國Qiagen）進行純化。濃縮的文庫經過緩沖液洗脫為最終的10 nmol/L濃度。最后利用Illumina Hiseq2000測序平臺對每個樣本進行50個循環的雙端測序。

1.4 數據處理與分析 RNA-Seq測序數據的處理與分析流程，見圖1。

1.4.1 實驗預處理及讀段定位 通過CASAVA（Illumina流程v1.38）提取出50 bp長的雙端讀段。樣本中滿足質量分數≥Q30的讀段用于生成完整的FASTQ格式文件，然后利用TopHat軟件[6]將讀段比對到UCSC大鼠參考基因組上，TopHat使用默認參數：每個讀段40次比對，每次比對最多允許2個錯配。序列比對文件（BAM）由RSeQC軟件包進行質量控制[7]。比對后的讀段由Cufflinks組裝成轉錄子并估計其相對豐度[8]。每個轉錄子的表達量用每百萬片段中來自于某基因每千堿基長度的片段數（fragments per kilobase of exon per million fragments，FPKM）來表示[9]。

1.4.2 基因注釋和差異表達分析 使用Ensembl轉錄子ID作為分析時的原始識別碼，并通過BioMart的Ensembl數據庫來獲取與原始識別碼相對應的注釋信息，最終檢索出39 550個轉錄子的注釋信息，其中22 919個轉錄子是有功能注釋的編碼蛋白質的基因。

Figure 1 The analysis process of RNA-Seq data圖1 RNA-Seq數據處理與分析流程

Figure 2 The distribution of probability density of original data圖2 原始數據的概率密度分布

Figure 3 The distribution of probability density of log transformed data圖3 對數轉換后數據的概率密度分布

首先對由FPKM組成的基因表達數據進行校正，忽略測量值接近于零的轉錄子（約占5%），在MATLAB平臺上自主開發程序考察數據分布，原始數據的曲線擬合分布無典型特點，見圖2。因此將數據進行底為2的對數轉換，轉換后的基因表達數據曲線擬合近似呈正態分布，見圖3。然后對數據進行異常值檢測[10]，并去除不符合條件的異常值。經過異常值檢測后，如果某個轉錄子在尼古丁處理組或鹽水對照組任何一組中的有效測量數小于6個，則認為該轉錄子的測量數據不準確，不做進一步分析。根據對應某個基因轉錄子數目的不同，考慮以下情況：（1）對應于某個基因有1個轉錄子，則使用對應的測量值進行后續分析。（2）對應于某個轉錄子有多個測量值（均滿足對比組中有效值高于6個），則使用測量值的均值作為1個記錄。（3）對應于某個基因有多個轉錄子，則將它們作為獨立的基因對待。最后以經過對數轉換和異常值去除的FPKM表達數據為對象，借鑒基因芯片的差異表達分析方法，在MATLAB平臺上使用雙端配對t檢驗進行差異表達分析[5]，并對P值進行FDR校正（false discovery rate）[11]。

1.4.3 GO功能類別及通路富集分析 GO（Gene Ontology）是一個標準化的基因功能分類體系，它通過一套動態更新的標準詞匯表來全面描述生物體中基因和基因產物的屬性。它采用分子功能（molecular function）、生物學過程（biological process）和細胞成分（cellular component）3個類別來描述每個基因，其中生物學過程描述的該基因參與的已知的有開始和結束狀態的分子生物學、生物化學或生理過程。KEGG（Kyoto Encyclopedia of Genes and Genomes）是一個基因組、代謝途徑及生物化學網絡的集合，其中的通路數據庫整合了已知的典型生化通路。實驗采用Onto-Express[12]和Onto Pathway-Express[13]對差異表達基因分別進行GO功能類別和通路富集分析。Onto-Express是一個用來挖掘功能注釋數據（主要是基于Gene Ontology），來幫助研究者找到相關生物通路的分子注釋工具[14]。Onto Pathway-Express是基于KEGG數據庫的一個通路分析工具[15]，它在通路富集分析時考慮了通路的拓撲結構和每個基因在通路中的位置，輸出結果包括通路信息以及對應的顯著性水平（P值），并進行了FDR校正[11]。

2 結果

2.1 RNA-Seq測序數據 每個樣本的初始測序數據大約有2.5 GB，共包括4 500～6 500萬個讀段，每個樣本的平均讀段數為5 400萬，其中超過95%的堿基質量分數≥Q30，平均質量分數為Q37。每個樣本有（40～50）×106個讀段（77%～83%）可以唯一地比對到大鼠參考基因組上，平均mRNA插入片段大小為（159±52）bp。Cufflinks將比對后的讀段進行轉錄子組裝、豐度估計以及標準化。所有樣本的平均轉錄子個數為50 071，經Ensembl數據庫注釋后鑒別出PFC區域有14 660個功能已知的轉錄子。

2.2 尼古丁作用下的差異表達基因 在MATLAB平臺上利用t檢驗篩選大鼠PFC區域尼古丁作用下的差異表達基因，選取P＜0.005（FDR＜0.20）且基因表達量差異倍數在1.2倍及以上的基因為差異表達基因。在此條件下與尼古丁處理相關的表達水平差異基因有885個。

2.3 差異表達基因的GO功能顯著性分析 經過Onto-Express軟件的處理后，PFC大腦區域中在尼古丁影響下差異表達基因的GO功能注釋集中在G蛋白偶聯受體信號通路、蛋白質泛素化、神經遞質吸收等GO生物學過程，見表1。

Table 1 Enrichment analysis of GO categories in the genes regulated by nicotine treatment(P＜0.05)表1 尼古丁影響下差異表達基因的GO富集分析（P＜0.05）

2.4 差異表達基因的通路顯著性富集分析 見表2。通路富集分析在Onto Pathway-express平臺上進行，以P＜0.05為篩選條件獲得了差異表達基因中顯著性富集的信號通路。尼古丁作用下差異表達基因參與的主要通路有胰島素信號通路、阿爾茨海默病及mTOR信號通路等。

Table 2 Pathways enriched in genes regulated by nicotine treatment表2 尼古丁作用下差異表達基因顯著性富集的通路

3 討論

已有的研究證明，在尼古丁作用下，多個大腦區域內許多基因的表達特性會發生改變，如轉錄速度、蛋白質合成速度或蛋白質轉運到特定亞細胞位置的速度等[10]。尼古丁進入人體內和乙酰膽堿的煙堿型受體（nAChRs）結合后，通過復雜的機制，能夠調控與多巴胺（Dopamine）、谷氨酸（Glutamate）、降腎上腺素（Norepinephrine）、5-羥色胺（Serotonin）和γ-氨基丁酸（GABA）等多種神經遞質的釋放及相關信號通路的功能[16]。近些年來，很多研究利用基因芯片、關聯分析和連鎖分析來分析尼古丁處理的影響[15]，Li等[10]利用基因芯片技術篩選出在PFC等多個區域內與尼古丁相關的多個差異表達基因，包括信號和轉錄因子，鈣離子結合蛋白，磷脂酶等。Wang等[15]通過文獻搜集和分析鑒別出多個與尼古丁相關的易感基因，這些基因涉及成癮、學習和記憶、能量代謝等多個生理過程。但是，目前對尼古丁對神經系統影響的分子機制的認識還很不充分。

本研究所確定的通路與神經系統功能或相關疾病有密切聯系。比如，與胰島素信號通路相關的差異表達基因有17個，其中包括多種磷酸化蛋白激酶如pdpk1、prkaa2等以及交換因子Sos2、鈣調蛋白calm2和生長因子受體束縛蛋白Grb2。與阿爾茨海默病相關的差異表達基因有18個，其中包括ATP連接酶atp5h、atp5g1，細胞色素C氧化酶cox6c、cox5a和多種NADH脫氫酶如ndufs4、ndufa2、ndufa7等,以及與胰島素通路共有的鈣調蛋白calm2。杜怡峰等[17]提出胰島素信號傳導通路的功能異常可能會介導Aβ引起神經元損害和學習基因功能障礙，從而引起阿爾茨海默病。本研究結果表明尼古丁處理可能調控胰島素信號通路和阿爾茨海默病信號通路，進而影響神經系統功能。此外本研究中參與mTOR信號通路的差異表達基因大部分也參與了胰島素信號通路，包括上調的3磷酸肌醇依賴的蛋白激酶pdpk1和腺苷酸活化蛋白激酶催化亞基prkaa2以及下調的促分裂原活化蛋白激酶mapk3，同時參與細胞凋亡的pik3r2和核糖體蛋白rps6。李艷琴等[18]提出mTOR可以作為篩選治療尼古丁成癮藥物的靶標，與本研究相符。

綜上，尼古丁對神經系統的影響涉及信號傳導、細胞生存、能量代謝等多個生物學過程。本研究表明尼古丁處理可能調節大腦PFC區域神經信號傳導和能量代謝等過程，并可能影響細胞凋亡等過程。

[1]Li MD,Burmeister M.New insights into the genetics of addiction[J]. Nat Rev Genet,2009,10(4):225-231.

[2]Anbarasi K,Vani G,Balakrishna K,et al.Effect of bacoside A on brain antioxidant status in cigarette smoke exposed rats[J].Life Sci, 2006,78(12):1378-1384.

[3]王曦,汪小我,王立坤,等.新一代高通量RNA測序數據的處理與分析[J].生物化學與生物物理學進展,2010,37(8):834-846.

[4]Vigorito M,LaShomb AL,Chang SL.Spatial learning and memory in HIV-1 transgenic rats[J].J Neuroimmune Pharmacol,2007,2(4): 319-328.

[5]Li MD,Cao J,Wang S,et al.Transcriptome sequencing of gene expression in the brain of the HIV-1 transgenic rat[J].PLoS One, 2013,8(3):e59582.

[6]Trapnell C,Williams BA,Pertea G,et al.Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation[J].Nat Biotechnol,2010,28 (5):511-515.

[7]Wang L,Wang S,Li W.RSeQC:quality control of RNA-seq experiments[J].Bioinformatics,2012,28(16):2184-2185.

[8]Trapnell C,Roberts A,Goff L,et al.Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks[J].Nat Protoc,2012,7(3):562-578.

[9]Mortazavi A,Williams BA,McCue K,et al.Mapping and quantifying mammalian transcriptomes by RNA-Seq[J].Nat Methods,2008, 5(7):621-628.

[10]Li MD,Kane JK,Wang J,et al.Time-dependent changes in transcriptional profiles within five rat brain regions in response to nicotine treatment[J].Brain Res Mol Brain Res,2004,132(2):168-180.

[11]Benjamini Y,Hochberg Y.Controlling the false discovery rate:a practical and p owerful approach to multiple testing[J].J R Stat Soci,1995,Series B(57):289-300.

[12]Draghici S,Khatri P,Martins RP,et al.Global functional profiling of gene expression[J].Genomics,2003,81(2):98-104.

[13]Draghici S,Khatri P,Tarca AL,et al.A systems biology approach for pathway level analysis[J].Genome Res,2007,17(10):1537-1545.

[14]Khatri P,Draghici S,Ostermeier GC,et al.Profiling gene expression using onto-express[J].Genomics,2002,79(2):266-270.

[15]Wang J,Li MD.Common and unique biological pathways associated with smoking initiation/progression,nicotine dependence,and smoking cessation[J].Neuropsychopharmacology,2010,35(3):702-719.

[16]Xiao C,Nashmi R,McKinney S,et al.Chronic nicotine selectively enhances α4β2＊nicotinic acetylcholine receptors in the nigrostriatal dopamine pathway[J].J Neurosci,2009,29(40):12428-12439.

[17]杜怡峰,閆鵬,韓曉.Alzheimer病胰島素信號傳導通路的研究[J].濟寧醫學院學報,2012,35(2):84-92.

[18]李艷琴,彭淑嫻,溫泉,等.mTOR作為標靶在篩選治療尼古丁成癮藥物中的應用:中國,201210216084[P].2012-10-03.

（2013-09-02收稿 2013-09-24修回）

（本文編輯 李國琪）

Analysis of Gene Expression through RNA Sequencing in Prefrontal Cortex of Nicotine-Treated Rats

ZHANG Lihua1,LI Mingding2,WANG Ju1,2
1 School of Biomedical Engineering,Tianjin Medical University,Tianjin 300070,China;2 School of Medicine, University of Virginia,VA 22911,USA

ObjectiveTo investigate the effects of nicotine on the gene expression profile in prefrontal cortex(PFC) region of rats.MethodsTwenty-four adult male F344 rats were randomly divided into treatment group and control group, with 12 animals in each group.Animals in treatment group were injected with nicotine solution while those in control group were given physiological saline for 17 days.At the end of drug administration,the tissue of PFC region was sliced from the brain of each animal.RNA sequencing(RNA-Seq)was utilized to analyze the gene expression profile in PFC of each animal. Bioinformatics tools including Bowtie,Tophat and Cufflinks were used to map the sequencing reads to the reference genome of rat,and to measure the relative expression level of each transcript.Genes whose expression levels were significantly different between the two groups were identified.The functional categories and the biochemical pathways associated with these genes were further analyzed.ResultsBy comparing the expression level of each transcript between two groups,885 significantly differentially expressed genes were identified.These genes were enriched in multiple Gene Ontology Biological Process categories(e.g.,G protein coupled receptor protein signaling pathway,protein ubiquitination,and neurotransmitter uptake)and biological pathways(e.g.,insulin signaling pathway,Alzheimer's disease,and long term potentiation).ConclusionNicotine treatment may regulate signaling transduction,energy metabolism and other biological processes in PFC of rat brain.

nicotine;nervous system;prefrontal cortex;sequence analysis,RNA;gene expression;signal transduction;energy metabolism;rats,inbred F344

R742，R332

A【DOI】10.3969/j.issn.0253-9896.2014.01.015

*國家自然科學基金資助項目（項目編號：31271411）

1天津醫科大學生物醫學工程學院（郵編300070）；2美國弗吉尼亞大學醫學院

△通訊作者 E-mail:wangju@tmu.edu.cn