Blood flow velocity analysis in cerebral perforating arteries on 7T 2D phase contrast MRI with an open-source software tool (SELMA)

doi:10.21203/rs.3.rs-5045336/v1

Blood flow velocity in the cerebral perforating arteries can be quantified in a two-dimensional plane with phase contrast magnetic imaging (2D PC-MRI). The velocity pulsatility index (PI) can inform on the stiffness of these perforating arteries, which is related to several cerebrovascular diseases. Currently, there is no open-source analysis tool for 2D PC-MRI data from these small vessels, impeding the usage of these measurements. In this study we present the Small vessEL MArker (SELMA) analysis software as a novel, user-friendly, open-source tool for velocity analysis in cerebral perforating arteries. The implementation of the analysis algorithm in SELMA was validated against previously published data with a Bland-Altman analysis. The inter-rater reliability of SELMA was assessed on PC-MRI data of sixty participants from three MRI vendors between eight different sites. The mean velocity (v_mean) and velocity PI of SELMA was very similar to the original results (v_mean: mean difference ± standard deviation: 0.1 ± 0.8 cm/s; velocity PI: mean difference ± standard deviation: 0.01 ± 0.1) despite the slightly higher number of detected vessels in SELMA (N_detected: mean difference ± standard deviation: 4 ± 9 vessels), which can be explained by the vessel selection paradigm of SELMA. The Dice Similarity Coefficient of drawn regions of interest between two operators using SELMA was 0.91 (range 0.69–0.95) and the overall intra-class coefficient for N_detected, v_mean, and velocity PI were 0.92, 0.84, and 0.85, respectively. The differences in the outcome measures was higher between sites than vendors, indicating the challenges in harmonizing the 2D PC-MRI sequence even across sites with the same vendor. We show that SELMA is a consistent and user-friendly analysis tool for small cerebral vessels.

Cerebral blood flow in large intracranial arteries can be quantified using phase contrast magnetic resonance imaging (PC-MRI). This technique estimates the blood flow velocity from the MR signal's phase accrual of moving spins as the velocity of the moving spins (i.e. blood) is proportional to this phase accrual. Recently, with the advantages of the increased sensitivity on high-field 7T MRI, granting higher spatial resolutions with adequate SNR, blood flow velocity of the brain's perforating arteries through a two-dimensional plane can be measured using 2D PC-MRI^{1, 2}.

Blood flow velocity measurements can inform on the condition of these arteries, notably also in disease states^3–6. In particular, the blood flow velocity's pulsatility index (PI) is an important indicator of arterial stiffness^7–9. Increased stiffness in the large intracranial arteries has been associated with damage to the brain parenchyma, such as microbleeds^{10, 11}, lacunar infarcts^{10, 12}, and white matter hyperintensities^{10, 13, 14}, while also being linked to cognitive impairment^{7, 15} and cognitive decline^{15, 16}. In the smaller arteries, information from these blood flow velocity measurements could be especially relevant in cerebral small vessel diseases (cSVDs), an important cause of stroke and dementia¹⁷. Recently, using 7T MRI, an increased velocity PI in the perforating arteries was found in patients with cSVDs and a lower blood flow velocity was linked to white matter damage increase/progression^{18, 19}. Age was also found to be a determinant of increased velocity PI in the perforating arteries of the basal ganglia in both cSVD and general population cohorts^{5, 20, 21}.

For the larger intracranial arteries, several publicly available image processing tools exist to analyze PC-MRI data^22–25. For the smaller arteries, however, there are no processing tools available yet. Reported analyses of PC-MRI data of the smaller arteries were performed with a collection of in-house developed code^{2, 26, 27}, which has several limitations. First, such code consists of scripts in which the parameters of the algorithms are hard-coded and not accessible as user-defined settings. Second, these scripts, including parameter settings, often evolve with time without explicit version tracking. Also, as part of the developments of the analysis methods, several variants of output parameters were computed in parallel, such as the PI from the mean normalized velocity traces from all arteries and the mean over the PIs of the individual velocity traces. This situation can potentially create ambiguities and unintentionally result in inconsistent application of output parameters in studies. Additionally, it impedes the reproducibility of previous analyses, even more so across centers. Finally, the lack of version-controlled software with a user-friendly user interface hampers the dissemination of the methods to other researchers, which limits the usage of blood flow velocity and PI measurements in the perforating arteries of the basal ganglia and semioval center on a larger scale. To support consistent, repeatable analysis and to promote the usage and possibly further development of the PC-MRI measurements in small arteries, an analysis tool should be established that is accessible (open-source), user-friendly (with a graphical user interface (GUI)), and facilitate repeatable analysis (logging of version and settings with results) of small artery PC-MRI data.

This work presents the Small vessEL MArker (SELMA) analysis software as a novel, open-source tool for cerebral small artery flow velocity analysis compatible with data from multiple MRI vendors. SELMA incorporates converged and updated analysis algorithms used in our previous publications^{2, 26, 28}. To validate the implementation of the algorithm in SELMA, we re-analyzed previously published data²⁹ and compared the results of SELMA with the original results obtained with a previously published iteration of the algorithm. Additionally, we assessed the inter-rater reliability of SELMA on data from three different 7T MRI vendors by two trained operators. A secondary aim was to assess the performance of measuring perforating artery velocity pulsatility at MRI scanners from different sites and vendors by comparing measurements between sites and testing for age and sex effects.

Algorithm

