張文華，王燕平

嗎啡急性處理小鼠的微陣列數(shù)據(jù)分析和生物標(biāo)記識(shí)別

張文華，王燕平△

目的通過生物信息學(xué)方法，識(shí)別嗎啡急性處理小鼠的差異表達(dá)基因及其富集的通路。方法從高通量基因表達(dá)數(shù)據(jù)庫(kù)（Gene Expression Omnibus）下載嗎啡急性處理小鼠的微陣列數(shù)據(jù)，并調(diào)用R語言3.1.0軟件的質(zhì)量控制包AffyQCReport 1.42.0對(duì)數(shù)據(jù)樣本進(jìn)行質(zhì)量控制分析，運(yùn)用R語言的微陣列數(shù)據(jù)線性模型識(shí)別嗎啡處理前后的差異表達(dá)基因，采用表達(dá)分析系統(tǒng)檢測(cè)算法進(jìn)行差異表達(dá)基因的通路富集分析。結(jié)果嗎啡急性處理小鼠的微陣列基因表達(dá)數(shù)據(jù)具有良好的均一性和陣列強(qiáng)度相似性，識(shí)別得到Gm11627、Zfand4、Zbtb16、Pkp2和Plin4等481個(gè)差異表達(dá)基因和癌癥、黑素原生成和絲裂原活化蛋白激酶信號(hào)等8條代謝通路。結(jié)論所識(shí)別的差異表達(dá)基因和代謝通路成為嗎啡急性處理小鼠潛在的生物標(biāo)記。

嗎啡；小鼠；微陣列分析；基因表達(dá)；生物學(xué)標(biāo)記；功能富集分析

嗎啡具有鎮(zhèn)痛、止咳、降低血壓等作用［1］。同時(shí)，嗎啡的過量吸入會(huì)形成神經(jīng)系統(tǒng)分子和細(xì)胞的適應(yīng)性，進(jìn)而導(dǎo)致對(duì)藥品的耐受、物理依賴和成癮［2］。因此，深入研究嗎啡對(duì)人體的生理作用機(jī)制，有助于對(duì)其高效地利用。生物信息學(xué)從核酸和蛋白質(zhì)序列出發(fā)，分析序列中包含的結(jié)構(gòu)功能的生物信息，可部分代替大量的分子生物學(xué)工作，具有快捷、高效、覆蓋面廣等優(yōu)勢(shì)［3］。本研究通過微陣列芯片數(shù)據(jù)的基因表達(dá)情況分析嗎啡對(duì)小鼠組織進(jìn)行急性處理前后基因表達(dá)的差異，并對(duì)差異表達(dá)基因進(jìn)行代謝通路的富集分析，為研究嗎啡對(duì)人體的作用機(jī)制提供參考。

1　材料與方法

1.1基因表達(dá)譜數(shù)據(jù)本文從高通量基因表達(dá)數(shù)據(jù)庫(kù)（Gene Expression Omnibus）［4］下載E-GEOD-7762［2］數(shù)據(jù)樣本表達(dá)微陣列數(shù)據(jù)，采用Mouse430_2平臺(tái)進(jìn)行檢測(cè)，分析流程見圖1。同時(shí)從數(shù)據(jù)庫(kù)下載原始raw文件以及該平臺(tái)探針注釋信息文件。實(shí)驗(yàn)選用來自美國(guó)緬因州巴爾港的Jackson實(shí)驗(yàn)室的4種不同種屬的小鼠（Mice），其種名分別為129P3/J、DBA/2J、C57BL/6J和SWR/J。該微陣列數(shù)據(jù)實(shí)驗(yàn)組和對(duì)照組分別選用12只小鼠為實(shí)驗(yàn)對(duì)象，每個(gè)種屬各3只，鼠齡8～10周，體質(zhì)量20～30 g，雌雄不限。嗎啡注射以前，2組小鼠均等量生理鹽水注射4 d，3次/d，食物與水正常喂養(yǎng)。實(shí)驗(yàn)組以20 mg/kg一次性皮下注射嗎啡，4 h以后取小鼠組織樣品進(jìn)行分析；對(duì)照組為等量生理鹽水注射，4 h后取小鼠組織樣品［2］。

Fig.1　The experimental protocol of microarray data on mice who were acutely administered with morphine圖1　嗎啡急性處理小鼠微陣列數(shù)據(jù)分析流程

