Visualization of Optic Nerve Structural Patterns in Papilledema Using Deep Learning Variational Autoencoders

Abstract

Purpose: To visualize and quantify structural patterns of optic nerve edema encountered in papilledema during treatment.

Methods: A novel bi-channel deep-learning variational autoencoder (biVAE) model was trained using 1498 optical coherence tomography (OCT) scans of 125 subjects over time from the Idiopathic Intracranial Hypertension Treatment Trial (IIHTT) and 791 OCT scans of 96 control subjects from the University of Iowa. An independent test dataset of 70 eyes from 70 papilledema subjects was used to evaluate the ability of the biVAE model to quantify and reconstruct the papilledema spatial patterns from input OCT scans using only two variables.

Results: The montage color maps of the retinal nerve fiber layer (RNFL) and total retinal thickness (TRT) produced by the biVAE model provided an organized visualization of the variety of morphological patterns of optic disc edema (including differing patterns at similar thickness levels). Treatment effects of acetazolamide versus placebo in the IIHTT were also demonstrated in the latent space. In image reconstruction, the mean signed peripapillary retinal nerve fiber layer thickness (pRNFLT) difference ± SD was −0.12 ± 17.34 µm, the absolute pRNFLT difference was 13.68 ± 10.65 µm, and the RNFL structural similarity index reached 0.91 ± 0.05.

Conclusions: A wide array of structural patterns of papilledema, integrating the magnitude of disc edema with underlying disc and retinal morphology, can be quantified by just two latent variables.

Translational Relevance: A biVAE model encodes structural patterns, as well as the correlation between channels, and may be applied to visualize individuals or populations with papilledema throughout treatment.

Introduction

Characterizing the swollen optic nerve head (ONH) is crucial for diagnosing and monitoring papilledema. Although Frisén grading is helpful,^1–3 its subjectivity limits interobserver repeatability, especially with subtle edema improvement within a single grade.^3–6 Commercial optical coherence tomography (OCT) devices can measure peripapillary retinal nerve fiber layer thickness (pRNFLT),^7–9 but this measurement is unstable with severe edema.^10,11 Our previously published three-dimensional (3D) graph-based approach enhances retinal segmentation stability, providing detailed maps of RNFL and total retinal thickness (TRT) maps that contain diverse spatial information about the swollen nerve.^12–17 This study aimed to quantify papilledema spatial patterns, visualize pattern changes over time, and relate the global summary measurements (i.e., pRNFLT and ONH volume [ONHV]) to the swelling spatial patterns. 

In the Idiopathic Intracranial Hypertension Treatment Trial (IIHTT) OCT substudy, we investigated OCT-related measures in relation to clinical parameters in IIH before and after treatment.^10,17–23 We found pRNFLT, peripapillary TR thickness (pTRT), and ONHV were correlated with Frisén grades at baseline but were not significantly correlated post-treatment.^5,20 Subjects in the treatment group showed significant improvement in the anterior displacement of Bruch's membrane over 6 months. The subjects in the placebo group did not show the same change despite reduced pRNFL thickness, implying that the translaminar pressure gradient across the optic nerve from elevated intracranial pressure (ICP) was still present in the placebo group,¹⁷ which suggests that pRNFLT and ONHV may not fully capture the important morphologic changes associated with papilledema. 

Deep-learning approaches^24–26 have rapidly developed and been increasingly applied to solve problems in neuro-ophthalmology,²⁷ including detection and quantification of papilledema severity from fundus photographs.^28–31 Variational autoencoders (VAEs)^32–34 are a well-known generative deep-learning architecture used to synthesize new data associated with the probability distribution of the training data and to enable visualization of smooth morphing between any two data points. A typical VAE concatenates an encoder and a decoder: The encoder is trained to deconstruct the input (e.g., images) into a succinct numeric representation, referred to as “latent variables,” and the decoder is simultaneously trained to reconstruct the input only based on these latent variables. VAEs are unsupervised, which means no labels (or ground truths) are needed during training. We previously proposed a basic VAE model to create a two-dimensional (2D) synthetic montage map (“latent space”) illustrating meaningful spatial pattern continuums of the retinal ganglion cell layer defects in glaucoma.³⁵ Use of a VAE latent space has also been applied to estimate the glaucomatous progression of the visual field and retinal thickness.^36–38 We believe the use of VAEs has strong potential to improve clinical interpretation based on the latent space. 

