Mitral Valve Regurgitation Murmurs—Insights from Hemoacoustic Computational Modeling

Abstract

Mitral regurgitation (MR) is the leakage of blood from the left ventricle into the left atrium during systole through a mitral valve that does not close fully. A systolic murmur is produced by MR and can be used to diagnose this disease. In the current study, we use hemoacoustic simulations to characterize the features of murmurs for a range of severities relevant to chronic MR. The incompressible Navier–Stokes equations are solved using an immersed boundary method to simulate the blood flow. The resultant pressure fluctuations on the lumen wall serve as the source for the murmur, and the murmur propagation through the thorax is modeled as a 3D elastic wave in a linear viscoelastic material. The resulting acceleration on the surface of the thorax is used as a surrogate for the measurement from a stethoscope, and these characteristics of the acceleration signal are examined in detail. We found that the intensity of the MR murmur is lower at the mitral point on the precordium, as compared with the aortic and pulmonic areas. This is somewhat counterintuitive but is supported by other studies in the past. We also found that the intensity of the murmur, as well as the break frequency, are well correlated with the severity of MR, and this information can be useful for automated auscultation and phonocardiographic applications.

1. Introduction

A healthy mitral valve facilitates unidirectional blood flow from the left atrium to the left ventricle. During ventricular systole, the left ventricle begins to contract to pump blood into the aorta, and the mitral valve closes to prevent backflow to the left atrium. Mitral regurgitation occurs when the mitral valve is unable to close properly, allowing flow back into the atrium from the ventricle. Mitral regurgitation (MR) is one of the most common valvular diseases in the world, with a significant annual incidence [1,2]. In the US alone, over 4 million people are affected by MR and the condition is also progressive [1,2]. Mitral regurgitation, if left untreated, can set off a chain of events that eventually lead to heart failure. However, almost half of the symptomatic patients with severe mitral regurgitation are unable to undergo surgery due to a high risk–benefit ratio [3,4,5]. Therefore, an early diagnosis of chronic mitral regurgitation is very important so that early intervention and treatments can be deployed.

However, early-stage MR might be asymptomatic [6]. This coupled with the relatively high diagnostic error limit early diagnosis of this condition, and MR is often discovered only incidentally during a routine examination [6]. Non-invasive diagnosis of mitral regurgitation relies mostly on cardiac auscultation and Doppler ultrasound. However, with regard to auscultation, studies show that with traditional stethoscopes, more than half of the patients with mild mitral regurgitation and one-third of patients with moderate or severe mitral regurgitation have no detectable murmur [7]. The recommended echocardiographic approaches to evaluate MR include the vena contracta (VC) method and the proximal isovelocity surface area (PISA) method. The vena contracta method determines the degree of valvular regurgitation by estimating the minimum diameter of the regurgitant blood jet. The measurement error caused by variation in orifice shape and other complex conditions can, however, lead to misclassification of the mitral regurgitation severity [8]. The proximal isovelocity surface area (PISA) method estimates the average regurgitant flow rate, effective regurgitant orifice area, and the regurgitant volume by using the color Doppler image of the converging flow [8]. The effective regurgitant orifice area is approximated by a hemisphere with the radius of the converging flow. The flow rate and volume are calculated by assuming a uniform velocity profile in space. The areal and volumetric estimations are based on a single instantaneous image frame, which inevitably induces error.

Our long-term goal is to diagnose various heart conditions such as mitral regurgitation using automated auscultation devices, e.g., the “StethoVest” [9]. This particular auscultation device has multiple electronic stethoscope sensors mounted on it to simultaneously measure the acceleration on the chest wall at multiple locations [9]. The potential advantage of such a device over conventional auscultation is that it can be used for longitudinal monitoring of heart conditions at home, without the need for trained medical personnel. Furthermore, such devices can also be integrated into a telemedicine framework, thereby improving the management of patients with chronic heart conditions. We are, therefore, interested in finding a correlation between the wall-normal fluctuation and the severity of MR that could be used in automated auscultation, and we achieve this through the use of computational modeling. We note that several previous studies have used computational modeling to study atrial and mitral valve hemodynamics. For instance, Vedula et al. (2015) examined the effect of atrial hemodynamics on ventricular filling patterns using simulations [10]. Collia et al. (2019) simulated the blood flow field in the left ventricle model including flow-driven dynamics of the mitral valve leaflets [11] and focused on the correlation between the regurgitant volume, the effective regurgitant orifice area, and the wash-out of the left ventricle. However, no previous study has explored murmurs originating from the atrium using computational modeling, and that is the focus of the current study.

