[摘要]"目的"利用基因共表達權(quán)重網(wǎng)絡（weighted"gene"co-expression"network"analysis，WGCNA）探索與尿道下裂發(fā)病相關的潛在基因。方法"利用WGCNA算法將數(shù)據(jù)集GSE35034處理并構(gòu)建基因共表達權(quán)重網(wǎng)絡，篩選出與尿道下裂相關性最高的模塊，通過基因本體論（gene"ontology，GO）、京都基因和基因組百科全書（Kyoto"Encyclopedia"of"Genes"and"Genomes，KEGG）對模塊中的基因進行富集分析。同時進行差異分析，篩選差異基因。將差異基因?qū)隨tring數(shù)據(jù)庫，利用Cytoscape軟件，找出網(wǎng)絡中樞紐基因。將以上3"種方法的結(jié)果系統(tǒng)分析，篩選出交集中的核心基因。利用外部數(shù)據(jù)集對核心基因進行mRNA表達變化及受試者操作特征（receiver"operating"characteristic，ROC）曲線診斷驗證。結(jié)果"基于WGCNA方法得到15個共表達模塊。通過差異分析獲得了93個滿足條件的共同差異基因。通過Cytoscape軟件分析最終得到10個核心基因。最后3種方法得到2個交集基因：FBXL16、SYNDIG1。ROC曲線驗證了交集基因在尿道下裂患者與正常人中表達存在差異。結(jié)論"本研究通過WGCNA獲得了2個與尿道下裂具有顯著相關性的關鍵基因，可用于尿道下裂發(fā)病及治療的探索。

[關鍵詞]"尿道下裂；基因共表達權(quán)重網(wǎng)絡；差異分析；富集分析