The original algorithm used to analyze small artery PC-MRI data as described in previous work^{2, 26, 28} was developed in MATLAB (version R2015b, the MathWorks, Natick, MA). These MATLAB scripts have been rewritten in Python 3.7 to develop SELMA for analysis of the cerebral perforating arteries at the level of the basal ganglia or the white matter at the semioval center (Fig. 1). SELMA was developed to be compatible with all 2D PC-MRI DICOM data, regardless of MRI vendor or field strength. The implementation of the algorithm in SELMA follows that of the original codes^{2, 26, 28}, but includes several updates and improvements, which will be detailed below.

Analysis in SELMA can be started by drawing a region of interest (ROI) on the 2D PC-MRI images using the software interface or loading such mask from the disk. The first step of the analysis is the estimation of the noise level in the phase and magnitude frames, which is computed on a voxel-wise basis using the standard deviation of the complex signal (constructed from the phase and magnitude data) over the cardiac cycle. The root mean square of the standard deviations of the imaginary and complex parts of the constructed complex signal are computed next. A median filter with a user-defined kernel size (default 10 mm) is then applied to the root mean square maps. The temporal signal-to-noise ratio of the magnitude frames (SNR_mag) is computed by dividing the magnitude frames with the median filtered root mean square maps. The phase maps are scaled to velocity maps with the pre-defined velocity encoding (v_enc) of the scan protocol parameters available from the DICOM header. The velocity frames are averaged over the cardiac cycle and median-filtered with a 10 mm kernel. The median filtered velocity map is subtracted from the raw velocity maps at each point in the cardiac cycle to remove the background phase of static tissue and to center the velocity around zero. Next, using SNR_mag, the standard deviation of the corrected velocity maps (σ_v) can be computed using the following formula³⁰:

$$\:{\sigma\:}_{v}=\frac{{v}_{enc}}{\pi\:}\frac{1}{{SNR}_{mag}}$$

The corrected velocity maps are divided by σ_v to obtain the temporal signal-to-noise ratio of the velocity frames (SNR_v). Based on the SNR_mag and SNR_v, voxels in the 2D PC-MRI image can be clustered into several profiles (Table 1). For the detection of vessels in the basal ganglia, clusters of voxels with an SNR_mag and SNR_v above a positive noise threshold (T_n) are detected as potential arteries. Clusters of voxels in the semioval center are identified as arteries only if SNR_v is lower than the negative T_n, regardless of SNR_mag. The sign for SNR_v in the semioval center is changed due to the opposite directionality of the blood flow of the perforating arteries at the semioval center compared to the basal ganglia. T_n can be changed by the user by changing the significance level in the settings (default 0.05, in which case Tn = 1.96). A new addition in SELMA is the scaling of T_n based on the number of acquired heart phases and the temporal resolution of the data. Estimation of the SNR in the magnitude and velocity frames are dependent on the heart rate, which affects the number of acquired heart phases, and on the number of k-space lines acquired per heartbeat and TR of the acquisition, which affects the temporal resolution.

After identification of the arteries that exceed T_n, operators of SELMA can select the options to automatically filter out arteries whose orientation is far from being perpendicular to the scanning plane and/or to deduplicate arteries that are too close to each other. The blood flow velocity in arteries that are not perpendicular to the scanning plane will be underestimated². SELMA considers the artery's roundness to determine whether the artery can be considered as perpendicular to the scanning plane, at least in first approximation, in which case it should be included in the analysis. An ellipse is fitted to the artery's circumference, and the axes ratio of the ellipse is computed, which discards the artery from further analysis if the axes ratio exceeds a user-set threshold (default 2). Artery detections that are too close to each other might be multiple detections on a single artery oriented parallel to the scanning plane or might be caused by ghosting artefacts. In such cases, SELMA discards all arteries except the one with the highest velocity if detected within a user-set distance (default 1.2 mm) from each other. These options can be applied both to basal ganglia and semioval center focussed PC-MRI scans. Additional options can be selected to erode the outer region of the ROI and/or filter out ghosting artefacts from the scans. Arteries in the outer edges of the ROI in the semioval center are more prone to motion artifacts. Users can specify how many voxels from the edges of the ROI can be eroded (default 80). Ghosting of larger arteries in the phase encoding direction can lead to erroneous artery detections. SELMA incorporates the automatic ghosting censoring method as previously described²⁶. First, SELMA identifies the large blood arteries by applying a relative intensity threshold on the magnitude images. Only clusters of voxels in the top magnitude percentile defined by the user (default 0.3%) and larger than a minimum size are included as potential large arteries. Depending on the size of the bright voxel cluster, an exclusion zone will be drawn around the cluster. Undesired detections in these zones will be discarded from further analysis. Ghosting artefacts are a larger issue in the semioval center when automatically segmented ROIs are used compared to the basal ganglia where they can be avoided when manually delineating a ROI. In the options, users can define custom values for the thresholds, voxel sizes of the arteries or length and width of the exclusion zones.

Another addition to SELMA is the ability for operators to censor artery detections in the basal ganglia manually. Selecting this option will override the in-plane artery censoring and deduplication steps. SELMA will visualize all detected arteries and guide the user systematically through every artery prompting the user if it should be discarded for analysis. The user can see the proximity of the artery to others in the user interface, and the axes ratio of the artery is provided to aid the user in this process. The user can opt to keep or discard the highlighted artery, and is then guided to the next one. After all detected arteries have been evaluated, SELMA prompts the user to either save their results, or to restart the manual vessel censoring process again from the beginning. All parts of the manual vessel censoring functionality can be directly performed in the GUI.

