Kinematic analysis of diastolic function using the freely available software Echo E-waves – feasibility and reproducibility

Background

Early diastolic left ventricular (LV) filling can be accurately described using the same methods used in classical mechanics to describe the motion of a loaded spring as it recoils, a validated method also referred to as the Parameterized Diastolic Filling (PDF) formalism. With this method, each E-wave recorded by pulsed wave (PW) Doppler can be mathematically described in terms of three constants: LV stiffness (k), viscoelasticity (c), and load (x₀). Also, additional parameters of physiological and diagnostic interest can be derived. An efficient software application for PDF analysis has not been available. We aim to describe the structure, feasibility, time efficiency and intra-and interobserver variability for use of such a solution, implemented in Echo E-waves, a freely available software application (www.echoewaves.org).

Results

An application was developed, with the ability to open DICOM files from different vendors, as well as rapid semi-automatic analysis and export of results. E-waves from 20 patients were analyzed by two investigators. Analysis time for a median of 34 (interquartile range (IQR) 29–42) E-waves per patient (representing 63 %, IQR 56–79 % of the recorded E-waves per patient) was 4.3 min (IQR 4.0–4.6 min). Intra-and intraobserver variability was good or excellent for 12 out of 14 parameters (coefficient of variation 2.5–18.7 %, intraclass correlation coefficient 0.80–0.99).

Conclusion

Kinematic analysis of diastolic function using the PDF method for Doppler echocardiography implemented in freely available semiautomatic software is highly feasible, time efficient, and has good to excellent intra-and interobserver variability.

Diastole, the time during which the ventricles of the heart fill with blood, is conventionally assessed by measuring various blood flow velocities, time intervals, ratios of velocities, and dimensions of chambers and walls of the heart [1]. The investigator then tries to fit these measurements with a theoretical and empirical body of knowledge constituting part of our understanding of left ventricular (LV) mechanics and hemodynamics. As such, the conventional method of assessing diastolic function lacks a unifying theory. Improving our ability to investigate diastolic function is of importance not only as it is vital to our understanding of cardiac function in health and disease, but also in order to better approach diagnostics and treatment options in the increasingly recognized condition of heart failure with preserved ejection fraction. We aim to describe an alternative strategy for assessing diastole, and a software implementation for this method that is ready for use in clinical research.

During the cardiac cycle, the LV contracts during systole, whereby the volume of blood in the LV decreases. During early diastole, the LV recoils, leading to an increase in the volume of blood in the LV, with blood entering the LV through the mitral valve. In the absence of aortic regurgitation or shunts, the velocity of blood flow through the mitral valve during early diastole will be the same as the velocity of volume expansion of the LV, and hence the velocity of LV recoil. The velocity of transmitral blood flow can readily be measured by pulsed wave (PW) Doppler echocardiography. It has been shown that the velocities of early diastolic transmitral blood flow, and thus LV recoil, can be accurately described as a case of damped harmonic motion [2]. Damped harmonic motion is a part of the analytical framework describing motion in classical mechanics, and can been used, among other things, to describe the recoil of a spring. From this framework one can derive a mathematical expression that describes velocity as a function of time in terms of three constants, namely stiffness (k), viscoelasticity, or energy loss, (c), and load (x₀). By curve fitting the velocity profile of the PW Doppler signal of the E-wave to the mathematical expression, the constants can be obtained. The theoretical framework describing LV filling in this manner has been termed the parameterized diastolic filling (PDF) formalism [2]. Using gold standard high fidelity invasive hemodynamic measurements, the stiffness constant k has been shown to have a close linear relationship with LV diastolic stiffness [3]. In the same manner, it has been shown that the influence of the damping constant c can be used to formulate an accurate estimate of the invasively determined time constant of isovolumic pressure decay, tau [4]. Furthermore, a load independent index of diastolic filling, called M, can be calculated by examining the relationship between the peak driving force and the peak resistive force of filling at different levels of load [5]. Furthermore, clinical studies have illustrated the utility of the PDF method in characterizing clinical disease including diabetes [6, 7], hypertension [8], as well as prognosis in congestive heart failure [9]. A brief overview of the PDF parameters is given in Table 1, and further details are given in the Appendix. In summary, this unified approach to studying early diastolic filling narrows the gap between clinical echocardiography and classical mechanics, and also promises new insights into different pathological states as well as basic cardiac physiology. However, a user-friendly solution enabling PDF analysis of echocardiographic DICOM images and output of results in a fashion optimized for clinical research use has not been available [10]. Therefore, we aim to describe the structure, feasibility, time efficiency and intra-and interobserver variability for use of such a solution, implemented in Echo E-waves, a freely available software application.