[中圖分類號]"R726.9""""""[文獻標識碼]"A""""""[DOI]"10.3969/j.issn.1673-9701.2024.14.002

Analysis"of"key"gene"related"to"hypospadias"based"on"gene"co-expression"weight"network"analysis

LIU"Xiaoya1，"CHANG"Mengmeng1，"MA"Wenyue1，"GAO"Hongjie2，"SUN"Fengyin1，3

1.Department"of"Pediatric"Surgery，"Qilu"Hospital"of"Shandong"University，"Jinan"250012，"Shandong，"China;"2.Department"of"Pediatrics，"Qilu"Hospital"of"Shandong"University，"Jinan"250012，"Shandong，"China;"3.Center"for"Post-Doctoral"Studies"of"Shandong"Qidu"Pharmaceutical"Co.，LTD，"Jinan"250012，"Shandong，"China

[Abstract]"Objective"To"explore"potential"genes"associated"with"the"pathogenesis"of"hypospadias"using"weighted"Gene"co-expression"network"analysis"（WGCNA）."Methods"The"WGCNA"algorithm"was"used"to"process"the"hypospadias-related"dataset"GSE35034，"and"then"a"gene"co-expression"weight"network"was"constructed"to"screen"the"modules"with"the"highest"correlation"with"hypospadias，"and"the"genes"in"the"modules"were"enriched"and"detected"by"gene"ontology"（GO），"Kyoto"Encyclopedia"of"Genes"and"Genomes"（KEGG）."Differential"analysis"was"also"performed"to"screen"out"differential"genes."The"differential"genes"were"imported"into"the"String"database."Using"Cytoscape"software，"the"hub"genes"in"the"network"were"identified."The"results"screened"by"the"above"threenbsp;methods"were"combined"and"analyzed，"and"the"core"genes"in"the"intersection"set"were"screened."Using"the"external"dataset"GSE121712，"the"core"genes"were"verified"by"mRNA"expression"changes"and"subject"work"characterization"receiver"operating"characteristic"（ROC）"curve"diagnosis."Results"Fifteen"co-expression"modules"were"obtained"based"on"the"WGCNA"method."Ninety-three"common"differential"genes""meeting"the"conditions"were"obtained"by"differential"analysis."Ten"core"genes"were"finally"obtained"by"Cytoscape"software"analysis."Finally"MEbrown"module，"differential"genes"and"the"10"core"genes"yielded"a"total"of"2"intersecting"genes："FBXL16"and"SYNDIG1."ROC"curves"verified"that"the"intersecting"genes"were"differentially"expressed"in"patients"with"hypospadias"versus"normal"subjects"."Conclusion"In"this"study，"two"key"genes"with"significant"correlation"with"hypospadias"were"obtained"by"WGCNA，"which"may"be"used"for"the"early"diagnosis"and"treatment"of"hypospadias"after"further"study.

[Key"words]"Hypospadias;"Weighted"gene"co-expression"network"analysis;"Differential"analysis;"Enrichment"analysis

尿道下裂是最常見的男性泌尿生殖系統(tǒng)畸形之一，但其發(fā)病機制仍不明確。近年來，隨著外科手術技術不斷進步，已出現(xiàn)多種手術方式，但治療效果仍不滿意；近幾十年來發(fā)病率有上升趨勢，因此人們對尿道下裂的發(fā)病原因也更加關注[1]。尿道下裂的遺傳易感性具有公認的家族性聚集性，患有尿道下裂的男孩的男性親屬更容易患這種疾病。當遺傳易感性與抗雄激素藥物暴露相結(jié)合時，基因和環(huán)境風險因素可能相互作用超過閾值，導致這種出生缺陷的發(fā)生[2]。因此，尿道下裂發(fā)育畸形和生殖疾病的發(fā)生，存在基因和環(huán)境的共同作用，并且攜帶遺傳易感基因合并環(huán)境內(nèi)分泌干擾物的暴露可能是發(fā)病的主要原因。

在基因組時代，隨著高通量測序技術發(fā)展和基因表達數(shù)據(jù)庫的完善，為尿道下裂的病因?qū)W研究提供了新的途徑。基因表達數(shù)據(jù)庫（gene"expression"omnibus，GEO）是當今最全面的基因表達數(shù)據(jù)庫，目前儲存了高達208萬樣本的表達譜數(shù)據(jù)。本研究從GEO上選擇了包含3種人類包皮成纖維細胞基因表達譜的芯片，進行加權(quán)基因共表達網(wǎng)絡分析。

WGCNA是一種基于基因表達數(shù)據(jù)集檢測，研究高度相關基因簇的常用方法，通過構(gòu)建基因表達網(wǎng)絡，篩選生物標志物或治療靶點[3-4]。WGCNA可根據(jù)基因表達模式將不同基因進行聚類，分析基因模塊與特定性狀之間的關系，篩選出有效生物標志物[5-6]。目前，WGCNA被廣泛用于疾病、生理學、藥物、進化和基因組的生物學研究[7]。

本研究擬通過對尿道下裂相關的基因表達數(shù)據(jù)分析整合后尋找出與尿道下裂發(fā)生相關的致病基因，為今后研究疾病的發(fā)病機制和靶向治療提供線索。

1""材料與方法

1.1""數(shù)據(jù)獲取

從GEO（https//www.ncbi.nlm.nih.gov/gds）數(shù)據(jù)庫中下載與尿道下裂有關的芯片包含表達數(shù)據(jù)及臨床信息（GSE35034）。其收錄了來自日本國家兒童健康與發(fā)展研究所的尿道下裂兒童患者（n=23，中位年齡2.3歲）接受外科手術治療的3種人類包皮成纖維細胞（human"foreskin"fibroblasts，hFFC）進行了全基因組篩選。用二甲基亞砜（dimethyl"sulfoxide，DMSO）作為對照，分別暴露于10nmol/L雙酚A（bisphenol"A，BPA）、0.01nmol/L17β-雌二醇（estradiol，E2）和1nmol/L"2，3，7，8-四氯二苯并對二英（tetrachlorodibenzo-p-dioxin，TCDD）24h后測量來自尿道下裂患兒的3種hFFC的基因表達譜。所有樣本均采用Agilent-028004"SurePrint"G3"Human"GE平臺檢測。

1.2""WGCNA共表達模塊分析與功能富集分析

WGCNA是一種系統(tǒng)的生物學方法，常用于描述不同樣本之間的基因關聯(lián)模式[8]。使用“WGCNA”"R包構(gòu)建共表達網(wǎng)絡[9]。選取合適軟閾值，并篩選出權(quán)重最大模塊。在WGCNA最佳模塊中選擇核心基因（Rgt;0.8，Plt;0.01），通過利用R的“clusterprofiler”包對模塊核心基因進行基因本體論（gene"ontology，GO）和基因和基因組百科全書（Kyoto"Encyclopedia"of"Genes"and"Genomes，KEGG）通路功能富集，篩選出具有顯著富集的功能和信號通路。

1.3""篩選差異基因與核心基因

利用R語言的“Limma”包對GSE35034所有轉(zhuǎn)錄數(shù)據(jù)進行差異分析，獲取與尿道下裂發(fā)病相關的差異基因[10]。對E2、TCDD、BPA組數(shù)據(jù)篩選得到的差異表達基因取交集，獲得表達一致的差異基因。將差異基因輸入到String數(shù)據(jù)庫中，篩選并繪制蛋白質(zhì)互作網(wǎng)絡（protein-protein"interaction"network，"PPI）圖[11]。使用Cytoscape軟件對網(wǎng)絡進行分析，選取連通性最高的前10個基因作為尿道下裂發(fā)病的樞紐基因。為使結(jié)果更為準確，將WGCNA核心模塊基因、共同差異基因及PPI網(wǎng)絡所篩選的基因，系統(tǒng)分析取其交集，篩選出重疊的關鍵基因。再對關鍵基因進行mRNA表達變化及受試者操作特征（receiver"operating"characteristic，ROC）曲線，計算ROC曲線下面積（area"under"curve，AUC），診斷驗證。

1.4""統(tǒng)計學方法

本研究所有數(shù)據(jù)采用R軟件分析，差異分析采用貝葉斯算法，相關性分析采用皮爾遜相關性分析，兩組間的數(shù)據(jù)分析采用雙t檢驗，適當?shù)臅r候用t檢驗，以Plt;0.05為差異有統(tǒng)計學意義。

2""結(jié)果

2.1""樣本預處理及基因共表達模塊的構(gòu)建

對原始數(shù)據(jù)進行預處理，共獲得109"834個基因的表達數(shù)據(jù)。對樣本進行離群性檢驗，然后進行聚類分析。首先對軟閾值進行篩選，當軟閾值β設置為8時，樣本獨立性較高，且可獲得較好的連接關系。將基因動態(tài)切割，得到基因聚類樹（圖1）。再連鎖分層聚類合并，最終得到15個模塊，各模塊基因個數(shù)從110～357不等。

2.2""基因模塊與臨床性狀關聯(lián)分析

計算模塊的特征基因與樣本來源的皮爾遜相關系數(shù)，并考慮到相關系數(shù)和P值，r=0.72，P=0.008作為具有臨床意義的模塊。然后對其進行基因重要性及基因表達水平分析，并繪制散點圖分析模塊內(nèi)基因是否符合線性（圖2）。

2.3""GO和KEGG富集分析

利用R的“clusterProfiler”包對選擇模塊基因進行GO和KEGG功能富集分析（圖3）。GO富集結(jié)

果顯示，模塊基因主要參與腎臟發(fā)育，生殖系統(tǒng)發(fā)育，促進脂質(zhì)生物合成等過程。在KEGG通路方面，主要富集于一些與細胞增殖分化等有關調(diào)節(jié)，如癌癥通路，PI3K-Akt信號通路，蛋白質(zhì)消化吸收，PPAR信號通路。

2.4""差異基因的PPI構(gòu)建

從GEO數(shù)據(jù)庫分別獲得的E2、TCDD、BPA組數(shù)據(jù)，經(jīng)過差異分析，得到滿足條件的差異基因，再將三組數(shù)據(jù)取交集共得到共同差異表達基因93個，其中上調(diào)基因"37個，下調(diào)基因56個，并用火山圖展示其表達的差異性（圖4）。

將篩選出的共同差異基因?qū)隨TRING，篩選條件設置為：置信度≥0.15，互作最大值等于0。之后把STRING的計算結(jié)果導入Cytoscape"軟件。應用cytoHubba插件篩選出節(jié)點打分排名前10的關鍵基因：CTNND2、SLITRK1、KCNF1、FAM107A、CTTNBP2、SCAMP5、TMEM132B、LAMP5、SYNDIG1、FBXL16。

2.5""關鍵基因表達驗證和ROC曲線分析

為進一步驗證核心分子準確性，對數(shù)據(jù)集中WGCNA核心模塊基因、共同差異基因及Cytoscape軟件所篩選的基因，整合分析，最后篩選出關鍵基因SYNDIG1和FBXL16為有效的生物標志物。

在數(shù)據(jù)集GSE121712中驗證所篩選的基因在尿道下裂患者與正常組表達存在差異（圖5），同時，ROC分析表明，發(fā)現(xiàn)SYNDIG1和FBXL16基因?qū)τ谀虻老铝训拇嬖诰哂休^好的預測能力（圖6）。

3""討論

尿道下裂是一種男性常見生殖系統(tǒng)畸形，目前手術治療效果不滿意且術后并發(fā)癥多，因此早期診斷預防具有重要意義。目前公認的尿道下裂發(fā)病可能與基因、環(huán)境等多因素相關，但其詳細機制仍不明確。攜帶易感基因并存在環(huán)境內(nèi)分泌干擾物暴露是其中一種較為合理的推論。本研究通過WGCNA方法，篩選與尿道下裂發(fā)病相關的基因，并進行了驗證分析。

通過基因共表達網(wǎng)絡，篩選出了與尿道下裂發(fā)病具有顯著相關性的棕色模塊。然后進行了功能富集分析。根據(jù)富集結(jié)果顯示，關鍵基因在生物學過程方面主要參與腎臟發(fā)育，泌尿系統(tǒng)發(fā)育和泌尿生殖系統(tǒng)發(fā)育；而在功能和信號通路方面，主要參與一些與細胞增殖分化等有關調(diào)節(jié)，如："PI3K-Akt信號通路，卵巢類固醇生成，PPAR信號通路。其中，PI3K-Akt信號通路在先天性尿道下裂中發(fā)揮重要作用，PPAR信號通路與尿道發(fā)育相關[12-13]。相關研究表明，PI3K/Akt信號廣泛參與調(diào)節(jié)各種細胞過程，包括增殖、分化、存活、遷移和凋亡等；PI3K/Akt的激活可以減少雄激素受體介導的轉(zhuǎn)錄[14]。

同時本研究進行差異分析和PPI網(wǎng)絡構(gòu)建，與棕色模塊聯(lián)合分析后，從中篩選出與尿道下裂生物發(fā)病過程最可能相關的兩個基因SYNDIG1和FBXL16。SYNDIG1即突觸分化誘導基因1，主要參與調(diào)控突觸后受體密度[15]。

FBXL16即F-box和富含亮氨酸重復蛋白16"，是一種F-box蛋白，參與幾乎所有生命活動的調(diào)節(jié)，如細胞周期、增殖和凋亡。有研究表明，F(xiàn)BXL16可能與泛素化、磷酸化有關，可以促進細胞增殖[16]。Sato等[17]的實驗也證明，F(xiàn)BXL16表達的缺失促進了細胞增殖。一項關于乳腺癌的研究中也發(fā)現(xiàn)FBXL16沉默抑制了細胞中p-AKT和p-mTOR的水平。FBXL16表達抑制乳腺癌上皮-間充質(zhì)轉(zhuǎn)化和血管生成也被證實。FBXL16還被證明是心血管祖細胞分化抑制因子之一[18-20]。在胚胎發(fā)育過程中泌尿系統(tǒng)及原始血管系統(tǒng)均起源于中胚層。原始血管的形成與血島邊緣細胞向血管內(nèi)皮細胞的分化有關。男性尿道則形成主要依賴于尿生殖竇的中下段，生殖結(jié)節(jié)演變?yōu)殛幥o，左、右尿生殖褶在中線閉合并形成尿道海綿體部。FBXL16已被證實可以抑制心血管祖細胞分化。而類似的機制也可能發(fā)生于尿道溝閉合、尿道形成時。其具體機制仍有待進一步探索和實驗驗證。

綜上所述，本研究通過構(gòu)建WGCNA構(gòu)建共表達網(wǎng)絡，分析處理后發(fā)現(xiàn)了SYNDIG1和FBXL16兩個基因與環(huán)境共同作用后在尿道下裂生物發(fā)病過程中可能起著重要作用。其機制可能與PI3K/Akt/"mTOR信號通路相關。這些發(fā)現(xiàn)對于研究尿道發(fā)病機制，以及早期診斷、靶向治療有重要參考價值。

利益沖突：所有作者均聲明不存在利益沖突。

[參考文獻]

[1] VAN"D"Z"L"F，"VAN"R"I"A，"FEITZ"W"F，"et"al."Aetiology"of"hypospadias："A"systematic"review"of"genes"and"environment[J]."Hum"Reprod"Update，"2012，"18（3）："260–283.

[2] MANSON"J"M，"CARR"M"C."Molecular"epidemiology"of"hypospadias："Review"of"genetic"and"environmental"risk"factors[J]."Birth"Defects"Res.，"Part"A，"2003，"67（10）："825–836.

[3] WEI"J，"ZHANG"J，"WEI"J，"et"al."Identification"of"AGXT2，"SHMT1，"and"ACO2"as"important"biomarkers"of"acute"kidney"injury"by"WGCNA[J]."PLoS"One，"2023，"18（2）："e0281439.

[4] ZHANG"T，"WONG"G."Gene"expression"data"analysis"using"hellinger"correlation"in"weighted"gene"co-expression"networks"（WGCNA）[J]."Comput"Struct"Biotechnol"J，"2022，"20："3851–3863.

[5] MASON"M"J，"FAN"G，"PLATH"K，"et"al."Signed"weighted"gene"co-expression"network"analysis"of"transcriptional"regulation"in"murine"embryonic"stem"cells[J]."BMC"Genomics，"2009，"10："327.

[6] LIANG"W，"SUN"F，"ZHAO"Y，"et"al."Identification"of"susceptibility"modules"and"genes"for"cardiovascular"disease"in"diabetic"patients"using"WGCNA"analysis[J]."J"Diabetes"Res，"2020，"2020："4178639.

[7] TIAN"Z，"HE"W，"TANG"J，"et"al."Identification"of"important"modules"and"biomarkers"in"breast"cancer"based"on"WGCNA[J]."Onco"Targets"Ther，"2020，"13："6805–6817.

[8] WANG"Y，"ZHOU"W，"CHEN"Y，"et"al."Identification"of"susceptibility"modules"and"hub"genes"of"osteoarthritis"by"WGCNA"analysis[J]."Front"Genet，"2022，"13："1036156.

[9] LANGFELDER"P，"HORVATH"S."WGCNA："An"R"package"for"weighted"correlation"network"analysis[J]."BMC"Bioinf，"2008，"9："559.

[10] WEN"J，"REN"T，"ZHENG"J，"et"al."Identification"and"verification"of"pivotal"genes"promoting"the"progression"of"atherosclerosis"based"on"WGCNA[J]."Artif"Cells，"Nanomed，"Biotechnol，"2023，"51（1）："276–285.

[11] SZKLARCZYK"D，"GABLE"A"L，"NASTOU"K"C，"et"al."The"STRING"database"in"2021："customizable"protein-"protein"networks，"and"functional"characterization"of"user-"uploaded"gene/measurement"sets[J]."Nucleic"Acids"Res，"2021，"49（D1）："D605–d612.

[12] CHAO"C，"JIAN"Y"Y，"ZHAO"X"Y，"et"al."The"involvement"of"hsa_circ_0000417"in"the"development"of"hypospadias"by"regulating"AR[J]."Differentiation，"2020，"116："9–15.

[13] MANSOUR"M"M，"GOYAL"H"O，"BRADEN"T"D，"et"al."Activation"of"penile"proadipogenic"peroxisome"proliferator-"activated"receptor"gamma"with"an"estrogen："interaction"with"estrogen"receptor"alpha"during"postnatal"development[J]."PPAR"Res，"2008，"2008："651419.

（下轉(zhuǎn)第18頁）