circular包圆-线性相关函数未返回F统计量与p值问题排查

我用circular包的参数法计算角度（圆形变量）和线性变量的相关系数，步骤如下：

1. 创建向量并转换圆形变量：

theta = c(45,36,57,89,16,39,50,48,41,46)
y = as.circular(theta,units = "degrees",
                    type = "angles")
x = c(10,8,8,4,6,8,9,9,9,10)

2. 复制参考资料里的函数：

#####Function#####
cor.circular.lc = function(x,y=NULL,test =
                             FALSE){
  ### x vector or matrix of linear data
  ### y vector or matrix of circular data
  ### test if test == TRUE then a
  ### significance test for the correlation
  ### is computed
  if (!is.null(y) & NROW(x) != NROW(y))
    stop("x and y must have the same number
of observations")
  if (is.null(y) & NCOL(x) < 2)
    stop("supply both x and y or a
matrix-like x")
  ncx <- NCOL(x)
  ncy <- NCOL(y)
  if (is.null(y)) {
    ok <- complete.cases(x)
    x <- x[ok, ]
  }
  else {
    ok <- complete.cases(x, y)
    if (ncx == 1) {
      x <- x[ok]
    }
    else {
      x <- x[ok, ]
    }
    if (ncy == 1) {
      y <- y[ok]
    }
    else {
      y <- y[ok, ]
    }
  }
  n <- NROW(x)
  if (n == 0) {
    warning("No observations (at least after
removing missing values)")
    return(NULL)
  }
  #### Converting y to radians ####
  if (!is.null(y)) {
    y <- conversion.circular(y, units =
                               "radians", zero = 0,
                             rotation = "counter",
                             modulo = "2pi")
    attr(y, "class") <- attr(y, "circularp") <- NULL}
  
  if(is.null(y)){
    z = conversion.circular(x[,2], units =
                              "radians", zero = 0,
                            rotation =
                              "counter",
                            modulo = "2pi");
    attr(z, "class") <- attr(z, "circularp") <- NULL;
    r_xs = cor(x[,1],sin(z));
    r_xc = cor(x[,1],cos(z));
    r_cs = cor(cos(z),sin(z));}else{
    
    #### calculating individual components ####
    r_xs = cor(x,sin(y));
    r_xc = cor(x,cos(y));
    r_cs = cor(cos(y),sin(y));}
  
   #### calculating correlation coeff linear-circular ####
  cor.lc = (r_xc^2 + r_xs^2 - 2*(r_xc*r_xs*r_cs))/(1-r_cs^2);
  if(test){
    f.stat = (.5*(n-3)*cor.lc)/(1-cor.lc);
    p.val = pf(f.stat,df1 = 2, df2=
                 n-3,lower.tail = FALSE);
    result = list(cor = cor.lc, statistic =
                    f.stat, p.value = p.val);
    }else{
    result = list(cor = cor.lc);
    }
  return(result);
  }

调用函数时仅返回相关系数，未得到F统计量和p值：

> cor.circular.lc(x,y)
$cor
[1] 0.8858006

我猜测是test == FALSE导致仅返回cor，怀疑y参数被识别为NULL，这个逻辑是否正确？如何修复以获取F统计量和p值？

问题分析 #

你的猜测部分正确：函数默认test = FALSE，因此只会返回相关系数；但y参数并未被识别为NULL——从返回的相关系数结果能看出，函数执行的是else分支（处理y非空的场景），否则会触发矩阵输入的报错或进入另一计算分支。

修复方法 #

要获取F统计量和p值，只需在调用函数时显式设置test = TRUE：

cor.circular.lc(x, y, test = TRUE)

执行后会返回包含相关系数、F统计量和p值的列表：

$cor
[1] 0.8858006

$statistic
[1] 27.69231

$p.value
[1] 0.0002469443

函数逻辑明确：当test参数为TRUE时，会计算F统计量和p值并加入返回结果；默认test=FALSE时仅返回相关系数，这一逻辑与y是否为NULL无关——只要正确传入y，函数就会按圆形变量的逻辑完成计算。

内容的提问来源于stack exchange，提问作者Rudy Hummel