The program, Echo E-waves, is freely available for download at www.echoewaves.org [11]. It is implemented in the MATLAB® release 2015a (Mathworks, Natick, Massachusetts, USA) environment for Windows (64-or 32-bit architecture). In the following section we describe the structure and functionality of the program. For analysis of inter-and intraobserver variability, sets of E-waves with visually moderate to good quality were selected from 20 patients undergoing echocardiography as part of ongoing clinical studies on adults with acute myocardial infarction with or without obstructive coronary artery disease, and excluding patients with significant valvular disease. The studies were approved by the Regional Ethical Review Board in Stockholm, Sweden, and all subjects provided written informed consent. All echocardiographic acquisitions were performed on Vivid E9 scanners (GE, Horten, Norway). The Echo E-waves program was run on an Intel® Core i5™ 2500 K processor with a NVIDIA® GeForce® GTX 670 graphics card and the 64-bit version of Windows 7.

Echocardiographic image acquisition

The program can analyze PW Doppler echocardiographic images of transmitral flow acquired in the same manner as is recommended by clinical guidelines [12]. However, some adjustments are suggested to ensure improved accuracy in analysis. A horizontal sweep speed of 100 mm/s is recommended, as lower sweep speeds lead to a substantial loss of temporal resolution. A sweep speed of 150 mm/s or 200 mm/s is also acceptable. However, in our experience these sweep speeds make the visual recognition of the shape of the E-wave somewhat more difficult, since they are less commonly used clinically. The program currently reads DICOM images from GE, Philips, Siemens and Acuson systems. When using a GE or Philips system, the user can export a cine recording of Doppler data stretching over several screen widths as a single DICOM multi frame file. A method for this is described in on the program website (www.echoewaves.org). From all the systems above, single frame DICOM image files can also be loaded. For the purposes of testing the intra-and interobserver variability, 20 echocardiographic exams collected during ongoing clinical trials were selected. The exams were performed by four dedicated sonographers, all with at least 8 years of experience with clinical echocardiography. A 3.5 mm sample volume was placed at the tips of the mitral leaflets, and PW Doppler recordings were performed both during free breathing and after cessation of the Valsalva maneuver.

Loading and visualization of Doppler data from DICOM images

Using DICOM metadata, the Doppler region of the image is automatically identified and displayed with correct time and velocity scaling. For ease of overview, the contents of the folder are displayed as thumbnails in the Heartbeat Browser, also enabling navigation of the data set using mouse or keyboard input. Extended (cine) recordings of Doppler registrations from GE and Philips scanners are automatically cropped to one image per cardiac cycle. In order to facilitate the detection of the part of the Doppler signal constituting the E-wave, the user can change color maps, gain settings and apply contrast stretching. Currently, the program does not load the contents of folders with a mix of single frame and cine DICOM files. The graphical user interface is described in Fig. 1.

Overview of the graphical user interface. 1: Display of Doppler data, velocity profile detection and fitted curve for a given heart beat. 2: Numerical display of results, per E-wave (list on left) and summarized (table on bottom right). 3: Overview of each E-waves in the current acquisition comprised of multiple E-waves. 4: Graphs of (top) the basis for calculation of the load independent index of filling, M, and (bottom) all accumulated PDF curves for each E-wave