In this study, we proposed a novel bi-channel VAE (biVAE) model that facilitates visualization and quantification of spatial patterns of neuroretinal thickening observed in papilledema. We linked the RNFLT and TRT maps as an input pair and created corresponding RNFLT and TRT latent space montage maps to display and stratify the IIHTT data points by swelling patterns over time and overlaying global measurements and Frisén grades. By taking the spatial pattern of the thickness into account, the latent variables encode more information about the discs than global thickness measurements alone, without manual demarcation of other disc features. This study extends our 2022 ARVO abstract (Wang J-K, et al. IOVS. 2022;63:ARVO E-Abstract 436) and will enable future work to compare spatial patterns as defined quantitatively by latent variables with specific morphometric disc features and clinical outcomes, as well as similar work with other optic neuropathies. 

Methods

Overview

To simultaneously quantify and visualize spatial patterns of optic disc edema in OCT, we proposed a biVAE model to convert OCT thickness maps into VAE latent spaces for a straightforward presentation. We also linked the latent spaces between the RNFLT and TRT maps (through adding a loss function term during deep-learning optimization to encourage correspondence of latent values), so the locations of the RNFL and TR data points in the latent spaces were correlated. Linking the RNFLT and TRT VAE models (as opposed to use of two completely independent VAE models) was motivated by the expected strong correspondence between the two thicknesses allowing for increased robustness by modeling them together; however, it is important to note that the biVAE also allows flexibility in deviating from the expected correspondence, as each thickness map still has its own latent variables. The design of the biVAE architecture in training is shown in Figure 1A, and the trained RNFL and TR latent space montage maps are shown in Figure 1B. The montage maps provide an intuitive view of optic disc edema increasing from practically no edema in the upper-right corner (i.e., blue region) to more severe edema (i.e., red and white regions) with various spatial patterns of swelling. It is also worth noting that the color scales in the RNFLT and TRT maps are adjusted differently to make the spatial patterns visually similar in both thickness maps. Please see more details in Supplementary Methods, Section 1. 

Figure 1.

 

(A) Flowchart of the proposed biVAE model with example paired input RNFLT and TRT maps and the reconstructed outcome images. The color fundus photograph is provided as a reference. (B) The RNFL and TR latent space montage maps were created by the corresponding VAE decoder with latent variable pairs ranging from −29 to 4 (both the x-axis and y-axis). The corresponding decoder creates each tile in the latent space montage map according to the coordinates of latent variables: (d_R1, d_R2) = (4, 4), (3, 4), …, (−28, −29), (−29, −29) and (d_T1, d_T2) = (4, 4), (3, 4), …, (−28, −29), (−29, −29). The example reconstructed images in (A) are highlighted in both montage maps in (B).

Data in Training and Test Datasets

Clinical data and findings from the IIHTT and OCT substudy were previously published, including the trial design, the treatment effects, and the OCT structural parameters over time.^{5,10,17–20,22,39} Overall, 165 patients (161 women and four men) were enrolled in the original trial; the mean age ± SD of the enrollees was 29.0 ± 7.4 years, with an age range between 18 and 52 years. The enrolled study eyes were defined by the worst eye of the IIH patients naïve to treatments with mild vision loss, a perimetric mean deviation between −2.0 and −7.0 dB using the Swedish Interactive Thresholding Algorithm (SITA) 24-2 pattern on a Humphrey Field Analyzer II with CIRRUS OCT (Carl Zeiss Meditec, Oberkochen, Germany). A subset of 125 subjects who underwent OCT imaging at each test date was included in the OCT substudy with OCT volumetric scans.¹⁰ 

In addition to the IIHTT OCT substudy data, normal subjects and those with papilledema were recruited from the Neuro-Ophthalmology Clinic at the University of Iowa and included in this study. Of the 70 papilledema subjects (primarily IIH), 64 were females (91%) and six were males (9%); the age range was from 18 to 59 years (mean ± SD, 31.6 ± 9.7 years). For the Iowa normal subjects, all 96 subjects were healthy male athletes ranging in age from 18 to 23 years (mean ± SD, 20.0 ± 1.4 years). 

