 Research article
 Open Access
 Published:
Exploring diurnal variation using piecewise linear splines: an example using blood pressure
Emerging Themes in Epidemiology volume 14, Article number: 1 (2017)
Abstract
Background
There are many examples of physiological processes that follow a circadian cycle and researchers are interested in alternative methods to illustrate and quantify this diurnal variation. Circadian blood pressure (BP) deserves additional attention given uncertainty relating to the prognostic significance of BP variability in relation to cardiovascular disease. However, the majority of studies exploring variability in ambulatory blood pressure monitoring (ABPM) collapse the data into single readings ignoring the temporal nature of the data. Advanced statistical techniques are required to explore complete variation over 24 h.
Methods
We use piecewise linear splines in a mixedeffects model with a constraint to ensure periodicity as a novel application for modelling daily blood pressure. Data from the Mitchelstown Study, a crosssectional study of Irish adults aged 47–73 years (n = 2047) was utilized. A subsample (1207) underwent 24h ABPM. We compared patterns between those with and without evidence of subclinical target organ damage (microalbuminuria).
Results
We were able to quantify the steepest rise and fall in SBP, which occurred just after waking (2.23 mmHg/30 min) and immediately after falling asleep (−1.93 mmHg/30 min) respectively. The variation about an individual’s trajectory over 24 h was 12.3 mmHg (standard deviation). On average those with microalbuminuria were found to have significantly higher SBP (7.6 mmHg, 95% CI 5.0–10.1) after adjustment for age, sex and BMI. Including an interaction term between each linear spline and microalbuminuria did not improve model fit.
Conclusion
We have introduced a practical method for the analysis of ABPM where we can determine the rate of increase or decrease for different periods of the day. This may be particularly useful in examining chronotherapy effects of antihypertensive medication. It offers new measures of shortterm BP variability as we can quantify the variation about an individual’s trajectory but also allows examination of the variation in slopes between individuals (randomeffects).
Background
There are many examples of physiological processes that follow a circadian cycle such as cortisol, intraocular pressure and body temperature where abnormalities in these patterns have been shown to be related to depression [1], glaucoma [2] and delayed sleepphase disorder [3]. The ability to analyse and capture features of these cycles remains a challenge but is necessary to get a deeper understanding of the mechanisms behind them. For example, the cardiovascular system shows clear circadian rhythmicity where researchers are interested in alternative methods to illustrate and quantify this diurnal variation [4]. Circadian blood pressure (BP) represents a situation where diurnal variation deserves additional attention given the uncertainty relating to the prognostic significance of BP variability (BPV) [5–8]. The benefits of using ambulatory blood pressure monitors (ABPM) in addition to clinic measurements in the diagnosis and management of hypertension are well established [9, 10]. As well as mean day, night and dip values, ABPM provides measures of shortterm BPV and individual profile patterns. The majority of studies examining shortterm BPV have focused on summary measures such as the standard deviation (SD) of ABPM readings over the day. These summary measures are easily obtained without the need for advanced statistical techniques [5–8, 11, 12] but ignore the temporal nature of the data. To date relatively little work has modelled 24 h ABPM profiles to exploit the full potential of ABPM data to capture shortterm BPV [13]. Moreover, there are a lack of studies exploring circadian patterns and specifically, studies examining differences in patterns among different groups of individuals.
Cosinor analysis which incorporates a sinusoidal function has been the most common approach to modelling 24 h blood pressure (BP) [14–17], while a similar method, Fourier analysis [18, 19], has also been implemented. These approaches have focused on betweengroup effects (fixedeffects) where typically inferences are based on estimated differences in model parameters between particular groups of patients, such as comparing the estimated amplitude or midline estimating statistic of rhythm (MESOR) between groups of individuals on different antihypertensive agents obtained in cosinor analysis [16]. The focus of fixedeffects is on population trajectories. However one of the main advantages of ABPM is that we obtain individual BP profiles and modelling subjectspecific trajectories involves incorporating subjectspecific effects (randomeffects). To model mean profiles Selwyn and Difranco [20] used a hierarchical model incorporating a 4th degree polynomial. Lambert et al. extended on this by incorporating restricted cubic splines to model the mean BP profiles [13]. More recently, Edwards and Simpson [21] utilised orthonormal polynomials in a linear mixed model in a group of hypertensive subjects. Both polynomials and cubic splines, by their nature, have the ability to produce wellfitting curves to the data but have the disadvantage that the corresponding coefficients are challenging to interpret directly.
As an alternative we propose using piecewise linear splines in a mixedeffects model as a different approach for modelling ABPM data. Although linear splines have been used to model BP change over years and gestational age [22], using them to explore daily patterns of BP represents a novel method for analysing ABPM. This approach has the advantage that coefficients represent something meaningful, in this case the slope of BP at different periods of the day. To date it is unclear if different underlying circadian BP patterns exist across various groups of the population. This method allows slopes at a group level (and individual level) to be easily compared. Furthermore, using randomeffects we want to predict and plot curves at an individual level and to explore BPV within each period of the day. Thus the aim of this study is twofold (1) to introduce and describe a mixedeffects piecewise linear model in relation to BP; (2) to apply our method to a middleaged population sample and explore their circadian BP patterns. We also introduce and present a constraint for our model that ensures periodicity, so that on the average BP is the same 24 h later. We are particularly interested in identifying distinct differences in the shape of mean curves at a group level. For purposes of illustration of the models at a group level, we will compare those with and without evidence of subclinical target organ damage (TOD), specifically microalbuminuria.
Methods
Study population
The analysis utilises data from the Mitchelstown Study, a cross sectional study of middleaged men and women, recruited in Ireland 2010–2011. A description of the study design is available from previous publications [23, 24]. The study recruited patients attending a single large primary care centre, the LivingHealth Clinic, in Mitchelstown. Participants completed a detailed health and lifestyle survey questionnaire, and attended for a physical examination including height, weight, blood pressure, fasting blood samples and urine samples. ABPM was offered to all participants. All participants provided written informed consent and ethical approval was obtained from the Research Ethics Committee of the Cork Teaching Hospitals.
Blood pressure measurements
Study BP was measured three times after 5 min of rest in a seated position by experienced research nurses using an OMRON M7 blood pressure monitor (OMRON Healthcare, The Netherlands). The average of the second and third measurements was used for analyses. Ambulatory BP was measured using dabl ABPM system (dabl ltd., Ireland) with the Meditech ABOM05 Monitor (Meditech LTD., Hungary). The monitors were programmed to obtain readings every 30 min and remained in place for 24h. Participants kept diaries of wake and sleep periods, which were used to calculate sleep and waking times. Only participants with a minimum of 20 measurements during the day and a minimum of 7 measurements during the night period were included in the analysis. Additionally, any participants with data lacking for more than two consecutive hourly intervals were excluded [25].
Target organ damage
Each participant provided an earlymorning spot urine sample on the day of their appointment. Laboratory analyses included analysis for albumin:creatinine ratio (ACR). Microalbuminuria is defined as ACR ≥1.1 mg mmol^{−1} [26].
Statistical analysis
Linear mixed model: linear splines
The linear mixed model [27, 28] is a wellrecognised tool in the analysis of longitudinal data and its ability to obtain both population (fixedeffects) and subjectspecific (randomeffects) trajectories makes it particularly appealing for the analysis of ABPM data. However, its use todate has focused on BP following a smooth curvature trajectory which results in spline and polynomial coefficients that are of no direct clinical relevance. Piecewise linear functions or linear splines offer an alternative. These involve segregating the data into different segments across time initially assuming the segments are the same for everyone. Within each partition, a linear spline is fitted and where these are connected are known as knot points. The corresponding coefficient of each spline represents the rate of increase or decrease of BP during each time period.
Knot selection
The position of the knot points were determined based on a number of factors. Firstly, to get a general sense of the shape of the data and determine regions of interest (how many knots were required), we plotted an average curve of BP including all participants to determine common knot points. In addition we incorporated prior known characteristics of BP. The period of awakening corresponds with an abrupt and steep acceleration of BP and for many the maximum value obtained during this morning period corresponds to their maximum BP reached throughout the day [29]. We also know that BP gradually falls throughout the day and usually dips to its lowest value during the sleeping period [9]. Since waking and sleeping times are clearly important in terms of changes in BP we decided that the use of these times as two additional subjectspecific knot points was appropriate. We were able to create these subjectspecific knots using the wake and sleep times reported by the participant. We had 49 readings for each individual, where the first reading was t_{1} (12 p.m.) and the final reading was t_{49} (12 p.m. the following day). Individual waking and sleeping times were included within this range (t_{1} − t_{49}). For each individual we created n linear splines, where the kth spline:
Incorporating these linear splines into a linear mixed effects model for BP we get:
where BP _{ ij } is the BP value for the jth measurement on the ith person, at time t _{ ij }, the β's are the fixed effects coefficients associated with the average intercept at β _{0} (BP at 12 p.m.) and the average slopes (β _{ k }’s) between knot points, b’s are the random effects associated with the average intercept (b _{0i }) and average slopes between knot points, and ε _{ ij } representing the individuallevel residuals from the model. The model is extended by incorporating the subjectspecific knots in the s _{ k } term. It is assumed the random effects (b _{ i }) have zero mean and an unstructured variance–covariance matrix Σ_{ b }. The individual level residuals have mean zero and variance–covariance matrix Σ_{ε}.
To expand Eq. (2) to include the restriction that on the average BP is the same 24 h later we define an equation which states average, subjectspecific change in BP over 24 h is zero:
where w _{ ki } is the width of the kth interval (and where those involving wake and sleep times are subjectspecific width intervals). Rewriting this in terms of β _{1} gives:
which implies:
where \(s_{ki}^{*} = s_{ki}  \frac{{w_{ki} }}{{w_{1i} }}s_{1i} ,\) which allows us rewrite (2) as
To explore a group effect the model can easily incorporate a variable of interest, in this case TOD (microalbuminuria), as a dichotomous covariate. We further extended the model allowing the shape of the trajectory to depend on TOD by including interactions between TOD and each linear spline slope. Comparing this model with one without any interactions allowed us to test if the overall trajectory of BP was different between the two groups across the day. Additionally we were able to test if slopes between the groups differed at specific locations throughout the day. We adjusted for confounders by adding them into the model as fixed effects. In additional models we tested the effect of allowing the residual variance to differ between those with and without microalbuminuria. Similarly we tested the impact of allowing the interaction terms of microalbuminuria with each linear spline to be random to determine if there was heterogeneity of variance between the groups at any period of the day. Although we used an unstructured covariance structure for our models, we assumed these interaction terms to be independent of the other randomeffects parameters. These interactions represented the difference in variation between the microalbuminuria groups within each segment. For all models explored we allowed all the linear spline terms to be random.
As individual ABPM readings taken close in time are likely to be correlated, a model with an independent residual correlation structure may not be appropriate. We compared this to a model with a firstorder autoregressive AR(1) structure and examined a plot of the autocorrelation function (ACF) to detect violations of the assumption of independence. Allowing for temporal correlation can potentially result in a large improvement in the precision of parameter estimates [30].
Models were compared formally by a likelihood ratio test (LRT) [28, 31]. The appropriate variance and residual function structures were also identified using a LRT in addition with an ACF plot. Rsquared (R ^{2}) statistic is often presented as a summary measure for linear models but due to theoretical or practical problems is rarely presented for mixedmodels. Nakagawa and Schielzeth discuss these issues and present a general but simple method for calculating an appropriate R ^{2} for random intercept mixedmodels [32]. Johnson extended this to include random slope models which we implement in our analysis [33].
The parameters for our final models were estimated using restricted maximum likelihood estimation (REML) as this method produces unbiased estimates unlike maximum likelihood (ML) estimation [30]. Subjectspecific trajectories were estimated using Empirical Best Linear Unbiased Predictors (EBLUPs) of the randomeffects [28]. Residual diagnostic plots were examined to verify model assumptions including normality of both randomeffects and residual errors and that the error terms had constant variance. In addition a visual predictive check (VPC) was performed in which the estimated mean and the 90% prediction interval from our model were plotted together with the observed BP values and the 90% interquantile range of the observations. The purpose of the VPC is to assess graphically if predictions from the fitted model reproduce the central trend and variability of BP in the observed data, when plotted against time. It is an internal validation method that assesses the goodnessoffit [34]. All analysis was completed for both SBP and DBP. All analysis were implemented in R [35] and parameter estimation for the mixedeffect model was carried out by means of the lme command in nlme package [36]. Sample R code is provided to ensure reproducibility and improve dissemination of methods (see Additional file 1).
Validation
Although difficult to interpret the coefficients, polynomial regression can still be a useful tool in the analysis of medical data to plot the trajectory of nonlinear or curvilinear relationships [37]. As a method of validation for our approach we additionally implemented a linear mixed model with orthogonal polynomials across time in both the fixed and random effects, similar to that of Edwards and Simpson [21]. We wanted to determine if linear splines were capable of capturing the circadian rhythm of BP. This was investigated by comparing the trajectories obtained from both methods to determine if they followed similar patterns. A similar process to that of the piecewise model was followed when fitting the polynomial model. As we were only concerned to know if the general shape could be captured by our piecewise approach we were not worried about overfitting the polynomial model. We implemented a model up to a 6th order polynomial allowing all the terms to be random.
Results
Of 3051 individuals invited to participate, 2047 (response rate: 67%) completed the questionnaire and physical examination component. ABPM was offered to all 2047 participants and it was completed by 1207 (response rate: 58%) people, of whom 1008 had a minimum of 20 day and 7 night measurements respectively. Of these 886 had no data missing for more than two consecutive hourly intervals, and the main clinical characteristics of these participants are presented in Table 1. Overall, participants had a mean age of 59.9 (5.5) and the majority were female (55%). Sixty percent were classified as hypertensive. Also presented in Table 1 are the characteristics of the full sample which shows the ABPM subsample follows a similar distribution in terms of age, sex, BMI and the presence of microalbuminuria. However, the proportion of those with hypertension was higher amongst those in ABPM sub group than in overall study population (60 vs 47%).
The plot of average SBP for the 886 subjects is presented in Fig. 1. Based on this plot we identified two common knot points where the trajectory of SBP changed notably at 6 p.m. and 4 a.m. In addition to these two points we were able to include two subjectspecific knot points for each participant based on the time an individual woke and went to sleep (Fig. 1). This meant each participant was assigned 4 knot points which in turn resulted in their SBP pattern being broken into 5 linear segments.
Figure 2 represents subjectspecific trajectories as a function of time only from a linear mixedeffects model using both orthogonal polynomials (6th order) (red line) and piecewise linear splines (blue lines). The plots suggest that the individual curves can be adequately captured by using piecewise linear splines.
We initially included the 5 linear splines as fixedeffects. A significant improvement in fit was observed when additionally including each term individually as a randomeffect, based on a LRT (all p < 0.001). As a consequence we included all the linear spline terms as random effects. To allow for temporal correlation we incorporated an AR1 structure which resulted in a significant improvement in fit (p < 0.001) (rho = 0.27). Examining the ACF plot indicated that the inclusion of an AR1 residual structure adequately accounted for the autocorrelation in the data. This unadjusted model, which only incorporates linear splines as a function of time was our base model, Table 2 (Model 1). Presented are the parameter estimates (fixedeffects, randomeffects correlation matrix, the autocorrelation decay ρ along with fit criteria values). With the exception of the slope for the period from 12.00 to 18.00 [0.02 (0.04) mmHg/30 min], all slopes differed significantly from zero (all p < 0.001). This suggests that on average this is the period during the day where average SBP remains constant. The largest rise and fall in SBP occurred between wake and 12.00 (2.23 mmHg/30 min) and, between sleep and 04.00 (−1.93 mmHg/30 min) respectively. These segments correspond to the period when an individual wakes up and the period immediately after they fall asleep. The variation in slopes was lowest from 12.00 to 18.00 where the variance was 0.51. The largest variation in slopes was observed between waking and 12.00 where the variance was 2.05 which is substantially larger in comparison to the rest of the day. The model R ^{2} value which illustrates the proportion of variance explained by both the fixed and random factors was quite high (0.67).
In subsequent models we adjusted for age, sex and BMI. We also included our variable of interest, microalbuminuria, to determine if it could help explain the larger variation in the period, wake to 12.00 (Model 2, Table 2). The residual variance, which represents the variation about an individual’s trajectory, was 12.3 mmHg. We additionally allowed the residual variance to vary between microalbuminuria groups (ratio of standard deviation of those with to without microalbuminuria was 1.09). On average, over the day, those with microalbuminuria were found to have significantly higher SBP (7.6 mmHg, 95% CI 5.0–10.1, p < 0.001). However, adjusting for age, sex, BMI and microalbuminuria had almost no effect on the model parameter estimates (except the intercept). To determine if slopes were different between groups at different times of the day we included an interaction between each linear spline and microalbuminuria (Model 3, Table 2). Although two of the interaction terms were marginally significant, a LRT suggested that including interaction terms did not improve the overall fit to the data (p = 0.12). Based on additional models (results not shown) we found no evidence that the variance of the randomeffects varied with microalbuminuria. With the inclusion of age, sex, BMI and microalbuminuria we concluded that Model 2 offered the best fit to the data (Model 2 vs Model 1, p < 0.01). Residual diagnostic plots of the models (residuals vs fitted values, histogram of randomeffects & residuals) showed no violation of assumptions (results not shown). The VPC plot showed the model was adequately predicting central trend and variability of SBP in the observed data, when plotted against time (see Additional file 2).
Figure 3 represents the average piecewise linear curve along with a 95% confidence interval for those with and without the presence of microalbuminuria using Model 2. The numbers on the plot correspond to the time periods presented in Table 2. It is clear that those with microalbuminuria have a higher average SBP throughout the day. For the purposes of this plot we have set the sleep and wake time knots at 23.00 and 08.00 respectively. A similar plot using model 3 can be found in the Additional file 2. Similar findings were found for all analysis when repeated using DBP (results not shown).
Discussion
In this large population based study we present an alternative method of modelling 24 h BP that can easily be applied to any physiological process that follows a circadian cycle. Our novel but simple approach utilising a piecewise linear randomeffects model, with an adjustment to ensure that the average level is the same at the beginning and end of each 24h period, offers a practical alternative to other methodological modelling techniques for researchers exploring circadian patterns. The flexible model has the ability to capture overall average, group and individual trajectories (in addition to being capable of examining slopes at different periods of the day).
Despite the large amounts of literature relating to BP, those specifically modelling 24 h ABPM remain sparse. Our method offers new measures of shortterm BP variability as we can quantify the variation about an individual’s trajectory but it also allows examination of the variation in slopes between individuals (randomeffects). Our results indicated that after adjustment for age, sex and BMI the sharpest fall in BP occurred just after an individual went to sleep and the steepest rise occurred just after waking. Although there was a significant difference on average between those with and without microalbuminuria we found there was no overall improvement in fit after including interaction effects with the spline terms. However interestingly we found that the variation after awaking, representing what is known as the morning surge was considerably larger than the other periods of the day.
It has been acknowledged there is not a generally accepted “standard” method of analysing 24h ABPM [21]. Cosinor analysis has been highlighted as the most common approach [14–17] while fourier analysis [18], has also been implemented which are both based on the idea that any time series can be described by a series of cosine (and sine) waves of various frequencies [38]. It has been suggested that these methods impose too many restrictions on the shape of the profile and have been shown to fit real profiles poorly [39]. Wang et al. [40] suggest problems with fitting a sinusoidal function to a circadian pattern include (1) that the pattern over time may not be symmetric; that is, the peak and nadir may not be separated by 12 h and/or the amplitude and width of the peak may differ from those of the nadir, (2) sometimes there are local minimum and maximum points. Additionally Wang et al. [41] suggests that the sinusoidal function is too restrictive and “rhythms with a shape closely approximating a cosine curve are uncommon” [42]. Alternative methods have examined restricted cubic splines and more recently orthonormal polynomials [13, 21]. As we highlighted previously these approaches may model the data quite well and their curvature nature may look graphically appealing but it is difficult to understand and compare their resulting coefficients.
Piecewise regression which allows separate slopes to be fitted to observations before and after a certain period or event (knot points) has been cited as a useful tool that should be implemented more often in the context of epidemiological studies [43] but has not, to the best of our knowledge been used with ABPM data or other physiological processes that are circadian. The benefit of this method as opposed to polynomials is that the regression coefficients represent something meaningful directly without the need for further manipulation of the results—in our context, the rate of increase or decrease of BP for a certain time of day. The position of the knot points can easily be altered depending on the requirements of a specific study. For example if we were examining the effect of dialysis on BP in haemodialysis patients we could fix knot points at the time their dialysis began and at period(s) a number of hours later.
The morning is recognised as the most important period in relation to cardiovascular diseases [44] and cardiovascular events occur more frequently in this period [44–46]. In our study we found that the steepest rise (slope) occurred during the period just after waking which is in line with the literature, thus verifying that our method is capturing known features of the data. It is suggested that the abrupt steep rise in BP may explain the link between cardiovascular events and the morning period [29]. In a review of morning surge with cardiovascular risk, three different definitions of morning surge were identified, all of which simply use BP differences where they subtracted some average night value minus an average of morning BP readings [45]. We argue that our method offers a more accurate estimate, as by definition of a slope we can specifically quantify the rate of “surge”. In fact, Parati et al. argue that a method that would be capable of capturing a slope similar to one purposed by our method would provide an accurate method of estimating the morning surge [47]. Considering that morning surge has been cited as a predictor of stroke and advanced target organ damage independent of ambulatory BP and nocturnal BP [44, 46], accurately quantifying it remains an important issue, particularly when we are assessing the benefits of antihypertensive medication in their ability to reduce this steep rise. This may not only have health implications but also financial benefits. A similar argument could be put forward for the dipping effect at night which is usually quantified just as a ratio of the mean BP between nightday periods. The slope at night obtained by our approach may represent a more accurate measure but further work would be needed to explore this.
Kario argues that the perfect 24 h BP control is not limited to reducing mean BP but includes restoring disrupted circadian BP rhythms and reducing exaggerated BP variability [44]. As highlighted previously most studies examining BPV have concentrated on summary measures of variability such as SD over 24 h or separated into day and night values [5–8, 11, 12, 24]. With the use of our mixedeffects model we were are able to obtain superior measures of BPV that take into account the temporal nature of the data. We were able to quantify the variation about an individual’s trajectory but also the variation in slopes between individuals. Our work highlighted that the largest variation between individuals occurred during the morning surge period. Adjusting for age, sex and BMI did not help explain this variation. Similarly the presence of microalbuminuria had little impact on the variation. Ideally we would have preferred to explore if the variation could in part predict CV events but as data is currently only available for wave one, we have been restricted to explore a surrogate marker in microalbuminuria and have acknowledged this as a limitation. Further work is warranted to include CVD endpoints but perhaps an underlining physiological phenomenon of BP is that it is most variable in the morning possibly because this period of the day has an abrupt rise. Although some of the knots are subjectspecific, others are at common fixed locations which may not represent the best position for a specific individual and this assumption is recognised as a limitation. In addition to the average plot, we have attempted to incorporate our knowledge of the underlying pattern of BP to help inform our knot positions as suggested by Howe et al. [48].
As debate remains in relation to how to correctly quantify shortterm BPV [24, 49] our approach offers new alternatives that utilise the full power of ABPM that is often lost when using summary measures such as SD because it only reflects the dispersion of measurements around a single value (mean) not accounting for the order in which BP measurements were obtained [12, 50]. As argued by Bamberger et al. there are a large and increasing variety of mathematical equations available which allow us to test specific hypotheses about, and estimate parameters of, growth or change (BP) with varying purposes such as summarizing the shape of change across a sample or determining patterns in variability [51]. Despite this however, Diehl et al. suggests that summary measures and multilevel models described here are two of the most common techniques used specifically for exploring BP changes [52]. They also suggest the use of dynamic developmental models with multiple time scales as another approach to address questions related to intraindividual variability. This is quite similar to the multilevel model framework but is slightly more flexible and can include more complex scenarios. These models are characterized by timedependent parametric changes occurring at different time scales, where each time scale defines a distinct level within a hierarchy of time scales [53]. In the context of our work this approach may be useful if we had additional data that could incorporate another time scale, such as BP over years e.g. the lifespan. The benefit from this would be that we could compare and relate BP fluctuations over 24 h to that over the lifespan. However, in this study we are only interested in 24 h variation where we only have recordings taken every 30 min. With only one time scale, we believe that the use of the traditional mixedeffects model is appropriate for our analysis. We argue that our approach, which has the ability to determine variation over specific periods of the day offers a novel measure of variability in the analysis of 24 h BP which may have benefits when attempting to determine the optimal timing of antihypertensive medication administration in future studies. Finally, as was briefly alluded to, the approach and discussion outlined is not restricted to the use of BP and can easily be implemented on any physiological process that demonstrates a circadian cycle. BP is not the only biological process where disruptions to circadian rhythms are clinical relevant. Wang et al. found that those with Cushing syndrome exhibited no circadian rhythm of cortisol, while those with depression showed a dampened rhythm compared to the normal group [41]. Liu et al. found that larger shortterm fluctuations in intraocular pressure are more common in glaucoma [54]. Similar to the morning BP surge, it was found that intraocular pressure was higher in the morning and more prevalent in those with glaucoma. This suggests that our approach may be beneficial to the exploration of other biological rhythms that have similar features to that of BP.
Conclusion
This study has introduced a novel but practical method for the analysis of ABPM data. Based on our work circadian BP patterns can be modelled using a mixedeffects model with piecewise linear splines. The main advantage of our method compared to other approaches is that the resulting regression coefficients have direct interpretation. We can determine the rate of increase or decrease at different periods of the day. In addition we can determine alternative measures of variability compared to classical BPV indices. Future research in this area should focus on the association between the measures obtained from this method to stronger clinical outcomes.
Abbreviations
 ABPM:

