Automated Volume Status Assessment Using Inferior Vena Cava Pulsatility

: Assessment of volume status is important to correctly plan the treatment of patients admitted and managed by cardiology, emergency and internal medicine departments. Non-invasive assessment of volume status by echography of the inferior vena cava (IVC) is a promising possibility, but its clinical use is limited by poor reproducibility of current standard procedures. We have developed new algorithms to extract reliable information from non-invasive IVC monitoring by ultrasound (US) imaging. Both long and short axis US B-mode video-clips were taken from 50 patients, in either hypo-, eu-, or hyper-volemic conditions. The video-clips were processed to extract static and dynamic indexes characterizing the IVC behaviour. Different binary tree models (BTM) were developed to identify patient conditions on the basis of those indexes. The best classiﬁer was a BTM using IVC pulsatility indexes as input features. Its accuracy (78.0% when tested with a leave-one-out approach) is superior to that achieved using indexes measured by the standard clinical method from M-mode US recordings. These results were obtained with patients in conditions of normal respiratory function and cardiac rhythm. Further studies are necessary to extend this approach to patients with more complex cardio-respiratory conditions.


Introduction
The intravascular volume status (i.e., the extent of vascular filling) is a relevant cardiovascular parameter related to the cardiac preload (i.e., the stretch of cardiac tissue in relaxed conditions), which in turn affects cardiac output and arterial blood pressure. Its assessment in critically ill patients is essential to establish and carefully balance the appropriate fluid therapy, whereby fluid supplementation may favor cardiac efficiency, but also increase the rate of complications and mortality [1][2][3]. Various pathological conditions are characterized by alteration of the volume status, e.g., heart failure causes overload while dehydration leads to volume depletion.
The non invasive evaluation of the volume status is very important, as it allows to inspect a patient in emergency conditions or during the follow-up. An approximate assessment of the volemic condition can be obtained from US imaging of the IVC [4,5]. In fact, the pulsatility of the IVC was found to correlate size and pulsatility have been automatically extracted and used to build a classifier able to discriminate patients in either hypo-, eu-, or hyper-volemic conditions.

US Video-Clip Processing
The physical dimension of a pixel was determined for each video-clip as a preliminary step, by scanning automatically a graduated length scale present in the frames. Then, the user selected some parameters (e.g., concerning the portion of frame to be processed, points to be tracked, etc.) needed by two algorithms (implemented in MATLAB R2018a, The Mathworks) that processed the B-mode US video-clips of the IVC in either long or short axis (notice that, as a preliminary interaction with the user is needed to process the US videos, we refer to our processing algorithms as semi-automated). Figure 1 shows an example of IVC border automated delineation obtained by those algorithms, described below. Once obtained the IVC borders in either of the two views, the mean diameter and pulsatility indexes (defined below) were estimated.

Identification of IVC Borders in Long Axis
The algorithm proposed in [20] (and already applied in [13,25,31]) was used. In the first frame of the clip, the user located the vein and the region of interest, indicating two reference points to be tracked to compensate for IVC movements and deformations, the leftmost and rightmost lines to be considered and the location of the borders of the vein in the leftmost line. In optimal conditions, the available tract was between the confluence of the hepatic veins into the IVC and the caudate lobe of the liver. Each frame was first pre-processed with a 2D median filter (neighborhoods of 9×9 pixels). Then, the software uniformly distributed 21 lines in the region of interest and identified the borders of the vein along these lines (as a jump of the US intensity along them). For each subsequent frame, the location and direction of those lines were updated based on the estimated movements of the reference points from the previous one.
Once obtained the superior and inferior borders of the vein, the software computed the IVC midline and distributed uniformly 5 points along it. For our specific application, the extension of the midline was considered only from the 20% to the 80% of its length, thus excluding the edges. Sections orthogonal to the IVC midline passing from each of these 5 points were considered and the IVC pulsatility was estimated for each of them.

