编程语言R中计算两个向量之间角度的最有效方法是什么？
问问题
43223 次
9 回答
47
根据本 PDF的第 5 页，是求向量和sum(a*b)的点积的R 命令，是求向量范数的 R 命令，是反余弦的 R 命令。由此可见，计算两个向量之间夹角的 R 代码是absqrt(sum(a * a))aacos(x)
theta <- acos( sum(a*b) / ( sqrt(sum(a * a)) * sqrt(sum(b * b)) ) )
于 2009-12-13T22:26:55.880 回答
24
我的回答由两部分组成。第 1 部分是数学 - 使线程的所有读者都清楚，并使后面的 R 代码易于理解。第 2 部分是 R 编程。
第 1 部分 - 数学
两个向量x和y的点积可以定义为：
在哪里 || x || 是向量x的欧几里得范数（也称为 L 2范数） 。
操作点积的定义，我们可以得到：
其中 theta 是以弧度表示的向量x和y之间的角度。请注意，theta 可以取一个位于从 0 到 pi 的闭合区间上的值。
求解 theta 本身，我们得到：
第 2 部分 - R 代码
要将数学转换为 R 代码，我们需要知道如何执行两个矩阵（向量）计算；点积和欧几里得范数（这是一种特定类型的范数，称为 L 2范数）。我们还需要知道反余弦函数 cos -1的 R 等效项。
从顶部开始。参考?"%*%"，可以使用算子计算点积（也称为内积）%*%。参考?norm，norm()函数（基本包）返回向量的范数。这里感兴趣的范数是 L 2范数，或者用 R 帮助文档的说法是“光谱”或“2”范数。这意味着函数的type参数norm()应该设置为等于"2"。最后，R中的反余弦函数由acos()函数表示。
解决方案
配备了数学和相关的 R 函数，原型函数（即非生产标准）可以放在一起 - 使用 Base 包函数 - 如下所示。如果上述信息有意义，那么后面的angle()功能应该很清楚，无需进一步评论。
angle <- function(x,y){
dot.prod <- x%*%y
norm.x <- norm(x,type="2")
norm.y <- norm(y,type="2")
theta <- acos(dot.prod / (norm.x * norm.y))
as.numeric(theta)
}
测试功能
验证该功能是否有效的测试。设x = (2,1) 和y = (1,2)。x和y之间的点积是 4。x 的欧几里得范数是sqrt(5)。y的欧几里得范数也是 sqrt(5)。因为θ = 4/5。Theta 约为 0.643 弧度。
x <- as.matrix(c(2,1))
y <- as.matrix(c(1,2))
angle(t(x),y) # Use of transpose to make vectors (matrices) conformable.
[1] 0.6435011
我希望这有帮助！
于 2014-07-28T16:29:41.590 回答
15
对于 2D 向量，在接受的答案和其他答案中给出的方式不考虑角度的方向（符号）（angle(M,N)与 相同angle(N,M)），它仅对 和 之间的角度返回正确0值pi。
使用该atan2函数获取定向角度和正确值（模数2pi）。
angle <- function(M,N){
acos( sum(M*N) / ( sqrt(sum(M*M)) * sqrt(sum(N*N)) ) )
}
angle2 <- function(M,N){
atan2(N[2],N[1]) - atan2(M[2],M[1])
}
检查是否angle2给出正确的值：
> theta <- seq(-2*pi, 2*pi, length.out=10)
> O <- c(1,0)
> test1 <- sapply(theta, function(theta) angle(M=O, N=c(cos(theta),sin(theta))))
> all.equal(test1 %% (2*pi), theta %% (2*pi))
[1] "Mean relative difference: 1"
> test2 <- sapply(theta, function(theta) angle2(M=O, N=c(cos(theta),sin(theta))))
> all.equal(test2 %% (2*pi), theta %% (2*pi))
[1] TRUE
于 2015-07-26T22:01:42.900 回答
7
您应该使用点积。假设你有V ₁ = ( x ₁, y ₁, z ₁) 和V ₂ = ( x ₂, y ₂, z ₂)，然后是点积，我将用V ₁·<em>V₂ 表示，计算为
V ₁·<em>V₂ = x ₁·<em>x₂ + y ₁·<em>y₂ + z ₁ ·<em>z₂ = | ₁ | · | Ⅴ2 | · cos( θ );
这意味着左边显示的和等于向量的绝对值乘以向量之间夹角的余弦的乘积。矢量V 1 和V 2的绝对值计算为
| ₁ | = √( x ₁² + y ₁² + z ₁²)，和
| Ⅴ2 | = √( x ²² + y ²² + z ²²),
所以，如果你重新排列上面的第一个方程，你会得到
cos( θ ) = ( x ₁·<em>x₂ + y ₁·<em>y₂ + z ₁ ·<em>z₂) ÷ (| V ₁|·| V ₂|)，
您只需要将 arccos 函数（或反余弦）应用于 cos( θ ) 即可获得角度。
根据您的 arccos 函数，角度可能以度或弧度为单位。
（对于二维向量，只需忘记z坐标并进行相同的计算。）
祝你好运，
约翰·多纳
于 2009-12-13T21:40:24.687 回答
6
另一种解决方案：两个向量之间的相关性等于两个向量之间夹角的余弦。
所以角度可以通过acos(cor(u,v))
# example u(1,2,0) ; v(0,2,1)
cor(c(1,2),c(2,1))
theta = acos(cor(c(1,2),c(2,1)))
于 2013-08-20T08:51:10.890 回答
5
如果您安装/上传库（matlib）：有一个名为 angle(x, y, degree = TRUE) 的函数，其中 x 和 y 是向量。注意：如果您有矩阵形式的 x 和 y，请使用 as.vector(x) 和 as.vector(y)：
library(matlib)
matA <- matrix(c(3, 1), nrow = 2) ##column vectors
matB <- matrix(c(5, 5), nrow = 2)
angle(as.vector(matA), as.vector(matB))
##default in degrees, use degree = FALSE for radians
于 2019-11-05T02:23:46.420 回答
1
我认为你需要的是内在产品。对于两个向量v,u（在R^n或任何其他内积空间中）<v,u>/|v||u|= cos(alpha)。（alpha是向量之间的角度）
有关更多详细信息，请参见：
于 2009-12-13T21:35:57.997 回答
1
如果要计算多个变量之间的角度，可以使用以下函数，它是@Graeme Walsh 提供的解决方案的扩展。
angles <- function(matrix){
## Calculate the cross-product of the matrix
cross.product <- t(matrix)%*%matrix
## the lower and the upper triangle of the cross-product is the dot products among vectors
dot.products<- cross.product[lower.tri(cross.product)]
## Calculate the L2 norms
temp <- suppressWarnings(diag(sqrt(cross.product)))
temp <- temp%*%t(temp)
L2.norms <- temp[lower.tri(temp)]
## Arccosine values for each pair of variables
lower.t <- acos(dot.products/L2.norms)
## Create an empty matrix to present the results
result.matrix <- matrix(NA,ncol = dim(matrix)[2],nrow=dim(matrix)[2])
## Fill the matrix with arccosine values and assign the diagonal values as zero “0”
result.matrix[lower.tri(result.matrix)] <- lower.t
diag(result.matrix) <- 0
result.matrix[upper.tri(result.matrix)] <- t(result.matrix)[upper.tri(t(result.matrix))]
## Get the result matrix
return(result.matrix)
}
此外，如果您以输入变量为中心并获得上面提供的结果矩阵的余弦值，您将获得变量的精确相关矩阵。
这是该功能的应用。
set.seed(123)
n <- 100
m <- 5
# Generate a set of random variables
mt <- matrix(rnorm(n*m),nrow = n,ncol = m)
# Mean-centered matrix
mt.c <- scale(mt,scale = F)
# Cosine angles
cosine.angles <- angles(matrix = mt)
> cosine.angles
[,1] [,2] [,3] [,4] [,5]
[1,] 0.000000 1.630819 1.686037 1.618119 1.751859
[2,] 1.630819 0.000000 1.554695 1.523353 1.712214
[3,] 1.686037 1.554695 0.000000 1.619723 1.581786
[4,] 1.618119 1.523353 1.619723 0.000000 1.593681
[5,] 1.751859 1.712214 1.581786 1.593681 0.000000
# Centered-data cosine angles
centered.cosine.angles <- angles(matrix = mt.c)
> centered.cosine.angles
[,1] [,2] [,3] [,4] [,5]
[1,] 0.000000 1.620349 1.700334 1.614890 1.764721
[2,] 1.620349 0.000000 1.540213 1.526950 1.701793
[3,] 1.700334 1.540213 0.000000 1.615677 1.595647
[4,] 1.614890 1.526950 1.615677 0.000000 1.590057
[5,] 1.764721 1.701793 1.595647 1.590057 0.000000
# This will give you correlation matrix
cos(angles(matrix = mt.c))
[,1] [,2] [,3] [,4] [,5]
[1,] 1.00000000 -0.04953215 -0.12917601 -0.04407900 -0.19271110
[2,] -0.04953215 1.00000000 0.03057903 0.04383271 -0.13062219
[3,] -0.12917601 0.03057903 1.00000000 -0.04486571 -0.02484838
[4,] -0.04407900 0.04383271 -0.04486571 1.00000000 -0.01925986
[5,] -0.19271110 -0.13062219 -0.02484838 -0.01925986 1.00000000
# Orginal correlation matrix
cor(mt)
[,1] [,2] [,3] [,4] [,5]
[1,] 1.00000000 -0.04953215 -0.12917601 -0.04407900 -0.19271110
[2,] -0.04953215 1.00000000 0.03057903 0.04383271 -0.13062219
[3,] -0.12917601 0.03057903 1.00000000 -0.04486571 -0.02484838
[4,] -0.04407900 0.04383271 -0.04486571 1.00000000 -0.01925986
[5,] -0.19271110 -0.13062219 -0.02484838 -0.01925986 1.00000000
# Check whether they are equal
all.equal(cos(angles(matrix = mt.c)),cor(mt))
[1] TRUE
于 2019-12-30T02:58:30.563 回答
0
获得两个向量之间的角度的传统方法（即acos(sum(a*b) / (sqrt(sum(a*a)) * sqrt(sum(b*b))))，如其他一些答案中所述）在几个极端情况下存在数值不稳定性。以下代码适用于n维和所有极端情况（它不检查零长度向量，但这很容易添加）。请参阅下面的注释。
# Get angle between two n-dimensional vectors
angle_btw <- function(v1, v2) {
signbit <- function(x) {
x < 0
}
u1 <- v1 / norm(v1, "2")
u2 <- v2 / norm(v2, "2")
y <- u1 - u2
x <- u1 + u2
a0 <- 2 * atan(norm(y, "2") / norm(x, "2"))
if (!(signbit(a0) || signbit(pi - a0))) {
a <- a0
} else if (signbit(a0)) {
a <- 0.0
} else {
a <- pi
}
a
}
此代码改编自 Jeffrey Sarnoff 的Julia 实现（MIT 许可），又基于W. Kahan 教授的这些注释（第 15 页）。
于 2022-02-03T00:39:12.183 回答