Multichannel seismocardiography: an imaging modality for investigating heart vibrations

Objective: Seismocardiography is the measurement of vibration waves caused by the beating heart with accelerometer(s) placed on the chest. Investigating the nature and the behavior of these vibration waves, by comparing measurements from multiple sites, would help to understand the heart’s mechanical contraction activity. Approach: Using newly designed multichannel seismocardiogram equipment, it was possible to investigate the vibration waves with 16 three-axis sensors. The equipment performed well with highly precise synchronization rate over 10 min, linear frequency response and high signal quality. The vibration waves were analyzed using the sagittal axis, a single cardiac cycle and focusing on four fiducial points. Two of the fiducial point where the negative and positive peaks associated with aorta valve opening, along with peaks associated with aorta valve closing. Main results: The respective average centers of mass of the four fiducial points in 13 subjects were at (frontal axis: 35 mm, vertical axis: 5 mm), (31, 6), (26, 24), and (4, −2), relative to the Xiphoid Process. Similar patterns among the subjects were identified for the propagation of the waves across the chest for the four fiducial points. Significance: The multichannel seismocardiogram equipment successfully revealed a general pattern present in chest surface vibration maps.


Introduction
Seismocardiography (SCG) is a method of measuring and analyzing the vibrations of the chest produced by the beating heart . The vibrations are related to cardiovascular motion and hemodynamics, including opening and closing of heart valves (2017( , Sørensen et al 2018. These cardiovascular events produce mechanical waves that propagate through the body and can be obtained as a SCG by measuring the vibrations on the body surface using accelerometers (Zanetti and Tavakolian 2013). Studies have shown the relationships between markers in the SCG and timing of specific cardiac events (Crow et al 1994). It has been demonstrated that these markers reflect pathological changes in patients with heart failure (Pouyan et al 2017), cardiac resynchronization therapy (Marcus et al 2007), and coronary artery diseases (Wilson et al 1993). Zanetti and Tavakolian (2013), proposed a multichannel SCG (mchSCG) as an imaging modality to improve our understanding of mechanical wave propagation and its effect on the SCG markers of cardiac function (Zanetti and Tavakolian 2013).
Multi-channel recording of heart sounds was reported as early as 1982 by (Okada 1982). Later Cozic et al (1998), proposed a modern microphone array method to investigate optimal auscultation positions at the chest (Cozic et al 1998). To investigate the relationships between blood flow and force at different sites of the chest and extremities, Rienzo et al first designed a mchSCG system in 2016 (Di Rienzo et al 2020). The system records from up to 12 SCG channels in combination with photoplethysmography. Lin et al (2018), reported that markers from multiple recording sites (channels) improved correlation to echocardiographic events (Lin et al 2018). Laser Doppler vibrometry (LDV) has also been used to record mechanical vibration of the chest (Morbiducci et al 2007). Heart beats where recorded asynchronously using LDV from multiple sites to obtain vibration maps, where the signals were synchronized using the R-peaks of simultaneously recorded ECG (Munck et al 2016). Shirkovskiy et al (2018), found that using airborne ultrasound waves, chest vibrations could be recorded at a high spatial resolution (Shirkovskiy et al 2018).
Limitations of these systems are either a low spatial resolution (⩽5 channels) (Di Rienzo et al 2016, Lin et al 2018, asynchronous sampling (Munck et al 2016), or low temporal resolution (⩽157 samples per second) (Shirkovskiy et al 2018). To gain better insight in cardiac-induced chest surface vibrations, we designed and validated an accelerometer-based mchSCG system with 16 three-axis accelerometers with a synchronized maximum sample rate at 1000 Hz.

System architecture
The system was designed to comprise 16 three-axis accelerometer sensor based sensorboards, a 3-lead electrocardiogram (ECG) and one additional channel for external equipment. The external equipment channel is practical for comparing the sensors to, e.g. other accelerometers or a respiration belt. The motherboard ensures synchronization of all sensor modalities. To ease sensor placement and to reduce noise sensitivity two sensorboard chains were designed to host up to eight sensorboards connected on a ribbon cable. Each sensorboard chain and the ECG/external equipment was connected to one of three microcontrollers. The microcontrollers were located on the motherboard and were responsible for collecting and transferring data to the computer. The computer acquisition software was designed in LabVIEW ® (National Instruments, US), to collect data in real-time. Collected data was saved to Matlab for data analysis. The system architecture is illustrated by figure 1.
The ADS1298 ® analog front-end (Texas Instrument, US) was used for the ECG and external equipment. Additionally, the ADS1298 supplied the 16 sensors with a synchronization signal. The ADS1298 transferred data to the microcontroller using the SPI communication protocol.
To obtain high data rates, the system used FT51 ® (FTDI Chip, US) microcontrollers, that can transfer directly using the USB protocol and has Direct Memory Access (DMA) for the SPI data transfer.

Sensors and communication
The accelerometer ADXL355 (Analog Devices, US) was used on the sensorboards. To record the small details of a SCG signal, the measurement had to have high resolution (3.9 µg bit −1 for the ADXL355), low noise (25 µg/√Hz), and linear frequency response at lower frequencies (0-1000 Hz) (Analog Devices 2016). In addition, the ADXL355 has a built-in synchronizable ADC, and uses SPI communication, which was important for the ribbon cable sensor chain design. The 10 wires ribbon cables included SPI communication, synchronization, power supply, and sensor addressing. Since the SPI communication protocol is not designed for off board communications, an isoSPI solution was implemented with the   LTC6820 communication interface (Linear Technology, US) (Linear Technology Corporation 2012). The isoSPI solution increases the signal stability and maximum cable length for the sensors (10 m), but, limits SPI bandwidth to 1 MHz, which limited the number of sensorboards on a chain to eight. The system can measure at three sampling frequencies, 250 Hz, 500 Hz, and 1 kHz. Sensorboards were placed inside 3D printed plastic casings as illustrated by figure 2.
The dimension of the finished sensor, 30 × 30 × 18 mm (length × width × height), resulted in 9 cm 2 sensor area. Taebi et al (2019) found that a larger sensor area reduces the SCG dependency of sensor placement.

System test
Two experiments were conducted in order to determine how the system performed with respect to synchronization and frequency response. In both experiments, all sensors were placed on a linear resonant actuator (LRA) along with a reference accelerometer (model 1521 (Silicon Designs Inc 2018)) connected to the external equipment channel of the mchSCG equipment. The experimental setup is shown in figure 3.
The first experiment investigated the sensors' synchronization performance by oscillating the LRA at 20 Hz for 10 min. The experiments were conducted for all three sampling frequencies and only the axis perpendicular to the LRA plate was used. The latency of the cross correlation between the reference accelerometer and the 16 sensors was used to determine the synchronization characteristics between the ECG channel and the SCG. The cross correlation was performed with 5 s windows and the mean and standard deviation latency of the windows was used as interpretation of the synchronization characteristics. The same procedure was used to compare the synchronization characteristics of any one sensor to the other 15 sensors.
The second experiment investigated the sensors' frequency response by running the LRA with a maximum magnitude response of ±2.5 g while doing a logarithmic sweep from 5 to 800 Hz with a sweep rate at 0.6 decade min −1 . Experiments were performed with sampling rates of 250, 500, and 1000 Hz. Since the LRA did not have a linear frequency response, the reference accelerometer measured at a 2 kHz sampling rate as reference. The frequency response was calculated based on a moving root mean square (mRMS) of the reference compared with a mean of all sensors mRMS.

Quantification of chest wall vibrations
Measurements were conducted on 13 male subjects with mean age 27 ± 4 years, and mean BMI of 24 ± 4 kg m −2 . The study was approved by the North Denmark Region Committee on Health Research Ethics: N-20 170 008. (Munck et al 2020).
The sensors were placed in a 4 by 4 grid formation on the volunteer's chest, along with a 3-lead ECG and a respiration belt; see figure 4. The sensor coordinate system is defined as in figure 5. The sensor grid was placed such that the sensor in location 3,2 was placed at the Xiphoid Process, and such that the centers of adjacent sensors were spaced at approximately 40 mm. Subjects were placed in a supine position and 3-min recordings were obtained with a 500 Hz sampling rate.

Pre-processing
Each SCG signal was filtered individually with zero-phase highpass and lowpass-filters of 0.05 Hz and 90 Hz, respectively (Sørensen et al 2018). For the current analyses a single cardiac cycle was extracted from the exhalation phase, where the latter was identified using the respiration belt signal.
To assess the signal quality, the signal to noise ratio (SNR) for each of the 16 SCG signals, measured from each of the 13 subjects, was calculated based on the extracted cardiac cycles. SNR was calculated by comparing the signals variance (σ) filtered to 90 Hz (f ) with the noise density (n) of the ADXL355: The ADXL355 has a noise density of 25 µg/√Hz.

Chest surface vibration colormaps
To demonstrate the potential of mchSCG as a tool for chest surface vibration analysis, colormaps of the z-axis (sagittal) vibrations were constructed at the times of occurrence of four high-amplitude fiducial points, two systolic and two diastolic fiducial points. Fiducial point was defined as the negative and positive peaks that align with aorta valve opening, for the systolic points, and with aorta valve closing for the diastolic points. The fiducial points were defined as AO hill , AO valley , AC hill , and AC valley . The systolic fiducial point AO hill corresponding to the opening of the aortic valve, was defined as the G s point in Sørensen et al (2018). The AO valley fiducial point, the F s point in Sørensen et al (2018), corresponding to the isovolumetric phase, was identified as the minimum point before the AO hill . The diastolic fiducial point AC valley was the C d point in Sørensen et al (2018), defined as the minimum of the first negative wave after the onset of aortic closure. AC hill was the maximum of the next positive wave, defined as D d in Sørensen et al (2018). All fiducial points were manually annotated at channel (3,2) the channel from the Xiphoid Process. Colormaps show the raw amplitudes of acceleration at the 16 sensor locations. Each set of colormaps belonging to an annotation has a color range of the maximum absolute acceleration of the given interval. Colormaps were composed with contours to represent gradients between measuring points and squares indicating sensor location on the maps. To trace the temporal development in chest surface vibrations, colormaps were constructed in a sequence starting and ending 8 ms before and after the fiducial point.

Quantification of colormaps
To qualify spatial movement of acceleration waves on the chest surface, the center of mass of each colormap was calculated and illustrated as a magenta point on the colormaps. The center of mass was calculated only Figure 6. Results from the frequency response test with the mchSCG equipment set to sample at the three sampling frequency settings. The response is a comparison between the accelerometers model 1521, from Silicon designs, and ADXL355, from Analog Devices.
for negative values for the AO valley and AC valley , and the positive values for AO hill and AC hill . The center of mass was calculated according to equations (2) and (3): where i is each of the N = 16 sensors, m is the amplitude, x and y are the sensor locations, and M is the mean amplitude. Before the center of mass calculations, a 3-dimensional cubic interpolation was performed on the data to increase the spatial resolution from 40 mm to 1 mm. Also, a 15% amplitude threshold was applied in order to include information about the hill and valley formations in the center of mass calculations only. For comparison the average colormap was calculated from normalized colormaps of all the subjects. This average colormap along with the center of mass was used for a cross subject interpretation of the vibrations, both in the temporal and spatial domain.

Results
The results from the System test and mchSCG in human validation methods are presented below in sections A and B, respectively.

System test
The latency between the ECG channel and the 16 sensors was 16 (SD: 0) ms, 8 (SD: 0) ms, and 4 (SD: 0) ms, for sampling frequencies 250 Hz, 500 Hz and 1 kHz, respectively. This means that there is a latency of 4 samples independent of the sampling frequency. The latency between the sensors was 0 (SD: 0) ms, for all sampling frequencies.
The use of 16 sensors meant that the system had a limited bandwidth of 48.000 samples per second or 1.2 Mbit s −1 .
The second experiment revealed a linear frequency response for frequencies below 20 Hz and a below −15 dB antialiasing attenuation for all sampling frequencies, as illustrated by figure 6. The stopband attenuation of −3 dB was found at 60 Hz, 116 Hz, and 174 Hz for the three sampling frequencies, respectively.
The signal quality assessment revealed an expected high SNR of 25.8 (SD: 3.0) dB for the 208 signals. The system conditions made it possible to analyze a single cardiac cycle without the use of other methods to enhance the signal quality such as segmentation.

Multichannel SCG in humans
Recordings where successfully collected from the 13 subjects. Examples of acquired signals after pre-processing are illustrated by figure 7, for one subject, and for all subjects in the supplementary material by figures S1-S13 (available online at https://stacks.iop.org/PMEA/41/115001/mmedia).
For each subject a set of colormaps was composed at each fiducial point. Colormaps of one subject are shown in figure 8, and colormaps for all subjects can be seen in the supplementary material in figures S14-S26.  Average colormaps of the 13 subjects and calculated center of mass, are illustrated in figure 9 and table 1. and includes the average location of center of mass from the 13 subjects.
A comparison between all the subjects was done based on average colormaps from the annotations and calculated center of mass, as illustrated by figure 8. Figure 9 and table 1 demonstrate how the location of center of mass is different among the four fiducial points. At the time of the fiducial point AO valley the average center of mass was 35 mm to the left of and 5 mm above the Xiphoid Process and at the time of AO hill the center of mass was at 31 mm to the left of and 6 mm  Average coordinates of center of mass in the 4 by 4 sensor grid, read as coordinate mean(±SD). All coordinates are in mm and relative to the Xiphoid Process at (0 mm, 0 mm). X-axis is distance in the transversal plane, with positive as left side, and Y-axis is distance in the sagittal plane, with positive above the Xiphoid Process.
above the Xiphoid Process which corresponds to the level of the 5th intercostal space at sternum. The AC valley center of mass clusters in two groups, one at the Xiphoid Process and one at a higher location near 40 mm to the left of and 40 mm above the Xiphoid Process, corresponding to the 3rd intercostal space at the left parasternal line. AC hill had an average center of mass at the Xiphoid Process. In addition, temporal trends can be observed: • AO valley : the negative wave starts at the left side of sternum and moves closer to the sternum at 0 ms and then disperses in all directions. • AO hill starts more diffuse on the chest, but at 0 ms the center of mass is primarily focused left of the Xiphoid Process, followed by the positive waves moving mainly down and to the left.

Discussion
In the current study we developed and validated multichannel SCG equipment for recording of human chest wall vibrations.

System validation
The system test showed that the mchSCG equipment performed as expected with a 0 ms variation in synchronization over the 10-min test period and with a linear frequency response below the Nyquist frequency. The frequency response showed that at sampling frequency of 500 and 1000 Hz includes most of the energy from the heart valve sounds, since most of the energy of heart sounds is below 150 Hz (Sakai et al 1971).
Human trials were conducted with a 500 Hz sampling rate because it was enough to obtain the desired bandwidth. The analysis revealed that the temporal resolution also was enough to investigate propagation of the vibration waves.
Additionally, a signal to noise evaluation of the human trial measurements showed a 25.8 (SD: 3.0) dB ratio for all channels and subjects. This high SNR was likely due to the high-end accelerometers and the short signal pathway, where the digital conversion is performed on the sensor.

Interpretation of systolic chest wall vibrations
The average colormaps suggest that AO valley and AO hill peak at the 5th intercostal space approximately 30-40 mm left of the midsternal line. This point is anterior to the left cardiac ventricle. Thereby it is likely that left ventricular activity is the driving force in the AO valley and AO hill events. The downward movement of the chest wall related to the AO valley can be explained by early systolic shortening where both the short axis and long axis diameter shorten (Remme et al 2008). Early systolic shortening is then interrupted by the closure of the mitral valve closing, which leads to the isovolumetric phase where the apex of the heart moves toward the base (Remme et al 2008). In that period the left ventricle becomes more circular and the lateral left ventricular wall will move outward toward the chest wall causing the AO hill event. After aortic opening the left ventricular wall will again continue its downward movement which corresponds to a downward wave post AO hill. This interpretation is in line with the experimentally determined LV diameter acceleration found by Remme et al (2008), but delayed compared to the findings by Sørensen et al (2018) who found that the onset and not the minimum of the AO valley coincide with the mitral valve closure and the AO hill maximum coincided with the aortic opening.
AO valley and AO hill occur during the isovolumetric contraction. During the isovolumetric contraction the left ventricle starts a clockwise rotation followed by a rapid counterclockwise rotation (Sengupta et al 2008). This rapid counterclockwise rotation could explain the propagation of AO valley and AO hill . Here the negative AO valley propagates opposite of the rotation, to the right, and the positive AO hill propagates in the same direction as the rotation, to the left.

Interpretation of diastolic chest wall vibrations
The average colormaps suggest that the AC valley was located at two spatially divided locations, one minimum located inferior to the Xiphoid Process and one minimum located more superior at the 2nd intercostal space at the left parasternal line. The AC valley at the inferior location could again be explained by LV movement. According to Remme et al (2008) lengthening of the myocardium starts before the aortic closure which will lead to an outward movement, but this lengthening is stopped or decelerated by the start of the isovolumetric relaxation period started by the aortic closure, leading to a downward acceleration, the AC valley . The superior AC valley could be explained by a pressure drop in the aortic arch, caused by the drop in left ventricular pressure and coinciding with backflow and aortic closure.
The AC hill peaks close to the Xiphoid Process. A first obvious explanation for an outward acceleration like the AC hill could be the re-initiation of the myocardial lengthening at the time of mitral valve opening, but this is not consisted with results demonstrating that the AC hill occurs much earlier (~50 ms) than the opening of the mitral valve (Sørensen et al 2018). Rather, AC hill could be the result of a water-hammer effect from the aortic closure.
AC valley and AC hill occur during the isovolumetric relaxation. During the isovolumetric relaxation the left ventricle turns from a rapid counter-clockwise rotation to a rapid clockwise rotation (Sengupta et al 2008). This rapid clockwise rotation could explain the propagation of AC valley and AC hill . Here the negative AO valley propagates opposite of the rotation, to the left, and the positive AC hill propagates in the same direction as the rotation, to the right.
Regarding the standard deviation, we mainly have smaller variations then the sensor size. However, for the fiducial points AO valley , the center of mass variates by ±38 mm at 0 ms on the y axis, which is due to the two earlier described epicentres of AO valley .

Limitations
The size of the sensors works as a spatial filter, where larger sensors dampen shorter wavelengths. Reducing the size of the sensors will make it possible to increase the spatial resolution of the map and the possible number of sensors on the chest. Shirkovskiy et al (2018) used a method with overlapping spatial measuring points along with a higher resolution and a similar measurement area. Studying these high spatial resolution videos will be useful in determining requirements for the spatial resolution. In addition, understanding how sensor size dampens shorter wavelengths and how it influences the vibration measurements will be necessary for further analysis of the vibration waves propagation. We know from seismology that the types of vibration waves can be described through the use of vibration measurements in all three-axis (Shearer 2009).
Understanding the types of vibration waves will help focus on the elements of the SCG that relate to certain cardiac events.
A reference accelerometer was used for the frequency response experiment because the LRA did not have a linear frequency response. However, the reference accelerometers linear frequency response was described as a less than 3 dB variation from 0 to 400 Hz (Silicon Designs Inc 2018), which accounts for the 3 dB variation seen above 20 Hz in the passband of the frequency response.

Future potential of technology
Further comparison of how spatial propagation of the waves relates to echocardiographic measures such as left ventricle rotation, could be used to estimate physiological variables such as preload, afterload, and contractility (Sengupta et al 2008). Furthermore, we believe that mchSCG would become a useful modality for optimizing cardiac resynchronization therapy and help determine who would benefit from the treatment.
For future reference through this study we have been made aware of multiple spatial sensor parameters regarding sensor size, spacing and placement. Understanding how these spatial sensors parameter influences the obtained information is vital if we are to better understand (SCG) as a modality. This should be done by comparing the mchSCG method with similar methods such as LDV and airborne ultrasound waves (Munck et al 2016, Shirkovskiy et al 2018.

Conclusions
Through the design, build and application of mchSCG equipment we demonstrated that cardiac events introduce complex vibration patterns on the chest wall and we showed their spatial locality. There is evidence showing a relationship between the patterns of the surface waves and individual cardiac events. We believe that continued investigation of these patterns will lead to further understanding of the main components in the SCG signal.