ambulatory blood pressure monitoring
 BP:

blood pressure
 BPV:

blood pressure variability
 SD:

standard deviation
 ACR:

albumin:creatinine ratio
 REML:

restricted maximum likelihood estimation
 EBLUP:

Empirical Best Linear Unbiased Predictors
 VPC:

visual predictive check
 AR(1):

autoregressive
References
 1.
SalgadoDelgado R, Tapia Osorio A, Saderi N, Escobar C. Disruption of circadian rhythms: a crucial factor in the etiology of depression. Depress Res Treat. 2011;2011:839743.
 2.
Sihota R, Saxena R, Gogoi M, Sood A, Gulati V, Pandey RM. A comparison of the circadian rhythm of intraocular pressure in primary phronic angle closure glaucoma, primary open angle glaucoma and normal eyes. Indian J Ophthalmol. 2005;53(4):243–7.
 3.
Zee PC, Attarian H, Videnovic A. Circadian rhythm abnormalities. Continuum (Minneap Minn). 2013;19(1 Sleep Disorders):132–47.
 4.
Chen L, Yang G. Recent advances in circadian rhythms in cardiovascular system. Front Pharmacol. 2015;6:71.
 5.
Hansen TW, Thijs L, Li Y, Boggia J, Kikuya M, BjorklundBodegard K, et al. Prognostic value of readingtoreading blood pressure variability over 24 hours in 8938 subjects from 11 populations. Hypertension. 2010;55(4):1049–57.
 6.
