This study was conducted as an exploratory observational cohort study. The analysis of ACC data from the motor active part of the study population is a reanalysis with the method developed in the present study. Data has previously been analyzed and reported using the DAR [15]. To make reporting as transparent as possible, the manuscript has been developed following the STROBE guidelines.

Patients were recruited from a single ward at a neurorehabilitation hospital. The ward specializes in the rehabilitation of patients with sABI as soon as they are weaned from mechanical ventilation and sedation in acute care. Patients were screened when they arrived and deputy consent was obtained before 21 days of data acquisition. Inclusion criteria were $\ge$18 years old, classified as $\le$7 on the Rancho Los Amigos Scale (RLAS) of consciousness and cognition [17], relatives able to give deputy consent, sABI, expected stay $>$3 weeks, no paroxysmal sympathetic hyperactivity, no spinal lesions, no terminal illness and no diagnosis of polyneuropathy.

Motor inactive patients were defined as patients with RLAS 3 or locked-in syndrome described in their hospital record as recommended by Zollman et al. [16]. Furthermore, patients with severe tetraparesis, defined as less motor activity than two times the standard deviation of the total activity level of the beforementioned inactive patients, were also classified as motor inactive. The remaining patients were classified as motor active.

Data from patients included in the study were recorded with two data recording modalities. To measure motor activity, patients wore an Actigraph wGTX3-BT device as soon as consent could be obtained after admission. The ACC device contains a 3-axis micro electro mechanical system (MEMS) accelerometer with a dynamic range of $\pm$ 8 G and uses a 12 bit analog to digital converter sampling in a 0.002 G resolution. The device was set to record at a 30 Hz sampling rate with idle sleep mode enabled to conserve battery. The accelerometer was fixed to the subject’s least affected arm with a watch strap. Data was aggregated into one-minute epochs of vector magnitude counts using Actilife vs. 6.13.3 software (Actigraph Inc., Pensacola, FL, USA). By using Actigraphs proprietary unit counts/minute, the Actigraph algorithm is applied to the data, and any non-human activity is filtered out. By keeping aggregation levels at one-minute epochs, the time series closely represent the continuous activity levels of the patients.

HR was collected using the Mindray TM80 telemetry three lead Electrocardiogram (ECG) device, with one lead placed below each clavicle and one on the left lateral distal side of the ribcage. This was done according to manufacturer instructions. Data was transmitted to the Central Monitoring System (CMS) via local Wi-Fi where data was stored for up to 11 days and then exported in “csv” file format for processing. HR trend data was exported as an aggregated measurement of HR every five seconds extrapolated from the raw ECG signal by the CMS. Data were merged, cleaned for duplicates, sorted, and averaged over each minute using the software package Matlab r2019b (MathWorks Inc., Natick, MA, USA) [18] before analysis.

To ensure synchronized data from ECG and ACC, we used the hospital’s local Network Time Protocol (NTP) server to update the time on both the Central Monitoring System (CMS) and the laptop installed with Actilife vs. 6.13.3 software (Actigraph Inc., Pensacola, FL, USA) also used for initialization of the wGTX3-BT devices. The CMS was specified to synchronize with the NTP server approximately three times per minute and to correct the time if discrepancies between the CMS and the NTP were more than three seconds. The computer installed with the Actigraph software was synchronized immediately before the initialization of individual accelerometers.

Due to rehabilitation activities, the ECG recording was occasionally interrupted, for example during bathing or examinations outside the range of the hospital WIFI. Additionally, it was necessary to stop the recording due to hygienic procedures and performance-related replacement of the electrode pads and battery pack each day. This implied that the ECG recordings were not completely seamless and a threshold for the allowed amount of missing data was needed. To allow for minor technical issues such as loose or disconnected electrodes, any missing data with a duration less than or equal to 30 minutes were allowed. If more than 85% of the data of the remaining 24 hr cycle were available, the cycle was included in the analysis.