After final artery selection, the number of perforating arteries (N_detected), their mean blood flow velocity (v_mean) and the velocity PI are assessed. The PI is generally defined as $\:\frac{{v}_{max}-{v}_{min}}{{v}_{mean}}$, where v_max, v_min, and v_mean, are the maximum, minimum, and mean of the velocity trace over the cardiac cycle in cm/s. We calculated the PI from the averaged velocity trace of all final included arteries, while the individual velocity traces were normalized before averaging (hence, v_mean in the PI formula was 1.0 by definition). During the development of SELMA, an approach using the median of the normalized velocity trace was also considered, however during internal testing, we found that the blood flow velocity distribution over all arteries for each time point in the cardiac cycle followed a normal distribution and that using a median estimator had a higher uncertainty than the mean estimator, which led to inflated values of the PI¹⁹. Thus, we opted to use the mean of the normalized velocity trace to compute the PI. During internal testing, we observed similar group differences in blood flow velocity and velocity PI using the mean versus median normalized velocity trace.

User interface

Analysis of small vessel 2D PC-MRI data has been made more user friendly with the addition of a GUI in SELMA (Fig. 1). 2D PC-MRI scans of either the basal ganglia or the semioval center in DICOM format can be loaded into SELMA via the menu in the top bar. ROIs that were either manually drawn in SELMA or segmented with external software, such as white matter masks, can also be loaded from the disk with this menu. Analysis can be started with the 'Analyse' option in the menu for single scans or, if masks are already provided, entire folders containing scans. The batch analysis option in SELMA allows for automatic analysis of multiple scans without any user input, drastically speeding up the analysis process. Two default voxel clustering settings have been defined, tailored to acquisitions at the basal ganglia and semioval center, respectively. Before the analysis, the actual slice location must be selected in the bottom right-hand corner to ensure that the corresponding predefined voxel clustering settings are applied in SELMA. Custom clustering options can also be selected in this menu for analysis of data that do not fall into the standard profiles for basal ganglia or semioval center scans. The remaining options in the top menu allow for changing the contrast or zoom of the image, which can also be directly manipulated in the user interface using the mouse. In the settings menu, the operator can change several options and thresholds as described in the 'Algorithm' section.

SELMA was developed specifically with the goal of an easy to maintain code base. To this end, the GUI and the processing parts of the code are intentionally separated in different classes that can only talk to each other using the Signal-and-Slot paradigm. This should enable future developers to easily understand the code, and make small changes to either the processing algorithms or the GUI without accidentally introducing unwanted behavior elsewhere in the program.

Test data

The implementation of the analysis algorithm in SELMA was validated with 7T MRI (Achieva 7T, Philips Medical Systems, Best, The Netherlands) data from previous work²⁹. These 2D PC-MRI data were acquired with a 32-channel receive head coil (Nova Medical, Wilmington, NC, United States) in 14 patients with coarctation of the aorta and 15 control subjects with no history of cardiovascular disease, neurological disease, or intellectual disability. The 2D PC sequence was planned on a 3D T1-weighted image at the level of the basal ganglia (Fig. 2) and was retrospectively gated using a peripheral pulse oximeter for triggering. The following scan parameters were used: 250 × 250 mm² field of view; acquired spatial resolution 0.3 × 0.3 × 2.0 mm³, reconstructed spatial resolution (through zero-filling in k-space) 0.2 × 0.2 × 2.0 mm³; TR/TE = 28/14.7–15.1 ms; flip angle = 50°; TFE factor = 2 (i.e. number of k-lines acquired per cardiac frame, per heart beat); velocity encoding = 20 cm/s; acquired temporal resolution = 112 ms; 13–15 reconstructed heart phases, depending on the heart rate; SENSE factor = 1.5; scan duration was about 5 minutes for a heart rate of 60 beats/min.

Inter-rater reliability was assessed with data from 8 different scanning sites comprising three different MRI vendors: Philips (Achieva 7T, Philips Healthcare, Best, The Netherlands), Siemens (MAGNETOM 7T, Siemens Healthineers, Erlangen, Germany) and General Electric (GE) (Discovery MR950, GE Healthcare, Chicago, Illinois, USA). All MRI data were acquired with a 32-channel receive head coil (Nova Medical, Wilmington, NC, United States), except for one site where a 24-channel receive head coil (Nova Medical, Wilmington, NC, United States) was used. All scanning sites participated in the European Ultrahigh-Field Imaging Network for Neurodegenerative Diseases (EUFIND)³¹. EUFIND comprises researchers from 22 7T MRI sites in Europe with the common aim of identifying opportunities and challenges of 7T MRI for clinical and research applications in neurodegeneration. Each site included healthy elderly subjects aged between 50 and 70 years with an equal gender distribution and attempted to include an equal number of volunteers between 50 and 60 years and between 60 and 70 years to allow for secondary analysis of age and gender effects on the outcome parameters. All institutions performed the scanning under the approval of the institution's local ethical review board, and all subjects provided written informed consent.

The 2D PC-MRI sequence was originally developed for the Philips 7T MRI and the protocol was harmonized across all sites with the parameters in Table 2. Harmonization was limited by the vendor-specific implementation of the 2D PC sequence and parameter ranges. The acquisition slice was planned in the basal ganglia, targeting the perforating lenticulostriate arteries that branch from the circle of Willis. Representative basal ganglia scans per vendor are shown in Fig. 2.