Identification of IVC Borders in Short Axis
The algorithm proposed in [28] (and already used in [29]) was employed. The user was asked to indicate the centre of the IVC and to draw a rectangle enclosing it in the first frame of the video-clip. Subsequent frames were cropped in a rectangular region with the same dimension, centered on the IVC estimated on the previous frame. The image was converted in gray-scale, contrast enhanced using histogram equalization and processed with a 2D median filter (neighborhoods of 11×11 pixels). The outline of the vein was then estimated by the algorithm. Twenty rays were defined, originating from the centre of the considered rectangular portion of image and sampling uniformly the directions around it. For each ray, the intensity of the image along it was estimated by cubic interpolation. The border of the vein was identified as an abrupt increase of the intensity (from the lumen to the outside tissues).
Once the 20 border points were found, their coordinates were low pass filtered (Butterworth non-causal, zero-phase IIR filter of order 4 with cut-off at 0.3) to get a smooth boundary of the vein. Furthermore, the maximum variation of the length of a ray was imposed to be 5 pixels; the rays which overcame such a threshold were removed and substituted by a quadratic interpolation of the 4 closest neighboring border points.

IVC Indexes
The mean diameter was estimated averaging both across different sections (i.e., 5 sections in long axis and 10 diameters corresponding to the 20 rays in short axis) and time (i.e., considering the frames of the video-clips).
Pulsatility was measured in terms of the CI where D indicates the dimension over the time variable t of IVC, expressed either as diameter or equivalent diameter (proportional to the square root of the area [28]), in the long and short axis, respectively, and max/min indicate local extrema. Local maxima and minima were computed for each respiratory cycle. A CI accounting for the overall pulsatility was obtained by averaging the estimations across different respiratory cycles and different sections (the latter, only in the case of the long axis approach). Additional indexes were also estimated by decomposing the time series reflecting IVC pulsations into low and high frequency components (below 0.4 Hz and above 0.8 Hz, respectively), assumed to reflect the stimulations induced by either respiration or heartbeats, respectively (both filters were 4 th order Butterworth, used twice, once with time reversed, to remove phase distortion and delay). From these filtered time series, applying again the definition of CI (1) on local maxima and minima, the respiratory caval index (RCI) and the cardiac caval index (CCI) were obtained. Stable estimations of both indexes were computed by averaging across either respiratory cycles or heartbeats (and on the 5 sections, in the case of the long axis). An example of time series extracted from a video-clip is given in Figure 1C.
IVC was also investigated by standard manual measurements, in both long and short axis, in M-mode. Stable estimations of the minimum and maximum IVC diameter were obtained by averaging across more measurements (up to 3). Then, the maximum and minimum diameters were used to compute the CI and the average IVC diameter (defined as the mean of the two diameters).

Experimental Data
Inclusion criteria were the presence of pathological conditions in the Emergency Department and in the Department of Medicine resulting in overload (heart failure) or volume depletion (dehydration or moderate bleeding). As a control group, patients without the previous conditions were selected. Exclusion criteria were chronic obstructive pulmonary disease, pulmonary hypertension, interstitial disease or thromboembolism, tension pneumothorax, cirrhosis and/or ascitic effusion, serum creatinine >3 mg/dl, constrictive pericarditis and cardiac tamponade. Fifty patients were included in the study. They were selected from a database of 69 patients (Table 1). On the basis of clinical considerations (based upon physical examination, laboratory data and imaging), each patient was associated to one of the following classes: 1.
US B-mode video-clips of about 15 s were recorded bedside in spontaneous breathing, with subxifoideal approach, using a MyLab Seven system (Esaote, Genova, Italy; frame rate 30 Hz, 256 gray levels) equipped with a convex 2-5 MHz probe. M-mode scans were also recorded to allow for standard manual measurements.
According to the Declaration of Helsinki, subjects provided written informed consent for the collection of data and subsequent analysis. The study was approved by the local Ethics Committee.
The data included here were only those for which both video-clips, recorded along either long or short axis, could be reliably processed. They were 20 from patients in overload, 19 controls and 11 patients with volume depletion (notice that video-clips of patients with volume depletion are more difficult to be processed, due to the small dimension of the IVC, which could even collapse in some frames, hindering proper processing). Figure 2 shows examples of patients in the 3 classes. Table 1. Number of patients included in different groups (with indication of the entire database and of the patients for which successful processing of both long and short axis ultrasound videos was achieved).

