Rで画像データをUMAPで次元圧縮する際のパラメータによる描画の差異を可視化する
塩漬けにしていた画像データを軽い気持ちでUMAPを使ってみようと思ったら10箇所くらい引っかかるポイントがあったので経緯をメモしておく。最初に感想だけ書いておくと、PCAなどではほとんど描画側がいじれるパラメータがなく、あっても1-2軸で見せるのか、1-3軸で見せるのか、はたまた2-3軸（殆ど見たことないけど）で見せるのか、それとも3次元で見せるのかくらいしか差がないのに、UMAPは様々にいじれてしまうパラメータがある（上に、そのパラメータの設定について現時点ではあまり何も批判的に言われていないように思う）ため、さまざまなパラメータを試してから自分の思い通りに、恣意的に選べる気がして危うい。
libraries
もはやどれが必要でどれが必要でなかったかわからなくなったが
library(umap)
library(tidyverse)
library(snedata) # install via remotes::install_github("jlmelville/snedata")
library(imager) # for images
library(purrr)
library(Rtsne) # https://github.com/jkrijthe/Rtsne
UMAPの練習
まずはUMAPが走るかを確認する。
iris_umap <- iris |>
select(where(is.numeric)) |>
umap()
iris_df <- iris_umap$layout |>
as_tibble() |>
transmute(x = V1, y = V2) |>
mutate(labels = iris$Species)
iris_df |>
ggplot(aes(x = x, y = y, color = labels)) +
geom_point()
ggsave("./output/iris_umap.png", height = 1000, width = 1500, unit = "px")
同様に、MNISTについてここを参考にやる。しかしmnistデータをsnedataから落としたものを直接UMAPしてもエラーになるので、as.numeric()でintegerをdoubleにする。。。これが原因なのかよくわからないが動くのでヨシ
また、Label列はクラスタリング用のラベル、いわば「正解」なのでほんとうは抜いてから計算すべきだがサンプルコードではそのようになっていない。
umap()はデフォルトで並列化されているらしいが、ちゃんとn_threadsを指定すると早くなる気がする…。verboseはTRUEにしておくと実行時間もわかるし、どこまで計算が進んでいるかわかって安心す。。
mnist <- snedata::download_mnist()
mnist_umap <- umap(mnist, pca = 100) # error... :(
mnist_df <- mnist |>
as_tibble() |>
mutate(across(everything(), as.numeric)) |>
select(-Label) # ラベル列を落とす
mnist_umap <- mnist_df |>
umap(
pca = 10,
verbose = TRUE,
n_threads = 6
)
mnist_umap_df <-
mnist_umap$layout |>
as_tibble() |>
transmute(x = V1, y = V2) |>
mutate(label = mnist$Label)
mnist_umap_df |>
ggplot(aes(x = x, y = y, colour = label)) +
geom_point(alpha = .5, size = .001) +
theme_minimal()
ggsave("./output/mnist_umap.png", height = 1000, width = 1500, unit = "px")
前処理–画像データを読み込み、ひとつの巨大なデータフレームにする
手元に数十枚程度の正方形のjpgデータがあるとする。{imager}を使うのだが、ここでもチュートリアルに従ってsystem.file('[path]', package = "imager")を使うとうまくいかないので、直接放り込んでしまう。input/フォルダに000.jpg, 001.jpg, 002.jpg... 041.jpgと42枚の正方形の画像があったとしよう。
実際の画像の例
他にファイルがなければlist.files("input/", full.names=T) |> でもよいのだが、他にもファイルがある場合は次のように書く。
df_image <-
paste0(
"input/",
formatC(0:41, width = 3, flag = "0"),
".jpg"
) |> # すべてのパスをstringのベクトルとして書き下す。フォルダ内に他ファイルがなければ`list.files()`でも可
map(~ load.image(file = .x)) |> # imagerをpurrr::mapに渡して一気にリストとして読み込む
map(~ resize(im = .x, size_x = 128, size_y = 128)) |> # 受け取った画像のリストをリサイズする。800*800くらいあったので…。
map(~ grayscale(im = .x)) |> # もともとグレイスケールの画像だったがファイルとしては3チャンネルあったので削減
map(~t(as.numeric(.))) |> # ここでまた数値に変換しないといけない。なぜ
map_dfr(~ as.data.frame(.)) # このままではマトリックスなのでデータフレームに変換してからそれらを縦にひっつける(`map_dfr`)
128 * 128 = 16384 + tの1列で16385列に対して、1行1画像なので42枚・行の非常に横幅の大きいデータフレームができる。
UMAP
先程のMNISTと違い、こちらは行数がものすごく少ないので一瞬で終わる。データが小さいのでPCAは必要ないと判断した。
image_umap <- df_image |>
umap(verbose = TRUE, n_threads = 6)
プロット
image_umap_df |>
ggplot(aes(x = x, y = y, colour = t)) +
geom_point(size = .3) +
geom_path() +
geom_text(aes(label = t), hjust = 0, nudge_x = .1, size = 2) +
scale_colour_viridis_c() +
theme_minimal() +
theme(
panel.grid = element_blank(),
axis.line = element_line())
ggsave("./output/image_umap.png", height = 1000, width = 1500, unit = "px")
tのleakがないバージョン
tの列を残したままやってしまったバージョン。こんなに違うのか…
t-SNEとの比較
R wrapper for Van der Maaten’s Barnes-Hut implementation of t-Distributed Stochastic Neighbor Embeddingを参考に、あとt-SNEのクラスターサイズは何も意味しない、距離もあんまり意味ないかもあたりを肝に銘じつつt-SNEもやってみる。
Perplexityはデータサイズの1/3が上限(p21)なので、デフォルトの30では動かない。
image_tsne <- df_image |>
as.matrix() |>
Rtsne(
perplexity = 12, # default 30, max datasize/3
n_threads = 4,
verbose = TRUE
)
image_tsne_df <- image_tsne$Y |>
as_tibble() |>
rename(x = "V1", y = "V2") |>
mutate(t = 0:41)
image_tsne_df |>
ggplot(aes(x = x, y = y, colour = t)) +
geom_point(size = .3) +
geom_path() +
geom_text(aes(label = t), hjust = 0, nudge_x = .1, size = 2) +
scale_colour_viridis_c() +
theme_minimal() +
theme(
panel.grid = element_blank(),
axis.line = element_line())
ggsave("./output/image_tsne.png", height = 1000, width = 1500, unit = "px")
LeakageがないUMAPと同じくらいのめちゃくちゃだが、なんとなくこちらのほうがすっきりしてる…？
UMAPのパラメータによる描画の差異を可視化する
このスライドのp43にあるmin_distとn_neighborsを変えた場合の描画の差をみて、「同じデータでこんなに違う印象与えられるUMAP、ちょっと危なすぎでは？」と思って自分でもやってみることにした。ひとつひとつマニュアルで描画してもよいのだが、せっかくなのでmap()などを多用して一気にやってしまおうと思ったらmatrix-columnが生成されて困ったが解決した。
まずmin_distとn_neighborsそれぞれについて描画したい変数の組み合わせをexpand_grid()で生成する。spreadは今回は2で固定した[1]。それぞれに対してumap()をmap2()を用いて行う。もし3つ以上のパラメタを同時に試したければ、たぶんpmap()を使うことになる。数秒かかる。warningsは…無視！
purrr inside mutateを参考にしつつ、
nested_df_grid_conditions <-
expand_grid(min_dist = c(2^(-4:-1)), n_neighbors = 2^(1:4)) |>
mutate(
umap = map2(
min_dist,
n_neighbors,
~ umap(
df_image,
spread = 2,
n_threads = 4,
min_dist = .x, # min_dist = min_distでも.min_distにも非ず
n_neighbors = .y # 同様
)
)
)
# A tibble: 16 × 3
min_dist n_neighbors umap
<dbl> <dbl> <list>
1 0.00781 4 <umap>
2 0.00781 8 <umap>
3 0.00781 16 <umap>
4 0.00781 32 <umap>
5 0.0312 4 <umap>
6 0.0312 8 <umap>
7 0.0312 16 <umap>
8 0.0312 32 <umap>
9 0.125 4 <umap>
10 0.125 8 <umap>
11 0.125 16 <umap>
12 0.125 32 <umap>
13 0.5 4 <umap>
14 0.5 8 <umap>
15 0.5 16 <umap>
16 0.5 32 <umap>
map2()はlistを返す。それぞれのumapの結果について、そのなかにある$layoutに格納されているmatrixだけに用があるので、取り出したい。rowwise()を使って1行ずつほどく。umap |> extract2(1)はumap[[1]]と同義。
df_unnested_matrix_column <-
nested_df_grid_conditions |>
rowwise() |>
mutate(umap_df = umap |> extract2(1) |> list()) |>
mutate(umap = NULL) |>
ungroup() |>
unnest_longer(umap_df, indices_to = "t")
このままではmatrix-columnだ。# A tibble: 672 × 4と出るのだが、
# A tibble: 672 × 4
min_dist n_neighbors umap_df[,1] [,2] t
<dbl> <dbl> <dbl> <dbl> <int>
1 0.00781 4 -2.32 -1.14 1
2 0.00781 4 -0.427 -1.95 2
3 0.00781 4 -0.784 -1.44 3
4 0.00781 4 -1.43 -1.26 4
5 0.00781 4 -0.0518 -2.62 5
6 0.00781 4 0.616 1.41 6
7 0.00781 4 -1.57 -2.19 7
8 0.00781 4 0.860 1.94 8
9 0.00781 4 2.20 2.68 9
10 0.00781 4 2.98 2.35 10
# … with 662 more rows
# ℹ Use `print(n = ...)` to see more rows
となり、umap_df[,1], [,2]が2列と思いきや1列で、それぞれを列としてアクセスすることができない（多分）。それはそれで有用なものらしいのだが、今回は単なるフラットなデータフレーム（というかtibble）にしたい。deframe()というものをr-wakalangで教えてもらったのだが、これは2列まででないとうまくいかない。例えば次のようなかんたんなtibbleがあったとする：
repl_df <- tibble(
a = c(1, 2, 3),
b = list(
matrix(data = c(1,2,3,4), nrow = 2, ncol = 2),
matrix(data = c(1,2,3,4), nrow = 2, ncol = 2),
matrix(data = c(1,2,3,4), nrow = 2, ncol = 2)
)
)
# A tibble: 3 × 2
a b
<dbl> <list>
1 1 <dbl [2 × 2]>
2 2 <dbl [2 × 2]>
3 3 <dbl [2 × 2]>
repl_df2 <- repl_df |>
unnest(b) # 非常に近いがよくみると3列にならない…？
# A tibble: 6 × 2
a b[,1] [,2]
<dbl> <dbl> <dbl>
1 1 1 3
2 1 2 4
3 2 1 3
4 2 2 4
5 3 1 3
6 3 2 4
repl_df3 <- repl_df |>
mutate(df = b |> as_tibble()) # error
repl_df3 <- repl_df |>
mutate(df = b |>
extract2(1) |>
as_tibble()) # error
repl_df3 <- repl_df |>
rowwise() |>
mutate(
df = b |> extract2(1) |> as_tibble(),
b = NULL
)
# A tibble: 3 × 2
# Rowwise:
a df$value
<dbl> <dbl>
1 1 1
2 2 1
3 3 1
そして最終的に欲しい結果は
tibble(a = c(1, 1, 2, 2, 3, 3), V1 = c(1,2,1,2,1,2), V2 = c(3,4,3,4,3,4))
# A tibble: 6 × 3
a V1 V2
<dbl> <dbl> <dbl>
1 1 1 3
2 1 2 4
3 2 1 3
4 2 2 4
5 3 1 3
6 3 2 4
だ。この例の場合であれば教えてもらった
repl_df %>%
unnest(b) %>%
deframe() %>%
as_tibble(rownames = "a")
という方法でうまくいくが、今回は列が多くうまくいかないので抽出して2ステップでなんとかする。
df_matrix_only <-
df_unnested_matrix_column |>
select(umap_df) |>
deframe() |>
as_tibble()
# A tibble: 672 × 2
V1 V2
<dbl> <dbl>
1 -2.32 -1.14
2 -0.427 -1.95
3 -0.784 -1.44
4 -1.43 -1.26
5 -0.0518 -2.62
6 0.616 1.41
7 -1.57 -2.19
8 0.860 1.94
9 2.20 2.68
10 2.98 2.35
# … with 662 more rows
# ℹ Use `print(n = ...)` to see more rows
bind_cols()でつなげる。
df_grid_conditions <-
bind_cols(
df_unnested_matrix_column,
df_matrix_only
) |>
select(-umap_df) |>
rename(x = "V1", y = "V2")
# A tibble: 672 × 5
min_dist n_neighbors t x y
<dbl> <dbl> <int> <dbl> <dbl>
1 0.00781 4 1 -2.32 -1.14
2 0.00781 4 2 -0.427 -1.95
3 0.00781 4 3 -0.784 -1.44
4 0.00781 4 4 -1.43 -1.26
5 0.00781 4 5 -0.0518 -2.62
6 0.00781 4 6 0.616 1.41
7 0.00781 4 7 -1.57 -2.19
8 0.00781 4 8 0.860 1.94
9 0.00781 4 9 2.20 2.68
10 0.00781 4 10 2.98 2.35
# … with 662 more rows
# ℹ Use `print(n = ...)` to see more rows
プロット
あとはfacet_grid()を使うだけ。
df_grid_conditions |> colnames()
df_grid_conditions |>
ggplot(aes(x = x, y = y, colour = t)) +
geom_point(size = .3) +
geom_path(size = .1) +
scale_colour_viridis_c() +
facet_grid(rows = vars(min_dist), cols = vars(n_neighbors)) +
theme_minimal() +
theme(
panel.grid = element_blank(),
axis.line = element_line())
ggsave("./output/image_umap_grid_conditions.png", height = 1000, width = 1500, unit = "px")
データが少なすぎて・イマイチすぎて差があまり出ていないが、きちんとしたデータであればかなり差があるはずだ…。
で、ここから都合のいい感じのみた感じでチェリーピックして素知らぬ顔で報告してもコードを公開しない限り誰も気づかない気がするのだが、いいのだろうか。
-
上記スライドでは
spreadが1になっているが、min_distよりも小さい値を設定できないはずなので、min_dist = 5は誤植かなと思う ↩︎
Discussion