1.2微陣列數(shù)據(jù)均一性分析調(diào)用R語言（R programming language）3.1.0軟件的質(zhì)量控制包AffyQCReport 1.42.0進(jìn)行微陣列數(shù)據(jù)的均一性分析。使用該質(zhì)量控制平臺(tái)讀取陣列數(shù)據(jù)文件，對(duì)微陣列數(shù)據(jù)進(jìn)行質(zhì)量評(píng)估，得到樣本數(shù)據(jù)的廂式圖、趨勢(shì)圖和質(zhì)量控制統(tǒng)計(jì)圖，以此衡量微陣列數(shù)據(jù)的均一性［5］。actin3/actin5基因的尺度因子＜3，且actin3/actin5內(nèi)參基因的尺度因子在-1和1之間，則為芯片數(shù)據(jù)均一性良好，以線條及端點(diǎn)表示［6-7］。

1.3微陣列數(shù)據(jù)陣列強(qiáng)度相似性分析調(diào)用R語言3.1.0軟件的質(zhì)量控制包AffyQCReport 1.42.0進(jìn)行微陣列數(shù)據(jù)的陣列強(qiáng)度相似性分析，相關(guān)系數(shù)大于0.90為相似性良好［8］。將微陣列數(shù)據(jù)輸入質(zhì)量控制平臺(tái)進(jìn)行質(zhì)量評(píng)估，由R語言輸出該樣本數(shù)據(jù)的陣列強(qiáng)度相似性報(bào)告。

1.4差異表達(dá)基因識(shí)別首先讀取全部芯片基因數(shù)據(jù)并進(jìn)行預(yù)處理。輸入原始數(shù)據(jù)，并對(duì)所有樣本表達(dá)數(shù)據(jù)進(jìn)行背景矯正。刪除沒有對(duì)應(yīng)基因的探針，同時(shí)調(diào)用R語言3.1.0的geneFilter包，運(yùn)行其中的featureFilter函數(shù)以過濾對(duì)應(yīng)有多個(gè)探針的基因，得到部分基因［9］。對(duì)上述預(yù)處理所得基因進(jìn)行進(jìn)一步篩選，運(yùn)用R語言的微陣列數(shù)據(jù)線性模型（linear models for microarray data，LIMMA）提取基因表達(dá)量差異倍數(shù)值大于2的基因［10］。并校正假陽性率（False discovery rate）＜5%和P＜0.05，得到相關(guān)的差異表達(dá)基因［11］。

1.5通路富集分析本文采用表達(dá)分析系統(tǒng)檢測(cè)算法（ex?pression analysis systematic explorer，EASE）進(jìn)行京都基因與基因組百科全書（kyoto encyclopedia of genes and genomes，KEGG）通路富集在線分析［12-13］?；贓ASE算法的KEGG通路富集分析通過生物功能信息學(xué)微陣列分析注釋（DAVID）網(wǎng)站［14］實(shí)現(xiàn)。在KEGG通路分析過程中，提交所得的481個(gè)差異表達(dá)基因，并設(shè)置P＜0.05和通路最小基因數(shù)=2作為差異表達(dá)基因的富集條件，以篩選差異表達(dá)基因富集程度較高的通路。

2　結(jié)果

2.1微陣列數(shù)據(jù)均一性分析2組全部數(shù)據(jù)的actin3/actin5的尺度因子均＜2，且內(nèi)參基因actin3/actin5的尺度因子在-1～1之間，芯片數(shù)據(jù)均一性良好，見圖2。

Fig.2　Intension similarity of microarray data of E-GEOD-7762圖2　E-GEOD-7762微陣列芯片數(shù)據(jù)陣列強(qiáng)度相似性

2.2微陣列芯片數(shù)據(jù)陣列強(qiáng)度相似性分析2組的數(shù)據(jù)相似性均＞0.90，整體陣列強(qiáng)度相似性良好，見圖3。

2.3差異表達(dá)基因識(shí)別本文共讀取45 101條芯片基因數(shù)據(jù)，其中10 313條探針基因，最終識(shí)別得到481個(gè)差異表達(dá)基因。差異表達(dá)基因的火山圖，見圖4。P值較小的前20個(gè)差異表達(dá)基因，見表1（僅列出NCBI官方基因簡(jiǎn)稱，以P值排序）。微陣列芯片數(shù)據(jù)與前20個(gè)差異表達(dá)基因的表達(dá)調(diào)節(jié)熱圖，見圖5。其中P值較小的Zbtb16、Plin4、Pkp2、Cebpa、Gm11627和Zfand4等差異表達(dá)基因在實(shí)驗(yàn)組的表達(dá)量顯著增大。

2.4通路富集分析嗎啡急性處理小鼠后，差異表達(dá)基因主要富集于癌癥、黑素原生成和絲裂原活化蛋白激酶信號(hào)等8條通路，以P值排序，見表2。