Desjardins et al. (1996) explored the correlation between murmur intensity and mitral regurgitation severity by using 2D Doppler imaging, confirming the hypothesis that the murmur intensity correlates well with the degree of mitral regurgitation [12]. However, the statistical relationship between them is unknown, due to substantial measurement error [12]. Lees et al. (1970) attributed the murmur generation to the turbulent pressure fluctuation on the lumen wall [13]. Consequently, it is believed that pressure fluctuations generated by the unsteady flow and leaflets induced by the regurgitant jet excite the heart chamber wall, resulting in the heart murmur. Thus, the current study takes the pressure on the heart wall simulated by computational fluid dynamics (CFD) as the source for the murmur.

2. Methods

The two-step, one-way coupled computational hemoacoustic (CHA) method consists of a hemodynamic simulation of the unsteady flow in the blood vessel in the first step, and the simulation of the propagation of the murmur propagation through the thorax in the subsequent step [14]. The blood flow field was computed by solving the incompressible Navier–Stokes equations, and the linear elastic wave equations were used to model the wavefield associated with the dispersion of the murmur through the thorax. The two steps were coupled by feeding the pressure fluctuation induced by blood flow on the vessel wall into the thorax model as the input source at the corresponding location.

2.1. Flow Simulation

The fluid was assumed to be incompressible and Newtonian. Thus, the regurgitant flow field was generated by incompressible Navier–Stokes equations as follows:

$$\nabla \cdot \overrightarrow{U}=0, \frac{\partial \overrightarrow{U}}{\partial t}+\left(\overrightarrow{U}\cdot \mathit{\nabla }\right)\overrightarrow{U}+\frac{1}{{\rho }_0}\mathit{\nabla }P={\nu }_0{\nabla }^2\overrightarrow{U},$$

(1)

where $\overrightarrow{U}=\left(u,v,w\right)$ is flow velocity; $P$ is pressure; ${\rho }_0$ and ${\nu }_0$ are the density and kinematic viscosity of blood.

Figure 1 displays the left atrium model from the CT scan that was used for regurgitant flow simulation. The left atrial volume of this patient-specific model was 94.6 mL, including the left atrial appendage volume of 4.3 mL. The lesions that induce mitral regurgitation, such as the mitral valve prolapse, annular dilatation, damaged leaflets, tissue chords, and muscles, all diminish the normal systolic coaptation between the anterior and posterior mitral leaflets [15]. Therefore, the mitral valve morphology was modeled as an ellipse, with an aspect ratio of three, positioned at the center along the black dashed line; a time-varying flow, as shown in Figure 2, was prescribed through this elliptic slit to model the regurgitant flow. The time interval between the systolic regurgitant flow and the diastole flow into the left atrium was so large that the flow in the left atrium was almost stagnant when the regurgitant jet entered. Therefore, only one systolic period was computed in this study. Neumann boundary conditions were applied to the outflow boundaries on the pulmonary veins, and Neumann pressure boundary conditions were used at all boundaries of the fluid domain.

