1

我正在尝试按照本出版物中描述的程序来确定组之间的分离是否具有统计学意义。该出版物的水平高于我的知识水平，但我正在尝试逐步接近它。

为了澄清和简单起见，以鸢尾花数据集为例，并在 R 中进行分析。正如 PCA 图所示，该方法应该使我能够确定组/物种之间的距离是否显着不同。
Iris Species PCA 图

据我了解，要获得此结果，该过程包括以下四个步骤：

  1. 距离计算：使用来自前两个主成分的组质心之间的马氏距离。
  2. 学生 t 检验：采用 Hotelling 的双样本 T^2 统计量来确定是否分离 btw。聚类具有统计学意义。
  3. 计算 F 统计量：将 T^2 统计量转换为 F 值并计算 F 检验以指示集群之间是否存在分离。
  4. 使用 F 统计量执行假设检验：如果 F 值大于临界 F 值，则可以拒绝假设组之间没有分离的原假设。

我被困在第一步和第二步之间。如何使用 Mahalanobis 距离计算的结果进行 Hotelling 的 T^2 检验。

MWE如下：

library(ICSNP)
library(ggbiplot)
data(iris)

# Mahalanobis Distance calculation Function from https://stackoverflow.com/a/34708113/5731401
D.sq <- function (g1, g2) {
    dbar <- as.vector(colMeans(g1) - colMeans(g2))
    S1 <- cov(g1)
    S2 <- cov(g2)
    n1 <- nrow(g1)
    n2 <- nrow(g2)
    V <- as.matrix((1/(n1 + n2 - 2)) * (((n1 - 1) * S1) + ((n2 - 1) * S2)))
    D.sq <- t(dbar) %*% solve(V) %*% dbar
    res <- list()
    res$D.sq <- D.sq
    res$V <- V
    res
}

iris.pca <- prcomp(iris[,-5], center = TRUE, scale. = TRUE)
str(iris)
# uncomment the next line for illustrative plot
# print(ggbiplot(iris.pca, obs.scale = 1, var.scale = 1, groups = iris$Species, ellipse = TRUE, circle = TRUE))
df.iris.x <- as.data.frame(iris.pca$x)
df.iris.x$Species <- iris$Species

split.data = split(df.iris.x[,-5],df.iris.x$Species)
S1 = split.data[['setosa']]
S2 = split.data[['versicolor']]
S3 = split.data[['virginica']]

# calculate mahalanobis distances for the first two principal components between the groups/species
d1 <- D.sq(S1[,1:2],S2[,1:2])
d2 <- D.sq(S1[,1:2],S3[,1:2])
d3 <- D.sq(S2[,1:2],S3[,1:2])

# T-test on the first two principal components 
HotellingsT2(S1[,1:2], S2[,1:2]) #btw setosa and versicolor
HotellingsT2(S1[,1:2], S3[,1:2]) #btw setosa and virginica
HotellingsT2(S2[,1:2], S3[,1:2]) #btw versicolor and virginica

马氏距离计算返回三个距离

d1$D.sq = 70.7239
d2$D.sq = 97.53259
d3$D.sq = 5.910554

就比例而言，它们似乎与 PCA 图所示的相当。对前两个主成分的Hotelling T2 检验对所有三个比较都返回显着结果。但我想知道如何使用参考出版物中所述的先前马氏距离计算进行 T 检验？

4

0 回答 0