Single-cell transcriptomic analysis revealing changes in retinal cell subpopulation levels and the pathways involved in diabetic retinopathy.

Background: Diabetic retinopathy (DR) is a common microvascular complication of diabetes and one of the most common causes of visual impairment and blindness. However, it is not yet known how abnormal retinal cell subpopulations contribute to disease progression.
Methods: In this study, we used the Gene Expression Omnibus database to construct a single-cell atlas of DR and healthy samples to explore changes in the abundance of different cell subpopulations, and identify the molecular pathways potentially involved in DR.
Results: Our results showed that DR was associated with significantly reduced numbers of bipolar cells, Müller glia, retinal pigment epithelial cells, and cone photoreceptors, but was also associated with significantly greater numbers of pericytes, rod photoreceptors, anaplastic cells, and microglia. Our results suggest that subpopulations of Müller glia, microglia, endothelial cells, and bipolar cells in DR tissues may be involved in various oxidative stress- and inflammation-related pathways. Conclusions: In summary, we showed that Müller glia, endothelial cells, microglia, and bipolar cells in DR tissues are involved in oxidative stress- and inflammation-related pathways, which may contribute to the progression of the disease and ultimately lead to visual impairment and blindness. This study will provide a theoretical basis for further exploring the specific mechanism of DR. 2022 Annals of Translational Medicine. All rights reserved.

Introduction

Diabetes mellitus refers to a group of metabolic diseases involving chronic hyperglycemia, resulting from multiple causes, and its incidence has been increasing dramatically worldwide (1). Diabetic retinopathy (DR) is a common microvascular complication of diabetes that affects 1/3 of diabetic patients and is associated with an increased risk of life-threatening systemic vascular complications, such as stroke, coronary artery disease, and heart failure (2). DR involves specific neuroglial and microvascular abnormalities, such as microaneurysms, soft or hard exudates, intraretinal microhemorrhages, beaded veins, and intraretinal microvascular abnormalities (3,4), which have been shown to progressively disrupt retinal function (5). Traditional DR treatments include laser therapy, vitrectomy, and conventional medications to control hyperglycemia, hypertension, and hyperlipidemia. These treatments inhibit retinal neovascularization, reduce retinal macular edema, and prevent the stretching of fibrous tissue during the proliferation of new vessels in the retina (6). However, laser therapy is often accompanied by vitreous hemorrhage, while 23G minimally invasive vitrectomy can lead to retinal rupture and postoperative cataract formation (7,8). As these treatments can have a number of adverse effects, such as intraoperative retinal fissures and postoperative retinal detachment (9), a more in-depth exploration of the disease is needed to develop safer treatment plans.

The pathogenesis of DR involves multiple biochemical signaling pathways, such as protein glycosylation and the increased activity of protein kinase C, which lead to the formation of advanced glycosylation end products (AGEs). These can promote intercellular interactions involving the vascular endothelial growth factor, leading to neointima formation in the anterior and posterior segments of the eye, increased vascular permeability, resulting in leakage, and the collapse of the internal blood-retinal barrier (10). DR progression appears to involve AGE accumulation, an increase in glucose flux into the polyol and hexosamine pathways, and the excessive activation of the plasma kininase-releasing kinase pathway (11).

DR also involves inflammation (12,13), which damages and permeabilizes the retinal vascular system and disrupts neointima formation, 2 of the leading causes of visual impairment in diabetes (14,15). The dysfunction of retinal glial cells, such as astrocytes, Müller glia, and microglia, may contribute to early retinal inflammatory processes in DR (16). Cellular levels change during disease progression. Thus, this study sought to explore how abnormalities at the level of subpopulations of retinal tissue cells play a role in the development of DR.