Full size image

Detection of the velocity profile

In the data structure of a DICOM image, Doppler data is represented as a two dimensional matrix of numbers, with each matrix index corresponding to a velocity at a given point in time. In these images, the horizontal axis of the image is time in seconds, and the vertical axis is velocity in meters per second. The value at a given matrix index indicates the brightness of the pixel at the corresponding position in the image, and the brightness corresponds to the amount of blood travelling at that given velocity at that point in time. For the purpose of detecting the velocity profile, the program searches each column (i.e. each time point) of the image from the top down for the first pixel having a brightness higher than or equal to a user defined threshold level, matching each time point with a velocity. However, as the Doppler recording is susceptible to noise, and due to the presence of non E-wave signal, the algorithm has been optimized by taking into account pre-existing knowledge of how the velocity profile of the E-wave can be differentiated from Doppler signal not originating from the E-wave. For instance, velocity at the moment when flow is about to begin must be zero. Furthermore, we can assume that there must be an upper limit to the velocity difference between two time points a few milliseconds apart. The algorithm accordingly discards a detected velocity that differs too much from the previous velocity. Discarded velocities are displayed for reference and transparency as to the behavior of the algorithm.

The user can adjust the input to the curve fitting algorithm with regards to the time points along the horizontal axis of the Doppler image at which the detection algorithm starts and ends, as well as the threshold of Doppler signal intensity at which detection occurs. Notably, it is of great importance that the starting point is set accurately at the onset of flow, as a small deviation can lead to substantial differences in the constants resulting from curve fitting. In many cases, the exact position of the onset of flow can be difficult to ascertain visually at first glance. However, the velocity profile during the acceleration time is often quite distinct, and if close attention is paid to how different starting points affect the visual correspondence of the PDF curve to the underlying Doppler data of the acceleration time and the upper, often well defined part of the E-wave, a convincing position for the starting point can be identified in effectively all cases. Positioning of the end point must take into account the various non-E-wave Doppler signals that frequently occur during the deceleration phase. A good curve fit will often be achieved if the end point is set at approximately 75 % of the total E-wave duration and typically at less than 50 % of the maximum E-wave velocity. In the case of partial overlap of the E-and the A-wave, the user must examine the E-wave to see if it is of sufficient duration to determine where it would have ended had there been no overlap between the E-and the A-wave. Typically, this is not possible when the crossover point between the E-and the A-wave occurs at a velocity is above 50 % of the maximum velocity of the E-wave. The program contains no particular feature to guide the user in this process, but using the alternative method to curve fitting described below can often make the task easier. We recommend excluding E-waves from analysis when there is an overlap with the A-wave at greater than 50 % of the maximum velocity of the E-wave. To further increase the agreement between the final fitted curve and the E-wave portion of the Doppler signal, the program also has two options that will make the deceleration time of the fitted curve shorter. The rationale for this is that the various artifacts and non-E-wave blood flow signals during the deceleration phase usually have higher velocities than the true E-wave, and thus will prolong the deceleration time of the fitted curve if they are included. The first approach to alleviate this is to simply replace the final two detected velocities with the value zero. The second approach is to first fit a curve to the detected velocities, as described below, and then calculate the tangent at the time point where the velocity has decreased to 70 % of the maximum velocity. Extrapolating this tangent to the velocity baseline will form a close approximation of the line drawn in clinical echocardiography to mark the deceleration time. This straight line of velocities is then substituted for the velocities of the original fitted curve during the corresponding time frame, and this new set of velocities is then fitted again, producing a final curve with a shorter deceleration time. On a modern computer, these calculations are performed virtually instantaneously. By default, both approaches for making the deceleration time agree more with clinical practice are enabled. It should be noted that these two alterations of the curve fitting process, although having justifications in physiology (if the E-wave ends, its velocity must be zero at the end) and the nature of the Doppler acquisition (noise and non-E-wave signal during the deceleration time and diastasis), they are based on empirically determined parameters. In our experience, using them will lead to a faster and more intuitive curve fitting workflow, in that the sensitivity to the problems mentioned above will be reduced. These fitting settings should be viewed as improvements to speed of effectively obtaining a visually acceptable fit, and not as changes to the PDF model per se. If the user so chooses, they can be disabled individually.

