An Explainable Machine Learning Model to Predict the Effects of Training and Match Load on Heart Rate Variability in Semi-Professional Basketball Players

1. Introduction

In team sports, where the competitive period spans several months, a well-structured approach to load quantification is essential [1,2]. During a competitive microcycle, sport scientists must monitor load and recovery to promote adaptation, minimize fatigue, and reduce injury risk in team-sports athletes [1,3]. Training load can be categorized into external and internal components, with external load referring to the work performed by the athlete during training or competition, and internal load referring primarily to the physiological response [2,4,5]. Basketball involves intermittent high-intensity actions and diverse player profiles across positions, each with distinct physical and tactical demands [6,7,8]. Combined with a congested schedule, this underscores the need for individual monitoring of internal and external load. A common method for load monitoring is heart rate (HR)-based measurements [9,10]. Metrics derived from HR monitoring, such as heart rate variability (HRV), can provide insights into the status of the autonomic nervous system and its potential relationship with fatigue [9,10]. HRV, sensitive to the prior day’s stimulus [4,11,12,13], can support assessment of physiological response and optimize the dose-response process throughout the competitive microcycle in basketball players.

HRV refers to the fluctuation in the time intervals between consecutive heartbeats [14]. It responds differently to high- and low-intensity stimuli [15], and exhibits individual variation [12,16,17], highlighting the need for individualized monitoring and analysis for each athlete. Among various HRV metrics, the Root Mean Square of Successive Differences (RMSSD) and its logarithmic transformation (LnRMSSD) have been widely used for fatigue monitoring in team sports [11,12,16,18,19]. However, individual responses to load must be considered, and individual LnRMSSD values have been used to observe training responses, highlighting the athlete heterogeneity in team sports [12,18]. HRV-based load modulation usually involves adjusting training based on the morning HRV reading [13], but this is often limited by time constraints in team settings. Therefore, predictive models of HRV responses to training loads could improve planning and decision-making during the competitive microcycle. Additionally, given the length of the competition period in team sports, considering HRV trends and load evolution across the entire season is also important. However, most studies have focused on partial periods such as pre-seasons or short cycles [11,17,18,20]. Likewise, there does not appear to be a clear consensus in the literature regarding the relationship between HRV and other load-related variables, as this relationship tends to be highly individual and vary depending on the monitoring method used [11,18,19,21,22]. Consequently, more research on predictive models is needed to better characterize these complex relationships and enhance proactive load management in basketball players.

In this context, advanced methods such as the SHapley Additive exPlanation (SHAP) technique offer a unique opportunity to improve the explainability of competitive black-box algorithms such as Gradient Boosting Models by identifying the contributing variables and enabling timely and informed adjustments [23]. This one can be used to relate a response variable, such as LnRMSSD with load variables that mostly influence on athlete’s fatigue reflected in HRV changes. SHAP technique provides values (SHAP values), which would allow to know the contribution (positive or negative) of each load variable to the LnRMSSD prediction given by the model for each athlete and, globally, the importance and contribution of each variable load in the model. This is particularly valuable in team sports, where the interaction of multiple complex variables can complicate data-driven decision-making.

Therefore, the aims of this study were: i) To assess the effect of the previous day’s stimulus on HRV; ii) To build a model that explains which load variables affect HRV in basketball players, based on monitoring them throughout an entire season.; iii) To identify individual profiles regarding the specific factors that influence HRV changes. We hypothesized that the training or match load from the previous day could influence the HRV response on the following day, and that specific load variables would explain this modification to a greater extent. Additionally, we expected to find individual differences among athletes in the factors influencing these responses.

2. Materials and Methods

2.1. Study Design

A longitudinal observational study with a repeated measures design was conducted, recording training load variables, rate of perceived exertion (RPE), training and match volume and HRV throughout a semi-professional basketball season in Spain (2021-2022). Daily recordings were performed, and only players who completed the full monitoring were included in the statistical analysis. Data were collected over a period of seven months (October to April) by the strength and conditioning coach as part of their professional development.

