Age-related changes in the fractal dimension of the retinal microvasculature, effects of cardiovascular risk factors and smoking behaviour.

PURPOSE: This cross-sectional study investigates the association between retinal vessel complexity and age and studies the effects of cardiovascular health determinants.
METHODS: Retinal vessel complexity was assessed by calculating the box-counting fractal dimension (Df ) from digital fundus photographs of 850 subjects (3-97 years). All photographs were labelled as 'non-pathological' by the treating ophthalmologist.
RESULTS: Statistical models showed a significantly decreasing relationship between age and Df (linear: R-squared = 0.1897, p < 0.0001; quadratic: R-squared = 0.2343, p < 0.0001; cubic: R-squared = 0.2721, p < 0.0001), with the cubic regression model offering the best compromise between accuracy and model simplicity. Multivariate cubic regression showed that age, spherical equivalent and smoking behaviour have an effect (p < 0.0001) on Df . A significantly increasing effect of the number of pack-years on Df was observed (effect: 0.0004, p = 0.0017), as well as a significantly decreasing effect of years since tobacco abstinence (effect: -0.0149, p < 0.0001).
CONCLUSION: We propose using a cubic trend with age, refractive error and smoking behaviour when interpreting retinal vessel complexity.

Introduction

The eye is a unique window into health. This concept is extensively used to assess physiological and pathological processes elsewhere in the body. Many conditions affecting the systemic and cerebral microvasculature, such as diabetes, arterial hypertension, stroke and dementia, are paralleled by changes in the retinal arborizing vessel network (Patton et al. 2007; Frost et al. 2013; Cheung et al. 2019; Lemmens et al. 2020; Zhu et al. 2020). A fundus examination allows an easy and non‐invasive investigation of the retinal microvasculature.

Digital fundus images can quantify changes in the retinal vessel network using fractal dimensions (Azemin et al. 2012, 2013; Wei et al. 2017; Cheung et al. 2018). The retinal vascular network is a transport network, permitting blood supply and thus delivering oxygen and nutrients in the retinal tissue. This network appears to strive for minimal energy consumption, reflected by its conformity to Murray’s Law of Minimal Work, describing an optimal relationship between the radii of mother and daughter branches in a network (Zamir et al. 1979; Zamir 2001; Masters 2004). This law leads to the typical arborization pattern of the retina, characterized by its complexity. This complexity is ideally at a certain optimum since a more or less complex network would function at a suboptimal level. We can loosely define a fractal object as a complex geometric structure that is self‐similar, meaning that it has a similar appearance across various spatial scales (Masters 2004). The FD can be determined in a monofractal or multifractal approach. The monofractal FD is a constant for all scales, whereas the multifractal FDq describes the multifractal behaviour in different scales. Mainster was the first to describe the retinal microvascular network as a fractal dimension (FD) (Mainster 1990). Because of its degree of spatial complexity, the retinal microvascular network may represent a composite of many monofractal dimensions, making a multifractal technique better suited to characterize such arrangement (Stosić & Stosić, 2006). According to Doubal et al., multifractal FD0 (capacity dimension) would be the most appropriate measure for the complexity of the retinal microvasculature because it appeared most sensitive to small vascular changes (Doubal et al. 2010).

Significant correlations between retinal FD and age (Cheung et al. 2009; Azemin et al. 2012, 2013; Zhu et al. 2014; Wei et al. 2017; Fan et al. 2018; Van Craenendonck et al. 2020), refractive error (Li et al. 2010, 2017; Al‐Sheikh et al. 2017), glaucoma (Kolář & Jan 2008; Wu et al. 2013), age‐related macular degeneration (Al‐Sheikh et al. 2018), cardiovascular health (Liew et al. 2008; Cheung et al. 2011; Van Craenendonck et al. 2020) and diabetes mellitus (Avakian et al. 2002; Cheung et al. 2009; Zhu et al. 2020) have been reported over the last two decades. However, findings remain inconclusive, with some studies reporting an increased FD in patients with diabetes mellitus (Cheung et al. 2009; Lim et al. 2017; Orlando et al. 2017), while others describe a decreased FD (Avakian et al.2002; Chen et al. 2018; Popovic et al. 2018; Zhu et al. 2020). In arterial hypertension, FD has been investigated less extensively, but most studies – such as those performed by Cheung & colleagues and Liew & colleagues – suggest a decreased FD in patients (Liew et al. 2008; Cheung et al. 2011). A significant decrease in FD with ageing has also been reported (Azemin et al. 2012, 2013; Zhu et al. 2014; Wei et al. 2017; Cheung et al. 2018; Fan et al. 2018). However, the exact association between FD and age remains unexplained and appears to be more complex than a simple linear relationship (Azemin et al. 2013). The wide variety in study methodologies can partly explain this inconclusiveness, but additional variables, such as lifestyle and environmental factors, can also change the retinal vessel dimensions (Serre & Sasongko 2012).

