Bathtub ECG as a Potential Alternative to Light Stress Test in Daily Life

: The exercise stress test (EST) is a common procedure to evaluate cardiovascular functions. However, the EST is not suitable for daily use, is sometimes risky, and even accompanies fatal incidents of myocardial infarction, arrhythmia, and sudden death during the test. The aim of this study was to evaluate heart rate variability (HRV) behaviors in the EST and during bathing, and to explore if daily bathing can serve as a potential alternative means of performing the EST. Electrocardiogram (ECG) signals were collected from 10 healthy subjects during the EST and bathing test (BT). The EST follows the modiﬁed Bruce protocol. ECG collection in the BT was conducted at ﬁve water temperatures ranging from 37 to 41 degrees Celsius ( ◦ C); each BT lasted 15 min. Twenty-three HRV features were used to group different bathing conditions corresponding to the EST stages using the Voronoi diagram method in terms of HRV behaviors. In all equivalent EST stages of BTs at the ﬁve water temperatures, the low stage, medium stage, and high stage account for 17.86%, 52.86%, and 29.29%, respectively. The results show that higher water temperatures and longer bathing durations in BT correspond to higher stages in the EST. The BT at the most severe condition of 41 ◦ C and 15 min corresponds to a high EST stage in terms of HRV behavior. The results suggest that daily bathing can serve not only for healthcare monitoring but also as a reference for an at-home alternative to the EST.


Introduction
The exercise stress test (EST) is a safe procedure as a common cardiological test that doctors use to diagnose coronary artery disease (CAD). The EST has been developed to apply widely so far, including evaluation of the anatomic and functional severity of CAD, evaluation of exercise-related symptoms, and assessment of the response to medical interventions [1]. In addition to the diagnostic value of screening cardiovascular disease, EST also has a significant prognostic value. Goldschlager et al. confirmed that EST can be used for diagnosing myocardial ischemia and CAD based on ESTs on 269 patients and 141 normal subjects [2]. Théroux et al. implemented a limited EST before hospital discharge of 210 patients after acute myocardial infarction and proved that this test is safe and can predict mortality in the subsequent year [3]. EST has also been used for valvular heart disease patients to quantify disability, reproduce exercise-induced symptoms, and assess responses to medical and surgical interventions [4]. There are multiple protocols for EST; the Bruce protocol is one of the most popular protocols [5]. It has been 59 years since the introduction of the Bruce protocol in 1963, and the protocol has developed from the first version of four stages to the modified Bruce protocol [6]. The EST is generally safe, but there have been reports of myocardial infarction, arrhythmia, and deaths during tests, expected to occur once every 2500 tests [7,8]. The EST should be implemented under the supervision of a physician [9]. Therefore, EST is not suitable for daily application.
In Japan, bathing is popular with many people in daily life and has become an almost routine activity [10]. Multiple studies on bathtub electrocardiograms (ECGs) exist. Kwatra et al. proposed a different method from the conventional bioelectric measurement for recording ECG signals from a home bathtub setup [11]. Tamura et al. designed a bathtub heart rate monitor and applied it in a fully automated monitoring system of health status at a pilot house [12,13]. Ogawa et al. completed identifying the bathtub ECG signal using a neural network with wavelet transform [14]. Mizukami et al. confirmed that a bathtub ECG is a suitable method to follow up those patients with a pacemaker implanted [15]. Xu et al. confirmed that heart rate variability (HRV) using a bathtub ECG is impacted by water temperature during bathing [16]. In our previous study, a stress index (SI) was calculated to evaluate daily stress using a bathtub ECG and using a convolutional neural network personal identification with a bathtub ECG [17,18]. A bathtub ECG is confirmed as convenient for collection and suitable for daily healthcare monitoring.
Exercise is one of many tolerable physiological stresses that can induce cardiovascular abnormalities that are not present at rest. In an EST, the most common exercise-induced stressor is an exercise on a treadmill or a bicycle ergometer. Sympathetic discharge is maximal and parasympathetic stimulation is withdrawn during strenuous exertion [19]. In bathing, the stressful load is from water pressure and thermal irritation. Gorman et al. identified sympathetic activation and parasympathetic withdrawal during thermal stress in a baboon [20]. There are potential shared mechanisms of action between thermal stress and exercise stress.
HRV is a reliable reflection of the many physiological factors modulating the cardiac rhythm [21]. Karthikeyan et al. used short-term ECG and HRV signals of the Stroop color word-based stress-inducing task to detect stress and achieved an overall average classification accuracy of 91.66% and 94.66% using the probabilistic neural network algorithm and k-Nearest Neighbor (KNN) algorithm classifiers [22]. Munla et al. used a support vector machine with radial basis kernel function to complete driver stress detection based on HRV analysis and achieved an accuracy of 83% [23]. Ferdinando et al. used KNN to implement emotion recognition based on new features from short ECG signals and HRV features and offered an approach to emotion recognition based on short ECG signals [24]. Orphanidou et al. proposed a quality assessment system based on wavelet entropy measurements of the HRV signal for wearable sensors [25]. Zarei et al. proposed an algorithm using new features extracted from HRV and ECG-derived respiration signals for the detection of obstructive sleep apnea [26]. In this study, HRV analysis is implemented to extract features of the EST and bathing test (BT) using the features to investigate the HRV behaviors of the BT and EST for similarity.
In this study, an alternative to the EST using HRV analysis of the ECG during bathing in daily life is introduced. This method achieves the same equivalent EST stages in terms of HRV behaviors by controlling the water temperature and bathing duration during bathing. Bathtub ECG as a light stress test is implemented during bathing by a natural way and without uncomfortableness in daily life. This paper is organized as follows: Section 2 describes the method for data collection, processing, analysis, and evaluation. The results are presented in Section 3 and discussed in Section 4. Finally, the conclusions are expounded on in Section 5.