Data analysis

The differences in the outcome measures between SELMA and the previously published results on the validation data were quantified on a group level with a Bland-Altman analysis. Inter-rater reliability of the ROIs drawn in SELMA was assessed using the Dice similarity coefficient (DSC) between ROIs drawn by two trained operators. Inter-rater reliability of the outcome measures using manual artery selection was assessed by calculating the Intra-class Correlation Coefficient (ICC) between two operators. An ICC above 0.75 was considered to be an excellent correlation between operators³². The coefficient of variation (CV) of all outcome measures and the SNR_v of the included arteries averaged over all subjects within each site were used to assess inter-site differences. The CV for all sites within a single vendor was computed and subsequently averaged over all vendors to assess the inter-vendor differences. The relation of age and sex with the outcome measures was assessed using a linear mixed model corrected with site and vendor as random effects. All statistical analyses were performed in MATLAB (version R2021a, the MathWorks, Natick, MA) except for the linear mixed modelling, which was performed in R (version 4.2.1, R Foundation for Statistical Computing, Vienna, Austria) with the 'lme4' package.

Validation of the implementation

Table 3 shows the results of the original measurements compared to SELMA on previously published data of 29 participants²⁴. In the same data, SELMA had a slightly higher N_detected compared to the previous version of the algorithm. The v_mean and velocity PI were similar between both approaches at a group level. Bland-Altman plots for the comparison between the original measurements and SELMA are shown in Figure 3.

Multi-vendor comparison

Data of 60 participants (mean age ± standard deviation: 59±6 years) were included (eight sites, three MRI vendors) for the multi-vendor analysis. Due to vendor-specific implementation of the 2D PC sequence, scan parameters slightly varied between subject, site and vendor. Table 4 provides an overview of the relevant parameters. The N_detected, v_mean, PI, and SNR_v are stratified per vendor in Table 5 and per site in Table 6. The CV for inter-site differences were larger for all measurements than the inter-vendor differences (Table 7).

Inter-rater analysis

The mean DSC for the manually drawn ROIs by the two operators was 0.91 (range 0.69-0.95) across all 60 participants (Table 4). All ICCs for the measurements exceeded 0.75 except for the velocity PI of vendor 1 and N_detectedof vendor 2. The ICC of N_detected, v_mean, and velocity PI for the entire dataset was 0.92, 0.84, and 0.85, respectively (Table 5).

Age and gender analyses

Age and gender were not significantly associated with N_detected, v_mean, and velocity PI over the entire population in multivariate models adjusted for site and vendor (Figure 4).

In this study, we developed SELMA, an open-source small cerebral artery analysis tool for 2D PC-MRI data. We have shown that the results using SELMA match well with the results of the original version of our algorithm using published data²⁹. Reproducible multi-vendor analysis of 2D PC-MRI data was also possible with SELMA.

With SELMA, we obtained a higher N_detected than in the previous analysis, which can be explained by the changes that were made in the deduplication steps. Previously, all perforating arteries within a 1.2 mm distance from each other were discarded in the analysis, which led to more false negatives during detection. In SELMA, only the artery with the highest v_mean is kept, and all other detections within a 1.2 mm distance are discarded. The increase in N_detected had no effect on the average v_mean and velocity PI, as these were similar to the original results over the entire group. This supports the original assumption that averaging over all available small arteries gives a representative number for the detected population of arteries². The Limits of Agreement (LoA) for the velocity PI found in the Bland-Altman analysis was similar to the LoA on a group-level for the original measurements on 7T (LoA original versus SELMA: 0.20 vs 0.16)²⁴, suggesting that the variability in velocity PI can mostly be explained by thermal or physiological noise during the measurement rather than a bias introduced by different analysis methods.

On multi-vendor data, analysis with SELMA showed that the CV in the outcome measures was higher between sites than between vendors. The CV for N_detected was higher than the other outcome measures in both comparisons. The high variance in N_detected can be explained by various factors causing the sensitivity differences for detecting arteries across sites. Differences in the reconstruction software of the scanners could likely partially explain this variation. The variation of N_detected between sites of the same vendor might reflect variation in patient handling by the local MRI technician, such as the usage of padding to limit motion, which was not standardized and may have varied between sites. Also, despite slice planning was standardized via an illustrated instruction, variations between technicians in the execution of this planning cannot be prevented and could effect the number of visible vessels. Another source of variation could be that the protocol could slightly vary between sites as gradient performance or the availability of (commercial) options in the software could vary and the protocol at each site was implemented as close to the requested parameters as possible. All these factors could have contributed to the variance of N_detected and these differences in sensitivity illustrate the challenges in harmonizing the 2D PC-MRI sequence across multiple sites and vendors, as shown in the differences in scan parameters between sites and vendors (Table 3) and SNR_v differences (Table 4).

In an early attempt by the EUFIND group, the 2D PC-MRI sequence was also attempted to harmonize for the perforating arteries in the semioval center. However, large sensitivity differences between vendors made the implementation of this sequence in these small (≤ 300µm) perforating arteries difficult³¹. These sensitivity differences might also affect the CVs of the outcome measures of the perforating arteries in the basal ganglia. The v_mean and velocity PI over the entire dataset were slightly higher than reported for control groups in other studies^{18, 20, 29}. This can likely be attributed to the high variance in N_detected due to the aforementioned sensitivity differences across sites. Smaller arteries with a blood flow velocity just above the noise threshold might not be detected in several sites, which could lead to an overestimation of v_mean and velocity PI in these participants. A similar effect has been observed when measuring these perforating arteries at 3T MRI²⁹.