The ageing‐complexity theory states that senescence associates with reduced complexity in a wide range of physiological processes and anatomic structures such as the vasculature (Lipsitz & Goldberger 1992; Goldberger et al. 2002), resulting in reduced functionality and ability to adapt to the environment (Kyriazis 2003; Vasto et al. 2010). In line with this theory, the retinal microvascular network becomes less complex with increasing age. Endothelial dysfunction is considered a key feature of vascular ageing (Vanhoutte et al. 2009). Exposure to toxic substances and unhealthy lifestyles cause a worsening of the cardiovascular risk profile, which can be reflected in structural changes in the retinal microvasculature. Interestingly, certain features of vascular dysfunction, such as oxidative stress and endothelial dysfunction, can be reversed by lifestyle interventions such as calorie restriction (Belsky et al. 2018; Most et al. 2018; Man et al. 2020), regular exercise (Seals et al. 2008; Durrant et al. 2009; Man et al. 2020), sodium restriction (Gates et al. 2004; Jablonski et al. 2013), weight loss (Pierce et al. 2008) and smoking cessation (Celermajer et al. 1993; Karatzi et al. 2007).

Many age‐related diseases have a vascular component for which retinal vessel changes can be a biomarker. However, innocuous retinal vessel changes related to physiological ageing, referred to as senescence, coincide (da Costa et al. 2016). Therefore, it is essential to differentiate between pathological and physiological age‐related retinal vessel changes. We have developed a statistical model in this study to describe the relationship between retinal FD and age over a wide age range. We investigated the relationship between cardiovascular risk factors and the retinal FD, particularly the effect of smoking and smoking cessation on the retinal FD.

Materials and methods

Study design and study population details

This single‐centre cross‐sectional study was executed at University Hospitals UZ Leuven Department of Ophthalmology (Leuven, Belgium) during summer 2018. Data were retrospectively included from individuals labelled in the hospital’s database as ‘non‐pathological’ by the treating ophthalmologist. In order to include 100 subjects per 10‐year age interval in the 0–99 years of age range, this cohort was enriched with participants recruited from the Department of Ophthalmology UZ Leuven. Exclusion criteria included insufficient clarity of optical media for retinal imaging; congenital pathology of the posterior segment; a medical history of retinal neovascularization, diabetic retinopathy or other retinal vascular pathology, retinal dystrophy or marked atrophy; malignancy; uncontrolled arterial hypertension, diabetes mellitus or hypercholesterolaemia; and participants unable or unwilling to give consent.

We got written informed consent from each volunteer before inclusion in the study in compliance with relevant regulations on clinical trials. All data were anonymized before analysis. The study was conducted in compliance with the principles of the European Union Directive on Clinical Trials (2001/20/EC) and all local/regional requirements required to conform to the Declaration of Helsinki (World Medical Association, Edinburgh, 2000). Approval was issued by the Ethics Committee of the University Hospitals Leuven before the study began (Reference numbers S59048, S59830 and S60649).

Data collection and image analysis

General and ophthalmological history, including retinography

Data on personal general history were collected from the electronic patient record for all eligible subjects, including their cardiovascular risk profile. In addition, information on personal ocular history was compiled. Digital, 45° field of view, optic disc centred retinal photographs with a spatial resolution of 5 megapixels, captured with the Visucam PRO NM (Carl Zeiss Meditec AG, Jena, Germany), were used for retinal image analysis. If only one of both eyes of a subject was labelled as non‐pathological, that eye was included. If both eyes were non‐pathological, a random eye was included in the study.