In this study, right-eye oriented (i.e., all the left eyes were horizontally flipped) ONH-centered CIRRUS OCT scans were used to create the RNFLT and TRT maps to train and test the proposed biVAE model. OCT scans that were significantly de-centered or affected by severe artifacts were excluded; each available OCT scan had 200 × 200 × 1024 voxels covering 6 × 6 × 2 mm³ and had a signal strength of ≥7. The training datasets included a total of 1498 papilledema scans from the IIHTT study (763 and 735 scans from the 62 treated and 63 placebo subjects’ both eyes over time, respectively) and 791 normal longitudinal scans from the 96 Iowa healthy subjects’ both eyes. OCT scans from both eyes, multiple visits, and repeated scans were used in the training dataset to make the VAE model access as many retinal spatial patterns as possible, because papilledema is often bilateral and asymmetric between the two eyes. To verify the ability of the proposed biVAE model to reconstruct the input thickness maps, we utilized 70 independent OCT scans (single eye, single visit) from 70 Iowa papilledema subjects. The overall data distributions are described in Tables 1 and 2. The study protocol was approved by the University of Iowa's Institutional Review Board and adhered to the tenets of the Declaration of Helsinki (ClinicalTrials.gov identifier NCT01003639). 

OCT Layer Segmentation and Retinal Layer Thickness Maps

For each of the available ONH-centered OCT scans in the datasets, the inner limiting membrane, RNFL, and the retinal pigment epithelium (RPE) complex were segmented using our previously published automated 3D graph-based approach.^10,14,16 In particular, for the lower RNFL and RPE surfaces, we applied a thin-plate spline (TPS)¹⁵ to extrapolate the surfaces at the region under the swelling and above the Bruch's membrane opening (BMO).^14,16 The TPS method ensures that the smoothness of these two surfaces provides stable reference planes for creating a thickness map according to the A-scan distance.^10,14,20 Figure 2 shows an example of the OCT layer segmentation and the RNFLT/TRT maps in 2D and 3D. Details regarding color scale adjustment and OCT thickness map alignment are described in the Supplementary Methods, Sections 1 and 2. 

Figure 2.

 

An example of retinal layer segmentation and the thickness maps of an ONH-centered OCT scan. (A) OCT B-scans and the segmented layer surfaces. (B, C) The corresponding RNFLT and TRT maps in three dimensions (using the same color scale). (D, E) The RNFLT and TRT maps in two dimensions (using adjusted color scales), which are the input for the proposed VAE model.

Bi-Channel Variational Autoencoders

In this study, we extended our previous single-channel VAE model (Wang J-K, et al. IOVS. 2022;63:ARVO E-Abstract 436) to a biVAE model^40,41 that simultaneously constrains the RNFLT and TRT maps in separate latent spaces through two decoders (Fig. 1A). This novel design allows the RNFLT and TRT maps to be trained end-to-end together, and the latent variables from both decoder channels can correspond to each other but not necessarily be the same. The encoder and two decoders were trained simultaneously but independently parameterized.³² During the training process, the goal of the encoder (plus the sampling units) was to learn how to effectively convert the stacked input of the RNFLT and TRT map into two pairs of parameters (i.e., [d_R1, d_R2] and [d_T1, d_T2])—the latent variables. Meanwhile, the two separate decoders were trained to reconstruct the RNFLT and TRT maps (respectively) as close as possible to the input thickness maps using only the corresponding latent variables. Overall, the encoder is trained to decompose the input images effectively and arranges the latent variable distributions that are encouraged to be smooth. The smooth latent variable distributions can be visualized by observing the edema spatial pattern of each tile morphing in the biVAE latent space montage maps in Figure 1B, in which a greater number of red/white latent space tiles in the montage maps indicates more severe swelling. Please see details regarding the neural network architecture and loss function design in Supplementary Methods, Section 3. 

Model Evaluations

In this study, we intentionally constrained the model to use just two latent variable pairs: (d_R1, d_R2) for the RNFL-channel decoder and (d_T1, d_T2) for the TR channel decoder. We used the mean signed and unsigned (absolute) errors to quantify the encoding of our biVAE model and its decoding ability. The pRNFLT and ONHV were computed based on the original (input) and VAE reconstructed (output) RNFLT/TRT maps. The signed difference was defined by the measure from the original thickness map minus the one from the reconstructed thickness map. The mean and absolute errors were then computed for the IIHTT, Iowa normal, or Iowa papilledema dataset. 

Because the mean square error is not always the best parameter to evaluate how similar two images are, we also introduced the structural similarity (SSIM) index^42,43 to compare the intensity, contrast, and structure between the input and VAE reconstructed thickness maps. Using the Python scikit-image package,⁴⁴ the structural similarity index ranges from 0 (no similarity between two images) to 1 (two identical images). 

Results

The biVAE was trained with a penalty for dissimilar latent variable values between the RNFL and TR channels, resulting in a new way to encode information. The correlation between features in RNFLT and TRT maps is represented by the agreement between the tiles in the latent space montage maps. Excellent agreement was expected in this study due to the strong correlation between the RNFLT and TRT maps (because TR can be approximately seen as RNFL plus a thickness offset; see Fig. 2 and Supplementary Methods, Section 1). Figure 1B shows that the latent space montage maps created through the RNFL and TR channels are similar, indicating successful training of the two channels. Although input image pairs are represented simultaneously in both channels, given the similarity and limited space we focus on the RNFL figures in the remainder of the Results section and provide TR-related results in the Supplementary Results. 

We demonstrate the RNFL latent space montage map with all of the training images represented by colored dots in the Cartesian location given by their two latent variables, d_R1 and d_R2 in Figure 3. All of the data points from the IIHTT (blue dots, 763 RNFLT maps from 62 subjects in the treatment group; red dots, 735 RNFLT maps from 63 subjects in the placebo group) and Iowa normal (green dots, 791 RNFLT maps from 96 normal subjects) datasets are plotted in this way, allowing visualization of the entire study population at once. It is worth noting that biVAE is unsupervised learning, so all the labels were not involved in the training process. The background montage map in Figure 3A is identical to the left panel image in Figure 1B, but the image contrast is adjusted for better visibility of the overlying points. Because there is a high density of data points concentrated in the upper-right corner, including all of the normal subjects, a zoomed-in montage map is shown which displays the spatial patterns from the normal and mild papilledema cases. The standard color scale (from 0 to 350 µm) from the Zeiss CIRRUS OCT report is reproduced here to provide a familiar view of these reconstructed thickness maps. 

Figure 3.

 

Distribution of subject eyes in the RNFL latent space montage map, where the green dots represent the data from the Iowa normal subjects, and the red dots and blue dots represent the data from the IIHTT placebo and treatment groups, respectively. The upper-right inset, with the dashed outline, is a zoomed-in view of the montage map; the black contours are the KDE plots based on the Iowa normal data points. Three types of severe edema cases are displayed in the upper-left, lower-left, and lower-right panels. (A) In the subpanels, (a) is the color fundus photograph, (b) is the OCT central B-scan with automated layer segmentation, and (c) and (d) are the encoder's input and decoder's output RNFLT maps with the peripapillary circle, respectively. (B) pRNFLT montage maps that correspond to each tile of the original montage maps in (A). (C) The RNFL spatial patterns differ (i.e., latent values vary) even when the pRNFLTs are similar.

In the top-right zoomed-in view, black contour lines were drawn using the kernel distribution estimation (KDE) plot (generated by the Python package Seaborn)⁴⁵ to highlight the distribution of the Iowa normal data points. The intersection of the KDE plots with the montage tiles of the latent space demonstrates the variety of normal nerves, with a clear difference in shape between those with lower and higher values of d_R2. The more negative d_R2 values appear to be associated with nasal thickening and smaller nerve area typical of angled nerves. The average ± SD latent variable values of the normal subjects are −0.59 ± 1.17 for d_R1 and 0.58 ± 1.02 for d_R2. Comparing the zoomed-in view to the larger latent space, we see that the IIHTT eyes segregate from the Iowa normals reasonably well, with some overlap primarily between the treatment (acetazolamide [ACZ]) group and the normals with lower values of d_R2. The subjects in the overlap area were mainly those IIH patients treated with ACZ with papilledema resolving. 