Regurgitant fraction (RF) and effective regurgitant orifice area (EROA) are the indices often used for grading the regurgitation severity. Regurgitant fraction is defined as the ratio of regurgitant blood volume per beat to the stroke volume. The effective regurgitant orifice area is estimated as the ratio of the regurgitant volume to the regurgitant jet time-velocity integral determined by color Doppler. In the present case, the effective regurgitant orifice area was simply the area of the ellipse prescribed with a uniform velocity profile on the mitral annulus plane. Some previous studies have shown that there is a logarithmic correlation between the regurgitant fraction and EROA [16]. The velocity profile of the regurgitant jet was obtained by multiplying the prescribed left ventricular outflow rate during systole, as shown in Figure 2, by the regurgitant fraction and divided by EROA. The prescribed flow rate from the left ventricle to the aorta was a generic profile based on previous studies [17,18]. To study the effect of different severities on chronic mitral regurgitation murmur, four cases were classified based on the ASE standard severity grading scheme [8]—namely, “mild*”, “mild”, “moderate”, and “severe”—and were simulated. The severity level was altered by changing RF and EROA, while the left atrium model geometry and other factors were kept the same. The “mild*” case and the “severe” case had the same velocity profiles, but their regurgitant fractions and EROAs were different. Maintaining the same anatomy for all the cases enabled us to determine the effect of MR severity without other confounding factors. The key parameters for each case are provided in

Table 1. Details of four simulation cases. $D_{jet}$ is the length of major axis. $U_{jet}$ and $\bar{U}$ are the peak velocity and mean velocity, respectively, of regurgitant jet. $R{e}_{\bar{U}}$ is the Reynolds number in terms of the mean velocity.

	RF (%)	$$\begin{gathered} {\mathit{D}}_{\mathit{j}\mathit{e}\mathit{t}} \\ \left(\mathbf{cm}\right) \end{gathered}$$	EROA $\left({\mathbf{mm}}^2\right)$	${\mathit{U}}_{\mathit{j}\mathit{e}\mathit{t}}$ (m/s)	$\bar{\mathit{U}}$ (m/s)	$\mathit{R}{\mathit{e}}_{\bar{\mathit{U}}}$
Mild*	6.8	0.412	4.4	7.692	1.768	1918
Mild	20	0.556	9.0	11.193	2.747	4016
Moderate	40	0.999	26.1	7.706	1.793	4713
Severe	60	1.702	75.8	3.980	0.927	4152

.

We employed direct numerical simulation (DNS) to solve for the flow field. To handle the complicated geometry, a sharp-interface immersed-boundary ViCar3D was employed, in which a virtual cut-cell method was applied through a Cartesian grid to improve accuracy, conservation property, and computational effectiveness [19,20]. The number of mesh nodes employed in the flow simulation was 194, 230, and 208 in the x, y, and z directions, respectively. The minimum grid size was around $1.2\times {10}^{-2} D$, with $D$ (=3.626 cm) as the major axis of the mitral annulus. The grid spacing was based on resolution requirements that were obtained in previous simulations of atrial [10] and ventricular flows [21,22], which employed the same flow solver and methodology. The time step sizes were $4.6\times {10}^{-2} D_{jet}/U_{jet}$, $5.0\times {10}^{-2} D_{jet}/U_{jet}$, $1.9\times {10}^{-2} D_{jet}/U_{jet}$, $5.8\times {10}^{-2 }{D}_{jet}/U_{jet}$ for “mild*”, “mild”, “moderate”, and “severe” cases, respectively. This provided 14,000 time-steps over the systolic duration, thus providing very high temporal resolution. The simulations were performed on the MARCC–Rockfish cluster, using 288 CPU cores, and required about 15 h wall time for each case. A hybrid second-order upwind and central difference scheme was used for the convection term, and a second-order central difference scheme was used for the diffusion terms. The second-order Adams–Bashforth and Crank–Nicolson methods were employed for the time integration of convection and diffusion terms, correspondingly.

2.2. Murmur Propagation

The stress and velocity fluctuations caused by the generation and propagation of elastic waves in biological tissue were modeled by the generalized Hooke’s law with a Kelvin–Voigt viscoelastic model [23,24,25],

$$\frac{\partial \sigma {{}^{\prime }}_{ij}}{\partial t}+\lambda \frac{\partial u{{}^{\prime }}_k}{\partial x_k}{\delta }_{ij}+\mu \left(\frac{\partial u{{}^{\prime }}_i}{\partial x_j}+\frac{\partial u{{}^{\prime }}_j}{\partial x_i}\right)=0$$

