An In Silica Model for RPE Loss Patterns in Choroideremia

Abstract

Purpose: To use empirical data to develop a model of cell loss in choroideremia that predicts the known exponential rate of RPE loss and central, scalloped preservation pattern seen in this disease.

Methods: A computational model of RPE loss was created in Python 3.7, which constructed an array of RPE cells clusters, binarized as either live or atrophic. Two rules were applied to this model: the background effect gave each cell a chance of dying defined by a background function, and the neighbor effect increased the chance of RPE cell death if a neighbor were dead. The known anatomic distribution of rods, RPE, choriocapillaris density, amacrine, ganglion, and cone cells were derived from the literature and applied to this model. Atrophy growth rates were measured over arbitrary time units and fit to the known exponential decay model. The main outcome measures: included topography of atrophy over time and fit of simulated residual RPE area to exponential decay.

Results: A background effect alone can simulate exponential decay, but does not simulate the central island preservation seen in choroideremia. An additive neighbor effect alone does not simulate exponential decay. When the neighbor effect multiplies the background effect using the rod density function, our model follows an exponential decay, similar to previous observations. Also, our model predicts a residual island of RPE that resembles the topographic distribution of residual RPE seen in choroideremia.

Conclusions: The pattern of RPE loss in choroideremia can be predicted by applying simple rules. The RPE preservation pattern typically seen in choroideremia may be related to the underlying pattern of rod density. Further studies are needed to validate these findings.

Choroideremia is a progressive retinal degeneration with an X-linked inheritance pattern.¹ It is characterized by early outer retinal atrophy in the mid periphery and peripapillary regions.^1–3 Loss of the RPE in choroideremia follows a characteristic geometric progression pattern that progresses from the periphery toward the center, leaving a central, scalloped, shrinking island of residual RPE late in the stages of the disease.^1,4,5 This pattern of preservation of a central scalloped island is characteristic of this disease, but there is no clear consensus on why patients with choroideremia develop this anatomic pattern.⁶ The absence of an animal model makes it difficult to discern the underlying pathophysiology.^1,4,5 

In addition to the scalloped pattern of RPE preservation, previous groups, including our own, have shown that the area of the preserved central island decays exponentially over time.^7,8 The presence of an exponential relationship for the residual RPE area in choroideremia implies that the loss rate of remaining RPE cells is proportional to the number of remaining cells. However, such a relationship does not explain the geometric, scalloped pattern of the residual island in choroideremia. 

The pathophysiology of choroideremia has been linked to the mutation of REP-1 (Rab escort protein-1) on chromosome Xq21.2, whose function is integral to intercellular vesicular transport.^1,6 REP-1 is ubiquitously expressed in retinal layers, with broad expression within the RPE, choriocapillaris, rods, cones, and other retinal cell types.⁶ There is a lack of consistency in the literature as to which outer retinal layer is the initial site of disease. Some studies suggest injury begins in the choroidal vasculature,^9,10 whereas others suggest that the RPE is the primary site of injury,^11,12 and a histopathologic case report of a female X-linked carrier suggests injury starts in the rods before propagating to the RPE and choriocapillaris.¹³ 

Recently, we have shown that the growth of atrophic lesions in AMD is governed linearly by the number of RPE cells exposed along the border of atrophy.¹⁴ 

The objective of this research is to provide a unifying explanation to these unique features of choroideremia and infer governing principles which drive this behavior. To this end, we developed a computational model of RPE loss for choroideremia. We inferred two different rules to model RPE loss at each time point using this model. From the principle of exponential decay,^7,8 and peripheral distribution of lesions in choroideremia,¹⁰ we infer a “background effect,” where each RPE cell has a random chance of becoming atrophic. From the principle that macular atrophic lesions grow in proportion to their perimeter,¹⁴ we infer a neighbor effect; an RPE cell is more likely to atrophy when a neighboring RPE is dead. We show the results of various combinations of these rules in approximating known findings in choroideremia. 

Methods

