Карта двоичной плотности в R - PullRequest
Новая
       22

Карта двоичной плотности в R

1 голос
/ 05 апреля 2020

У меня большой набор данных с парами координат x и y мест проживания домохозяйств.

Я хочу создать бинарный график плотности, который различает плотные и не плотные области (точный порог будет определен эвристически). Поэтому я хочу сказать, что если пара координат находится в определенной области, она считается плотной. Я не уверен, что подходящие подходы в этом случае - кригинг или обратное взвешивание.

Мой подход на данный момент следующий: 

library(ggplot2)

m <- ggplot(df, aes(x = df$V1, y = df$V2)) +
  geom_point() 
m + geom_density_2d()

К сожалению, это не позволяет мне извлечь binary map, но отображает только плотности. Мыслящим решением является третий столбец, который сортирует пары координат по плотным и не плотным категориям.

Ваша помощь очень ценится!

данные выдержка из данных (всего 162'000 строк): 

df <- structure(list(V1 = c(2690503, 2689797, 2690685, 2690685, 2689409, 
2689409, 2689141, 2689141, 2689141, 2689141, 2690515, 2690515, 
2689474, 2689148, 2689148, 2689148, 2689148, 2689148, 2689148, 
2690199, 2689760, 2689473, 2689198, 2689198, 2689558, 2690020, 
2690020, 2689612, 2689132, 2689208, 2689208, 2690244, 2689614, 
2689614, 2690264, 2690264, 2689625, 2690315, 2689720, 2689720, 
2689720, 2691037, 2691037, 2691037, 2691037, 2689433, 2689433, 
2690715, 2690715, 2690715, 2689560, 2689655, 2689563, 2690240, 
2690240, 2678781, 2689498, 2689206, 2690521, 2689662, 2689662, 
2690027, 2690027, 2689383, 2689383, 2690588, 2688999, 2689397, 
2689725, 2689725, 2689100, 2689100, 2689100, 2689100, 2689906, 
2690529, 2690529, 2690199, 2690164, 2689638, 2689638, 2689498, 
2689093, 2689093, 2690502, 2689740, 2689153, 2690536, 2689027, 
2689950, 2690959, 2690959, 2690486, 2690486, 2689794, 2689307, 
2690010, 2690010, 2689599, 2689599, 2689599, 2689464, 2689464, 
2689464, 2689464, 2689711, 2689711, 2690222, 2690602, 2690602, 
2689790, 2689790, 2689404, 2689428, 2689428, 2689387, 2688960, 
2688960, 2689258, 2689258, 2689258, 2689355, 2690795, 2689521, 
2689954, 2689954, 2688926, 2689543, 2689363, 2689363, 2689186, 
2689186, 2689571, 2689571, 2689970, 2689970, 2689675, 2689498, 
2690941, 2690941, 2689060, 2689442, 2690122, 2690725, 2690725, 
2689419, 2689366, 2689366, 2689097, 2689332, 2690123, 2690123, 
2690994, 2690994, 2690180, 2690180, 2689706, 2689706, 2689612, 
2690074, 2688961, 2688961, 2689692, 2689260, 2689419, 2689419, 
2689146, 2690083, 2690625, 2690167, 2690167, 2689540, 2689540, 
2689512, 2689512, 2690469, 2689720, 2689711, 2690874, 2690072, 
2690072, 2690072, 2688946, 2689502, 2689431, 2689531, 2689131, 
2689131, 2690257, 2690001, 2689608, 2689843, 2689502, 2689773, 
2689773, 2689507, 2690060, 2678781, 2689500, 2689260), V2 = c(1254816, 
1254916, 1254061, 1254061, 1255542, 1255542, 1255220, 1255220, 
1255220, 1255220, 1254872, 1254872, 1255561, 1255199, 1255199, 
1255199, 1255199, 1255199, 1255199, 1255390, 1255667, 1255233, 
1255830, 1255830, 1255029, 1254812, 1254812, 1255297, 1255391, 
1255728, 1255728, 1254961, 1255385, 1255385, 1255149, 1255149, 
1255704, 1255312, 1254949, 1254949, 1254949, 1253836, 1253836, 
1253836, 1253836, 1255130, 1255130, 1253886, 1253886, 1253886, 
1255124, 1254928, 1255858, 1255267, 1255267, 1237314, 1255231, 
1255426, 1254796, 1255315, 1255315, 1255231, 1255231, 1255065, 
1255065, 1254882, 1255504, 1255493, 1255279, 1255279, 1256005, 
1256005, 1256005, 1256005, 1255418, 1254909, 1254909, 1255390, 
1255233, 1255716, 1255716, 1255231, 1255787, 1255787, 1253745, 
1255672, 1255827, 1254775, 1255813, 1255187, 1254105, 1254105, 
1255155, 1255155, 1255128, 1255623, 1255448, 1255448, 1255397, 
1255397, 1255397, 1255353, 1255353, 1255353, 1255353, 1255306, 
1255306, 1254824, 1254771, 1254771, 1255170, 1255170, 1255380, 
1255919, 1255919, 1255204, 1255885, 1255885, 1256001, 1256001, 
1256001, 1255552, 1254091, 1255334, 1255052, 1255052, 1255609, 
1254960, 1255090, 1255090, 1255426, 1255426, 1255140, 1255140, 
1254886, 1254886, 1255464, 1255231, 1254052, 1254052, 1255454, 
1255219, 1255265, 1254818, 1254818, 1255145, 1255407, 1255407, 
1255180, 1255187, 1255179, 1255179, 1253875, 1253875, 1255001, 
1255001, 1255573, 1255573, 1255232, 1255023, 1255458, 1255458, 
1255305, 1255453, 1255145, 1255145, 1255275, 1255209, 1254931, 
1255447, 1255447, 1255003, 1255003, 1255760, 1255760, 1254749, 
1254949, 1255306, 1253851, 1255286, 1255286, 1255286, 1255599, 
1255183, 1255580, 1255320, 1255460, 1255460, 1254804, 1255360, 
1255731, 1255470, 1255183, 1255131, 1255131, 1255081, 1255245, 
1237314, 1255586, 1255419)), row.names = c(NA, 200L), class = "data.frame")

1 Ответ

2 голосов
/ 05 апреля 2020

Аналогичный вопрос был задан в этой теме . Я изменяю полученный ответ.

Я использую MASS:kde2d, чтобы создать фрейм данных, который вы на самом деле также упоминали. Обратите внимание, что эта функция даст матрицу z, которая - как обычно - упорядочена по столбцу. Поэтому вам нужно повторить координаты y на each, а x на times 

library(tidyverse)

# Manually estimate the density for each coordinate value with MASS::kde2d
# I scale the data first and remove this weird outlier in your data 
scaled_df <- df %>% scale() %>% as.data.frame() %>% filter (V1 > -7)

kde2d_est <- MASS::kde2d(scaled_df$V1, scaled_df$V2)

res <- data.frame(x = rep(kde2d_est$x, times = dim(kde2d_est$z)[1]),
                  y = rep(kde2d_est$y, each = dim(kde2d_est$z)[2]),
                  density = as.numeric(kde2d_est$z))

ggplot() +
  geom_point(data = scaled_df, aes(V1, V2)) +
  geom_tile(data = res, aes(x = x, y = y, fill = density), alpha = 0.5)

Это выглядит уже неплохо

Теперь давайте сделаем "двоичную плотность", с произвольной отсечкой. Вам не нужны сокращения, как в ответе, на который я ссылаюсь. Просто используйте логическое утверждение. 

arbitrary_cut <- 1
res_bin <- res %>% mutate(binary_dens = density < arbitrary_cut)

ggplot() +
  geom_point(data = scaled_df, aes(V1, V2)) +
  geom_tile(data = res_bin, aes(x = x, y = y, fill = binary_dens), alpha = 0.5)

Создано в 2020-04-05 пакетом Представить (v0.3.0)

...