White matter alterations in glaucoma and monocular blindness differ outside the visual system

The degree to which glaucoma has effects in the brain beyond the eye and the visual pathways is unclear. To clarify this, we investigated white matter microstructure (WMM) in 37 tracts of patients with glaucoma, monocular blindness, and controls. We used brainlife.io for reproducibility. White matter tracts were subdivided into seven categories ranging from those primarily involved in vision (the visual white matter) to those primarily involved in cognition and motor control. In the vision tracts, WMM was decreased as measured by fractional anisotropy in both glaucoma and monocular blind subjects compared to controls, suggesting neurodegeneration due to reduced sensory inputs. A test–retest approach was used to validate these results. The pattern of results was different in monocular blind subjects, where WMM properties increased outside the visual white matter as compared to controls. This pattern of results suggests that whereas in the monocular blind loss of visual input might promote white matter reorganization outside of the early visual system, such reorganization might be reduced or absent in glaucoma. The results provide indirect evidence that in glaucoma unknown factors might limit the reorganization as seen in other patient groups following visual loss.

www.nature.com/scientificreports/ age-related diseases and cause RGC degeneration and irreversible neuronal cell loss. Furthermore, research suggests a genetic link between glaucoma and AD [27][28][29] . More recently, dysfunction of the glymphatic system (i.e. a cerebrospinal fluid transport system that facilitates clearance of neurotoxic molecules) was proposed as a novel and missing link between AD and glaucoma [30][31][32][33][34][35] . At present, the most parsimonious explanation for changes to the brain tissue in glaucoma remains that these are a consequence of the decreased visual input due to retinal and optic nerve damage 19,36 . Yet, a complete picture of the alterations of the WM tissue in glaucoma is missing. Understanding the degree to which alterations in the brain tissue in glaucoma go beyond the early visual system, and the degree to which the effect of glaucoma outside of the visual system is similar to that of other diseases also reducing visual input, can help to clarify the pathophysiology of glaucoma. Answering this question is critical to advance our understanding of the disease and may, ultimately, inform future therapeutic approaches.
Our current study tested whether WM microstructure (WMM) alterations are unique to the glaucomatous brain or depend primarily on decreased visual input. To do so, the study used reproducible cloud computing methods to measure WMM in a much larger set of tracts 37,38 than previously considered. We compared WMM between two independent sets of glaucoma patients (using two independently acquired and demographically different glaucoma datasets to replicate the pattern of results) and one monocular blind patients group. We note that blindness was defined as: "Light-perception negative for at least five years, due to perforation, enucleation, evisceration or ablatio retinae leading to blindness to the affected eye". For the purpose of our study, monocular blindness acquired later in life was considered a reference model for the effect of reduced visual input to the brain. This is because due to unilateral sensory deafferentation half of the normal visual input is lost in these patients 39 .

Results
The goal of this study was twofold. First, we wanted to determine the extent of the WMM alterations in the glaucomatous brain outside of the visual WM 40 . Second, we wanted to establish whether the pattern of alterations across a large set of WM tracts was unique to glaucoma or was instead similar to those observed in a model of reduced visual input (i.e. monocular blindness). A total of 102 subjects participated in our study, including subjects with glaucoma (GL), monocular blindness (MBL) and healthy controls (HC). To validate our results in glaucoma, we used two independent datasets with major differences in the demographics and clinical characteristics of the subjects (one group from the Netherlands and another from Japan, see Methods: subjects). This resulted in a total of two glaucoma patients groups (GL1 and GL2) and three control groups (MBL, HC1 and HC2). The participant groups did not differ in mean age (GL1 vs HC1: p = 0.80, MBL vs HC1: p = 0.65, GL2 vs HC2: p = 0.90). The gender distribution did not differ in the NL subjects (GL1 vs HC1: p = 0.52, MBL vs HC1: p = 0.73), but differed in the JP subjects (GL2 and HC2 p = 0.04). The duration of the defect and the age at diagnosis/acquisition differed between the MBL and GL1 subjects (p = 0.02). All available demographics and clinical characteristics are summarized for the first cohort in Table 1 and the second cohort in Table 2.
Diffusion-weighted magnetic resonance imaging data was processed using an automated pipeline implemented as reproducible web-services on www. brain life. io (see Table 3 for the full processing pipeline 38,42,43 ). Using this pipeline, 37 major WM tracts were identified (Fig. 1a right-hand side; see also segmentation of all tracts, we computed WMM fractional anisotropy (FA) and mean diffusivity (MD) along the tract-length ( Fig. 1a left-hand side; see also Methods: Tract profile generation) 41,45 . The 37 WM tracts were also assigned to 7 categories varying from tracts highly related to vision and tract more clearly related to cognition and motor control ( Fig. 1b; see also Table 4). Categories were defined as follows: (1) Early vision: the optic radiation, OR 40 , (2) Ventral visual stream: inferior lateral fasciculus (ILF) and inferior fronto-occipital fasciculus (IFOF) 46  The average fractional anisotropy (FA) and mean diffusivity (MD) for each tract in each participant was estimated and used to compute Hedges' g effect size (independently for FA and MD) 49 . To compute g the following participant groups were compared: GL1 and HC1, MBL and HC1 and GL2 and HC2. The average (and standard error) effect size across all 37 WM tracts is reported in Supplementary Table S1 (see also Methods: Statistical analysis). This procedure returned three sets of 37 estimated g. The variable group consisted of three different levels (GL1 vs HC1, MBL vs HC1 and GL2 vs HC2) and the variable category consisted of seven levels (early vision, ventral, occipital, dorsal, vertical, callosal and other tracts ( Fig. 1 and Table 4). For all the subsequent analyses Two-Way Factorial ANOVAs were conducted independently for FA and MD. The ANOVA compared the main effects of the subjects' group comparisons (see above) and tract category and the interactions between group comparisons and tract category using the Hedges' g.
Which white matter tracts show the strongest difference between group comparisons? To gain insight into which WM tracts were most different in our individual group comparisons, and before discussing the specific statistical tests, we provide a visualization of the Hedges' g, effect sizes for all tracts ranked by magnitude (Fig. 2). We further describe the strongest effect sizes (g <>|0.7|) based on FA and MD (see also  Table S2). For the first glaucoma cohort (Fig. 2a), we found large effect sizes were observed for the right VOF (g = − 0.90), bilateral ILF (L g = − 0.72; R g = − 0.72) and bilateral OR (L g = − 0.72; R g = − 0.70). When comparing MBL vs HC1, the left IFOF (g = 0.77) revealed a large effect size. For the second glaucoma cohort, we found large effect sizes for the left OR (g = − 0.91), bilateral ILF (L g = − 0.76; R g = − 0.89) and bilateral IFOF (L g = − 0.82 R g = − 0.81). A similar analysis was performed for MD to examine effect sizes for each tract (Fig. 2b). For our first glaucoma cohort, MD did not reveal large effect sizes. On the other hand, when comparing MBL vs HC1, we found a large effect for the left IFOF (g = − 0.83) and left CST (g = − 0.71). Our second glaucoma cohort revealed large effect sizes for MD in the parietal CC (g = 0.75) and right OR (g = 0.73). As can be seen in Fig. 2a, when examining FA g estimates, tracts involved in visual processing (e.g. OR, ILF, VOF, IFOF) were negatively affected in all group comparisons suggesting decreased FA compared to controls. For monocular blindness, tracts less involved in vision, however, revealed a positive effect size indicating increased FA in comparison to controls. This pattern is uniquely different from glaucoma group comparisons: both glaucoma groups reveal a negative effect on the WMM for the majority of WM tracts, indicating FA is decreased in glaucoma compared to controls. For MD g estimates (Fig. 2b), the pattern is less apparent, but can still be observed. For monocular blindness, the majority of tracts revealed a negative effect size (i.e., decreased MD in comparison to controls). In opposition, the majority of tracts in both glaucoma cohorts revealed a positive effect size (i.e., increased MD). This suggests that in monocular blindness a decreased FA in vision tracts is accompanied by an increased FA in non-vision tracts (this could be interpreted as a plasticity and compensatory process). In glaucoma instead no plasticity was observed in non-vision tracts.  www.nature.com/scientificreports/ (g = − 0.08, SE = 0.04) and MBL vs HC1 (g = − 0.37, SE = 0.04). To determine which group comparisons differed from each other, post hoc comparisons were made using the Tukey HSD test. These tests indicated that the effect size of the glaucoma patients differed from that of MBL group (FA p < 0.0001). We repeated the same analysis for our second cohort of glaucoma patients and found that the effect size of the glaucoma patients differed from the MBL group as well (FA p < 0.0001). Furthermore, we compared effect size for FA and MD between the two glaucoma cohorts and found no difference overall for (FA: p = 0.8638 and a small but significant difference for MD: p = 0.0043). This suggests that the glaucomatous-and visual-deprived brain respond differently to the loss of visual input.
Are white matter tract categories similarly affected? Next, we were interested in establishing whether the WM tract categories were similarly or different in the group comparisons used to estimate g. For FA, a significant main effect of category was found F(6,90) = 16.14, p < 0.0001, indicating a difference between tract mean g for each WM tract category -early These results show that both for FA and MD there was a significant difference in effect size across WM tract categories. This difference in effect size across tract categories suggests that the impact of each disease to tracts beyond the early visual WM is not uniform. were interested in determining which WM tract categories differed between group comparisons. We first plotted the Hedges' g organized by tract category and group comparison for FA ( Fig. 3a) and MD (Fig. 3b). We found a significant interaction effect between group comparisons and the tract categories (F(12,90) = 3.23, p < 0.0007). To determine which tract categories differed between group comparisons based on the g estimates of FA, post hoc analysis was performed using Tukey HSD. Our post hoc analysis revealed that both glaucoma group comparisons differed from monocular blindness in tracts of the ventral and dorsal stream and callosal, vertical and other tracts (all p < 0.0024; Fig. 3a and Table 5). The early vision and occipital WM structures did not reveal differences between any of the group comparisons. None of the categories differed between the two glaucoma group comparisons. The g estimates of MD did not reveal a significant interaction effect between group comparisons and tract categories (Fig. 3b). When focussing on FA, Hedges' g was strongest and negative for the OR (early vision), with a magnitude similar across group comparisons ( Fig. 3a; green, orange and blue). This indicated that FA for the OR was reduced for all patient groups. Hedges' g for the rest of the tract categories was strikingly different between group comparisons. Notably, whereas the effect size was primarily positive for the comparison between the monocular blind patients and controls ( Fig. 3a orange), g was primarily negative for the comparisons between the two glaucoma groups and matched controls ( Fig. 3a green and blue). A positive Hedges' g for the monocular blind subjects can be interpreted as potential indication for WM reorganization processes beyond the visual system as result of prolonged visual input reduction. The opposing pattern of results in the two independent glaucoma group comparisons can be interpreted as an indication of a possible lack of reorganization in the two glaucoma patient groups. The pattern of results for MD is similar -but opposite in the sign of g. Namely, whereas the two glaucoma groups have similar and positive g values in the majority of the tract categories, the monocular blinds subjects have negative g for the majority of the tract categories. As in the case of FA also for MD we found no difference between the Hedges' g across the three group comparisons of the early visual WM (i.e., the OR).
We performed an additional analysis to explore whether the results are robust to the patient's selection criteria (see Supplementary Analysis 1). A subset of patients in the GL2 cohort (i.e those with the largest visual field defect, VFMD > − 2 dB in the better eye) were compared versus HC2. A similar pattern of results to that reported in Fig. 3 was observed indicating that our analyses were not biased by our selection criteria (see Supplementary  Fig. S1 and Supplementary Table S4).
Is there an association between white matter tracts microstructure and patients' clinical outcomes? We performed additional analyses to evidence any association between WM tract properties and a Table 4. White matter tracts of interest and their classification. C Co ol lo or r T Tr ra ac ct t N N= =( (3 37 7) ) C Ca at te eg go or ry y  www.nature.com/scientificreports/ series of clinical outcomes available (Visual Field Mean Deviation, VFMD, in the worse eye, age of diagnosis/ acquisition, duration of the disease, and retinal nerve fiber layer measured by NFI or OCT) in the three cohorts (GL1, GL2 and MBL). For all these analyses, we first selected the tracts with the strongest effect size in the previous analysis ( Fig. 2; g <>|0.7|) because these tracts demonstrated the most meaningful group differences. This approach returned 5 tracts for the GL1 cohort (bilateral OR, bilateral ILF and right VOF), 5 tracts for GL2 cohort (left OR, bilateral ILF and bilateral IFOF) and 1 tract for MBL cohort (left IFOF). We then used a Spearman's rank correlation coefficient (ρ) to associate individual tracts microstructure and the clinical outcome across subjects. In general, the associations were small (with the highest ρ's ranging between |0.5| and |0.57|) and only a few of them passed statistical significance at p < 0.05. Below, we report specific details of these additional analyses that returned significant results. For the first patient cohort (GL1), we examined the correlation between VFMD in the worse eye, retinal nerve fiber layer measured by NFI in the worse eye, age at diagnosis and duration of the disease with the FA values of the bilateral OR, bilateral ILF and right VOF. We found a positive correlation between the VFMD of the worse eye and the right ILF (ρ = 0.57, p = 0.02; See Supplementary Analysis 2 and Supplementary Fig. S2), a negative correlation between the age at diagnosis and the right ILF (ρ = − 0.57, p = 0.02), left ILF (ρ = − 0.55, p = 0.02), and left OR (ρ = − 0.52, p = 0.05). For the second patient cohort (GL2), we examined the correlation between VFMD in the worse eye and retinal nerve fiber layer measured by OCT in the worse eye with the FA values of the left OR,

Discussion
This study set out to investigate the effect of glaucoma on the human WM as it is an ongoing debate whether glaucoma is purely an eye disease or whether it affects the whole brain. To address this, we characterized whether the effect of glaucoma on the brain WM is comparable for WM tracts supporting communication within the visual system and for WM tracts serving brain structures beyond the visual system. We used two independent datasets to validate our results. We used two independent glaucoma datasets and the overall pattern of results was consistent between the two datasets. The effect of glaucoma on WM was strongest in the early visual pathways and it decreased as the tracts of interest were categorized as less involved in serving the visual system. We also used an additional control group of monocular blind individuals. This group was used as a model of reduced visual input to the brain. We compared the pattern of results in the three participant cohorts. We estimated which tracts and categories were most affected in the glaucoma subjects and monocular blind patients to test whether any change in the WMM in the glaucoma brain could be explained primarily by a reduction of visual input. We found that whereas the pattern of effects within the two glaucoma patients group was very similar, the effect was strikingly different from that measured in the monocular blind subjects. The results indicated that in monocular blind patients processes of reorganization of the large set of tracts WMM: FA was increased and MD decreased for WM tracts outside of the visual system. This can be cautiously interpreted as plasticity and reorganization following long-term visual deprivation to possibly aid behavior (even though we have no direct evidence that the pattern of alterations in FA and MD relate to behavior). In the two glaucoma groups we found no evidence for a similar reorganization process in the WM. We found that in both glaucoma groups, FA was reduced and MD increases across all tract categories as compared to controls. Overall, these results show that glaucoma shows a clear difference in the majority of the WM tracts with our model of reduced visual input (i.e. MBL). Below, we will discuss the current findings, the relevance for understanding the aetiology of glaucoma and its clinical relevance. As it becomes apparent from our results from two independent cohorts − collected at different sites with different parameters − glaucoma consistently shows that WMM is affected strongest in early vision WM structures which then gradually decreases towards cognition structures. The involvement of the early vision structures is consistent with previous studies [5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20] . Our novel findings show a graded alteration of the WMM. The most severe alterations were found in early vision structures, followed by less severe alterations in tracts categorized as motor-, or cognitive-related (Fig. 3). This pattern of graded increasing alteration was specific to the glaucoma patients and differed in the monocular blind patients.
The pattern of results reported here indirectly support the network degeneration hypothesis (NDH). NDH suggests that neurodegenerative diseases primarily spread through distinct brain networks, An interesting result to emerge from the data is that glaucomatous WM alterations are significantly different from non-glaucomatous alterations of the WM. This suggests that dissimilarities between glaucoma and decreased visual input alterations likely reflect different underlying mechanisms. In both glaucoma group comparisons, on average, a negative effect size for FA and a positive effect size for MD was observed indicating increased and decreased values in glaucoma patients compared to healthy controls respectively, while the opposite effect is observed for monocular blindness. Glaucoma groups did not differ from each other, despite the fact that glaucoma groups had different pathology (NTG vs monocular POAG). This strengthens the results of comparable patterns of alterations. As little is known about mechanisms that play a role, the discrepancies between glaucoma and monocular blindness groups may favor the hypothesis that WM alterations are not caused by decreased visual input alone, likely reflecting different underlying mechanisms. Although the interpretation of diffusion measures is not straightforward 51,52 , FA is assumed to be highly sensitive for microstructural changes in general, other scalars are more specific for the type of microstructural changes. MD is thought to be sensitive to cellularity, necrosis and edema [53][54][55][56] . We speculate that on one hand, decreased FA and increased MD in glaucoma may reflect axonal degeneration and demyelination, respectively. On the other hand 57 , in monocular blindness the decreased MD and increased FA may reflect WM maturation and high myelination, respectively [54][55][56] . In support of this latter notion, decreases in MD and increases in FA have been observed after training, supporting the www.nature.com/scientificreports/ idea that this reflects neuroplasticity or reorganization 57 . Our findings corroborate a meta-analysis of 10 dMRI studies examining diffusion values (i.e., FA and MD) in patients with glaucoma, revealing FA decreases and MD increases in the early vision structures. Therefore, we hypothesize that glaucomatous changes are a consequence of a disease progress that involves neurodegeneration, whereas in monocular blindness WM alterations are a consequence of neuroplasticity and reorganization following loss of input. Several studies have previously reported that glaucoma worsens with age 10 . An alternative proposal has been to consider glaucoma as a neurodegenerative disease that accelerates brain aging and also involves damage to the optic nerve 18,58,59 . Previous studies have suggested that visual WM structures such as the ILF or IFOF are particularly prone to age-related decline, showing the highest rate of decline as compared to other WM tracts 60 . Because our subjects were age-matched, we were not in the position to directly analyze the clinical outcomes, or WMM, as a function of the patients' age. Instead, we looked at the age at which glaucoma was first diagnosed (Supplementary Fig. S2). In line with the proposal, we found a negative association between the age at which glaucoma was diagnosed and WMM in two visual pathways-the OR and ILF. More specifically, the older the age at which glaucoma was diagnosed, the lower the FA in cohort GL1 (Supplementary Fig. S2). The ILF FA was also positively associated with VFMD, indicating the greater the visual field defect the worse the WMM. These results are especially interesting when contrasted with those in the MBL cohort. Indeed, in contrast to the glaucoma patients, results with the MBL patients did not show an association between age of acquisition or duration of the visual defect and FA. Finally, several issues stand between a full interpretation of the results above and suggest caution. In most if not all glaucoma cases, visual field loss remains unnoticed until it becomes severe, so age of diagnosis is only a marginal proxy to look at a degenerative etiology of the disease. Given the difficulty to establish disease onset, we believe that future studies will be necessary to address how age of disease onset or progression rate affects WMM, VFMD and their association.
Our research advances the scientific understanding of glaucoma by using open science methods and data sharing. In addition to improving understanding WM alterations and their potential underlying mechanisms in glaucoma, we have used the open cloud computing platform brainlife.io for secure neuroscience data analysis. The brainlife.io platform allows publishing processed data and associated analyses steps (Apps) integrated into a single record referenced by a digital-object-identifier (DOI) 38 . This is a major advancement in supporting reproducibility efforts because it removes the complexities involved with modifying the computational environment for new software, increasing the ability of researchers to easily replicate or add to prior work. All of our data and methods are fully available to the wider scientific community on https:// doi. org/ 10. 25663/ brain life. pub. 18.
Our research has several clinical implications. First, affected WM integrity shows that damage extends beyond the visual pathways in glaucoma. In addition, we show that monocular glaucoma and monocular blindness affects the WM differently. Determining the full extent of the alternation to the brain WM due to glaucoma could guide the future developments, for example suggesting the span of brain tissues to best target with neuroprotective treatments [61][62][63][64] . For glaucoma, such neuroprotective medication is currently under development (Doozandeh & Yazdani 2016). Neuroprotection may either prevent healthy RGCs from dying or "slow down" the process of already sick RGCs. Combining neuroprotection with IOP reduction-which is clinically proven to prevent RGC loss in glaucoma-may therefore optimize treatment.
As for all dMRI studies, the interpretation of differences between groups in diffusion values is not straightforward 51,52 . Differences may reflect WM abnormalities, but can also be affected by the anatomical features of a tract, such as its curvature, by partial volume effects with neighbouring fiber tracts due to crossing, and by merging or kissing fibers. Fibers originating from the occipital lobe can merge, kiss and cross with fibers from neighbouring tracts 65 . For example, WM alterations in the OR may consequently have caused (part of) the diffusion changes in neighbouring tracts such as the forcep major as well. However, changes in merging or crossing of fibers could also be part of clinical manifestation itself. Further, for our first cohort (NL), the T1-weighted image had a relatively low contrast ratio between grey-and WM compared to the second cohort (JP). The T1-weighted image is used for segmentation of ROIs and to create a GMWMI mask to constrain tractography. We resolved this by using different WM intensity level settings for each dataset. Although we visually checked all segmentations and masks, a lower contrast ratio between grey-and WM in the first cohort could result in misclassifications in WM segmentation. Additionally, susceptibility-induced distortions can only be corrected using a second diffusion scan acquired with opposing phase encoding directions. As common in clinical studies, the JP cohort did not include a second scan. As a result distortion correction for susceptibility artifacts was not applied on this dataset. Furthermore, our results may be influenced by differences in clinical outcomes for each patient. For example our subjects showed differences in disease duration, IOP treatment, glaucoma form (high-versus normal-tension glaucoma) and visual field defects (binocular or monocular glaucoma). Future studies would be necessary to either measure groups with more homogeneous clinical measures or with larger groups of patients so to allow measuring explicitly the relation between each clinical measure the the individual variability in the groups. We make all our data readily available online at DOI. We believe this is important for evaluating and reproducing our results but also for allowing future studies to build upon our current dataset. This will ultimately improve our understanding of the heterogeneous character of glaucoma.
In conclusion, glaucoma is associated with widespread WM reduction in fractional anisotropy and increase in mean diffusivity in both early vision as well as non-vision related WM tracts. Although the early vision and ventral tracts are affected in monocular blindness as well, the WM changes in glaucoma show a different and opposing pattern. Finding such widespread abnormalities in glaucoma suggests the presence of degeneration beyond what can be explained on the basis of decreased visual input. Together with degeneration of the RGCs, the degeneration of visual pathways and ventral WM tracts may be part of the manifestation of glaucoma.

Materials and methods
Data sources and App processing records. All data are de-identified and published online using the brainlife.io platform. Data and Apps are interlinked and each dataset is preserved with associated provenance information to allow other researchers to download the minimally processed or processed data and to reuse the Apps on the cloud platform on their own data. Data and Apps associated with the present project are published and archived at https:// doi. org/ 10. 25663/ brain life. pub. 18.

Study subjects.
Data from two sites were included in the current study. Data were collected among the ophthalmologic patient populations of the University Medical Center Groningen (Groningen, the Netherlands (NL)) and Jikei University hospital (Tokyo, Japan (JP)). This resulted in a total of two patient groups and three control groups. The first cohort included 17 NL subjects with primary open angle glaucoma (POAG) with a monocular visual field defect (GL1), 20 healthy age-matched NL subjects (HC1) and 14 NL age-matched subjects with non-glaucomatous monocular blindness (MBL). A second glaucoma cohort was included to validate our results. For the second cohort, we included 30 JP subjects with glaucoma with a high prevalence of normal tension glaucoma (NTG) (GL2) and 21 healthy age-matched JP subjects (HC2). The first cohort has been previously studied using morphometry of structural MRI data (T1-weighted images) 66 , whereas the second cohort has been previously studied using Tract-Based Spatial Statistics (TBSS) 18 . The inclusion criteria can be found in the Supplementary Methods. Characteristics of the NL subjects can be found in Table 1 and of the JP subjects in Table 2. In addition to Tables 1 and 2 of the average clinical characteristics of our cohorts, all individual subject data and their disease characteristics are openly available online and accessible on https:// doi. org/ 10. 25663/ brain life. pub. 18. This study conformed to the tenets of the Declaration of Helsinki. The study on the NL subjects was approved by the medical ethical committee (METC) of the University Medical Center Groningen while that on the JP subjects was approved by the METC of the Jikei University School of Medicine. All methods were performed in accordance with the relevant guidelines and regulations. All subjects gave their informed written consent prior to participation.
Clinical data acquisition. Visual acuity was measured for all NL subjects with a Snellen chart with optimal correction for the viewing distance.
Visual fields of the GL1, GL2 and HC2 subjects were assessed using a Humphrey Field Analyser (HFA; Carl Zeiss Meditec, Dublin, CA, USA) 30-2 Swedish Interactive Threshold Algorithm SITA) fast. An abnormal visual field was defined as a glaucoma hemifield test 'outside normal limits' for at least two consecutive fields; defects had to be compatible with glaucoma and without any other explanation; for a normal visual field the glaucoma hemifield test had to be within normal limits. Visual field loss is reported in terms of mean deviation (VFMD). The Frequency Doubling Technology perimeter (FDT; C20-1 screening mode) was used in the MBL and HC1 subjects. An abnormal test result was defined as at least 1 reproducibly abnormal test location at p < 0.01. Both eyes of the HC1 subjects and the contralateral eye of the MBL subjects had to be normal.
Retinal nerve fiber layer (RNFL) was assessed for the GL1 patients using laser polarimetry (GDx; Carl Zeiss Meditec, Jena, Germany); outcome measure was the Nerve Fiber Indicator (NFI) (a summary value that indicates the likelihood of glaucomatous RNFL loss). RNFL thickness of the GL2 patients was measured by means of Optical Coherence Tomography (OCT; Stratus OCT 3000; Carl Zeiss Meditec, Dublin, CA, USA).
Cohort 2 (collection site: Japan) JP subjects were imaged using a 3.0 T Siemens MAGNETOM Trio A Tim System scanner (Siemens, Erlangen, Germany) with a 32 matrix head coil. Diffusion-weighted magnetic resonance imaging (dMRI) data were collected. The following parameters were used for the dMRI pulse sequence using single shot echo planar sequence. The scanning parameters were: 64 diffusion-weighting directions including an additional b0 image, TR = 8800 ms, TE = 87 ms, field of view = 210 ~ 230, matrix = 140, slice distance factor = 5%, voxel size = 1.5 × 1.5 × 1.575 mm, b value = 1000 s/mm2 and acceleration factor = 2 (GRAPPA). One T1w anatomical image was acquired for each participant using the following T1W 3D MPRAGE sequence using the following parameters: flip angle = 9 degrees, TR = 2300 ms, TE = 2.98 ms, matrix size = 256 × 256, field of view = 256, number of slices = 176 slices and voxel size = 1 × 1x1mm3. The brain images were corrected for intensity inhomogeneity using the built-in console.
Data processing using open cloud computing services on brainlife.io. All data was processed using the reproducible, open cloud computing services available at brainlife.io 38,67,68 . The workflow used consisted of a series of 20 brainlife.io Apps (  Table 3 reports the DOIs to reproducible web services on brainlife.io and the URLs to the source code on GitHub.com. Anatomical data (T1w) processing. The T1w anatomical images were preprocessed using the FMRIB Software Library (FSL) function fsl_anat (brainlife.app.273). Using this app, the T1w anatomical images were reoriented and linearly-aligned to the MNI152 1 mm atlas using FLIRT, and non-linearly aligned to the same atlas using FNIRT. Subsequently, the anatomical images were segmented into multiple tissue types (e.g. cortical grey matter, subcortical grey matter, WM, cerebrospinal fluid and pathological tissue) using MRtrix3's 5ttgen function 69 (brainlife.app.239). A mask representing all of the voxels of grey matter-white matter interface (GMWMI) was created and used to seed tractography (see Methods: White matter microstructure modelling). Additionally, the aligned T1w anatomical image was used for segmentation and surface generation using Freesurfer's recon-all function 70 (bl.app.0). For the NL dataset, we used a WM intensity level of 105 using the flag -seg-wlo. Next, the bilateral lateral geniculate nuclei were segmented using Freesurfer's developmental seg-ment_thalamic_nuclei function 71 (brainlife.app.222). To be able to reconstruct regions of interest (ROIs) for tractography of the OR, population receptive field (pRF) maps were fit to the cortical surfaces generated from Freesurfer using a retinotopy atlas (brainlife.app.187). This app performs a retinotopic mapping and segmentation of visual areas V1, V2, and V3, as well as higher-order visual areas, from from the surface topology of the T1w anatomical images using an open source library (github.com/noahbenson/neuropythy) 72 .
Diffusion data (dMRI) processing. The DWI images were preprocessed using the procedure dwiprepoc provided by MRtrix3 (bl.app.68). The app is largely based on the DESIGNER pipeline (Diffusion parameter EStImation with Gibbs and NoisE Removal), developed to identify and minimize common sources of methodological variability. The steps consisted of performing PCA denoising, correction on Gibbs ringing artifact, motion and eddy current correction with FSL, bias correction with ANTs N4 and removal of Rician background noise. If AP and PA diffusion scans are fed into this App, FSL's TOPUP 73,74 is additionally run as well. This was the case for the NL data and as such, AP and PA diffusion scans were combined into a single corrected one by estimating the susceptibility-induced off-resonance field using TOPUP. This step was omitted for the JP data automatically by the brainlife.io App because only one diffusion scan was provided. As a result susceptibility-induced distortions were not corrected in the JP data. DWI slices whose average intensity is at least four standard deviations lower than the expected intensity were considered as outliers and replaced with predictions made by the Gaussian Process. The DWI data and gradients were aligned to the ACPC aligned T1-weighted anatomical scan using FSL's BBR. The final isotropic voxel dimension in mm was set to 1.875 and 1.5 for the NL and JP data, respectively.

Region-of-interest (ROI) generation.
A multiple region of interest (ROIs) approach was used to be able to reconstruct the OR using the roiGenerator app (brainlife.app.223). This app uses AFNI's 3dROIMaker functionality 75 to generate ROIs in dMRI space by taking the Freesurfer's parcellation volume, thalamic nuclei segmentation and pRF visual area segmentation, to generate exclusion ROIs, the lateral geniculate nuclei (LGN) and V1, respectively. Exclusions ROIs used from the Freesurfer's aseg segmentation were as follows: the bilateral cerebral WM, cerebral cortex, lateral ventricle, cerebellum WM and cerebellum cortex.
White matter microstructure modelling (tractography). Probabilistic whole-brain fiber tracking was performed using the MRtrix3 app 76 (bl.app.101) which uses MRtrix3's Anatomically Constrained Tractography (ACT) algorithm 69 . Constrained spherical deconvolution (CSD) was used to estimate local fiber orientation distributions (FOD) (Tournier et al. 2007). CSD estimates were generated with a maximum spherical harmonic order of Lmax = 8. The GMWMI mask was used as a seed mask to constrain fiber tracking. The following parameters were used: a maximum curvature angle of 45° between successive steps, a step size of 0.2 mm, a maximum track length of 250 mm and a minimum length of 29 mm. A total of 1,500,000 streamlines per parameter combination was produced.
White matter microstructure modelling (tracking the human optic radiation (THORA)). The OR was reconstructed using the THORA App (brainlife.app.347). This app will generate ORs using MRtrix3's iFOD2 algorithm with the capability of using the ACT framework and performs backtracking between the bilateral LGN and bilateral V1 as an end ROI. Previously generated exclusion ROIs were used to exclude fibers that are not part of the OR. The following parameters were used: a maximum spherical harmonic order of Lmax = 8, number of repetitions, minimum and maximum length of streamlines of 10 and 180 mm, respectively, a step size of 0.5 mm, a streamline count of 1000, a minimum FOD amplitude of 0.05 as a cutoff and a maximum curvature angle of 35°.
White matter microstructure modelling (segmentation). To examine WM differences in early vision, ventral stream, dorsal stream, occipital, vertical, callosal and other tracts, we identified and segmented a total of 37 WM tracts (see Table 4) using the app WM anatomy segmentation 37  www.nature.com/scientificreports/ than 3 standard deviations (SDs) from the core of the tract (i.e. maximum distance) as well as fibers that are more than 4 SDs above the mean length of the tract (i.e. maximum length) were removed.
Tract profile generation. Tract profile generation is performed after the reconstruction of tracts via tractography and subsequent segmentation. For each of the segmented tract, plots of diffusion metrics (i.e. FA, MD, RD, AD) were created, known as Tract Profiles 41 ; bl.app.43). The app utilizes the Compute_FA_AlongFGcommand provided by Vistasoft. For each tract the core representation is estimated by applying a gaussian weight to each streamline based on distance away from the "core". Each tract is resampled into 100 equidistant nodes along their trajectories and the weighted average metric at each node for each participant can be extracted.
Quality assurance. All raw data was manually inspected by two expert observers. Processed images were evaluated for quality of processing and data quality using the streamlined method available on brainlife.io. To check the quality of our analyses, we used a series of Quality Assurance (QA) Apps provided by brainlife.io that allows evaluating the quality of the initial data and the processing steps. Using these Apps, we generated images of each subject of the T1-weighted scan, tissue type masks, DWI scans and DWI scan overlaid on T1-weighted scan. We further plotted the response function, computed SNR and calculated statistics associated with average streamline characteristics of each segmented tract (i.e. count, volume occupied, average length and length distribution). Table 3 provides an overview of the Apps used. The output of these QA Apps was visually inspected for each participant by two expert observers. The output is openly available on https:// doi. org/ 10. 25663/ brain life. pub. 18 (see tag "QA_image") and readers can appreciate it by visualizing the data from the platform to judge the quality of the processing achieved with the pipeline. Additional information about how to visualize data using brainlife.io can be found here: https:// brain life. io/ docs/ user/ tutor ial/# visua lize-datas et.

Statistical analyses. Extraction of white matter properties of tracts
For each WM tract, we computed diffusion along the tract by resampling each tract into 100 equidistant nodes. At each node, FA MD were calculated. Individual tracts were excluded from the analysis if they had less than 50 streamlines (see Supplementary  Table S3 for overview of excluded tracts per group). Additionally, the first and last 10 nodes of each tract were excluded from the analyses since the ends of the tract may contain ambiguous voxels of the grey-white matter junction. The average FA and MD of each tract was calculated. Effect size calculation To examine the effect of glaucoma and monocular blindness on the WMM compared to healthy controls, we calculated the effect size for both FA and MD separately between groups using the average tract values. We compared subjects of GL1 vs HC1 and MBL versus HC1. To validate results, we additionally compared GL2 vs HC2 in our second cohort. Due to heterogeneity of our groups (neuroimaging parameters and disease characteristics, effect sizes were calculated using within-cohort group comparisons only. This is an approach similar to that used by ENIGMA, the largest multi-site neuroimaging study of clinical populations. The approach has been benchmarked by the consortium and it has been shown to allow comparison of very heterogeneous measurements and clinical populations 77 . The effect sizes were computed using Hedges' g for each tract using The Measures of Effect Size (MES) Toolbox 49 (github.com/hhentschke/measures-of-effect-size-toolbox). Hedges' g was calculated using the following equation (1): Hedges ′ g =ȳ 1 −ȳ 2 s p , with ȳ 1 ,ȳ 2 and s p denoting the mean of sample 1, the mean of sample 2, and the pooled standard deviation, respectively. The equation for the pooled standard deviation is s p = √ (n 1 −1)s 2 1 +(n 2 −1)s 2 2 n 1 +n 2 −2 , with s 1 and n 1 denoting the standard deviation and number of observations for sample 1, respectively, and s 2 and s 2 denoting the standard deviation and number of observations for sample 2, respectively. Additional subsets of GL2 were also compared versus HC2 and are included in the Supplementary Analysis 1.
Two-Way Factorial ANOVAs were conducted independently for FA and MD to compare the main effects of participant group and tract category and the interactions between group and category on the g estimates. The variable group consisted of three different levels (GL1 vs HC1, MBL vs HC1 and GL2 vs HC2) and the variable category consisted of seven levels (early vision, ventral, occipital, dorsal, vertical, callosal and other tracts. What are the similarities between glaucomatous WM alterations and those caused by decreased visual input? To examine whether the effect of glaucomatous alterations is similar to decreased input (as measured by monocular blindness), we investigated the main effect of group for Hedges' g FA and MD, separately, using a Two-Way Factorial ANOVA. The main effect of group consisted of three levels (GL1 vs HC1, MBL vs HC1, GL2 vs HC).
Are alterations in the visual WM similar to those in the rest of the WM? To examine whether alterations in visual WM are similar to the rest of the WM, all WM tracts were categorized based on their involvement in vision. Categorizations were made as follows: (1) early vision tracts are tracts part of the visual pathways 40 , (2) ventral stream tracts are associated with ventral stream processing 46 , (3) dorsal stream tracts are associated with dorsal stream processing 46 , (4) occipital tracts are defined as tracts originating from or projecting to the occipital lobe 47 , (5) vertical tracts are defined as vertically-oriented tracts 37 , (6) callosal tracts are tracts projecting through the body or genu of the corpus callosum 48 ) and (7) other tracts are defined as WM tracts that are not categorized in previous categories. An overview of these categorizations of WM tracts can be found in Table 4 and Fig. 1. Next, the Hedges' g value for FA and MD of each tract was averaged per tract category to create an average effect size measure of WM property changes in glaucoma groups and non-glaucomatous blindness compared to healthy controls. To examine whether the effect of glaucomatous alterations in visual WM is similar to other WM tracts), we investigated the main effect of WM tract categories for Hedges' g FA and MD, separately, using a Two-Way Factorial ANOVA. www.nature.com/scientificreports/ Interaction effect group and category To examine whether WM tract categories are differently affected between groups, we examined the interaction effect between group and categories using Two-Way Factorial ANOVA. If significant at p < 0.05, then post hocs using Tukey HSD were performed to examine which categories differed between the groups. If significant at p < 0.05, then post hocs were performed to examine which categories differed within the group. Post hoc results were Bonferroni corrected for multiple comparisons and significance was set at p < 0.05/21 (category*group) = p < 0.0024.
Investigating the relation between clinical outcomes and the white matter tracts microstructure. We investigated the association between WM tract microstructure properties and a series of clinical outcomes available for the difference datasets. More specifically, VFMD in the worse eye, age of diagnosis/acquisition, disease duration, and RNFL measured by NFI or OCT were used. For each of these analyses, we first selected the tracts with the strongest effect size in the previous analysis ( Fig. 2; g <>|0.7|) because these tracts demonstrated the most meaningful group differences. Five tracts passed the threshold of g <>|0.7|, more specifically for the GL1 cohort: bilateral OR, bilateral ILF and right VOF. Five tracts passed the threshold of g <>|0.7| for the GL2 cohort: left OR, bilateral ILF and bilateral IFOF. Only one tract passed the threshold of g <>|0.7| for the MBL cohort: the left IFOF. Spearman's rank correlation coefficient was used to associate individual tracts microstructure and the clinical outcome across subjects. More specifically for the first cohort we correlated VFMD in the worse eye, age of diagnosis/acquisition, duration of the disease, and RNFL measured by NFI with bilateral OR, bilateral ILF and right VOF. For the second cohort we correlated VFMD in the worse eye and RNFL measured by OCT in the worse eye with left OR, bilateral ILF and bilateral IFOF and for the last cohort we correlated age of acquisition and disease duration and with the left IFOF.
Demographics Group differences in age, duration of defect, age at diagnosis/acquisition were tested using two-sample t-tests. Gender differences were tested using Fisher's exact tests.