Curve fitting and adjustment of the fitting algorithm input

The paired time-velocity data obtained during detection is used as input to a native curve fitting function of MATLAB, which employs the Levenberg-Marquardt algorithm, a standard algorithm for solving non-linear least square problems. After curve fitting, the constants c, k, and x₀ are obtained with values that, when used as input to an expression of velocity as a function of time, would yield the closest match to the originally detected velocity profile. This calculated curve is displayed as an overlay to the original Doppler image, along with the result of the velocity profile detection algorithm for visual confirmation. After display of the fitted curve, the user can adjust the start and end points of the velocity profile detection time interval, as well as the threshold for detection, resulting in an automatic update of the fitted curve. These adjustments can be made by clicking and dragging with the mouse, or by using keyboard shortcuts. No measures of goodness of fit for the E-wave are displayed. Mathematically, the curve fitting process will produce an excellent fit in nearly all cases, and differences in e.g. sum of squared errors cannot be used to decide if an adjustment to a fit will result in a physiologically more accurate fit. Rather, we have found that visual assessment is the most reliable and useful method for assessing the appropriateness of a given fit.

Automatic detection and curve fitting

In order to facilitate curve fitting of data sets comprised of large numbers of E-waves, a semi-automatic algorithm for curve fitting several E-waves has been implemented. The user must first set the starting point and adjust the settings so that the first E-wave is fit in an acceptable manner. The program stores the part of the image corresponding to the start of this user-defined E-wave as a template. The template is then used to find the starting point of the E-wave on the subsequent images using normalized cross correlation. The velocity profile detection algorithm described above is then run from this starting time point to either 1/3 of a cardiac cycle, or to the time point where velocity has fallen below 35 % of the peak velocity of the E-wave. The detected velocities are used for curve fitting as described above. The resulting parameters are compared to the mean of the previously obtained parameters, and if any parameter of the current E-wave deviates more than 40 % from the mean of the previous, the program discards the fit and moves on to the next E-wave. Usage of this automated method can speed up analysis, however, it requires close scrutiny of the obtained results.

Quality control

In order to review produced results, and also facilitate the identification of anomalies, the results are displayed simultaneously in multiple fashions. The numerical values of the constants c, k, and x₀ are displayed in a list, and the mean of the constants and the derived parameters are also displayed. Furthermore, the calculated curve for each E-wave is displayed in a graph, where each curve is selectable. Each E-wave is also displayed as a selectable marker in a plot of peak driving force vs. peak resistive force. Selecting an E-wave in any of these displays retrieves and displays the corresponding E-wave with its fitted curve, which can then be adjusted or deleted if the user so chooses.

Alternative method to curve fitting

It has been shown that the curves yielded by applying the mathematics of damped harmonic motion to Doppler E-wave velocity profiles can be fully described in terms of unique combinations of their acceleration and deceleration times and peak velocities [13]. In other words, each unique combination of the constants c, k, and x₀ also has a unique combination of acceleration and deceleration time and peak velocity. As a consequence of this, the constants for an E-wave can be calculated if these two time intervals and the velocity are known. As an alternative to fitting the full data set of velocities, it is thus possible to use just these three measurements. In the program, this is implemented both as a primary way of producing a curve, and as a method for adjusting a previously obtained curve. The user can produce a curve by marking the start, peak velocity and end points. For any curve, regardless of if it was produced by curve fitting or this alternative method, the start, peak and end points are manually adjustable by dragging them with the mouse, enabling an easy user interface for adjusting the shape of the curve, which is updated in real time.

Export of data

