使用dplyr和broom在训练集和测试集上计算kmeans

Question

我正在使用dplyr和broom计算我的数据的kmeans。我的数据包含一组测试和训练集的X和Y坐标，并且根据某个参数值（在这种情况下为lambda）进行分组：

mds.test = data.frame()
for(l in seq(0.1, 0.9, by=0.2)) {
  new.dist <- run.distance.model(x, y, lambda=l)
  mds <- preform.mds(new.dist, ndim=2)
  mds.test <- rbind(mds.test, cbind(mds$space, design[,c(1,3,4,5)], lambda=rep(l, nrow(mds$space)), data="test"))
}

> head(mds.test)
                        Comp1       Comp2 Transcripts Genes Timepoint Run lambda data
7A_0_AAGCCTAGCGAC -0.06690476 -0.25519106       68125  9324     Day 0  7A    0.1 test
7A_0_AAATGACTGGCC -0.15292848  0.04310200       28443  6746     Day 0  7A    0.1 test
7A_0_CATCTCGTTCTA -0.12529445  0.13022908       27360  6318     Day 0  7A    0.1 test
7A_0_ACCGGCACATTC -0.33015913  0.14647857       23038  5709     Day 0  7A    0.1 test
7A_0_TATGTCGGAATG -0.25826098  0.05424976       22414  5878     Day 0  7A    0.1 test
7A_0_GAAAAAGGTGAT -0.24349387  0.08071162       21907  6766     Day 0  7A    0.1 test

我已经准备好了以上所述的测试数据集，但我还有一个名为mds.train的数据集，其中包含我的训练数据坐标。我的最终目标是按lambda分组对两个数据集运行k-means，然后计算测试数据在训练中心处的内部平方和（within.ss）、外部平方和（between.ss）和总平方和（total.ss）。感谢broom的一个很棒的资源，我可以通过简单地执行以下操作来为每个lambda在测试集上运行kmeans：

test.kclusts  = mds.test %>% 
  group_by(lambda) %>% 
  do(kclust=kmeans(cbind(.$Comp1, .$Comp2), centers=length(unique(design$Timepoint))))

那么我可以计算每个 $\lambda$ 中每个聚类的数据中心：

test.clusters = test.kclusts %>% 
  group_by(lambda) %>%  
  do(tidy(.$kclust[[1]])) 

我现在遇到了瓶颈。如何像参考页面上所示计算特征分配（例如 kclusts %>% group_by(k) %>% do(augment(.$kclust[[1]], points.matrix))），其中我的 points.matrix 是 mds.test，它是一个数据框，行数应该是 length(unique(mds.test$lambda)) 倍。是否有一种方法可以使用来自训练集的中心来计算基于测试分配的glance()统计信息？

非常感谢任何帮助！谢谢！

编辑：更新进展。我已经找出了如何聚合测试/训练分配，但仍然无法尝试从两个集合（基于测试中心的训练分配和基于训练中心的测试分配）计算k均值统计信息。下面是更新后的代码：

test.kclusts  = mds.test %>% group_by(lambda) %>% do(kclust=kmeans(cbind(.$Comp1, .$Comp2), centers=length(unique(design$Timepoint))))
test.clusters = test.kclusts %>% group_by(lambda) %>%  do(tidy(.$kclust[[1]])) 
test.clusterings = test.kclusts %>% group_by(lambda) %>% do(glance(.$kclust[[1]]))
test.assignments = left_join(test.kclusts, mds.test) %>% group_by(lambda) %>% do(augment(.$kclust[[1]], cbind(.$Comp1, .$Comp2)))

train.kclusts  = mds.train %>% group_by(lambda) %>% do(kclust=kmeans(cbind(.$Comp1, .$Comp2), centers=length(unique(design$Timepoint))))
train.clusters = train.kclusts %>% group_by(lambda) %>%  do(tidy(.$kclust[[1]])) 
train.clusterings = train.kclusts %>% group_by(lambda) %>% do(glance(.$kclust[[1]]))
train.assignments = left_join(train.kclusts, mds.train) %>% group_by(lambda) %>% do(augment(.$kclust[[1]], cbind(.$Comp1, .$Comp2)))

test.assignments$data = "test"
train.assignments$data = "train"
merge.assignments = rbind(test.assignments, train.assignments)
merge.assignments %>% filter(., data=='test') %>% group_by(lambda) ... ? 

我在下面附上了一张图表，展示了我到目前为止的进展。再次强调，我想计算kmeans统计量（组内平方和、总平方和和组间平方和）应用于测试分配/坐标的训练数据中心（中心看起来不准确的那些图表）：

- user2117258

---

看起来有一些关于制作tidymodels配方或模型规范的讨论，例如recipes::step_kmeans()（请参见此处的对话以及从那里的链接：https://github.com/tidymodels/embed/issues/77#issuecomment-801300749）。但是尚未完成。还有widyr::widely_kmeans()，但这也没有一个易于应用于测试集的方法。 - Bryan Shalloway