The following model was conducted in a Python 3.7 environment, importing NumPy and Scikit-learn libraries for computational functions, and the MatPlotLib library for visualization.^15–17 All simulations were performed on a Dell Inspiron 14 7490 with 8 gb of RAM and using an Intel i7-10510u processor at 1.80 gHz processing speed. The RPE was modeled as a 300 × 300 dimensional array, wherein each array element represents a cluster of cells, and the array was made circular by excluding array elements farther than 150 units from the center (see Fig. 1). The size of this array was modeled as being equivalent to 55° of retinal eccentricity, centered on the fovea, as this was the field limit used in the comparative meta-analysis.⁸ As a result, each array element or “cell cluster” represents approximately 205 to 260 cells.¹⁸ Each array element was assigned one of two states; either a “1” indicating the cluster of cells was healthy, or “0” indicating the cell cluster was dead. This assignment was represented by either dark or light clusters, respectively, on a diagram to visualize this effect (see Fig. 1). In the initial condition, each cell cluster was assumed to contain only live cells, which gives the array a homogeneous dark grey color (Fig. 1, column 1). When a cell dies, the dark grey cluster turns to light grey (Fig. 1, column 2). Then, an algorithm was implemented to simulate rules that may govern cell death, namely, either a background effect, neighbor effect, or a combination of the two over a course of time steps. 

Figure 1.

 

Simulation of cell loss in an RPE monolayer over time. Each circle has a diameter of 55°, which is roughly the region captured on wide angle retinal photographs. Dark gray indicates healthy RPE and light gray indicates areas of RPE death or atrophy. Optic disc and vascular arcades are depicted for scale only. Top row: Background effect only, in which the probability of any RPE cell dying is the same. Middle row: Neighbor effect added to the background effect. In this model, there is background loss of RPE similar to the top row, but the presence of an adjacent area of atrophy will increase the risk of atrophy for an individual RPE. Bottom row: Background effect × Neighbor effect. In this model, the presence of an atrophic neighbor will multiply the size of the background effect. The background function used is rod density (see Fig. 2), which leaves peripheral RPE atrophy with a preserved central island.

Background Effect

For the background effect, each cell has a chance of dying as defined by a background function, which governed risk of death by distance from the foveal center. Computationally, the background effect function was written so each element in the RPE array was cycled through using for loops; if the element was within bounds and not already atrophy, a random float value between 0 and 1 was generated using the random() function. If this value was less than the background function at that position, then that array element's value was changed to 0. These background functions were then linearly scaled in each model for residual area to approach −0.5 log units at approximately 75 time steps, where each time step is to approximate 1 year, to best match with scaling as presented in the compared meta-analysis.⁸ This practice affects horizontal scaling of data, but not pattern or rate of area loss. 

Initially, we used a uniform background function where each array element had the same chance of atrophy. As an initial modification of this process, we repeated the analysis using an RPE density function from the literature, because it is known that the RPE density is greatest centrally and decreases in the periphery.¹⁸ We then hypothesized that the pattern that was seen in RPE cell loss might be affected by the target cell of injury in choroideremia. For example, if choroideremia starts in cone cells one might see central loss initially with peripheral preservation, versus a different a pattern if the initial site of injury for choroideremia is in either rods, RPE, choriocapillaris, or other neural cellular elements of the retina. We determined these background density functions, in addition to ganglion cell and amacrine cell density, by review of the published literature^18–21 searched using PubMed and Google Scholar. For amacrine and ganglion cell density, only a single article was identified, which presented cell density with respect to the fovea. When multiple literature sources were found for a cell type (rod and cone cell density), literature was chosen which used an in vivo technique (adaptive optics scanning laser ophthalmoscopy), in contrast with histologic techniques, to best approximate live cell density distributions. Background functions of these respective layer densities with respect to distance from the foveal center were extracted from this literature using Data Thief 3²² (Supplementary Material). If a cellular density was not measured as far out as 55° of retinal eccentricity, data representing the broadest range of values for cellular density and the data from the most peripheral 5° of eccentricity were extrapolated outward in a linear fashion as an approximation. Dimensionless units were used for all measurements, with all functions normalized from 0 to the maximum value of 1 to make it unitless, with all other values being proportional to that maximum. 

This function modified the relative chance to atrophy as distance from the foveal center.¹⁷ 

Neighbor Effect

The neighbor effect was defined such that a cell is more likely to undergo death if it is adjacent to a dead cell than if it is surrounded by live cells. In other words, the chance of cell death increased per neighbor atrophic cell cluster (including diagonally adjacent cell clusters in this abstraction). For this analysis, we modeled the neighbor effect under two distinct conditions: (1) a linear neighbor effect model, with a fixed increase in the chance of atrophy if a neighboring cell is atrophic, independent of the background effect (Background + Neighbor), that is:  

\begin{eqnarray*}Chance\;of\;Atrophy = Fixed\;Constant\;Probability\end{eqnarray*}

 