All constants and derived parameters can be saved and exported to a spreadsheet. The Appendix details the mathematical calculations employed in the software for subsequent calculation of derived parameters. If Excel® is installed on the computer on which Echo E-waves is running, the results can be saved in.xlsx format, and otherwise in a comma separated values.txt format. Saving of images and graphs is also possible. It is also possible to copy the results to the computer’s memory clipboard and paste them into a program of choice. The program also stores the part of the DICOM images containing the Doppler data as well as any results of curve fitting and user interaction, so that images and/or curve fits for a given set of DICOM files can be retrieved, reviewed and/or edited at a later time. Saving in this manner also enables full optional anonymization, leaving the saved file without any trace of patient or source file identification. These data files are created in the native .mat format of MATLAB. When saving the results of an analyzed case, it also possible to designate a results folder, which will save all results to a common space, and also append each case to a file (.xlsx or.txt) with the means and standard deviations of all constants for each case, effectively creating a study database.

Statistical analysis

All statistical analyses were performed in MATLAB release 2015a (Mathworks, Natick, Massachusetts, USA). The basic PDF constants and the derived parameters were analyzed as mean values for each patient. The coefficient of variation (CV) for the difference between measurements, intraclass correlation coefficient, and percentage difference with standard deviation were calculated and presented as needed. The CV was calculated as the standard deviation of the mean differences between analyses divided by the mean and expressed as a percent. Durations of analysis and numbers of analyzed E-waves are reported as median and interquartile range (IQR).

Representative results of the performance of the velocity profile detection algorithm and the corresponding fitted curves are shown in Fig. 2.

Representative examples of Doppler velocity profile edge detection and the resulting fitted curves for two different patients (top and bottom, respectively). Left panels: original Doppler images. Mid panels: detected edges of the velocity profiles in white. Right panels: fitted curves (white). PDF parameters for the top row: c 17.3 g/s, k 135 g/s² x₀ 11.4 cm; bottom row: c 17.7 g/s, k 226 g/s², x₀ 16.0 cm

Full size image

For each patient, 34 (29–42) E-waves were analyzed. Thus, 63 (56–79) % of the recorded E-waves per patient were deemed to be of sufficient quality to perform curve fitting with visually acceptable results. The time for analyzing one patient case was 4.3 (4.0–4.6) minutes. These timings were calculated by letting the automatic fitting algorithm analyze all E-waves, after which the investigator reviewed the whole set, making adjustments as necessary. The median, range and interquartile range of the analyzed parameters are shown in Table 2. For the peak driving force, the per patient mean range was 8 mN.

Inter-and intraobserver repeatability

The results for observer variability are presented in Table 3. Overall, inter-and intraobserver reliability was good to excellent. The investigators were free to discard E-waves for which a visually convincing fit was not possible to obtain. Limiting analysis to only those E-waves for which both investigators had accepted fits yielded similar results. The notable exception to this was the parameters M and B. For M, interobserver variability was high when all fitted E-waves were considered. Discarding those E-waves for which one of the investigators had not provided a result yielded a substantial improvement. By comparison, the interobserver variability for B was substantial using both approaches.

