I am working on a phenology study of equatorial plants. I have a predictor variable (Year of observation) which is non-circular. And a response variable (day of observation) which has to be considered circular, since the first of January comes right after December 31.

The goal is to fit a model on the data in order to predict a LS regression that can “cross” the arbitrary boundary of the new year. My approach was to use the predict function of the fitted model but this does not seems to be compatible with lm.circular from the circular package.

I used the following code:

install.packages("circular")
library(circular)

observation_Year <- c(1983, 1965, 1968, 1968, 1970, 1970, 1970, 1970, 1970, 1967, 1973, 1977, 1983, 1983, 1986, 1988, 1994, 1994, 1994, 1994, 1994, 1998, 1997, 1998, 1998, 2001, 2005, 2005, 2005, 2008, 2011, 2011, 2017, 2021)
day_of_observation <- c(228, 166, 78, 83, 110, 92, 109, 147, 197, 44, 86, 238, 228, 228, 43, 75, 213, 212, 214, 217, 231, 19, 12, 20, 19, 25, 140, 141, 141, 52, 17, 19, 67, 210)

day_of_observation_rad <- day_of_observation * 2 * pi / 365
day_of_observation_rad <- circular(day_of_observation_rad)

#fit circular model
model <- lm.circular(y = day_of_observation_rad,
                     x = observation_Year,
                     type = "c-l",
                     init = 0)

#create predictions
Plants <- data.frame(
  StartDate = seq(from = min(observation_Year),
                       to = max(observation_Year)))


Plants$DayOfYearPredRad <- predict(model, newdata = Plants)

I assumed that predict was compatible with lm.circular which is not the case.
I am looking for viable alternatives on order to predict a model on this data.