Construction for Long Non-Coding RNA (lncRNA)-Associated Competing Endogenous RNA (ceRNA) Network in Human Retinal Detachment (RD) with Proliferative Vitreoretinopathy (PVR).

BACKGROUND The aim of this study was to analyze the long non-coding RNA (lncRNA)-associated competing endogenous RNA (ceRNA) network in human retinal tissues following detachment with proliferative vitreoretinopathy (PVR). MATERIAL AND METHODS Expression data of 19 human detached retinas with PVR and 19 normal retinas from postmortem donors were downloaded from Gene Expression Omnibust (GEO) database (GSE28133). The R package "limma" was utilized to discriminate the dysregulated lncRNA and mRNA profiles. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses of differentially expressed mRNAs were performed using R packages "Clusterprofiler." The ceRNA network of dysregulated genes was constructed by using mircode, miRDB, miRTarBase and TargetScan databases, and was visualized by Cytoscape v3.6.1. RESULTS A total of 23 lncRNAs and 994 mRNAs were identified significantly expressed between the human detached retinas with PVR and the normal retina tissues, with thresholds of |log₂FoldChange| >1.0 and adjusted P-value <0.05. The constructed ceRNA network (lncRNA-miRNA-mRNA regulatory axis) included 9 PVR-specific lncRNAs, as well as 27 miRNAs and 73 mRNAs. CONCLUSIONS We demonstrated the differential lncRNA expression profile and constructed a lncRNA-associated ceRNA network in human detached retinas with PVR. This may ferret out an unknown ceRNA regulatory network in human retinal detachment with PVR.

Background

Retinal detachment (RD) occurs when the neurosensory retina layer separates from the retinal pigment epithelium. RD involving the foveal center can lead to profound loss of vision [1]. In most cases, RD occurs following a full thickness retinal break, allowing the ingress of fluid from the vitreous cavity to the subretinal space, which can result in retinal separation and is so-called “rhegmatogenous” [1,2]. A systematic review showed that rhegmatogenous RD incidence varied between 6.3 and 17.9 per 100 000 population, and was strongly associated with myopia, increasing age and certain vitreoretinal degenerations [3]. In addition to the loss of photoreceptors following RD, an inflammatory response also develops. RD triggers cell migration and proliferation and production of extracellular matrix proteins, which in turn results in the accumulation of vitreal and periretinal membranes, both hallmarks of proliferative vitreoretinopathy (PVR) [4]. In total, PVR occurs in 5–10% RD cases [5]. PVR can be detected during the late presentation of RD and complicate the post-operative period after a surgery for RD, leading to a re-detachment or limiting visual recovery [6,7]. Clinical studies using adjuvant therapy for the treatment of PVR have been conducted, including anti-inflammatory agents, anti-neoplastic/anti-proliferative agents, anti-growth factor pathway inhibitors and antioxidants, etc. [6,8-11]. However, the results are often contradictory or inconclusive with only limited success.

Long non-coding RNA (lncRNA) is a group of non-coding transcripts 200–10 000 bp in length but lacks significant protein-coding capacity. They interact in a regulatory manner before, during and after transcription [12]. Salmena et al. made a competing endogenous RNA (ceRNA) hypothesis that lncRNA, mRNA and other RNAs might act as microRNA (miRNA) sponges to inhibit miRNA function through sharing the miRNA response elements (MREs), which is a complicated post-transcriptional regulatory network [13]. There is much evidence to support this hypothesis [14-16]. On one hand, miRNAs can combine with their target mRNAs and inhibit their expression in the ceRNA network. On the other hand, lncRNA can regulate gene encoding protein level and following cell biology by competing with miRNAs [15,16]. Besides, each miRNA can influence up to hundreds of expressions of transcription, while each RNA transcription with different MREs can be targeted by multiple miRNAs [16,17].

