Youri Duchene
email: youri.duchene@univ-lorraine.fr

Gérome C Gauchard

Guillaume Mornieux

Youri Duchene Orcid

Influence of sidestepping expertise and core stability on knee joint loading during change of direction

Keywords: trunk lateral lean, trunk rotation, pelvis, neuromuscular control, ACL

The aims of this study were twofold: firstly, to compare core stability and knee joint loading between sidestepping experts and nonexperts; secondly, to determine core predictors of knee joint loading. Thirteen handball male players (experts) and fourteen karatekas (nonexperts) performed six unanticipated 45° sidestepping manoeuvers, while trunk and pelvis 3D kinematics as well as ground reaction forces were measured, and peak knee abduction moment (PKAM) was determined. Student t-tests enabled a comparison of both groups and a linear mixed model approach was used to identify PKAM predictors. Sidestepping experts demonstrated significantly lower pelvis rotation towards the new movement direction at the initial contact than nonexperts (4.9° vs. 10.8°) and higher PKAM (0.539 vs. 0.321 Nm/kg-bwt). Trunk medial lean, trunk axial rotation and pelvis anterior tilt at the initial contact predicted PKAM, while trunk axial rotation, pelvis medial lean and posterior ground reaction force predicted PKAM during the weight acceptance phase. Despite higher PKAM, handball players might not be at a higher risk of anterior cruciate ligament injury as the knee joint loading remained at a relatively low level during this sidestepping task. Core stability, in its three dimensions, is a key determinant of knee joint loading.

Introduction

Anterior cruciate ligament (ACL) rupture is a major injury occurring mostly during noncontact changes of direction or landings [START_REF] Boden | Mechanisms of anterior cruciate ligament injury[END_REF]. Different biomechanical risk factors at the knee level [START_REF] Alentorn-Geli | Prevention of non-contact anterior cruciate ligament injuries in soccer players. Part 1 : Mechanisms of injury and underlying risk factors[END_REF] but also related to the trunk position and control [START_REF] Hughes | A review of recent perspectives on biomechanical risk factors associated with anterior cruciate ligament injury[END_REF], have been identified in the literature.

Indeed, an elevated knee joint loading possibly stems from higher lateral trunk motion [START_REF] Hewett | The mechanistic connection between the trunk, hip, knee, and anterior cruciate ligament injury[END_REF][START_REF] Jamison | Randomized controlled trial of the effects of a trunk stabilization program on trunk control and knee loading[END_REF] and is also associated with trunk rotation [START_REF] Critchley | The effects of mid-flight whole-body and trunk rotation on landing mechanics : Implications for anterior cruciate ligament injuries[END_REF][START_REF] Dempsey | The effect of technique change on knee loads during sidestep cutting[END_REF][START_REF] Dempsey | Whole body kinematics and knee moments that occur during an overhead catch and landing task in sport[END_REF]. Moreover, [START_REF] Zazulak | Deficits in neuromuscular control of the trunk predict knee injury risk[END_REF] highlighted a potential connection between core stability and ACL injury risk in female athletes. Recent studies provided an analysis of trunk neuromuscular control to further understand trunk motion during landings [START_REF] Haddas | Volitional spine stabilization during a drop vertical jump from different landing heights : Implications for anterior cruciate ligament injury[END_REF][START_REF] Lessi | Effects of fatigue on lower limb, pelvis and trunk kinematics and muscle activation : Gender differences[END_REF] and changes of direction [START_REF] Jamison | The effects of core muscle activation on dynamic trunk position and knee abduction moments : Implications for ACL injury[END_REF]Oliveira et al., 2013a;[START_REF] Oliveira | Effects of perturbations to balance on neuromechanics of fast changes in direction during locomotion[END_REF]. Accordingly, future research on trunk control during changes of direction and its consequences on knee joint injury risk should not only determine trunk kinematics but also analyze the related neuromuscular control, for instance based on trunk muscle co-contractions [START_REF] Jamison | The effects of core muscle activation on dynamic trunk position and knee abduction moments : Implications for ACL injury[END_REF].

ACL injury risk is dependent on different factors, including gender or the type of sport [START_REF] Beynnon | The effects of level of competition, sport, and sex on the incidence of first-time noncontact anterior cruciate ligament injury[END_REF][START_REF] Renstrom | Non-contact ACL injuries in female athletes : An international olympic committee current concepts statement[END_REF]. Moreover, an influence of the skill level on knee joint loading could be expected as high-level athletes are supposed to perform sidestepping changes of direction more efficiently [START_REF] Fujii | Superior reaction to changing directions for skilled basketball defenders but not linked with specialised anticipation[END_REF][START_REF] Lockie | Planned and reactive agility performance in semiprofessional and amateur basketball players[END_REF]. For instance, Sigward and Powers (2006b) reported higher knee joint moments for elite athletes, whereas greater knee abduction angle was found for novices but without any abduction moment difference in another study [START_REF] Kipp | Decision making and experience level influence frontal plane knee joint biomechanics during a cutting maneuver[END_REF]. Likewise, the influence of athletic background on knee joint biomechanics during various tasks associated with high knee joint loading remains controversial in the literature [START_REF] Cowley | Differences in neuromuscular strategies between landing and cutting tasks in female basketball and soccer athletes[END_REF][START_REF] Herrington | Knee valgus angle during landing tasks in female volleyball and basketball players[END_REF][START_REF] Orishimo | Comparison of landing biomechanics between male and female dancers and athletes, part 1 : influence of sex on risk of anterior cruciate ligament injury[END_REF].

In addition, it has been demonstrated that athletic background has an impact on trunk control during unstable sitting tasks and trunk perturbations [START_REF] Barbado | Sports-related testing protocols are required to reveal trunk stability adaptations in highlevel athletes[END_REF][START_REF] Glofcheskie | Athletic background is related to superior trunk proprioceptive ability, postural control, and neuromuscular responses to sudden perturbations[END_REF]. Differences between high-level sport athletes, i.e. judokas vs. kayakers [START_REF] Barbado | Sports-related testing protocols are required to reveal trunk stability adaptations in highlevel athletes[END_REF] and runners vs. golfers [START_REF] Glofcheskie | Athletic background is related to superior trunk proprioceptive ability, postural control, and neuromuscular responses to sudden perturbations[END_REF], were only noticed during specific tasks that were similar to those realized during their practice. Therefore, the comparison of task experts and nonexperts might give a more accurate insight of an adequate trunk control while performing changes of direction.