Kikuya M, Hozawa A, Ohokubo T, Tsuji I, Michimata M, Matsubara M, et al. Prognostic significance of blood pressure and heart rate variabilities: the Ohasama study. Hypertension. 2000;36(5):901–6.
 7.
Rothwell PM, Howard SC, Dolan E, O’Brien E, Dobson JE, Dahlof B, et al. Prognostic significance of visittovisit variability, maximum systolic blood pressure, and episodic hypertension. Lancet. 2010;375(9718):895–905.
 8.
Schillaci G, Verdecchia P, Borgioni C, Ciucci A, Porcellati C. Lack of association between blood pressure variability and left ventricular mass in essential hypertension. Am J Hypertens. 1998;11(5):515–22.
 9.
Mancia G, Fagard R, Narkiewicz K, Redon J, Zanchetti A, Bohm M, et al. 2013 ESH/ESC Guidelines for the management of arterial hypertension: the Task Force for the management of arterial hypertension of the European Society of Hypertension (ESH) and of the European Society of Cardiology (ESC). J Hypertens. 2013;31(7):1281–357.
 10.
O’Brien E, Parati G, Stergiou G, Asmar R, Beilin L, Bilo G, et al. European Society of Hypertension position paper on ambulatory blood pressure monitoring. J Hypertens. 2013;31(9):1731–68.
 11.