In recent years, a growing number of studies have proven that the lncRNA-miRNA-mRNA regulation network is involved in the course of many diseases, including the tumors, Alzheimer disease and dermatitis etc. [12,16,18,19]. However, whether the ceRNA network takes part in human detached retinas with PVR has not been studied. Here, we downloaded expression data from Gene Expression Omnibust (GEO) database (GSE28133). It contained 19 specimens from patients for severe retinal detachment with PVR and 19 normal control retina specimens from postmortem donors. Bioinformatics analysis and the ceRNA network were conducted to explore the pathological mechanism and potential therapeutic targets of retinal detachment with PVR.

Material and Methods

Tissue samples from the GEO database and bioinformatics analysis

The RNA expression data of tissue samples were obtained from GEO database. The dataset (GSE28133) contained 19 specimens from patients for severe retinal detachment with PVR and 19 normal control retina specimens from postmortem donors. Clinical and pathological features of patients were described [2]. All the samples were analyzed using Affymetrix Human Genome U133 Plus 2.0 Array.

Differentially expressed gene analysis

We used the “limma” package in R software to identify the differentially expressed genes (lncRNAs and mRNAs) between human detached retinas with PVR and the control group, with |log2FoldChange (FC)|>1.0 and adjusted P-value <0.05. Besides, mRNA and lncRNA annotation was performed using the Encyclopedia of DNA Elements (ENCODE) with ENSEMBL.

GO and KEGG functional enrichment analysis

To understand the potential biological functions and processes of differentially expressed genes, Gene Ontology database (GO, ) and Kyoto Encyclopedia of Genes and Genomes (KEGG, ), as well as the “clusterProfilerGO” package and “clusterProfilerKEGG” in R software, were utilized to perform GO and KEGG pathway analysis. Records with P-value <0.05 and enrichment >2.0 were preserved.

Construction of the lncRNA-miRNA-mRNA ceRNA and protein–protein interaction (PPI) networks

Based on the hypothesis that lncRNA could sponge the common miRNA and thus prevent miRNA from binding to their target genes [20], a ceRNA network was constructed. The mircode database () was used for lncRNA to predict their targeted miRNAs. The miRNA-mRNA interactions were predicted by the miRTarBase (), miRDB () and TargetScan (). Finally, the predicted mRNAs were cross matched with the differentially expressed ones, and the mRNA with no negatively regulated lncRNA or miRNA were abandoned. The lncRNA, miRNA and mRNA with |log2FC| >1.0 and adjusted P-value <0.05 were preserved. Construction and visualization of the lncRNA-miRNA-mRNA ceRNA network were conducted using Cytoscape v3.6.1. Besides, PPI network of mRNAs involved in the ceRNA network was set up with high confidence (0.700) by String ().

Results

Differentially expressed lncRNAs and mRNAs between the human detached retinas with PVR and normal controls

Clinical characteristics of patients were described in the previous study [2]. In total, 994 genes were significantly changed between human detached retinas with PVR and normal human retinas with thresholds of |log2FC| >1.0 and adjusted P-value <0.05, as is shown in the volcano map (Figure 1A). Among them, 23 genes were lncRNA, with 16 genes downregulated and 7 genes upregulated. The aberrantly expressed lncRNAs can be viewed in the heatmap (Figure 1B) and their gene IDs and log2FC are listed (Table 1), in which the 9 lncRNAs taking part in the ceRNA network were marked with asterisk keys. Besides, 971 mRNAs were identified as differentially expressed in which 655 (67.46%) were upregulated and 316 (32.54%) were downregulated. The differentially expressed mRNA levels were also presented in the heatmap (Figure 1C).

Figure 1

Differentially expressed genes between the human detached retinas with PVR and normal controls. (A) The volcano plot of the differently expressed genes. The x axis is -log10(adjusted P-value) while the y axis is log2FC. Black dots indicated non-differently expressed genes. The red dots indicate upregulated genes while the green dots indicate downregulated ones. (B) Heatmap of differently expressed lncRNAs. (C) Heatmap of differently expressed mRNAs. (B, C) The x axis is for the sample serial number and the y axis is for the differentially expressed gene. Red blocks represent upregulation and green blocks represent downregulation.