More generally, researchers seeking to further understand the knee joint injury risk and therefore predict peak knee abduction moment (PKAM) have recently provided models putting different trunk and lower limb variables together [START_REF] Frank | Trunk and hip biomechanics influence anterior cruciate loading mechanisms in physically active participants[END_REF][START_REF] Jones | Technique determinants of knee joint loads during cutting in female soccer players[END_REF][START_REF] Kristianslund | Sidestep cutting technique and knee abduction loading : Implications for ACL prevention exercises[END_REF][START_REF] Staynor | By failing to prepare, you are preparing your anterior cruciate ligament to fail[END_REF][START_REF] Weir | A reliable video-based ACL injury screening tool for female team sport athletes[END_REF]. While trunk lateral lean and trunk rotation were predictive variables for PKAM, little is known about the influence of the pelvis, despite its central role within the core stability. Also, while [START_REF] Jamison | The effects of core muscle activation on dynamic trunk position and knee abduction moments : Implications for ACL injury[END_REF] underlined the effect of erector spinae muscles on knee joint loading and trunk flexion, the whole comprehension of the influence of core neuromuscular control (including glutei muscles, erector spinae muscles and abdominal muscles) on trunk and pelvis 3D motion during sidestepping maneuvers is lacking. Such a modeling approach, using data from a wide range of level of expertise in changes of direction, might provide a better understanding of the effect of the core kinematics and neuromuscular control on injury predictive variables during sidestepping tasks and help improve prevention programs to decrease the knee injury risk.

Thus, the purpose of this study was i) to determine the difference in core stability and knee joint loading between task experts and nonexperts during sidestepping maneuvers, ii) to understand core parameters predicting peak knee abduction moment, and iii) to determine neuromuscular parameters controlling core kinematics. We hypothesized that experts will have their trunk more orientated in the targeted direction, a higher neuromuscular co-contraction ratio towards this direction, and higher PKAM than nonexperts. Moreover, we suggested that PKAM will be predicted by trunk lateral lean and trunk rotation, and that co-contraction ratio will explain core kinematics.

Methods

Participants

The sample size was estimated based upon a power analysis performed using a power of 80% and the effect sizes reported for trunk kinematics [START_REF] Glofcheskie | Athletic background is related to superior trunk proprioceptive ability, postural control, and neuromuscular responses to sudden perturbations[END_REF] and PKAM [START_REF] Sigward | Predictors of frontal plane knee moments during side-step cutting to 45° and 110° men and women : Implications for ACL injury[END_REF], Sigward et al., 2006b) differences between male and female as well as experts and novice (alpha level=0.05). Also, based on [START_REF] Staynor | By failing to prepare, you are preparing your anterior cruciate ligament to fail[END_REF] and [START_REF] Sigward | Predictors of frontal plane knee moments during side-step cutting to 45° and 110° men and women : Implications for ACL injury[END_REF], 62 trials and 22 participants are necessary to achieve an effect size of 0.3 with a power of 80% with maximum 9 predictors of PKAM in the regression models. Accordingly, at least 12 subjects performing 5 trials would be necessary in each group (120 trials in aggregate) and in line with other studies.

The task-expert group was composed of handball players, as sidestepping maneuver is an element of performance for this population. Karatekas were recruited for the nonexpert group.

The rationale for the choice of this second population was the rareness of sidestepping maneuvers but a high-level of trunk control in their practice. Comparing experts and nonexperts of the task would allow determining the practice-based behaviors to enhance performance and reduce knee joint loading. Thus, the best trunk and pelvis kinematics could be identified and be set as targets for future training programs. Thirteen male handball players (age: 21.5 ± 2.4 years old; height: 1.81 ± 0.06 m; mass: 74.6 ± 7.9 kg) and 14 male karatekas (age: 26.9 ± 8.2 years old; height: 1.77 ± 0.06 m; mass: 72.1 ± 6.3 kg) were included in the study. All participants had at least 10 years of experience in their respective sport. Handball players were playing at the nonprofessional national level and karatekas were black belt.

Participants were wearing their own inside sport shoes. No participant had reported an ACL injury or any lower limb or trunk injury during the past 2 months before the experimentation.

Prior to testing, all participants were informed about possible risks and gave written informed consent. All procedures were approved by the ethics committee Sud Mediterranee III (approval reference 2018.07.03 bis) and conformed to requirements of the Declaration of Helsinki.

Protocol

Participants were asked to perform three different changes of direction tasks on a force plate in a randomized order, including a crossover to -20° to the left, a straightforward deceleration and a sidestepping maneuver to 45° to the right after a dynamic two-step approach. Movement direction was indicated by a light signal triggered at the end of the two-step approach to create an unanticipated change of direction paradigm. Eighteen randomized trials were carried out (6 in each direction) with a 1-minute rest between trials. A trial was defined as successful if the left foot fully hit the force plate and the change of direction was conducted towards the light stimulus directly on the first step after force plate contact.

Moreover, in order to normalize electromyography (EMG) recordings, participants performed a maximal broad jump.

Material

The visual light stimulus was triggered automatically via Optojump ™ (Microgate, Bolzano, Italy) on the last step before the force plate contact (Bertec Corp, Columbus, Ohio). Force plate sampling frequency was of 1000 Hz. Light stimulus randomization was computed on Labview (version 15.0.1f1).

Three-dimensional trunk, pelvis and left-leg kinematics were recorded using reflective skin markers (14 mm) placed on anatomical landmarks with self-adhesive tape. Specifically, markers were placed on the upper sternum and the xiphoid process, the C7 vertebra, the anterior superior iliac spines, mid posterior superior iliac spines, the great trochanter, medial and lateral epicondyle of the knee, at the midway of the peroneus, the medial and lateral malleolus, the first and fifth metatarsal, and calcaneus. This marker placement was adapted from the Lyon Whole Body Model [START_REF] Dumas | Adjustments to McConville et al. and Young et al. body segment inertial parameters[END_REF]. Markers were captured with a 6-camera motion analysis system (SIMI Reality Motion Systems, Germany) with a sampling frequency of 100 Hz.

Surface electromyography recordings of core and hip muscles, sampled at 1000 Hz and synchronized with the motion analysis system, were obtained from the rectus abdominis (RAB), the external oblique (EOB), the erector spinae (ESP), the gluteus maximus (GMax) and the gluteus medius (GMed) of the right and left sides. The skin was first shaved, sanded and cleaned with alcohol to reduce impedance, then wireless surface electrodes (Trigno™, Delsys, Natick, MA, USA) were attached parallel to the muscle-fiber orientation, following SENIAM guidelines for sensor locations.

Data analysis