Roman MJ, Pickering TG, Schwartz JE, Pini R, Devereux RB. Relation of blood pressure variability to carotid atherosclerosis and carotid artery and left ventricular hypertrophy. Arterioscler Thromb Vasc Biol. 2001;21(9):1507–11.
 12.
Pierdomenico SD, Di Nicola M, Esposito AL, Di Mascio R, Ballone E, Lapenna D, et al. Prognostic value of different indices of blood pressure variability in hypertensive patients. Am J Hypertens. 2009;22(8):842–7.
 13.
Lambert PC, Abrams KR, Jones DR, Halligan AW, Shennan A. Analysis of ambulatory blood pressure monitor data using a hierarchical model incorporating restricted cubic splines and heterogeneous withinsubject variances. Stat Med. 2001;20(24):3789–805.
 14.
Portaluppi F, Montanari L. Consistency of circadian blood pressure pattern assessed by noninvasive monitoring and cosinor analysis in hospitalized hypertensive patients. Acta Cardiol. 1988;43(5):605–13.
 15.
Portaluppi F, Bagni B, degli Uberti E, Montanari L, Cavallini R, Trasforini G, et al. Circadian rhythms of atrial natriuretic peptide, renin, aldosterone, cortisol, blood pressure and heart rate in normal and hypertensive subjects. J Hypertens. 1990;8(1):85–95.
 16.
