Модель поширення інфекційних захворювань
Раніше я писала про модель Шеллінга. А зараз мені захотілось придумати щось своє. Тож тримайте: модель, яка пояснює процес поширення інфекційних захворювань. Маємо двовимірну решітку; в моїй моделі її розмір буде фіксований, становитиме 30х30. Частина клітинок зайнята агентами, частина – порожні. Порожні клітинки мають сірий колір, їхня кількість теж фіксована і дорівнює 90. Агенти, що знаходяться у клітинках, які мають хоча б одну точку дотику, є сусідами. Маємо три групи агентів:
- Інфіковані. Їх буде позначено червоним кольором. Кількість інфікованих агентів у моделі фіксована, становить 1.
- Неінфіковані та вакциновані (далі – вакциновані). Вони мають блакитний колір.
- Неінфіковані та невакциновані (далі – невакциновані). Вони будуть фіолетового кольору.
Початково всі агенти можуть бути розташовані на решітці: 1) випадковим чином; 2) так, щоб навколо інфікованого агента були лише вакциновані та/або вільні клітинки. Якщо у вакцинованих агентів є хоча б один інфікований сусід, вони стають інфікованими з імовірністю piv, при цьому 0 <= piv <= 1. Якщо ж у невакцинованих агентів наявний хоча б один інфікований сусід, вони стають інфікованими з імовірністю pi, 0 <= pi <= 1. Процес інфікування відбувається в декілька кроків (дискретних часових проміжків) і завершується тоді, коли кількість інфікованих агентів на даному і попередньому кроці є однаковою. Загальна кількість неінфікованих агентів (є сумою кількостей вакцинованих і невакцинованих агентів) на першому кроці фіксована, дорівнює 809.
Код моделі наведено наприкінці. Бо мені здається, що так текст краще сприйматиметься.
Перейдемо до проведення комп’ютерних експериментів з моделлю. Розглянемо ситуацію, коли агенти розташовані на решітці випадковим чином.
library(kableExtra)
library(dplyr)
kable(rand, align = "c") %>%
kable_styling(bootstrap_options = "hover", font_size = 13) %>%
row_spec(0, extra_css = "vertical-align:middle;") %>%
add_header_above(c(" ", "Відсоток інфікованих на останньому кроці" = 9, " "), font_size = 14) %>%
footnote(general = "У стовчиках 3-9 наведено інтервали. Може виникнути питання, куди, наприклад, входить 25%: в інтервал Від 10% до 25% чи інтервал Від 25% до 50%? Тож я наведу для всіх інтервалів таку інтерпретацію: (0%; 1%), [1%; 10%), [10%; 25%), [25%; 50%), [50%; 75%), [75%; 90%), [90%; 99%), [99%; 100%]. Кругла дужка -- відсоток не входить в інтервал, квадратна -- входить. Відсотки розраховані по відношенню до кількості неінфікованих агентів на першому кроці", footnote_as_chunk = T)| Кількість вакцинованих на першому кроці | Немає | Від 0% до 1% | Від 1% до 10% | Від 10% до 25% | Від 25% до 50% | Від 50% до 75% | Від 75% до 90% | Від 90% до 99% | Від 99% до 100% | Всього експериментів |
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 2 | 1 | 0 | 0 | 0 | 0 | 0 | 0 | 47 | 50 |
| 1 | 1 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 49 | 50 |
| 8 | 1 | 0 | 0 | 0 | 0 | 0 | 0 | 2 | 47 | 50 |
| 81 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 50 | 0 | 50 |
| 202 | 5 | 0 | 0 | 0 | 0 | 0 | 0 | 45 | 0 | 50 |
| 405 | 5 | 6 | 2 | 0 | 0 | 0 | 2 | 35 | 0 | 50 |
| 607 | 14 | 23 | 5 | 0 | 0 | 0 | 0 | 8 | 0 | 50 |
| 728 | 26 | 20 | 3 | 0 | 0 | 0 | 0 | 1 | 0 | 50 |
| Note: У стовчиках 3-9 наведено інтервали. Може виникнути питання, куди, наприклад, входить 25%: в інтервал Від 10% до 25% чи інтервал Від 25% до 50%? Тож я наведу для всіх інтервалів таку інтерпретацію: (0%; 1%), [1%; 10%), [10%; 25%), [25%; 50%), [50%; 75%), [75%; 90%), [90%; 99%), [99%; 100%]. Кругла дужка – відсоток не входить в інтервал, квадратна – входить. Відсотки розраховані по відношенню до кількості неінфікованих агентів на першому кроці |
Експериментуючи з моделлю, варіювали значення параметру, що позначає кількість вакцинованих на першому кроці. Всього було обрано 8 значень цього параметру: 0, 1, 8 (1% від кількості неінфікованих агентів на першому кроці), 81 (10% від кількості неінфікованих агентів на першому кроці), 202 (25% від кількості неінфікованих агентів на першому кроці), 405 (50% від кількості неінфікованих агентів на першому кроці), 607 (75% від кількості неінфікованих агентів на першому кроці), 728 (90% від кількості неінфікованих агентів на першому кроці). Для кожного з 8 значень параметру моделі було проведено 50 експериментів. Імовірність зараження для невакцинованих (pi) становила 0,5, а для вакцинованих (piv) – 0,05. pi та piv є параметрами моделі, тож, якщо буде бажання, можете поекспериментувати з іншими їх значеннями.
Перейдемо до інтерпретації результатів експериментів (див. таблицю вище). Сума значень у стовпчиках 2-9 таблиці становить 50, тобто дорівнює кількості експериментів для одного значення параметру, що позначає кількість вакцинованих на першому кроці. Якщо серед агентів не було вакцинованих, у переважній більшості випадків (а саме у 47 випадках із 501) інфікований агент міг заразити від 99% до 100% неінфікованих агентів. Також коли кількість вакцинованих становила 1 та 8, у переважній більшості випадків інфікований агент міг заразити від 99% до 100% неінфікованих агентів. Якщо кількість вакцинованих агентів дорівнювала 81 та 202, інфікований агент, як правило, заражав від 90% до 99% неінфікованих агентів. У разі, коли кількість вакцинованих становила 405, інфікований агент теж у більшості випадків міг заразити від 90% до 99% неінфікованих агентів. Лише коли кількість вакцинованих досягала 607 та 728, інфікований агент заражав, в основному, до 1% неінфікованих агентів.
Як приклад можна переглянути перший і останній крок процесу інфікування для трьох випадків: 1) коли кількість вакцинованих становить 0; 2) коли кількість вакцинованих дорівнює 8; 3) коли кількість вакцинованих досягає 728.
Model_c(vaccinated = 0, pos = "random", mode = "static")
Model_c(vaccinated = 8, pos = "random", mode = "static")
Model_c(vaccinated = 728, pos = "random", mode = "static")
Тепер розглянемо ситуацію, коли агенти розташовані на решітці так, щоб навколо інфікованого агента були вакциновані та/або вільні клітинки.
kable(arnd, align = "c") %>%
kable_styling(bootstrap_options = "hover", font_size = 13) %>%
row_spec(0, extra_css = "vertical-align:middle;") %>%
add_header_above(c(" ", "Відсоток інфікованих на останньому кроці" = 9, " "), font_size = 14)| Кількість вакцинованих на першому кроці | Немає | Від 0% до 1% | Від 1% до 10% | Від 10% до 25% | Від 25% до 50% | Від 50% до 75% | Від 75% до 90% | Від 90% до 99% | Від 99% до 100% | Всього експериментів |
|---|---|---|---|---|---|---|---|---|---|---|
| 8 | 29 | 4 | 0 | 0 | 0 | 0 | 0 | 0 | 17 | 50 |
| 81 | 36 | 3 | 0 | 0 | 0 | 0 | 0 | 11 | 0 | 50 |
| 202 | 37 | 1 | 0 | 0 | 0 | 0 | 0 | 12 | 0 | 50 |
| 405 | 38 | 4 | 0 | 0 | 0 | 0 | 1 | 7 | 0 | 50 |
| 607 | 35 | 9 | 3 | 0 | 0 | 0 | 1 | 2 | 0 | 50 |
| 728 | 34 | 12 | 2 | 0 | 0 | 0 | 1 | 1 | 0 | 50 |
Експериментуючи з моделлю, обрали 6 значень параметру, що позначає кількість вакцинованих агентів на першому кроці: 8, 81, 202, 405, 607, 728. Для кожного з 6 значень було проведено 50 експериментів. Імовірність зараження для невакцинованих становила 0,5, а для вакцинованих – 0,05.
Щодо результатів експериментів маємо наступне (див. таблицю). Якщо кількість вакцинованих агентів становить 8, інфікований агент у 29 випадках із 50 не зможе нікого заразити. Разом з тим, у 17 випадках із 50 інфікований агент заражає від 99% до 100% неінфікованих агентів. Якщо кількість вакцинованих становить 81, 202 та 405, інфікований агент в більшості випадків нікого не зможе заразити. Щоправда, в деяких випадках він цілком може заразити від 90% до 99% неінфікованих агентів. Якщо кількість вакцинованих становить 607 і 728, інфікований агент, знову ж таки в більшості випадків нікого не зможе заразити. В деяких випадках він може заразити від 0% до 1% неінфікованих агентів.
А ось і приклади візуалізації для двох випадків: коли кількість вакцинованих дорівнює 8 і коли вона становить 728.
Model_c(vaccinated = 8, pos = "around", mode = "static")
Model_c(vaccinated = 728, pos = "around", mode = "static")
Які ж висновки можна зробити з усього наведеного вище? По-перше, навіть дуже невелика кількість інфікованих агентів може заразити значну кількість агентів, що не є інфікованими. По-друге, важлива не стільки кількість вакцинованих агентів, скільки їхнє розташування по відношенню до інфікованих агентів. Якщо вакциновані знаходяться навколо інфікованого агента, це знижує його шанси заразити інших.
А ось код моделі, може хтось захоче погратись:
# Model_c - модель, що дозволяє пояснити процес поширення інфекційних захворювань
# Маємо двовимірну решітку 30х30. Кількість порожніх клітинок фіксована, становить 90. Кількість інфікованих теж фіксована, дорівнює 1
# Параметри:
## vaccinated - кількість вакцинованих на першому кроці, є цілим числом, не може перевищувати 809
## pi - ймовірність зараження для невакцинованих; 0 <= pi <= 1
## piv - ймовірність зараження для вакцинованих; 0 <= piv <= 1
## pos - набуває значень random (агенти розташовані на решітці випадковим чином) та around (навколо інфікованого агента розташовуються вакциновані агенти та/або порожні клітинки)
## mode - візуалізація: static (візуалізуються перший та останній часові проміжки), dynamic (візуалізуються всі часові проміжки), none (немає візуалізації)
## res - результати у вигляді числових значень. Якщо res = T, результати виводяться, якщо res = F - ні
Model_c <- function(vaccinated,
pi = 0.5,
piv = 0.05,
pos = c("random", "around"),
mode = c("static", "dynamic", "none"),
res = F) {
# Перший крок
count <- 1
# Розподіл інфікованих (2), вакцинованих(3) та невакцинованих(4) на першому кроці. 1 - порожні клітинки
inf_vac <- c(rep(1, 90),
rep(2, 1),
rep(3, vaccinated),
rep(4, 900 - 90 - 1 - vaccinated))
if (pos == "random") {
if (vaccinated < 0) {
} else {
# Рандомізуємо
grid <- matrix(sample(inf_vac, 900), 30, 30)
}
}
if (pos == "around") {
if (vaccinated < 8) {
} else {
# Спочатку рандомізуємо
grid <- matrix(sample(inf_vac, 900), 30, 30)
# Визначаємо положення інфікованого агента
i_num <- which(grid == 2, arr.ind = T)
# Визначаємо положення всіх невакцинованих агентів
nv_num <- which(grid == 4, arr.ind = T)
# Визначаємо невакцинованих сусідів інфікованого агента
diff <- cbind(abs(nv_num[,1] - i_num[1]),
abs(nv_num[,2] - i_num[2]))
nhbrs_nv <- subset(nv_num, diff[,1] <= 1 & diff[,2] <= 1)
# Визначаємо положення всіх вакцинованих агентів
v_num <- which(grid == 3, arr.ind = T)
# Визначаємо, хто з вакцинованих агентів не є сусідом інфікованого агента
diff_2 <- cbind(abs(v_num[,1] - i_num[1]),
abs(v_num[,2] - i_num[2]))
nnhbrs_v <- subset(v_num, (diff_2[,1] > 1 & diff_2[,2] > 1) |
(diff_2[,1] <= 1 & diff_2[,2] > 1) |
(diff_2[,1] > 1 & diff_2[,2] <= 1))
# Відбираємо (випадковим чином) вакцинованих агентів, що не є сусідами інфікованого агента. Їхня кількість має дорівнювати кількості невакцинованих сусідів інфікованого агента
v_num2 <- sample(1:dim(nnhbrs_v)[1], dim(nhbrs_nv)[1])
# На місце невакцинованих сусідів інфікованого агента поміщаємо вакцинованих агентів
grid[nhbrs_nv] <- 3
# А на місце вакцинованих агентів поміщаємо невакцинованих
grid[nnhbrs_v[v_num2,]] <- 4
}
}
# Візуалізація першого кроку (а в подальшому й останнього)
if (mode == "static") {
par(mfrow = c(1, 2))
grid_col <- as.character(grid)
grid_col[grid_col == "1"] <- "gray80"
grid_col[grid_col == "2"] <- "brown1"
grid_col[grid_col == "3"] <- "deepskyblue2"
grid_col[grid_col == "4"] <- "blueviolet"
colors <- names(table(grid_col))[order(match(names(table(grid_col)),
c("gray80", "brown1", "deepskyblue2", "blueviolet")))]
image(t(apply(grid, 2, rev)),
col = colors,
oldstyle = T,
axes = F,
asp = 1)
}
# Візуалізація першого кроку (і в подальшому всіх інших кроків)
if (mode == "dynamic") {
par(mfrow = c(1, 1))
grid_col <- as.character(grid)
grid_col[grid_col == "1"] <- "gray80"
grid_col[grid_col == "2"] <- "brown1"
grid_col[grid_col == "3"] <- "deepskyblue2"
grid_col[grid_col == "4"] <- "blueviolet"
colors <- names(table(grid_col))[order(match(names(table(grid_col)),
c("gray80", "brown1", "deepskyblue2", "blueviolet")))]
image(t(apply(grid, 2, rev)),
col = colors,
oldstyle = T,
axes = F,
asp = 1)
}
df <- NULL
repeat {
# Кількість інфікованих на попередньому кроці
inf <- length(grid[grid == 2])
# Встановлюємо сусідів для кожної з клітинок
w <- cbind(rep(1:30, 30),
sort(rep(1:30, 30)))
d <- as.matrix(dist(w, method = "maximum"))
a <- apply(d, 1, function(i) grid[i == 1])
names(a) <- grid
# Розраховуємо кількість інфікованих сусідів для кожного агента
i_n <- sapply(1:900,
function(x) table(factor(a[[x]], levels = 1:4)))
i_n <- i_n[2,]
names(i_n) <- grid
i_n <- ifelse(names(i_n) == "1", NA, i_n)
# Результати: кількість інфікованих, вакцинованих і невакцинованих агентів на кожному кроці
dt <- c(length(grid[grid == 2]),
length(grid[grid == 3]),
length(grid[grid == 4]))
df <- rbind(df, dt)
# Якщо є хоча б один інфікований сусід, вакцинований агент інфікується з імовірністю piv
grid[grid == 3 & i_n > 0] <- sample(c(2, 3),
length(grid[grid == 3 & i_n > 0]),
replace = T,
prob = c(piv, 1 - piv))
# Якщо є хоча б один інфікований сусід, невакцинований агент інфікується з імовірністю pi
grid[grid == 4 & i_n > 0] <- sample(c(2, 4),
length(grid[grid == 4 & i_n > 0]),
replace = T,
prob = c(pi, 1 - pi))
# Наступний крок
count <- count + 1
# Кількість інфікованих на наступному кроці
inf_next <- length(grid[grid == 2])
# Візуалізація для кожного з кроків (окрім першого)
if (mode == "dynamic") {
grid_col <- as.character(grid)
grid_col[grid_col == "1"] <- "gray80"
grid_col[grid_col == "2"] <- "brown1"
grid_col[grid_col == "3"] <- "deepskyblue2"
grid_col[grid_col == "4"] <- "blueviolet"
colors <- names(table(grid_col))[order(match(names(table(grid_col)),
c("gray80", "brown1", "deepskyblue2", "blueviolet")))]
image(t(apply(grid, 2, rev)),
col = colors,
oldstyle = T,
axes = F,
asp = 1)
}
# Якщо кількість інфікованих на попередньому і наступному кроці співпадає, процес інфікування зупиняється
if(inf_next == inf) break
}
# Візуалізація на останньому кроці
if (mode == "static") {
grid_col <- as.character(grid)
grid_col[grid_col == "1"] <- "gray80"
grid_col[grid_col == "2"] <- "brown1"
grid_col[grid_col == "3"] <- "deepskyblue2"
grid_col[grid_col == "4"] <- "blueviolet"
colors <- names(table(grid_col))[order(match(names(table(grid_col)),
c("gray80", "brown1", "deepskyblue2", "blueviolet")))]
image(t(apply(grid, 2, rev)),
col = colors,
oldstyle = T,
axes = F,
asp = 1)
}
# Виведення результатів у вигляді числових значень
if(res == T) {
df <- as.data.frame(df)
colnames(df) <- c("Кількість інфікованих",
"Кількість вакцинованих",
"Кількість невакцинованих")
rownames(df) <- 1:dim(df)[1]
df
}
}Тобто результати 47 експериментів із 50 показали, що на останньому кроці буде від 99% до 100% інфікованих агентів (без урахування агента, що був інфікований на першому кроці). Відсоток інфікованих агентів на останньому кроці розраховувався по відношенню до кількості неінфікованих агентів на першому кроці.↩︎