The mean DSC for the manually drawn ROIs on all scans across all sites and vendors (mean: 0.91; range: 0.69–0.95) was comparable to earlier results on 3T 2D PC-MRI data which was single-vendor only (mean: 0.90; range: 0.85–0.95)³³. The option to manually censor the results of the automatic artery detection in the basal ganglia was added after the observation that the automatic detection could not always successfully remove in-plane arteries (i.e. arteries not perpendicular to the imaging plane). The lenticulostriate arteries in the basal ganglia are quite tortuous³⁴, which can be a challenge for the automatic in-plane artery removal and deduplication steps in SELMA. The arteries in the semioval center are smaller than in the basal ganglia and are censored better with the automatic censoring steps. The operator can freely choose which arteries to include or exclude with manual censoring. The ICCs after manual artery censoring exceeded 0.75 for all outcome measurements over the entire dataset, indicating that despite the differences in acquisition parameters in the entire dataset, manual censoring yields good to excellent inter-rater reliability. The overall ICC for the velocity PI was higher than the reported ICCs for test-retest reliability by Schnerr et al²¹., suggesting that the SELMA results are probably less dependent on the operator (manual ROI delineation and vessel censoring) than on the measurement noise (scan-to-scan variability).

We found no relation between age or gender and the outcome measures when corrected for sites and vendors using a linear mixed model. In literature, velocity PI of the lenticulostriate arteries has been found to be positively associated with age^{5, 20, 21}, and males had a higher velocity PI in the internal carotid arteries³⁵ but not in the lenticulostriate arteries²⁰ compared to females. The lack of an age effect on the velocity PI in our dataset could be explained by both the differences in acquisition parameters and sensitivity between sites and vendors and the relatively small age range of the study participants in our study (mean age: 59; range: 50–70; standard deviation: 6 years) when compared to the other studies in the literature^{5, 20, 21}. In an additional analysis, we found velocity PI to be associated with N_detected in our dataset (Supplementary Table 1). As a result of the differences in acquisition between sites, velocity PI seems to be dependent on N_detected, which could also explain the lack of an age effect on the velocity PI. The association between N_detected and PI was negative, which suggests that a lower N_detected means that only the relatively larger vessels (with accordingly higher PI) are detected.

Analysis of cerebral perforating arteries of the basal ganglia and semioval center on 2D PC-MRI data has been performed with several other methods using similar algorithmic steps that have been incorporated in SELMA. These methods are tailored towards the analysis of one single artery at a time, require multiple voxels across the diameter of an artery, or require the need for manual delineation of the perforating arteries^36–38. One method exists that uses convolutional neural networks to automatically segment perforating arteries³⁹. SELMA has been developed specifically for the automatic analysis of small arteries with subvoxel size and pooling those to obtain metrics from them. Blood flow velocity in the smaller distal branches of the lenticulostrate arteries included in SELMA is more attenuated, which could explain the lower observed velocities (4.0-4.7 cm/s versus 7.0 cm/s³⁸ and 11.8 cm/s³⁶) and PI (0.27–0.49 versus 0.79³⁸) compared to those obtained with methods that (manually) delineated a larger singular branch of the lenticulostriate arteries.

Arterial pulsatility is an important functional measure in neurodegenerative diseases such as cSVDs and can be quantified in a 2D plane, such as done in this study, or with a 4D flow method. The latter approach has been used in the larger intracranial arteries in cSVDs where associations between pulsatility and SVD markers were found^{4, 40}. In the perforating cerebral arteries, 4D flow is more difficult to implement due to the small vessel diameters, thus blood flow velocity in these smaller arteries are computed in a 2D plane and pulsatility is derived from the blood flow velocity alone. Especially in cSVDs, we found that the pulsatility of these small vessels are important markers for disease burden^{18, 19}.

SELMA and its source code were made publically available and can be installed from: https://github.com/TNI-UMCU/SELMA/ and is provided with step-by-step installation instructions, including a user manual. The tool is supported on Windows and Linux.

This study does have some limitations that have to be addressed. First, the results of our multi-vendor comparison and inter-rater reliability analyses are potentially affected by the differences in hardware, acquisition parameters and reconstruction settings across all sites and vendors. It is unclear how much this could have affected the manual ROI delineation in the basal ganglia or the manual artery censoring. However, these differences are assumed to affect the inter-site and inter-vendor differences in the outcome measures more than the inter-rater reliability. Future work on harmonization of the 2D PC-MRI sequence across multiple vendors, by e.g. using a vendor neutral MR pulse sequence program^{41, 42}, should address these differences to further stimulate the usage of these measurements of hemodynamic function in the perforating arteries in multicenter studies. Second, the SELMA approach to include subvoxel arteries in size might lead to underestimating the blood flow velocity due to partial volume effects². It is assumed that the effect of underestimation is comparable between subjects and patients, thus still allowing the study of disease or physiology. Third, the relatively low number of included participants per site could also have affected the inter-site and inter-vendor comparisons for all outcome measures. We assume, however, that the contribution of the study population and sample size was small compared to the contribution of the sensitivity differences with regard to the variability in the outcome measures we found.