Answer 1

一个方法是通过broom提取您的聚类表（建立在训练集上）中指定的质心。

1. 使用fuzzyjoin包，计算测试集中每个点与使用训练集构建的每个聚类质心之间的距离。
2. 测试点与其最短欧几里得距离的聚类质心代表其分配的聚类。
3. 从那里，您可以计算任何感兴趣的度量。

请参见下面使用从tidymodels中提取的聚类example的简化数据集。

library(tidyverse)
library(rsample)
library(broom)
library(fuzzyjoin)

# data and train / test set-up
set.seed(27)
centers <- tibble(
  cluster = factor(1:3), 
  num_points = c(100, 150, 50),  # number points in each cluster
  x1 = c(5, 0, -3),              # x1 coordinate of cluster center
  x2 = c(-1, 1, -2)              # x2 coordinate of cluster center
)

labelled_points <- 
  centers %>%
  mutate(
    x1 = map2(num_points, x1, rnorm),
    x2 = map2(num_points, x2, rnorm)
  ) %>% 
  select(-num_points) %>% 
  unnest(cols = c(x1, x2))

points <- 
  labelled_points %>% 
  select(-cluster)

set.seed(1234)

split <- rsample::initial_split(points)
train <- rsample::training(split)
test <- rsample::testing(split)

# Fit kmeans on train then assign clusters to test
kclust <- kmeans(train, centers = 3)

clust_centers <- kclust %>% 
  tidy() %>% 
  select(-c(size, withinss))

test_clusts <- fuzzyjoin::distance_join(mutate(test, index = row_number()), 
                         clust_centers,
                         max_dist = Inf,
                         method = "euclidean",
                         distance_col = "dist") %>% 
  group_by(index) %>% 
  filter(dist == min(dist)) %>% 
  ungroup()
#> Joining by: c("x1", "x2")

# resulting table
test_clusts
#> # A tibble: 75 x 7
#>     x1.x    x2.x index  x1.y  x2.y cluster  dist
#>    <dbl>   <dbl> <int> <dbl> <dbl> <fct>   <dbl>
#>  1  4.24 -0.946      1  5.07 -1.10 3       0.847
#>  2  3.54  0.287      2  5.07 -1.10 3       2.06 
#>  3  3.71 -1.67       3  5.07 -1.10 3       1.47 
#>  4  5.03 -0.788      4  5.07 -1.10 3       0.317
#>  5  6.57 -2.49       5  5.07 -1.10 3       2.04 
#>  6  4.97  0.233      6  5.07 -1.10 3       1.34 
#>  7  4.43 -1.89       7  5.07 -1.10 3       1.01 
#>  8  5.34 -0.0705     8  5.07 -1.10 3       1.07 
#>  9  4.60  0.196      9  5.07 -1.10 3       1.38 
#> 10  5.68 -1.55      10  5.07 -1.10 3       0.758
#> # ... with 65 more rows

# calc within clusts SS on test
test_clusts %>% 
  group_by(cluster) %>% 
  summarise(size = n(),
            withinss = sum(dist^2),
            withinss_avg = withinss / size)
#> # A tibble: 3 x 4
#>   cluster  size withinss withinss_avg
#>   <fct>   <int>    <dbl>        <dbl>
#> 1 1          11     32.7         2.97
#> 2 2          35     78.9         2.26
#> 3 3          29     62.0         2.14

# compare to on train
tidy(kclust) %>% 
  mutate(withinss_avg = withinss / size)
#> # A tibble: 3 x 6
#>        x1    x2  size withinss cluster withinss_avg
#>     <dbl> <dbl> <int>    <dbl> <fct>          <dbl>
#> 1 -3.22   -1.91    40     76.8 1               1.92
#> 2  0.0993  1.06   113    220.  2               1.95
#> 3  5.07   -1.10    72    182.  3               2.53

# plot of test and train points
test_clusts %>% 
  select(x1 = x1.x, x2 = x2.x, cluster) %>% 
  mutate(type = "test") %>% 
  bind_rows(
    augment(kclust, train) %>% 
      mutate(type = "train") %>% 
      rename(cluster = .cluster)
    ) %>% 
  ggplot(aes(x = x1, 
             y = x2, 
             color = as.factor(cluster)))+
  geom_point()+
  facet_wrap(~fct_rev(as.factor(type)))+
  coord_fixed()+
  labs(title = "Cluster Assignment on Training and Holdout Datasets",
       color = "Cluster")+
  theme_bw()

^{2021年08月19日由reprex包 (v2.0.0)创建}

(有关在tidymodels中更轻松地进行此操作的对话链接，请参见OP上的评论。)

- Bryan Shalloway