Munakata M, Imai Y, Hashimoto J, Sakuma H, Sekino H, Abe K, et al. The influence of antihypertensive agents on circadian rhythms of blood pressure and heart rate in patients with essential hypertension. Tohoku J Exp Med. 1992;166(2):217–27.
 17.
Ayala DE, Hermida RC, Mojon A, Fernandez JR, Iglesias M. Circadian blood pressure variability in healthy and complicated pregnancies. Hypertension. 1997;30(3 Pt 2):603–10.
 18.
Somes GW, Harshfield GA, Arheart KL, Miller ST. A Fourier series approach for comparing groups of subjects on ambulatory blood pressure patterns. Stat Med. 1994;13(12):1201–10.
 19.
Staessen JA, Fagard R, Thijs L, Amery A. Fourier analysis of blood pressure profiles. Am J Hypertens. 1993;6(6 Pt 2):184s–7s.
 20.
Selwyn MR, Difranco DM. The application of large Gaussian mixed models to the analysis of 24 hour ambulatory blood pressure monitoring data in clinical trials. Stat Med. 1993;12(18):1665–82.
 21.
Edwards LJ, Simpson SL. An analysis of 24h ambulatory blood pressure monitoring data using orthonormal polynomials in the linear mixed model. Blood Press Monit. 2014;19(3):153–63.
 22.
