Functional Linear Modelling (FLM)

Import required libraries:

import circstudio as cs
import plotly.graph_objs as go
import os
import numpy as np
import plotly.io as pio
pio.templates.default = "plotly_white"
pio.renderers.default = "notebook"

Open a sample file:

# Create file path for sample files within the circStudio package
fpath = os.path.join(os.path.dirname(cs.__file__), 'data')

raw = cs.io.read_awd(
    os.path.join(fpath, 'example_01.AWD'), 
    start_time='1918-01-24 08:00:00',
    period='9 days'
)

Overview

Actigraphy time series are often noisy: brief spikes, non-wear artifacts, and rapid transitions can obscure the underlying 24-hour structure. The FLM class in circStudio provides two complementary smoothing workflows:

1. Basis function expansion: build an average 24-hour profile and represent it as a smooth function using either a Fourier basis or B-splines.

2. Gaussian kernel smoothing - using either the average daily profile or the full time series.

Basis function expansion

The goal of functional linear modelling (FLM) is to represent a noisy activity signal as a smooth function constructed from simpler building blocks called basis functions.

Formally, we represent the observed signal x(t_i) as a weighted sum of basis functions: $$x(t_i) \approx \beta_1 \phi_1(t_i) + \beta_2 \phi_2(t_i) + … + \beta_n \phi_n(t_i) $$

Here:

• φ_n(t) are known functions (the basis functions)

• β_k are coefficients estimated from the data

Intuitively, this means we are rebuilding the activity signal as a combination of smooth, structured components, instead of describing it point-by-point.

Different choices of basis functions lead to different types of smooth representations. In circStudio, you can find two:

• Fourier basis functions (periodic waves).

• B-splines (flexible smooth curves).

Using a Fourier basis expansion

For circadian rhythms, Fourier functions are a natural choice because they are inherently periodic. Instead of describing the activity signal point-by-point, we represent it as a weighted sum of smooth oscillatory components: $$x(t) = \beta_0 \phi_0(t) + \sum_{n=1}^{N}[\beta_{2n-1} \sin(n \omega t) + \beta_{2n}cos(n \omega t)]$$

Where:

• φ₀(t) = 1 (constant term)

• ω = 2π / T

• T is the period (for daily activity profiles, T = 24 hours).

• β_k are coefficients estimated from the data.

These basis functions arise from Fourier series theory, which states that any periodic signal with period T (e.g., actigraphy-derived time series) can be expressed as a sum of sine and cosine waves, with frequency nω = 2πn / T.

Geometric intuition behind sine and consine in Fourier basis expansion

A useful way to build intuition is to think of sine and cosine as coordinates on a unit circle (a circle of radius one). If we move around a unit circle at angular speed ω = 2π / T, cos(ωt) gives the horizontal coordinate, whereas sin(ωt) gives the vertical coordinate. These two functions are called the fundamental (first harmonic) components of the 24-hour rhythm. Together, they complete one full rotation every T hours. In other words, a 24-hour oscillation corresponds to one full rotation per day.

Higher harmonics sin(nωt), cos(nωt) simply rotate around the circle n times within the same 24-hour interval. Thus, these additional harmonics correspond to faster oscillatory components of the same daily rhythm, and they allow the model to capture wake-sleep transitions and secondary, rhythmic structure.

Applying Fourier FLM in circStudio

First, define the sampling frequency (freq) and the number of harmonics (max_order):

freq = '1min'
max_order = 9

Now instantiate and fit the model:

flm = cs.FLM(
    basis="fourier",
    sampling_freq=freq,
    max_order=max_order
)
flm.fit(raw.activity)

Finally, reconstruct the smooth 24-hour profile:

daily_curve = flm.evaluate()

The result is a smooth, periodic representation of the average daily activity pattern, capturing its circadian structure while filtering out high-frequency noise. You can visualize and compare it to the raw average daily profile:

# Compute average daily profile for activity
daily_avg = cs.daily_profile(raw.activity)

# Create figure layout
labels = ['00:00', '06:00', '12:00', '18:00']
tickvals = ['00:00:00', '06:00:00', '12:00:00', '18:00:00']

layout = go.Layout(
    autosize=False, width=900, height=600,
    title="Daily profile",
    xaxis=dict(
        title="Time of day (HH:MM)",
        ticktext=labels,
        tickvals=tickvals),
    yaxis=dict(title="Activity"),
    shapes=[], showlegend=True)