Table 1

Dysregulated lncRNAs between human detached retinas with PVR and normal controls.

lncRNA	Gene ID	Expression change	logFC (PVR/N)	Adj. P-value
HCP5*	10866	Up-regulation	1.947938297	1.53E-09
LINC00623	728855	Up-regulation	1.230129889	5.99E-08
LINC01094	100505702	Up-regulation	1.217325671	5.36E-06
PSMB8-AS1	100507463	Up-regulation	1.210172661	1.66E-16
MIAT*	440823	Up-regulation	1.209215343	1.06E-06
CRNDE*	643911	Up-regulation	1.055524554	1.51E-08
BLACAT1	101669762	Up-regulation	1.048582535	9.62E-07
ZNF571-AS1	100507433	Down-regulation	−1.031576283	1.01E-06
AP000462*	---	Down-regulation	−1.066771719	0.000173
SH3BP5-AS1	100505696	Down-regulation	−1.08329906	4.33E-07
JPX*	554203	Down-regulation	−1.203614666	1.33E-07
PGM5-AS1	572558	Down-regulation	−1.248531249	6.67E-06
OTX2-AS1	100309464	Down-regulation	−1.26206851	1.89E-05
RNF139-AS1	101927612	Down-regulation	−1.383351463	1.96E-10
TDRG1*	732253	Down-regulation	−1.399211813	2.00E-05
ELOVL2-AS1	100506409	Down-regulation	−1.472180825	6.24E-09
LINC01137	728431	Down-regulation	−1.613551017	1.38E-09
FAM13A-AS1*	285512	Down-regulation	−1.689264402	5.73E-09
SPAG5-AS1*	100506436	Down-regulation	−1.774851965	1.87E-11
MIR124-2HG	100130155	Down-regulation	−1.894825134	3.78E-10
AC005592*	101926975	Down-regulation	−2.035310985	1.55E-07
RBFADN	100506070	Down-regulation	−2.175539179	4.42E-13
LINC00844	100507008	Down-regulation	−2.290846454	1.01E-08

PVR – proliferative vitreoretinopathy; N – normal human retinas.

Represents the lncRNA taking part in the ceRNA network.

Functional enrichment analysis of differentially expressed mRNAs

GO and KEGG pathway analyses allow for the molecular function and signal pathway annotation of dysregulated mRNAs. Records with P-value <0.05 and enrichment >2.0 were preserved. The complete list of the 74 terms of GO analysis is presented in Supplementary Table 1. The top 20 GO terms were visualized in the bubble diagram and the most enriched pathways were the cell adhesion molecule binding (GO: 0050839), actin binding (GO: 0003779) and enzyme inhibitor activity (GO: 0004857), which contained 60, 49, and 44 genes respectively. As shown in Figure 2A), 13 GO terms of the top 20 took part in the binding function of the cell.

Figure 2

Functional enrichment analysis of differentially expressed mRNAs in the human detached retinas with PVR. (A) Top 20 most significant biological processes in GO analysis. (B) Top 20 most significant KEGG pathways.

The KEGG pathway analysis indicated that 39 pathways were associated with differentially expressed mRNAs (Supplementary Table 2). The top 20 KEGG pathways enriched are visualized in the bubble diagram (Figure 2B). It shows that infection probably played an important role in retinal detachment with PVR, including human papillomavirus, Epstein-Barr virus, human cytomegalovirus, as well as tuberculosis, etc. In addition, phagosome, focal adhesion, and complement and coagulation cascade were also suggested as participate in PVR progress, containing 32, 29, and 22 genes respectively.

Construction of a ceRNA network and PPI network in the human detached retinas with PVR