We present SELMA, a novel small cerebral artery blood flow velocity analysis tool of 2D PC-MRI data. We achieved good inter-rater reliability of the analysis on data acquired from MRI systems with different hardware, vendor, and scan parameters. Consistent and user-friendly analysis of small cerebral arteries is possible using SELMA. Differences in the implementation of 2D PC-MRI across vendors currently hampers multicenter studies on the hemodynamic function of the perforating arteries in the brain. Still, the tool can be used for single-center studies regardless of vendor or site in which 2D PC-MRI is performed in both patients and controls with the same system and protocol.

Author Contribution

S.D.T.P.: Conceptualization; Tool development; Formal analysis; Writing original draft; Final approval.C.C.: Tool development; Writing review & editing; Final approval.J.T.v.V.: Formal analysis; Writing review & editing; Final approval.R.J.v.T.: Formal analysis; Writing review & editing; Final approval.M.B., M.C., L.d.R., O.K., M.E.L., K.P., I.R., J.C.W.S., M.T., A.V.,: Data acquisition, Writing review & editing; Final approval.G.J.B.: Writing review & editing; Final approval.J.J.Z.: Conceptualization; Data acquisition; Formal analysis; Writing review & editing; Final approval.

Acknowledgement

This work received funding from the European Research Council under the European Union’s Seventh Framework Program (FP7/2007-2013)/ERC grant agreements nr. 337333 (SmallVesselMRI) and nr. 841865 (SELMA)S.D.T.P. and J.C.W.S. are funded by the UMCU Brain Center Young Talent Fellowship 2019.

Data Availability

The data that support the findings of this study are available, upon reasonable request, from the corresponding author.