Materials and Method
The scheme of the implementation procedure of this study is shown in Figure 1. The ECG signals are collected during the EST and BT. Extraction of HRV features is implemented after raw signal denoising and R wave peak detection on preprocessing. The principal component analysis (PCA) is used for data dimensionality reduction of HRV features. The data after dimensionality reduction via the Voronoi diagram are used to evaluate the equivalent EST stages in terms of HRV behaviors of the BT. EST Stages 1-7 are grouped to low stage (Stages 1 and 2), medium stage (Stages 3 and 4), and high stage (Stages 5-7) for investigating and evaluating the equivalent EST stage of the BT.

Data Collection
Data collection consists of two parts in this study, including the BT and the EST. Ten healthy male subjects were enrolled for this data collection. All the subjects were university students in their twenties to thirties without cardiac disease. The data collection procedure was explained to subjects, and written informed consent was obtained from all subjects before data collection. Basic information of subjects, including age, gender, body weight, height, and health status, was recorded in the experiment diary before data collection.

Data Collection
Data collection consists of two parts in this study, including the BT and the EST. Ten healthy male subjects were enrolled for this data collection. All the subjects were university students in their twenties to thirties without cardiac disease. The data collection procedure was explained to subjects, and written informed consent was obtained from all subjects before data collection. Basic information of subjects, including age, gender, body weight, height, and health status, was recorded in the experiment diary before data collection.

