Applied Example with Bison Movement Data

This vignette illustrates a compact workflow using the bison data included in CircularRegression. The response y.dir is the movement direction at the current time. The covariate y.prec is the previous movement direction, and x.meadow:z.meadow represents the direction toward the nearest meadow weighted by the associated non-negative distance variable.

library(CircularRegression)

data(bison)

d <- bison[seq_len(200), c("y.dir", "y.prec", "x.meadow", "z.meadow")]
d <- na.omit(d)
str(d)
#> Classes 'tbl_df', 'tbl' and 'data.frame':    182 obs. of  4 variables:
#>  $ y.dir   : num  -2.27 -2.27 -2.27 -2.27 -2.26 ...
#>  $ y.prec  : num  -2.27 -2.27 -2.27 -2.27 -2.27 ...
#>  $ x.meadow: num  -2.22 -2.22 -2.22 -2.22 -2.22 ...
#>  $ z.meadow: num  600 590 580 570 560 ...
#>  - attr(*, "na.action")= 'omit' Named int [1:18] 1 2 66 67 74 75 99 100 121 122 ...
#>   ..- attr(*, "names")= chr [1:18] "1" "2" "66" "67" ...

Fit the default two-step model

fit <- circular_regression(
  y.dir ~ y.prec + x.meadow:z.meadow,
  data = d
)

fit
#> Call:
#> circular_regression(formula = y.dir ~ y.prec + x.meadow:z.meadow, 
#>     data = d)
#> 
#> Method: two_step 
#> Number of observations: 182 
#> 
#> Call:
#> angular(formula = formula, data = data, reference = c("name", 
#>     ref_name), initbeta = initbeta, control = control_angular, 
#>     na.action = na.action)
#> 
#> Reference angle: y.prec 
#> Maximum cosine: 0.9813 
#> Concentration parameter: 27.05 
#> 
#> Parameters (non-reference terms):
#>                     Estimate Std. Error z value P(|z|>.)   
#> x.meadow:z.meadow 0.00051640 0.00017127  3.0151 0.002569 **
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The default fit uses the consensus model to select a reference direction and initialize the homogeneous angular fit.

summary(fit)
#> Circular regression summary
#> Method: two_step 
#> Number of observations: 182 
#> 
#> Call:
#> angular(formula = formula, data = data, reference = c("name", 
#>     ref_name), initbeta = initbeta, control = control_angular, 
#>     na.action = na.action)
#> 
#> Reference angle: y.prec 
#> Maximum cosine: 0.9813 
#> Concentration parameter: 27.05 
#> 
#> Residuals:
#>       Min.    1st Qu.     Median       Mean    3rd Qu.       Max. 
#> -8.859e-01 -6.634e-02 -1.000e-08 -6.543e-03  7.021e-02  8.253e-01 
#> 
#> Parameters:
#>                     estimate     stderr    z value P(|z|>.)
#> x.meadow:z.meadow 0.00051640 0.00017127 3.01506345   0.0026
#> 
#> Number of observations: 182
coef(fit)
#> [1] 0.0005163958

Predictions and residuals

Predicted values are angles in radians. Residuals are wrapped circular differences between the observed and fitted directions.

pred <- predict(fit)
res <- residuals(fit)

head(data.frame(observed = d$y.dir, fitted = pred, residual = res))
#>    observed    fitted     residual
#> 1 -2.268730 -2.257728 -0.011002591
#> 2 -2.268787 -2.257641 -0.011146246
#> 3 -2.268723 -2.257640 -0.011082586
#> 4 -2.265709 -2.257558 -0.008150492
#> 5 -2.264139 -2.255194 -0.008945436
#> 6 -2.264147 -2.253935 -0.010211883
summary(res)
#>       Min.    1st Qu.     Median       Mean    3rd Qu.       Max. 
#> -8.859e-01 -6.634e-02 -1.000e-08 -6.543e-03  7.021e-02  8.253e-01

A simple base-R diagnostic plot can be useful before calling the package plot method.

plot(
  pred,
  res,
  xlab = "Fitted direction",
  ylab = "Circular residual",
  pch = 19,
  col = "gray30"
)
abline(h = 0, lty = 2)

Direct consensus fit

The consensus model can be fitted directly when the interest is in the consensus coefficients.

fit_cons <- consensus(
  y.dir ~ y.prec + x.meadow:z.meadow,
  data = d
)

coef(fit_cons)
#>            y.prec x.meadow:z.meadow 
#>       2.115901258       0.004961799
head(predict(fit_cons))
#> [1] -2.241449 -2.241157 -2.240906 -2.240630 -2.239088 -2.238146

This applied example is intentionally small so that it can be built quickly as part of package checks. Larger analyses should keep the same structure: inspect the angular variables, fit the model, examine circular residuals, and verify that conclusions are stable under reasonable changes to the formula.