A single patient with a fully matured, functioning radio-cephalic AVF (RCAVF) was consented for this study.

Inclusion Criteria: Fully matured upper limb AVF being used for regular haemodialysis (HD) with 2 needles. 18 years of age or older and able to provide informed consent.

Exclusion Criteria: Previous surgical or endovascular salvage of the AVF. Previous central venous catheter on the side of the AVF. Repeated episodes of elevated venous pressures or high recirculation rates on HD in the previous 3 months.

At the time of recruitment, the patient had been using the AVF for regular HD for a period of 3 months. Following an initial ultrasound examination and magnetic resonance imaging (MRI) scan, the patient was followed up with ultrasound surveillance over the next 8 months allowing for structural and functional changes in the vascular anatomy to be monitored and compared with the initial ultrasound derived turbulence intensity ratio (USTIR) and shear stress measurements.

The ethical considerations surrounding the use of radiological contrast agents in patients with end stage renal failure (ESRF) meant that axial imaging was limited to a single timepoint. Computed tomography angiography (CTA) contrast agents are nephrotoxic and risk depleting any residual kidney function. It was decided that imaging would be obtained from contrast enhanced magnetic resonance angiography (CEMRA). The rapid change in flow direction at the anastomosis means that resolution of non-contrast magnetic resonance angiography (MRA) can be limited in this region. Due to risks of nephrogenic systemic fibrosis (NSF), gadolinium based contrast agents (GBCA) should be used with caution in patients with ESRF [7], however following a literature review and full risk assessment we were unable to identify any reported cases of NSF resulting from an isolated administration of Dotarem® (gadoteric acid).

The original timeline for the planned surveillance scans was adapted to match the patient’s clinical requirements. Table 1 shows the time points used for data collection.

We have previously developed a standardised methodology, using a diagnostic ultrasound scanner, to quantify the degree of turbulence at a selected point within an AVF, and have demonstrated a correlation between elevated USTIR and NIH development in the efferent vein of maturing AVF [8].

Cardiac gated Doppler frequency spectra were recorded as audio data in a 24 bit .wav format, with a sample rate of 44.1 kHz. Recordings of 25 cardiac cycles were obtained from each acquisition zone, Fourier transforms were performed, and an ensemble averaging technique was used to isolate the non-cyclical components of the resulting spectrograms. The isolated components of the frequency spectra, assumed to correspond to random fluctuations due to turbulence, were averaged over the cardiac cycle and then normalised by the average frequency shift at each location to obtain the USTIR. We have previously detailed the data collection and analysis protocol, and validated this technique in patients with surgically created AVF [8].

USTIR was calculated at 1 cm intervals along the radial artery and the cephalic vein. The complete AVF circuit was also assessed for clinical functionality and the presence of any hamodynamically significant stenotic disease. A moderate stenosis was identified in the radial artery, which was also evident on the CEMRA and the CFD mesh, however it did not appear to significantly alter flow waveforms at the level of the anastomosis. Despite the irregularities in the luminal diameter of the cephalic vein, and the presence of mild NIH at the anastomosis there were no flow limiting lesions identified throughout the efferent veins (anastomosis to subclavian vein).

Follow up ultrasound scans were conducted to assess for structural remodelling of the vessels and NIH development, along with changes in access flow over time. Where detected on ultrasound, NIH was graded as minor, moderate or significant, based on sonographic appearance and peak systolic flow velocities (PSV), as defined in Table 2.

The geometry for the simulation was based on the CEMRA of the upper limb. The data were exported in Digital Imaging and Communications in Medicine (DICOM) format and axial image slices were loaded into Simpleware${}^{\text{TM}}$ Scan IP (Synopsys Inc, Irvine, CA, USA).