(2) We also simulated it as a modifying factor of the background effect (Background × Neighbor), in which a “coupling constant” multiplied the chance of atrophy from the background effect for each atrophic neighboring cell cluster, that is:  

\begin{eqnarray*} \begin{array}{@{}l@{}}Chance\;of\;Atrophy = Background\;Function\;*\;Coupling\;Constant\end{array}\end{eqnarray*}

 

In other words, the coupling constant is the factor by which a neighboring dead RPE cell increases the chance of cell death compared with the background effect; in the case of using rod cells as a background function, this was scaled by comparing the peak RPE cell density in the macula with the peak rod cell density.^18,20 This coupling constant was varied for each background function by running the simulation ten times at each coupling constant value in increments of ten to optimize the exponential fit of the residual area (discussed elsewhere in this article). 

Residual Area Calculation

For each simulation, the residual area was calculated at each time point in a stepwise fashion. Simulations were repeated in an iterative fashion until the residual area was reduced below an equivalent size of 316 microns; this value was used because 316 microns, or −0.5 log mm², was the minimum area observed in the compared meta-analysis.⁸ Because prior studies have shown that the residual area in choroideremia follows exponential decay,^7,8 we plotted the residual area on a logarithmic plot (Fig. 2) in the form:  

\begin{eqnarray*}log\left( {Area} \right) = mX + b\end{eqnarray*}

 

Figure 2.

 

(A) Residual RPE area as a function of time for all three models on a logarithmic scale. Top row: Background effect only, in which the probability of any RPE cell dying is the same. Middle row: Neighbor effect added to the background effect. In this model, there is background loss of RPE like the top row, but the presence of an adjacent area of atrophy will add to the risk of atrophy for an individual RPE. Bottom row: Background effect × Neighbor effect. In this model, the presence of an atrophic neighbor will multiply the size of the background effect. Both the background effect (top row) and the Background × Neighbor effect (bottom row) demonstrate exponential decay. The bottom row shows an r² = 0.992. (B) Empirically derived data to demonstrate the tight exponential distribution over 68 individual eyes analyzed in a previously published meta-analysis.⁸

We then performed linear regression on these area calculations to determine whether the data indeed fit followed exponential decay on a logarithmic plot using the Sci-kit-learn linear regression function (“sklearn.linear_model.LinearRegression()”) to calculate the linear regression and an R² value.¹⁶ 

No human subjects were used in this study, and this methodology followed the Declaration of Helsinki. 

Results

The implementation of a background effect alone accurately simulated exponential decay, but did not lead to a residual island of RPE that resembled the central scalloping of choroideremia; rather, the pattern resembles a random “raindrop” pattern (Fig. 1, row 1). When the neighbor effect was applied in addition to the background effect, a residual island of RPE centrally was not observed, which suggest that this could not explain the anatomy seen in choroideremia (Fig. 1, row 2). When the neighbor effect was applied multiplicatively with the background effect, a central scalloped area of residual RPE was observed if the rod density was used as the background function (Fig. 1, row 3). We then determined the residual RPE area as a function of time, to see which model led to the best expected fit for exponential decay on a logarithmic plot. In this model, the best fit was observed for both the background effect alone (Fig. 2A, row 1) and Background × Neighbor effect (Fig. 2A, row 3). The fit to exponential decline was poor when we assumed that background effect had a superimposed linear neighbor effect (Fig. 2A, row 2). These models used the rod density function as the background function; other background functions were also explored, as described elsewhere in this article. 

We then determined the effects of background cell density on the residual area of RPE in choroideremia, using the Background × Neighbor effect model. The density of various cell types (RPE, amacrine cells, ganglion cells, cones, and 1/cones determined from a review of the literature) was used in the simulation. When we used the rod density function, we were able to generate a residual island of RPE centrally with progressive loss of peripheral RPE, like what is seen in choroideremia. Note the scalloped RPE appearance at the latest time point (Fig. 3, row 1, right). Using the choriocapillaris and RPE density function background shows essentially random, growing areas of RPE atrophy that do not result in a preserved central island (Fig. 3, rows 2 and 3). Applying the amacrine cell density function background shows progressive loss of the RPE layer with a residual circumferential zone of RPE from the equator peripherally, and a very small central island (Fig. 3, row 4). Using the ganglion cell density function background shows progressive loss of the RPE layer, with small patches of residual RPE centrally but larger patches from about 30° to 55° (Fig. 3, row 5), in a pattern that is not characteristic of RPE sparing in choroideremia. Neither the cone density function (Fig. 3, row 6) nor the 1/cone density function (Fig. 3, row 7) showed a pattern of RPE sparing in choroideremia. 