To further understand how lncRNA regulates mRNA through combining with miRNA in human detached retinas with PVR, a lncRNA-miRNA-mRNA (ceRNA) network was constructed. Consequently, we found that 9 lncRNAs interacted with the 27 miRNAs in the ceRNA network using the miRcode database (Table 2). A total of 1269 miRNA-targeted mRNAs based on the 27 miRNAs was predicted through miRDB, miRTarBase, and TargetScan databases. MiRNA-targeted mRNAs which were not contained in the differentially expressed mRNAs, which were discarded, and the 73 mutual mRNAs which were preserved (Figure 3A, Table 3). Finally, 9 lncRNAs, 27 miRNAs, and 73 mRNAs constituted the ceRNA network (Figure 4). Besides, KEEG pathways enriched by mRNAs in ceRNA network showed that most of them were bound with the PI3K-Akt signaling pathway and human cytomegalovirus infection (Figure 3B). The PPI network constructed for the 73 mRNAs furnished 21 genes with high confidence 0.700 (Figure 3C).

Table 2

lncRNAs and specific targeted miRNAs in ceRNA network.

lncRNA	miRNA
TDRG1	miR-17-5p, miR-20b-5p, miR-27a-3p, miR-125a-5p, miR-125b-5p
HCP5	miR-137, miR-139-5p, miR-140-5p, miR-17-5p, miR-20b-5p, miR-216b-5p, miR-22-3p, miR-23b-3p, miR-24-3p, miR-363-3p, miR-1297, miR-27a-3p, miR-107, miR-425-5p, miR-125a-5p, miR-125b-5p, miR-10a-5p
JPX	miR-301b-3p, miR-140-5p, miR-193a-3p, miR-216b-5p, miR-23b-3p, miR-24-3p, miR-363-3p, miR-449c-5p, miR-129-5p
MIAT	miR-301b-3p, miR-139-5p, miR-140-5p, miR-17-5p, miR-20b-5p, miR-206, miR-613, miR-216b-5p, miR-22-3p, miR-23b-3p, miR-24-3p, miR-363-3p, miR-27a-3p, miR-107, miR-449c-5p, miR-125a-5p, miR-125b-5p, miR-10a-5p,miR-129-5p
SPAG5-AS1	miR-17-5p, miR-20b-5p, miR-429, miR-217, miR-24-3p, miR-107, miR-425-5p
AC005592	miR-139-5p, miR-17-5p, miR-20b-5p, miR-206, miR-613, miR-216b-5p, miR-107, miR-425-5p, miR-125a-5p, miR-125b-5p, miR-137, miR-140-5p, miR-429, miR-23b-3p, miR-363-3p, miR-33a-3p, miR-10a-5p, miR-129-5p
CRNDE	miR-140-5p, miR-142-3p, miR-193a-3p, miR-216b-5p, miR-217, miR-22-3p, miR-23b-3p, miR-363-3p, miR-1297, miR-27a-3p, miR-129-5p
FAM13A-AS1	miR-137, miR-139-5p, miR-142-3p, miR-17-5p, miR-20b-5p, miR-217, miR-22-3p, miR-23b-3p, miR-24-3p, miR-363-3p, miR-107, miR-449c-5p, miR-125a-5p, miR-125b-5p, miR-129-5p
AP000462	miR-17-5p, miR-20b-5p, miR-1297, miR-206, miR-613, miR-22-3p, miR-24-3p, miR-449c-5p, miR-1297

Figure 3

mRNAs in ceRNA network. (A) Venn diagram of differentially expressed mRNAs in ceRNA network. (B) KEEG pathways enriched by mRNAs in ceRNA network. (C) The protein-protein interaction network constructed for mRNA involved in ceRNA network.

Table 3

miRNAs and specific targeted mRNAs in ceRNA network.