(2)

$${\rho }_s\frac{\partial u{{}^{\prime }}_i}{\partial t}+\frac{\partial \sigma {{}^{\prime }}_{ij}}{\partial x_j}=\eta \frac{\partial }{\partial x_j}\left(\frac{\partial u{{}^{\prime }}_i}{\partial x_j}+\frac{\partial u{{}^{\prime }}_j}{\partial x_i}\right)$$

(3)

where $\sigma {{}^{\prime }}_{ij}$ is the fluctuating elastic stress; $u{{}^{\prime }}_i$ is the fluctuating velocity; $\lambda$, $\mu$, ${\rho }_s$ and $\eta$ are the tissue properties: the first and second Lame’s constants, the density and the viscosity, and ${\delta }_{ij}$ is the Kronecker delta. The first equation describes the normal and shear stresses generated by tissue deformation. The second equation describes the propagation of bulk (compression) and shear waves along with the dissipation due to viscosity.

The model of the human thorax was based on the CT-scan-based image available from the Visible Human Project [26] and is displayed in Figure 3. The upper and lower boundaries of the thoracic cage were assumed to be fixed, and the following boundary conditions were applied to these boundaries:

$$\frac{\partial \sigma {{}^{\prime }}_{ij}}{\partial n}=0, u_i=0.$$

(4)

The tissue–air interface is a traction-free boundary that is expressed as

$$\sigma {{}^{\prime }}_{ij}{n}_j=0, [\lambda \frac{\partial u{{}^{\prime }}_k}{\partial x_k}{\delta }_{ij}+\mu (\frac{\partial u{{}^{\prime }}_i}{\partial x_j}+\frac{\partial u{{}^{\prime }}_j}{\partial x_i})]n_j=0.$$

(5)

The normal stress on the blood vessel wall arising from hemodynamics is the major contributor to the force between the blood flow field and the surrounding tissue. By assuming that the shear stress is negligible, the boundary condition on the wall of the atrium is specified as

$$\sigma {{}^{\prime }}_{ij}{n}_j=P{}^{\prime }{n}_i, [\lambda \frac{\partial u{{}^{\prime }}_k}{\partial x_k}{\delta }_{ij}+\mu (\frac{\partial u{{}^{\prime }}_i}{\partial x_j}+\frac{\partial u{{}^{\prime }}_j}{\partial x_i})]n_j=-\frac{\partial P{}^{\prime }}{\partial t}{n}_i.$$

(6)

where $n_j$ is the normal vector, and $P{}^{\prime }=P-\bar{P}$ is the pressure fluctuation on the vessel wall from the unsteady blood flow. The pressure source was placed within the thorax model at the left atrium wall that would be consistent with anatomy (Figure 3), and the stress-matching boundary conditions were transferred to the thorax model in a manner consistent with this placement.

The above equations for murmur propagation were implemented in a separate immersed-boundary solver [14]. The linear elastic wave equations were discretized by a sixth-order central compact finite difference scheme [27] in space, and time integration was implemented by the four-stage Runge–Kutta (RK4) method [28]. The murmur simulation was carried out on an $82\times 98\times 66$ grid; the minimum grid size was around $0.10 D$, and the time step was roughly $8.0\times {10}^{-2} D/{\mathrm{c}}_0$, where ${\mathrm{c}}_0$ is the speed of sound (compression wave) in blood. The grid cell size and the time step size for these murmur simulations were based on previous heart murmur simulations that were carried out by Seo et al. (2017) [14] and validated against experimental measurements. The murmur simulation required around 10 h wall time on the MARCC–Rockfish cluster using 288 CPU cores. Further details of the linear elastic wave solver, relevant verification, and validation against experimental results can be found in Seo et al. (2017) [14].

3. Results and Discussion