All data were only recorded during sidestepping maneuvers to 45°. Marker trajectory and ground reaction force signals were both filtered with a low-pass Butterworth filter (4th order, 15 Hz cutoff frequency) prior to calculating external joint moments with a standard inverse dynamics approach (custom development realized by Moveck Solution Inc.). Three-dimension trunk and pelvis kinematics were defined relative to the global coordinate system. Core kinematics was positive in the direction of the cut, i.e. a rotation towards the right in transverse plane, an anterior tilt in sagittal plane, a lean towards the right in frontal plane (medial lean). A lean towards the left in the frontal plane is therefore defined as a lateral lean. Kinematics and kinetics of left leg were analyzed to determine knee joint loading, with a focus on peak knee abduction moment (PKAM), as this latter parameter is associated with knee injury [START_REF] Hewett | Biomechanical measures of neuromuscular control and valgus loading of the knee predict anterior cruciate ligament injury risk in female athletes : A prospective study[END_REF]. All EMG data were band-pass filtered (10 Hz-500 Hz). Muscle activation of the broad jump was then low-pass filtered at 30 Hz prior to determining a reference peak value during the jump. Heart muscle electrical activity was filtered with a Butterworth high-pass filter at 30 Hz on the RAB and the EOB. Root mean square (RMS) values were determined during the preactivation phase (100 ms prior to the initial contact with the force plate) and during the weight acceptance (WA) phase, which represented the time between IC and the first trough in the vertical force. The activation of the different muscles was normalized to their peak filtered value recorded during the maximal broad jump. The simultaneous action of agonist and antagonist muscles was analyzed via directed co-contraction ratios (DCCR) during the preactivation phase and WA [START_REF] Donnelly | Changes in muscle activation following balance and technique training and a season of Australian football[END_REF]. DCCR were defined according to the anatomical function of the muscles with respect to motion during the change of direction.

Therefore, a positive ratio in trunk flexion/extension would indicate a higher co-contraction of RAB and EOB (agonists) than ESP (antagonists). In the frontal plane, right EOB and ESP were defined as agonists and left EOB and ESP were antagonists. According to the change of direction to the right, left EOB was defined as the agonist and right EOB was the antagonist with respect to trunk rotation. Although fewer muscles were recorded at the pelvic level, the anterior/posterior pelvic tilt could also be analyzed using directed co-contraction ratio, where a positive ratio would indicate a higher co-contraction of ESP (agonists) than RAB, EOB and GMax (antagonists). DCCR values were between -1 (towards antagonists) and 1 (towards agonists), while a zero would indicate equal activation of agonist and antagonist muscle groups.

Statistics

All variables were averaged across the six sidestepping trials for each participant. Results are presented as group mean (standard deviation, SD). Normal distribution and variance homogeneity were respectively verified via Lilliefors' and Levene's test. The influence of the sidestepping expertise (Handball vs. Karate) on the dependent variables during changes of direction was analyzed using independent Student's t-tests. The level of significance was set at p < 0.05. Moreover, for each significant variable, effect sizes (ES) were calculated with Cohen's d (0.2, 0.5 and 0.8 correspond to small, medium and large effects) and confidence intervals (CI) of 95% were reported.

A total of 162 trials were used to perform linear mixed model analyses. Such models take account of within-participant effect, allowing an analysis with individual trials as data points.

First, two linear mixed models were used to determine the association between PKAM (dependent variable) and independent variables. At IC, the predictors were 3D core kinematics (model no.1) and during WA, 3D GRF were added to 3D core kinematics (model no.2).

Sidestepping expertise was included as fixed effect covariates and subject effect as a random effect. Backwards elimination of independent variables (fixed effect) was performed by sequentially removing nonsignificant predictors until all predictors were significant. Then, kinematic predictors of PKAM were defined as dependent variables of linear mixed models at IC and during WA. Trunk DCCR values were the independent variables for trunk kinematics and RMS of the gluteal muscles and pelvis DCCR were the predictors for pelvis kinematics, during the preactivation phase(model no.1) and WA (model no.2). The same backwards elimination procedure was processed.

Standardized effects of each predictor were computed by multiplying the predictor standard deviation with its estimate. Estimate values were coefficients of independent variables in a linear mixed model. Standardized effects were expressed with respect to the mean of the dependent variable, i.e., the one SD change of the dependent variable when all other predictors were set to their average value [START_REF] Staynor | By failing to prepare, you are preparing your anterior cruciate ligament to fail[END_REF].

Results

Kinetics and kinematics

At IC, none of the trunk kinematical variables differed between handball players and karatekas (Table 1). Pelvis rotation to the right was greater for karatekas than handball players (p<0.01).

Its mean difference was -5.9° (CI = [-10.1; -1.7]; ES = 1.11).

Handball

Karate

Pooled subjects

Trunk medial lean -3.8 (2.6) -3.1 (3.4) -3.4 (4.4)

Trunk axial rotation -4.3 (6.3) -6.9 (5.5) -5.6 (7.0)

Trunk anterior tilt 12.9 (5.2) 16.0 (6.0) 14.5 (6.2)

Pelvis

Neuromuscular control

Neither DCCR values nor glutei RMS values were significantly different between handball players (HB) and karatekas (K) (Figure 1). During the preactivation phase, pelvis DCCR supported pelvis posterior tilt (HB: -0.15 ± 0.49; K: -0.22 ± 0.31) and during WA phase, pelvis DCCR remained negative (HB: -0.26 ± 0.37; K: -0.03 ± 0.36).

Knee joint loading and core kinematics predictors

Predictors of PKAM and predictors of significant kinematic variables are reported in Table 3.

In the first model, at IC, PKAM was negatively associated with trunk medial lean and increased with greater trunk axial rotation (towards the new direction) and pelvis anterior tilt. In the second model, during WA, PKAM increased with larger trunk axial rotation (towards the new direction) and peak posterior GRF, but was reduced when pelvis medial lean increased.

Moreover, the different kinematical variables could be predicted by trunk and pelvis neuromuscular control.

Discussion and Implications

The present study's objectives were to determine the differences of core stability and knee joint loading between experts and nonexperts of sidestepping maneuvers and the influence of core stability on a knee joint injury predictor. The main findings were: 1) that there was no evidence of major biomechanical differences or injury risk during a sidestepping maneuver between handball players and karateka, despite a higher PKAM and lower pelvis rotation towards the new movement direction for handball players; 2) that trunk and pelvis kinematics predicted PKAM; and 3) that directed co-contraction ratio could help in the understanding of trunk and pelvis kinematics.

Core stability and knee joint loading differences between experts and nonexperts