Accelerometers were preset to record for 21 days. When recording had started, it was not possible to turn off the devices during data acquisition. In cases where the accelerometer needed to be removed, e.g., during bathing, the device was kept together with the patient’s belongings. These events were not documented, and therefore, could not be separated from genuine patient activity in the data. Bathing frequency and duration was similar for all patients, and usually took place in the mornings and evenings but not daily. In cases where patients had to leave the hospital for shorter or longer durations, non-use was documented. On rare occasions, the accelerometer was placed on the wrong arm after having been removed. When this was noticed, the duration of the mistake was recorded as starting from the last bathing session or documented nonuse and deemed as nonuse. These periods were validated by visual inspection of individual actograms. Staff was instructed to check placement each day. If any part of a 24 hr cycle was missing due to nonuse, the entire day was excluded from the analysis.

The model for the circadian rhythm is exemplified in Fig. 1. The code for the model was developed and executed in R Version 3.6.2 (R Core Team Inc., Vienna, Austria) [19]. For this model, a circadian rhythm was conceptualized as a persistent dip of the signal for at least five hours during a 24 hr cycle from 12:00 h to 12:00 h. The circadian rhythm was modeled by a partially linear model. For each fixed starting and ending point (t${}_1$ and t${}_2$), the three unknown levels, annotated pre-dip level, dip-level, and post-dip level of the measure in question were fitted by the least-squares method to the log(1 + x) transformed ACC or log-transformed HR data. This allowed only a drop in the signal at the start of the dip and an increase at the end of the dip. These were annotated the dip-drop and the dip-rise. From the timing of the dip-drop and the dip-rise, the dip-duration could be estimated. The least-squares criterion was then minimized by iterating over all eligible values of t${}_1$ and t${}_2$ in steps of one minute. Log and log(1 + x) transformation was used to avoid excessive influence of extreme values of HR and ACC data, respectively, and at the same time allowing zero values of ACC. A back transformation of model estimates was performed before further statistical analysis. The lower limit of the dip duration of five hours was determined based on what could represent a meaningful night for this patient group. It was assumed that a conceptualized nighttime would only be reasonable if its duration was longer than that of the daytime resting periods. At the recruitment ward, 80% of the daytime resting periods have earlier been documented being between 30 minutes to two and a half hours. On occasions, adjacent resting periods are interrupted by short activities, creating possible false nights during the day from a model perspective. Someren et al. [20] used the 10 most active hours versus the five least active hours as a means of documenting rest-activity cycles with ACC in Alzheimer’s patients. The purpose of our model was similar although more nuanced since no time periods were discarded. All data recorded from the 24 hr cycle were used for analysis if sufficient data were available for relevant analysis. Based on the pre-dip level, dip level and post-dip level, the phase shift expressing the difference between the expected midpoint of dip and the actual midpoint of the dip was calculated. The expected midpoint of dip was set to 02:00 h based on the rhythm of activity at the rehabilitation hospital and a normal night time from 22:00 h to 06:00 h. The phase shift was calculated from the actual midpoint of the dip and the expected midpoint of the dip. The phase shift expresses whether the timing of the dip is adequate relative to the hospital’s rhythm. Lastly, the absolute level difference between pre-dip and dip-level, and post-dip and dip-level was calculated. These shifts in level express the intensity of the difference between the dip and the rest of the 24 hr cycle. The mathematical expression for the model is presented in Eqn. 1:

(1) $\log (\mathrm{measure}(\mathrm{t}))=$

Fig. 1.

Visual representation of the circadian dip model used in the present study.

The mathematical model will attempt to estimate the dip-duration, dip-drop, and dip-rise regardless of their magnitude. Therefore, a set of assumptions and restrictions was defined for respectively HR data and ACC data to classify the presence of circadian rhythm. The restrictions were set depending on the completeness and quality of the recorded data. The absence of circadian rhythm was defined for 24 hr cycles in which any of the following three conditions were satisfied:

- The estimated dip-duration did not exceed the imposed minimal limit of five hours.

- No least-squares regression could be fitted.

- A least-squares regression could be fitted, but the estimates did not fulfill the set of thresholds defined below.