Bathing Test
A bathtub ECG monitoring system was built for data collection of the BT, as shown in Figure 2. This system consists of three components: bathtub electrodes, an ECG monitoring device (OpenBCI Ganglion system, OpenBCI, Inc., Brooklyn, NY, USA), and a computer (MacBook Pro, Apple, Cupertino, CA, USA). There are four electrodes attached to the inner wall of the bathtub. The electrodes are hemispheres (stainless steel, Dragonmarts Co., Ltd., Hong Kong, S.A.R. of China) with diameters of 25 mm and placed on both sides at the same level as limbs. The electrodes are connected to the ECG monitoring device with shielded cables. OpenBCI Ganglion Board is the biosensing hardware for providing imaging and recording of electromyography, ECG, and electroencephalography signals. The board communicated wirelessly to a computer via Bluetooth. The board's sample rate is limited by the Bluetooth bandwidth. Bluetooth transmission rate is limited to 200 Hz. Therefore, data were sampled at 200 Hz on each of the four channels and recorded on a computer. ECG signals were collected during subjects' bathing in a relaxed position, as shown in Figure 2. In previous studies, the impact of different water temperatures on HRV during bathing was examined and confirmed. With the increase in water temperature, some HRV features displayed different degrees of impact, such as decrease or increase [16,27,28]. In daily practice, the water temperature mostly is set between 38 • C and 42 • C during bathing [29]. Some subjects could not complete taking a bathing at a water temperature of 42 • C during rehearsal. Therefore, each subject completed five BTs at different water temperatures: 37 ± 0.5 • C, 38 ± 0.5 • C, 39 ± 0.5 • C, 40 ± 0.5 • C, and 41 ± 0.5 • C. The initial water temperature was set before data collection. After the subjects entered the bathtub, a thermal insulation film was used to cover the bathtub to reduce the loss of water temperature during bathing. The duration of each BT was controlled to 15 min. The lead II ECG signals were collected and used in this study.  Figure 2. This system consists of three components: bathtub electrodes, an ECG monitoring device (OpenBCI Ganglion system, OpenBCI, Inc., Brooklyn, NY, USA), and a computer (MacBook Pro, Apple, Cupertino, CA, USA). There are four electrodes attached to the inner wall of the bathtub. The electrodes are hemispheres (stainless steel, Dragonmarts Co. Ltd., Hong Kong, S.A.R. of China) with diameters of 25 mm and placed on both sides at the same level as limbs. The electrodes are connected to the ECG monitoring device with shielded cables. OpenBCI Ganglion Board is the biosensing hardware for providing imaging and recording of electromyography, ECG, and electroencephalography signals. The board communicated wirelessly to a computer via Bluetooth. The board's sample rate is limited by the Bluetooth bandwidth. Bluetooth transmission rate is limited to 200 Hz. Therefore, data were sampled at 200 Hz on each of the four channels and recorded on a computer. ECG signals were collected during subjects' bathing in a relaxed position, as shown in Figure 2. In previous studies, the impact of different water temperatures on HRV during bathing was examined and confirmed. With the increase in water temperature, some HRV features displayed different degrees of impact, such as decrease or increase [16,27,28]. In daily practice, the water temperature mostly is set between 38 °C and 42 °C during bathing [29]. Some subjects could not complete taking a bathing at a water temperature of 42 °C during rehearsal. Therefore, each subject completed five BTs at different water temperatures: 37 ± 0.5 °C, 38 ± 0.5 °C, 39 ± 0.5 °C, 40 ± 0.5 °C, and 41 ± 0.5 °C. The initial water temperature was set before data collection. After the subjects entered the bathtub, a thermal insulation film was used to cover the bathtub to reduce the loss of water temperature during bathing. The duration of each BT was controlled to 15 min. The lead II ECG signals were collected and used in this study.

Exercise Stress Test
To collect the EST data, a treadmill ECG monitoring system was used. A schematic image is shown in Figure 3. This system consists of two components: an ECG monitoring device (BIOPAC MP36, BIOPAC Systems, Inc., Goleta, CA, USA) and a computer (Think-Pad X1 Carbon 7th, Lenovo, Hong Kong, S.A.R. of China). The exercise-induced stressor

Exercise Stress Test
To collect the EST data, a treadmill ECG monitoring system was used. A schematic image is shown in Figure 3. This system consists of two components: an ECG monitoring device (BIOPAC MP36, BIOPAC Systems, Inc., Goleta, CA, USA) and a computer (ThinkPad X1 Carbon 7th, Lenovo, Hong Kong, S.A.R. of China). The exercise-induced stressor was an exercise on a treadmill (MATRIX T5x, Johnson Health Tech Japan Co. Ltd., Minato-ku, Tokyo, Japan). Three adhesion-type electrodes (BlueSensor SP-00-S, Ambu A/S, Ballerup, Denmark) were attached to the subject's right chest, left subcostal, and right subcostal during data collection, as shown in Figure 3. The electrodes were connected to the ECG monitoring device. The ECG monitoring device sampled the ECG signal at 1 kHz on one channel during the EST. was an exercise on a treadmill (MATRIX T5x, Johnson Health Tech Japan Co. Ltd., Minatoku, Tokyo, Japan). Three adhesion-type electrodes (BlueSensor SP-00-S, Ambu A/S, Ballerup, Denmark) were attached to the subject's right chest, left subcostal, and right subcostal during data collection, as shown in Figure 3. The electrodes were connected to the ECG monitoring device. The ECG monitoring device sampled the ECG signal at 1 kHz on one channel during the EST. Data collection of the EST included the resting stage (R1) before the EST, the EST, the recovery stage (R2), and the relax stage (R3) after the EST. Data collection started at R1 and completed at R3. The EST was terminated before completion at the request of the subject or in case of poor signs, e.g., fatigue, shortness of breath, wheezing, claudication, leg cramps, or chest pain, being observed.
The modified Bruce protocol was used in this EST. As one of the most popular protocols utilized in exercise laboratories, the Bruce protocol has also been considered before data collection of the EST. Due to the heavy workload of the Bruce protocol for subjects, as illustrated in Table 1, no subject completed it in its entirety during rehearsal. Therefore, the modified Bruce protocol that was utilized had a lower workload than the standard test, as illustrated in Table 1 [9]. In R1, R2, R3, and the EST, the duration of each stage was 3 min, and the total duration of data collection was 30 min. Each subject completed data collection of the EST twice. The lead II ECG signals were collected for processing and analysis in this study. Data collection of the EST included the resting stage (R1) before the EST, the EST, the recovery stage (R2), and the relax stage (R3) after the EST. Data collection started at R1 and completed at R3. The EST was terminated before completion at the request of the subject or in case of poor signs, e.g., fatigue, shortness of breath, wheezing, claudication, leg cramps, or chest pain, being observed.
The modified Bruce protocol was used in this EST. As one of the most popular protocols utilized in exercise laboratories, the Bruce protocol has also been considered before data collection of the EST. Due to the heavy workload of the Bruce protocol for subjects, as illustrated in Table 1, no subject completed it in its entirety during rehearsal. Therefore, the modified Bruce protocol that was utilized had a lower workload than the standard test, as illustrated in Table 1 [9]. In R1, R2, R3, and the EST, the duration of each stage was 3 min, and the total duration of data collection was 30 min. Each subject completed data collection of the EST twice. The lead II ECG signals were collected for processing and analysis in this study.