Single-cell RNA sequencing (scRNA-seq) technology can reveal unique changes in each cell and enable precise classification of different cell types, as well as help us understand the changes in tumor cells during metastasis to target organs (17). It provides new and precise insights into our understanding of genetic heterogeneity during tumor evolution and tumor cell metastasis. In this study, we used single-cell transcriptome data from DR and healthy tissues to construct a single-cell atlas of DR and explore changes in the levels of different cell subpopulations in the disease. A functional enrichment analysis was also conducted to identify the signaling pathways that may be involved in the disease. In contrast to others’ studies, we will also explore intercellular interactions to further analyze their role in the development of DR. Our findings may help clarify the disease mechanism and inform the development of effective treatment strategies. We present the following article in accordance with the MDAR reporting checklist (available at https://atm.amegroups.com/article/view/10.21037/atm-22-1546/rc).

Methods

Data source

Single-cell transcriptomic data were obtained from the Gene Expression Omnibus database (18). The GSE178121 data set was based on the GPL24247 platform and included retinal tissue samples from healthy C57BL/6 mice (n=1) and C57BL/6 mice with streptozotocin-induced diabetes (n=1) (a total of 2 samples). In GSE178121, the mice kept high blood glucose concentration (>300 mg/dL) for 25 weeks before sacrificed to prepare retinal cell suspension for scRNA-Seq, all cells were coming from mice, thus there no cell lines.

Construction of the single-cell atlas of DR and cell-type annotation

Cell clustering was performed using the Seurat package in R with default parameters (19). The clustering results were uniformly downscaled and visualized using a principal component analysis and t-distributed stochastic neighbor embedding (20), and then projected in a 2-dimensional image defined as the single-cell atlas. The FindAllMarkers function of the Seurat package was also used to identify the marker genes that were highly expressed in each cell cluster. Cell-type annotation was performed based on known retinal cell markers (18) and cell markers of the Human Proteome Project (21).

Single-cell component analysis

A single-cell component analysis was performed based on the markers of different cell clusters to investigate differences in the numbers of retinal cell subpopulations between the DR and healthy samples and to identify cells dysregulated in DR.

Analysis of protein-protein interaction network and receiver operating characteristic (ROC) curve

Moreover, the differentially expressed genes with adjusted P<0.05 was conducted to constructed the protein-protein interaction (PPI) network based on STRING database (22). Then, the interaction between cells was plot with iTalk package (23).

To observe the diagnostic potential of differentially expressed genes (adjusted P<0.05, average log2 fold change >1), area under the receiver operating characteristic (ROC) curve was applied to the analysis using pROC package (24). While the closer the area under the curve (AUC) is to zero, the higher the diagnostic potential.

Cell subcluster analysis

The subcluster analysis of the identified cells was performed with the Seurat package, and the FindAllMarkers function was used to identify marker genes expressed in each cell subcluster. Cell subclusters were then classified based on which marker genes were expressed the most abundantly.

Functional enrichment analysis

A Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was conducted using the clusterProfiler package (25) to identify the potential functions of molecular pathways occurring in each cell subpopulation. Pathways were considered significantly involved in marker genes when they had an adjusted P value <0.05. The adjusted P value was obtained using the Benjamin-Hochberg method (26), which can correct the P value and control the false/true positive ratio within a certain range.

Data and statistical analyses

All statistical analyses were performed using R (https://www.r-project.org/). We analyzed the expression levels of genes using unpaired t-tests. P<0.05 was consistent with statistical significance. And the analyses were performed using the online Platform Bioinforcloud platform (http://www.bioinforcloud.org.cn) and Biomarker Technologies Corporation (Beijing, China).

Results

Single-cell atlas of DR

A single-cell atlas comprising 55 cell clusters () was constructed using a total of 21,588 cells, of which 13,336 came from the DR mice and 8,252 came from the healthy controls (). The obtained clusters were classified into the following 10 cell types: retinal pigment epithelial cells, pericytes, endothelial cells, astrocytes, microglia, anaplastic cells, cone photoreceptors, rod photoreceptors, bipolar cells, and Müller cells (). The following marker genes were the most strongly expressed in the various cell clusters: cone photoreceptors (Cone), Guca1; Müller cells, Rlbp1; bipolar cells (BCs), Trpm1; rod photoreceptors (Rod), Nrl; pericytes (Per), Sparcl1; endothelial cells (ECs), Pecam1; anaplastic cells (ACs), Pax6; retinal pigment epithelial cells (RPEs), Rpe65; microglia (Mic), Aldoc; astrocytes (Ast), and Gfap (). We explored the differences in the numbers of various cell types between the DR and healthy control samples (). The results indicated that Müller glial cells were significantly decreased in the DR samples. The analysis revealed a wide range of intercellular interactions, suggesting that DR development is closely related to intercellular interactions (). In the present study, Cirbp was found to have potential (AUC =0.795) as a diagnostic marker for DR ().

Figure 1

Single-cell mapping of DR. (A) Single-cell mapping of 21,588 cells from DR and healthy control samples. Each color indicates a different cell cluster (n=55). (B) Single-cell atlas based on the sample source. (C) Classification of cell clusters into cell types. Each color represents a different cell type. (D) Expression of the marker genes in each cell cluster. (E) Levels of various cellular components in DR and healthy control samples. (F) The interaction between cells. (G) The ROC of marker gene. ACs, anaplastic cells; Ast, astrocytes; BCs, bipolar cells; Cone, cone photoreceptors; DR, diabetic retinopathy; ECs, endothelial cells; Mic, microglia; Müller, Müller cells; Per, pericytes; Rod, rod photoreceptors; RPEs, retinal pigment epithelial cells; tSNE, t-distributed stochastic neighbor embedding; ROC, receiver operating characteristic; AUC, area under the curve.

Müller glial cell function in DR

The Müller cells were divided into 78 subpopulations (), and a single-cell atlas was constructed according to the sample source (). The subpopulations Müller_Mt2 and Müller_Mir124-2hg were significantly larger in the DR tissues than the control tissues, while the other 5 subpopulations were significantly smaller (). The expression of the marker genes in each Müller cell subpopulation was determined by a cell subcluster analysis (). The KEGG pathway enrichment analysis showed that the subpopulation genes were involved in several inflammation-related pathways, such as synaptic vesicle cycling, interleukin 17 (IL-17) signaling, Toll-like receptor signaling, and the lysosomal pathways ().

Figure 2

Müller glial cell function in DR. (A) Single-cell mapping of Müller cell subpopulations. Each color indicates a different subpopulation. (B) Single-cell atlas of Müller cells based on the sample source. (C) Differences in the abundance of the Müller cell subpopulations between the DR and healthy control samples. (D) Expression of the marker genes in each subpopulation. (E) Signaling pathways enriched in Müller cell subpopulations. DR, diabetic retinopathy; KEGG, Kyoto Encyclopedia of Genes and Genomes; tSNE, t-distributed stochastic neighbor embedding.

Endothelial cell function in DR

The endothelial cells were divided into 4 subpopulations () and a single-cell atlas was constructed based on the sample source (). The levels of ECs_Mfsd2a were significantly lower in the DR group than in the healthy control group, while the other subpopulations were significantly increased (). The expression of the marker genes in each endothelial cell subpopulation was determined by a cell subcluster analysis (). The enrichment analysis showed that the subpopulation genes were involved in the cell adhesion-related signaling pathways, such as extracellular matrix-receptor interactions and gap junctions, and the chemokine signaling pathway (chemokine is known to play an important role in inflammation) (27) ().

Figure 3

Endothelial cell function in DR. (A) Single-cell mapping of endothelial cell subpopulations. Each color indicates a different subpopulation. (B) Single-cell atlas of endothelial cells based on the sample source. (C) Differences in the abundance of the endothelial cell subpopulations between the DR and healthy control samples. (D) Expression of the marker genes in each subpopulation. (E) Signaling pathways enriched in endothelial cell subpopulations. DR, diabetic retinopathy; KEGG, Kyoto Encyclopedia of Genes and Genomes; tSNE, t-distributed stochastic neighbor embedding.

Microglial function in DR

The microglia were divided into 3 microglia subpopulations () and were used to construct a single-cell atlas depending on the sample source (). The abundance of Mic_Fcer1g and Mic_Slco2b1 was lower in the DR group than the healthy control group, while that of the Mic_Bex2 was significantly higher in the DR group than the healthy control group (). The expression of the marker genes in each microglia subpopulation was determined by a cell subcluster analysis (). Consistent with previous findings (28), the subpopulation genes were found to be involved in several inflammation-related pathways, such as synaptic vesicle cycling, IL17 signaling, and Toll-like receptor signaling ().

Figure 4

Microglial function in DR. (A) Single-cell mapping of microglia subpopulations. Each color indicates a different subpopulation. (B) Single-cell atlas of microglia based on the sample source. (C) Differences in the abundance of the microglial subpopulations between the DR and healthy control samples. (D) Expression of the marker genes in each subpopulation. (E) Signaling pathways enriched in microglial subpopulations. DR, diabetic retinopathy; KEGG, Kyoto Encyclopedia of Genes and Genomes; tSNE, t-distributed stochastic neighbor embedding.

Bipolar cell function in DR

The bipolar cells were divided into 6 subpopulations () and a single-cell atlas was constructed based on the sample source (). Our results showed that only 2 of the 6 subpopulations were larger in the DR tissues than the control tissues (). The expression of the marker genes in each bipolar cell subpopulation was determined by a cell subcluster analysis (). The enrichment analysis showed that the subpopulation genes were involved in inflammation-related pathways, including chemokine signaling, advanced glycosylation end product (AGE)-receptor for advanced glycation end product (RAGE) signaling in diabetic complications, and calcium signaling ().

Figure 5

Bipolar cell function in DR. (A) Single-cell mapping of bipolar cell subpopulations. Each color indicates a different subpopulation. (B) Single-cell atlas of bipolar cells based on the sample source. (C) Differences in the abundance of the bipolar cell subpopulations between the DR and healthy control samples. (D) Expression of the marker genes in each subpopulation. (E) Signaling pathways enriched in bipolar cell subpopulations. DR, diabetic retinopathy; KEGG, Kyoto Encyclopedia of Genes and Genomes; tSNE, t-distributed stochastic neighbor embedding.

Discussion

In this study, we used DR-related single-cell transcriptome data to construct a single-cell atlas of the DR and healthy control samples. We also performed a cell subpopulation analysis to explore changes in the abundance of different cell subpopulations in DR and identify pathways involved in disease progression. Our results showed that DR was associated with significantly reduced numbers of bipolar cells, Müller glia, retinal pigment epithelial cells, and cone photoreceptors, but was also associated with significantly greater numbers of pericytes, rod photoreceptors, anaplastic cells, and microglia. The enrichment analysis suggests that the retinal cell subpopulations altered in DR are involved in pathways related to oxidative stress and inflammation.

Müller glial cells are the predominant cell type in the mammalian retina and play an important role in maintaining homeostasis and retinal integrity (29,30), as they can alter their morphology to respond to retinal damage and disease (31). However, retinal damage can reduce their abundance, thereby limiting the healing response and promoting the development of DR (32). Müller cells have also been found to produce high levels of pro-inflammatory molecules under hyperglycemic conditions (33). Our study showed that Müller cell subpopulation genes are involved in lysosomal pathways, and lysosomal membrane permeabilization, and the leakage of lysosomal contents are known to lead to cell death, which may explain their decreased abundance in DR. The genes are also involved in several inflammation-related pathways. These findings imply that restoring the response of Müller cells in DR lesions and lowering blood glucose levels may inhibit the secretion of pro-inflammatory factors and mitigate DR.

Endothelial cells are one of the main components of the retinal microvasculature, which is essential for retinal function (34). They play an important role in both inflammation and angiogenesis in DR (15) and are critical in attracting immune cells to sites of inflammation (35). Inflammatory cells produce angiogenic cytokines and growth factors, which in turn activate the microvascular endothelium, which expresses pro-inflammatory molecules (36-38). In our study, we showed that endothelial cell subpopulations are involved in similar inflammation-related pathways.

Microglia are a key component of the retina and are known to enhance neurodegeneration in ocular disease by producing pro-inflammatory neurotoxic cytokines and phagocytosing neurons (28). The activation of microglia may help initiate neuroinflammation in DR (14). Activated microglia produce pro-inflammatory mediators that promote neuronal cell death but that may also act as neuroprotective agents (39-41). In the present study, the abundance of microglia were increased, which might drive DR progression, and we found that microglia subpopulations in DR participate in inflammation-related pathways.

In the retina, bipolar cells receive input from rod and cone photoreceptors, ganglion cells receive input from bipolar cells, and amacrine cells regulate these interactions (42). Bipolar cells are a fundamental component of vision, as they serve as a link between the 2 synaptic layers of the retina (43,44). An earlier study using the GSE178121 data set showed that the bipolar cell population is not significantly affected in DR (18); however, we found that the abundance of bipolar cells was significantly lower in DR tissues than healthy tissues. We also found that the BCs_Gngt1 subpopulation participates in the calcium signaling pathway. It may be that the abnormal levels of bipolar cells in DR are the result of dysfunctional calcium homeostasis, as decreased calcium signaling inhibits neurotransmission (45), particularly in patients with diabetes (46-48).

Inflammation plays a key role in the development of DR (12,13,49). In this study, we found that Müller_Slc2a1, Müller_Nupr1, and Mic_Slco2b1 subpopulation genes are involved in the Toll-like receptor signaling pathway, which promotes the secretion of angiogenic growth factors in arterial endothelial cells and leads to choroidal neovascularization (50,51). We also found that the ECs_Cirbp subpopulation is involved in the chemokine and nucleotide-binding oligomerization domain (NOD)-like signaling pathways. Chemokines are a large class of chemotactic cytokines that act on G protein-coupled receptors and regulate various biological processes, and they can help drive diseases involving pathological inflammation (52,53). Conversely, NOD-like receptor (NLR) proteins are pattern recognition receptors that mediate the innate immune response to cellular injury and stress (54). NLRs can also form inflammatory vesicles, which mediate the activation of caspase-1 and the secretion of the pro-inflammatory cytokines IL-1β and IL-18 (55,56). Endothelial cells contain abundant danger signals, such as reactive oxygen species, which activate the NLRP3 inflammasome, exacerbating inflammatory signaling, and cellular damage (57). These results suggest that endothelial cell abnormalities in DR may be associated with the abnormal activation of these pathways.

The results of the enrichment analysis of the above subgroups revealed that the “AGE-RAGE signaling pathway in diabetic complications” is significantly enriched in the Müller_Slc2a1, Müller_Vegfa, BCs_Sebox, BCs_Ndnf, and BCs_Cnmd subpopulations. AGEs in diabetes are formed by the non-enzymatic reaction between glucose and long-lived proteins, such as collagen (58), and the AGE-RAGE axis is strongly associated with the progression of diabetes-related macrovascular complications, arterial damage, diabetic nephropathy, and retinopathy (59). AGEs have been found to induce inflammation and immunosuppression by binding to RAGE, leading to the production of pro-inflammatory cytokines, reactive oxygen species, and reactive nitrogen intermediates that alter innate and adaptive immune responses (60). These molecules also promote the pathogenesis of diabetic vascular complications mediated by oxidative stress (61-63), and activate nicotinamide adenine dinucleotide phosphate oxidase and the nuclear factor kappa light chain enhancer of activated B cells by binding to the specific RAGE on the cell surface (64), thus triggering a vicious cycle of oxidative stress and inflammation (62,65,66). Based on these results, we hypothesize that the AGE-RAGE signaling pathway in diabetic complications may be associated with abnormalities in bipolar cells and Müller glia. Further, we found that Müller glia, endothelial cell, microglia, and bipolar cell subpopulation genes are significantly involved in the phosphatidylinositol-3-kinase and protein kinase B signaling pathway, which accelerates the endothelial-mesenchymal transition in DR (67).

This study was done based on bioinformatics analysis and single-cell sequencing of cells, and the lack of clinical experimental validation is one of the drawbacks of this study. Since human single cell data on DR are scarce, this study used data from mice, whose samples are also small, and in future plans will expand the samples and collect data from humans as much as possible.

In summary, we showed that Müller glia, endothelial cells, microglia, and bipolar cells in DR tissues are involved in oxidative stress- and inflammation-related pathways, which may contribute to the progression of the disease and ultimately lead to visual impairment and blindness. If our findings can be confirmed, they suggest that the early detection of changes in specific retinal cell subpopulations may lead to earlier and more effective treatments.

The article’s supplementary files as