Fig.4　Volcano plot for P-value versus ratio for each gene samples including subjects and normal samples圖4　2組所有基因數(shù)據(jù)的P值與表達(dá)差異值對(duì)應(yīng)關(guān)系火山圖

Tab.1　Top 20 differentially expressed genes identified from the microarray datasets of E-GEOD-7762表1　前20個(gè)E-GEOD-7762微陣列數(shù)據(jù)識(shí)別的差異表達(dá)基因

3　討論

Korostynski等［2］對(duì)小鼠紋狀體轉(zhuǎn)錄組進(jìn)行分析，研究嗎啡急性處理后4只不同種屬的小鼠之間轉(zhuǎn)錄應(yīng)答的區(qū)別，以及嗎啡慢性處理后的基因表達(dá)情況。本文對(duì)嗎啡急性處理小鼠的微陣列數(shù)據(jù)進(jìn)行質(zhì)量控制分析，識(shí)別其差異表達(dá)基因，并進(jìn)行差異表達(dá)基因的通路富集分析。質(zhì)量控制報(bào)告顯示芯片數(shù)據(jù)均一性和陣列強(qiáng)度相似性良好。本文共識(shí)別到Gm11627、Zfand4、Zbtb16、Pkp2和Plin4等481個(gè)差異表達(dá)基因。表達(dá)差異最為顯著的前20個(gè)基因在嗎啡處理后的表達(dá)量變化顯著，差異表達(dá)基因富集于癌癥通路、黑素原生成和絲裂原活化蛋白激酶信號(hào)通路等8條代謝通路之中，且差異表達(dá)基因的富集程度較高，通路最小基因數(shù)均≥5，具有統(tǒng)計(jì)學(xué)意義。

Tab.2　Pathway enrichment analysis based on KEGG表2　基于KEGG的通路富集分析

本文所識(shí)別的Gm11627、Zfand4、Zbtb16、Pkp2 和Plin4等差異表達(dá)基因在嗎啡對(duì)小鼠急性處理前后顯著差異表達(dá)，并與該過程導(dǎo)致的生理變化有一定關(guān)系。已有研究證明，uGTZB7基因序列的多態(tài)性影響嗎啡的體內(nèi)代謝［15］。大鼠JNK3和GRK基因的表達(dá)分別與腦內(nèi)神經(jīng)細(xì)胞損傷和嗎啡成癮，以及嗎啡急性處理過程中的耐受有關(guān)［16］。差異表達(dá)基因Zbtb16作為一種轉(zhuǎn)錄抑制劑，與嗎啡成癮現(xiàn)象和小鼠腦紋狀體的細(xì)胞分化功能密切相關(guān)［17］。此外，Pkp2基因是導(dǎo)致心律失常性右心室心肌病的重要致病基因［18］。本研究KEGG通路富集分析結(jié)果顯示，包括Pkp2在內(nèi)的6個(gè)差異表達(dá)基因富集在致心律失常性右心室心肌病代謝通路中，進(jìn)一步印證了Pkp2在心肌病表達(dá)方面的功能，提示嗎啡急性處理后Zbtb16、Pkp2等基因的差異表達(dá)可能影響細(xì)胞分化，并可能導(dǎo)致嗎啡成癮和心肌病。

已有研究顯示，絲裂原活化蛋白激酶信號(hào)廣泛存在于神經(jīng)元中［19］，通過影響神經(jīng)突觸的可塑性產(chǎn)生對(duì)嗎啡的病理性鎮(zhèn)痛耐受作用［20］。絲裂原活化蛋白激酶信號(hào)通路在前額葉皮質(zhì)磷酸化過程中導(dǎo)致細(xì)胞外調(diào)節(jié)蛋白激酶表達(dá)降低，進(jìn)而參與嗎啡依賴的形成和戒斷反應(yīng)，造成嗎啡成癮［21］。此外，絲裂原活化蛋白激酶信號(hào)通路具有信號(hào)轉(zhuǎn)導(dǎo)作用，因此其在細(xì)胞的增殖、分化及凋亡等生化反應(yīng)過程中起到重要作用。而鞘脂類在脂筏的形成過程中扮演著重要角色，其代謝過程與細(xì)胞的信號(hào)轉(zhuǎn)導(dǎo)有關(guān)［22］。同時(shí)，焦點(diǎn)黏附通常在焦點(diǎn)黏附激酶的作用下發(fā)生，該過程可以間接反映細(xì)胞與細(xì)胞外基質(zhì)之間的信號(hào)傳遞，與細(xì)胞遷移、侵襲、存活等生物過程密切相關(guān)［23］。代謝通路的富集反映出，嗎啡急性處理后可能發(fā)生絲裂原活化蛋白激酶信號(hào)轉(zhuǎn)導(dǎo)、鞘脂的代謝和焦點(diǎn)黏附等生理過程。