We provide three severe swelling examples (Fig. 3A) to highlight three types of swelling patterns located in the VAE montage map. We colored each tile of the montage map according to its pRNFLT (Fig. 3B). The upper-right corner of the map contains dark blue shades corresponding to normal pRNFL thickness. More negative values of d_R1 (moving left across the graph) correspond to greater pRNFLT. This provides a straightforward visualization of the global (average) pRNFLT for individuals and populations in the latent space. 

In addition to thickness, the latent space stratifies nerves based on other morphological features. To highlight this stratification, we defined six arbitrary latent variable (d_R1, d_R2) pairs along the yellow–orange band (having similar pRNFLT) in Figure 3B (i.e., samples 1–6). Next, these six latent variable pairs were sent to the decoder to reconstruct the RNFLT maps (Fig. 3C). The biVAE model demonstrates an outstanding ability to define the continuum of RNFL spatial patterns (i.e., from a typical inversed “C-shape” [sample 1] to a more congested dome shape [sample 6]). This task is not easy for other types of interpolation approaches, demonstrating one of the major strengths of the biVAE model. 

To compare the treatment effect on IIHTT subjects over time, we investigated the study eye (considered the eye with the worse perimetric mean deviation)¹⁹ for each subject who had the OCT images available at the baseline, 3-month, and 6-month visits. This selection resulted in 42 and 37 subjects in the treatment and placebo groups having 126 and 111 OCT scans, respectively. We then created the KDE plots to help visualize the data distributions between the treatment and placebo groups over time (Fig. 4A). It is visually apparent that the treatment group (i.e., the blue contours) migrates to the upper-right corner closer to the normal subject group (i.e., the green contours) much faster than the placebo group does (i.e., the red contours) over time. Next, we created the box plots, separating the latent variables d_R1 and d_R2 and comparing them with treatment effect over time (Fig. 4B). Tukey's honestly significant difference (HSD) test showed that both latent variables d_R1 and d_R2 show a significant mean value change from baseline to 3 and to 6 months in the treatment group (using a Python package: statemodels⁴⁶; familywise error rate α was 0.05). Whereas the IIHTT OCT substudy found a significant reduction in ONHV and pRNFLT in the placebo group at 6 months,^17,20 the changes in mean values of d_R1 and d_R2 did not reach significance in this study (which is consistent with the Bruch's membrane shape measure in a previous IIHTT report).¹⁷ 

Figure 4.

 

IIHTT data distribution changes over time between the treatment (ACZ) and placebo groups. (A) KDE plots of data points in the RNFL latent space montage maps help visualize how the treatment group (blue contours) migrates to the upper-right corner close to the normal group (green contours) faster than the placebo group (red contours) does. (B) Box plots of latent variable d_R1 and d_R2 in the different treatment groups over time. The red stars indicate significant mean differences based on Tukey's HSD test (α = 0.05).

In the IIHTT dataset, 45 and 38 subjects had study-eye OCT scans with good-quality color fundus photographs for the Frisén grade assessment at both the baseline and 6-month visits. The data points by Frisén grades in the RNFL latent space montage map are shown in Figure 5A. The background montage map is color coded to indicate the pRNFLT at the location associated with each montage map tile. Figure 5B displays the trajectories of the same data points (from Fig. 5A) over time with treatment, migrating from the baseline visit to the 6-month visit. It is visually apparent that more data points in the treatment group moved closer to the normal subject distribution in 6 months than the ones in the placebo group, whereas a subset of the data points in the placebo group shifted to the right without a significant change in d_R2. This observation is consistent with a similar finding in the IIHTT OCT substudy in which the Bruch's membrane shape did not significantly improve in the IIHTT placebo group in 6 months, implying an alternative mechanism for improving the disc volume and reducing the ICP differential.¹⁷ Although we do not recommend interpreting latent variables separately due to the high nonlinearity, details of the complex relationships among the biVAE latent variables, OCT measurements, age effects, and Frisén grades are described in Supplementary Results, Sections 1 to 3. 

Figure 5.

 

IIHTT data (baseline and 6 months) relationships among the Frisén grades, pRNFLT measurements, and RNFL latent variable d_R1 and d_R2 pairs. (A) Data points with labels of Frisén grades superimposed on the RNFL latent space montage map, in which the color scale indicates the pRNFLT. (B) Data point trajectories (per eye) from the baseline to 6-month visits.

