如何计算 R 中一系列样本中两个线向量之间的角度列表?

How to calculate a list of angles between two line vectors in a series of specimens in R?

我正在研究具有二维笛卡尔坐标(即地标)的生物标本的几何形态数据集。我对这个数据集的部分兴趣是生成一个数据框,其中包含重要地标之间的距离和地标向量之间的角度,所有样本都具有同源坐标,以便进一步比较。点之间的距离不是问题,但是当我尝试计算每个样本的角度时,我 运行 遇到了问题。

我正在尝试找到一种方法来一次计算代表一系列标本上相同角度的线向量之间的角度(这些向量描述了一个具有生物学意义的特征,该特征存在于我的数据集中的所有标本中)。我用来操作形态测量数据的 R 包将地标输出为 class array() 的三维对象,其中行是地标,列是坐标的维度,第三维是标本.这是此处提供的样本数据中三个标本的地标坐标。

, , Specimen1

         X         Y
1  -0.2411  0.060183
2  -0.0677 -0.029954
3   0.1147  0.012111
4   0.0085  0.000462

, , Specimen2

          X       Y
1  -0.22509  0.0764
2  -0.09437 -0.0202
3   0.09135 -0.0182
4   0.00367 -0.0045

, , Specimen3

         X        Y
1  -0.2223  0.06122
2  -0.1001 -0.02366
3   0.0553  0.00577
4  -0.0617 -0.01557

下面是一个可读数组,它复制了我作为前三个样本的输出获得的数据。

data<-array(data=c(-0.2411,-0.0677,0.1147,0.0085,0.060183,-0.029954,0.012111,0.000462,-0.22509,-0.09437,0.09135,0.00367,0.0764,-0.0202,-0.0182,-0.0045,-0.2223,-0.1001,0.0553,-0.0617,0.06122,-0.02366,0.00577,-0.01557),c(4,2,3))

我想找出点 1 和点 2 与点 3 和点 4 之间的线形成的角度。为此,我通过从一个点的坐标中减去另一个点的坐标来计算这两条线的向量。当我这样做时,我为每个向量得到一个 2*N 矩阵,其中行对应于 X 和 Y 轴上的向量,列对应于样本。

我想找到样本 1(两个矩阵的第 1 列)、2、3 等的这两个向量之间的角度。理想情况下,我希望将输出转换为 R 报告每个单独样本的角度的格式,如下所示:

     Angle1
[,1] 146
[,2] 152
[,3] 135

然后以某种方式将其添加到数据框中,如下所示:

          Distance Angle1 Angle2
Specimen1 100      146    100
Specimen2 100      152    100
Specimen3 100      135    100

Angle2 的距离和值只是填充数字,以显示我正在尝试实现的格式。

我尝试使用 mathlib 包中的 angle() 函数和 Morpho 包中的 angle.calc() 函数] 包裹。这是我使用 angle() 函数计算角度的代码。

vector1<-data[1,,]-data[2,,]
vector2<-data[3,,]-data[4,,]
angle(as.vector(vector1),as.vector(vector2),degree=TRUE)

这些函数将计算出看似准确的角度,但如果另有说明,它们只会对第一个样本或单个样本进行计算,并且不会生成所有样本的角度列表。我还在 Stack Overflow 上检查了 previous ,但给出的答案中 none 似乎适用于包含多个样本的数据集,两个矩阵上的每一列代表一个不同的样本(当我最好地应用它们时他们似乎有同样的问题,他们只有 return 第一个样本的角度)。这里只有三个样本,我在实际数据集中有大量样本和每个样本的多个角度,不可能为每个样本单独测量每个角度。

有什么方法可以编写公式来生成包含每个样本各自角度的列表吗?

首先,将 2x3 矩阵命名为 vector1 有点令人困惑。除此之外,您可以使用 matlib::angle 来计算角度。

library(matlib)
phi <- setNames(mapply(
    function(x1, x2) angle(x1, x2),
    as.data.frame(vector1), as.data.frame(vector2)),
    paste0("Specimen", 1:ncol(vector1)))
phi
#Specimen1 Specimen2 Specimen3
# 146.2674  152.4439  134.8753

解释:使用mapply同时遍历vector1vector2的列,已转换为data.frame。使用 stack 将命名向量转换为 data.frame:

stack(phi)
#    values       ind
#1 146.2674 Specimen1
#2 152.4439 Specimen2
#3 134.8753 Specimen3

一般来说,您可能会发现首先以 long/tidy 格式获取数据更简单。这允许更容易的分组和聚合,这是 vectorized/mapped 操作真正发挥作用的地方。

这是使用您提供的数据执行此操作的一种方法:

library(matlib)
library(tidyverse)

n_specimen <- ncol(vector1)

data.frame(t(vector1)) %>%
  rename_all(~str_replace(., "X", "v1_")) %>%
  bind_cols(
    data.frame(t(vector2)) %>%
      rename_all(~str_replace(., "X", "v2_"))
    ) %>%
  mutate(specimen = paste0("Specimen", 1:n_specimen)) %>%
  pivot_longer(-specimen, names_to = "vector", values_to = "value") %>%
  separate(vector, into = c("vector", "ix"), sep = "_", remove = TRUE) %>%
  group_by(specimen) %>%
  summarise(angle = angle(value[vector == "v1"], value[vector == "v2"])) 

  specimen  angle
  <chr>     <dbl>
1 Specimen1  146.
2 Specimen2  152.
3 Specimen3  135.

作为解释,下面是数据框在 separate 步骤后的样子:

# A tibble: 12 x 4
   specimen  vector ix      value
   <chr>     <chr>  <chr>   <dbl>
 1 Specimen1 v1     1     -0.173 
 2 Specimen1 v1     2      0.0901
 3 Specimen1 v2     1      0.106 
 4 Specimen1 v2     2      0.0116
 5 Specimen2 v1     1     -0.131 
 6 Specimen2 v1     2      0.0966
 7 Specimen2 v2     1      0.0877
 8 Specimen2 v2     2     -0.0137
 9 Specimen3 v1     1     -0.122 
10 Specimen3 v1     2      0.0849
11 Specimen3 v2     1      0.117 
12 Specimen3 v2     2      0.0213

每一列都有独特的信息级别,您现在可以分组并汇总到所需的具体级别。