Retinal image analysis

Technical factors can affect fractal dimension calculations. We standardized the protocol for image capturing in the clinical setting (location, darkened room, same camera type and camera settings) to avoid this. The remaining variation is cancelled maximally by using image pre‐processing and post‐processing. These steps ensure minimal impact of variations in brightness, focus and contrast. An element in our processing pipeline is the selection and inversion of the green channel, after which the application of a Gaussian filter reduces the noise. Photographs were analysed by a trained grader using the MONA REVA software (version 2.1.1) developed by VITO NV (Mol, Belgium; http:\\mona.health – Luyten et al. 2020). This software package allows the calculation of the monofractal dimension (Df) of the retinal vascular network using the box‐counting method in a predetermined region of interest (ROI). The region of interest was set to the annulus between 1.5 and 4 times the radius of the optic disc, centred on the optic disc (see Fig. S1). This outer radius was the largest zone that could be demarcated in all images. In each image, the optic disc was automatically located with the possibility of manual correction. Before cropping the images to this region, a computer algorithm automatically segmented the retinal vessels. The segmentation algorithm is based on a multiscale line filtering algorithm inspired by Nguyen et al. (2013) and post‐processing steps such as double thresholding, blob extraction, removal of small connected regions and filling holes. The graders (CL, RP, AS, JV) were first trained by SL, after which they analysed each about one‐fourth of the images, verified and corrected vessel segmentation with the MONA REVA vessel editing toolbox using an agreed protocol (cfr. Supplement 1). SL checked the quality and the consistency of the results randomly. The agreement between the different graders was very good, as evidenced by an intraclass correlation coefficient of 0.87, based on a random set of 50 fundus images. Images that did not allow to apply the protocol correctly were excluded.

The monofractal dimension Df was computed using the sliding box‐counting method on the segmented vessel image. Boxes with side length δ slide across the image, and for each δ the mean of the pixel count (N) for the segmented image is recorded. The fractal dimension Df is expressed as follows (Ţǎlu 2013):

Theory/calculation

Multiple imputation

A multiple imputation (MI) technique (Rubin 1976, 1987) was used to handle the incomplete data structure (Table 1). This approach is the current golden standard to deal with missing data and has been accepted by the European Medicines Agency for application in clinical trials (O’Neill & Temple 2012; Dziura et al. 2013). Briefly, three steps are involved in the MI process. In the first step, each missing value is replaced with M imputations, sampled from the predictive distribution of what is missing, given what is observed. As a result, M completed datasets are obtained. In the second phase, properly chosen statistical inferences are made on each dataset, implying M estimates of the model parameters. In the last step, Rubin’s rules (Rubin 1987) combine these M estimates into a single set of parameter and precision estimates. Multiple imputations with M = 100 are applied using a wholly conditional specification approach (Bartlett et al. 2015).

Table 1

Baseline characteristics of the study cohort.

Characteristics, n = 850		Missing data, n (%)
Monofractal Df, range (median)	1.307–1.611 (1.395)
Age (year range), n (%)	3–97 (51)
0–9	71 (8)
10–19	49 (6)
20–29	101 (12)
30–39	91 (11)
40–49	107 (13)
50–59	101 (12)
60–69	114 (13)
70–79	110 (13)
80–89	94 (11)
90–99	12 (1)
Male, n (%)	398 (47)
Right eye, n (%)	464 (55)
Spherical equivalent – Diopter, mean ± standard deviation (missing)	−0.43 ± 2.65	191 (22)
Macular drusen, n (%)	13 (2)
Smoking status, n (%) (missing)	604 (71)	246 (29)
Current smoker	50 (7)	97 (11)
Pack years, range (median)	0.5 ‐ 90.5 (20.0)	14 (28)
Years since tobacco abstinence, range (median)	0.0 – 0.0 (0.0)
Former smoker	143 (24)	246 (29)
Pack years, range (median)	0.0 ‐ 60.0 (10.0)	20 (14)
Years since tobacco abstinence, range (median)	0.5 – 66.0 (25.0)	10 (7)
Never smoked	411 (68)	246 (29)
Pack years, range (median)	0.0 – 0.0 (0.0)
Years since tobacco abstinence, range (median)	3.0 – 93.0 (34.0)
Diabetes, n (%) (missing)	61 (7)	14 (2)
Arterial hypertension, n (%) (missing)	209 (25)	17 (2)
Hypercholesterolaemia, n (%) (missing)	196 (24)	19 (2)
BMI ‐ kg/m², mean ± standard deviation (missing)	25.10 ± 4.03	278 (33)

BMI = body mass index.

Descriptive statistics

Descriptive statistics for continuous and categorical variables in the original cohort are reported with means, standard deviations and proportions respectively. All statistical inferences were conducted with SAS® (version 9.4; SAS Institute Inc., Cary, NC, USA). We compared several descriptive statistical models with increasing levels of complexity: from linear, quadratic and cubic regression, over penalized B‐spline up to locally estimated scatterplot smoothing (LOESS) to identify the best approach to describe the relationship between retinal Df and age. The influence of binary systemic variables on retinal Df was assessed using linear regression on the baseline variables to apply the imputation rule properly. We have performed normality tests. Correlations between continuous systemic variables and retinal Df were evaluated with Fisher‐Z tests. Statistical significance was considered when the two‐sided p‐value < 0.05. Potential confounders were considered in a multivariate cubic regression model, selecting significantly contributing variables (p < 0.05).

Backward elimination was used as the variable selection method. Initially, a complete model including all variables was considered. Variables were then deleted and refitted one by one until all the remaining variables were considered to have some significant contribution (p < 0.05) to the outcome of interest.

Multiple linear regression was applied to assess the effect of pack‐years and smoke‐stop on the complexity of the retinal vascular network. A pack‐year is the equivalent of smoking 20 cigarettes (one pack) per day for 1 year. The analysis was performed in the subpopulation aged 25–75 years, and the model was corrected for the individual’s age.

Results

Study population characteristics

The final cohort consisted of 850 subjects and comprised 47% men. The median age was 51 years, with a range from 3 to 97 years. Detailed baseline characteristics are listed in Table 1.

All subjects with diabetes, arterial hypertension, and/or hypercholesterolaemia received medication, which resulted in disease control. The amount of missing data is indicated in Table 1, and these concern spherical equivalent, smoking status and body mass index (BMI).

Retinal fractal dimension in association with age

Parametric (linear, quadratic, cubic, penalized B‐splines) and non‐parametric (LOESS) regression techniques were performed to assess the Df‐age relationship. All parametric analyses show a decreasing significant relationship (e.g. linear: R‐squared = 0.1897, p < 0.0001; quadratic: R‐squared = 0.2343, p < 0.0001; cubic: R‐squared = 0.2721, p < 0.0001). Similar monotonic decreasing trends are observed between cubic, penalized b‐spline & LOESS regression. Of all simple applied parametric methods of analysis, the cubic regression model appears to offer the best compromise between accuracy and model simplicity over the age range, compared to linear and quadratic regression (linear: LR = 91.1896; df = 2; p < 0.0001; quadratic: LR = 43.077; df = 1; p < 0.0001). The benefit of a more complex model is negligible. A comparison between the different techniques is presented in Fig. 1.

Fig. 1

Monofractal (Df) as a function of age, as well as the five statistical models describing the relationship between FD and age. Note the similarity between the results obtained with the cubic regression model and the more complex penalized B‐spline and LOESS techniques.

Fractal dimension (Df) in association with cardiovascular risk factors

Table 2 summarizes the results of the univariate analysis, showing a significant (p < 0.05) result for the relationship between the Df metric and most variables. Of note, Df did not differ significantly between right and left eyes (p = 0.4266). A significant correlation was noted between Df and spherical equivalent (Pearson correlation 0.168053, p < 0.0001, average increase in Df of 0.001398 per diopter increase in spherical equivalent), and in subjects with macular drusen, Df was on average 0.031255 lower compared to subjects without drusen (p = 0.0028). Drusen are retinal pigment epithelium deformation or thickening that may form irregularities and undulations.