MacdonaldWallis C, Tilling K, Fraser A, Nelson SM, Lawlor DA. Associations of blood pressure change in pregnancy with fetal growth and gestational age at delivery: findings from a prospective cohort. Hypertension. 2014;64(1):36–44.
 23.
Kearney PM, Harrington JM, Mc Carthy VJ, Fitzgerald AP, Perry IJ. Cohort profile: the Cork and Kerry diabetes and heart disease study. Int J Epidemiol. 2013;42(5):1253–62.
 24.
Madden JM, O’Flynn AM, Dolan E, Fitzgerald AP, Kearney PM. Shortterm blood pressure variability over 24 h and target organ damage in middleaged men and women. J Hum Hypertens. 2015;29:719–25.
 25.
Hermida RC, Smolensky MH, Ayala DE, Portaluppi F. 2013 ambulatory blood pressure monitoring recommendations for the diagnosis of adult hypertension, assessment of cardiovascular and other hypertensionassociated risk, and attainment of therapeutic goals. Chronobiol Int. 2013;30(3):355–410.
 26.
Matsushita K, van der Velde M, Astor BC, Woodward M, Levey AS, de Jong PE, et al. Association of estimated glomerular filtration rate and albuminuria with allcause and cardiovascular mortality in general population cohorts: a collaborative metaanalysis. Lancet. 2010;375(9731):2073–81.
 27.