2.2. Participants

Five semi-professional basketball players (age: 23±5 years; height: 191 ± 7 cm; body mass: 90 ±11 kg) voluntarily participated in this study. Based on playing positions, they were grouped as: 2 point guards, 2 wings, and 1 center. Each player contributed an average of 212.6 ± 23.9 HR and HRV recordings and 257 RPE and volume recordings. The typical weekly schedule during the season included four training sessions and one match. On average, each player had 98 ± 9.77 training session recordings and 26.6 ± 2.19 match recordings. All participants were fully informed about the study procedures and gave their consent for the use of their data. Additionally, written informed consent was obtained from the club to analyze the data collected throughout the season. As data used in this study were collected as part of routine player monitoring, no ethics committee approval was required [24].

2.3. Procedures

2.3.1. Training Days Classification (MD-TD-NTD)

During the season, days were classified into three categories based on their content: match day (MD), training day (TD), or non-training day (NTD). A match day was defined as any day the player participated in at least one minute of an official game. A training day was recorded when there was a team training session in which the player participated fully. A non-training day was recorded when the player did not participate in either a game or training session that day. The typical in-season microcycle was structured as shown in Figure 1.

2.3.2. Training Load

The primary load metrics recorded were training session duration or match minutes played, and subjective effort assessment via RPE. From these two metrics, session-RPE (sRPE) was calculated by multiplying volume and RPE. RPE is a valid tool for recording training load[25]. All metrics used can be found in Table A1. RPE data were collected using the TrainingFeel app (https://trainingfeel.com/), which enables automated reporteing of Borg’s 1-10 scale by each player on their own mobile device. Data collection was preceded by a two-week familiarization phase and conducted individually. Each player answered the question “How intense did the session feel?” on their mobile device 15-30 minutes post-session. Session duration for each player was monitored by the team’s strength and conditioning coach, while match minutes were obtained from the official records on the Spanish Basketball Federation’s website (https://www.feb.es/).

For each of these variables, the moving average, weighted average, and exponentially weighted moving average (EWMA) were calculated[26]. A more detailed explanation of the EWMA can be found in Table A1. The average was calculated as the mean of the last few days, with the number of days included indicated in the variable name (e.g., sRPE_avg4). The weighted average was calculated by weighting the same average, giving a higher value to the most recent day (4 in the case of a four-day period) and reducing by one for each preceding day. The EWMA was calculated using the formula previously applied in other studies on load monitoring in sports [27].

2.3.3. Heart Rate Variability

To monitor load responses, HR data were collected during training sessions and matches [10]. From these, HR, RMSSD, and LnRMSSD were derived. HRV monitoring, particularly through LnRMSSD, has been used by various authors as a method to assess the response to training load [9,16].

For HR data collection, participants downloaded the validated HRV4Training app [28] (http://www.hrv4training.com/) and consented to data recording and processing for professional development and research purposes via the app itself. Each morning, upon waking, HRV was measured in a one-minute supine recording to avoid external stressors [29]. A measurement was considered valid if confirmed by the app’s algorithm. Ultrashort recordings have shown validity and are sensitive to the prior day’s stimulus [30]. The measurement was conducted using the photoplethysmography (PPG) method. PPG was used as a non-invasive and validated method to assess autonomic nervous system status[30], with participants placing their finger over the mobile device’s camera and flash. Data were uploaded automatically to a cloud server and accessed by the strength and conditioning coach.

2.4. Statistical Analyses

Descriptive data are presented as mean ± standard deviation. Firstly, the linear mixed-effects model analysis was conducted to examine the effects of different types of days (MD, TD or NTD) on HRV, specifically LnRMSSD. The day category was set as a fixed effect factor, whereas the intercept was allowed to vary randomly across participants to account for individual differences in LnRMSSD for the reference condition (i.e., NTD). Post-hoc contrast were performed by paired t-tests with Bonferroni’s correction. Normality of the residuals was visually inspected (i.e. histogram and Q-Q plots) whereas kurtosis and skewness statistics were also considered.

On the other hand, we aimed to adjust a model that allowed establishing the dependency relationships between the LnRMSSD variable and load indicators of basketball players. Although there are quite popular and interpretable techniques such as decision trees and regression models, the fact that the relationships involved in the current study can be especially complex led us to opt for more sophisticated models. Among them, we used an Extreme Gradient Boosting (XGBoost), which consists of the sequential adjustment of hundreds of decision trees, where each one assigns greater weight to the patterns worst predicted by the previous trees. This is why it is considered as a black-box type model, whose results are competitive, from a predictive point of view. In this context, the SHAP technique is applied to quantify the contribution of each variable in the model (i.e., load variables; see Table A1) to the predicted outcome (LnRMSSD) through the computation of SHAP values. Furthermore, this interpretation can be done both at a global level, to identify the most relevant variables to characterize the LnRMSSD, and at a local level, that is, at the athlete level. The latter is decisive to be able to make personalized decisions about each of them. Thus, one variable may be the one that generally contributes the most to an increase or decrease in the LnRMSSD, and yet another variable may be the one that is most relevant for a specific athlete. In this study normal level of LnRMSSD was defined as the average of NTD ± 0.5 SD within the last 30 days and therefore being updated daily.

LnRMSSD was estimated using an Extreme Gradient Boosting Model. Regarding the model configuration, hyperparameters were tuned using grid search five-fold cross-validation minimizing the root mean squared error (RMSE). The final model included 400 trees, a learning rate of 0.01, a maximum depth of 4 and a subsample of 0.6 as hyperparameters. To obtain explanations of the features that drive players-dates estimates, we used a SHAP algorithm. SHAP plots provided a clear and intuitive way to visualize the contributions of various features to our model’s estimates, proving to be an invaluable tool in understanding the factors that underline predictions related to LnRMSSD. The application of SHAP plots allowed us to analyze the impact of individual features on LnRMSSD estimates, encompassing both the direction and magnitude of their influence. In these plots, the magnitude of the value of the load variables is represented by colors, while the influence over the predicted given by the model is expressed by a Shapley value. In fact, these values represent the contribution (positive or negative) of each load variable to the LnRMSSD prediction given by the model.

By closely examining SHAP values, this analytical approach facilitated the identification of the input variables in the estimate of LnRMSSD. Thus, once a variable and an observation have been fixed, it will be represented by a dot. This dot will be associated with a more intense red color, the higher the value of the variable on that observation, and a more intense blue color the lower the value. In addition, the higher the Shapley value (i.e., the farther to the right the observation is represented), the greater the contribution of the variable to the increase in the estimate given by the model for that observation. Conversely, the lower the Shapley value (i.e., the further to the left the observation is represented), the greater the contribution of the variable to the decrease in the estimate given by the model for that observation. In consequence, the sign informs whether the variable contributes positively or negatively to the estimate made for the observation. According to their global influence on the model, the variables will be presented from greater to lesser importance on the Y axis.

On the other hand, linear mixed model analysis was carried out using SPSS software v.28 (IBM Corp., Armonk, NY, USA), with a significance level set at p < 0.05, whereas the explanatory framework model was performed using the Python package XGBoost in Python 3.11.3 (Python Software Foundation, Oregon, USA)

3. Results

Values of LnRMSSD were analyzed according to the categorization of days using mixed models. The estimated marginal means for each type of day were: 4.711 ± 0.146 ln (ms) for NTD, 4.653 ± 0.146 ln (ms) for TD, and 4.408 ± 0.149 ln (ms) for MD.

Pairwise comparisons (Figure 2) revealed significant differences between NTD and MD ( mean difference = 0.304, p < 0.001) and between TD and MD (mean difference = 0.245, p < 0.001), whereas non-significant differences were observed between NTD and TD (mean difference = 0.058, p = 0.578).

The SHAP analysis shows the contribution of different variables to the modification of LnRMSSD (Figure 3). The most impactful variables include RPE, the number of days since the last match (DaysLastMatch), minutes played in the last match (Volume_LastMatch), the average sRPE of the last four days, and the RPE of the last match (RPE_LastMatch), among others. Additionally, variables such as pre-training or pre-match RMSSD and HR also influence the modification of LnRMSSD on the following day.

Figure 4 displays the individual analyses, highlighting inter-individual differences in the order, weight, and impact of the variables associated with HRV changes for each player.

4. Discussion

The main findings of this study are: i) HRV is influenced by the training load from the previous day, with matches being the stimulus that induces the greatest change; ii) intensity metrics appear to be one of the most important variables explaining HRV variation on the following day; iii) individual differences exist regarding which variables most influence HRV modification in each player.