Figure 4 shows the vortex structures obtained from our simulations for the “mild” and the “moderate” cases. As can be seen, the flow is transitional in nature and very complex. As the vortex structures from the regurgitant jet hit the lumen of the atrium, they generate pressure fluctuations that become the source for the murmur that propagates through the thorax.

We now describe the characteristics of the murmur-associated elastic wavefield through the thorax. The magnitude of instantaneous velocity fluctuation ($U^{{}^{\prime }}=\sqrt{u{{}^{\prime }}^2+v{{}^{\prime }}^2+w{{}^{\prime }}^2}$), average normal stress ($P{}^{\prime }=({\sigma }_{xx}+{\sigma }_{yy}+{\sigma }_{zz})/3$), and shear stress (${\sigma }_{xy}$) in the moderate MR case are plotted in Figure 5, on a transverse plane crossing the point at which the root-mean-squared (rms) surface pressure fluctuation is maximum on the left atrium wall. The velocity fluctuation ($U^{{}^{\prime }}$) represents the tissue vibration due to the elastic waves generated by the mitral regurgitation jet, and it is computed relative to the equilibrium state when there are no elastic waves. In Figure 5, not only does the left atrium reveal notable velocity fluctuation, as is expected in the left lung, some weaker fluctuations are found in the upper lobe of the right lung as well. An isosurface of relatively acute velocity fluctuation and part of the skeletal system colored by the magnitude of velocity fluctuation is presented in Figure 6. By referring to the tissue’s physical properties and anatomical position of the fluctuating conduit extending from the heart to the right lung, the four main pulmonary veins are observed as the major media for the wave propagating into both the left and the right lungs. Due to high density and compression/shear modulus, the blood vessel wall and the blood flowing inside them fluctuate more intensely than the lungs filled by air. In general, the normal stress and shear stress are concentrated in the heart, aorta, and pulmonary veins, as well as superior vena cava. They propagate the wave to the sternum and vertebral column, and then to the costal cartilages, ribs, muscles, and skin of the chest wall. It is also worth noting that the pressure fluctuation is distributed asymmetrically between the left and right lungs, with the magnitude of pressure fluctuation in the left lung being considerably larger than that in the right lung. In addition, the upper lobe and more than half of the lower lobe of the left lung have positive pressure fluctuation when negative pressure fluctuation dominates in the right lung. Due to the high Poisson’s ratio ($\upsilon$) of the lung [29], Young’s modulus (E) can be related to the shear stress (G) by $E=2G\left(1+\upsilon \right)$, assuming the lungs are isotropic and homogeneous tissue [30], which could explain why the magnitude of shear stress in the lungs is negligible, compared with the magnitude of the normal stress.

The instantaneous x-direction (which is normal to the chest precordium) velocity fluctuation in the moderate case is presented in Figure 7. Notably, the fluctuations are observed to propagate through the front wall of the thorax intensively on each side, resulting in two apparent independent regions of sharp fluctuation in x-direction on the chest wall due to the complex composition of various anatomical structures.

Figure 8a displays the root-mean-squared (rms) wall-normal velocity fluctuation magnitude over the systolic period on the anterior surface of the chest for the moderate case. Additionally, included for reference are the five standard points on the precordium where cardiac auscultation measurement is usually made: the aortic, pulmonic, ERB’s point, tricuspid, and mitral areas. Figure 8b illustrates the wall-normal velocity over time at the five auscultation sites. We should note that there are two local maxima of fluctuation at points 1 and 3, and this is consistent with Figure 7 since the wall-normal magnitude is dominated by the x-direction component. Although points 2 and 4 are local maxima as well, they have a lower fluctuating kinetic energy compared with the other two points. The simulation results for different severities exhibit very similar patterns. This would seem to contradict the well-established notion in clinical practice that mitral regurgitation murmurs are best heard near the mitral point. However, previous studies actually support our simulation results that the overall intensity of MR murmur is not the highest at the mitral point. Sh-Hussain et al. (2013) used an electronic stethoscope, which is more sensitive, especially to high-pitched regurgitation murmur [31], than a traditional device to collect training and testing data from each auscultation site [32]. Their results showed that the mitral area is the worst auscultation point for mitral and tricuspid regurgitation based on a comparison of identification accuracy between the aortic, pulmonic, tricuspid, and mitral areas [32]. Jost et al. (2000) found that only 5 (18%) of the 22 patients with mitral regurgitation have measured intensities of the murmur that are the largest at the mitral point [33]. However, in practice, the problem encountered by clinicians is that measurement of MR-induced murmurs at the aortic or pulmonic points may be contaminated by murmurs generated by the aortic or pulmonary valves [34,35,36]. In contrast, even though the overall murmur intensity of the MR murmur is lower at the mitral point, it is likely less contaminated by murmurs from the other valves. Thus, in clinical practice, the mitral point is favored as a point on the precordium for detecting mitral-valve-related murmurs.