miRNA	mRNA
miR-107	VCAN, TGFBR3, GABRB1, ITGA2, LCOR, PLEKHF2,PCSK5
miR-10a-5p	ELOVL2
miR-125a-5p	TMEM136, EIF4EBP1, MTUS1, TMEM136, STAT3, LIPA, EIF4EBP1
miR-129-5p	SPRY4, C1S
miR-1297	ADM, CKS2, MAN2A1
miR-137	GLIPR1
miR-139-5p	LCOR, ZBTB34
miR-140-5p	LAMC1
miR-142-3p	EGR2, ZNF217, LCOR, MTUS1
miR-17-5p	TMEM123, CDKN1A, PLS1, NETO2, FJX1, CCND1, FAM57A, FAM129A, RAPGEF4, SLC16A9, CYBRD1, PPP3R1, PLEKHO2, CADM2, BTN3A1, STAT3, WEE1, TXNIP, SSX2IP, LCOR, TMEM138
miR-193a-3p	LAMC1
miR-206	GJA1, KCNJ2, CERS2, BDNF, SFRP1, WEE1
miR-20b-5p	TXNIP, CADM2, RAPGEF4, CCND1, TMEM123, SLC16A9, PLEKHO2, FJX1, SSX2IP, STAT3, CYBRD1, PPP3R1, CDKN1A, PLS1, BAMBI, NETO2, FAM129A
miR-216b-5p	SMAD1
miR-217	NR4A2
miR-22-3p	RGS2, CSF1R
miR-23b-3p	GJA1
miR-24-3p	GBA2, SCML1, MLEC, ZNF217, ADD1, FSCN1
miR-27a-3p	PPIF, TGFBR3, LPCAT1, PHLPP2, PLXND1, ADD1, WEE1, ADORA2B
miR-301b-3p	GRB10, PRUNE2, TRPC3, IRF1
miR-33a-3p	PPP3R1
miR-363-3p	CNNM4, PDPN, GFPT2, PHLPP2, LHFPL2
miR-425-5p	LCOR, RAB31
miR-429	TPD52L1
miR-449c-5p	MYC
miR-613	CERS2, WEE1

Figure 4

CeRNA network in human detached retinas with PVR. Diamonds represent lncRNAs, triangles represent miRNAs and rounds represent protein-coding genes. Gray edges indicate lncRNA-miRNA-mRNA interactions.

Discussion

PVR is known to complicate the post-operative period after a surgery for RD, leading to a re-detachment or limiting visual recovery. Traditional adjuvant therapy included anti-proliferative agents, anti-inflammatory agents, and anti-growth factor pathway inhibitors. But so far, studies for the treatment of PVR report limited success [6,8-11]. The previous study detected expression profile of 19 PVR patients and controls, but only focused on part of differentially expressed mRNAs [2]. Neither GO and KEGG functional analysis nor construction of ceRNA network were established. Thus, we re-analyzed the published microarray data of GEO database (GSE28133).

In our analysis, we found 971 mRNAs significantly differentially expressed in the human detached retinas with PVR compared to normal human retinas. The GO pathways enriched by the differentially expressed mRNA showed that main pathways of the top 20 took part in the binding function of the cell, such as cell adhesion molecule binding which contained CLIC1 (logFC=2.76), S100A11 (logFC=4.31) and ICAM1 (logFC=3.06), etc. Chloride intracellular channel 1 (CLIC1), and S100 calcium-binding protein A11 (S100A11) were positively correlated with cell proliferation, invasion, and migration and angiogenesis [21,22]. Intercellular adhesion molecule-1 (ICAM1) is involved in the adhesion of leukocytes to the blood vessel wall [23]. Such significantly dysregulated genes involved in cell adhesion and binding function may have an important role in PVR pathogenesis and deserve more exploration. The KEGG pathways showed that infection and Phosphoinositide 3-kinase (PI3K)-Akt signaling pathway come out in front in RD tissues with PVR. Cytomegalovirus, chlamydia trachomatis, and human immunodeficiency virus have been confirmed in many RD cases [24-27]. These are in conformance with our analysis and remind us of necessary etiological detection for tissues in vitreous cavity and possible targeted treatment. PI3K plays a crucial role as a mediator of growth factor signaling, cell proliferation, cell survival, and apoptotic inhibition. As the major component of the extracellular matrix in PVR, type I collagen was found to be regulated by the PI3K/Akt pathway in human retinal pigment epithelial cells [28]. Treatment targeting PI3K/Akt pathway to prevent PVR after surgery is worth being studied.