Previous research has shown that HRV is influenced by training load in team-sport athletes[12,19]. Consistent with these findings, our results indicate that LnRMSSD is affected by training load in semiprofessional basketball players. However, those studies focused on short-term periods [11,12,18] such as a single competitive microcycle or preaseason. In contrast, the present study monitored players across a full competitive season. Specifically, our results show a significant difference between MD and both TD and NTD, while no differences were found between TD and NTD. This may be due to the variability of training sessions, which ranged from low to high intensity. Nakamura et al. [31] reported a decrease in LnRMSSD the morning following matches during a beach volleyball tournament, though they did not compare this response to training sessions, leaving uncertainty about HRV responses to different stimuli. Our findings appear to support structuring the microcycle around the match and its proximity, as the match appears to elicit the largest suppression of LnRMSSD. One of the limitations of previous studies is the lack of direct comparisons between physiological responses to matches and training. Future research should address this to determine whether matches consistently cause greater HRV reductions and to clarify the implications of this response.

Through SHAP analysis, RPE and days since the last match were identified as the two most influential variables in HRV changes. Consistent with previous research [15,30], perceived intensity through RPE appears to be the main factor influencing HRV [32]. In a previous study, O’Connor et al. [19] reported a relationship between sprint volume and a reduction in RMSSD upon waking the following day. Similarly, Stanley et al. [15], found that while HRV tends to recover within 24 hours following low-intensity exercise, high-intensity sessions suppress HRV for up to 48 hours post-exercise. Our findings, highlighting RPE as the main variable, align with these authors, reinforcing the impact of exercise intensity on next-day HRV suppression. The number of days since the last match ranks second, suggesting that having a match the previous day negatively affects HRV, which highlights the importance of organizing the days following the match due to their high-intensity nature, and the need for at least 48 hours for full parasympathetic cardiac reactivation [15]. These findings suggest that match intensity and recovery window should be key considerations in basketball microcycle planning. Minutes played in the match also significantly influenced HRV, with higher values of minutes played linked to decreased LnRMSSD the following day. We also identified a third key metric related to the match: the reported RPE metric that refers to the RPE reported in the last match played. Therefore, these results seem to highlight the importance of the match in HRV response, along with minutes played and match RPE as a key metrics for predicting players’ HRV responses. Since HRV can also be influenced by psychological variables[33], match demands and outcomes may play a role in the magnitude of these changes. The average sRPE of the last 4 days showed an inverse relationship with the next-day LnRMSSD. Since it links both volume and perceived intensity metrics, the response may be primarily influenced by the volume, as it takes higher values than RPE. This effect was clearer when analyzing the 7-day average volume, which appeared to positively influence HRV. This suggests that training exposure may help sustain favorable HRV values [21], although acute responses depend on session intensity [15,19]. Pre-stimulus variables such as HR and pre-match RMSSD also play a significant role in the response the following day. Higher resting HR or lower RMSSD values prior to training seem to influence a greater reduction in LnRMSSD the next day. Taken together, these results suggest the need to plan training loads during the microcycle by taking into account not only the match and the intensity of training sessions, but also the athlete’s physiological status prior to each session. Considering the clear impact that match-related variables and perceived intensity have on HRV modulation, practitioners should incorporate these findings into the weekly planning process to optimize the balance between stimulus and recovery and ultimately support performance and adaptation.