Following clinical practice, we focus the rest of our discussion and analysis on point 4, i.e., the mitral point. The frequency spectra of the acceleration data for different severities are presented in Figure 9. The amplitude increases slowly and then decreases gradually in all cases, and the spectra can, therefore, be approximated by fitting the spectral data to two piecewise linear functions. Overall, we found that the slopes of the ascending portions have an average slope of 5/3, whereas the descending portions have a slope that is nearly equal to −4. It is important to note that the maximum amplitude and the overall amplitude increase with the severity and regurgitant fraction rather than the Reynolds number since the regurgitant fraction is not linearly correlated with the Reynolds number. The frequency at which the slope of trend lines intersect is commonly called the “break frequency”, which has been proposed as a diagnostic metric for the severity of valve disease when using electronic auscultation [37]. The Strouhal number ($St=f{D}_{jet}/U_{jet}$) in terms of break frequency in mild cases fall in the range of 0.5 to 1. The moderate and severe cases break at Strouhal numbers around 2 and 7, respectively. Figure 10 indicates that the break frequency Strouhal number grows rapidly with the regurgitant fraction, and we can fit an exponential curve to this date with an R-squared value of 0.9937. Thus, the current hemoacoustic simulations indicate a clear quantitative correlation between the spectral characteristics of the murmur that would be detected at the chest surface and the severity of the mitral regurgitation. Correlations such as these could be used to diagnose the severity of mitral regurgitation-induced murmurs.

It is worthwhile to point to some of the key limitations of the current study. First, the left atrium model does not include a detailed representation of the mitral valve leaflets. In practice, the leaflets can appear in many different prolapsed conformations, which could interact with the regurgitant jet and, in turn, affect the location of jet impact on the left atrium lumen and the resultant murmur. The inclusion of the leaflets would, however, likely require a fluid–structure interaction model, which would significantly complicate the modeling procedure and analysis. The second limitation worth noting is that the left atrium wall was assumed to be rigid, and therefore, the interplay between the deformable atrial wall and the hemodynamics of the sound was neglected. Such simulation with fluid–structure interaction model has been conducted by Feng et al. (2019) [38]. The third limitation is the use of a single patient-specific model in our study, which did not allow us to include subject-to-subject variabilities in anatomy. Future studies can incorporate other factors not considered here.

4. Conclusions

A computational hemoacoustic modeling approach was used to simulate the murmurs generated by mitral regurgitation with different degrees of regurgitation severity. The model employed unsteady Navier–Stokes model-based simulations of the regurgitant jet, and this was coupled with a high-fidelity model of the murmur propagation through a thorax model based on imaging data. To our knowledge, this is the first such model of mitral regurgitation murmur. Two key conclusions are drawn from the current study: first, the murmur intensity for the MR murmur was not the highest at the mitral point, which is typically used for the detection of this condition by clinicians. Murmur intensity was highest at the aortic and pulmonic points, and with the use of electronic stethoscopes and more sophisticated signal analysis techniques, it might become possible to disambiguate the MR murmur signal at these points from other murmurs and provide more accurate detection of MR. Second, the simulations showed that there is a clear quantitative correlation between the severity of chronic MR and the break frequency, and this might again provide a target for automated detection and stratification of chronic MR.