Figure 3.

 

Effects of background cell density on the residual area of RPE in choroideremia, using the Background × Neighbor effect model. The density of various cell types was determined from a review of the literature, and this cell density, unique for each cell population, was used in the simulation. Columns represent four different time points as the model for cell loss evolves. Top row: Rod density function background shows progressive loss of the RPE layer with a residual scalloped island of RPE centrally similar to what is seen in choroideremia. Note that scalloped RPE appearance at the latest time point (top row, right column). Second row: Choriocapillaris density function background shows progressive loss of the RPE layer small patches of residual RPE centrally (last column) but larger patches from about 30° to 55° (second row, right column) at the latest time point. This pattern is not characteristic of RPE sparing in choroideremia. Third row: RPE density function background shows progressive loss of the RPE layer with residual patches peripherally, but not centrally (third row, right column) at the latest time point. Fourth row: Amacrine density function background shows progressive loss of the RPE layer with a residual circumferential zone of RPE from the equator peripherally (fourth row, right column) at the latest time point. Fifth row: Ganglion cell density function background shows progressive loss of the RPE layer small patches of residual RPE centrally (last column) but larger patches from about 30° to 55° (second row, right column) at the latest time point. This pattern is not characteristic of RPE sparing in choroideremia. Sixth row: Cone cell density function background shows progressive loss of the RPE layer small, round area of residual RPE centrally (last column) at the latest time point. This pattern is not characteristic of RPE sparing in choroideremia. Bottom row: Inverse of cone cells: density function background shows progressive loss of the RPE layer small patches of residual RPE centrally (last column), but larger patches from about 30° to 55° (second row, right column) at the latest time point. This pattern is not characteristic of RPE sparing in choroideremia.

We then determined the residual RPE area as a function of time, to see which model led to the best expected logarithmic fit. For this, the background function was modified based on known cell topographies (Fig. 4).^18–20 Topologically, a central scalloped island was only found when we used data on known rod densities as a background function. The cone background function was inverted to test the hypothesis that cones may be the cause of a preserved central scalloped island, which was not found. The amacrine density also left a central island. 

Figure 4.

 

Effects of background cell density on the residual area of RPE in choroideremia, using the Background × Neighbor effect model. The logarithm of the residual area was plotted versus time assuming different background cell density function. The density of various cell types was determined from a review of the literature, and this cell density, unique for each cell population, was used in the simulation.

Figure 5.

 

Model fit as coupling constant is changed for choroideremia, using the rod density background function. Coupling constant is defined as the factor by which the neighbor effect multiplicatively increases the background effect, after normalizing for rod density and RPE density.^18,20 Gray shading indicates 95% confidence intervals. The coupling constant plateaus beyond values of approximately 5000.

These background functions were also tested for whether they demonstrate an exponential decay in the residual area (Fig. 4). Rod and amacrine cell functions approximated an exponential function, with R² of 0.992 and 0.986, respectively. 

We then determined the ideal coupling constant for choroideremia, using the rod density background function, wherein the coupling constant is defined as the factor by which the neighbor effect enhances the background effect (Fig. 5). Gray shading indicates 95% confidence intervals. This analysis shows that the R² value plateaus when the coupling constant is approximately 5000. 

Discussion

In this study, we developed an in silica computational model for cell loss in choroideremia that accurately models two known attributes of the residual central RPE island, namely, the exponential loss of residual RPE area as a function of time, and the late presence of a residual central scalloped island of RPE. This matches the empiric behavior of choroideremia as reported in a meta-analysis previously published by our group examining RPE atrophy expansion rates in 68 individual eyes (Fig. 2B).⁸ 