Individual differences were observed in the training load metrics associated with HRV responses. Similarly, Flatt et al. [16] reported in American football, differences among players. This highlights the importance of individualizing the monitoring process based on the specific characteristics of each athlete. Each player exhibits a unique profile in which specific metrics have varying degrees of importance and magnitude in HRV modification. Differences are observed between Player 1 (wing) and Player 5 (point guard) regarding the position and effect of the RMSSD_pre variable. For instance, in Player 1, RMSSD_pre ranks eighth in importance for HRV response, while in Player 5, it ranks third. Furthermore, in Player 5, there is a clear tendency where low RMSSD_pre values significantly and negatively influence HRV, whereas in Player 1, the effect is less pronounced, and it appears that high RMSSD_pre values have a negative influence. Differences are also observed in Player 3 (point guard), where the HR_pre metric holds high importance, ranking third. It shows a clear influence, with low values having a positive effect, while high or medium values have little to no effect, or even a negative influence. In the other players, this metric appears with lower importance, and its interpretation is less clear, with mixed influences from both high and low values. These findings suggest that individual profiles should be considered when managing training load and monitoring potential fatigue. The SHAP analysis provides an easy-to-interpret individual and global interpretation, allowing us to translate it into a practical application for coaches and strength and conditioning professionals.

To the best of the authors’ knowledge, this is the first study to apply SHAP analysis to model athlete fatigue, considering different types of variables recorded throughout an entire competitive season. This technique enables the development of both general and individual predictive models, rather than relying on the retrospective decision-making approach provided by other analytical methods. Our results showed that RPE and the number of days since the last match were the most predictive variables of player condition, as assessed through HRV. However, distinct individual profiles were also identified. Machine learning techniques, combined with interpretability methods such as SHAP, emerge as valuable tools for identifying predictive models of fatigue and performance in team sports, allowing for the anticipation of training load adjustments.