Signal Processing
All the signal processing procedures were implemented using MATLAB (R2020b, The MathWorks, Inc., Natick, MA, USA). The ECG signals were collected using different devices at different sampling rates. The ECG signals of the EST were sampled using BIOPAC MP36 at 1000 Hz. The ECG signals of the BT were sampled using OpenBCI Ganglion at 200 Hz and were upsampled from 200 Hz to 1 kHz for coordination with the former.
First, the ECG signals of the BT and EST were decomposed into seven levels by multilevel one-dimensional wavelet decomposition using one of Daubechies wavelets, "db8", and the final approximation coefficient was taken as the baseline drift and subtracted from the original signal. Second, a Butterworth bandpass filter was implemented to remove power-line noise and high-frequency distortions. Then, R wave peaks detection was implemented after noise removal. The R-R interval (RRI) is defined as the sequence of the time intervals occurring between each pair of consecutive R wave peaks. The RRIs were calculated based on the ECG signal R wave peaks detected. Finally, the normal-to-normal RRI (NNI) was calculated based on the RRI where unreliable RRI were excluded. The NNI dataset of the EST was calculated based on 3 min ECG signals of each stage. The NNI dataset of the BT was calculated based on each 30 s segmentation ECG signal. Figure 4 shows an example of fully BT and EST ECG signal preprocessing in 10 s, including signal detrending, signal denoising, and R wave peaks detection. In Figure 4, (a) shows the signal processing of a BT ECG and (b) shows the signal processing of an EST ECG signal. The 10 s BT and EST ECG signals are extracted from the ECG signals of subject V's BT and EST. The collection condition of BT ECG signals is at the water temperature of 39 • C. The EST ECG signals are 10 s of stage 5 of test I.

HRV Features' Extraction
In this study, the obtained NNI data after signal preprocessing including denoising and R wave peak detection were used to extract HRV features. To reflect the heart functionalities comprehensively during the EST and BT, HRV features were extracted in the time domain, frequency domain, and nonlinear domain. All extracted HRV features are common features in use. In the nonlinear domain, features alpha 1 and alpha 2 of detrended fluctuation analysis (DFA) describe short-term and long-term fluctuation, respectively. Therefore, DFA alpha 1 was extracted. The extracted HRV features are summarized in Table 2 [30][31][32].

HRV Features' Extraction
In this study, the obtained NNI data after signal preprocessing including denoising and R wave peak detection were used to extract HRV features. To reflect the heart functionalities comprehensively during the EST and BT, HRV features were extracted in the time domain, frequency domain, and nonlinear domain. All extracted HRV features are common features in use. In the nonlinear domain, features alpha 1 and alpha 2 of detrended fluctuation analysis (DFA) describe short-term and long-term fluctuation, respectively. Therefore, DFA alpha 1 was extracted. The extracted HRV features are summarized in Table 2 [30][31][32].