go.Figure(data=[
    go.Scatter(x=daily_avg.index.astype(str),y=daily_avg,name='Raw activity'),
    go.Scatter(x=daily_avg.index.astype(str),y=daily_curve,name='Fourier expansion (9th order)')
],layout=layout)

Using B-splines

While Fourier basis functions are ideal when we want to enforce strict 24-hour periodicity, B-splines provide a more flexible alternative. Instead of representing the daily rhythm as a combination of sine and cosine waves, splines approximate the signal using local smooth polynomial pieces that are stitched together.

Resample the data

Before fitting splines, we resample the activity signal to a 30-minute resolution and then compute the average daily profile:

# Resample activity trace to 30 min
raw.apply_filters(new_freq='30min')

# Compute average daily profile for resampled activity
daily_avg = cs.daily_profile(raw.activity)

Fit the spline model

We now define a FLM with a spline basis:

flm_spline = cs.FLM(
    basis='spline',
    sampling_freq='30min',
    max_order=3
)
flm_spline.fit(raw.activity)

Evaluate the smooth curve

The spline can be evaluated at a higher resolution than the original data using the parameter r:

r = 10
spline_est = flm_spline.evaluate(r=r)

The parameter r controls how finely we sample the smooth curve relative to the original number of points:

• If r = 1, the curve is evaluated at the same resolution as the original daily profile.

• If r = 10, the curve is evaluated at 10 times more points, producing a visually smoother curve.

Importantly, this does not change the fitted model, it only increases the resolution at which we visualize it.

Visualization

Because the evaluated spline contains r times more points, we need two different time axes:

t = np.linspace(0,daily_avg.index.size,daily_avg.index.size,endpoint=True)
t_10 = np.linspace(0,daily_avg.index.size,daily_avg.index.size*r,endpoint=True)

We can then plot both the original daily average and the spline estimate:

data = [go.Scatter(x=t,y=daily_avg,name='Raw activity'),
        go.Scatter(x=t_10,y=spline_est,name='B-spline')
       ]
go.Figure(data=data, layout=layout)

The result is a smooth curve that follows the general shape of the daily rhythm while reducing short-scale fluctuations.

Gaussian kernel smoothing

Basis expansions (Fourier or B-splines) smooth the signal by restricting it to a particular family of functions (periodic waves or piecewise polynomials). A complementary approach is kernel smoothing, which smooths the signal locally by averaging each time point with its neighbors.

In circStudio, Gaussian kernel smoothing is implemented by convolving the signal with a Gaussian kernel. Intuitively, this replaces each sample by a weighted average of nearby samples, where points closer in time receive higher weights.

Smoothing the average daily profile

We begin by computing the average 24-hour activity profile and then smoothing it using several kernel bandwidth choices.

raw.reset_filters()

daily_avg = cs.daily_profile(raw.activity)

names = ['Raw activity', 'Scott', 'Silverman', 'Bandwith: 20']
daily_avg_smoothed = []

daily_avg_smoothed.append(flm.smooth_daily_profile(raw.activity, method='scott', verbose=True))
daily_avg_smoothed.append(flm.smooth_daily_profile(raw.activity, method='silverman', verbose=True))
daily_avg_smoothed.append(flm.smooth_daily_profile(raw.activity, method=20))

Kernel size used to smooth the data: 39.16586876177629
Kernel size used to smooth the data: 33.28544087403084

The key parameter is the Gaussian kernel width. A small bandwidth corresponds to miniomal smoothing (more detail retained, but more noise remains). A large bandwidth, corresponds to stronger smoothing (noise reduced, but sharp features may be flattened). circStudio supports automatic bandwidth selection rules that are commonly used in smoothing (Scott and Silverman):

data = [go.Scatter(x=daily_avg.index.astype(str),y=daily_avg_smoothed[n], name=names[n+1]) for n in range(0,len(daily_avg_smoothed))]
data.insert(0,go.Scatter(x=daily_avg.index.astype(str),y=daily_avg,name=names[0]))
go.Figure(data=data,layout=layout)

Smoothing the full time series

Sometimes you do not want to smooth the average day, but rather the entire continuous recording (e.g., before fitting downstream models or visualizing long recordings). For this, use smooth_timeseries():

ts_smoothed = flm.smooth_timeseries(raw.activity, method='scott', verbose=True)

Kernel size used to smooth the data: 28.587203569652328

Visualize the result:

data = [go.Scatter(
    x=raw.activity.index.astype(str),
    y=ts_smoothed,
    name='scott')]

data.insert(0,go.Scatter(
    x=raw.activity.index.astype(str),
    y=raw.activity,
    name='Raw activity'))

go.Figure(data=data,layout=layout)