The miRNA is an extensive class of endogenous, noncoding and single-strand RNAs with 18–24 nucleotides that negatively regulates gene expression through interacting with the 30-untranslated regions (30UTR) of their target mRNAs [29]. Thus, miRNAs have essential roles in homeostasis and pathogenesis [29]. In the eye, various miRNAs could act on the retina and have an important role in neuroprotection and angiogenesis [30-32]. In our analysis, 27 miRNAs were predicted by the 9 lncRNAs in the ceRNA network using the miRcode database. MiR-107, miR-125a-5p, miR-17-5p, miR-20b-5p, and miR-27a-3p interrelated with 7 or more protein-coding genes, which may play a more important role in ceRNA network. These aberrantly expressed miRNAs also played key roles in multiple biological processes of various diseases [33-35]. For example, miR-107 was reported to inhibit cell migration and invasion by modulating Notch2 expression and regulate autophagy and apoptosis by targeting TRAF3 [36,37]. However, the influence of these miRNAs on PVR is rarely explained.

Compared with protein-coding genes and miRNAs, lncRNAs have significant advantages as prognostic biomarkers or therapeutic targets [38,39]. LncRNA regulates gene encoding protein level and participates in the regulation of cell biology through competing with miRNAs in the ceRNA network. There were 23 lncRNAs that were detected to be significantly differentially expressed in the human detached retinas with PVR and 9 of them were involved in the ceRNA network. The most upregulated lncRNA in the ceRNA network is reported to be HCP5. HCP5 was reported to express mainly in the immune system and often considered to be associated with herpes zoster and cancers [40,41]. HCP5 can promote cancer via the PI3K/AKT pathway, which was found to be the most enriched in the KEGG pathways by mRNAs in the ceRNA network [40]. In our analysis, HCP5 was also predicted to interact with 16 miRNAs. Other notable lncRNAs in the ceRNA network reported have included TDRG1 and MIAT [42,43]. TDRG was reported to promote the proliferation and progression of cells through PI3K/Akt/mTOR signaling [42]. Targeting MIAT was found to protect against myocardial hypoxia/reoxygenation injury, but its function in RD and PVR has not been studied [43]. In addition, MALAT1 was also found to be significantly upregulated in the fibrovascular membranes and the peripheral blood samples of PVR patients. In vitro studies revealed a critical role of MALAT1 in RPE proliferation and migration [44]. However, studies that focus on the function of lncRNA as miRNA sponges in human RD and PVR are in deficiency.

Our analysis revealed how specific lncRNAs interact with miRNAs and mRNAs through the successful construction of a lncRNA-miRNA-mRNA network in human detached retinas with PVR. This may reveal unknown pathological mechanisms of RD with PVR, which may help to provide potential therapeutic targets. However, we did not conduct quantitative real-time polymerase chain reaction analysis to validate the results of our bioinformatics analysis. The results would be more convincing with such verification experiments.

Conclusions

We demonstrated the differential expression profiles of lncRNAs, and we constructed a lncRNA-associated ceRNA network in human retinal detachment with PVR. Our analysis may contribute to increased understanding of the pathogenesis of human retinal detachment with PVR and provide novel lncRNAs as potential therapeutic targets.

Supplementary Table 1

GO pathways enriched by the differentially expressed coding genes.