Principal Component Analysis
Twenty-three HRV features in the time domain, frequency domain, and nonlinear domain were extracted, some of which have strong correlations. The high-dimensionality limits exploration of the data, therefore, it is necessary to extract the important information from the 23 HRV features. PCA is a multivariate technique that analyzes the data and describes several intercorrelated quantitative dependent variables. It is used to extract the important information from the data and represent the data as new variables called principal components. It has become one of the most useful tools for data modeling, compression, and visualization [33][34][35]. PCA was implemented for each of the 23 HRV features data of the BT, and the sums of percentages of the first two principal components (PCs) are above 95%, as shown in Figure 5.

Principal Component Analysis
Twenty-three HRV features in the time domain, frequency domain, and nonlinear domain were extracted, some of which have strong correlations. The high-dimensionality limits exploration of the data, therefore, it is necessary to extract the important information from the 23 HRV features. PCA is a multivariate technique that analyzes the data and describes several intercorrelated quantitative dependent variables. It is used to extract the important information from the data and represent the data as new variables called principal components. It has become one of the most useful tools for data modeling, compression, and visualization [33][34][35]. PCA was implemented for each of the 23 HRV features data of the BT, and the sums of percentages of the first two principal components (PCs) are above 95%, as shown in Figure 5. Therefore, the first two PCs as the two main PCA components were used to approximate the 23 HRV features data. The 23 HRV features data of the BT were implemented in dimensionality reduction to two-dimensional (2-D) HRV data of the BT. The new 2-D HRV data of the EST were calculated using the PC coefficients of the BT HRV features data. Therefore, the first two PCs as the two main PCA components were used to approximate the 23 HRV features data. The 23 HRV features data of the BT were implemented in dimensionality reduction to two-dimensional (2-D) HRV data of the BT. The new 2-D HRV data of the EST were calculated using the PC coefficients of the BT HRV features data.

Correlation Analysis of the Bathing Test and Exercise Stress Test
To explore the potential correlation between the BT and EST, the Voronoi diagram was used to cluster the BT data based on the Voronoi diagram generated from the EST data. In the 2-D HRV data of the BT and EST after PCA, two main PCA components (PC1 and PC2) were selected. PC1 and PC2 of the EST were used to generate a 2-D Voronoi diagram of EST Stages 1-7. For the 2-D HRV data of each subject's two ESTs, the mean of PC1 and mean of PC2 of each stage were calculated to be the site point in the Voronoi diagram of each stage in the Voronoi diagram. PC1 is the abscissa value and PC2 is the ordinate value. The site point of each stage was used to generate a Voronoi region on the Voronoi diagram. As shown in Figure 6, St. 1-7 are the site points of Stages 1-7, and the seven Voronoi regions correspond to the seven stages of the EST. The 2-D HRV data after PCA of 30 s duration before and after per minute during the BT are calculated as the mean of the data for the corresponding duration. In each duration of the BT, the mean of PC1 and mean of PC2 are the feature point's abscissa value and ordinate value, respectively. To evaluate the equivalent EST stage of each duration of the BT, the EST site point closest to the BT feature points was queried. The HRV behaviors of the EST stage correspond to the closest site point, and this BT feature point has the highest similarities. The EST stage corresponding to the closest site point serves as the equivalent EST stage of this BT feature point. Euclidean distance is used to query the EST site point closest to the BT feature points. The Euclidean distance can be calculated using the following equation: seven Voronoi regions correspond to the seven stages of the EST. The 2-D HRV data after PCA of 30 s duration before and after per minute during the BT are calculated as the mean of the data for the corresponding duration. In each duration of the BT, the mean of PC1 and mean of PC2 are the feature point's abscissa value and ordinate value, respectively. To evaluate the equivalent EST stage of each duration of the BT, the EST site point closest to the BT feature points was queried. The HRV behaviors of the EST stage correspond to the closest site point, and this BT feature point has the highest similarities. The EST stage corresponding to the closest site point serves as the equivalent EST stage of this BT feature point. Euclidean distance is used to query the EST site point closest to the BT feature points. The Euclidean distance can be calculated using the following equation: where is the abscissa value of the BT feature point, is the ordinate value of the BT feature point, is the abscissa value of the EST site point, is the ordinate value of the EST site point, and D is the Euclidean distance between the BT feature point and EST site point. The EST site point with the smallest distance from the BT feature point is the site point with the most similar features. The EST stage represented by this site point is the stage to which the highest probability of the BT at the present bathing conditions represented by this feature point is equivalent. The equivalent EST stage in terms of HRV behaviors of the BT feature point was determined based on the corresponding EST stage of the closest EST site point.  As shown in Figure 6, the BT feature points of 37 • C water temperature are mostly located in the Voronoi region of Stage 2. The BT feature points of 38 • C water temperature are mostly located in the Voronoi region of Stages 3 and 4. The BT feature points of 39, 40, and 41 • C water temperature are mostly located in the Voronoi region of Stages 4 and 5.
To explain the similarity of HRV behaviors between the BT and EST more clearly, the original seven EST stages are grouped into three Ex. stages, namely, low stage, medium stage, and high stage. Stages 1 and 2 in the EST are merged as the low stage. The medium stage consists of Stages 3 and 4 in the EST. The high stage contains Stages 5, 6, and 7 in the EST. The mean of Stages 1 and 2 of two ESTs is used as the site point for the low stage. The site point of the medium stage is generated by the mean of Stages 3 and 4. The mean of Stages 5, 6, and 7 of two ESTs is used as the site point for the high stage. A 2-D Voronoi diagram of the three regions for the low, medium, and high stages was generated using the same method. The regions where the BT feature points are located are confirmed in the Voronoi diagram with the same method. As shown in Figure 7, for BTs at the five water temperatures of one subject, the BT feature points are located on a 2-D Voronoi diagram. The Voronoi diagram is generated and three Voronoi regions are generated by the site points of low, medium, and high stages.
of Stages 5, 6, and 7 of two ESTs is used as the site point for the high stage. A 2-D Voronoi diagram of the three regions for the low, medium, and high stages was generated using the same method. The regions where the BT feature points are located are confirmed in the Voronoi diagram with the same method. As shown in Figure 7, for BTs at the five water temperatures of one subject, the BT feature points are located on a 2-D Voronoi diagram. The Voronoi diagram is generated and three Voronoi regions are generated by the site points of low, medium, and high stages. For all subjects, the equivalent EST stages of BTs at the five water temperatures were calculated. The means of the corresponding stage of all subjects under the same water temperature and bathing duration during the BT were calculated. With water temperature as the x-axis, bathing duration as the y-axis, and corresponding stage as the z-axis, the fitting surface was implemented using cubic spline interpolation and polynomials.