A 3D volume was generated from the segmented images and a smoothing algorithm was applied to compensate for the image artefact introduced by the discrete slice thickness of each image. The smoothed volumetric model was then imported into Ansys MESH (Ansys Inc, Canonsburg, PA, USA) for initial meshing and refinement. Prior to meshing comparisons were made between luminal diameters of the original CEMRA images and the imported fluid volume at selected regions (anastomosis, max diameters and min diameters), and were found to be comparable within the limits of the resolution of the original scan. The only exception to this was the moderate narrowing in the inflow artery, which was slightly exaggerated on the smoothed model, however given that this region of the model was not our main area of interest it was decided that the improvements gained by smoothing the model justified the introduction of this minor discrepancy in the geometry. An unstructured mesh consisting of tetrahedral elements was created, and wedge shaped elements were used to generate inflation layers to ensure optimum resolution in our main region of interest, adjacent to the vessel walls.

The mesh was refined until successive versions showed less than 10% variation in the average oscillatory shear index (OSI) values for each vessel segment, ensuring good agreement on both distribution and magnitude of oscillatory shear between successive simulations. In order to achieve this, considerable refinement was required and the final mesh consisted of approximately 3.5 million elements, and incorporated 10 inflation layers with a minimum height of 0.025 mm and a total thickness of 0.65 mm.

The magnitude and distribution of a range of haemodynamic parameters were obtained from the results of the CFD simulation. Our primary interest was in identifying regions of oscillatory flow, quantified using the OSI, which has previously been linked to NIH formation in vascular access AVF [13], and intuitively was the most likely metric to map to regions of high USTIR.

We know that when an AVF is formed the act of bypassing the distal capillary bed results in a sudden decrease in distal resistance and a corresponding increase in flow. This leads to an elevation in wall shear stress (WSS), triggering a vasodilatory response, which plays a key role in the vascular remodelling observed following AVF formation. Due to the pulsatile nature of blood flow WSS is continuously variable over the cardiac cycle, for this reason time averaged values of WSS (TAWSS) are often quoted rather than absolute maximum or minimum values.

$$TAWSS=\frac{1}{T}{\int }_0^T|WSS|𝑑t$$

The disturbed nature of the flow around the anastomosis of a native AVF means that the direction and magnitude of WSS is continuously changing throughout the cardiac cycle. Ku et al. [14] devised the OSI to quantify the deviation of WSS from the average direction.

$$OSI=0.5\left[1-\frac{\left|{\int }_0^{T}WSS𝑑t\right|}{{\int }_0^T|WSS|𝑑t}\right]$$

It has been previously shown that both low TAWSS and high OSI may be contributing factors to NIH formation and resultant vascular access dysfunction [13, 15, 16].

Another quantitative metric of shear forces, commonly encountered in literature in relation to AVF failure, is relative residence time (RRT), which combines the TAWSS and OSI, and is defined as [17]

$$RRT={[(1-2\cdot OSI)\cdot TAWSS]}^{-1}$$

A variation of RRT is Highly Oscillatory Low MagnitudE Shear (HOLMES), which is equal to half of the reciprocal of RRT. HOLMES has the advantage of providing a linear index, which is more conceptually intuitive, providing a value corresponding to TAWSS, which continues to drop in the presence of high OSI [18]. Holmes has previously been shown to be a predictor of NIH growth in venous bypass conduits exposed to arterial flow in the lower limbs [19].

$$HOLMES=TAWSS\cdot (0.5-OSI)$$

The CFD model was divided into segments, 1 cm in length, corresponding to the USTIR observation sites (Fig. 1). The shear stress parameters were averaged over the vessel walls for each of the segments allowing for magnitude and distribution of OSI and HOLMES to be compared with our ultrasound-based predictions of turbulent flow regions. We also recorded the maximum OSI and minimum HOLMES values for each vessel segment.

Fig. 1.

Observation sites used for analysis. Radial Artery: A, B, C, D. Cephalic Vein: a, b, c, d, e, f.

Different fluid properties and turbulence models were compared, whilst the patient specific mesh and boundary conditions were identical for all 3 simulations. Fig. 2 shows the average OSI in the pre-defined regions of the cephalic vein, as produced by each of the 3 simulations.

Fig. 2.