Table 2

Results of univariate linear regression and correlations between fractal dimension Df of the retinal vasculature and cardiovascular risk parameters.

Variable	Effect on Df (for binary variables) Correlation with Df (for continuous variables)	p
Sex	In females 0.000096 higher	0.9704
Diabetes	In (treated) diabetes 0.001224 lower	0.0133
Arterial hypertension	In (treated) arterial hypertension 0.018003 lower	<0.0001
Hypercholesterolaemia/dyslipidaemia	In (treated) hypercholesterolaemia/dyslipidaemia 0.014290 lower	<0.0001
Body Mass Index (BMI; kg/m²)	Pearson correlation: −0.229527	<0.0001
Smoker (current or former)	Meaningless	0.0048
Current smoker	In current smokers 0.008691 higher (compared to former and never smokers)	0.1136
Former smoker	In former smokers 0.012350 lower (compared to current and never smokers)	0.0002
Pack years	Pearson correlation: −0.041230	0.3017
Years since tobacco abstinence	Pearson correlation: −0.314599	<0.0001

Significant p‐values (p < 0.05) in bold.

Multivariate cubic regression was performed to consider potentially confounding effects. Backward and forward selection of the most significant effects have been used and yielded an equal selection of variables. The final multivariate model describing FD is presented: with

Df = fractal dimension as monofractal.

SE = spherical equivalent.

I1Smoker = 1 for current smokers and 0 otherwise.

I2Smoker = 1 for former smokers and 0 otherwise.

After correction of age, spherical equivalent and smoking behaviour are significantly affecting Df. Variables that did not make the model selection were diabetes (p = 0.8947), hypercholesterolaemia/dyslipidaemia (p = 0.6577), years since tobacco abstinence (p = 0.5575), arterial hypertension (p = 0.5425), AMD (p = 0.2318) and BMI (p = 0.1175).

To appraise the association between smoking behaviour and retinal Df, we focused on a subset of the cohort for which the Df‐age relationship is stable (see Fig. 1) and which contains most smokers (96% of current smokers and 70% of former smokers). By doing this, the confounding effect of age was minimized. The results of the multiple linear regression showed a significant positive effect with an increase in the number of pack‐years on the fractal dimension (effect 0.0004, p‐value 0.0017), and a significant negative effect on the fractal dimension with an increased weighted proportion of years since tobacco abstinence (effect −0.0149, p‐value < 0.0001). This observation means that smoking behaviour significantly affects the complexity of the retinal vascular network, with an increase in complexity proportional to the number of pack‐years. In contrast, the complexity of this network significantly decreases proportionally to the number of years since smoke‐stop.

Discussion

Following the ageing‐complexity theory, previous studies report a linear attenuation of the retinal vascular network with increasing age (Azemin et al. 2012; Zhu et al. 2014; Wei et al. 2017; Van Craenendonck et al. 2020). However, only one publication addressed the non‐linear relationship between age and FD, concluding with a quadratic model based on analysis of fundus photographs of 215 subjects aged between 10 and 73 years (Azemin et al. 2013). The current study comprises a much larger cohort of 850 subjects aged between 3 and 97 years. Our study identifies a cubic model as the most appropriate to describe the relationship between age and the box‐counting FD (Df). Figure 1 shows a better performance of the cubic model than the quadratic, mainly below 40 years and above 80 years. This observation emphasizes the importance of adequate representation of a broad age spectrum.