Laird NM, Ware JH. Randomeffects models for longitudinal data. Biometrics. 1982;38(4):963–74.
 28.
Fitzmaurice G, Laird NM, Ware JH. Applied longitudinal analysis. Hoboken: Wiley; 2011.
 29.
Parati G, Bilo G, Redon J. Morning and smooth 24h ambulatory blood pressure control is not achieved in general practice: results from the SURGE observational study. J Hypertens. 2013;31(3):616–23 (discussion 23).
 30.
Zuur A, Ieno EN, Walker NJ, Saveliev AA, Smith GM. Mixed effects models and extensions in ecology with R. Berlin: Springer; 2009.
 31.
Pinheiro JCD, Bates DM. Mixedeffects models in S and SPLUS. New York: Springer; 2000.
 32.
Nakagawa S, Schielzeth H. A general and simple method for obtaining R2 from generalized linear mixedeffects models. Methods Ecol Evolution. 2013;4(2):133–42.
 33.
Johnson PC. Extension of Nakagawa & Schielzeth’s to random slopes models. Methods Ecol Evol. 2014;5(9):944–6.
 34.
Bergstrand M, Hooker AC, Wallin JE, Karlsson MO. Predictioncorrected visual predictive checks for diagnosing nonlinear mixedeffects models. AAPS J. 2011;13(2):143–51.
 35.
Team; RC. R: a language and environment for statistical computing. Vienna: R Foundation for Statistical Computing; 2014. http://www.Rproject.org/.
 36.
Pinheiro JBD, DebRoy S, Sarkar D, R Core Team. nlme: Linear and nonlinear mixed effects models. R package version 31119. 2015.
 37.
