通过构建长链非编码RNA（lncRNA）、microRNA（miRNA）、mRNA参与的竞争性内源RNA（ceRNA）网络，进一步探究lncRNA在非梗阻性无精子症（NOA）发生发展过程中的调控机制。
收集2017年6月至2018年6月南方医科大学南方医院生殖中心收治的无精子症患者的睾丸组织。本课题组前期使用NOA患者（试验组）和生精正常的梗阻性无精子症（OA）患者（对照组）的睾丸组织进行基因芯片检测，筛选出差异表达的lncRNA和mRNA。本研究采用在NOA患者中差异表达的lncRNA和mRNA构建NOA ceRNA调控网络，并通过GeneMANIA数据库对ceRNA中编码蛋白的mRNA进行蛋白质-蛋白质互作（PPI）网络构建。然后，采用WebGestalt对PPI中编码蛋白的基因进行功能富集分析。最后，通过实时荧光定量聚合酶链式反应（qRT-PCR）和双荧光素酶报告基因试验进行验证。
NOA的ceRNA调控网络共包括21个节点和26条边，其中4个lncRNA节点，13个miRNA节点以及4个mRNA节点。PPI网络分析发现有19个蛋白与ceRNA网络中mRNA编码蛋白有相互作用。功能富集分析结果显示这些编码蛋白显著性富集于戊糖代谢相关的生物学功能和信号通路。进一步的试验验证，证实lncRNA ANXA2P3通过与miR-613和miR-206竞争性结合，进而抑制mRNA TKT的表达水平。
lncRNA通过ceRNA调控网络影响NOA的发生发展，为NOA的诊治提供了新靶点和治疗方案。
版权归中华医学会所有。
未经授权，不得转载、摘编本刊文章，不得使用本刊的版式设计。
除非特别声明，本刊刊出的所有文章不代表中华医学会和本刊编委会的观点。
非梗阻性无精子症（non-obstructive azoospermia，NOA）是男性生殖系统疾病中最为严重的情况之一，精子发生障碍是其主要病因[1]，由于具体的发病机制不详，故目前针对NOA缺乏明确的治疗靶点和有效的治疗措施。长链非编码RNA（long non-coding RNA，lncRNA）是一类长度>200 nt的非编码RNA分子，与人类众多疾病的发生发展都有着密切的联系[2]。LncRNA可通过多种作用机制发挥其重要的生物学功能，其中竞争性内源RNA（competing endogenous RNA，ceRNA）假说[3]揭示了lncRNA作用的新机制。近年来大量研究结果表明，lncRNA在睾丸组织中呈现高度特异性表达，并且可作为重要的调控分子参与精子发生过程[4,5,6,7,8]。然而关于lncRNA在NOA中的作用机制研究却鲜有报道，且lncRNA在NOA中的ceRNA调控机制也尚未明确。因此，本课题组前期使用NOA患者和生精正常的梗阻性无精子症（OA）患者的睾丸组织进行基因芯片检测，筛选出差异表达的lncRNA和mRNA，初步探究了lncRNA在NOA发生发展中的作用机制。本研究拟通过构建ceRNA调控网络，蛋白互作网络及功能富集分析，并通过实时荧光定量聚合酶链式反应（qRT-PCR）验证、双荧光素酶报告基因试验等进一步探究lncRNA作为ceRNA在NOA发生发展中的作用，以期揭示NOA相关的发病机制。
本项目由南方医科大学南方医院伦理委员会授批立项（批准号：NFEC-2017-102），并于试验开展前对每例参与者进行知情告知，同时获取所有研究对象自愿签署的同意书。根据世界卫生组织（WHO）《人类精液检查与处理实验手册》（第5版）推荐的诊断标准，本课题组将2017年6月至2018年6月于南方医科大学南方医院生殖中心就诊的无精子症患者纳入研究。试验开展前，对所有研究对象进行全面的男性科学测试，包括病史采集、体格检查、精液分析、激素水平检测、染色体核型分析及Y染色微缺失检测及组织病理学检查等，排除：（1）既往有隐睾、小睾丸、睾丸炎、睾丸创伤、腮腺炎等病史；（2）有射线暴露、毒物接触史等；（3）染色体异常、Y染色体微缺失；（4）激素水平异常的患者；（5）家族遗传病史、重大疾病史及手术史等，通过睾丸穿刺取精术（TESA）抽取每例患者约5 mg的睾丸组织用冰冷的磷酸盐缓冲液（PBS）清洗2~3遍后放入去RNA酶的冻存管并立即投入液氮保存。依据病理组织学检查以鉴别诊断NOA组和CON组。所有组织标本均由本院2位病理科医师双盲法确认。本课题组前期研究使用NOA患者（病例组，NOA组）和生精正常的OA患者（对照组，CON组）各4例的睾丸组织进行基因芯片检测。
本研究使用前期通过基因芯片筛选出的差异表达lncRNAs和mRNAs构建NOA ceRNA调控网络。利用miRcode（http://mircode.org/）数据库预测lncRNA与miRNA的相互作用，并使用TargetScan（http://www.targetscan.org）、miRDB（http://mirdb.org/）等数据库对miRNA-mRNA的靶向关系进行预测。而后对mRNA-miRNA和miRNA-lncRNA的靶向关系进行超几何检验，从而得到mRNA-miRNA-lncRNA的ceRNA调控网络。最后，使用Cytoscape（v3.2.0）对NOA中的ceRNA调控网络进行可视化。
PPI网络采用Gene MANIA数据库[9]进行构建，通过对ceRNA调控网络中编码蛋白的mRNA进行PPI网络构建，获得与ceRNA调控网络中mRNA编码蛋白相互作用的蛋白。
基因功能和信号通路富集分析采用WebGestalt[10]工具进行分析，分别选择Gene Oncology[11]和KEGG信号通路[12]数据库对PPI网络中的编码基因进行功能和信号通路富集分析。
根据ceRNA网络构建结果，对相关基因的表达水平进行qRT-PCR检测。NOA组和CON组各取10例睾丸组织（具体同方法1），组织裂解后，使用TRIzol法提取组织样本中的总RNA，并将总RNA逆转录为cDNA。使用Bio-Rad公司的CFX96 Real time PCR进行qRT-PCR检测，将β-肌动蛋白（H）设为内参基因，qRT-PCR结果采用2-∆∆Ct表示和计算。按照primer 6.0设计引物，目标基因引物序列见表1。
目标基因引物序列
目标基因引物序列
| 基因名称 | 双向引物序列 | 产物长度（bp） |
|---|---|---|
| β-肌动蛋白（H） | F:5 ′ GTGGCCGAGGACTTTGATTG3 ′ | 73 |
| R:5′CCTGTAACAACGCATCTCATATT3′ | ||
| ZNF674-AS1 | F:5′GTTCGCTCTAAGGGGAGAAGGA3′ | 113 |
| R:5′GGGTCTGGCTGTATTGAGTTTG3′ | ||
| KLF7-IT1 | F:5′ TGACCTCAGTTTCCTCATTTCT 3′ | 117 |
| R:5′ AATTCACACCAATTCCAAATAA 3′ | ||
| KGFLP1 | F:5′ CAGTGCCCTTCAGTTAAAGAC 3′ | 119 |
| R:5′ ATCCTGTCAAAATGATGAGAAG 3′ | ||
| ANXA2P3 | F:5′ ATGGGTCAGTCAAAGCCTACA 3′ | 139 |
| R:5′ TCCTCTCTCCGAGCATTCCT 3′ | ||
| KPNA4 | F:5′ CTATGAGAAGACAACGAAATGA 3′ | 222 |
| R:5′ GCTTCCTAGCAGCTTGAACT 3′ | ||
| EMB | F:5′ ATTCGCCTTTTACAAGTCCAC 3′ | 207 |
| R:5′ CATCTTTTTTCCAAGTCACATT 3′ | ||
| TAOK1 | F:5′ CTTGAAACTCAGCGTAACAAT 3′ | 95 |
| R:5′ TTGGACATCACTTTAGCCTCT 3′ | ||
| TKT | F:5′ CTATTGCTTGCTGGGAGACG 3′ | 263 |
| R:5′ CTGTTGGCTGGTGCTTGG 3′ |
分别构建miRNA过表达对照质粒、hsa-mir-613的过表达质粒、hsa-mir-206的过表达质粒和TKT-3′UTR野生型和突变型双荧光素酶质粒、ANXA2P3野生型和突变型双荧光素酶质粒、对照双荧光素酶质粒，质粒构建及具体试验分组见表2，表3。将目的基因3′ UTR区域构建至报告基因荧光素酶后面，通过比较过表达miRNA后荧光素酶的活性变化，来定量反映miRNA对目的基因的抑制作用。简言之，将质粒转染到293T细胞系，每组3复孔，不同组间的荧光素酶活性检测通过Promega公司的GloMax生物发光检测仪检测。
双荧光素酶报告基因检测质粒构建表
双荧光素酶报告基因检测质粒构建表
| 质粒名称 | 质粒描述 |
|---|---|
| miRNA-NC | microRNA空载质粒，作为目的microRNA质粒阴性对照 |
| miRNA（hsa-mir-206） | microRNA载体质粒，表达目的microRNA[hsa-mir-206(15268-1)] |
| miRNA（hsa-mir-613） | microRNA载体质粒，表达目的microRNA[hsa-mir-613(50769-11)] |
| 3′ UTR-NC | 3′ UTR空载质粒，作为靶基因3′ UTR质粒阴性对照 |
| 3′ UTR（NM_001258028） | 靶基因3′ UTR质粒[TKT(50767-1)] |
| 3 ′ UTR-M（NM_001258028-mut） | 靶基因3′ UTR突变体质粒[TKT(50766-1)] |
| 3′ UTR（NR_001446） | 靶基因3′ UTR质粒[ANXA2P3(50765-1)] |
| 3 ′ UTR-M（NR_001446-mut） | 靶基因3′ UTR突变体质粒[ANXA2P3(50764-1)] |
| 阳参3′ UTR | TRAF6基因3′ UTR质粒 |
| 阳参miRNA | hsa-mir-146b载体质粒 |
双荧光素酶报告基因检测试验分组
双荧光素酶报告基因检测试验分组
| 组别 | 类型 |
|---|---|
| 试验组1 | 3′ UTR-NC+miRNA（hsa-mir-206） |
| 试验组2 | 3′ UTR（NM_001258028）+miRNA（hsa-mir-206） |
| 试验组3 | 3′ UTR-M（NM_001258028-mut）+miRNA（hsa-mir-206） |
| 试验组4 | 3′ UTR（NR_001446）+miRNA（hsa-mir-206） |
| 试验组5 | 3′ UTR-M（NR_001446-mut）+miRNA（hsa-mir-206） |
| 试验组6 | 3′ UTR-NC+miRNA（hsa-mir-613） |
| 试验组7 | 3′ UTR（NM_001258028）+miRNA（hsa-mir-613） |
| 试验组8 | 3′ UTR-M（NM_001258028-mut）+miRNA（hsa-mir-613） |
| 试验组9 | 3′ UTR（NR_001446）+miRNA（hsa-mir-613） |
| 试验组10 | 3′ UTR-M（NR_001446-mut）+miRNA（hsa-mir-613） |
| 阳参miRNA NC组 | 阳参3′ UTR+miRNA-NC |
| 阳参miRNA组 | 阳参3′ UTR+阳参miRNA |
数据采用SPSS 20.0软件进行统计学分析，计量数据采用
±s进行表示。两组差异统计学意义比较采用t检验，多组间比较采用单因素方差分析（One-way ANOVA），以P<0.05为差异有统计学意义。
本课题组前期研究通过基因芯片筛选出差异表达的lncRNA和mRNA（差异倍数fold chang≥2，P<0.05），初步探究了lncRNA在NOA中的作用。为了更深入地研究NOA中差异表达的lncRNA的作用机制，本研究对差异表达的lncRNA和mRNA进行lncRNA-miRNA- mRNA的ceRNA调控网络构建。lncRNA-miRNA- mRNA ceRNA调控网络是根据miRcode数据库中lncRNA-miRNA之间的关系，以及TargetScan、miRDB等数据库预测出的miRNA-mRNA的靶向关系，并利用差异表达的lncRNA和mRNA而建立。识别出具有miRNA结合位点的差异表达的lncRNA和mRNA后，通过建立成对的lncRNA-miRNA- mRNA（表4）和关系匹配对，以lncRNA（诱饵），mRNA（靶标）和miRNA（中心），构建lncRNA- miRNA-mRNA的ceRNA调控网络（图1）。ceRNA网络共包括21个节点和26条边。其中，有4个lncRNA节点，13个miRNA节点以及4个mRNA节点。ceRNA调控网络可从整体上直观的反映lncRNA-miRNA-mRNA之间的相互作用关系，快速定位在NOA调控网络中发挥作用的关键分子。
mRNA-miRNA-lncRNA匹配对统计结果
mRNA-miRNA-lncRNA匹配对统计结果
| mRNA | miRNA | lncRNA | Cor（mRNA-lncRNA） | P-val（mRNA-lncRNA） |
|---|---|---|---|---|
| KPNA4 | hsa-miR-544a | ZNF674-AS1 | 0.99409569 | 5.18E-08 |
| KPNA4 | hsa-miR-23b-3p | ZNF674-AS1 | 0.99409569 | 5.18E-08 |
| KPNA4 | hsa-miR-23c | ZNF674-AS1 | 0.99409569 | 5.18E-08 |
| KPNA4 | hsa-miR-23a-3p | ZNF674-AS1 | 0.99409569 | 5.18E-08 |
| EMB | hsa-miR-455-5p | KLF7-IT1 | 0.991881433 | 1.58E-07 |
| TAOK1 | hsa-miR-15b-5p | KGFLP1 | 0.990476539 | 2.75E-07 |
| TAOK1 | hsa-miR-15a-5p | KGFLP1 | 0.990476539 | 2.75E-07 |
| TAOK1 | hsa-miR-497-5p | KGFLP1 | 0.990476539 | 2.75E-07 |
| TAOK1 | hsa-miR-424-5p | KGFLP1 | 0.990476539 | 2.75E-07 |
| TAOK1 | hsa-miR-16-5p | KGFLP1 | 0.990476539 | 2.75E-07 |
| TAOK1 | hsa-miR-195-5p | KGFLP1 | 0.990476539 | 2.75E-07 |
| TKT | hsa-miR-613 | ANXA2P3 | 0.993815821 | 6.09E-08 |
| TKT | hsa-miR-206 | ANXA2P3 | 0.993815821 | 6.09E-08 |
通过对NOA ceRNA调控网络中编码蛋白的mRNA进行PPI网络构建（图2），本研究发现19个蛋白与ceRNA调控网络中mRNA编码蛋白有相互作用。其中，绝大多数的相互作用为物理间的相互作用（67.64%），其次是共表达（13.50%）。这揭示了NOA ceRNA网络中的编码蛋白主要通过物理间的相互作用和共表达机制与其他19个蛋白协同调控NOA的发生和发展。
为了进一步研究PPI网络中的蛋白所具有的生物学功能以及参与的信号通路，本研究对PPI网络中的编码蛋白基因进行GO和KEGG信号通路富集分析。结果发现这些基因主要富集于戊糖代谢过程、磷酸戊糖旁路、甘油醛-3-磷酸代谢过程等生物学过程（图3）；KEGG富集显示这些基因主要富集于磷酸戊糖信号通路、碳代谢信号通路等（图4），这说明NOA的发病机制与戊糖代谢的失调有着密切的关系。
NOA组和CON组各取10例睾丸组织，qRT-PCR法检验lncRNA ZNF674- AS1，KLF7-IT1，KGFLP1，ANXA2P3（图5A）和mRNA KPNA4，EMB，TAOK1，TKT（图5B）的表达。结果发现，相较于CON组，NOA组中lncRNA ANXA2P3和KGFLP1的表达显著性上调，差异有统计学意义（P<0.05），而lncRNA KGF7-IT1和ZNF674-AS1表达差异无统计学意义；mRNA EMB和TKT的表达显著上调，差异有统计学意义，而mRNA KPNA4和TAOK1的表达差异无统计学意义（图5）。这从试验角度证实了lncRNA ANXA2P3和KGFLP1以及mRNA EMB和TKT的失调在NOA的发生发展过程中可能发挥着重要的作用。
根据qRT-PCR结果选取lncRNA ANXA2P3-TKT作为后续研究对象。通过ceRNA调控网络和生物信息学结合位点预测，发现lncRNA ANXA2P3通过与has-miR-613以及has-miR-206相互作用，进而影响TKT的表达水平。has-miR-613和has-miR-206与TKT 3′ UTR的结合位点见图6，同时分析发现ANXA2P3与hsa-miR-613和hsa-miR-206的结合位点与TKT为同一段序列。这提示，lncRNA ANXA2P3可能通过has-miR-613和has-miR-206调节TKT表达水平，进而在NOA中发挥作用。
根生物信息学lncRNA-miRNA以及miRNA-mRNA结合位点预测结果，进一步采用双荧光素酶报告基因实验进行验证，检测结果显示（图7）：相比于阳性参照miRNA NC组，阳性参照miRNA组荧光素酶的活性有显著性下降（P<0.05），说明整个转染检测体系没有问题。试验组1和2相比、试验组1和4相比、试验组6和7相比以及试验组6和9相比，荧光素酶活性均发生显著性降低（P<0.05）。因此，证实了lncRNA ANXA2P3通过与miR-206和miR-613竞争性结合TKT的3′UTR，从而影响TKT的表达水平。
近年来不孕不育在全球的发生率呈逐年升高趋势。而NOA作为男性生殖系统疾病中最为严重的情况之一，迫切需要了解其发病机制，寻找明确的治疗靶点和有效的治疗措施。
随着生物信息学的发展，lncRNA对生长发育的影响逐渐彰显，可通过多种调控机制参与疾病的发生发展。目前，在男性生殖方面已有许多研究证实lncRNA与精子发生过程相关，如Wen等[5]通过构建果蝇模型，系统地鉴定出了睾丸特异性表达的lncRNA，并证实其在精子发生过程中起重要作用；Hu等[13]发现lncRNA Gm2044在精母细胞中高度表达，且通过与mRNA Utf1相互作用从而抑制其翻译表达；Wichman等[14]揭示了哺乳动物特异性lncRNA的失调是导致精子数目减少或造成不育的一种新机制。然而，lncRNA与NOA的相关性研究却鲜有报道，lncRNA在NOA中是否发挥着ceRNA作用及其机制也尚未明确。
我们前期通过对微阵列基因芯片分析，鉴定出了在NOA患者中差异表达的lncRNA，初步探究了lncRNA在NOA发生发展中的作用机制。为了进一步探究其作用机制，本研究对这些差异表达的lncRNA构建ceRNA调控网络，并通过PPI网络、功能富集分析、qRT-PCR验证、双荧光素酶报告基因试验等进一步探究了lncRNA参与调控的ceRNA网络在NOA发生发展中的作用及潜在机制。结果发现NOA的ceRNA调控网络共包括21个节点和26条边。其中，有4个lncRNA节点，13个miRNA节点以及4个mRNA节点。通过PPI网络预测发现，19个与其相互作用的蛋白。这些基因主要富集在与戊糖代谢相关的生物学功能和信号通路。进一步的试验验证，证实lncRNA ANXA2P3与TKT在NOA中显著上调，生物信息学分析结合双荧光素酶试验证实lncRNA ANXA2P3可竞争性结合miR-613和miR-206，进而消除miR-613和miR-206对TKT的沉默作用，从而上调TKT的表达。
LncRNA ANXA2P3（Annexin A2 pseudogene 3）是膜联蛋白A2假基因3，可作为信号通路的调节因子参与多种细胞功能的调节。肝细胞体外试验显示，过表达ANXA2P3后能显著的抑制细胞的增殖并且促进细胞的凋亡[15]。TKT（Transketolase）是转酮醇酶的编码基因，转酮醇酶在戊糖磷酸循环以及光合成的还原型戊糖磷酸循环中起着重要作用。此酶的作用是把磷酸酮糖上的乙酮醇基转移给磷酸醛糖。TKT在戊糖磷酸途径中的糖酵解中起着重要的作用。研究显示，TKT活性与癌症的发生相关。转酮醇酶活性在多种肿瘤中被检测到升高，利用转酮醇酶抑制剂处理肿瘤细胞后发现细胞增殖明显减少，而转酮醇酶激活剂的应用则可以刺激肿瘤生长[16]，这提示TKT参与细胞的生长与再生。在NOA中，lncRNA ANXA2P3通过与miR-613和miR-206相互结合，调控miR-613/miR-206对TKT表达水平的抑制作用，从而引起戊糖代谢功能和信号通路的失调，促进了NOA的发展。
综上所述，本研究通过构建NOA中差异表达的lncRNA参与调控的ceRNA网络，发现了lncRNA ANXA2P3通过与miR-613和miR-206相互竞争性结合而调控TKT的表达从而影响NOA的作用机制，这为揭示NOA的潜在病理机制提供了理论依据和新思路。
利益冲突 所有作者均声明不存在利益冲突