The thresholds were intended to ensure that the estimated model is an expression of a period in time resembling rest as compared to the remaining hours of the 24 hr cycle. The ACC signal is quite dichotomous, i.e., active or passive, Bigué et al. [21] advised using low thresholds when estimating sleep-wake cycles, and Duclos et al. [11, 14] used thresholds of, respectively, 10 counts and 20 counts per epoch (minute) in patients with TBI admitted to intensive care when evaluating rest-activity cycles. Based on the above, the thresholds for this model were based on the pre-dip, during-dip, and post-dip levels of ACC. Both the estimated pre-dip and post-dip were required to be greater than 10 counts and the during dip level to be less than or equal to 10 counts to classify as a circadian rhythm. The HR measure is less dichotomic in its nature as compared to ACC and fluctuates between resting HR and the individual maximum HR. For this reason, a lower limit was set at five bpm in both dip-drop and dip-rise, ensuring an absolute minimum difference in HR level between the pre-and post-dip levels and the dip level. The resting HR was not used since it is likely to be abnormal for these patients and fluctuant across days.

Evaluation of distributional characteristics was done by visual inspection of histograms. For both HR and ACC data, estimates of mean or median dip duration and phase shift were calculated for all subjects taking into account the non-independent sampling of 24 hr cycles by using the sandwich covariance estimator [22]. All statistics reported in the paper have been calculated using Stata IC 16.1/MP with an $\alpha$-level of 5% [23]. The sample size was implied by the chosen trial period and the consecutive patient flow during data collection.

In total, 609 cycles were planned to be recorded. Data recordings were available in 80% of all cycles for HR and 76% for ACC. In total, 420 cycles were planned to be recorded on the 20 motor active patients. Because of various factors, such as removal of equipment before transfer to an external unit or technical or human errors, etc. as described under the sections describing the preprocessing, forty-nine 24 hr cycles had insufficient data for analysis of both HR and ACC recordings. In 27 cycles, there were insufficient data for analysis in ACC recordings, but sufficient in HR recordings. In 62 cycles, ACC recording was sufficient, but HR recording was not. In total, this left 344 cycles for analysis on ACC and 309 cycles for analysis on HR data in the motor active patients.

One hundred and eighty-nine 24 hr cycles were planned to be recorded in the nine motor inactive patients. Again, several cycles were excluded based on the preprocessing procedures. In 33 cycles, there were insufficient data in both HR and ACC recordings. In two cycles insufficient data were available on ACC but sufficient data was available in HR. In cases where HR data were unavailable, ACC data were available in 17 cycles. In total, this left 154 cycles available in ACC and 139 on HR data in the motor inactive patients. Results of the preprocessing procedures can also be inspected in Table 2.

The examples in Fig. 3 show typical 24 hr cycles of patients with sABI, and both motor active and motor inactive patients are represented. Furthermore, it visually depicts how ACC signals of patients with sABI struggle to fit the cosinor model, and how our model fits and identifies a relevant shift between day and night and estimates dip duration and phase shift using an ACC signal, but also using an HR signal in motor inactive patients.

Fig. 3.

Figures exemplifying the circadian dip model and the cosinor model on heart rate (HR) and accelerometry (ACC) data. (A) is an example of the circadian dip model performance. It identifies what is perceived as the circadian dip from the data. (C) is the equivalent cosinor model fitted to the same data. (B) and (D) depict a typical patient case where motor activity is not sufficient to identify a rhythm. In (B) only limited activity is broadly spread over most of the 24 hr cycle and there is only a small indication of a difference between night and day. Moreover, a relevant phase shift or dip duration is not estimated adequately since model thresholds are not reached in ACC data. In Fig. 3D, a distinct difference between night and day is relatively clear in the HR signal, estimating a dip duration of 9:18 hours starting at 21:46 h and ending at 07:03 h and a phase shift of 24.5 minutes delay as compared to the hospital fixed rhythm.