Zwanenburg JJM, van Osch MJP. Targeting cerebral small vessel disease with mri. Stroke. 2017;48:3175-3182
Bouvy WH, Geurts LJ, Kuijf HJ, Luijten PR, Kappelle LJ, Biessels GJ, et al. Assessment of blood flow velocity and pulsatility in cerebral perforating arteries with 7-t quantitative flow mri. NMR Biomed. 2016;29:1295-1304
Mitchell GF, Van Buchem MA, Sigurdsson S, Gotal JD, Jonsdottir MK, Kjartansson Ó, et al. Arterial stiffness, pressure and flow pulsatility and brain structure and function: The age, gene/environment susceptibility – reykjavik study. Brain. 2011;134:3398-3407
O'Rourke MF, Hashimoto J. Mechanical factors in arterial aging: A clinical perspective. J Am Coll Cardiol. 2007;50:1-13
Webb AJS, Simoni M, Mazzucco S, Kuker W, Schulz U, Rothwell PM. Increased cerebral arterial pulsatility in patients with leukoaraiosis. Stroke. 2012;43:2631-2636
Van Sloten TT, Protogerou AD, Henry RMA, Schram MT, Launer LJ, Stehouwer CDA. Association between arterial stiffness, cerebral small vessel disease and cognitive impairment: A systematic review and meta-analysis. Neuroscience & Biobehavioral Reviews. 2015;53:121-130
Zhai F-F, Ye Y-C, Chen S-Y, Ding F-M, Han F, Yang X-L, et al. Arterial stiffness and cerebral small vessel disease. Frontiers in Neurology. 2018;9
Chuang S-Y, Cheng H-M, Bai C-H, Yeh W-T, Chen J-R, Pan W-H. Blood pressure, carotid flow pulsatility, and the risk of stroke. Stroke. 2016;47:2262-2268
Aribisala BS, Morris Z, Eadie E, Thomas A, Gow A, Valdés Hernández MC, et al. Blood pressure, internal carotid artery flow parameters, and age-related white matter hyperintensities. Hypertension. 2014;63:1011-1018
Poels MMF, Zaccai K, Verwoert GC, Vernooij MW, Hofman A, Van Der Lugt A, et al. Arterial stiffness and cerebral small vessel disease. Stroke. 2012;43:2637-2642
Singer J, Trollor JN, Baune BT, Sachdev PS, Smith E. Arterial stiffness, the brain and cognition: A systematic review. Ageing Res Rev. 2014;15:16-27
Zeki Al Hazzouri A, Newman AB, Simonsick E, Sink KM, Sutton Tyrrell K, Watson N, et al. Pulse wave velocity and cognitive decline in elders. Stroke. 2013;44:388-393
Wardlaw JM, Smith C, Dichgans M. Small vessel disease: Mechanisms and clinical implications. Lancet Neurol. 2019;18:684-696
van den Brink H, Kopczak A, Arts T, Onkenhout L, Siero JCW, Zwanenburg JJM, et al. Cadasil affects multiple aspects of cerebral small vessel function on 7t-mri. Ann Neurol. 2023;93:29-39
Van Den Brink H, Pham S, Siero JC, Arts T, Onkenhout L, Kuijf H, et al. Assessment of small vessel function using 7t mri in patients with sporadic cerebral small vessel disease. Neurology. 2024;102
Perosa V, Arts T, Assmann A, Mattern H, Speck O, Oltmer J, et al. Pulsatility index in the basal ganglia arteries increases with age in elderly with and without cerebral small vessel disease. American Journal of Neuroradiology. 2022;43:540-546
Vikner T, Nyberg L, Holmgren M, Malm J, Eklund A, Wahlin A. Characterizing pulsatility in distal cerebral arteries using 4d flow mri. J Cereb Blood Flow Metab. 2020;40:2429-2440
Schnerr RS, Jansen JFA, Uludag K, Hofman PAM, Wildberger JE, Van Oostenbrugge RJ, et al. Pulsatility of lenticulostriate arteries assessed by 7 tesla flow mri—measurement, reproducibility, and applicability to aging effect. Frontiers in Physiology. 2017;8
van Hespen KM, Kuijf HJ, Hendrikse J, Luijten PR, Zwanenburg JJM. Blood flow velocity pulsatility and arterial diameter pulsatility measurements of the intracranial arteries using 4d pc-mri. Neuroinformatics. 2022;20:317-326
Köhler B, Born S, Van Pelt RFP, Hennemuth A, Preim U, Preim B. A survey of cardiac 4d pc‐mri data processing. Computer Graphics Forum. 2017;36:5-35
Köhler B, Grothoff M, Gutberlet M, Preim B. Bloodline: A system for the guided analysis of cardiac 4d pc-mri data. Computers & Graphics. 2019;82:32-43
Arts T, Siero JCW, Biessels GJ, Zwanenburg JJM. Automated assessment of cerebral arterial perforator function on 7t mri. J Magn Reson Imaging. 2021;53:234-241
Van Den Kerkhof M, Van Der Thiel MM, Van Oostenbrugge RJ, Postma AA, Kroon AA, Backes WH, et al. Impaired damping of cerebral blood flow velocity pulsatility is associated with the number of perivascular spaces as measured with 7t mri. Journal of Cerebral Blood Flow & Metabolism. 2023;43:937-946
Geurts L, Biessels GJ, Luijten P, Zwanenburg J. Better and faster velocity pulsatility assessment in cerebral white matter perforating arteries with 7t quantitative flow mri through improved slice profile, acquisition scheme, and postprocessing. Magn Reson Med. 2018;79:1473-1482
Arts T, Meijs TA, Grotenhuis H, Voskuil M, Siero J, Biessels GJ, et al. Velocity and pulsatility measures in the perforating arteries of the basal ganglia at 3t mri in reference to 7t mri. Front Neurosci. 2021;15:665480
Düzel E, Acosta‐Cabronero J, Berron D, Biessels GJ, Björkman‐Burtscher I, Bottlaender M, et al. European ultrahigh‐field imaging network for neurodegenerative diseases (eufind). Alzheimer's & Dementia: Diagnosis, Assessment & Disease Monitoring. 2019;11:538-549
Cicchetti DV. Guidelines, criteria, and rules of thumb for evaluating normed and standardized assessment instruments in psychology. Psychological Assessment. 1994;6:284-290
Van Tuijl RJ, Pham SDT, Ruigrok YM, Biessels GJ, Velthuis BK, Zwanenburg JJM. Reliability of velocity pulsatility in small vessels on 3tesla mri in the basal ganglia: A test–retest study. Magnetic Resonance Materials in Physics, Biology and Medicine. 2022;36:15-23
Ma SJ, Sarabi MS, Yan L, Shao X, Chen Y, Yang Q, et al. Characterization of lenticulostriate arteries with high resolution black-blood t1-weighted turbo spin echo with variable flip angles at 3 and 7 tesla. NeuroImage. 2019;199:184-193
van Tuijl RJ, Ruigrok YM, Velthuis BK, van der Schaaf IC, Rinkel GJE, Zwanenburg JJM. Velocity pulsatility and arterial distensibility along the internal carotid artery. J Am Heart Assoc. 2020;9:e016883
Chengyue S, Yue W, Chen L, Zhiying X, Yunchuang S, Zhihao X, et al. Reduced blood flow velocity in lenticulostriate arteries of patients with cadasil assessed by pc-mra at 7t. Journal of Neurology, Neurosurgery &amp; Psychiatry. 2022;93:451
Zong X, Lin W. Quantitative phase contrast mri of penetrating arteries in centrum semiovale at 7t. NeuroImage. 2019;195:463-474
van den Kerkhof M, Jansen JFA, van Oostenbrugge RJ, Backes WH. 1d versus 3d blood flow velocity and pulsatility measurements of lenticulostriate arteries at 7t mri. Magnetic Resonance Imaging. 2023;96:144-150
Moore J, Jimenez J, Lin W, Powers W, Zong X. Prospective motion correction and automatic segmentation of penetrating arteries in phase contrast <scp>mri</scp> at 7 t. Magnetic Resonance in Medicine. 2022;88:2088-2100

Table 1: Standard clustering profiles for basal ganglia and semioval center 2D PC-MRI scans in SELMA

	SNR_vel < - T_n	- T_n < SNR_vel < T_n	SNR_vel > T_n
Basal ganglia
SNR_mag > T_n	-	-	+
SNR_mag < T_n	-	-	-
Centrum semioval
SNR_mag > T_n	+	-	-
SNR_mag < T_n	+	-	-

SNR_mag= magnitude signal-to-noise ratio; SNR_vel = velocity signal-to-noise ratio; T_n= noise threshold

Table 2: Desired scan parameters for the two dimensional phase contrast sequence for all three MRI vendors

FOV RL x AP (mm)	250x250
Foldover direction	AP
Acquired matrix	832x833
Reconstructed voxel size (mm)	0.2 x 0.2
Parallel imaging	none
Slice thickness (mm)	2.0
TFE factor (nr. phase enc. lines per cardiac time point)	2
TR/TE (ms)	27/16
Flip angle	60°
Cardiac synch.	Retrospective triggering
Reconstructed nr. of heart phases	14
Acquired temporal resolution	107 ms
Venc (cm/s)	20
Scan duration	4 minutes with heart rate of 60 bpm

Table 3: Results of the groupwise comparison between SELMA and the originally published results.

	Original	SELMA