We tested the reconstruction ability of the proposed biVAE model using an independent dataset with 70 papilledema eyes (one OCT scan per subject) from the University of Iowa. After training, 70 OCT RNFLT and TRT map pairs were first sent to the encoder to compute the latent variable (i.e., [d_R1, d_R2] and [d_T1, d_T2]) pairs, then the decoder reconstructed the corresponding thickness maps back based only on these latent variables (Fig. 1A). For the 70 test subjects, the mean signed and absolute pRNFLT errors were −0.12 ± 17.34 and 13.68 ± 10.65 µm, respectively. The structural similarity indices for the original and reconstructed thickness maps through the RNFL and TR channels were 0.91 and 0.95, respectively. Table 3 lists all of the evaluation outcomes, including both training and test datasets. Figure 6 shows a scatterplot of the pRNFLT measures for the original RNFLT maps versus the biVAE reconstructed ones from the test 70 subjects; the Pearson correlation coefficient was 0.992 (P < 0.01). Three examples are shown in the scatterplot. The biVAE model demonstrated excellent reconstruction ability using just two latent variables, preserving many local details and various global patterns of the input image in the output image. More qualitative results regarding reconstruction are shown in Supplementary Results, Section 4. Note that, for the TR channel, the Pearson correlation coefficient of the volumes was 0.995 (P < 0.01). 

Table 1.

 

Distribution of ONH-Centered OCT Scans in the Training Dataset Derived From the IIHTT OCT Substudy (1498 Scans from 125 Subjects)

 

Table 2.

 

Distribution of OCT Scans in the Training and Test Datasets

 

Table 3.

 

Reconstruction Errors of pRNFLT, ONHV, and Structural Similarity Based on Input and biVAE Reconstructed Thickness Maps

 

Figure 6.

 

A scatterplot showing a high correlation of pRNFLT measurements derived from the original and reconstructed thickness maps based on 70 independent subjects in the test dataset. Qualitative examples show that the VAE model can catch prominent, local spatial features in the input images and create closely reconstructed outputs.

Discussion

The proposed biVAE model is an unsupervised method that quantifies the ONH RNFLT and TRT maps by latent variables and visualizes the integrated information in the corresponding latent space montage maps. Because the ONH spatial patterns of the disc edema may vary among nerves with similar severity, the computed latent variables can simultaneously stratify the data population by the amount of swelling and its spatial patterns, indicating factors such as the chronicity of the edema and the underlying architecture of the nerve. The clinical utility is based on the idea that the data points clustering in the latent space maps show similar ONH spatial patterns that imply similar risks of vision loss. 

As the biVAE latent variables are nonlinear indices representing the variability in the training dataset, we sought to influence the training so that normal patterns collected in one area and the swollen patterns fell on a continuum from mild to severe that was easy to see in the 2D latent space (Figs. 1B, 3–5). The treatment effect of acetazolamide (blue KDE plot in Fig. 4A) can be visualized as the collapse of the 2D distribution of the population toward the upper right, where the normal nerves are located. In Figures 4 and 5, we see that, although all eyes moved left to right as the pRNFL and d_R1 decreased, the blue treatment group had much greater normalization of d_R2 (top-right corner), but the placebo group continued to have abnormal d_R2 levels as their nerves adopted the morphology more typical of chronic papilledema (more dome-shaped). Previous analysis of the IIHTT OCT data showed that the placebo group continued to have abnormal bowing forward of Bruch's membrane even at 6 months, when pRNFLT and ONHV showed improvement, and this feature has been associated with ongoing elevated ICP.¹⁷ 

Most of the subjects’ trajectory lines in Figure 5B show movement toward the normal region (the upper-right corner). It is also intuitive and easy to spot subjects that do not follow the trend and may need special medical attention. The wide range of the structural distributions of Frisén grades in Figure 5A shows how ordinal visual cues are somewhat dissociated from structural measures and demonstrates how continuous quantitative measures such as latent variables may be better suited to demonstrate changes in papilledema during longitudinal monitoring. 