We have described how the freely available software application Echo E-waves can be used to analyze transmitral Doppler recordings in the PDF framework. The application can be used to obtain results in a timely fashion, and makes it feasible to investigate further usage of this approach to measuring and assessing diastolic function. Furthermore, intra -and interobserver reliability was good, with the exception of the load independent index of filling, M, and the intercept B. There may be several reasons for these findings. As it is defined, M is the slope of a regression line describing the relationship between driving and resistive forces under different loading conditions. As for all linear regression slope values, the accuracy of M can be disproportionately influenced by outliers, and the decision to include a particular E-wave can have a sizable impact on M. In order to alleviate the disproportionate impact of potential outliers, one would prefer a data set with a wide variation in load, with a roughly equal spread of E-waves with regards to differences in loading conditions. In our population, using the load variation after the Valsalva maneuver, the per patient mean range of peak driving force was only 8 mN. By comparison, in the study by Shmuylovich et al, the range was 30–40 mN, with load variation induced by varying the body position of the subjects using a tilt board and body positions ranging from head up 90° to head down 90° [5]. It seems likely that this larger load variation could lead to a more accurate assessment of M, however, the present study is the first to analyze the inter-and intraobserver reliability regardless of loading conditions. It should also be noted that the optimal cut off values and ranges for using M for clinical or scientific work have yet to be determined. Other methods for inducing variations in load are worthy of exploration. Furthermore, the optimal number of E-waves necessary for reliable analysis is also worthy of further investigation. It should also be noted that we did not examine the test-retest variability. Future studies incorporating multiple examinations per subject are motivated in order to address the effects of user dependency in echocardiographic image acquisition on quantification of PDF measures.

To our knowledge, there is only one other publicly available software application for PDF analysis [10]. The advantages of Echo E-waves are the ability to analyze DICOM files directly without pre-processing or file conversion, general ease of use, ability to review, compare, quality control, and revise E-wave delineations in the program environment, semi-automatic fitting and efficient export of results in various formats. Furthermore, Echo E-waves does not require installation of other software requiring licensing costs.

The software program Echo E-waves is freely available and can be used to quantify and investigate diastolic function and dysfunction using insights from classical mechanics in the PDF framework, producing results rapidly and reliably enough for application of this method for use in the cardiac imaging community.

DICOM:

Digital imaging and communications in medicine

KFEI:

Kinematic filling efficiency index

LV:

Left ventricle

PDF:

Parameterized diastolic filling

PW:

Pulsed Wave (Doppler)

The authors thank and acknowledge the following individuals for discussions, bug reports and/or other input into the development of the software: Sándor J Kovács PhD MD, Peter A Cain MBBS PhD, Ythan H Goldberg MD, Brett D Atwater MD, Emelie Mörtsell MD, and medical students Fahim Khan, David Megyessi, Mischa Flam, and Madeleine Turesson.

Funding

The studies from which the patients were selected were supported by grants from the Swedish Research Council, the Swedish Heart and Lung Foundation, Stockholm County Council ALF grant and Karolinska Institutet/Stockholm County Council Strategic Cardiovascular Programme. The trials were performed without any relationships to industry.

Availability of data and materials

Project name: Echo E-waves

Project home page: www.echoewaves.org

Operating system: Windows 32-/64-bit

Programming language: MATLAB, Java

Other requirements: MATLAB Runtime version 8.5

Restrictions: All publications containing data produced using Echo E-waves must clearly state this in the description of the study methods, and also reference the present software article.

Due to patient privacy protection, the clinical dataset of the study is only available upon individual request directed to the corresponding author.

Authors’ contributions

MS participated in study design, designed and programmed the software, performed PDF analysis, statistical analyses and drafted the manuscript. KS performed PDF analysis for interobserver variability testing and participated in testing of the software. PT was responsible for the clinical studies from which patients were included into the study, and participated in study design and drafting of the manuscript. MU designed the study, participated in design and testing of the software, and drafted the manuscript. All authors have read and approved of the final manuscript.

Competing interests

The authors declare that they have no competing interests.

Consent for publication

Not applicable.

Ethics approval and consent to participate

The studies were approved by the Regional Ethical Review Board in Stockholm, Sweden, and all subjects provided written informed consent.

Affiliations

1. Department of Clinical Science and Education, Södersjukhuset, Karolinska Institutet, and Cardiology Clinic, Södersjukhuset, SE-188 83, Stockholm, Sweden
   • Martin G. Sundqvist
   •  & Per Tornvall
2. Department of Clinical Physiology, Karolinska Institutet, and Karolinska University Hospital, SE-171 76, Stockholm, Sweden
   • Katrin Salman
   •  & Martin Ugander

Corresponding author

Correspondence to Martin Ugander.

Reprints and Permissions