Comparison of OSI in the cephalic vein based on 3 different simulations. OSI, oscillatory shear index.

The Newtonian and non-Newtonian fluid models provide very similar results for the low OSI regions, however there is some degree of discrepancy between the 2 models in zones d, e and f, which correspond the regions with the highest OSI. Given that the assumption of approximated Newtonian behaviour of blood, is only valid for laminar flow in large vessels, these results suggest that the addition of the non-Newtonian model is necessary to accurately quantify OSI in regions of disturbed flow. It should be noted however, that despite the quantitative differences between these 2 models the OSI distribution patterns remain similar, with both simulations successfully differentiating between regions of high and low OSI.

The addition of the k-$\omega$ SST model resulted in a reduction in the OSI in all regions of the cephalic vein, except for zone a; the juxta-anastomotic segment, where calculations of OSI was significantly higher than the 2 simulations based on the laminar flow model. The k-$\omega$ SST model is known to over-estimate turbulence in regions of rapidly accelerating, pulsatile flow, which may explain this discrepancy. All other zones showed very similar distributions of OSI to the non-Newtonian, laminar flow model, however the magnitude was reduced across the entire vessel.

Due to the limitations in validating the k-$\omega$ SST model, it was decided that the non-Newtonian, laminar flow simulation would be used for the remainder of the analysis. Variations in the magnitude of OSI on the laminar and non-laminar simulations are not of great significance as we do not have a range of normal, or expected values for comparison, so for the purposes of this case study it is the distribution of OSI that is of interest.

Figs. 3,4 show the distribution of OSI and HOLMES respectively, across 3 different projections of the 3D geometry. The colour map used in Fig. 3 to represent the values of OSI throughout the modelled anatomy provides useful visual guide to both localised OSI maximums, as well as regions with globally elevated oscillatory shear. In the artery slight elevations are seen distal to regions of vessel curvature, and also at the distal edge of the anastomosis, where retrograde flow from the hand mixes with the antegrade flow, which provides the majority of the inflow.

Fig. 3.

Distribution of oscillatory shear index (OSI) as calculated by the computational fluid dynamics (CFD) simulation.

Fig. 4.

Distribution of Highly Oscillatory Low MagnitudE Shear (HOLMES) as calculated by the computational fluid dynamics (CFD) simulation.

Distribution patterns in the vein appear far more complex, with elevated OSI spread over larger regions with greater maximum values. Interestingly the juxta-anastomotic segment of the vein (zones a and b on Fig. 1), show relatively low OSI, with values increasingly as flow travels along the vein into the dilated segment (zones d, e and f on Fig. 1).

The colour map used in Fig. 4 is less intuitive due to the wide scale required to capture localised regions of high HOLMES. This is problematic as it means that the regions of low HOLMES, in which we are interested are poorly represented by the chosen scale. Although the scale is easily adjusted to allow us to better highlight these regions, the presented images provide a useful insight into why metrics which incorporate TAWSS may prove less reliable in the high flow environment typical of AVF. It is clear from Fig. 4, that the distribution of HOLMES in both the venous and arterial sections of the model, is heavily dependent on the vessel diameter. We know that under laminar flow conditions WSS is proportional to flow velocity, and therefore vessel diameter. AVF typically exhibit high velocity flow, and commonly have irregular luminal diameters, resulting in large differences in TAWSS throughout the flow circuit.

If we compare the measurements of average HOLMES obtained from the three flow simulations (Fig. 5), whilst we still see variations in magnitude between the 3 flow models, they converge in regions of low HOLMES, and the overall distribution patterns are very similar. This suggests that this metric is less sensitive in detecting subtle differences in the flow patterns, however the limitation is not with the metric itself, but likely results from our chosen methodology of averaging values over relatively large vessel segments.

Fig. 5.

Comparison of Highly Oscillatory Low MagnitudE Shear (HOLMES) in the cephalic vein based on 3 different simulations.

Despite the fact that this AVF was considered fully matured at the start of this study, the cephalic vein underwent significant structural remodelling over the course of the follow up period.

The follow up was stopped after 13 months due to the need for a fistulaplasty to restore diminished flow volumes through the AVF. We decided that following this procedure, further comparisons to our CFD simulation would be potentially misleading, and in the absence of post plasty axial imaging an updated CFD model was not possible.

At the 10-month ultrasound scan, 2 ectatic segments had developed, corresponding to the regular needling sites in the upper half of the forearm. Moderate NIH was observed in the juxta-anastomotic segment of the cephalic vein, becoming increasingly prevalent in the segment of the vein immediately distal to the needling zones. No significant, focal stenotic lesions were identified, but flow volumes were slightly reduced when compared with the baseline ultrasound scan.

The 13-month ultrasound showed progression of NIH in the cephalic vein, resulting in a significant stenosis, extending from 22 mm above the anastomosis up to the level of the ectatic, needling zone. Flow volumes in the feeding artery had decreased from 700 mL/min to approximately 250 mL/min. High recirculation rates were reported on dialysis and the patient was referred for venoplasty of the compromised vessel segment.

Figs. 3,4 show the distribution of OSI and Holmes respectively, whilst the data presented in Table 3 shows the maximum, minimum and average values across the pre-defined regions of the upper limb anatomy (Fig. 1). Fig. 6 shows regions of NIH development, as visualized on the 13-month ultrasound.

Fig. 6.

Ultrasound images showing development of neointimal hyperplasia in the efferent vein.

In the venous portion of the AVF the segments with moderate or significant NIH development corresponded with the regions of lowest average HOLMES and average WSS values, as well as the highest average OSI. The absolute maximum OSI and minimum HOLMES recorded from each venous segment did not map to the regions of NIH development, but this is not surprising given the non-uniform distribution of WSS throughout the AVF circuit.

No obvious relationship was identified between any of the WSS indices and development of NIH in the arterial segments of the AVF, however given the lack of significant structural remodelling which occurred in the artery, meaningful interpretation of the data obtained from this region of the model is very limited.

USTIR is highly elevated at the anastomosis despite the presence of only minor NIH in this region on the 13-month ultrasound, and the comparatively low OSI, and high HOLMES values. This is likely due to the inertial forces resulting from the sudden change in flow direction at the anastomosis leading to an elevated USTIR in the absence of true turbulence. The streamline plot in Fig. 7, shows a helical flow pattern through the anastomotic region, at the point where fluid from the 2 independent inlets mixes.

Fig. 7.

Streamline plots generated from the CFD simulation. Note: The side branch is not visible on these images as only 2% of the total flow was diverted into this vein. The geometry used to generate these plots was identical to that represented in Figs. 1,3,4. CFD, computational fluid dynamics.

Beyond the anastomosis USTIR in the vein corresponds well with both the WSS indices and the observed structural remodelling and NIH development.

The extent of the inward luminal remodelling in this particular case meant that localised regions of elevated OSI could not be directly correlated to regions of smooth muscle cell proliferation. Future studies of this type would benefit from more regular surveillance of the vessel walls to better establish the relationship between elevated OSI and the onset of NIH formation. Although the ethical considerations surrounding most multi-planar image acquisition are a limiting factor in such work, ultrasound imaging at regular intervals may provide a degree of insight into the natural history of NIH progression.

Although care was taken to make the CFD model as robust as possible, a number of assumptions and simplifications were made in order to balance physiological accuracy with computational expense, including rigid vessel walls and laminar flow conditions. It is also important to recognise that the haemodynamic scenario modelled in our simulation corresponds to resting flow conditions at a single time point. During dialysis the position and angle of the needles, along with the pump speed will have a significant impact on the flow through the circuit, resulting in a very different set of haemodynamic conditions, with localised regions of both high and low WSS combined with oscillatory flow patterns [20]. Further to this, flow will likely vary in the periods between dialysis sessions due to filtration induced fluctuations in viscosity, and blood pressure [21].