ID	Description	Count	P-value
GO: 0005201	Extracellular matrix structural constituent	22	2.10E-11
GO: 0005178	Integrin binding	26	1.92E-10
GO: 0050839	Cell adhesion molecule binding	60	2.62E-10
GO: 0019838	Growth factor binding	26	2.98E-09
GO: 0003779	Actin binding	49	2.78E-08
GO: 0042605	Peptide antigen binding	11	5.64E-08
GO: 0061134	Peptidase regulator activity	31	2.34E-07
GO: 0004857	Enzyme inhibitor activity	44	1.17E-06
GO: 0061135	Endopeptidase regulator activity	25	4.70E-06
GO: 0048407	Platelet-derived growth factor binding	6	6.50E-06
GO: 0042277	Peptide binding	32	1.06E-05
GO: 0033218	Amide binding	34	2.05E-05
GO: 0004866	Endopeptidase inhibitor activity	23	2.51E-05
GO: 0051015	Actin filament binding	23	2.76E-05
GO: 0044548	S100 protein binding	6	3.70E-05
GO: 0005539	glycosaminoglycan binding	26	4.10E-05
GO: 0002020	Protease binding	19	4.63E-05
GO: 0008201	Heparin binding	21	5.60E-05
GO: 0030414	Peptidase inhibitor activity	23	5.65E-05
GO: 0008191	Metalloendopeptidase inhibitor activity	6	5.90E-05
GO: 0005518	Collagen binding	12	7.07E-05
GO: 0023026	MHC class II protein complex binding	6	9.03E-05
GO: 0031406	Carboxylic acid binding	23	0.000101
GO: 0043177	Organic acid binding	23	0.000119
GO: 0019864	IgG binding	5	0.000123
GO: 0016641	Oxidoreductase activity, acting on the CH-NH2 group of donors, oxygen as acceptor	6	0.000133
GO: 0004859	Phospholipase inhibitor activity	5	0.000202
GO: 0019955	Cytokine binding	15	0.000229
GO: 0001968	Fibronectin binding	7	0.000249
GO: 0030507	Spectrin binding	7	0.000249
GO: 0050840	Extracellular matrix binding	10	0.000259
GO: 0008020	G-protein coupled photoreceptor activity	5	0.000315
GO: 0023023	MHC protein complex binding	6	0.000366
GO: 0005504	Fatty acid binding	7	0.00041
GO: 0030674	Protein binding, bridging	20	0.000439
GO: 0016504	Peptidase activator activity	8	0.000464
GO: 0005520	Insulin-like growth factor binding	7	0.000517
GO: 0005546	Phosphatidylinositol-4,5-bisphosphate binding	11	0.000544
GO: 0030246	Carbohydrate binding	27	0.000547
GO: 0016638	Oxidoreductase activity, acting on the CH-NH2 group of donors	6	0.000646
GO: 0019865	Immunoglobulin binding	6	0.000646
GO: 0009881	Photoreceptor activity	5	0.000674
GO: 0070492	Oligosaccharide binding	5	0.00094
GO: 0033293	Monocarboxylic acid binding	10	0.000973
GO: 0042578	Phosphoric ester hydrolase activity	33	0.001054
GO: 0046982	Protein heterodimerization activity	42	0.001109
GO: 0019203	Carbohydrate phosphatase activity	4	0.001111
GO: 0035325	Toll-like receptor binding	4	0.001111
GO: 0050308	Sugar-phosphatase activity	4	0.001111
GO: 0008329	Signaling pattern recognition receptor activity	5	0.001276
GO: 0038187	Pattern recognition receptor activity	5	0.001276
GO: 0055102	Lipase inhibitor activity	5	0.001276
GO: 0019956	Chemokine binding	6	0.001345
GO: 0015298	Solute: cation antiporter activity	7	0.001425
GO: 0004197	Cysteine-type endopeptidase activity	12	0.001529
GO: 0071949	FAD binding	6	0.001673
GO: 0050786	RAGE receptor binding	4	0.001676
GO: 0043531	ADP binding	7	0.001704
GO: 0060090	Molecular adaptor activity	20	0.001883
GO: 0005507	Copper ion binding	9	0.00199
GO: 0003785	Actin monomer binding	6	0.00206
GO: 1901681	Sulfur compound binding	23	0.002196
GO: 0098641	Cadherin binding involved in cell-cell adhesion	5	0.002201
GO: 0001540	Amyloid-beta binding	9	0.002258
GO: 0036041	Long-chain fatty acid binding	4	0.002413
GO: 0045296	Cadherin binding	29	0.002488
GO: 0043394	Proteoglycan binding	7	0.002799
GO: 0031994	Insulin-like growth factor i binding	4	0.003346
GO: 0032561	Guanyl ribonucleotide binding	33	0.003348
GO: 0019001	Guanyl nucleotide binding	33	0.003483
GO: 0001848	Complement binding	5	0.003539
GO: 0098631	Cell adhesion mediator activity	7	0.003784
GO: 0015297	Antiporter activity	11	0.003925
GO: 0008236	Serine-type peptidase activity	25	0.004086