Biometric Information of Subjects and Completion Status of Data Collection
All the biometric information of subjects was obtained, and the completion status of the EST and BT is illustrated in Table 3. Fifty BTs were implemented completely, and the end-stages in 7 of 20 ESTs are Stage 7. Three subjects reached Stage 7 at the end-stage in tests I and II. According to the experiment diary, subjects who completed Stage 7 exercise frequently. To enhance the repeatability and improve the reliability, two ESTs were implemented. As illustrated in Table 3, the number values of tests I and II indicate the endstage in which a subject reached. Decimals indicate that the test of the stage was started For all subjects, the equivalent EST stages of BTs at the five water temperatures were calculated. The means of the corresponding stage of all subjects under the same water temperature and bathing duration during the BT were calculated. With water temperature as the x-axis, bathing duration as the y-axis, and corresponding stage as the z-axis, the fitting surface was implemented using cubic spline interpolation and polynomials.

Biometric Information of Subjects and Completion Status of Data Collection
All the biometric information of subjects was obtained, and the completion status of the EST and BT is illustrated in Table 3. Fifty BTs were implemented completely, and the end-stages in 7 of 20 ESTs are Stage 7. Three subjects reached Stage 7 at the end-stage in tests I and II. According to the experiment diary, subjects who completed Stage 7 exercise frequently. To enhance the repeatability and improve the reliability, two ESTs were implemented. As illustrated in Table 3, the number values of tests I and II indicate the end-stage in which a subject reached. Decimals indicate that the test of the stage was started but not fully implemented. Further, 3, 5, and 7 after the decimal point represent 1 min, 1.5 min, and 2 min, respectively. For example, 5.3 indicates Stage 5 was completed and Stage 6 was completed for 1 min; 6.5 represents Stage 6 was implemented and Stage 7 was implemented for 1.5 min; 6.7 indicates Stage 6 was finished and Stage 7 was stopped at 2 min. Note: all subjects were healthy males. * Decimals indicate the test of the stage was started but not fully implemented. 6.7 indicates Stage 6 was completed and Stage 7 was implemented for 2 min. 3 and 5 after the decimal point represent 1 min and 1.5 min, respectively.