Despite significantly higher PKAM for handball players, which validated our hypothesis, it could not be considered as a higher risk of ACL injury, as its mean value remained relatively low and was in line with the literature [START_REF] Sigward | Predictors of frontal plane knee moments during side-step cutting to 45° and 110° men and women : Implications for ACL injury[END_REF]. Moreover, [START_REF] Hewett | Biomechanical measures of neuromuscular control and valgus loading of the knee predict anterior cruciate ligament injury risk in female athletes : A prospective study[END_REF] reported PKAM values on ACL injured females of around 0.73 Nm/kg-bwt during drop jump while uninjured females were at 0.31 Nm/kg-bwt. Furthermore, using a biomechanical model, [START_REF] Lin | Stochastic biomechanical modeling of the risk and risk factors for non-contact ACL injury in a stop-jump task[END_REF] estimated ACL injury risk for males to increase significantly with a PKAM of around 1.25 Nm/kg-bwt. Thus, handball players maintained knee joint loading at a relatively low level, unlikely to be associated with ACL injury. Otherwise, lower PKAM values for karatekas might find an explanation in their absence of expertise in the sidestepping task, yielding more challenging changes of direction execution for this population. In order to keep their joints safe while experiencing a new task, karatekas might have adopted some protective behavior. This would be in line with Sigward and Powers (2006b), who reported higher PKAM values for experts than nonexperts in the same sport. One potential adaptation could be the pelvis orientation towards the new direction already at IC. As the pelvis is central for core stability [START_REF] Houck | Comparison of frontal plane trunk kinematics and hip and knee moments during anticipated and unanticipated walking and side step cutting tasks[END_REF][START_REF] Weltin | Effect of gender on trunk and pelvis control during lateral movements with perturbed landing[END_REF], it could be suggested that better orienteering of the pelvis in the new sidestepping direction might stabilize the whole body and be a protective strategy for the knee. However, the sidestepping task's objective is to eliminate defender in handball. An obvious rotation toward the direction would give early information to the opponent, explaining why handball players presented a more neutral position. Greater lateral and vertical GRF were observed in handball players, whereas knee joint loading was kept at a safe level, suggesting a more performance-orientated execution by handball players, which is in line with the potential more aggressive behavior hypothesized by Sigward and Powers (2006b). Despite pelvic kinematics changes, no trunk motion difference was observed. Thus, contrary to our hypothesis, handball players did not better orientate their trunk in the targeted direction, despite differences in trunk stability among athletes with different disciplines and levels of practice reported in the literature [START_REF] Barbado | Sports-related testing protocols are required to reveal trunk stability adaptations in highlevel athletes[END_REF]Gloscheskie & Brown, 2017).

However, these studies revealed that athletic background doesn't necessarily impact every component of the core stability (e.g. core proprioception evaluated during seated postural control; Gloscheskie & Brown, 2017). Therefore, our results show that handball players and karatekas might have quite similar core stability and thus karatekas could transfer their ability to control their trunk to this new sidestepping task. In addition, the realization of unanticipated sidestepping might imply too high constraints, which subjects would handle almost in the same manner whatever the level of expertise or core stability.

Kinematic predictors of injury risk

In the frontal plane, this study has confirmed previous observations of an increase of trunk lateral lean at IC increasing PKAM [START_REF] Dempsey | Whole body kinematics and knee moments that occur during an overhead catch and landing task in sport[END_REF][START_REF] Jamison | Randomized controlled trial of the effects of a trunk stabilization program on trunk control and knee loading[END_REF][START_REF] Jones | Technique determinants of knee joint loads during cutting in female soccer players[END_REF][START_REF] Kristianslund | Sidestep cutting technique and knee abduction loading : Implications for ACL prevention exercises[END_REF][START_REF] Mornieux | Anticipatory postural adjustments during cutting manoeuvres in football and their consequences for knee injury risk[END_REF][START_REF] Staynor | By failing to prepare, you are preparing your anterior cruciate ligament to fail[END_REF]. Moreover, pelvis lateral lean during WA phase increased PKAM. Trunk lateral lean had the highest one SD percentage of change meaning that this variable explained more variation of PKAM than any others in the model. Thus, training interventions seeking to reduce ACL injury risk should primarily target trunk control in the frontal plane.

In the transverse plane, larger rotation towards the new direction was associated with larger PKAM. If trunk rotation towards the stance limb seems to be a predictor of PKAM during landings [START_REF] Critchley | The effects of mid-flight whole-body and trunk rotation on landing mechanics : Implications for anterior cruciate ligament injuries[END_REF][START_REF] Dempsey | Whole body kinematics and knee moments that occur during an overhead catch and landing task in sport[END_REF] and planned sidestepping maneuvers [START_REF] Frank | Trunk and hip biomechanics influence anterior cruciate loading mechanisms in physically active participants[END_REF], our results are in line with previous studies on unplanned sidestepping maneuvers [START_REF] Kristianslund | Sidestep cutting technique and knee abduction loading : Implications for ACL prevention exercises[END_REF][START_REF] Staynor | By failing to prepare, you are preparing your anterior cruciate ligament to fail[END_REF]. Therefore, trunk rotation during changes of direction should be taken into account with caution. While increased trunk rotation towards the new movement direction could be seen as a strategy to reorient the whole body and improve performance [START_REF] Marshall | Biomechanical factors associated with time to complete a change of direction cutting maneuver[END_REF], we found that this could also be associated to some extend with increased knee joint injury risk.

In the sagittal plane, only the pelvis anterior tilt at IC appeared to predict PKAM. From an anatomical point of view, pelvis anterior tilt flexes and rotates the hip internally, which changes gluteal muscles moment arms and weakens the hamstrings, therefore possibly impacting PKAM [START_REF] Alentorn-Geli | Prevention of non-contact anterior cruciate ligament injuries in soccer players. Part 1 : Mechanisms of injury and underlying risk factors[END_REF]. Overall, our results confirm the importance of an adapted core stability to limit ACL injury risk [START_REF] Larwa | Stiff landings, core stability, and dynamic knee valgus: a systematic review on documented anterior cruciate ligament ruptures in male and female athletes[END_REF].

Neuromuscular predictors of core kinematics and knee injury risk

As no specific trunk kinematics differed between our two populations, seeing no difference in directed co-contraction ratio was not surprising. However, trunk and pelvis neuromuscular control predicted core stability. DCCR in the frontal plane explained trunk lateral lean at IC, suggesting that higher contractions towards the stance leg during the flight-phase supported trunk lateral lean. Later on, during WA, lateral pelvis lean was predicted by both GMed and the right GMax. According to [START_REF] Neumann | Kinesiology of the hip : A focus on muscular actions[END_REF], this might underline the role of the left GMed (stance leg) as being associated with the pelvis contralateral drop, while right Gmed and GMax muscles probably abducted the contralateral lower limb. Surprisingly, DCCR in rotation did not predict trunk rotation, while DCCR in flexion did. DCCR in rotation was computed, based on only external oblique muscles, whose 3D actions might have blurred their impact on trunk rotation.

Pelvis anterior tilt at IC was predicted by anterior-posterior DCCR, meaning that core muscles acted to limit excessive pelvis anterior tilt. However, as DCCR explained only 13% of the pelvis tilt, the role of each muscle integrated in this ratio might not be fully understood. For instance, gluteus maximus probably also acted to control hip flexion as antagonists [START_REF] Alentorn-Geli | Prevention of non-contact anterior cruciate ligament injuries in soccer players. Part 1 : Mechanisms of injury and underlying risk factors[END_REF] and maybe had a larger role in force absorption than tilting the pelvis [START_REF] Powers | The influence of abnormal hip mechanics on knee injury : a biomechanical perspective[END_REF].

Limitations