Automated Identification of the Volemic Status
Three different classification approaches were fit to our dataset, as a preliminary step: the error-correcting output codes (ECOC) model, using support vector machines (SVM [32]) for binary one-to-one classifications [33,34]; the Naive Bayes classifier (estimating data distributions using smoothed densities with normal kernel) [35]; the BTM [35]. For each approach, different models were fit to our data, considering all possible combinations of input features (detailed below). The performances of different classifiers were compared in terms of a 10-fold cross-validation test, which allowed to select the best input features and classification approach. Then, the selected classifier was tested by a leave-one-out approach and, finally, trained on the entire dataset, to provide an ultimate prediction model. In the following, we will focus only on the BTM, as best results were obtained using this approach.
Different BTMs were fit to our multi-class classification problem (including 3 classes), selecting the simplest one (i.e., with minimum dimension) with best performances. A BTM iteratively splits the dataset in two groups, after comparing an index with a threshold (Gini's diversity index was used as splitting criterion). Thus, it is built by choosing the optimal number of splittings, the specific index to be considered for each binary separation and selecting the threshold value for each splitting. Different BTMs were developed considering all possible combinations of input indexes (exhaustive search): all possible choices of a single index, all pairs, triplets, ... until using all indexes. Examples of data from patients in either hypo-, eu-or hyper-volemic conditions. The first frames of the long and short axis scans are shown (left and right, respectively), together with the IVC boundaries identified by the algorithm. Time series are also shown for the diameters in 5 sections of the IVC (in gray, with superimposed the mean diameter in black) and for the IVC area estimated from the long and short axis scans, respectively. In the case of long axis scans, pulsatility indexes were computed as averages of estimations from each of the 5 sections; in the case of short axis scans, they were computed from the equivalent diameter, proportional to the square root of the IVC cross-section area. D m : mean diameter; A m : mean area; CI: caval index; RCI: respiratory caval index; CCI: cardiac caval index. Different sets of indexes were used, considering the possibility of either employing the semi-automated processing or not (so that in the latter case only manual measurements were considered).
The set of indexes obtained by semi-automated video processing was the following: For each set of indexes, different BTMs were developed using all possible combinations of features taken from it and the one with highest performance was selected (thus, they were 255 and 15, for the first and second features set, respectively). Specifically, the best categorical predictor split was chosen from all possible combinations of choices. As mentioned above, the models were cross-validated considering 10 folds. The order of the data was random, so that the three categories of patients had a similar representation in each fold (however, they could not be equally represented; this problem is emphasized by the small size of our dataset). The one providing minimum average root mean squared regression error (or loss) on the validation sets was then selected. This specific model was then tested by a leave-one-out approach, to reduce the bias in error estimation (considering our small dataset) [36].

Results
Indexes characterizing the IVC were extracted from long and short axis views by either semi-automated processing or manual estimation (performed in M-mode). Then, they were used to classify patients. As using indexes extracted with the automated processing resulted in better performances, figures and tables shown below refer to those data, indicating in the text some performance indexes of the best BTM developed using the set of indexes obtained by standard manual measurements. Figure 3 shows the BTM selected as the classifier with best performances on our dataset. The shown BTM was trained on the entire dataset, including the best input features selected by the cross-validation test (described in Section 2.3), where minimum loss was obtained (equal to 0.26; the loss of the best classifiers using either ECOC or Naive Bayes models was 0.28).
Two pulsatility indexes are included: CCI in long axis and CI in short axis. The same loss was obtained by other 4 BTMs: the one with minimum number of input features was selected. The CCI in long axis was included in 4 of these BTMs with minimum loss; the CI in short axis was included in 2 of them. Another feature which was often included was the RCI in short axis, which was used in 3 among the 5 BTMs with minimum loss. In the case in which standard manual measurements were employed, the best BTM was unique, it had a loss of 0.28 and included two indexes: IVC diameter estimated in long axis and CI in short axis.
Distributions of the indexes are shown in Figure 4. The mean Fisher ratios (FR, considering all 3 binary comparisons) of the indexes selected by the best BTM between those estimated by semi-automated processing are among the highest. However, they have not the highest FRs: indeed, the best discrimination in terms of average FR is provided by the mean diameter estimated from the long axis view. This indicates that the selected indexes are those that are both informative and not much redundant, allowing a peak in performance of the classifier using them as inputs. Notice also that the FR is an index of linear discrimination, whereas the adopted classifier allows for nonlinear separation.
It is interesting to see that the indexes estimated manually have even higher FRs, indicating a better linear discrimination of the patients. The two indexes with highest FRs are those selected by the best BTM using only indexes measured manually. However, the semi-automated processing allows to extract additional information: specifically, the two pulsatility indexes RCI and CCI reflect the effect of different stimulations (respiration and heartbeat, respectively). This further information (and specifically that coming from the CCI) allows the BTM from automated processing to get better performances than the one developed on the basis of the set of manually estimated indexes.
The confusion matrix of the best BTM shown in Figure 3 is given in Table 2. Notice that all hypo-volemic patients were correctly identified. A few eu-volemic and hyper-volemic subjects were misclassified. No hyper-volemic patient was confused as hypo-volemic or vice-versa. Common performance indexes are the followings: mean sensitivity 90.0% (86.0% for the BTM built using the manually estimated indexes); mean specificity 95.0% (91.9% with manual indexes); positive predictive value 90.0% (86.2% with manual indexes); negative predictive value 94.2% (91.8% with manual indexes); mean accuracy 92.9% (89.8% with manual indexes).  Notice that these performances were obtained using the entire dataset to train our model. As some misclassifications were obtained, we deduce that some information is still missing and/or the features extracted by our processing contain some residual noise. To get a more faithful indication of performances, a leave-one-out test was performed (i.e., the best features selected before were kept, but each sample was excluded in turn from the training set and used for testing). The confusion matrix in Table 3 was obtained. Some degradation of the performance can be observed, especially in the discrimination of the control and hyper-volemic groups. The following performance indexes were achieved: mean sensitivity 70.0% (66.0% for the BTM built using the manually estimated indexes and tested by a leave-one-out approach); mean specificity 83.2% (80.4% with manual indexes); positive predictive value 70.0% (65.1% with manual indexes); negative predictive value 82.1% (80.5% with manual indexes); mean accuracy 78.0% (75.3% with manual indexes). Table 2. Confusion matrix of the best binary tree model classifying the volemic status, shown in Figure 3 (for comparison, the best error-correcting output codes and Naive Bayes classifiers trained on the entire dataset show a predictive value of 78% and 86%, respectively).

Predicted
Target Score Predictive

Discussion
The accurate assessment of the volume status is of relevance for a high percentage of patients who either access the emergency room or enter the medical wards. The development of new standardized clinical procedures and the support of automatic algorithms can help to correctly plan the treatment and monitor the follow-up.
We have developed two algorithms that help standardizing the assessment of IVC pulsatility. They process B-mode video-clips, allowing to compensate for either longitudinal or transverse respirophasic movements and to delineate the vessel edges in an entire region (either a tract of longitudinal section or a cross-section). As detailed in the Methods section, the IVC indexes extracted are not related to a single diameter along the IVC, but reflect the average size and pulsatility, calculated over the whole portion of considered IVC length (in long axis) or over the entire IVC cross-section (in short axis). Thus, the considered IVC indexes reflect the overall behaviour of the investigated regions (both static and dynamic behaviour, reflected by the size and pulsatility of the IVC, respectively). We have already documented for the long axis scans that this approach is more reliable and repeatable than standard clinical assessment [13,19,25]. Moreover, we have shown that the IVC in a short axis view can pulsate differently along different directions [28], so that an average indication of cross-sectional pulsatility is preferable and less subjective than referring to an arbitrarily chosen diameter.
Here, we have built BTMs using those indexes estimated by our algorithms to assess automatically the volemic conditions of patients. The indexes used by the best BTM reflect IVC pulsatility. Referring to Figure 3, the joint integration of information from CCI from long axis US scans and CI in short axis allows to identify the different conditions, with an accuracy of 78% in a leave-one-out test (larger than what could be achieved with the best classifier using only manual indexes). Notice that CCI is an index that was introduced recently [13,25,37,38] and whose estimation is expected to be stable, as the heartbeats are much less variable than respiratory cycles, mainly affecting the measurements of CI and RCI. Other 4 BTMs achieved the same loss in cross-validation as the best one (which was chosen because it had the smallest dimension). This can be interpreted as a consequence of the redundancy included in the pulsatility indexes, whereby one index may be obtained from a combination of the others. This result may also descend from the small sample size, which does not allow to appreciate fine differences in performance among models with high classification rates. Hence, this should be considered as a pilot study. Augmenting the numerosity of the sample would be important to get a more stable estimation of the classification model.
We have compared the classification performances of the above mentioned fully automated method (based on indexes extracted by processing US B-mode video-clips), with a BTM using indexes measured manually, with M-mode scans along an US ray selected either from a longitudinal or a transverse view of the IVC. It is interesting to notice that the indexes measured manually allowed in general to get a better linear discrimination of the volemic conditions (measured in terms of the average Fisher ratios comparing all pairs of groups). However, the automated processing allowed to extract more indexes describing IVC and the final best classifier showed better performances than that obtained using manual measurements. In particular, additional information on IVC pulsatility induced by either respiratory cycles or heartbeats was available and CCI (from long axis view) was selected by the best BTM. We deduce that this index includes additional/not redundant information that, together with other characterizations of IVC pulsatility (provided by the CI in short axis, in the best BTM), can be useful to disentangle the complex/nonlinear relation between IVC dynamics and volume status of the patient.
It must be underlined that the present results were obtained from a selected group of patients in which pathologies specifically affecting the respiratory system were excluded. Moreover, we expect that the selection of CCI as optimal feature depends also on the regular hearth rhythm shown by the patients included in our dataset; in the case of arrhythmia (typically due to atrial fibrillation, not shown by our data sample), a reliable estimation of this important parameter would be hindered. Thus, the application of our classification approach to different patient populations could result in different selections of parameters and thresholds. Nevertheless, as an effort to overcome the subjectivity of the measurement, our approach (but, probably, not the classification model) remains valid and worth to be investigated and extended to other patients groups.
Another limitation of the method is the need to rely on good quality imaging. Indeed, only 87% and 77% of long and short axis video-clips were properly processed, respectively, so that only 72% of our patients could be included in this study (as the processing of both recordings was required). Improvements could be obtained by adopting higher level US machines or by more effective image processing. We are currently trying to optimize our algorithms in order to process US recordings in real time, providing a feedback to the operator. We expect that this could help in getting successful processing in more US video-clips. Indeed, our present offline approach requires that the operator acquires data blindly, i.e., without knowing if the recorded video-clip will be adequate for processing. Instead, a real time software could guide the acquisition and indicate to the operator if there are problems in processing the data, in which case the operator could work at improving the quality of the imaging. This is exactly what happens in manual measurements: the operator may try different approaches and strategies to improve image quality until he is satisfied with the result and IVC measurement is made possible.
In summary, we have shown the joint application of long and short axis US views of the IVC, to assess the volume status of patients. The US videos have been automatically processed by multi-section and multi-directional algorithms, which track IVC movements and compute its size and pulsatility either over a longitudinal portion of the vessel or a cross-section, respectively. The IVC pulsations have been also split into two contributions, reflecting either the respiratory cycles or the heartbeats. The algorithms have been widely tested on healthy subjects in laboratory conditions in the past [13,19,20,28,29] and in a single clinical study, aimed at estimating right atrial pressure based on the analysis of IVC pulsatility [25,31]. Here, the different indexes were jointly applied in a clinical setting and used to solve the multiclass problem of discriminating patients with different volume status, showing better performance than when using manually measured indexes. Pulsatility indexes estimated from both long and short axis have been included in the best classification model, which supports the concept that they convey complementary information. Even considering the preliminary nature of these results (given the small sample size), the approach appears to be very promising. Extending the dataset and improving the processing algorithms (e.g., allowing real time interaction with the operator) may prospectively lead to obtain efficient systems for diagnostic support and follow-up.

Conclusions
The identification of the volume condition is important for the clinical management of many patients in the emergency room or in wards of general medicine and cardiology. We propose an automated approach for the classification of the volemic status, based on the processing of B-mode US video-clips of the IVC and on the extraction of pulsatility features. The presented results suggest that this approach may be useful to get more reliable clinical indication from the US monitoring of IVC. Investigation over a larger dataset will however be necessary to test the actual effectiveness of the proposed method. Moreover, our results hold true in conditions of normal respiratory function and cardiac rhythm. It is reasonable that our classifier will not apply to patients with more complex cardio-respiratory conditions; however, the same approach could be applied to develop models fitting their conditions.

Patents
An instrument implementing the algorithms for IVC delineation used in this paper was patented by Politecnico di Torino and Universitá di Torino (WO 2018/134726). Funding: This research was carried out as part of the project "Method and apparatus to characterise non-invasively images containing venous blood vessels", funded through the PoC Instrument initiative, implemented by LINKS, with the support of LIFFT, with funds from Campagnia di San Paolo.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: