Shape‐Invariant Modeling of Circadian Rhythms with Random Effects and Smoothing Spline ANOVA Decompositions

Fibroblasts as an in vitro model of circadian genetic and genomic studies

Bipolar disorder (BD) is a heritable disorder characterized by shifts in mood that manifest in manic or depressive episodes. Clinical studies have identified abnormalities of the circadian system in BD patients as a hallmark of underlying pathophysiology. Fibroblasts are a well-established in vitro model for measuring circadian patterns. We set out to examine the underlying genetic architecture of circadian rhythm in fibroblasts, with the goal to assess its contribution to the polygenic nature of BD disease risk. We collected, from primary cell lines of 6 healthy individuals, temporal genomic features over a 48 h period from transcriptomic data (RNA-seq) and open chromatin data (ATAC-seq). The RNA-seq data showed that only a limited number of genes, primarily the known core clock genes such as ARNTL, CRY1, PER3, NR1D2 and TEF display circadian patterns of expression consistently across cell cultures. The ATAC-seq data identified that distinct transcription factor families, like those with the basic helix-loop-helix motif, were associated with regions that were increasing in accessibility over time. Whereas known glucocorticoid receptor target motifs were identified in those regions that were decreasing in accessibility. Further evaluation of these regions using stratified linkage disequilibrium score regression analysis failed to identify a significant presence of them in the known genetic architecture of BD, and other psychiatric disorders or neurobehavioral traits in which the circadian rhythm is affected. In this study, we characterize the biological pathways that are activated in this in vitro circadian model, evaluating the relevance of these processes in the context of the genetic architecture of BD and other disorders, highlighting its limitations and future applications for circadian genomic studies.

Fibroblasts as an in vitro model of circadian genetic and genomic studies: A temporal analysis

Bipolar disorder (BD) is a heritable disorder characterized by shifts in mood that manifest in manic or depressive episodes. Clinical studies have identified abnormalities of the circadian system in BD patients as a hallmark of underlying pathophysiology. Fibroblasts are a well-established in vitro model for measuring circadian patterns. We set out to examine the underlying genetic architecture of circadian rhythm in fibroblasts, with the goal to assess its contribution to the polygenic nature of BD disease risk. We collected, from primary cell lines of 6 healthy individuals, temporal genomic features over a 48 hour period from transcriptomic data (RNA-seq) and open chromatin data (ATAC-seq). The RNA-seq data showed that only a limited number of genes, primarily the known core clock genes such as ARNTL, CRY1, PER3, NR1D2 and TEF display circadian patterns of expression consistently across cell cultures. The ATAC-seq data identified that distinct transcription factor families, like those with the basic helix-loop-helix motif, were associated with regions that were increasing in accessibility over time. Whereas known glucocorticoid receptor target motifs were identified in those regions that were decreasing in accessibility. Further evaluation of these regions using stratified linkage disequilibrium score regression (sLDSC) analysis failed to identify a significant presence of them in the known genetic architecture of BD, and other psychiatric disorders or neurobehavioral traits in which the circadian rhythm is affected. In this study, we characterize the biological pathways that are activated in this in vitro circadian model, evaluating the relevance of these processes in the context of the genetic architecture of BD and other disorders, highlighting its limitations and future applications for circadian genomic studies.

LimoRhyde2: Genomic analysis of biological rhythms based on effect sizes

Genome-scale data have revealed daily rhythms in various species and tissues. However, current methods to assess rhythmicity largely restrict their focus to quantifying statistical significance, which may not reflect biological relevance. To address this limitation, we developed a method called LimoRhyde2 (the successor to our method LimoRhyde), which focuses instead on rhythm-related effect sizes and their uncertainty. For each genomic feature, LimoRhyde2 fits a curve using a series of linear models based on periodic splines, moderates the fits using an Empirical Bayes approach called multivariate adaptive shrinkage (Mash), then uses the moderated fits to calculate rhythm statistics such as peak-to-trough amplitude. The periodic splines capture non-sinusoidal rhythmicity, while Mash uses patterns in the data to account for different fits having different levels of noise. To demonstrate LimoRhyde2's utility, we applied it to multiple circadian transcriptome datasets. Overall, LimoRhyde2 prioritized genes having high-amplitude rhythms in expression, whereas a prior method (BooteJTK) prioritized "statistically significant" genes whose amplitudes could be relatively small. Thus, quantifying effect sizes using approaches such as LimoRhyde2 has the potential to transform interpretation of genomic data related to biological rhythms.

Meta-analysis of Diurnal Transcriptomics in Mouse Liver Reveals Low Repeatability of Rhythm Analyses

To assess the consistency of biological rhythms across studies, 57 public mouse liver tissue timeseries totaling 1096 RNA-seq samples were obtained and analyzed. Only the control groups of each study were included, to create comparable data. Technical factors in RNA-seq library preparation were the largest contributors to transcriptome-level differences, beyond biological or experiment-specific factors such as lighting conditions. Core clock genes were remarkably consistent in phase across all studies. Overlap of genes identified as rhythmic across studies was generally low, with no pair of studies having over 60% overlap. Distributions of phases of significant genes were remarkably inconsistent across studies, but the genes that consistently identified as rhythmic had acrophase clustering near ZT0 and ZT12. Despite the discrepancies between single-study analyses, cross-study analyses found substantial consistency. Running compareRhythms on each pair of studies identified a median of only 11% of the identified rhythmic genes as rhythmic in only 1 of the 2 studies. Data were integrated across studies in a joint and individual variance estimate (JIVE) analysis, which showed that the top 2 components of joint within-study variation are determined by time of day. A shape-invariant model with random effects was fit to the genes to identify the underlying shape of the rhythms, consistent across all studies, including identifying 72 genes with consistently multiple peaks.

Uncovering circadian rhythms in metabolic longitudinal data: A Bayesian latent class modeling approach

Researchers in biology and medicine have increasingly focused on characterizing circadian rhythms and their potential impact on disease. Understanding circadian variation in metabolomics, the study of chemical processes involving metabolites may provide insight into important aspects of biological mechanism. Of scientific importance is developing a statistical rigorous approach for characterizing different types of 24-hour patterns among high dimensional longitudinal metabolites. We develop a latent class approach to incorporate variation in 24-hour patterns across metabolites where profiles are modeled with finite mixtures of distinct shape-invariant circadian curves that themselves incorporate variation in amplitude and phase across metabolites. An efficient Markov chain Monte Carlo sampling is used to carry out Bayesian posterior computation. When the model was fit separately by individual to the data from a small group of participants, two distinct 24-hour rhythms were identified, with one being sinusoidal and the other being more complex with multiple peaks. Interestingly, the latent pattern associated with circadian variation (simple sinusoidal curve) had a similar phase across the three participants, while the more complex latent pattern reflecting diurnal variation differed across individual. The results suggested that this modeling framework can be used to separate 24-hour rhythms into an endogenous circadian and one or more exogenous diurnal patterns in describing human metabolism.

LimoRhyde2: genomic analysis of biological rhythms based on effect sizes

Genome-scale data have revealed daily rhythms in various species and tissues. However, current methods to assess rhythmicity largely restrict their focus to quantifying statistical significance, which may not reflect biological relevance. To address this limitation, we developed a method called LimoRhyde2 (the successor to our method LimoRhyde), which focuses instead on rhythm-related effect sizes and their uncertainty. For each genomic feature, LimoRhyde2 fits a curve using a series of linear models based on periodic splines, moderates the fits using an Empirical Bayes approach called multivariate adaptive shrinkage (Mash), then uses the moderated fits to calculate rhythm statistics such as peak-to-trough amplitude. The periodic splines capture non-sinusoidal rhythmicity, while Mash uses patterns in the data to account for different fits having different levels of noise. To demonstrate LimoRhyde2's utility, we applied it to multiple circadian transcriptome datasets. Overall, LimoRhyde2 prioritized genes having high-amplitude rhythms in expression, whereas a prior method (BooteJTK) prioritized "statistically significant" genes whose amplitudes could be relatively small. Thus, quantifying effect sizes using approaches such as LimoRhyde2 has the potential to transform interpretation of genomic data related to biological rhythms.

Core and peripheral site measurement of body temperature in short wool sheep

Understanding circadian rhythms of body temperature is important for the interpretation of single body temperature measurements and the assessment of the physiological state of an animal. The ability to measure body temperature at peripheral locations may also be important in the development of minimally invasive tools for remote temperature measurement in livestock. This study aimed to investigate how well body temperature measured at peripheral sites reflected a commonly used core measurement (vaginal temperature) and the circadian rhythmicity of the body temperature of sheep with a view to practical application in extensive sheep production systems. Eleven crossbred ewes were implanted with peripheral temperature sensing microchips (LifeChip®) which were positioned transversely in the sternocleidomastoid (neck) muscle and subcutaneously under the tail. iButton® temperature loggers were placed intravaginally to record core body temperature measurements (Tv). The body temperature measurements observed at the peripheral sites in the neck (Tn) and tail (Tt) differed significantly to those measured at the core site, Tv (P < 0.05), with Tn lower than Tv and Tt lower than both Tv and Tn. Similarities in circadian rhythm patterns were observed across the day between Tv, Tn and Tt in repeated measures analysis, with a short period of difference between Tv and Tn (from 1400 to 1600 h) and a long period of difference between Tv and Tt (from 1000 to 2100 h) (P < 0.05). These results suggest that neck muscle temperature measurements may have utility in detecting circadian rhythm patterns in core temperature in sheep, but may not accurately reflect absolute core temperatures. Peripheral measures may require adjustment or correction to more accurately reflect absolute core temperature with respect to determining accurate clinical thresholds relative to the expected normal temperature for the time of day observed. Further investigation into the utility and application of peripheral measurement of body temperature is warranted.

Latent Variable Poisson Models for Assessing the Regularity of Circadian Patterns over Time

Many researchers in biology and medicine have focused on trying to understand biological rhythms and their potential impact on disease. A common biological rhythm is circadian, where the cycle repeats itself every 24 hours. However, a disturbance of the circadian pattern may be indicative of future disease. In this article, we develop new statistical methodology for assessing the degree of disturbance or irregularity in a circadian pattern for count sequences that are observed over time in a population of individuals. We develop a latent variable Poisson modeling approach with both circadian and stochastic short-term trend (autoregressive latent process) components that allow for individual variation in the degree of each component. A parameterization is proposed for modeling covariate dependence on the proportion of these two model components across individuals. In addition, we incorporate covariate dependence in the overall mean, the magnitude of the trend, and the phase-shift of the circadian pattern. Innovative Markov chain Monte Carlo sampling is used to carry out Bayesian posterior computation. Several variations of the proposed models are considered and compared using the deviance information criterion. We illustrate this methodology with longitudinal physical activity count data measured in a longitudinal cohort of adolescents.

Morning surge in blood pressure using a random-effects multiple-component cosinor model

Blood pressure (BP) fluctuates throughout the day. The pattern it follows represents one of the most important circadian rhythms in the human body. For example, morning BP surge has been suggested as a potential risk factor for cardiovascular events occurring in the morning, but the accurate quantification of this phenomenon remains a challenge. Here, we outline a novel method to quantify morning surge. We demonstrate how the most commonly used method to model 24-hour BP, the single cosinor approach, can be extended to a multiple-component cosinor random-effects model. We outline how this model can be used to obtain a measure of morning BP surge by obtaining derivatives of the model fit. The model is compared with a functional principal component analysis that determines the main components of variability in the data. Data from the Mitchelstown Study, a population-based study of Irish adults (n = 2047), were used where a subsample (1207) underwent 24-hour ambulatory blood pressure monitoring. We demonstrate that our 2-component model provided a significant improvement in fit compared with a single model and a similar fit to a more complex model captured by b-splines using functional principal component analysis. The estimate of the average maximum slope was 2.857 mmHg/30 min (bootstrap estimates; 95% CI: 2.855-2.858 mmHg/30 min). Simulation results allowed us to quantify the between-individual SD in maximum slopes, which was 1.02 mmHg/30 min. By obtaining derivatives we have demonstrated a novel approach to quantify morning BP surge and its variation between individuals. This is the first demonstration of cosinor approach to obtain a measure of morning surge.

Exploring diurnal variation using piecewise linear splines: an example using blood pressure

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 mixed-effects model with a constraint to ensure periodicity as a novel application for modelling daily blood pressure. Data from the Mitchelstown Study, a cross-sectional study of Irish adults aged 47–73 years (n = 2047) was utilized. A subsample (1207) underwent 24-h 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 short-term BP variability as we can quantify the variation about an individual’s trajectory but also allows examination of the variation in slopes between individuals (random-effects).

Electronic supplementary material

The online version of this article (doi:10.1186/s12982-017-0055-5) contains supplementary material, which is available to authorized users.

Modeling diurnal hormone profiles by hierarchical state space models

Adrenocorticotropic hormone (ACTH) diurnal patterns contain both smooth circadian rhythms and pulsatile activities. How to evaluate and compare them between different groups is a challenging statistical task. In particular, we are interested in testing (1) whether the smooth ACTH circadian rhythms in chronic fatigue syndrome and fibromyalgia patients differ from those in healthy controls and (2) whether the patterns of pulsatile activities are different. In this paper, a hierarchical state space model is proposed to extract these signals from noisy observations. The smooth circadian rhythms shared by a group of subjects are modeled by periodic smoothing splines. The subject level pulsatile activities are modeled by autoregressive processes. A functional random effect is adopted at the pair level to account for the matched pair design. Parameters are estimated by maximizing the marginal likelihood. Signals are extracted as posterior means. Computationally efficient Kalman filter algorithms are adopted for implementation. Application of the proposed model reveals that the smooth circadian rhythms are similar in the two groups but the pulsatile activities in patients are weaker than those in the healthy controls.

Joint analysis of longitudinal and survival data measured on nested timescales by using shared parameter models: an application to fecundity data

We consider the joint modelling, analysis and prediction of a longitudinal binary process and a discrete time-to-event outcome. We consider data from a prospective pregnancy study, which provides day level information regarding the behaviour of couples attempting to conceive. Reproductive epidemiologists are particularly interested in developing a model for individualized predictions of time to pregnancy (TTP). A couple's intercourse behaviour should be an integral part of such a model and is one of the main focuses of the paper. In our motivating data, the intercourse observations are a long series of binary data with a periodic probability of success and the amount of available intercourse data is a function of both the menstrual cycle length and TTP. Moreover, these variables are dependent and observed on different, and nested, timescales (TTP is measured in menstrual cycles whereas intercourse is measured on days within a menstrual cycle) further complicating its analysis. Here, we propose a semiparametric shared parameter model for the joint modelling of the binary longitudinal data (intercourse behaviour) and the discrete survival outcome (TTP). Further, we develop couple-based dynamic predictions for the intercourse profiles, which in turn are used to assess the risk for subfertility (i.e. TTP longer than six menstrual cycles).

Estimating long-term multivariate progression from short-term data

MOTIVATION

Diseases that progress slowly are often studied by observing cohorts at different stages of disease for short periods of time. The Alzheimer's Disease Neuroimaging Initiative (ADNI) follows elders with various degrees of cognitive impairment, from normal to impaired. The study includes a rich panel of novel cognitive tests, biomarkers, and brain images collected every 6 months for as long as 6 years. The relative timing of the observations with respect to disease pathology is unknown. We propose a general semiparametric model and iterative estimation procedure to estimate simultaneously the pathological timing and long-term growth curves. The resulting estimates of long-term progression are fine-tuned using cognitive trajectories derived from the long-term "Personnes Agées Quid" study.

RESULTS

We demonstrate with simulations that the method can recover long-term disease trends from short-term observations. The method also estimates temporal ordering of individuals with respect to disease pathology, providing subject-specific prognostic estimates of the time until onset of symptoms. When the method is applied to ADNI data, the estimated growth curves are in general agreement with prevailing theories of the Alzheimer's disease cascade. Other data sets with common outcome measures can be combined using the proposed algorithm.

AVAILABILITY

Software to fit the model and reproduce results with the statistical software R is available as the grace package. ADNI data can be downloaded from the Laboratory of NeuroImaging.

Direct regression models for longitudinal rates of change

Comparing rates of growth, or rates of change, across covariate-defined subgroups is a primary objective for many longitudinal studies. In the special case of a linear trend over time, the interaction between a covariate and time will characterize differences in longitudinal rates of change. However, in the presence of a non-linear longitudinal trajectory, the standard mean regression approach does not permit parsimonious description or inference regarding differences in rates of change. Therefore, we propose regression methodology for longitudinal data that allows a direct, structured comparison of rates across subgroups even in the presence of a non-linear trend over time. Our basic longitudinal rate regression method assumes a proportional difference across covariate groups in the rate of change across time, but this assumption can be relaxed. Rates are compared relative to a generally specified time trend for which we discuss both parametric and non-parametric estimating approaches. We develop mixed model longitudinal methodology that explicitly characterizes subject-to-subject variation in rates, as well as a marginal estimating equation-based method. In addition, we detail a score test to detect violations of the proportionality assumption, and we allow time-varying rate effects as a natural generalization. Simulation results demonstrate potential gains in power for the longitudinal rate regression model relative to a linear mixed effects model in the presence of a non-linear trend in time. We apply our method to a study of growth among infants born to HIV infected mothers and conclude with a discussion of possible extensions for our methods.

Modeling the salivary cortisol profile in population research: the multi-ethnic study of atherosclerosis

In many studies, it has been hypothesized that stress and its biologic consequences may contribute to disparities in rates of cardiovascular disease. However, understanding of the most appropriate statistical methods to analyze biologic markers of stress, such as salivary cortisol, remains limited. The authors explore the utility of various statistical methods in modeling daily cortisol profiles in population-based studies. They demonstrate that the proposed methods allow additional insight into the cortisol profile compared with commonly used summaries of the profiles based on raw data. For instance, one can gain insights regarding the shape of the population average curve, characterize the types of individual-level departures from the average curve, and better understand the relation between covariates and attained cortisol levels or slopes at various points of the day, in addition to drawing inferences regarding common features of the cortisol profile, such as the cortisol awakening response and the area under the curve. The authors compare the inference and interpretations drawn from these methods and use data collected as part of the Multi-Ethnic Study of Atherosclerosis to illustrate them.

Summer activity patterns among teenage girls: harmonic shape invariant modeling to estimate circadian cycles

BACKGROUND

Physical activity as measured by activity counts over short time intervals across a 24 h period are often used to assess circadian variation. We are interested in characterizing circadian patterns in activity among adolescents and examining how these patterns vary by obesity status. New statistical approaches are needed to examine how factors affect different features of the circadian pattern and to make appropriate covariate adjustments when the outcomes are longitudinal count data.

METHODS

We develop a statistical model for longitudinal or repeated activity count data that is used to examine differences in the overall activity level, amplitude (defined as the difference between the lowest and highest activity level over a 24 hour period), and phase shift. Using seven days of continuous activity monitoring, we characterize the circadian patterns and compare them between obese and non-obese adolescent girls.

RESULTS

We find a statistically significant phase delay in adolescent girls who were obese compared with their non-obese counterparts. After the appropriate adjustment for measured potential confounders, we did not find differences in mean activity level between the two groups.

CONCLUSION

New statistical methodology was developed to identify a phase delay in obese compared with non-obese adolescents. This new approach for analyzing longitudinal circadian rhythm count data provides a useful statistical technique to add to the repertoire for those analyzing circadian rhythm data.

Salivary cortisol in ambulatory assessment--some dos, some don'ts, and some open questions

The impact of stress on health and disease is an important research topic in psychosomatic medicine. Because research on hypothalamic-pituitary-adrenal (HPA) axis regulation under controlled laboratory studies lacks ecological validity, it needs to be complemented by a research program that includes momentary ambulatory assessment. The measurement of salivary cortisol offers the possibility to trace the free steroid hormone concentrations in ambulant settings. Therefore, in this article, we first discuss the role of salivary cortisol in ambulatory monitoring. We start with a brief description of HPA axis regulation, and we then consider cortisol assessments in other organic materials, followed by a presentation of common salivary markers of HPA axis regulation suitable for ambulatory assessment. We further provide an overview on assessment designs and sources of variability within and between subjects (intervening variables), acknowledge the issue of (non)compliance, and address statistical aspects. We further give an overview of associations with psychosocial and health-related variables relevant for ambulatory assessment. Finally, we deal with preanalytical aspects of laboratory salivary cortisol analysis. The relative simplicity of salivary cortisol assessment protocols may lead to an overoptimistic view of the robustness of this method. We thus discuss several important issues related to the collection and storage of saliva samples and present empirical data on the stability of salivary cortisol measurements over time.

Lack of fit in self modeling regression: application to pulse waveforms

Self modeling regression (SEMOR) is an approach for modeling sets of observed curves that have a common shape (or sequence of features) but have variability in the amplitude (y-axis) and/or timing (x-axis) of the features across curves. SEMOR assumes the x and y axes for each observed curve can be separately transformed in a parametric manner so that the features across curves are aligned with the common shape, usually represented by non-parametric function. We show that when the common shape is modeled with a regression spline and the transformational parameters are modeled as random with the traditional distribution (normal with mean zero), the SEMOR model may surprisingly suffer from lack of fit and the variance components may be over-estimated. A random effects distribution that restricts the predicted random transformational parameters to have mean zero or the inclusion of a fixed transformational parameter improves estimation. Our work is motivated by arterial pulse pressure waveform data where one of the variance components is a novel measure of short-term variability in blood pressure.

Analyzing biological rhythms in clinical trials

BACKGROUND

The human body exhibits a variety of biological rhythms. There are patterns that correspond, among others, to the daily wake / sleep cycle, a yearly seasonal cycle and, in women, the menstrual cycle. Sine/cosine functions are often used to model biological patterns for continuous data, but this model is not appropriate for analysis of biological rhythms in failure time data.

METHODS

We consider a method appropriate for analysis of biological rhythms in clinical trials. We present a method to provide an estimate and confidence interval of the time when the minimum hazard is achieved. A motivating example from a clinical trial of adjuvant of pre-menopausal breast cancer patients provides an important illustration of the methodology in practice.

RESULTS

Adapting the Cosinor method to the Weibull proportional hazards model is proposed as useful way of modeling the biological rhythm data. It presents a method to estimate the time that achieves the minimum hazard along with its associated confidence interval. The application of this technique to the breast cancer data revealed that the optimal day for pre-resection incisional or excisional biopsy of 28-day cycle (i.e. the day associated with the lowest recurrence rate) is day 8 with 95% CI 5-10. We found that older age, fewer positive nodes, smaller tumor size, and experimental treatment are important prognostic factors of longer relapse-free survival.

CONCLUSIONS

The analysis of biological/circadian rhythms is usually handled by Cosinor rhythmometry method. However, in FTD this is simply not possible. In this case, we propose to adapt the Cosinor method to the Weibull proportional hazard model. The advantage of the proposed method is its ability to model survival data. This method is not limited to breast cancer data, and may be applied to any biological rhythms linked to right censored data.

Procedures for numerical analysis of circadian rhythms

This article reviews various procedures used in the analysis of circadian rhythms at the populational, organismal, cellular and molecular levels. The procedures range from visual inspection of time plots and actograms to several mathematical methods of time series analysis. Computational steps are described in some detail, and additional bibliographic resources and computer programs are listed.