綜上所述，嗎啡急性處理后可能發(fā)生細(xì)胞分化、絲裂原活化蛋白激酶信號(hào)轉(zhuǎn)導(dǎo)、鞘脂的代謝和焦點(diǎn)黏附等生理過程，并可能導(dǎo)致成癮和心肌病。上述差異表達(dá)基因及其富集的代謝通路可成為嗎啡對(duì)小鼠急性處理潛在的生物標(biāo)記，為進(jìn)一步研究嗎啡對(duì)人體的作用機(jī)制提供了依據(jù)。目前生物信息學(xué)分析主要集中在多元數(shù)據(jù)網(wǎng)絡(luò)的分析及可視化，蛋白質(zhì)相互作用的評(píng)估和預(yù)測(cè)，網(wǎng)絡(luò)聚類分析和相關(guān)子網(wǎng)絡(luò)的構(gòu)建，以及生物功能和代謝通路分析等方面［24］。因此，微陣列數(shù)據(jù)的擴(kuò)展、差異表達(dá)基因的蛋白質(zhì)相互作用網(wǎng)絡(luò)分析、基因本體分析等，都可能成為今后關(guān)于嗎啡對(duì)人體作用機(jī)制的研究方向。

（圖3、5見插頁）

［1］趙倩華.嗎啡給藥途徑及其新劑型的研究進(jìn)展［J］.醫(yī)藥論壇雜志，2013，34（1）:143-145.

［2］Korostynski M，Piechota M，Kaminska D，et al.Morphine effects on striatal transcriptome in mice［J］.Genome Biol，2007，8（6）:R128.

［3］徐娟.miRNA-miRNA協(xié)同調(diào)控網(wǎng)絡(luò)：構(gòu)建、疾病miRNA拓?fù)涮卣骷靶蛄泻徒Y(jié)構(gòu)相似性分析［D］.哈爾濱:哈爾濱醫(yī)科大學(xué)，2011.

［4］Li C，Shen W，Shen S，et al.Gene expression patterns combined with bioinformatics analysis identify genes associated with cholan?giocarcinoma［J］.Comput Biol Chem，2013，47:192-197.

［5］Gautier L，Cope L，Bolstad BM，et al.Affy--analysis of Affymetrix GeneChip data at the probe level［J］.Bioinformatics，2004，20（3）: 307-315.

［6］Wilson CL，Miller CJ.Simpleaffy:a BioConductor package for Af?fymetrix Quality Control and data analysis［J］.Bioinformatics，2005，21（18）:3683-3685.

［7］Raman T，O'Connor TP，Hackett NR，et al.Quality control in micro?array assessment of gene expression in human airway epithelium［J］. BMC Genomics，2009，10（1）:493.

［8］DUNDAS JB.Affymetrix Quality control［EB/OL］.（2011-3-7）［2014-07-26］http://bioinformatics.gr/tutorials/microarray-technol?ogy-basics/maqc/95-affyqc.html.

［9］Kim Y，Doan BQ，Duggal P，et al.Normalization of microarray ex?pression data using within-pedigree pool and its effect on linkage analysis［J］.BMC proceedings BioMed Central Ltd，2007，1（Suppl 1）：S152.

［10］Smyth GK.Linear models and empirical bayes methods for assess?ing differential expression in microarray experiments［J］.Stat Appl Genet Mol Biol，2004，3（1）:3.

［11］Reiner A，Yekutieli D，Benjamini Y.Identifying differentially ex?pressed genes using false discovery rate controlling procedures［J］. Bioinformatics，2003，19（3）:368-375.

［12］Huang DW，Sherman BT，Lempicki RA.Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources ［J］.Nature protocols，2009，4（1）:44-57.

［13］Ford G，Xu Z，Gates A，et al.Expression Analysis Systematic Ex?plorer（EASE）analysis reveals differential gene expression in per?manent and transient focal stroke rat models［J］.Brain Res，2006，1071（1）:226-236.

［14］Huang DW，Sherman BT，Lempicki RA.Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources ［J］.Nat Protoc，2009，4:44-57.