Armitage P, Berry G, Matthews JNS. Statistical methods in medical research. Oxford: Blackwell; 2009.
 38.
Refinetti R, Lissen GC, Halberg F. Procedures for numerical analysis of circadian rhythms. Biol Rhythm Res. 2007;38(4):275–325.
 39.
Streitberg B, MeyerSabellek W, Baumgart P. Statistical analysis of circadian blood pressure recordings in controlled clinical trials. J Hypertens Suppl. 1989;7(3):S11–7.
 40.
Wang Y, Brown MB. A flexible model for human circadian rhythms. Biometrics. 1996;52(2):588–96.
 41.
Wang Y, Ke C, Brown MB. Shapeinvariant modeling of circadian rhythms with random effects and smoothing spline ANOVA decompositions. Biometrics. 2003;59(4):804–12.
 42.
Reilly T, Atkinson G, Waterhouse J. Biological rhythms and exercise. New York: Oxford; 1997.
 43.
Hardy R. Commentary: are piecewise mixed effects models useful in epidemiology? Int J Epidemiol. 2001;30(6):1341–2.
 44.
Kario K. Essential manual of 24 hour blood pressure management: from morning to nocturnal hypertension. New York: WileyBlackwell; 2015.
 45.
Kario K. Morning surge in blood pressure and cardiovascular risk: evidence and perspectives. Hypertension. 2010;56(5):765–73.
 46.
Kario K, Pickering TG, Umeda Y, Hoshide S, Hoshide Y, Morinari M, et al. Morning surge in blood pressure as a predictor of silent and clinical cerebrovascular disease in elderly hypertensives: a prospective study. Circulation. 2003;107(10):1401–6.
 47.
Parati G, Vrijens B, Vincze G. Analysis and interpretation of 24h blood pressure profiles: appropriate mathematical models may yield deeper understanding. Am J Hypertens. 2008;21(2):123–5 (discussion 7–9).
 48.
Howe LD, Tilling K, Matijasevich A, Petherick ES, Santos AC, Fairley L, et al. Linear spline multilevel models for summarising childhood growth trajectories: a guide to their application using examples from five birth cohorts. Stat Methods Med Res. 2016;25(5):1854–74.
 49.
Taylor KS, Heneghan CJ, Stevens RJ, Adams EC, Nunan D, Ward A. Heterogeneity of prognostic studies of 24hour blood pressure variability: systematic review and metaanalysis. PLoS ONE. 2015;10(5):e0126375.
 50.
Bilo G, Giglio A, Styczkiewicz K, Caldara G, Maronati A, KaweckaJaszcz K, et al. A new method for assessing 24h blood pressure variability after excluding the contribution of nocturnal blood pressure fall. J Hypertens. 2007;25(10):2058–66.
 51.
Bamberger KT. The application of intensive longitudinal methods to investigate change: stimulating the field of applied family research. Clin Child Fam Psychol Rev. 2016;19(1):21–38.
 52.
Diehl M, Hooker K, Sliwinski MJ. Handbook of intraindividual variability across the life span. New York, NY: Routledge, Taylor & Francis Group; 2015.
 53.
Boker SM, Molenaar PC, Nesselroade JR. Issues in intraindividual variability: individual differences in equilibria and dynamics over multiple time scales. Psychol Aging. 2009;24(4):858–62.
 54.
Liu JH, Zhang X, Kripke DF, Weinreb RN. Twentyfourhour intraocular pressure pattern associated with early glaucomatous changes. Invest Ophthalmol Vis Sci. 2003;44(4):1586–90.
Authors’ contributions
JM drafted the initial manuscript and performed analysis, JM, AF and XL conceptualized the study and designed the methodology, AF and PK designed the original study, AF, XL, KT and PK critically reviewed the manuscript. All authors read and approved the final manuscript.
Acknowledgements
We especially thank all participants in the study, the study nurses and administrators and the staff at the LivingHealth Clinic.
Competing interests
The authors declare that they have no competing interests.
Ethics approval and consent
All participants provided written informed consent and ethical approval was obtained from the Research Ethics Committee of the Cork Teaching Hospitals.
Funding
This research was funded by the Health Research Board PhD/2007/16. The Cork and Kerry Diabetes and Heart Disease Study was funded by a research grant from the Irish Health Research Board (Ref. HRC 2007/13).
Financial disclosure statement
The authors have no financial relationships relevant to this article to disclose.
Author information
Affiliations
Corresponding author
Additional files
Additional file 1.
Example of statistical code in R.
Additional file 2.
VPC plot and predicted average piecewise linear trajectory for Model 3 (interaction model).
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Madden, J.M., Li, X., Kearney, P.M. et al. Exploring diurnal variation using piecewise linear splines: an example using blood pressure. Emerg Themes Epidemiol 14, 1 (2017). https://doi.org/10.1186/s1298201700555
Received:
Accepted:
Published:
Keywords
 Circadian modeling
 Mixedeffects models
 Biostatistics
 Blood pressure variability
 Blood pressure patterns