Modeling of hormone secretion-generating mechanisms with splines: a pseudo-likelihood approach

A flexible and robust approach is proposed for the investigation of underlying hormone secretion-generating mechanisms. Characterizing hormone time series is a difficult task as most hormones are secreted in a pulsatile manner and pulses are often masked by slow decay. We model hormone concentration as a filtered counting process where the intensity function of the counting process is modeled nonparametrically using periodic splines. The intensity function and parameters are estimated using a combination of weighted least squares and pseudo-likelihood based on the first two moments. Our method uses concentration measurements directly, which avoids the difficult task of estimating pulse numbers and locations. Both simulations and applications suggest that our method performs well for estimating the intensity function of the pulse-generating counting processes.

Infant growth modelling using a shape invariant model with random effects

Models for infant growth have usually been based on parametric forms, commonly an exponential or similar model, which have been shown to fit poorly especially during the first year of life. An alternative approach is to use a non-parametric model, based on a shape invariant model (SIM), where a single function is transformed by shifting and scaling to fit each subject. In the model a regression spline is used as the function, with log transformation of the data and a simplification of the SIM, obtained from the relationship with the exponential model. All subjects are fitted as a nonlinear mixed effects model, allowing the variation in the parameters between subjects to be determined. Methods for the inclusion of covariates in growth models based on SIM are developed, with parameters for time independent covariates included in the model by varying either the shape, the size parameter or the growth parameter and time-dependent co-variates included by transforming the time axis, to either increase or decrease the growth rate dependent on the co-variate, similar to methods used for accelerated failure-time models. The model is used to fit weight data for 602 infants, measured from 0 to 2 years as part of the Childhood Asthma Prevention Study (CAPS) trial, and to determine the effect of breastfeeding on infant weight.

Absence of circadian salivary cortisol rhythm in women with anorexia nervosa

STUDY OBJECTIVE

To compare the cortisol levels and 24 hour salivary cortisol rhythm in patients with anorexia nervosa (AN) and normal controls.

DESIGN

Prospective transversal controlled study.

SETTING

Tertiary-referral University Hospital.

PARTICIPANTS

Twenty-five patients aged 15 to 35 years, 13 of them with regular ovulatory cycles, and 12 with diagnosis of AN.

INTERVENTIONS

Salivary and blood collection for cortisol 24-hour rhythm determination.

MAIN OUTCOME

Salivary cortisol was determined at 9 am, 5 pm, and 11 pm. Seric follicle-stimulating hormone, luteinizing hormone (LH), prolactin, estradiol (E2), progesterone, dehydroepiandrosterone-S (DHEA-S), and cortisol were sampled together with the 9 am salivary sample.

RESULTS

LH, E2, and DHEA-S levels were reduced in patients with AN. A correlation between salivary and serum cortisol levels was observed in the 9 am sample only in controls (r = 0.67, P = 0.01; AN: r = 0.48, P = 0.12). Cortisol rhythm was present in all control subjects, whereas it was absent in one third of AN patients. The area under the curve for the AN group with preserved rhythm was significantly higher than for the control group (Me = 6811 ng/dl/24h vs 3708 ng/dl/24 h; P = 0.034).

CONCLUSION

Patients with AN have higher salivary cortisol levels when compared to normal women and some of them do not present circadian rhythm.

Modeling biological rhythms in failure time data

Background

The human body exhibits a variety of biological rhythms. There are patterns that correspond, among others, to the daily wake/sleep cycle, a yearly seasonal cycle and, in women, the menstrual cycle. Sine/cosine functions are often used to model biological patterns for continuous data, but this model is not appropriate for analysis of biological rhythms in failure time data.

Methods

We adapt the cosinor method to the proportional hazards model and present a method to provide an estimate and confidence interval of the time when the minimum hazard is achieved. We then apply this model to data taken from a clinical trial of adjuvant of pre-menopausal breast cancer patients.