When interpreting latent variables (d_R1, d_R2) and (d_T1, d_T2), it is important to consider their locations in the montage maps and compare them with the cohort. Due to the high nonlinearity of the latent variables, we can only state that d1 shows strong correlations with the pRNFLT and ONHV (Supplementary Results, Section 1) and d2 is visually associated with the edema patterns. Although a certain level of ambiguity exists in the VAE model, the nonlinearity allows for a data-driven interpretation/visualization in the latent space without rigidly defining each latent variable. 

Our proposed biVAE model uses only two latent variables for each decoder channel (a significant constraint chosen to depict the main structural patterns of papilledema over a wide range of severity). Although this makes it easy for clinicians and viewers to quickly and effortlessly identify a given patient's eye structure over time in the 2D latent space montage maps, it may also limit the reconstruction ability of decoders. For the thickness maps in papilledema, our current design demonstrated promising encoder decomposition and decoder reconstruction ability (Fig. 6; see also Supplementary Results, Section 4), but the number of latent variables may have to be increased if the input images contain complex detailed features (e.g., color fundus or OCT en face images). 

Our segmentation approach also has a limitation in using a TPS to extrapolate the RNFL and TR lower bounding surfaces over the center of the ONH rather than strictly adhering to anatomical structures.^10,14–16 Although BMO and lamina cribrosa might offer more anatomical fidelity, accurately segmenting these structures when there is moderate to severe swelling can be challenging. Our 3D graph-based method is a workable approach to consider volumetric contextual information to segment the retinal surfaces, which can be stably sampled at the peripapillary region (which is less affected by edema). The TPS method then utilized these sample points to minimize bending energy, creating smoothly extrapolated surfaces beneath the central region of the ONH, all artificial but preserving a natural and continuous surface topology. We highlight this point to aid interpretation. 

Another limitation is that papilledema due to IIH, by its nature, predominantly affects females more frequently than males, and this fact has influenced the sex distribution of the latent space montage maps, which limits the ability of VAEs to characterize potential sex-specific differences. Additionally, our VAE model has only been trained by papilledema due to IIH and control subjects, which may restrict the generalizability of the current model to other types of optic disc swelling (e.g., ischemic optic neuropathy, optic neuritis, pseudopapilledema). We are currently training a VAE model with other causes of optic disc edema. Finally, our proposed biVAE model did not consider the information regarding each subject's racial background, and this represents another limitation, as racial disparities in papilledema manifestations may exist and could be represented in the latent space to determine any differences. 

Overall, the proposed biVAE model provides an efficient visualization and quantification for categorizing the spatial patterns of papilledema for the severity range commonly encountered in clinical practice. VAEs also provide a useful means of visualizing and quantifying changes in swelling within an eye over time to aid the clinical assessment of treatment. Because the inputs of the proposed biVAE model are flexible, we expect that adding more channels (e.g., Bruch's membrane displacement) can potentially help the VAE model differentiate true papilledema from forms of pseudopapilledema and causes of edema. VAEs also provide a new perspective to identify features of spatial patterns of edema that may have a higher risk for poor visual outcome so that more aggressive treatment can be instituted in such cases early in their clinical course. 

Acknowledgments

The authors appreciate Michael Wall, MD and Matthew J. Thurtell, MBBS providing the Frisén grades assessment of the IIHTT data at the baseline and 6-month visits. We thank Young Kwon, MD, PhD and Andrew Pouw, MD, as well as Johannes Ledolter, PhD for their efforts in providing comprehensive discussion related to the overall applications of VAEs. We also thank Louis Pasquale, MD, Tobias Elze, PhD, Brian Woods, MD, Joseph Branco, MD, and David Szanto for their comments in our research meetings. A part of this study is extended from ARVO 2022 (Wang J.-K., et al. IOVS. 2022;63:ARVO E-Abstract 436). 

Supported, in part, by grants from the Department of Veteran Affairs Center for the Prevention and Treatment of Visual Loss, Rehabilitation Research and Development (I50RX003002, I01RX003797, I01RX001786); National Institutes of Health (R01EY031544, R01EY023279); and the New York Eye and Ear Infirmary Foundation (NEI EY032522). 

Disclosure: J.-K. (Ray) Wang, None; E.F. Linton, None; B.A. Johnson, None; M.J. Kupersmith, None; M.K. Garvin, University of Iowa (P); R.H. Kardon, FaceX LLC (F)