We confirm associations between cardiovascular risk factors and the retinal vascular complexity as assessed by Df. Although univariate analysis showed a significant (p < 0.05) relationship between Df and most variables, the association with diabetes, arterial hypertension, hypercholesterolaemia/dyslipidaemia and BMI no longer hold in the final multivariate model. However, all subjects in our cohort suffering from diabetes mellitus, arterial hypertension and/or hypercholesterolaemia/dyslipidaemia were under treatment with adequate disease control. Of note, the absence of BMI in the multivariate cubic regression model could be the consequence of a correlation between BMI and age and/or smoking status. This could explain why the multivariate cubic regression yielded only age, refractive error and smoking behaviour as significant variables. Consistent with previous reports (Li et al. 2010, 2017; Lim et al. 2011, 2017; Azemin et al. 2014; Yang et al. 2016; Al‐Sheikh et al. 2017; Tai et al. 2017), an effect of refractive error on the geometry of the retinal vascular network was observed in the current study, showing an increase in Df of 0.001398 per increase in the spherical equivalent of 1 diopter. Previous reports on the effect of smoking behaviour on retinal vessel complexity considered smoking a simple binomial snapshot variable, without considering past smoking, the number of pack‐years, or the years that have passed since smoke‐stop. The current study shows a significant effect of pack‐years and years since smoke‐stop on the retinal Df. Whereas more pack‐years were associated with a higher Df, a more extended period of tobacco abstinence was associated with a lower Df. This finding aligns with tobacco's known pathological angiogenic effect, mediated by multiple growth factors, including vascular endothelial growth factor (Heeschen et al. 2001; Lee & Cooke 2012). A more complex vascular network may compensate for smokers' decreased blood oxygen‐carrying capacity and decreased blood‐O2 affinity to meet the retinal high oxygen demand.

We should interpret the present study within the context of its strengths and limitations. On the one hand, the large cohort, with a fair representation of all age ranges, is one of the significant strengths of this study. It is also the first study to describe the relationship between age and retinal Df across a large age span. It is the first to assess the relationship between smoking behaviour and retinal Df in more detail, considering the complexity of this variable. One limitation of this study is missing data, at least partly compensated for by the thorough statistical approach applying the multiple imputation technique. The number of current smokers is a small subset of our cohort. Another limitation is the exclusive use of a monofractal analysis to describe retinal vessel complexity. The software package available for this study allows for the calculation of the monofractal FD by use of the box counting method. The retinal vessel network can also be quantified using multifractal metrics, providing additional information about the retinal microvasculature. However, monofractal analysis is widely used, but a systematic assessment of the effect of age and smoking behaviour is absent. This investigation was this study’s primary aim. A combined mono‐ and multifractal analysis could have contributed to the robustness of our findings, and we acknowledge the importance of such a follow‐up study. Using only monofractals may explain why diabetes, arterial hypertension and hypercholesterolaemia/ hyperlipidaemia did not make the model selection for Df. However, another explanation is that all patients diagnosed with a clinical condition were treated and in follow‐up. We can assume that the conditions are under control and therefore have a negligible effect on the status of the retinal vasculature.

Of note, there is a lack of standardization across the many publications related to retinal vascular complexity in terms of the acquisition method and settings and the algorithms used for calculating fractal dimensions from fundus images. A detailed discussion of this observation is beyond the current scope but should be considered when comparing different studies' results.

Our findings imply that age, refractive error, past and current smoking behaviour are relevant when interpreting retinal vascular fractal dimension in public health and clinic research. This aspect is essential when studying retinal diseases or systemic diseases such as Alzheimer’s disease, stroke or even COVID‐19, all characterized by a vasculopathogenic component.

Conclusions

Ageing comes with a decrease in the complexity of the retinal vascular network. We can assess this complexity from digital fundus images by calculating the monofractal dimension (Df). A cubic regression model can adequately describe the relationship between Df and age. Besides age, spherical equivalent and smoking status significantly affect the retinal Df. Smoking leads to a compensatory increase in vascular complexity proportional to the number of pack‐years, whereas this complexity significantly decreases proportionally to the number of years since smoke‐stop.

Figure S1. Indication of the annulus representing the region of interest (shaded area) and of retinal vessels up to the third degree (arrows and colours).

Click here for additional data file.

Supplementary Material. Protocol applied for vessel editing.

Click here for additional data file.

Effect	Parameter	Estimate	p‐Value
	β_0	1.477263	<0.0001
Age	β_1	‐0.004728	<0.0001
Age^2	β_2	0.000086353	<0.0001
Age^3	β_3	‐0.000000524	<0.0001
SE	β_4	0.001398	0.0025
I1_Smoker	β_5	0.013358	0.0066
I2_Smoker	β_6	0.001533	0.6356
	ε	∼N0,σ^2