N_detected	24 ± 6	28 ± 6
v_mean	3.9 ± 0.6	4.0 ± 0.5
PI	0.28 ± 0.08	0.27 ± 0.07

Values are given as mean ± SD. N_detected = amount of detected arteries; v_mean = mean blood flow velocity of the perforating arteries given in cm/s; PI = pulsatility index

Table 4: Scan parameters (mean ± SD) for the two-dimensional phase-contrast sequence for all three MRI vendors.

Scan parameter	Vendor 1	Vendor 2	Vendor 3
Head coil	32-channel receiver Circular polarized transmit coil	32-channel receiver Circular polarized transmit coil	24-channel / 32-channel receiver
FOV, mm	242 ± 4 x 242 ± 4	250 x 250	154 ± 0.4 x 154 ± 0.4
Reconstructed voxel size, mm	0.23 x 0.23	0.2 x 0.2	0.3 x 0.3
Flip angle, °	60	60	54-60
TR/TE, ms	32.4 ± 1.8 / 11.3 ± 0.5	26.3 ± 0.4 / 15.9 ± 0.4	28.1 ± 1.1 / 17.2 ± 0.4
Bandwidth, Hz/pixel	34.3 ± 1.6	59	60.8 ± 4.2
Echo Train Length	2	2 - 3	2
Acquired temporal resolution, ms	129.6 ± 7.2	128.7 ± 27.3	112.1 ± 4.3
Time points	12 ± 1	14 ± 1	14
Scan time, seconds	384 ± 72	242 ± 73	246 ± 39
Average heartrate/min	68 ± 16	60 ± 11	66 ± 10

Table 5: Outcome measures and inter-rater results per vendor

	Vendor 1 (n=10)	Vendor 2 (n=20)	Vendor 3 (n=30)	All (n=60)
Outcome measures
N_detected (n)	4 ± 3	18 ± 5	7 ± 5	10 ± 7
v_mean (cm/s)	3.6 ± 1.3	3.8 ± 0.7	5.7 ± 1.4	4.7 ± 1.5
PI	0.87 ± 0.27	0.36 ± 0.11	0.46 ± 0.24	0.49 ± 0.27
SNR_v	5.9 ± 2.5	6.7 ± 2.2	5.3 ± 1.5	5.9 ± 2.0
Inter-rater analysis
DSC	0.88 [0.81 - 0.93]	0.93 [0.90 - 0.95]	0.90 [0.69 - 0.95]	0.91 [0.69 - 0.95]
ICC N_detected	0.75 [0.22 - 0.94]	0.72 [0.43 - 0.88]	0.88 [0.76 - 0.94]	0.92 [0.87 - 0.95]
ICC v_mean	0.80 [0.21 - 0.96]	0.75 [0.47 - 0.89]	0.76 [0.54 - 0.88]	0.84 [0.74 - 0.90]
ICC PI	0.32 [-0.50 - 0.84]	0.79 [0.54 - 0.91]	0.93 [0.86 - 0.97]	0.85 [0.76 - 0.91]

Values are given in mean ± SD for the outcome measures and mean with 95% confidence interval for the inter-rater analysis. N_detected = amount of detected arteries; v_mean = mean blood flow velocity of the perforating arteries given in cm/s; PI = pulsatility index; ICC = intra-class coefficient; DSC = Dice similarity coefficient

Table 6: Outcome measures per site

Measurement	Site 1 (n=10)	Site 2 (n=2)	Site 3 (n=9)	Site 4 (n=9)	Site 5 (n=5)	Site 6 (n=6)	Site 7 (n=10)	Site 8 (n=9)
N_detected (n)	4 ± 3	17 ± 4	16 ± 5	21 ± 3	3 ± 2	4 ± 4	11 ± 3	7 ± 4
v_{mean (}cm/s)	3.6 ± 1.3	3.0 ± 0.3	4.0 ± 0.7	3.8 ± 0.5	6.6 ± 1.6	5.4 ± 1.2	4.9 ± 0.8	6.1 ± 1.6
PI	0.87 ± 0.27	0.43 ± 0.27	0.36 ± 0.11	0.34 ± 0.07	0.59 ± 0.56	0.46 ± 0.05	0.42 ± 0.14	0.45 ± 0.10
SNR_v	5.9 ± 2.5	4.5 ± 2.0	7.6 ± 2.7	6.3 ± 1.2	5.1 ± 1.4	4.1 ± 0.8	6.2 ± 1.7	5.1 ± 1.2

Values are given as mean ± SD. N_detected = amount of detected arteries; v_mean = mean blood flow velocity of the perforating arteries given in cm/s and PI = pulsatility index

Table 7: Coefficient of variation for inter-site and inter-vendor measurements.

Measurement	Inter-site	Inter-vendor
N_detected(n)	62	31
v_mean (cm/s)	21	12
PI	39	12

N_detected = amount of detected arteries; v_mean = mean blood flow velocity of the perforating arteries given in cm/s; PI = pulsatility index; CV = coefficient of variation. Coefficients of variation are given as percentages.

No competing interests reported.

Blood flow velocity analysis in cerebral perforating arteries on 7T 2D phase contrast MRI with an open-source software tool (SELMA)

Status:

Version 1

Abstract

Figures

Introduction

Methods

Algorithm

User interface

Test data

Data analysis

Results

Discussion

Conclusion

Declarations

Author Contribution

Acknowledgement

Data Availability

References

Tables

Additional Declarations

Status:

Version 1