Results

The application of this technique to the breast cancer data revealed that the optimal day for pre-resection incisional or excisional biopsy of 28-day cycle (i. e. the day associated with the lowest recurrence rate) is day 8 with 95% confidence interval of 4–12 days. We found that older age, fewer positive nodes, smaller tumor size, and experimental treatment were predictive of longer relapse-free survival.

Conclusion

In this paper we have described a method for modeling failure time data with an underlying biological rhythm. The advantage of adapting a cosinor model to proportional hazards model is its ability to model right censored data. We have presented a method to provide an estimate and confidence interval of the day in the menstrual cycle where the minimum hazard is achieved. This method is not limited to breast cancer data, and may be applied to any biological rhythms linked to right censored data.

Comparing tongue shapes from ultrasound imaging using smoothing spline analysis of variance

Ultrasound imaging of the tongue is increasingly common in speech production research. However, there has been little standardization regarding the quantification and statistical analysis of ultrasound data. In linguistic studies, researchers may want to determine whether the tongue shape for an articulation under two different conditions (e.g., consonants in word-final versus word-medial position) is the same or different. This paper demonstrates how the smoothing spline ANOVA (SS ANOVA) can be applied to the comparison of tongue curves [Gu, Smoothing Spline ANOVA Models (Springer, New York, 2002)]. The SS ANOVA is a technique for determining whether or not there are significant differences between the smoothing splines that are the best fits for two data sets being compared. If the interaction term of the SS ANOVA model is statistically significant, then the groups have different shapes. Since the interaction may be significant even if only a small section of the curves are different (i.e., the tongue root is the same, but the tip of one group is raised), Bayesian confidence intervals are used to determine which sections of the curves are statistically different. SS ANOVAs are illustrated with some data comparing obstruents produced in word-final and word-medial coda position.

On analyzing circadian rhythms data using nonlinear mixed models with harmonic terms

Wang, Ke, and Brown (2003, Biometrics59, 804-812) developed a smoothing-based approach for modeling circadian rhythms with random effects. Their approach is flexible in that fixed and random covariates can affect both the amplitude and phase shift of a nonparametrically smoothed periodic function. In motivating their approach, Wang et al. stated that a simple sinusoidal function is too restrictive. In addition, they stated that "although adding harmonics can improve the fit, it is difficult to decide how many harmonics to include in the model, and the results are difficult to interpret." We disagree with the notion that harmonic models cannot be a useful tool in modeling longitudinal circadian rhythm data. In this note, we show how nonlinear mixed models with harmonic terms allow for a simple and flexible alternative to Wang et al.'s approach. We show how to choose the number of harmonics using penalized likelihood to flexibly model circadian rhythms and to estimate the effect of covariates on the rhythms. We fit harmonic models to the cortisol circadian rhythm data presented by Wang et al. to illustrate our approach. Furthermore, we evaluate the properties of our procedure with a small simulation study. The proposed parametric approach provides an alternative to Wang et al.'s semiparametric approach and has the added advantage of being easy to implement in most statistical software packages.

Spline-based non-parametric regression for periodic functions and its application to directional tuning of neurons

The activity of neurons in the brain often varies systematically with some quantitative feature of a stimulus or action. A well-known example is the tendency of the firing rates of neurons in the primary motor cortex to vary with the direction of a subject's arm or wrist movement. When this movement is constrained to vary in only two dimensions, the direction of movement may be characterized by an angle, and the neuronal firing rate can be written as a function of this angle. The firing rate function has traditionally been fit with a cosine, but recent evidence suggests that departures from cosine tuning occur frequently. We report here a new non-parametric regression method for fitting periodic functions and demonstrate its application to the fitting of neuronal data. The method is an extension of Bayesian adaptive regression splines (BARS) and applies both to normal and non-normal data, including Poisson data, which commonly arise in neuronal applications. We compare the new method to a periodic version of smoothing splines and some parametric alternatives and find the new method to be especially valuable when the smoothness of the periodic function varies unevenly across its domain.