3
$\begingroup$

I'm currently facing a problem which I thought simpler in the beginning. Basically I have outputs of wind direction from a model and observed wind direction from my sensors (both in degrees). What I would like to do is the equivalent of linear regression and correlation coefficient computation that is used when comparing linear data.

I have both MATLAB and R available and here's something that I've attempted:

  • I am aware of the existence of the Circular Statistics Toolbox for MATLAB, but, while there are functions for correlations there's nothing for regression.

  • I've tried to use the package circular for R, but, for both regression and correlation data needs to be "circular objects" and I do not know how to correctly set the function as.circular for wind direction data: I understand that the unit flag must be set to degree and the rotation to "clock"...but I'm not entirely sure about the zero and modulo flag since the function reasons in radians. Would converting the degrees to radians work? While their modulo will became $2\pi$ what about the zero flag?

$\endgroup$
2
$\begingroup$

You may use the circular package. It's probably best to convert your data to radians by multiplying by $\pi/180$, or using xcirc <- rad(x).

The type, unit, zero, modulo flags are mostly there for plotting and storing information about a dataset. For your analysis, they don't matter.

Here is some example analysis, adapted from example(lm.circular)

# Example data
x <- runif(n, 0, 2*pi)
y <- atan2(0.15*cos(x) + 0.25*sin(x), 0.35*sin(x)) + rnorm(n, 0, 1)

# Compute the model, get a circular correlation.
circ.lm <- lm.circular(y, x, order=1, type = "c-c")
cor.circular(y=y, x=x)

# Plot the results.
plot.default(x, y)
circ.lm$fitted[circ.lm$fitted>pi] <- circ.lm$fitted[circ.lm$fitted>pi] - 2*pi 
points.default(x[order(x)], circ.lm$fitted[order(x)], type='l')
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.