Supplementary Table 2

KEGG pathways enriched by the differentially expressed coding genes.

ID	Description	Count	P-value
hsa04744	Phototransduction	17	1.47E-14
hsa05150	Staphylococcus aureus infection	21	1.87E-10
hsa04145	Phagosome	32	2.46E-10
hsa04610	Complement and coagulation cascades	22	6.39E-10
hsa05169	Epstein-Barr virus infection	34	2.67E-08
hsa05152	Tuberculosis	31	6.42E-08
hsa05145	Toxoplasmosis	23	1.69E-07
hsa05133	Pertussis	18	3.54E-07
hsa04612	Antigen processing and presentation	18	4.37E-07
hsa05140	Leishmaniasis	17	1.19E-06
hsa05167	Kaposi sarcoma-associated herpesvirus infection	29	1.68E-06
hsa05146	Amoebiasis	19	3.14E-06
hsa04512	ECM-receptor interaction	17	5.38E-06
hsa04510	Focal adhesion	29	6.69E-06
hsa04933	AGE-RAGE signaling pathway in diabetic complications	18	2.26E-05
hsa04940	Type I diabetes mellitus	11	3.20E-05
hsa05165	Human papillomavirus infection	38	6.93E-05
hsa05163	Human cytomegalovirus infection	29	7.05E-05
hsa04974	Protein digestion and absorption	16	7.51E-05
hsa05134	Legionellosis	12	7.74E-05
hsa05164	Influenza A	24	7.76E-05
hsa05132	Salmonella infection	15	0.000157
hsa05166	Human T-cell leukemia virus 1 infection	27	0.000263
hsa04514	Cell adhesion molecules (CAMs)	20	0.00035
hsa05144	Malaria	10	0.00054
hsa04151	PI3K-Akt signaling pathway	37	0.000605
hsa05170	Human immunodeficiency virus 1 infection	25	0.000858
hsa04210	Apoptosis	18	0.001211
hsa04380	Osteoclast differentiation	17	0.00158
hsa04672	Intestinal immune network for IgA production	9	0.002223
hsa04670	Leukocyte transendothelial migration	15	0.002691
hsa05014	Amyotrophic lateral sclerosis (ALS)	9	0.00296
hsa04010	MAPK signaling pathway	30	0.002984
hsa04218	Cellular senescence	19	0.003215
hsa05130	Pathogenic Escherichia coli infection	9	0.005004
hsa04650	Natural killer cell mediated cytotoxicity	16	0.005015
hsa05160	Hepatitis C	18	0.005171
hsa00532	Glycosaminoglycan biosynthesis – chondroitin sulfate/dermatan sulfate	5	0.005557
hsa04142	Lysosome	15	0.006575