The dynamic two-step approach used in the present study was similar to handball practice but differed from most the setups used in previous changes of direction studies with a classical run-up at a given approach speed. Nevertheless, PKAM values and trunk kinematics were in the range of the results found in the literature with a classic 4 m.s -1 running speed [START_REF] Weir | A reliable video-based ACL injury screening tool for female team sport athletes[END_REF][START_REF] Mornieux | Is there a sex difference in trunk neuromuscular control among recreational athletes during cutting maneuvers[END_REF] or a 3 steps approach [START_REF] Jamison | Randomized controlled trial of the effects of a trunk stabilization program on trunk control and knee loading[END_REF]. Also, the sampling frequency used for kinematics recording (i.e. 100 Hz) was lower than the 200 Hz typically reported in the literature. However, given the maximal frequency of the kinematic signals during changes of direction (up to 20 Hz), the present sampling frequency should be enough to gather the full information of the kinematic signals. Finally, trunk neuromuscular control could not be fully evaluated thanks to surface EMG, and future research might want to add an evaluation of deeper muscle activity.

Conclusion

Nonexperts in changes of direction reduced ground reaction forces and knee joint abduction moment by potentially rotating their pelvis more towards the new movement direction during sidestepping. However, the higher level of knee joint loading during the sidestepping task in experts might not be seen as a higher risk of ACL injury. Moreover, this study confirms that core stability (especially the trunk in the frontal plane) predicts knee joint loading and therefore possibly ACL injury risk. Finally, our results might guide neuromuscular prevention programs to focus on trunk lateral flexors' eccentric action in order to limit trunk lateral lean and on gluteal muscles' ability to stabilize pelvic lateral lean. This research could also influence interventions on kinematic changes of athletes' techniques while sidestepping. Further research could explore the balance between performance and injury during changes of direction tasks.

Figure 1 .

 1 Figure 1. Directed Co-Contraction Ratio (DCCR) during the preactivation phase (A) and the weight acceptance (WA) phase (B) for karatekas (grey bars) and handball players (black bars). Positive values indicate a co-contraction ratio towards trunk flexion, trunk medial lean and trunk right rotation.

Table 1

 1

	medial lean	1.1 (2.6)	1.7 (2.9)	1.4 (3.6)
	Pelvis axial rotation	4.9 (6.0)	10.8 (4.5)*	7.9 (6.6)
	Pelvis anterior tilt	10.0 (4.1)	11.0 (3.8)	10.5 (4.4)
	Knee abduction angle	2.4 (1.4)	0.21 (4.1)	1.3 (3.4)
	was 0.218 Nm/kg-bwt (CI = [0.008; 0.427]; ES = 0.84). Finally, the mean difference of peak
	medial GRF was 2.2 N/kg-bwt (CI = [0.58; 3.80]; ES = 1.03) and 3.6 N/kg-bwt (CI = [0.84;
	6.34]; ES = 1.01) for peak vertical GRF.

. Comparison of mean (SD) kinematic variables

During WA (Table

2

), PKAM was significantly greater for handball players than for karatekas (p = 0.04), as well as peak medial and vertical GRF (p = 0.01). The mean difference of PKAM

Table 2

 2

. Comparison of mean (SD) kinematic variables (°), peak knee abduction moment (Nm/kg-bwt) and GRF (N/kg-bwt) between handball players and karatekas during weight acceptance phase. *: p<0.05

Table 3

 3

		Predicted variable	Independent variable Estimate	CI	One SD change (%)
		PKAM (model no. 1)	Trunk medial lean Trunk axial rotation Pelvis anterior tilt	-33.9 13.1 20.4	[-51.9; -15.9] [1.6; 24.6] [2.3; 38.5]	35 22 21
	IC	Trunk medial lean	Medial DCCR	3.11	[1.38; 4.85]	35
		Trunk axial rotation	Flexion DCCR	-3.54	[-5.56; -1.51]	33
		Pelvis anterior tilt	Pelvis DCCR	2.93	[1.49; 4.38]	13
		PKAM (model no. 2)	Trunk axial rotation Peak posterior GRF Pelvis medial lean	11.8 40.9 -23.2	[3.0; 20.6] [10.0; 71.7] [-43.7; -2.7]	26 24 22
	WA	Trunk axial rotation	Flexion DCCR	-7.00	[-9.64; -4.36]	98
			Left Gluteus Medius	0.030 [0.011; 0.049]	164
		Pelvis medial lean	Right Gluteus Medius Right Gluteus Maximus	-0.026 [-0.044; -0.007] -0.081 [-0.156; -0.006]	149 95

. Kinematic and kinetic variables at initial contact (IC) and during weight acceptance (WA) predicting peak knee abduction moment (PKAM). EMG variables (DCCR or RMS) predicting core kinematics significant PKAM predictors. CI: Confidence interval; DCCR: directed co-contraction ratio. Gluteus Medius and Gluteus Maximus independent variables refer to their RMS activity.

Acknowledgments

The authors would like to thank all the participants of this study and Kristin L. Sainani for the statistical help.

Acknowledgments

The authors would like to all the participants of this study and Kristin L. Sainani for the statistical help.

Disclosure statement

No potential conflict of interest was reported by the authors.