Evaluation of the Relationship between Bathing Conditions and EST Stages
EST stages are grouped into three stages (low, medium, and high) with a Voronoi diagram using the 2-D HRV data of the EST after PCA analysis of each subject. The EST stages with the highest probability corresponding to the BT conditions at different water temperatures and bathing durations were calculated. In all equivalent EST stages of the BTs, the low stage accounts for 17.86%, medium stage accounts for 52.86%, and high stage accounts for 29.29%. The means and standard deviations of water temperatures of equivalent low, medium, and high stages are 38.21 ± 1.28, 38.93 ± 1.38, and 39.61 ± 1.28 • C, respectively. Of water temperatures corresponding to the low stage, 37 • C accounts for 40.80% and 38 • C accounts for 24.00%. Of water temperatures corresponding to the medium stage, 37, 38, 39, 40, and 41 • C account for 19.19%, 23.51%, 20.54%, 18.92%, and 17.84%, respectively. In addition, 39, 40, and 41 • C account for 21.84%, 25.73%, and 32.04%, respectively, of the water temperatures corresponding to the high stage in terms of HRV behaviors. The details of investigating the relationship between water temperature and the equivalent EST stage (Bruce and grouped) in terms of HRV behaviors are illustrated in Tables 4 and 5, respectively. As shown in Table 4, in the BT at 37 • C water temperature, the equivalent EST Bruce stages are mainly concentrated in Stages 2 (23.57%) and 3 (38.57%).
In the BT at 41 • C, the equivalent EST Bruce stages are mainly distributed in Stages 4 (37.14%) and 5 (26.43%). As the water temperature increases from 37 • C to 41 • C, the main distribution of equivalent EST Bruce stages is from Stages 2 and 3 up to Stages 4 and 5. The relationship between water temperature and the equivalent EST grouped stages is more obviously illustrated in Table 5. As the water temperature increased from 37 • C to 41 • C, the proportion of the low stage of the equivalent EST grouped stage decreased from 36.43% to 5.71%, and the proportion of the high stage increased from 12.86% to 47.14%. The means of the corresponding stage of all subjects under the same water temperature and bathing duration during the BT were calculated. The 3-D relationship diagram between water temperature, bathing duration, and EST stage is shown in Figure 8. Cubic splines were used to interpolate in Figure 8. The means of the corresponding stage of all subjects under the same water temperature and bathing duration during the BT were calculated. The 3-D relationship diagram between water temperature, bathing duration, and EST stage is shown in Figure 8. Cubic splines were used to interpolate in Figure 8. To show the relationship more intuitively between water temperature, bathing duration, and equivalent EST stage, surface fitting was implemented using polynomials. Figure 9 shows the relationship surfaces fitted by the first order (a), second order (b), third order (c), and fourth order (d) polynomial surfaces. At the water temperature of 41 °C, the To show the relationship more intuitively between water temperature, bathing duration, and equivalent EST stage, surface fitting was implemented using polynomials. Figure 9 shows the relationship surfaces fitted by the first order (a), second order (b), third order (c), and fourth order (d) polynomial surfaces. At the water temperature of 41 • C, the equivalent EST stage after 14 min of bathing is the highest. As the water temperature and bathing duration increase, so does the equivalent stage in terms of HRV behaviors.
The approximation equations of first to fourth order polynomials are as follows: where WT is the water temperature of the BT, BD is the bathing duration of the BT, and ST is the equivalent EST stage.
Electronics 2022, 11, x FOR PEER REVIEW 13 of 17 equivalent EST stage after 14 min of bathing is the highest. As the water temperature and bathing duration increase, so does the equivalent stage in terms of HRV behaviors.
where is the water temperature of the BT, is the bathing duration of the BT, and is the equivalent EST stage.
Equations (2) to (5) are approximation equations of first to fourth order polynomials of the relationship fitting surface (a) to (d) in Figure 9, respectively. Equation (2) is the bivariate first order polynomials and the relationship fitting surface is flat. Water temperature and bathing duration are proportional to the equivalent EST grouped stage. Equations (3) and (4)   where ST m andŜT m represent the original and fitted values of the equivalent EST grouped stage at the same water temperature and bathing duration, respectively, and M is the number of sample sets. The closer SSE is to 0, the better the fitting effect. The coefficient of determination (R 2 ) is calculated to evaluate the fitting effect of the fitting function. When the R 2 result is between 0 and 1 and the closer the result is to 1, the better the fitting effect is. The equation is shown in (7). As illustrated in Table 6, comparing the SSE and R 2 results of first to fourth order polynomials, the fourth order polynomial has the best fitness performance, i.e., SSE is closest to 0, and R 2 is closest to 1.

