Модель Шеллінга: випадкове розміщення vs розміщення у шаховому порядку

З'явилося у мене бажання щось написати про модель Шеллінга. І я вирішила втілити його у життя. Не буду довго мучити і перейду одразу до теорії. Модель Шеллінга - одна з перших агентних моделей. Вона дозволяє пояснити, як бажання агентів проживати серед представників своєї групи може призводити до виникнення сегрегації. Агенти розташовані на двовимірній решітці і належать до однієї з двох груп. Решітка не повністю заповнена агентами, там є і вільні клітинки. Якщо агент не задоволений своїм розташуванням, він переходить на вільну клітинку. В іншому разі він залишається на місці. Агент задоволений своїм розташуванням тоді, коли частка сусідів, які є представниками тієї ж групи, не менша за поріг. Поріг характеризує рівень толерантності до представників іншої групи, його значення знаходяться в межах від 0 до 1. Чим менше значення порогу, тим вищий рівень толерантності і навпаки. Так, значення порогу 0.2 вказує на те, що агентові комфортно проживати щонайменше з 20% сусідів, що належать до тієї ж групи. Решта (максимум 80% сусідів) можуть належати до іншої групи. Рух агентів завершується тоді, коли всі вони будуть задоволені своїм розташуванням.

Модель Шеллінга передбачає, що початково агенти розташовуються на решітці випадковим чином. Але я ще вирішила розглянути варіант, коли агенти розташовуються у шаховому порядку (частина клітинок будуть вільними).