All procedures were approved by the ethics committee Sud Mediterranee III (approval reference 2018.07.

Appendices


Yuchao Hua
email: yuchao.hua@univ-nantes.fr

Lingai Luo
email: lingai.luo@univ-nantes.fr

Online Learning Framework for Nonlinear Dynamic Thermal Modeling

Keywords: Thermal modeling, Online learning, Nonlinear, ensemble Kalman filter, Gaussian process 0

In order to obtain highly-efficient management schemes for thermal components and systems, the thermal dynamic modeling must be robust and accurate. Resistancecapacitance (R-C) thermal models involving unknown coefficients have been widely employed to realize the dynamic modeling following the gray-box methodology. In most of cases, R-C thermal models are amused to be linear, neglecting all the possible dependence of heat transfer coefficients on the node temperatures. Those linear models may become inapplicable, when nonlinearity effects are non-negligible. In addition, the training process of R-C models to determine the unknown coefficients is usually offline currently, which indicates the users have to re-train the models once the heat transfer features change during the operation process. Here, an ensemble Kalman filter (EnKF)-based framework integrated with Gaussian process (GP) is adopted to train nonlinear R-C thermal models to realize the self-adaptive dynamic modeling. The machine learning method, GP, guarantees the high-degree self-adaption of the framework, while the EnKF achieves the online training and the integration of nonlinear R-C thermal model and GP. Furthermore, the performance of our algorithm is evaluated on the synthetic datasets, which well proves its feasibility.

I. INTRODUCTION

Robust and accurate modeling of dynamic thermal behavior is the precondition to achieve highly-efficient thermal management from buildings [START_REF] Luo | Innovative systems for storage of thermal solar energy in buildings[END_REF] [START_REF] Hui | Evaluation of a seasonal storage system of solar energy for house heating using different absorption couples[END_REF] to electronics [START_REF] Rencz | Non-linearity issues in the dynamic compact model generation[END_REF][4] [START_REF] Hadeed | Development of a compact thermal model for electronic package A[END_REF] [START_REF] Shen | Spectral thermal spreading resistance of wide-bandgap semiconductors in ballistic-diffusive regime[END_REF]. Gray-box methodology is the most widely utilized approach for modeling various thermal problems [START_REF] Rao | Modeling of room temperature dynamics for efficient building energy management[END_REF][8] [START_REF] Quinten | Dynamic modelling of multidimensional thermal bridges in building envelopes: Review of existing methods, application and new mixed method[END_REF]. It employs lumped models that reflect some basic physical laws with some unknown coefficients to be determined by training according to experimentally-or numerically-derived datasets [START_REF] Deconinck | A maximum likelihood estimation of the thermal resistance of a cavity wall from on-site measurements[END_REF]. As the combination of white-and black-box methodologies, it has advantages of both sides: some physical information can be captured by the lumped models, while the corrections in terms of datasets guarantee the accuracy for any specific problem.

The flexibility and capability of lumped models can influence the performance of thermal modeling based on the gray-box approach [START_REF] Quinten | Dynamic modelling of multidimensional thermal bridges in building envelopes: Review of existing methods, application and new mixed method[END_REF]. Several types of lumped models have been developed, such as the structure factors method [START_REF] Kossecka | Multi-dimensional heat transfer through complex building envelope assemblies in hourly energy simulation programs[END_REF], the matrix of transfer functions method [START_REF] Quinten | Dynamic modelling of multidimensional thermal bridges in building envelopes: Review of existing methods, application and new mixed method[END_REF], and the resistancecapacitance (R-C) thermal model [START_REF] Alshatshati | Estimating building envelope thermal characteristics from single-point-in-time thermal images[END_REF] [START_REF] Górecki | Compact thermal modeling of power semiconductor devices with the influence of atmospheric pressure[END_REF], etc. Among them, the R-C thermal model has been widely adopted [START_REF] Alshatshati | Estimating building envelope thermal characteristics from single-point-in-time thermal images[END_REF][12] [START_REF] Torzewicz | Parametric compact thermal modeling of power LEDs[END_REF], where the dynamic thermal behavior is described using the model of resistance-capacitance circuits through the thermoselectric analogy. Most of R-C thermal models are linear, that is to say, the heat transfer coefficients (including conduction and convection) are assumed to be constant [START_REF] Wang | Development of RC model for thermal dynamic analysis of buildings through model structure simplification[END_REF]. This may make them inapplicable in the cases where nonlinearity has significant effects. Moreover, the training process of R-C models is frequently offline [START_REF] Brastein | Parameter estimation for grey-box models of building thermal behaviour[END_REF], which means the users have to re-train the models once the heat transfer properties change during the operation of a thermal device or system. Recently, although some researchers used the Kalman filters (KFs) to conduct the online estimations of heat transfer coefficients in thermal models [START_REF] Rao | Modeling of room temperature dynamics for efficient building energy management[END_REF] [START_REF] Elizabeth | Online building thermal parameter estimation via unscented Kalman filtering[END_REF], they are assumed to be constant or described by some simple functions of which the exact expressions should be given in advance. Therefore, their self-adaption ability could be limited.

In the present work, regarding these issues above, an online learning framework integrated with Gaussian process (GP) is utilized to train the nonlinear R-C thermal model, based on the ensemble Kalman filter (EnKF), with the numerical experiments to demonstrate the feasibility of our algorithm.

II. METHODOLOGY

A. Nonlinear R-C Thermal Model with GP-predicted Unknown Coefficients

In general, R-C thermal model is given by,

𝐶𝐶 𝑘𝑘 𝜕𝜕𝑇𝑇 𝑘𝑘 𝜕𝜕𝜕𝜕 = ∑ 𝑇𝑇 𝑖𝑖 -𝑇𝑇 𝑘𝑘 𝑅𝑅 𝑖𝑖𝑘𝑘 (𝑇𝑇 𝑖𝑖 ,𝑇𝑇 𝑘𝑘 ) 𝑖𝑖 + ∑ 𝑤𝑤 𝑗𝑗 𝑄𝑄 𝑗𝑗 𝑗𝑗 , (1)
where Tk(K) is temperature of node k, Ti(K) is temperature of node i, Ck (J/K) is heat capacity of node k, Rik (K/W) is thermal resistance between nodes i and k, Qj(W) is j th heat (or cold) input into node k, and wj is weighting factor. Some assumptions are adopted to build the R-C models: 1) the heat transfer between the components (nodes) is assumed to be one-dimensional; 2) the heating or cooling is uniform within each region corresponding to each node.

Different from the conventional R-C thermal models of which coefficients are constants, the present model is extended to nonlinear case that assumes them dependent on node temperatures. The GP predictors are employed to describe these coefficients. The GP is a representative method in the "kernel machines", which is characterized by the covariance functions constructed using kernels [START_REF] Seeger | [END_REF]. It is not needed to set a specific form of the unknown functions in advance, and thus the regression based on GP will not be restricted by the form of unknown functions, which can guarantee the high-degree selfadaption. Furthermore, the GP only relies on a few parameters to conduct predictions, in contrast to many of the state-of-the-art ML algorithms. Fewer parameters indicate the less training cost and data amount required, which could be beneficial to the efficiency of online learning for the dynamic thermal modeling problems where the amount of training data may be limited in practice.

The generic form of a GP model is given by,

𝑝𝑝 �� 𝐠𝐠 ob 𝐠𝐠 * �� ~ 𝑵𝑵 �𝟎𝟎, � 𝐾𝐾 cov (𝐗𝐗, 𝐗𝐗) + 𝜎𝜎 𝑔𝑔 2 𝑰𝑰 𝐾𝐾 cov (𝐗𝐗, 𝐗𝐗 * ) 𝐾𝐾 cov (𝐗𝐗 * , 𝐗𝐗) 𝐾𝐾 cov (𝐗𝐗 * , 𝐗𝐗 * ) �� (2)
in which 𝐠𝐠 ob is the assemble of observations for the target function 𝑓𝑓(𝐗𝐗) at the input variables X (that is called the grid points), 𝐠𝐠 * is the predictions at the newly-input variables 𝐗𝐗 * , 𝐾𝐾 cov is the covariance function, 𝜎𝜎 𝑔𝑔 2 is the variance corresponding to the Gaussian noise embedded in the observations, and N refers to a multivariate normal distribution with the mean vector equal to zero. The key prediction expression in the GP regression for 𝐠𝐠 * is derived from Eq. ( 2),

