- 1Laboratory of Animal Husbandry Resources, Division of Applied Biosciences, Graduate School of Agriculture, Kyoto University, Kyoto, Japan
- 2Division of Grassland Farming, Institute of Livestock and Grassland Science, National Agriculture and Food Research Organization, Tochigi, Japan
- 3Department of Animal and Grassland Sciences, Faculty of Agriculture, University of Miyazaki, Miyazaki, Japan
- 4Department of Bioinformatics, College of Life Sciences, Ritsumeikan University, Shiga, Japan
Heart rate variability (HRV) is the heart beat-to-beat variation under control of the cardiovascular function of animals. Under stressed conditions, cardiac activity is generally regulated with an upregulated sympathetic tone and withdrawal of vagal tone; thus, HRV monitoring can be a non-invasive technique to assess stress level in animals especially related to animal welfare. Among several stress-induced factors, heat stress is one of the most serious causes of physiological damage to animals. The aim of this study was to assess the effects of heat stress on HRV in small ruminants under free-moving conditions. In three experimental periods (June, August, and October), inter-beat intervals in sheep and goats (three for each) in two consecutive days were measured. HRV parameters were calculated from the inter-beat interval data by three types of analyses: time domain, frequency domain, and non-linear analyses. The temperature–humidity index (THI) was used as an indicator of heat stress, and vectorial dynamic body acceleration (VeDBA) was calculated to quantify the physical activity of the animals tested. First, we investigated correlations of THI and VeDBA with HRV parameters; subsequently, THI was divided into five categories according to the values obtained (≤ 65, 65–70, 70–75, 75–80, and >80), and the effects of the THI categories on HRV parameters were investigated with and without correcting for the effects of physical activity based on the VeDBA. The results indicated that HRV significantly decreased with increasing THI and VeDBA. For non-linear HRV parameters that were corrected for the effects of physical activity, it was suggested that there would be a threshold of THI around 80 that strongly affected HRV; high heat stress can affect the autonomic balance of animals non-linearly by inducing the sympathetic nervous system. In conclusion, to assess psychophysiological conditions of unrestrained animals by HRV analysis, the confounding effect of physical activity on HRV should be minimized for a more precise interpretation of the results.
Introduction
Heart rate variability (HRV) has been used as a sensitive indicator of the functional regulatory characteristics of the autonomic nervous system (1, 2). Inter-beat interval fluctuation expressed as HRV reflects the sympathetic and parasympathetic activity of the autonomic nervous system, and healthy cardiac function is characterized by irregular time intervals between consecutive heart beats (2). A decrease in HRV can be caused by an increase in sympathetic activity and/or withdrawal of parasympathetic activity, and HRV is a particularly good indicator for the non-invasive assessment of autonomic nervous system activity in response to various internal and external stressors.
Several studies have investigated HRV in animals: farm animals (2), companion animals (3, 4), monkeys (5), seabirds (6), etc. In particular, monitoring HRV in animals has recently gained attention as a non-invasive technique of assessing stress levels of animals related to animal welfare (7). For example, numerous studies have investigated internal and external psychological and physiological factors that affect the HRV of dairy cattle: temperament and reactivity to humans (8–10), bell noise (11), palpation (12, 13), insect harassment (14), seasonality (15), milking system (16, 17), pregnancy and calving (18–21), and diseases such as lameness (22), diarrhea (14), and bovine spongiform encephalopathy (23). However, physical activity (movement) influences HRV, and it can cloud the regulation linked to cognitive, emotional, social, and health process (24, 25). In fact, HRV is greatly affected by the behavior of animals, and some studies have focused on the effects of behavior on HRV, such as body posture (9), step behavior (17), and behavioral types (26). When the effects of psychophysiological changes on HRV of animals have to be evaluated, particularly under free-moving conditions, the quantified physical activity level of the target animals should also be considered as a key element that influences HRV (27).
Heat stress is a serious cause of physiological damage on animals and consequently affects their production (28–31). Many previous studies on heat stress for animals have focused on its negative effects on productivity (32–36) and physiological reactions related to hormonal and immune responses (37–40). For example, in dairy cattle, as milk production increases, metabolic heat production rises with the metabolism of large amounts of nutrients, which makes the high-producing cows more vulnerable to high ambient temperatures and humidity than animals that are less active metabolically (28). Also, heat stress compromises productivity in small ruminants, increasing maintenance energy requirement (41). For investigating such effects of heat stress on animals, thermoregulatory function traits such as respiration rate, rectal temperature, and heart rate were widely evaluated as physiological indicators of heat stress (34). However, although the effects of shading on HRV in cattle under extreme heat loads have recently been investigated (42, 43), no studies have analyzed the effects of heat stress on the autonomic nervous system of animals in detail using a variety of HRV analyses, in particular by considering the effect of physical activity. Regarding the variety of HRV analyses, not only conventional time domain and frequency domain indices but also non-linear indices have been used as reliable markers of sympathetic and vagal activation in HRV analysis (1). Quantifications in different domains of HRV analysis can focus on different characteristics of HRV: quantity of variance (time domain), periodic processes due to autonomic regulation (frequency domain), and chaotic phenomena in the regulation of cardiac activity (non-linear); and HRV parameters derived from such different domains can complement each other (16). Therefore, evaluation of different HRV parameters can be adequate to judge different quantitative and qualitative stress loads in animals (14). In order to successfully manage healthy animals, HRV assessment using a variety of parameters in multiple domains appears to be required, which provides quantitative and objective assessments of the effect of heat stress on the autonomic nervous system activity under free-moving conditions.
The aim of the present study was to assess the effects of heat stress on HRV in small ruminants using three types of analyses: time domain, frequency domain, and non-linear analyses. To elucidate the physiological interrelationships between heat stress and HRV under free-moving conditions, the effects of heat stress on HRV were evaluated with and without correcting for the confounding effect of physical activity on HRV by simultaneously measuring the dynamic body acceleration of the animals.
Materials and Methods
Animals and Experimental Periods
All of the data were obtained from experiments conducted at Kyoto University, Japan (35°02′N, 135°47′E). Three castrated Corriedale sheep that are 2.5 years old and three castrated Japanese Saanen goats that are 9 years old were used. The experiments were conducted in three periods: June (early summer), August (midsummer), and October (early autumn) in 2015. Data collection was carried out for all sheep and goats simultaneously, and the data were collected consecutively for about 2 days per animal in each period. The body weights (kg) of the animals were measured on the first day of each respective period: 34.3 ± 4.9, 34.9 ± 5.4, and 35.0 ± 2.5 kg for the sheep and 73.9 ± 3.9, 76.1 ± 2.2, and 72.8 ± 4.4 kg for the goats. The sheep and goats were separately managed into two pens (each species per pen: ~10 m2 per pen) in an open-sided, metal-roofed animal shelter. The sheep was shorn in the spring following conventional management. The animals were fed twice a day (9:30 and 15:30), and the amount of feed offered was 2% of their body weight per day (Italian ryegrass: alfalfa hay cube: concentrate = 10:3:3) as the amount required for maintenance. Mineral blocks and water were freely available. All of the animal experiments were approved by the Animal Experiment Committee of Kyoto University (permit number: 27-56). All of the procedures for equipping the animals with data loggers were performed as quickly as possible to minimize the animals' discomfort. All animals were housed for use in future researches.
Data Collection
Inter-beat interval data (ms) were obtained using heart rate monitors (RS800CX and H2 heart rate sensors, Polar Electro Ltd., Finland). Each heart rate monitor consisted of a transmitter with two electrodes and a logger. The electrodes were placed on the animal's right shoulder and left anterior thorax and attached with a homemade chest belt (44). The electrode sites of animals were roughly shaved and covered with a conductive gel to optimize electrode contact. In addition, three-dimensional accelerations were simultaneously recorded using acceleration data loggers (USB Accelerometer X6-1A and X6-2A, Gulf Coast Data Concepts, Waveland, MS, USA). The acceleration logger was attached to the back of the animals when the heart rate monitor was attached. Accelerations were recorded at 10 Hz, with a 16-bit resolution. The positions of the heart rate monitor and acceleration logger are shown in Figure 1.
Figure 1. Positions of heart rate (HR) transmitter, electrodes for the HR monitor, and the accelerometer.
During each experiment, ambient temperature and relative humidity per minute were measured using a thermo-hydro data logger (TR-72wf-H, T&D Corporation, Japan). The logger was placed in the animal shelter.
In addition, respiration rate and rectal temperature, traditional physiological indices, were measured on the last day of each experimental period. Respiration rate was measured visually by counting thoracoabdominal movements for 30 s twice by each of three observers on the last day of each experimental period at a distance of ~1 m between animals and observers. Then, the doubled value of the average of six counting values by the three observers was taken as the respiration rate per minute for each animal in each period. The observation was carried out for the target animals in standing position. Rectal temperature was measured once for each animal at the end of each experiment (before evening feeding) by a digital thermometer (MC-170, Omron Healthcare Co., Ltd., Japan).
Temperature–Humidity Index
The effect of heat is aggravated when heat stress is accompanied with high ambient humidity (34). Therefore, the THI was used as an indicator of heat stress in the present study, which was calculated from the average values of ambient temperature and relative humidity in 5 min using the following equation (45):
Effect of Physical Activity on HRV
A change in physical activity levels is thought to be one of the main modulators of HRV (46). Because cardiac activity is affected by variations in the respiration rate via sympathetic and parasympathetic nerves (47), HRV is greatly influenced by the change in respiration caused by physical activity. Therefore, to improve the precision of HRV evaluation under free-moving conditions, the effect of physical activity on HRV should be considered. One of the possible methods to consider the effect is to delete the sections of the heart inter-beat interval data with excessive movement before the HRV analysis (48–50), but there have been no general ways for separating the effect of changes in quantified physical activity levels on HRV from the influence of other regulatory processes (24). Therefore, in this study, the vectorial dynamic body acceleration (VeDBA), which can be a proxy for the activity-specific energy requirements of animals (25, 51, 52), was used to remove the confounding effect of physical activity on HRV based on the method by Oishi et al. (27). To calculate VeDBA, triaxial acceleration values obtained by the acceleration loggers were used in the following equation:
where Ax, Ay, and Az are triaxial dynamic body accelerations calculated by subtracting static accelerations from raw accelerations (53). Initially, the VeDBA values were transformed into natural logarithms for normalization, which were then averaged in 5 min to correspond with the respective HRV data and used for the analysis.
HRV Analysis
Kubios HRV 2.2 software (Kubios Ltd., Finland) (54) was used to calculate HRV parameters from inter-beat interval data. All of the HRV parameters for both animal species were calculated by 5-min windows according to the Task Force of the European Society of Cardiology, North American Society of Pacing and Electrophysiology (55). A total of 9,963 adjacent 5-min windows were obtained for the analysis. Artifacts were removed from the data using a threshold-based artifact correction algorithm built-in in this software with a medium correction level (54). Moreover, data that included 5-min-averaged inter-beat interval values outside the 3-sigma range of individual datasets were also removed. The proportion of data removed was 1.13%. A total of 9,850 5-min intervals, with an average of 547.2 intervals (45.6 h) per animal in each period, were included in the analysis. For each animal, the intervals were almost evenly distributed across each period.
In the present study, the mean heart rate per 5 min (HR) and the following six HRV parameters were calculated using three domain analyses: time domain parameters (SDNN and RMSSD) from time domain analysis, frequency domain parameters (HF and LF/HF) from frequency domain analysis, and non-linear domain parameters (Lmax and %DET) from non-linear domain analysis. These parameters are defined in Table 1.
Time domain analysis is the simplest form of HRV analysis. Time domain parameters reflect various aspects of the statistical variability of inter-beat intervals and are often used to interpret HRV characteristics. The standard deviation of inter-beat intervals (SDNN) is a good predictor of overall variability influenced by both sympathetic and parasympathetic activity, and the root mean square of successive inter-beat interval differences (RMSSD) is the primary time domain measure used to estimate the high-frequency beat-to-beat variation that represents vagal tone activity (2).
Frequency domain analysis is the procedure of decomposing a waveform of inter-beat intervals by range of the frequency using fast Fourier transformation (FFT). In this analysis, the waveform of inter-beat intervals is regarded as a synthesized waveform, and the power spectra of the waveform are calculated to estimate the HRV. An increase in high-frequency (HF) components is generally caused by increasing HRV, and HF has been used to describe the function of vagal tone (56). In addition, the ratio of the power spectra of low-frequency (LF) components to HF (LF/HF) is utilized for describing the sympathovagal balance of the autonomic nervous system (57). The ranges of the HF and LF components for the tested animals were set to be 0.20–0.40 and 0.04–0.20 Hz, respectively, as recommended by von Borell et al. (2).
Non-linear domain analysis elucidates the chaotic behavior of HRV using non-linearity indicators. A great deal of information can be extracted from physiological signals by describing their dynamic behavior, and non-linearity is the representative indicator of such complex dynamical systems (58). Recurrence quantification analysis (RQA) is a method of non-linear HRV analysis that was developed by Eckmann et al. (59), which has been used to detect hidden and complex characteristics of HRV (60). In the present study, the length of the longest line of recurrence points (Lmax) and the percentage of determinism (%DET) were used as the quantitative parameters of RQA. These parameters reflect the richness of the deterministic structure of inter-beat intervals in a time series, and an increase in these parameters indicates a decrease in HRV (61). The evaluation of non-linear parameters with time domain and frequency domain parameters is a useful method of analyzing HRV for the accurate interpretation of the autonomic nervous system function (14, 16).
Statistical Analysis
Ambient temperature, relative humidity, and THI were analyzed by one-way analysis of variance (ANOVA), with period (June, August, and October) included as a fixed effect. Respiration rate and rectal temperature were analyzed by two-way ANOVA, with period (June, August, and October) and species (sheep and goats) included as fixed effects.
Regarding HRV, Pearson's correlation coefficients of mean HR and HRV parameters with THI and VeDBA were firstly calculated. In addition, the THI was divided into five categories (≤ 65, 65–70, 70–75, 75–80, and >80), and the effects of the THI category on mean HR and HRV parameters were analyzed using the following linear mixed model:
where Yijkl is the value of the mean HR and HRV parameters (SDNN, RMSSD, HF, LF/HF, Lmax, and %DET) per 5 min, μ is the overall mean, catTHIi (i = 1–5) is the fixed effect of the THI category, Speciesj (j = 1 or 2) is the fixed effect of species, Animalk(j) (k = 1–6) is the random effect of individual animals nested within species, and eijkl is the error. Furthermore, based on the method by Oishi et al. (27), the above statistical model was transformed into the following model in order to include the effect of quantified physical activity on mean HR and HRV parameters:
where β (VeDBA)ijkl is the covariate effect of VeDBA per 5 min.
Differences were analyzed using the least squares means with the Tukey–Kramer post-hoc test (62) and were considered significant at P < 0.05. Correlation coefficients were calculated using PROC CORR, and the other analyses were performed using PROC MIXED in SAS 9.3 (SAS Institute) (63).
Results
Respiration Rate and Rectal Temperature and Changes in Environmental Conditions
The least square means of ambient temperature, relative humidity, and THI during each experimental period are shown in Table 2. The ambient temperature and THI in August were the highest, and those in October were the lowest (P < 0.05). The relative humidity in June was higher than that in August or October (P < 0.05).
The least square means of respiration rate and rectal temperature are shown in Table 3. The respiration rate in August was the highest (P < 0.05), and the respiration rate in sheep tended to be higher than that in goats (P = 0.053). The rectal temperature of sheep was higher than that of goats (P < 0.05), but no significant effect of period was found.
Effects of Heat Stress on HRV
Pearson's correlation coefficients between HRV parameters and THI and between HRV parameters and VeDBA are shown in Table 4. Both the THI and VeDBA were significantly correlated with the mean HR and all HRV parameters (P < 0.05). From the results of the coefficients, both the THI and VeDBA were negatively correlated with HRV, i.e., positively correlated with LF/HF and the non-linear parameters and negatively correlated with the other HRV parameters. For most of the HRV parameters, HRV was more strongly correlated with VeDBA than with THI.
Results of the changes in the mean HR and HRV parameters in the THI categories by analyzing the models with and without the effect of VeDBA are illustrated in Figure 2 (mean HR and time domain parameters) and Figure 3 (frequency domain and non-linear parameters). Regardless of whether the effect of physical activity was included in the model, the fixed effect of species was not significant for the mean HR and HRV parameters, except for the two frequency domain parameters; HF in goats was significantly higher than that in sheep, and LF/HF in goats was significantly lower than that in sheep (P < 0.05). The fixed effect of the THI category was significant for the mean HR and all HRV parameters (P < 0.05). When the effect of physical activity was included, the covariate of VeDBA was significant for the mean HR and all HRV parameters (P < 0.05).
Figure 2. Changes in mean HR and time domain HRV parameters (SDNN and RMSSD) by THI category, (A) without and (B) with the effect of physical activity quantified by VeDBA. The bars represent least square means, and error bars show standard errors. Values with different letters differ significantly (P < 0.05). HR, heart rate (bpm); HRV, heart rate variability; SDNN, standard deviation of inter-beat intervals (ms); RMSSD, square root of the mean squared differences of successive inter-beat intervals (ms); THI, temperature–humidity index; VeDBA, natural logarithmically transformed vectorial dynamic body acceleration (g).
Figure 3. Changes in frequency and non-linear domain parameters by THI category, (A) without and (B) with the effect of physical activity quantified by VeDBA. HF and LF/HF are frequency domain parameters (gray bars), and Lmax and %DET are non-linear parameters (white bars). The bars represent least square means, and error bars show standard errors. Values with different letters differ significantly (P < 0.05). HRV, heart rate variability; HF, normalized power of the high-frequency band (n.u.); LF/HF, ratio of the normalized power of the low-frequency (LF) band to HF; Lmax, the length of the longest line of recurrent points (beats); %DET, percentage of recurrent points that appear in sequence (%); THI, temperature–humidity index; VeDBA, natural logarithmically transformed vectorial dynamic body acceleration (g).
Regarding the mean HR and time domain parameters, when the effect of physical activity was not included, the mean HR gradually increased when the THI was ≤ 65 to >80, the SDNN decreased when the THI was ≤ 65 to 70–75 (but did not significantly change from 70–75 to >80), and the RMSSD decreased step by step from THI ≤ 65 to >80 (Figure 2A). Similar results were obtained when the effect of physical activity was included, although differences between 70–75 and 75–80 in mean HR, SDNN, and RMSSD were not significant (Figure 2B). Regarding the frequency domain parameters, when the effect of physical activity was not included, HF decreased and LF/HF increased from 65–70 to >80 (Figure 3A). However, when including VeDBA, HF and LF/HF had maximum and minimum values, respectively, in 65–70 (Figure 3B). With regard to the non-linear HRV parameters, the two non-linear parameters increased from THI ≤ 65 to 65–70, slightly changed from 65–70 to 75–80, and increased again from 75–80 to >80 (Figure 3A). Surprisingly, when the effect of physical activity was included, the non-linear HRV parameters only increased from 75–80 to >80 (Figure 3B), indicating that there was a threshold level between the physiological states of the two THI categories.
Discussion
Monitoring HRV has been used as a non-invasive method of investigating characteristics of the autonomic nervous system. In the present study, HRV was characterized using three domain parameters in order to reveal the general features of HRV in ruminants under heat-stress conditions.
Correlations of Heat Stress and Physical Activity With HRV
The two time domain HRV parameters (SDNN and RMSSD) were significantly, negatively correlated with the THI (Table 4). As for the frequency domain parameters, HF was negatively correlated with the THI, whereas LF/HF was positively correlated with it. Therefore, HRV parameters which represent vagal tone function in the two domains (RMSSD and HF) decreased under high-THI conditions, which might be in agreement with the previous study (43) showing that RMSSD of calves decreased under highly heat-stressed daytime. As for the non-linear parameters, Lmax and %DET were significantly, positively correlated with the THI, which reflected a change toward a decrease in HRV by increasing sympathetic nervous system activity and a more periodic heart rate under stress. Hence, all of the HRV parameters showed that HRV decreased with increasing THI. However, it is noticed that the correlation between THI and HRV was significant but considerably weak, which might be due to the existence of thermoneutral range for tested animals during the experimental periods. In addition, despite the fact that the animals were housed in an animal shelter and their physical activity levels were low, VeDBA was significantly correlated with all HRV parameters (Table 4), which was in accordance with our previous study (27). Furthermore, the absolute values of the correlation coefficients between VeDBA and the HRV parameters were mostly higher than those between the THI and the HRV parameters. These results suggest that animals' physical activity influences HRV and that the relationship between physical activity and HRV may be stronger than that between the THI and HRV. Therefore, we can conclude that correcting for the confounding effect of physical activity is necessary to evaluate the effect of heat stress on HRV more precisely.
Effects of Heat Stress and Physical Activity on the Three Domains of HRV Parameters
It was supposed that heat stress can interfere non-linearly with the physiological function of animals, although a weak linear correlation between the THI and HRV was found. In the present study, therefore, the THI was divided into five categories (levels) in the statistical model, and the effect of the THI categories on HRV was analyzed with and without correcting for the effect of physical activity (Figures 2, 3). First, regardless of whether VeDBA was included as a covariate in the statistical models, the effect of species was significant only for HF and LF/HF. Frequency domain parameters can be modified when the autonomic nervous system responds to changes in the respiration rate (64). In the present study, the respiration rate of sheep tended to be higher than that of goats (Table 3). Machando et al. (65) suggested a physiological susceptibility to heat in sheep with higher respiration rates when compared with goats. However, Johnson (66) reported that sheep and goats showed similar changes in respiration rate when the animals were shorn to the same hair length before the experiment. Hence, the result of differences in respiration rate in the present study might be due the effect of the regrown hair of sheep during the three experimental periods. Besides, the younger age and smaller body weight of sheep compared with goats also might be causes of this result. Thus, with the inclusion of such several differences between the two small ruminant species, the difference in respiration rate between the species was expressed as the effect of species on the frequency domain parameters of HRV. As for the effect of THI categories on the HRV parameters, when the effect of physical activity on HRV was not included, most of the HRV parameters showed that HRV decreased with increasing THI. However, when the effect of physical activity on HRV was included, we found specific changes with changes in THI category, which highlighted the characteristics of the three domains of analysis.
The HRV parameters in the time domain analysis (SDNN and RMSSD) decreased with increases in the THI categories in both analyses, with and without the inclusion of VeDBA. This finding is in accordance with those from our previous study, indicating that these time domain parameters do not strongly correspond with short-term changes in physical activity (27). In contrast, changes in the frequency domain parameters (HF and LF/HF) differed by the inclusion of the effect of physical activity; they formed curved patterns. It is possible that changes in respiration rate caused by an increase in THI affected the frequency domain parameters, even after correcting for the effect of physical activity. As already indicated, frequency domain HRV parameters can be strongly affected by respiration rate, and it is therefore crucial to control breathing in order to accurately interpret HRV when frequency domain parameters are used (67, 68). In fact, the respiratory frequency measured in the present study ranged widely between 0.24 and 1.41 Hz, while the setting HF range for the frequency domain analysis was limited to 0.20–0.40 Hz as generally suggested for sheep and goats (2). For the HF to have a bell-shaped distribution in the frequency domain analysis, any deviation of the respiratory frequency from the HF range simply masked the effect of respiration on the frequency domain HRV parameters, particularly at lower and higher HR conditions that corresponded to lower and higher THI categories. It is though that HF is one of the major parameters of cardiac vagal activity which is strongly linked with a range of self-control such as cognitive performance, emotion and stress regulation, health, and social interactions (69, 70). However, since controlling the breathing of animals is almost impossible, the frequency domain parameters are of limited use for evaluating autonomic regulation in freely-moving animals.
In comparison with the other two domains of HRV analysis, the non-linear domain parameters (Lmax and %DET) showed characteristic changes in the present study; they increased with the THI categories without correcting for the effect of physical activity, but only increased from 75–80 to 80 < when including the effect of physical activity. This result indicated that there was a THI threshold level of around 80 affecting the HRV parameters. Although the insusceptibility of non-linear HRV parameters to respiration is still in debate, it has been suggested that the non-linear complexity and determinism of HRV do not arise as a consequence of a respiration input into the cardiovascular oscillator (71, 72). The results of the present study suggest that, as long as animals breathe merely for gas exchange required for the metabolism, respiration does not induce non-linear HRV. However, if the breathing pattern is changed for some reason other than basal metabolism, the respiratory pattern may interfere with the central cardiovascular oscillator. Since one possible predictive marker of alterations in central autonomic regulation that may precede metabolic stress is a non-linear domain component of HRV, as suggested by Hoffmann et al. (73), we can conclude that it is beneficial to use non-linear HRV parameters as physiological heat-stress indicators.
It should be noted that the present study used only six individuals from two small ruminant species, which was one of the limitations of the present study. However, the fixed effects of THI category and species were properly analyzed with considering the effect of individuals in each species, since the effect of individual differences was treated as a random effect nested within the species in the statistical model with enough data for each animal. The results showed that the effect of species was not significant for HRV parameters except for frequency domain parameters, which indicated that the finding in particular for non-linear HRV parameters could be correctly evaluated. However, in order to strengthen the validity of the analysis and the conclusions for the effect of THI on HRV parameters and also to clarify the difference in the effect between the two small ruminant species in more detail, further studies using experimental designs with more animals would be required.
Heat Dissipation and Heart Rate Regulation
Many of the physiological responses under heat load are evoked to maintain the core temperature constant. For this purpose, heat dissipation from animal body is promoted by vasodilation, panting, and sweating. In the present study, moderate tachycardia was observed as the THI increased, which might be in accordance with the interpretation that increases in the heart rate under heat load is evoked due to increases in blood circulation to transfer heat from the core to the periphery (34).
The initial response that maintains the core temperature under heat stress is vasodilation (74), which results from the withdrawal of sympathetic nervous activity that governs the vessel tone (sensible heat dissipation). As the heat load increases, latent heat dissipation is gradually induced. Cholinergic sympathetic activity that regulates activity of sweat glands also acts as an active local vasodilator (75) that synergistically dissipates heat. Both responses to external heat stress, vasodilation and sweating, are regulated by the thermoregulatory center in the preoptic area. However, in the case of panting, the modified breathing pattern interferes with systemic circulation, and directly changes the heart rate via cardio-pulmonary baroreceptors. This mechanism may make the regulation of heart rate more complex and chaotic.
Ambient temperature and humidity are critical factors that affect the efficacy of both sensible and latent heat dissipation, and the THI is a useful index that indicates the heat-stress level in homeothermic animals. If body temperature cannot be maintained by the responses discussed so far because of rigid heat load, heat production is suppressed by reducing feed intake, which results in reduced production (76). However, in the present study, the rectal temperatures of the animals tested did not significantly differ among the three experimental periods (Table 3). Moreover, a decrease in feed intake was not observed for the tested animals. In dairy cows, Kadzere et al. (32) reported that THI values > 78 might cause extreme distress, with lactating cows being unable to maintain their thermoregulatory mechanisms or normal body temperatures. Lemerle and Goddard (77) reported that homeostatic mechanisms could prevent a rise in rectal temperature until the THI reaches 80, which might also be applied for the small ruminants in the present study. Therefore, although some degree of disorder of thermoregulatory control shown as the change in non-linear HRV parameters might occur over 80 of THI, the tested animals could mostly cope with heat stress by promoting heat dissipation until the THI reached the threshold value.
Conclusions
The present study revealed the effects of heat stress on HRV of sheep and goats. Under high THI conditions, HRV of the animals was decreased which might be due to the increase of sympathetic nervous system activity on heart rate regulation. From the evaluation of non-linear HRV parameters with correction for the effect of physical activity, this study could suggest the existence of a threshold value of THI around 80 for HRV. The threshold value might indicate a limit that the external stress imposes physiological non-linear heart rate regulation for the heat dissipation in order to maintain core temperature constant.
In recent years, increased interest has been paid to HRV measurement as a non-invasive technique to assess stress in animals in particular related to animal welfare. The results of the present study indicated both heat stress and physical activity levels affected HRV. Therefore, in order to investigate the effects of psychophysiological stress factors on HRV of unrestrained animals, the confounding effect of physical activity on HRV should be taken into consideration. Moreover, the use of multiple domains of HRV parameters, in particular non-linear parameters, should be recommended for investigating different characteristics of the effect of stressors on HRV.
Data Availability Statement
The datasets generated for this study are available on request to the corresponding author.
Ethics Statement
All of the animal experiments were approved by the Animal Experiment Committee of Kyoto University (Permit number: 27-56).
Author Contributions
KK and KO designed the study, collected and analyzed the data, and wrote the manuscript. MM and HA collected the data and helped write the manuscript. AS and YY collected the data. YH wrote the manuscript. HK and HH helped write the manuscript. All authors discussed the results, contributed to manuscript revisions and provided approval for publication.
Funding
This study was supported by Grant-in-Aid for Scientific Research (C) (Grant Numbers 15K07705 and 19K06352) of Japan's Society for the Promotion of Science.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
References
1. Acharya UR, Joseph KP, Kannathal N, Lim CM, Suri JS. Heart rate variability: a review. Med Biol Eng Comput. (2006) 44:1031–51. doi: 10.1007/s11517-006-0119-0
2. von Borell E, Langbein J, Després G, Hansen S, Leterrier C, Forde JM, et al. Heart rate variability as a measure of autonomic regulation of cardiac activity for assessing stress and welfare in farm animals – a review. Physiol Behav. (2007) 92:293–316. doi: 10.1016/j.physbeh.2007.01.007
3. Khor KH, Shiels IA, Campbell FE, Greer RM, Rose A, Mills PC. Evaluation of a technique to measure heart rate variability in anaesthetised cats. Vet J. (2014) 199:229–35. doi: 10.1016/j.tvjl.2013.11.006
4. Wormald D, Lawrence AJ, Carter G, Fisher AD. Reduced heart rate variability in pet dogs affected by anxiety-related behaviour problems. Physiol Behav. (2017) 168:122–7. doi: 10.1016/j.physbeh.2016.11.003
5. Grandi LC, Ishida H. The physiological effect of human grooming on the heart rate and the heart rate variability of laboratory non-human primates: a pilot study in male rhesus monkeys. Front Vet Sci. (2015) 2:50. doi: 10.3389/fvets.2015.00050
6. Carravieri A, Müller MS, Yoda K, Hayama S, Yamamoto M. Dominant parasympathetic modulation of heart rate and heart rate variability in a wild-caught seabird. Physiol Biochem Zool. (2016) 89:263–276. doi: 10.1086/686894
7. Kovács L, Jurkovich V, Bakony M, Póti P, Szenci O, Tözsér J. Welfare assessment in dairy cattle by heart rate and heart rate variability- literature review and implications for future research. Animal. (2014) 8:316–30. doi: 10.1017/S1751731113002140
8. Sutherland MA, Rogers AR, Verkerk GA. The effect of temperament and responsiveness towards humans on the behavior, physiology and milk production of multi-parous dairy cows in a familiar and novel milking environment. Physiol Behav. (2012) 107:329–37. doi: 10.1016/j.physbeh.2012.07.013
9. Frondelius L, Järvenranta K, Koponen T, Mononen J. The effect of body posture and temperament on heart rate variability in dairy cows. Physiol Behav. (2015) 139:427–41. doi: 10.1016/j.physbeh.2014.12.002
10. Kovács L, Kézér FL, Tozsér J, Szenci O, Póti P, Pajor F. Heart rate and heart rate variability in dairy cows with different temperament and behavioral reactivity to humans. PLoS ONE. (2015) 10:e0136294. doi: 10.1371/journal.pone.0136294
11. Johns J, Patt A, Hillmann E. Do bells affect behavior and heart rate variability in grazing dairy cows? PLoS ONE. (2015) 10:e0131632. doi: 10.1371/journal.pone.0131632
12. Kovács L, Tozsér J, Szenci O, Póti P, Kézér FL, Ruff F, et al. Cardiac responses to palpation per rectum in lactating and nonlactating dairy cows. J Dairy Sci. (2014) 97:6955–63. doi: 10.3168/jds.2014-8327
13. Kovács L, Kézér FL, Kulcsár-Huszenicza M, Ruff F, Szenci O, Jurkovich V. Hypothalamic-pituitary-adrenal and cardiac autonomic responses to transrectal examination differ with behavioral reactivity in dairy cows. J Dairy Sci. (2016) 99:7444–57. doi: 10.3168/jds.2015-10454
14. Mohr E, Langbein J, Nürnberg G. Heart rate variability: a noninvasive approach to measure stress in calves and cows. Physiol Behav. (2002) 75:251–9. doi: 10.1016/S0031-9384(01)00651-5
15. Kovács L, Kézér FL, Ruff F, Szenci O. Cardiac autonomic activity has a circadian rhythm in summer but not in winter in non-lactating pregnant dairy cows. Physiol Behav. (2016) 155:56–65. doi: 10.1016/j.physbeh.2015.11.031
16. Hagen K, Langbein J, Schmind C, Lexer D, Waiblinger S. Heart rate variability in dairy cows – influences of breed and milking system. Physiol Behav. (2005) 85:195–204. doi: 10.1016/j.physbeh.2005.03.019
17. Kézér FL, Kovács L, Tozsér J. Step behavior and autonomic nervous system activity in multiparous dairy cows during milking in a herringbone milking system. Animal. (2015) 9:1393–6. doi: 10.1017/S1751731115000130
18. Kovács L, Tozsér J, Kézér FL, Ruff F, Aubin-Wodala M, Albert E, et al. Heart rate and heart rate variability in multiparous dairy cows with unassisted calvings in the periparturient period. Physiol Behav. (2015) 139:281–9. doi: 10.1016/j.physbeh.2014.11.039
19. Trenk L, Kuhl J, Aurich J, Aurich C, Nagel C. Heart rate and heart rate variability in pregnant dairy cows and their fetuses determined by fetomaternal electrocardiography. Theriogenol. (2015) 84:1405–10. doi: 10.1016/j.theriogenology.2015.07.027
20. Kovács L, Kézér FL, Ruff F, Szenci O. Timing of obstetrical assistance affects peripartal cardiac autonomic function and early maternal behavior of dairy cows. Physiol Behav. (2016) 165:202–10. doi: 10.1016/j.physbeh.2016.08.001
21. Nagel C, Trenk L, Aurich C, Ille N, Pichler M, Drillich M, et al. Sympathoadrenal balance and physiological stress response in cattle at spontaneous and PGF2α-induced calving. Theriogenology. (2016) 85:979–85. doi: 10.1016/j.theriogenology.2015.11.009
22. Kovács L, Kézér FL, Jurkovich V, Kulcsár-Huszenicza M, Tozsér J. Heart rate variability as an indicator of chronic stress caused by lameness in dairy cows. PLoS ONE. (2015) 10:e0134792. doi: 10.1371/journal.pone.0134792
23. Konold T, Bone GE, Simmons MM. Time and frequency domain analysis of heart rate variability in cattle affected by bovine spongiform encephalopathy. BMC Res Notes. (2011) 4:259. doi: 10.1186/1756-0500-4-259
24. Laborde S, Mosley E, Thayer JF. Heart rate variability and cardiac vagal tone in psychophysiological research–recommendations for experiment planning, data analysis, and data reporting. Front Psychol. (2017) 8:213. doi: 10.3389/fpsyg.2017.00213
25. Brouwer AM, van Dam E, van Erp JBF, Spangler DP, Brooks JR. Improving real-life estimates of emotion based on heart rate: a perspective on taking metabolic heart rate into account. Front Hum Neurosci. (2018) 12:284. doi: 10.3389/fnhum.2018.00284
26. Kovács L, Kézér FL, Bakony M, Hufnágel L, Tozsér J, Jurkovich V. Associations between heart rate variability parameters and housing- and individual- related variables in dairy cows using canonical correspondence analysis. PLoS ONE. (2015) 10:e0145313. doi: 10.1371/journal.pone.0145313
27. Oishi K, Himeno Y, Miwa M, Anzai H, Kitajima K, Yasunaka Y, et al. Correcting the activity-specific component of heart rate variability using dynamic body acceleration under free-moving conditions. Front Physiol. (2018) 9:1063. doi: 10.3389/fphys.2018.01063
28. Ravagnolo O, Misztal I. Effect of heat stress on nonreturn rate in Holsteins: fixed-model analysis. J Dairy Sci. (2002) 85:3101–6. doi: 10.3168/jds.S0022-0302(02)74397-X
29. Finocchiaro R, van Kaam JBCHM, Portolano B, Misztal I. Effect of heat stress on production of Mediterranean dairy sheep. J Dairy Sci. (2005) 88:1855–64. doi: 10.3168/jds.S0022-0302(05)72860-5
30. Garcia AB, Angeli N, Machado L, de Cardoso FC, Gonzalez F. Relationships between heat stress and metabolic and milk parameters in dairy cows in southern Brazil. Trop Anim Health Prod. (2015) 47:889–94. doi: 10.1007/s11250-015-0804-9
31. Wheelock JB, Rhoads RP, VanBaale MJ, Sanders SR, Baumgard LH. Effects of heat stress on energetic metabolism in lactating Holstein cows. J Dairy Sci. (2010) 93:644–55. doi: 10.3168/jds.2009-2295
32. Kadzere CT, Murphy MR, Slanikove N, Maltz E. Heat stress in lactating dairy cows: a review. Livest Prod Sci. (2002) 77:59–91. doi: 10.1016/S0301-6226(01)00330-X
33. West JW. Effects of heat-stress on production in dairy cattle. J Dairy Sci. (2003) 86:2131–44. doi: 10.3168/jds.S0022-0302(03)73803-X
34. Marai IFM, El-Darawany AA, Fadiel A, Abdel-Hafez MAM. Physiological traits as affected by heat stress in sheep – a review. Small Rumin Res. (2007) 71:1–12. doi: 10.1016/j.smallrumres.2006.10.003
35. Rhoads ML, Rhoads RP, VanBaale MJ, Collier RJ, Sanders SR, Weber WJ, et al. Effects of heat stress and plane of nutrition on lactating Holstein cows: I. Production, metabolism, and aspects of circulating somatotropin. J Dairy Sci. (2009) 92:1986–97. doi: 10.3168/jds.2008-1641
36. Bernabucci U, Biffani S, Buggiotti L, Vitali A, Lacetera N, Nardone A. The effect of heat stress in Italian Holstein dairy cattle. J Dairy Sci. (2014) 97:471–86. doi: 10.3168/jds.2013-6611
37. Silanikove N. Effects of heat stress on the welfare of extensively managed domestic ruminants. Livest Prod Sci. (2000) 67:1–18. doi: 10.1016/S0301-6226(00)00162-7
38. Sevi A, Annicchiarico G, Albenzio M, Taibi L, Muscio A, Dell'Aquila S. Effects of solar radiation and feeding time on behavior, immune response and production of lactating ewes under high ambient temperature. J Dairy Sci. (2001) 84:629–40. doi: 10.3168/jds.S0022-0302(01)74518-3
39. Amaral BC, Connor EE, Tao S, Hayen MJ, Bubolz JW, Dahl GE. Heat stress abatement during the dry period influences metabolic gene expression and improves immune status in the transition period of dairy cows. J Dairy Sci. (2011) 94:86–96. doi: 10.3168/jds.2009-3004
40. Dahl GE, Tao S, Laporta J. Heat stress impacts immune status in cows across the life cycle. Front Vet Sci. (2020) 7:116. doi: 10.3389/fvets.2020.00116
41. Joy A, Dunshea FR, Leury BJ, Clarke I, DiGiacomo K, Chauhan SS. Resilience of small ruminants to climate change and increased environmental temperature: a review. Animals. (2020) 10:867. doi: 10.3390/ani10050867
42. Bun C, Watanabe Y, Uenoyama Y, Inoue N, Ieda N, Matsuda F, et al. Evaluation of heat stress response in crossbred dairy cows under tropical climate by analysis of heart rate variability. J Vet Med Sci. (2018) 80:181–5. doi: 10.1292/jvms.17-0368
43. Kovács L, Kézér FL, Ruff F, Jurkovich V, Szenci O. Heart rate, cardiac vagal tone, respiratory rate, and rectal temperature in dairy calves exposed to heat stress in a continental region. Int J Biometeorol. (2018) 62:1791–7. doi: 10.1007/s00484-018-1581-8
44. Miwa M, Oishi K, Nakagawa Y, Maeno H, Anzai H, Kumagai H, et al. Application of overall dynamic body acceleration as a proxy for estimating the energy expenditure of grazing farm animals: relationship with heart rate. PLoS ONE. (2015) 10:e0128042. doi: 10.1371/journal.pone.0128042
45. Mader TL, Davis MS, Brown-Brandl T. Environmental factors influencing heat stress in feedlot cattle. J Anim Sci. (2006) 84:712–9. doi: 10.2527/2006.843712x
46. Bernardi L, Valle F, Coco M, Calciati A, Sleight P. Physical activity influences heart rate variability and very-low-frequency components in Holter electrocardiograms. Cardiovasc Res. (1996) 32:234–7. doi: 10.1016/0008-6363(96)00081-8
47. Rainville P, Bechara A, Naqvi N, Damasio AR. Basic emotions are associated with distinct patterns of cardiorespiratory activity. Int J Psychophysiol. (2006) 61:5–18. doi: 10.1016/j.ijpsycho.2005.10.024
48. Hansen AL, Johnson BH, Thayer JF. Vagal influence on working memory and attention. Int J Psychophysiol. (2003) 48:263–74. doi: 10.1016/S0167-8760(03)00073-4
49. Johnson BH, Thayer JF, Laberg JC, Wormnes B, Raadal M, Skaret E, et al. Attentional and physiological characteristics of patients with dental anxiety. J Anxiety Disord. (2003) 17:75–87. doi: 10.1016/S0887-6185(02)00178-0
50. Verkuil B, Brosschot JF, Tollenaar MS, Lane RD, Thayer JF. Prolonged non-metabolic heart rate variability reduction as a physiological marker of psychological stress in daily life. Ann Behav Med. (2016) 50:704–14. doi: 10.1007/s12160-016-9795-7
51. Halsey LG, Shepard ELC, Quintana F, Laich AG, Green JA, Wilson RP. The relationship between oxygen consumption and body acceleration in a range of species. Comp Biochem Physiol A Mol Integr Physiol. (2009) 152:197–202. doi: 10.1016/j.cbpa.2008.09.021
52. Miwa M, Oishi K, Anzai H, Kumagai H, Ieiri S, Hirroka H. Estimation of the energy expenditure of grazing ruminants by incorporating dynamic body acceleration into a conventional energy requirement system. J Anim Sci. (2017) 95:901–9. doi: 10.2527/jas2016.0749
53. Qasem L, Cardew A, Wilson A, Griffiths I, Halsey LG, Shepard ELC, et al. Tri-axial dynamic acceleration as a proxy for animal energy expenditure; should we be summing values or calculation the vector? PLoS ONE. (2012) 7:e31187. doi: 10.1371/journal.pone.0031187
54. Tarvainen MP, Niskanen JP, Lipponen JA, Ranta-aho PO, Karjalainen PA. Kubios HRV-Heart rate variability analysis software. Comput Methods Programs Biomed. (2014) 113:210–20. doi: 10.1016/j.cmpb.2013.07.024
55. Task Force of the European Society of Cardiology North American Society of Pacing and Electrophysiology. Heart rate variability: standards of measurement, physiological interpretation, and clinical use. Circulation. (1996) 93:1043–65. doi: 10.1161/01.CIR.93.5.1043
56. Akselrod S, Gordon D, Ubel FA, Shannon DC, Barger AC, Cohen RJ. Power spectrum analysis of heart rate fluctuation: a quantitative probe of beat-to-beat cardiovascular control. Science. (1981) 213:220–2. doi: 10.1126/science.6166045
57. Malliani A, Pagani M, Lombardi F, Cerutti S. Cardiovascular neural regulation explored in the frequency domain. Circulation. (1991) 84:482–92. doi: 10.1161/01.CIR.84.2.482
58. Filligoi GC, Felici F. Detection of hidden rhythms in surface EMG signals with a non-linear time-series tool. Med Eng Phys. (1999) 21:439–48. doi: 10.1016/S1350-4533(99)00073-9
59. Eckmann JP, Kamphorst SO, Ruelle D. Recurrence plots of dynamical systems. Europhys Lett. (1987) 4:973–7. doi: 10.1209/0295-5075/4/9/004
60. Zbilut JP, Thomasson N, Webber CL. Recurrence quantification analysis as a tool for nonlinear exploration of nonstationary cardiac signals. Med Eng Phys. (2002) 24:53–60. doi: 10.1016/S1350-4533(01)00112-6
61. Dabiré H, Mestivier D, Jarnet J, Safar ME, Chau NP. Quantification of sympathetic and parasympathetic tones by nonlinear indexes in normotensive rats. Am J Physiol Heart Circ Physiol. (1998) 275:1290–7. doi: 10.1152/ajpheart.1998.275.4.H1290
62. Kramer CY. Extension of multiple range test to group means with unequal numbers of replications. Biometrics. (1956) 12:307–10. doi: 10.2307/3001469
63. SAS Institute Inc. SAS/STAT User's Guide: Version 9. 2. North Carolina: SAS Institute Inc. (2008).
64. Pagani M, Lombardi F, Guzzetti S, Rimoldi O, Furlan R, Pizzinelli P, et al. Power spectral analysis of heart rate and arterial pressure variabilities as a marker of sympatho-vagal interaction in man and conscious dog. Circ Res. (1986) 59:178–93. doi: 10.1161/01.RES.59.2.178
65. Machando NAF, Barbosa-Filho JAD, Oliveira KPL, Parente MOM, Siqueira JC, Pereira AM, et al. Biological rhythm of goats and sheep in response to heat stress. Biol. Rhythm Res. (2020) 51:1044–52. doi: 10.1080/09291016.2019.1573459
66. Johnson KG. Body temperature lability in sheep and goats during short-term exposures to heat and cold. J Agric Sci Camb. (1971) 77:267–72. doi: 10.1017/S0021859600024412
67. Penttilä J, Helminen A, Jartti T, Kuusela T, Huikuri HV, Tulppo MP, et al. Time domain, geometrical and frequency domain analysis of cardiac vagal outflow: effects of various respiratory patterns. Clin Physiol Funct Imaging. (2001) 21:365–76. doi: 10.1046/j.1365-2281.2001.00337.x
68. Billman GE. The LF/HF ratio does not accurately measure cardiac sympatho-vagal balance. Front Physiol. (2013) 4:26. doi: 10.3389/fphys.2013.00026
69. Laborde S, Mosley E. Commentary: heart rate variability and self-control – a meta-analysis. Front Psychol. (2016) 7:653. doi: 10.3389/fpsyg.2016.00653
70. Laborde S, Mosley E, Mertgen A. A unifying conceptual framework of factors associated to cardiac vagal control. Heliyon. (2018) 4:e01002. doi: 10.1016/j.heliyon.2018.e01002
71. Kanters JK, Højgaard MV, Agner E, Holstein-Rathlou NH. Influence of forced respiration on nonlinear dynamics in heart rate variability. Am J Physiol Reg Int Comp Physiol. (1997) 272:R1149–54. doi: 10.1152/ajpregu.1997.272.4.R1149
72. Young H, Benton D. We should be using nonlinear indices when relating heart-rate dynamics to cognition and mood. Sci Rep. (2015) 5:16619. doi: 10.1038/srep16619
73. Hoffmann G, Herbut P, Pinto S, Heinicke J, Kuhla B, Amon T. Animal-related, non-invasive indicators determining heat stress in dairy cows. Biosyst Eng. (2019) 199:83–96. doi: 10.1016/j.biosystemseng.2019.10.017
74. Tanoue A, Nasa Y, Koshimizu T, Shinoura H, Oshikawa S, Kawai T, et al. The α1D-adrenergic receptor directly regulates arterial blood pressure via vasoconstriction. J Clin Invest. (2002) 109:765–75. doi: 10.1172/JCI200214001
75. Tansey EA, Johnson CD. Recent advances in thermoregulation. Adv Physiol Educ. (2015) 39:139–48. doi: 10.1152/advan.00126.2014
76. Thompson VA, Barioni LG, Oltjen JW, Rumsey T, Fadel JG, Sainz RD. Development of a heat balance model for cattle under hot conditions. In: Sauvant D, van Milgen J, Faverdin P, Friggens N, editors. Modelling Nutrient Digestion and Utilization in Farm Animals. Wageningen: Wageningen Academic Publishers (2011). p. 243–51.
Keywords: physical activity, non-linear analysis, heat stress, heart rate variability, dynamic body acceleration
Citation: Kitajima K, Oishi K, Miwa M, Anzai H, Setoguchi A, Yasunaka Y, Himeno Y, Kumagai H and Hirooka H (2021) Effects of Heat Stress on Heart Rate Variability in Free-Moving Sheep and Goats Assessed With Correction for Physical Activity. Front. Vet. Sci. 8:658763. doi: 10.3389/fvets.2021.658763
Received: 26 January 2021; Accepted: 27 April 2021;
Published: 01 June 2021.
Edited by:
Jan Langbein, Leibniz Institute for Farm Animal Biology (FBN), GermanyReviewed by:
Viktor Jurkovich, University of Veterinary Medicine Budapest, HungaryAnnika Krause, Leibniz Institute for Farm Animal Biology (FBN), Germany
Copyright © 2021 Kitajima, Oishi, Miwa, Anzai, Setoguchi, Yasunaka, Himeno, Kumagai and Hirooka. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Kazato Oishi, oishi.kazato.3m@kyoto-u.ac.jp