# Маємо двовимірну решітку 32 Х 32
# r - поріг, змінюється від 0 до 1
# chs - початкове розташування агентів
# *якщо chs = F, агенти розміщені випадковим чином
# *якщо chs = T, агенти розміщені у шаховому порядку
# g1 - кількість агентів у 1-ій групі, ціле число > 0 
# *якщо chs = T, параметр g1 не має сенсу 
# n - кількість порожніх клітинок, ціле число > 0
# *n + g1 має бути меншим за 1024
# vz - візуалізація
# *якщо vz = "none", візуалізація не виводиться
# *якщо vz = "static", виводимо лише початкове і кінцеве розміщення
# *якщо vz = "dynamic", спостерігаємо рух агентів на кожному кроці
# res - потрібно виводити результати симуляцій чи ні
# *якщо res = F, результати не виводяться
# *якщо res = T, результати виводяться
Schelling_model <- function(r, chs = F, g1 = 462, n = 100, vz = c("none", "static", "dynamic"), res = F) {
  # Перший крок (початкове розміщення)
  count <- 1
  if (g1 <= 0 | n <= 0 | g1 + n >= 1024 | r < 0 | r > 1 | is.logical(chs) == F | vz != "none" & vz != "static" & vz != "dynamic" | is.logical(res) == F | is.na(res) == T) {
  } else {
    # Початкове розміщення агентів
    if (chs == F) {
      group <- c(rep(0, n), rep(1, g1), rep(2, 1024 - g1 - n))
      grid <- matrix(sample(group, 1024), 32, 32)
    } else {
      grid <- matrix(rep(c(rep(c(1, 2), 16), rep(c(2, 1), 16)), 16), 32, 32)
      k <- sample(c(floor(n/2), ceiling(n/2)), 1)
      grid[sample(which(grid == 1), k)] <- 0
      grid[sample(which(grid == 2), n - k)] <- 0
    }
  }
  # Візуалізація початкового розміщення агентів
  # Gray - вільні клітинки
  # Dark green - агенти,  представники 1-ої групи
  # Red - представники 2-ої групи
  if (vz == "static") {
    par(mfrow = c(1,2))
    image(t(apply(grid, 2, rev)), col = c("gray", "dark green", "red"), axes = F, asp = 1)
  }
  if (vz == "dynamic") {
    par(mfrow = c(1,1))
    image(t(apply(grid, 2, rev)), col = c("gray", "dark green", "red"), axes = F, asp = 1)
  }
  # Крім візуалізації ще будуть результати симуляцій
  df <- NULL
  repeat {
    # Встановлюємо сусідів для кожної з клітинок
    w <- cbind(rep(c(1:dim(grid)[1]), dim(grid)[2]), sort(rep(c(1:dim(grid)[1]), dim(grid)[2])))
    d <- as.matrix(dist(w, method = "maximum"))
    a <- apply(d, 1, function(i) grid[i == 1])
    names(a) <- grid
    # Розраховуємо загальну кількість сусідів для кожного агента
    t_ns <- sapply(1:length(grid), function(x) sum(table(a[[x]], exclude = 0)))
    names(t_ns) <- grid
    t_ns <- ifelse(names(t_ns) == "0", NA, t_ns)
    # Розраховуємо кількість сусідів агента, що належать до тієї ж групи
    g_s <- as.data.frame(sapply(1:length(grid), function(x) table(factor(a[[x]], levels = c("0", "1", "2")))))
    colnames(g_s) <- grid
    g_ns <- sapply(1:1024, function(x) g_s[match(colnames(g_s)[x], rownames(g_s)), x])
    names(g_ns) <- grid
    g_ns <- ifelse(names(g_ns) == "0", NA, g_ns)
    # Розраховуємо частку сусідів "своєї" групи
    p <- ifelse(t_ns == 0, 0, g_ns/t_ns)
    names(p) <- grid
    # Визначимо середню частку сусідів "своєї" групи
    p_sg <- round(mean(p,na.rm = T), 3)
    # Визначаємо агентів, задоволених і незадоволених своїм розташуванням
    satisfied <- p
    satisfied[satisfied < r] <- 0
    satisfied[satisfied >= r] <- 1
    # Розрахуємо частку агентів, задоволених своїм розташуванням
    p_st <- round(mean(satisfied, na.rm = T), 3)
    # Формування результатів симуляцій
    df <- rbind(df, c(p_sg, p_st))
    # Незадоволені агенти переходять на вільну клітинку
    if (length(which(is.na(satisfied) == T)) <= length(which(satisfied == 0 & is.na(satisfied) == F))) {
      ls1 <- sample(which(satisfied == 0 & is.na(satisfied) == F), length(which(is.na(satisfied) == T)))
      names(satisfied)[which(is.na(satisfied) == T)] <- names(satisfied)[ls1]
      names(satisfied)[ls1] <- "0"
    } else {
      ls2 <- sample(which(is.na(satisfied) == T), length(which(satisfied == 0 & is.na(satisfied) == F)))
      names(satisfied)[ls2] <- names(satisfied)[which(satisfied == 0 & is.na(satisfied) == F)] 
      names(satisfied)[which(satisfied == 0 & is.na(satisfied) == F)] <- "0"
    }
    # Наступний крок
    count <- count + 1
    # Нове розміщення
    grid <- matrix(as.numeric(names(satisfied)), 32, 32)
    # Візуалізація розміщення агентів на кожному кроці
    if (vz == "dynamic") {
      image(t(apply(grid, 2, rev)), col = c("gray", "dark green", "red"), axes = F, asp = 1)
    }
    # Умови для завершення руху агентів
    if (mean(satisfied, na.rm = T) == 1 | count == 300) break
  }
  # Візуалізація кінцевого розміщення агентів
  if (vz == "static") {
    image(t(apply(grid, 2, rev)), col = c("gray", "dark green", "red"), axes = F, asp = 1)
  }
  # Результати симуляцій
  if (res == T) {
    df <- as.data.frame(df)
    colnames(df) <- c("Середня частка своїх", "Частка задоволених")
    df
  }
}