𝐠𝐠 * = 𝐾𝐾 cov (𝐗𝐗 * , 𝐗𝐗)�𝐾𝐾 cov (𝐗𝐗, 𝐗𝐗) + 𝜎𝜎 𝑔𝑔 2 𝑰𝑰� -1 𝐠𝐠 ob . ( 3
)
Once the covariance function 𝐾𝐾 cov is known, we can make the predictions, 𝐠𝐠 * , for the target function at the newly-input variables 𝐗𝐗 * , given the existing observations 𝐠𝐠 ob at X. The covariance function is built using the "kernel trick" [START_REF] Seeger | [END_REF]. Some commonly-used kernels include the squared-exponential, the Matérn, the dot product covariance functions, and their combinations. They generally hold several free parameters, i.e. hyperparameters, to be determined through a learning (training) process.

It is noted that the learning process for the conventional GP regression is offline: a group of data points will be collected in advance; afterwards the learning process is conducted to determine the covariance function. Moreover, with the increasing size of dataset (i.e. the number of the grid points), the computational time scales cubically, and the storage memory scales quadratically [START_REF] Kuzin | Ensemble Kalman filtering for online Gaussian process regression and learning[END_REF]. This may be a major obstacle for the utilization of GP regression for online learning, since the observed data points will be generated in a temporal sequence and thus the size of dataset will keep increasing with time. A popular means to resolve this issue is the sparse GP [START_REF] Bauer | Understanding probabilistic sparse Gaussian process approximations[END_REF]: the targeted function can be modelled "sufficiently well" by some "inducing points", 𝐠𝐠 in , at a small subgroup 𝐗𝐗 in of X; and thus the number of grid points is fixed to be the size of 𝐗𝐗 in and will not increase with the increasing number of observations; thus, the prediction expression becomes slightly different from that in the conventional GP,

𝐠𝐠 * = 𝐾𝐾 cov (𝐗𝐗 * , 𝐗𝐗 in )�𝐾𝐾 cov (𝐗𝐗 in , 𝐗𝐗 in ) + 𝜎𝜎 𝑔𝑔 2 𝑰𝑰� -1 𝐠𝐠 in , (4)
Besides the hyperparameters in the kernel function, 𝐠𝐠 in will also be taken as the parameters and optimized during the learning process.

More details about GP regression can be found in Refs. [17][20]. An online learning process will be employed to train the GP predictors for these heat transfer coefficients using the streaming data generated during operation process. In this sense, once heat transfer features change due to the aging effect or incidence, the online learning framework will adjust the model automatically.

B. Algorithm

The overall framework is designed based on a modified EnKF. The detailed explanations of its modules and mathematical clarification based on the theory of Bayesian estimation can be found in our recent paper [START_REF] Hua | An online learning framework for self-adaptive dynamic thermal modeling of building envelopes[END_REF]. Here, we just specify the basic procedure of the algorithm, as given in Fig. 1.

a) Initialization

IN1 is for initializing the GP [START_REF] Hua | An online learning framework for self-adaptive dynamic thermal modeling of building envelopes[END_REF]. Build the covariance function using the selected kernels; specify the hyperparameters within the kernel; choose the coordinates of inducing points.

IN2 is for initializing the EnKF. Set the size of ensemble N. Initialize the state ensemble and the parameter ensembles [START_REF] Hua | An online learning framework for self-adaptive dynamic thermal modeling of building envelopes[END_REF].

b) Prediction

P1 is for predicting the state ensemble from time t -1 to t. The prediction function in P1 is the R-C thermal model involving the unknown coefficients calculated by the GP model using the hyperparameters and the values of inducing points at time t -1.

P2 is for predicting the parameter ensembles from time t -1 to t in terms of the Liu-West filter [START_REF] Hua | An online learning framework for self-adaptive dynamic thermal modeling of building envelopes[END_REF] [START_REF] Yang | Combined parameter and state estimation in particle filtering[END_REF].

c) Update

UP1 is for updating the state ensemble at time t in terms of the observations. The measurement function depends on the problem investigated.

This work is financially supported by Région Pays de la Loire France within the NExT2Talents program TOP-OPTIM project (998UMR6607 EOTP NEXINTERTALENTHUA).

UP2 and UP3 are for updating the hyperparameters and the values of inducing points, respectively.

d) Resampling

The simple random resampling is performed every Nrs time steps to avoid the degeneracy problem [START_REF] Labbe | Kalman and Bayesian Filters in Python[END_REF].

e) Output

In our framework, the nonlinear R-C thermal model will be updated once a newly-observed data point is input to capture the timely-varying thermal features of a specific object, like the surface aging, which can avoid the data collection and offline retraining. As one period of observation ends, we will obtain a predictive model capable of predicting the thermal response of next period. When a new period of observation starts, the online learning will continue to further refine the R-C thermal model to guarantee its self-adaptive ability.

III. NUMERICAL EXPERIMENTS

Our algorithm will be tested on the synthetic datasets derived from a pre-defined benchmark model. Here, a 2R2C R-C model with temperature-dependent thermal resistances, R12(T1, T2) and R23(T2, T3), is chosen to generate the benchmark data points.

For convenience, we set some reference quantities, 𝑇𝑇 min (K), 𝑇𝑇 max (K) , and tperiod(hour), to convert the benchmark model above to be dimensionless,

𝜕𝜕𝑇𝑇 1 𝜕𝜕𝜏𝜏 = 𝑟𝑟 12 (𝑇𝑇 2 -𝑇𝑇 1 ) + 𝑞𝑞, ( 7
) 𝜕𝜕𝑇𝑇 2 𝜕𝜕𝜏𝜏 = 𝑟𝑟 23 (𝑇𝑇 3 -𝑇𝑇 2 ) + 𝑟𝑟 12 𝑟𝑟 𝐶𝐶12 (𝑇𝑇 1 -𝑇𝑇 2 ) (8)
with

𝑇𝑇 = 𝑇𝑇 -𝑇𝑇 min 𝑇𝑇 max -𝑇𝑇 min , 𝜏𝜏 = 𝑡𝑡 𝑡𝑡 period , 𝑟𝑟 12 = 𝑡𝑡 period 𝑅𝑅 12 𝐶𝐶 1 , 𝑟𝑟 23 = 𝑡𝑡 period 𝑅𝑅 23 𝐶𝐶 2 , 𝑟𝑟 𝐶𝐶12 = 𝐶𝐶 1 𝐶𝐶 2 , 𝑞𝑞 = 𝑄𝑄 1 𝑡𝑡 period 𝐶𝐶 1 .
tperiod is the duration of one training period that is set as 24 hours, i.e. 1day. Moreover, the relations between the dimensionless thermal resistances and the dimensionless temperature are assumed to be nonlinear,