［15］龔妍，曹迪，李智平.UGT2B7代謝酶，阿片μ受體及ABCB1基因多態(tài)性對(duì)小兒?jiǎn)岱孺?zhèn)痛效能影響的研究進(jìn)展.第二十三屆全國(guó)兒科藥學(xué)學(xué)術(shù)會(huì)議論文集［C］.石家莊:中國(guó)藥房雜志社，2012.

［16］范學(xué)良.嗎啡對(duì)GRKZ、GRKS和JNK3基因在大鼠腦內(nèi)表達(dá)調(diào)控的研究［D］.上海:復(fù)旦大學(xué)，2002.

［17］Piechota M，Korostynski M，Solecki W，et al.The dissection of tran?scriptional modules regulated by various drugs of abuse in the mouse striatum［J］.Genome Biol，2010，11（5）:R48.

［18］張明昌，張寧生，童奇，等.PKP2基因突變與青壯年猝死綜合征.全國(guó)第九次法醫(yī)學(xué)術(shù)交流會(huì)論文集［C］.北京:中國(guó)法醫(yī)學(xué)會(huì)，2013.

［19］Shen N，Guo RX，Hu F，et al.Spinal Cx43 mediates chronic mor?phine antinociceptive tolerance through JNK pathway in rats［J］.CHI?NESE PHARMACOLOGICAL BULLETIN，2013，29（03）372-376.

［20］張樂，朱永平.絲裂原活化蛋白激酶信號(hào)通路在嗎啡依賴中的作用［J］.中國(guó)藥物依賴性雜志，2012，21（2）:100-104.

［21］閻春霞，李陽，張宏星，等.腦前額葉皮層MAPK通路在嗎啡成癮記憶中的作用.全國(guó)第九次法醫(yī)學(xué)術(shù)交流會(huì)論文集［C］.北京:中國(guó)法醫(yī)學(xué)會(huì)，2013.

［22］劉圣紅，茍萍.鞘脂類的研究進(jìn)展［J］.生物技術(shù)，2009，19（2）:96-98.

［23］袁周:焦點(diǎn)粘附激酶在肝癌中表達(dá)及其對(duì)人肝癌細(xì)胞MHCC97-H作用的研究［D］.上海:中山醫(yī)院，2006.

［24］陳剛.生物網(wǎng)絡(luò)分析及其在復(fù)雜疾病研究中的應(yīng)用［D］.長(zhǎng)沙:中南大學(xué)，2012.

（2014-07-28收稿2014-08-15修回）

（本文編輯陸榮展）

Microarray datasets and biomarker in mice who were acutely administered with morphine

ZHANG Wenhua1，WANG Yanping2
1 Department of Anesthesiology，Third Hospital Affiliated to Qiqihaer Medical College，Qiqihaer，Heilongjiang 161000，China；2 Evidence-Based Medicine Technology Development Center in Jinan
△Corresponding AuthorE-mail:1102207555@qq.com

ObjectiveTo identify the differentially expressed genes and functional enrichment pathways using bioin?formatics technology in mice who were acutely administered with Morphine.MethodsFirst，we downloaded microarray da?tasets of mice which were acutely administered with Morphine from Gene Expression Omnibus，and assess the quality control parameters of the microarray datasets using AffyQCReport 1.42.0 package by R programming language 3.1.0 software.Subse?quently，the differentially expressed genes of the microarray datasets were identified using the linear models for microarray data of R language.Finally，functional pathways of the differentially expressed genes were enriched based on expression anal?ysis and systematic explorer.ResultsThe microarray datasets showed preferable uniformity and intension similarity.A to?tal of 481 differentially expressed genes including Gm11627，Zfand4，Zbtb16，Pkp2 and Plin4 were identified.While 8 func?tional enrichment pathways，including pathways in cancer，melanogenesis and mitogen-activated protein kinase signaling were revealed.ConclusionThe differentially expressed genes and functional pathways were the underlying biomarkers of mice who were acutely administered with Morphine.

morphine；mice；microarray analysis；gene expression；biological markers；enrichment of functional analysis

R318.04

ADOI：10.11958/j.issn.0253-9896.2015.02.013

2012年國(guó)家醫(yī)學(xué)教育發(fā)展中心課題項(xiàng)目（2012-03-07-142）

1黑龍江齊齊哈爾，齊齊哈爾醫(yī)學(xué)院附屬第三醫(yī)院麻醉科（郵編161000）；2濟(jì)南，濟(jì)南循證醫(yī)藥科技開發(fā)中心科研

張文華（1972），女，碩士，副主任醫(yī)師，主要從事麻醉技術(shù)及麻醉用藥研究

△通迅作者E-mail:1102207555@qq.com