The presence or absence of circadian rhythm and 24 hr characteristics were estimated with the circadian dip model. The median proportion (IQR) (range) of cycles defined as reflecting a circadian rhythm in motor active patients of respectively ACC and HR was 57% (19–76) (0–95) for ACC and 50% (21–74) (0–90) for HR. Individual time series plots of model estimations per motor active patient per cycle can be found in the supplementary material.

The overall range of dip-durations in cycles with rhythm estimated by the model varied from 302 minutes to 1349 minutes in ACC data and from 303 minutes to 1368 minutes in HR. Fig. 4 (middle and left column) presents examples of the estimated dip-duration across admission days for HR and ACC in two motor active patients. Dip-durations reaching the lower time limit of an accepted dip-duration of 5 hrs are marked by plusses in Fig. 4. The estimated dip-duration was evaluated to be non-normal distributed in both HR and ACC data and was log-transformed before analysis, and back-transformed as median (95% CI) upon reporting (Table 3).

Fig. 4.

Examples of dip duration across days post-admission. The upper row represents examples of a circadian rhythm being estimated as overall present and the lower row presents cases in which the circadian rhythm is estimated as overall not present. Motor inactive (MI), Motor active (MA). ID24 is a 47-year-old patient who suffered a traumatic brain injury 64 days before the first day of recording. The subject had a high consciousness and cognitive level at baseline (Rancho Los Amigos Scale (RLAS) = 6 and Coma Recovery Scale-Revised (CRS-R) = 23). ID29 is a 74-year-old patient who had suffered a cerebral infection 274 days before the start of recording. This subject was in a confusional state and had a very low functional level at baseline (RLAS = 5 and Functional independence measure (FIM) = 19). ID9 is a 68-year-old patient who suffered a hemorrhagic stroke 157 days before the start of recording. This patient was in a state of prolonged disorder of consciousness at baseline (RLAS = 2 and CRS-R = 5). ID28 is a 60-year-old patient who suffered an ischemic stroke 66 days before the start of recording. The patient ID28 was in the process of emerging from a minimal consciousness state at baseline and had a very low functional level (RLAS = 3, CRS-R = 19, and FIM = 18). All patients displayed here were medicated with melatonin administered between 8 pm and 10 pm every evening as part of their medical care. ID28 had signs of infection (temperature $\ge$38 ${}^{\circ }$C, the temperature was measured at least once per day) during the first two recording days, and ID29 had signs of infection in 19/21 recording days, scattered from the second day to the last day.

The motor active subject ID2 presented in Fig. 5 represents the borderline between the relevance of the ACC and the HR analysis. Only 2 of 21 24 hr cycles with complete ACC data can be defined as circadian rhythm. On the opposite, the HR analysis estimates 17 of 18 recorded 24 hr cycles as having a circadian rhythm. In the ACC analysis, the dominant reason for the analysis to decline the presence of circadian rhythm is the fact that the model thresholds of a minimum shift and dip duration level are rejected.

Fig. 5.

Visualisation of the estimated circadian rhythm, dip duration, and phase shift of the patient with ID2, as an example of the lower limit of the relevance of the accelerometer (ACC) analysis, and the potential relevance of the heart rate (HR) analysis. ID2 is a patient who suffered encephalopathic brain injury more specifically an anoxic brain injury 38 days before the data recording started. At baseline, the patient had a low but not undetectable functional level while in a confused state (Functional independence measure (FIM) = 26 and Rancho Los Amigos Scale (RLAS) = 5). MA, Motor active.

Phase shift, the deviance of the actual midpoint from the expected midpoint, of the estimated model in ACC data ranged from –664 minutes to 334 minutes. In HR, similar results were found with a range from –621 minutes to 428 minutes. The estimated phase shifts for both HR and ACC were found to be normally distributed, hence the mean (95% CI) is reported (Table 2). In the above analyses, only estimates from cycles with rhythm detected are included. Fig. 6 (middle and left column) presents examples of the estimated phase shift across admission days for respectively ACC and HR in motor active patients.

Fig. 6.