𝑟𝑟 12 = 𝑐𝑐 12 �1 + 𝑎𝑎 1 𝑇𝑇 1 𝑎𝑎 2 + 𝑎𝑎 3 𝑇𝑇 2 𝑎𝑎 4 �, (9)
𝑟𝑟 23 = 𝑐𝑐 23 �1 + 𝑎𝑎 5 𝑇𝑇 2 𝑎𝑎 6 + 𝑎𝑎 7 𝑇𝑇 3 𝑎𝑎 8 �, (10)
with the constants{𝑐𝑐 Days #1-4: {2.0, 2.2, 1.1, 1.0, 0.8, 4.0, 2.0, 0.8, 1.0, 1.5} Days #5-8: {2.5, 1.8, 1.5, 0.8, 1.0, 2.5, 3.0, 2.0, 0.8, 2.0}

Furthermore, the prediction function in our framework is given,

𝑇𝑇 1(𝜏𝜏) = 𝑇𝑇 1(𝜏𝜏-1) + 𝑟𝑟 12_GP �𝑇𝑇 2(𝜏𝜏-1) -𝑇𝑇 1(𝜏𝜏-1) �∆𝜏𝜏 + 𝑞𝑞∆𝜏𝜏 (11) 𝑇𝑇 2(𝜏𝜏) = 𝑇𝑇 2(𝜏𝜏-1) + 𝑟𝑟 23_GP �𝑇𝑇 3(𝜏𝜏-1) -𝑇𝑇 2(𝜏𝜏-1) �∆𝜏𝜏 + 𝑟𝑟 12_GP 𝑟𝑟 𝐶𝐶12 �𝑇𝑇 1(𝜏𝜏-1) -𝑇𝑇 2(𝜏𝜏-1) �∆𝜏𝜏 (12)
where ∆𝜏𝜏 = 0.01 is the time interval, and 𝑟𝑟 12_GP & 𝑟𝑟 23_GP are predictions by the GP predictor that will be gradually refined during the learning process.

Here, the total observation period is set as 8 days. Each day will sever as a sub-period. The temperature evolutions of a subperiod will be estimated using the R-C model within the GP predictor trained during the sub-periods before. As for Day #1, the initial GP model will be utilized. The temperature of node 3, T3, and the dimensionless heat input will be taken as input variables, of which evolutions are illustrated in Figs. 3(a) and (b), respectively. Moreover, in order to quantitatively evaluate the prediction accuracy, the normalized mean squared error (NMSE) is calculated,

NMSE = 1 𝑀𝑀 ∑ � �𝑦𝑦 𝑚𝑚 test -𝑓𝑓�𝜏𝜏 𝑚𝑚 test �� 2 �𝑦𝑦 𝑚𝑚 test � 𝑀𝑀 𝑚𝑚=1 , (13)
where 𝑦𝑦 𝑚𝑚 test is the benchmark data point (without perturbation) at the test time 𝜏𝜏 𝑚𝑚 test , and 𝑓𝑓(𝜏𝜏 𝑚𝑚 test )is the prediction at the same time. According to Figs.4 &5, the nonlinear R-C model with the initial GP model gives a poor prediction of Day#1's temperature evolutions: the relative deviations from the benchmark data, i.e. NMSEs reach more than 0.2, which is much larger than the artificial relative measurement error. With the online learning process beginning and continuing, the model is self-adapted and able to give better predictions: on Day#4, the NMSEs become less than 0.05.

In order to test the self-adaptive ability of our framework, the dependence of dimensionless resistances on temperature is changed after Day#4 to simulate an incidental variation of heat transfer characteristics. In this case, the NMSEs increase considerably at Day#5, since the model is trained on the data before Day#5 that does not characterize the influence of that incidental change. Nevertheless, that errors will be corrected by the online learning process without artificial intervention. After Day #5, the NMSEs are largely reduced. Therefore, our algorithm may also serve as a means of fault-detecting for the thermal devices by evaluating the NMSE values.

IV. CONCLUSIONS

The gray-box approach with R-C models has been extensively utilized to handle dynamic thermal modeling problems from buildings to electronic devices. Nevertheless, currently most of R-C thermal models are assumed to be linear, which neglect the possible dependence of coefficients on the node temperature. Meanwhile, the training process of R-C models is frequently offline and lacks of the self-adaption ability.

Here, we propose that a modified EnKF-based framework can be employed to train nonlinear R-C thermal models involving GP predictors online. The R-C thermal model defines the physical constraints during the online learning process; the node-temperature-dependent coefficients are characterized by the GP, without restricting to some specific expressions, which guarantees the high-degree self-adaption of the framework. The algorithm is tested on the synthetic datasets obtained from the 2R2C benchmark model. The numerical experiments well demonstrate that the GP-embedded nonlinear R-C thermal model can gradually improve itself during the learning process and hold the strong self-adaptive ability that can handle the incidental variations of heat transfer characteristics.

Our framework enables training the nonlinear R-C thermal models capable of self-improvement and self-adaption, which will be much helpful for the development of the optimal control of various thermal systems and components.

Fig. 1 .

 1 Fig. 1. Flowchart of the algorithm [20].

Fig. 2 .

 2 Fig. 2. Schematic for the 2R2C R-C model.

Fig. 2 2 𝑅𝑅 12 (

 2212 Fig.2shows the schematic of 2R2C model. It is usually employed for modeling the energy consumption and temperature evolutions of buildings[START_REF] Wang | Development of RC model for thermal dynamic analysis of buildings through model structure simplification[END_REF], and it is given by𝐶𝐶 1 𝜕𝜕𝑇𝑇 1 𝜕𝜕𝜕𝜕 = 𝑇𝑇 2 -𝑇𝑇 1 𝑅𝑅 12 (𝑇𝑇 1 ,𝑇𝑇 2 ) + 𝑄𝑄 1 ,(5)

Fig. 3 .

 3 Fig. 3. The input variables: (a) Node temperature T3; (b) Dimensionless heat input q.

Fig. 4 .Fig. 5 .

 45 Fig. 4. Temperature evolutions from Day #1 to Day #8.

  12 , 𝑎𝑎 1 , 𝑎𝑎 2 , 𝑎𝑎 3 , 𝑎𝑎 4 , 𝑐𝑐 23 , 𝑎𝑎 5 , 𝑎𝑎 6 , 𝑎𝑎 7 , 𝑎𝑎 8 }. Additionally, an artificial relative error εre is added to the benchmark dimensionless temperature solutions to generate the observations. The setting values of the benchmark model are given in TABLE I.

TABLE I

 I 𝑎𝑎 1 , 𝑎𝑎 2 , 𝑎𝑎 3 , 𝑎𝑎 4 , 𝑐𝑐 23 , 𝑎𝑎 5 , 𝑎𝑎 6 , 𝑎𝑎 7 , 𝑎𝑎 8 }

	.	SETTING VALUES OF BENCHMARK MODEL
	Symbol	Setting value
	𝑇𝑇 min (K)	0
	𝑇𝑇 max (K)	30
	tperiod (hour)	24
	𝑟𝑟 𝐶𝐶12	2
	εre	0.02
	{𝑐𝑐 12 ,

ACKNOWLEDGMENT

This work is financially supported by Région Pays de la Loire France within the NExT2Talents program TOP-OPTIM project (998UMR6607 EOTP NEXINTERTALENTHUA).

AUTHORS' BACKGROUND