I am trying to model circular data that has a slight shift every cycle. The aim is to model actogram data, the daily activity of animals. In these experiments one of the most obvious phenotypes is the shift in the periodic activity. I should probably treat this as a time series, rather than doing this with circular distributions, but I don’t see why this shouldn’t work .
Fitting circular data with the von Mises distribution worked fine but when I added a shift in each cycle, things started to fall apart. Here is the simulated data I used.
library(circular)
library(tidyverse)
library(rstan)
options(mc.cores =parallel::detectCores())
# simulate circular data with shift
##############################################################################################
n_cycle <- 8
df <- tibble(cycle=1:n_cycle) %>%
mutate(y=map(cycle, ~rvonmises(100,
mu=circular(pi-2*pi/18*.x),
kappa=2*pi))) %>%
unnest(y)
ggplot(df, aes(x=y)) + geom_histogram() + facet_grid(cycle~.)
model_data <- list(
n = nrow(df),
n_cycle = n_cycle,
cycle = df$cycle,
y = df$y
)
model_name <- "von mises with shift.stan"
model <- stan(
file = model_name,
data = model_data, iter=2000, chains=10, warmup=1500
)
summary(model)
and here is the Stan model
data{
int<lower=1> n;
int<lower=1> n_cycle;
int<lower=1, upper=n_cycle> cycle[n];
real<lower=0, upper=2*pi()> y[n];
}
parameters{
real<lower=0, upper=2*pi()> mu;
real<lower=0> kappa;
// shift should be between -pi and pi
real<lower=-pi()/4, upper=pi()/4> shift;
}
model{
shift ~ normal(0,pi()/12);
for (i in 1:n){
if (kappa < 100)
y[i] ~ von_mises(mu+shift*(cycle[i]-1), kappa);
else
y[i] ~ normal(mu+shift*(cycle[i]-1), sqrt(1 / kappa));
}
}
The problem is that some chains seem to converge on the proper values, but others don’t. There seems to be a trade off between kappa and the amount of shift. I thought that giving the shift a prior that favours smaller values would fix this, but even a very strong prior on the shift does not work. The pairs plot look like polka dots, some chains seem to find the right answer, and others converge to a very small value of kappa. In the von Moses distribution Kappa approaches 0 as the distribution flattens (kappa =0 is uniform distribution). I thought maybe there wasn’t enough data points, but increasing the amount of data doesn’t change things, and fitting a similar model without the shift with the same amount of data seems to work just fine.
Am I going on about the wrong way? I know that what I am really trying to estimate is the period of the circular data, rather than mu and kappa, but to me this should be equivalent to estimating the amount of “shift”, assuming this shift is constant.
I feel like I am missing something, probably some embarrassing mistake, but can’t figure out what, and I need to know! Any help or insight would be much appreciated.
