“显着”不同的问题总是以数据的统计模型为前提。 该答案提出了与问题中提供的最少信息一致的最通用模型之一。简而言之,它适用于各种情况,但它可能并不总是检测差异的最有效方法。
数据的三个方面真正重要:点所占据空间的形状;该空间内的点分布;以及由具有“条件”的点对形成的图形——我将其称为“治疗”组。“图表”是指治疗组中点对所暗示的点和互连模式。例如,图的十个点对(“边”)可能涉及多达 20 个不同的点或少至五个点。在前一种情况下,没有两条边共享一个公共点,而在后一种情况下,边由五个点之间的所有可能对组成。
为了确定治疗组中边缘之间的平均距离是否“显着”,我们可以考虑一个随机过程,其中所有 n = 3000点由排列随机排列σ. 这也会置换边缘:边缘(v一世,vj)被替换为(vσ(一),vσ( j )). 零假设是边缘的治疗组作为其中之一出现3000 !≈1021024排列。如果是这样,它的平均距离应该与这些排列中出现的平均距离相当。我们可以通过对所有这些排列中的几千个进行抽样来相当容易地估计这些随机平均距离的分布。
(值得注意的是,这种方法只需稍作修改即可适用于与每个可能的点对相关的任何距离或任何数量。它也适用于距离的任何总结,而不仅仅是平均值。)
为了说明,这里有两种情况涉及n = 100点和28治疗组的边缘。在第一行中,每条边的第一个点是从100点,然后每条边的第二个点是独立且随机地从100−1点与他们的第一点不同。全部一起39这些点涉及到28边缘。
在最下面一行,八个100点是随机选择的。这28边由所有可能的对组成。
右侧的直方图显示了样本分布10000配置的随机排列。数据的实际平均距离用垂直的红色虚线标记。两种方法都与采样分布一致:既不位于右侧也不位于左侧。
采样分布不同:尽管平均距离相同,但在第二种情况下,由于边缘之间的图形相互依赖性,平均距离的变化更大。 这是不能使用中心极限定理的简单版本的原因之一:计算该分布的标准差很困难。
以下是与问题中描述的数据相当的结果:n=3000点大致均匀分布在一个正方形内,并且1500他们的对在治疗组。 计算只用了几秒钟,证明了它们的实用性。
顶行中的对再次被随机选择。在底行中,处理组中的所有边仅使用56离左下角最近的点。它们的平均距离远小于采样分布,因此可以认为具有统计学意义。
一般来说,模拟和治疗组的平均距离等于或大于治疗组平均距离的比例可以作为这个非参数置换检验的p值。
这是R
用于创建插图的代码。
n.vectors <- 3000
n.condition <- 1500
d <- 2 # Dimension of the space
n.sim <- 1e4 # Number of iterations
set.seed(17)
par(mfrow=c(2, 2))
#
# Construct a dataset like the actual one.
#
# `m` indexes the pairs of vectors with a "condition."
# `x` contains the coordinates of all vectors.
x <- matrix(runif(d*n.vectors), nrow=d)
x <- x[, order(x[1, ]+x[2, ])]
#
# Create two kinds of conditions and analyze each.
#
for (independent in c(TRUE, FALSE)) {
if (independent) {
i <- sample.int(n.vectors, n.condition)
j <- sample.int(n.vectors-1, n.condition)
j <- (i + j - 1) %% n.condition + 1
m <- cbind(i,j)
} else {
u <- floor(sqrt(2*n.condition))
v <- ceiling(2*n.condition/u)
m <- as.matrix(expand.grid(1:u, 1:v))
m <- m[m[,1] < m[,2], ]
}
#
# Plot the configuration.
#
plot(t(x), pch=19, cex=0.5, col="Gray", asp=1, bty="n",
main="The Data", xlab="X", ylab="Y",
sub=paste(length(unique(as.vector(m))), "points"))
invisible(apply(m, 1, function(i) lines(t(x[, i]), col="#80000040")))
points(t(x[, unique(as.vector(m))]), pch=16, col="Red", cex=0.6)
#
# Precompute all distances between all points.
#
distances <- sapply(1:n.vectors, function(i) sqrt(colSums((x-x[,i])^2)))
#
# Compute the mean distance in any set of pairs.
#
mean.distance <- function(m, distances)
mean(distances[m])
#
# Sample from the points using the same *pattern* in the "condition."
# `m` is a two-column array pairing indexes between 1 and `n` inclusive.
sample.graph <- function(m, n) {
n.permuted <- sample.int(n, n)
cbind(n.permuted[m[,1]], n.permuted[m[,2]])
}
#
# Simulate the sampling distribution of mean distances for randomly chosen
# subsets of a specified size.
#
system.time(
sim <- replicate(n.sim, mean.distance(sample.graph(m, n.vectors), distances))
stat <- mean.distance(m, distances)
p.value <- 2 * min(mean(c(sim, stat) <= stat), mean(c(sim, stat) >= stat))
hist(sim, freq=FALSE,
sub=paste("p-value:", signif(p.value, ceiling(log10(length(sim))/2)+1)),
main="Histogram of mean distances", xlab="Distance")
abline(v = stat, lwd=2, lty=3, col="Red")
}