Examples of phase shift across days post-admission. The top three tiles represent examples of circadian rhythm being estimated as present and the bottom three tiles indicate cases in which circadian rhythm is not estimated as present. Motor inactive (MI), Motor active (MA). ID24 is a 47-year-old patient who suffered a traumatic brain injury 64 days before the first day of recording. The subject had a high consciousness and cognitive level at baseline (Rancho Los Amigos Scale (RLAS) = 6 and Coma Recovery Scale-Revised (CRS-R) = 23). ID29 is a 74-year-old patient who had suffered a cerebral infection 274 days before the start of recording. This subject was in a confusional state and had a very low functional level at baseline (RLAS = 5 and Functional independence measure (FIM) = 19). ID9 is a 68-year-old patient who suffered a hemorrhagic stroke 157 days before the start of recording. This patient was in a state of prolonged disorder of consciousness at baseline (RLAS = 2 and CRS-R = 5). ID28 is a 60-year-old patient who suffered an ischemic stroke 66 days before the start of recording. The patient ID28 was in the process of emerging from a minimal consciousness state at baseline and had a very low functional level (RLAS = 3, CRS-R = 19, and FIM = 18). All patients displayed here were medicated with melatonin administered between 8:00 h and 22:00 h every evening. ID28 had signs of infection during the first two recording days and ID29 had signs of infection on 19/21 recording days, scattered from the second day to the last day ($\ge$38 ${}^{\circ }$C body temperature measured at least once per day).

For the motor inactive patients, the heterogeneity between days and between individuals is similar to the motor active patients. The proportions of 24 hr cycles with a circadian rhythm as measured by HR for motor inactive patients ranges from 0% to 86% with a median (IQR) of 29% (24–52). Individual time series plots of model estimations per motor inactive patient per cycle can be found in the supplementary material.

Meaningful estimation of dip duration and phase shift was possible in motor inactive patients using heart rate measurements. Dip duration ranged from 303 to 855 minutes. Phase shift was estimated to have a range from –678 to 293 minutes. The estimated dip duration was evaluated to be log-normally distributed and therefore log-transformed before analysis and reported as median (95% CI) while phase shift was evaluated to be normally distributed and therefore reported as mean (95% CI) (Table 3). The two outer right tiles in Figs. 4 & 6 represent examples of motor inactive patients and respectively dip duration and phase shift across days.

In this study, we presented a pragmatically developed model for the identification of circadian rhythm and estimation of dip duration and phase shift for individual 24 hr cycles. The conception of the model is primarily based on ACC recordings but it may also be useful in HR recordings. We have applied the model to recordings from patients with sABI. The rationale behind this approach was that established models are not suitable for the immensely fluctuating and variable signals of HR and ACC in sABI patients. We have attempted to fit the classic cosinor model as described by Cornelissen [10], and also the alternative of the extended cosinor as presented by Marler [27], which could not be fitted to our data (estimation algorithms failed to converge). The extended cosinor is a mathematical attempt to comply with the more dichotomous nature of the ACC signal, often resembling a square increase in activity during daytime and a square decrease in activity during nighttime. We suggest the presented model as a more direct and adequate alternative for this kind of data compared to existing models. Furthermore, the model thresholds can be transferred into clinically relevant measurements which could be adjusted individually from patient to patient. The model can be used as a simple and coarse way of detecting circadian rhythm on a day-to-day basis when clear signs of sleep are not yet present or consciousness is still fluctuating. If a circadian rhythm was detected, we were able to estimate nighttime, and the timing of the night relative to the expected midpoint of night (The phase shift). In a clinical implementation, it could potentially supply clinicians with valuable patient monitoring and information about circadian rhythm in conjunction with clinical knowledge, which in turn could be translated into relevant and timely intervention focusing on improvement of circadian disruptions.