Discussion
In this study, we collected 10 subjects' ECG signals during the BT and EST and explored whether daily bathing can serve as a potential alternative means of EST based on HRV behaviors evaluated in three analysis domains. Stressful loads by different water temperatures and bathing durations of the BT correspond to different exercise stages in terms of HRV behaviors. As shown in Figure 9, in general, the higher the water temperatures during bathing and the longer the bathing durations, the higher the equivalent stages. More specifically, in the high stage (Stages 5-7), 32.04% equivalent stages in terms of HRV behaviors correspond to bathing at 41 • C water temperatures, and 25.73% correspond to bathing at 40 • C. In all equivalent EST stages of BTs, the medium stage accounts for 52.86%. Therefore, bathing can be an alternative to light EST in daily life.
In our previous study, we defined an SI to quantify the stress level [17]. The SI mainly explains psychological stress and is not applicable to this study. In this study, a new procedure is explored for an alternative to EST. HRV features were extracted to reflect the heart behavior during the EST and BT. Different levels of load are generated by bathing with different water temperatures and bathing durations. Different bathing loads correspond to different exercise stages in terms of HRV behaviors, as shown in Figure 9. In all equivalent EST stages of BTs, the medium stage accounts for 52.86%. The equivalent EST stages in terms of HRV behaviors of BTs are mostly distributed in the low and medium stages. It is hard to correspond Stage 7 to the BT with normal and nonextreme conditions. If a certain equivalent EST stage is necessary to implement during daily health monitoring, it is required to bathe at the corresponding water temperature for the corresponding bathing duration.
The EST is a routine procedure in the clinic to evaluate cardiovascular functions, but is commonly too demanding and sometimes risky for some subjects [7,8]. As illustrated in Table 3, even with the modified Bruce protocol during ESTs, we confirmed that most subjects did not complete Stage 7 in the end-stage. In addition, there are some contraindications to the EST, such as acute myocardial infarction within 2 days, ongoing unstable angina, active endocarditis, decompensated heart failure, acute myocarditis or pericarditis, etc. [1]. In Japan, bathing is popular with many people and has become a daily activity [10]. Bathing may not only be safe for patients with heart failure but may even briefly improve hemodynamics with a short bathing duration [36]. It was confirmed that warm bathing, sauna bathing, and spring bathing improve cardiac function and have beneficial effects with heart failure patients [37][38][39][40]. Lee et al. proved that forest bathing has positive effects on physical and mental health [41]. Bathing can not only correspond to a light stress stage but also benefits cardiovascular function and physical and mental health. The BT is safer, more comfortable, and more convenient.
In this study, the relationship between water temperature, bathing duration, and EST grouped stages was investigated. A detailed analysis of the relationship with Stages 1-7 of the modified Bruce protocol was not implemented yet. Due to the limited data volume, complex classification methods are not used. All subjects were male students. There is great individual variability among subjects. More data collection is in progress. Future work will involve more data from different age and gender groups, as well as grouping subjects for more precise results. More methods will be explored to reveal more outcomes.

Conclusions
This study revealed the relationship between BT and EST by investigating the HRV behaviors in two tests. As the water temperature and bathing duration increase during bathing, the equivalent EST stages also increase. The equivalent EST grouped stages are mainly concentrated in the medium stage (52.86%). We confirmed that an alternative procedure to light EST is able to be provided by daily bathing. In the BT, the water temperature and bathing duration can be chosen to implement the corresponding exercise stage. This can serve as a tool for daily monitoring of the health condition or as a reference for clinical application.  Institutional Review Board Statement: The ECG recording procedures were approved by the Public University Corporation, the University of Aizu Research Ethics Committee.
Informed Consent Statement: Written informed consent was obtained from each subject before data collection.
Data Availability Statement: Data are available from the biomedical information engineering laboratory (contact via e-mail: wenxi@u-aizu.ac.jp) for researchers who meet the criteria for access to confidential data.