Трохи детальніше про результати симуляцій, які дозволяє вивести функція (якщо параметр res = T). Розраховується середня частка сусідів своєї групи (середня частка своїх) та частка агентів, задоволених своїм розташуванням (частка задоволених) на кожному кроці. На першому кроці маємо справу з початковим розміщенням агентів. Середня частка сусідів своєї групи розраховується так. Спочатку для кожного агента визначається частка сусідів, які належать до тієї ж групи. Потім знаходимо середнє арифметичне цих часток. Середня частка сусідів своєї групи набуває значень від 0 до 1. Наприклад, якщо значення цього показника становить 0.3, це означає, що в середньому агенти мають 30% сусідів, що належать до тієї ж групи; решта 70% сусідів належать до іншої групи. Також середня частка сусідів своєї групи вказує на те, наскільки сильною є сегрегація. Чим вищий цей показник, тим більш виражена сегрегація. А частка агентів, задоволених своїм розташуванням, обчислюється так. Спочатку визначаємо агентів, задоволених та незадоволених своїм розташуванням. А потім ділимо кількість агентів, задоволених своїм розташуванням на загальну кількість агентів. Частка агентів, задоволених своїм розташуванням теж набуває значень від 0 до 1. За високих значень порогу - орієнтовно 0.8 і більше - рух агентів може тривати дуууже довго (а то й нескінченно). Тож було встановлено значення кроку (300), досягнення якого означає завершення руху агентів.

Також деякі зі значень параметрів задані за замовчуванням. Так, кількість вільних клітинок дорівнює 100: це приблизно 10% від загальної кількості клітинок. А кількість агентів, що належить першої групи становить 462; такою ж буде кількість представників другої групи (за умови, якщо кількість вільних клітинок становить 100). Я буду мати справу лише зі значеннями параметрів n, g1 та res, заданими за замовчуванням. А ви, якщо буде бажання, можете експериментувати з іншими значеннями вказаних вище параметрів:) Мало не забула: всі агенти мають однаковий поріг.

Фух, перейдемо до симуляцій з моделлю. Розпочнемо з випадку, коли агенти розташовані випадковим чином. Вже за значення порогу 0.3 відмінності між початковим і кінцевим розташуванням агентів візуально помітні. Після завершення руху агентів можемо помітити утворення своєрідних "кварталів", де проживають представники однієї групи. Про всяк випадок уточню: зліва - початкове розташування агентів, справа - кінцеве розташування, яке можна спостерігати після завершення руху агентів.

Schelling_model(0.3, vz = "static")

Тепер поглянемо на початкове і кінцеве розміщення агентів за значень порогів 0.4, 0.5 та 0.7:

Schelling_model(0.4, vz = "static")

Schelling_model(0.5, vz = "static")

Schelling_model(0.7, vz = "static")

До речі, у випадку з порогом 0.7 спостерігаємо значну сегрегацію. Вільні клітинки створюють "паркан", який розділяє обидві групи агентів.

Тепер розглянемо випадок, коли агенти розміщені у шаховому порядку. Візьмемо ті ж значення порогів, які розглядалися вище: 0.3, 0.4, 0.5 та 0.7.

Schelling_model(0.3, chs = T, vz = "static")

Schelling_model(0.4, chs = T, vz = "static")

Schelling_model(0.5, chs = T, vz = "static")

Schelling_model(0.7, chs = T, vz = "static")

За значення порогу 0.3 різниця між початковим та кінцевим розташуванням агентів візуально майже не помітна. За значення порогу 0.4 вже можемо спостерігати утворення "кварталів", де проживають представники однієї групи. А за значення порогу 0.7 спостерігаємо значну сегрегацію: групи відділені одне від одного "парканом" із вільних клітинок.

Тож сегрегація може виникати і тоді, коли агенти початково розміщені випадковим чином, і тоді, коли вони розташовані у шаховому порядку. Так, за значення порогу 0.7 в обох випадках спостерігаємо значну сегрегацію. Однак є і відмінності. Якщо агенти розташовані у шаховому порядку, за значення порогу 0.3 різниці між початковим та кінцевим розміщенням візуально майже не видно. Якщо агенти розташовані випадковим чином, ця різниця (за значення порогу 0.3) візуально більш помітна.