We found the circadian rhythm of patients with sABI to be fluctuating from day-to-day. Furthermore, the rhythm does not necessarily develop exclusively in a positive manner across time. The circadian rhythm was even nonexistent per our definition on some days. Our analysis also indicates high between-subject variations possibly linked to the heterogeneity of the study population, despite all patients being diagnosed with sABI. These variations limit the meaningfulness of aggregating data for all individuals and estimating a common slope of change in dip duration or phase shift. The heterogeneity is underpinned by the fact that the analysis estimates circadian rhythm by definition for respectively ACC and HR recording to be as little as 0% for both ACC and HR and as many as respectively 95% and 90% of 24 hr cycles per subject regardless of motor activity level. It is unlikely that the observed fluctuations in the estimated circadian rhythm through ACC and HR are caused by changes in the endogenous circadian rhythm generated within the hypothalamus in the suprachiasmatic nucleus. It is more likely that these fluctuations are the result of environmental modifiers and the overall neurological and physical impact of brain injury. This is backed by the findings in patients with traumatic brain injury from Duclos et al. 2020 [2]. In this study, the authors found the endogenous circadian rhythm as measured with urinary melatonin excretion to be intact in patients with disrupted sleep-wake cycles measured with ACC [2].

Our results suggest that HR in some cases can reflect circadian rhythm, even in severe brain-injured people who are almost completely motor inactive. In one specific subject, the relevance of the two modalities seemed especially clear. In ID2, a per definition motor active subject HR indicated circadian rhythm while ACC did not. HR could be an alternative to ACC as a long-term, non-invasive, and easy data source for the detection of circadian rhythm in immobilized patients.

We have presented data on circadian rhythm along with information about the type of injury, presence of infection, and use of melatonin in individual patients. We propose that the combination of these kinds of data could be a valuable tool for clinicians in the context of neurorehabilitation of patients with sABI when considering the choice of interventions and medication. To reach a clinically applicable solution, the analysis of circadian rhythm needs to be introduced as part of the monitoring systems in clinical practice. This goal requires collaboration between clinicians, circadian rhythm scientists, and commercial companies. For such clinical applications to be implemented in a clinical setting, it is imperative that solutions are simple to use, and do not require expert knowledge of circadian analysis. Therefore the circadian dip model is a good option since it is fairly simple to understand, and relies on clinically interpretable model limits and not the goodness of fit or test of statistical significance.

We initially sought to estimate circadian rhythms in patients using both ACC and HR, and only in hindsight, we became aware that none of the existing methods would apply to the data collected from the study population. For this reason, no validation of the models has been undertaken yet. A limitation of this study is the lack of validation of the model implemented. However, validation could be hard to accomplish since the model we have developed does not attempt to infer anything about sleep quantity or quality, which most established tools do. It estimates whether a minimum difference between two distinct periods during each 24-hour cycle is present, how long the dip is, and when in time it is situated. Validation may be possible through subjective staff reports of night duration and temporal placement of the night. However, to our knowledge, there is no golden standard that is suitable for the type of data in the present study. The thresholds of ACC and HR in the definitions of a circadian rhythm applied in our model are based on similar thresholds in similar patients but are subject to further scrutiny. We did not have any alternative objective or subjective data on the presence or absence of circadian rhythm and its characteristics, and so it was impossible to test the validity of these thresholds. Furthermore, the thresholds might even be individual from patient to patient. Further research is needed to elucidate this. To counteract information bias, we tried to reduce the possibility of estimating circadian rhythm when none was present by setting the thresholds of the ACC and HR model rather conservatively. A further limitation in the reporting of the results of the present study is the possible selection of good data. Data recording was paused if patients were transferred for external examinations or treatment. This means that conclusions inferred can only be generalized to patients admitted for in-hospital subacute neurorehabilitation. It is likely, that a temporary transfer to an external unit could briefly influence circadian rhythm, but the directionality of the influence is yet unknown and possibly dependent on the patient characteristics such as the location and severity of their brain injury and also the reason for the transfer. It presents as a new research question in need of further investigation whether transfers affect the circadian rhythm. Lastly, the sample size of the present study is relatively small, reducing the overall generalizability. Nevertheless, the sample size is quite acceptable for a study investigating a group of severely brain injured patients with complex neurological and medical problems.