In at least two causes of macular atrophy, AMD and Stargardt disease, our group and others have found that the growth rate of geographic atrophy is proportional to the perimeter of the atrophic lesion, implying it is at the edge of geographic atrophy where injury occurs.^14,23 In essence, the lesions expand as the RPE cells “at risk” along the lesion border undergo cell death, and these cells are lost at a higher rate than RPE elsewhere in the monolayer that have healthy neighbors. We hypothesized that this principle could be generalized to disease with a well-defined peripheral chorioretinal atrophy, such as choroideremia. To test this hypothesis, we developed a model of cell loss in which there is random loss of RPE anywhere within the monolayer (background effect), in addition to loss of RPE at a higher rate if they were adjacent to dead RPE (neighbor effect). Importantly, neither the background effect nor the neighbor effect alone was enough to simulate the anatomy observed in the late stages of choroideremia. However, when the two effects were combined multiplicatively, we demonstrated that a residual island would be present late in the stages of the disease. This model demonstrates that the residual island of central RPE has an area that decreases exponentially with time, similar to what has been seen in clinical studies. We then used the cell density of different background functions for neural cells, RPE, and choriocapillaris from the literature and used these density function to determine which cell type might be involved early in choroideremia. Interestingly, a central scalloped island was only seen when we used the rod cell density as a background function. The use of the rod density function results in the preservation of this central island because the density of rod cells decreases as one approaches the foveal center, thus decreasing the proportional rate of RPE loss. Eventually, the central island is lost, as seen in real-world examples.⁴ 

Based on this result, we propose that the underlying rod cell density may generate the residual central RPE island seen in choroideremia, that is, the density distribution of rods correlates with the pattern of choroideremia. This finding suggests that rod cells may be the primary cells involved in choroideremia, but other explanations also exist. For example, an alternative hypothesis is a failure of transport of an important gene product from the RPE or choriocapillaris, and that the rods are the most susceptible to their absence. Interestingly, some histopathology has suggested choroideremia is a rod-based disease. For example, a case examining a female carrier of the CHM mutation showed outer retinal atrophy only where there was significant rod atrophy.¹¹ 

This model provides a mechanism that can unify the mechanism of cell loss seen in choroideremia and very different diseases, such as geographic atrophy in AMD; namely, that cell death is more common in cells along a border that neighbor dead cells rather than healthy cells. In cell biology, it is known that cell loss can have a profound influence on the behavior and fate of adjacent cells. For example, cell–cell interactions are crucial to maintaining the integrity of a cellular monolayer, including the corneal endothelium, RPE, and epithelia and endothelial layers.^24–26 Another model to explain the tendency of RPE atrophy to occur along a border with dead cells is that RPE cell loss may be repaired if the neighboring RPE cells are healthy, where atrophy may propagate when the ability of neighboring cells to migrate and expand fails.^18,26,27 We propose this computational method may be useful in modeling other retinal degenerations, whether they be peripheral or macular diseases. Additional study is required to fully understand the pathogenesis of this important disease. 

There are several limitations to this study. For computational purposes, many features were abstracted in this model; individual RPE cells were not modeled, but rather clusters of RPE that were modeled to behave en bloc. The behavior of the choriocapillaris and photoreceptors were not modeled directly, except by applying their cell or vessel densities to the model. Further, the state of these clusters was simplified to a binary state of live or atrophic; variations along the spectrum toward cell death were not simulated. The presence of the optic nerve was not simulated. We only modeled cell density as being responsible for topographic variations in rates of cell death. However, empiric data for other features that may vary over topography were not available. Further, this model only shows that the mechanism that drives RPE atrophy in choroideremia may follow the form of the rod density function; this function could have a cause unrelated to rod density or may be modulated by a complex combination of functions which happens to approximate the rod density function. It is also possible entirely other hypotheses may be able to explain both the rate of loss and topography of choroideremia using more complex mechanisms. However, we propose that the model presented is the simplest possible explanation for the strongly consistent patterns of choroideremia.^7,8 Last, other features to verify this model, such as actual topography of choroideremia lesions over time, were not available. Although this model is driven by data from a meta-analysis, and qualitative topography noted in choroideremia, it did not otherwise have in vivo or histologic data to test predictions made by this model. Examples of such predictions include lesions in the midperiphery being larger than average than those in the far periphery or macula, and the possibility of being able to predict longitudinal changes in real-world atrophic lesions, sufficient clinical longitudinal images of choroideremia were not available to test this prediction. Testing these hypotheses would be worthy of further study. 

Acknowledgments

Supported by P30 EY026878 from the National Eye Institute (NEI) (Recipient: Yale Vision Science Core). 

Disclosure: B.K. Young, None; L.L. Shen, Consultant - Boehringer Ingelheim; L.V. Del Priore, Consultant - Boehringer Ingelheim. Scientific advisory board - Tissue Regeneration Sciences; LambdaVision, CavtheRx; Scientific advisor - Seeing Medicines; Research through university - Stealth Pharmaceuticals