There are some limitations in the present study that should be acknowledged. Firstly, the sample size limits the study’s ability to draw broader conclusions, although a substantial amount of data per player was used to train the mathematical model. Secondly, we analyzed only a single season, whereas a larger longitudinal study would be beneficial to confirm the findings. Nevertheless, to our knowledge, this is the first study to use data from an entire competitive season to develop a predictive model of fatigue in team sports. Finally, training load was assessed using RPE, while integrating additional objective measures, such as accelerometers or other sensors, would be valuable for future research. However, RPE and sRPE have been shown to reflect internal load in team sports [8,34].

5. Conclusions

The results of the present study showed that HRV is significantly influenced by the training load from the previous day, with matches being the stimulus that induces the greatest change. RPE and the number of days since the last match, were identified as the main factors affecting LnRMSSD modification. Additionally, playing time in the last match and the physiological state before training or competition (HRV and HR) also influenced HRV responses on the following day. In addition, the individual analyses revealed differences among players regarding the most influential variables and their impact on HRV modification. These findings highlight the importance of considering individual profiles when managing training load and monitoring fatigue status in basketball players.

To the best of the authors’ knowledge, this study is the first to apply SHAP analysis to model athlete fatigue, as measured by HRV, over an entire competitive season in a team sport, enabling the identification of both general and individual response patterns to training load. The integration of machine learning techniques with interpretable methods such as SHAP provides valuable tools for anticipating training load adjustments, optimizing performance, and enhancing recovery in team sports.