Минулого разу мова йшла про алгоритм округлення чисел (він же спеціальний алгоритм; так його і будемо обзивати). Він дозволяє зменшити розбіжності між генеральною та вибірковою сукупностями. Звісно, повної відповідності не досягти, однак мінімізувати розбіжності можна. А сьогодні поговоримо про те, як покращити спеціальний алгоритм. Та спочатку згадаємо той самий вигаданий приклад неокруглених квот (див. першу частину: https://rpubs.com/ruslana/708189).
# Приклад неокруглених квот. Вигаданий
q <- t(matrix(rep(c(1.9, 1.4, 0.9, 2.3, 2.1, 1.4), 11), 6, 11))
colnames(q) <- c("Чоловіки_18-30", "Чоловіки_31-54", "Чоловіки_55+","Жінки_18-30", "Жінки_31-54", "Жінки_55+")
rownames(q) <- 1:11
# Додаємо суму за рядками і стовпчиками
q_total <- cbind(rbind(q, apply(q, 2, sum)), apply(rbind(q, apply(q, 2, sum)), 1, sum))
colnames(q_total)[7] <- "Всього"
rownames(q_total)[12] <- "Всього"
# Виводимо неокруглений варіант
# Всього 110 респондентів, по 10 у кожному селі (номери 1-11 - назви сіл)
q_total %>%
kbl(align = "c") %>%
kable_classic() %>%
column_spec(1, bold = T) %>%
row_spec(0, bold = T) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = F)| Чоловіки_18-30 | Чоловіки_31-54 | Чоловіки_55+ | Жінки_18-30 | Жінки_31-54 | Жінки_55+ | Всього | |
|---|---|---|---|---|---|---|---|
| 1 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 2 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 3 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 4 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 5 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 6 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 7 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 8 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 9 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 10 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| 11 | 1.9 | 1.4 | 0.9 | 2.3 | 2.1 | 1.4 | 10 |
| Всього | 20.9 | 15.4 | 9.9 | 25.3 | 23.1 | 15.4 | 110 |
А тепер подивимося на результати округлення за спеціальним алгоритмом (код функції наведено у першій частині, тут його немає).
# Округлюємо за спеціальним алгоритмом
q_rrnd <- r_round(q)
q_rrnd_total <- cbind(rbind(q_rrnd, apply(q_rrnd, 2, sum)), apply(rbind(q_rrnd, apply(q_rrnd, 2, sum)), 1, sum))
colnames(q_rrnd_total)[7] <- "Всього"
rownames(q_rrnd_total)[12] <- "Всього"
# Результати округлення за спеціальним алгоритмом
q_rrnd_total %>%
kbl(align = "c") %>%
kable_classic() %>%
column_spec(1, bold = T) %>%
row_spec(0, bold = T) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = F)| Чоловіки_18-30 | Чоловіки_31-54 | Чоловіки_55+ | Жінки_18-30 | Жінки_31-54 | Жінки_55+ | Всього | |
|---|---|---|---|---|---|---|---|
| 1 | 1 | 1 | 1 | 2 | 2 | 2 | 9 |
| 2 | 2 | 2 | 1 | 3 | 2 | 1 | 11 |
| 3 | 2 | 1 | 1 | 2 | 2 | 1 | 9 |
| 4 | 2 | 1 | 1 | 2 | 2 | 1 | 9 |
| 5 | 2 | 1 | 1 | 2 | 2 | 1 | 9 |
| 6 | 2 | 1 | 1 | 2 | 2 | 2 | 10 |
| 7 | 2 | 2 | 1 | 3 | 2 | 2 | 12 |
| 8 | 2 | 2 | 1 | 3 | 2 | 1 | 11 |
| 9 | 2 | 1 | 1 | 3 | 2 | 2 | 11 |
| 10 | 2 | 1 | 1 | 2 | 2 | 1 | 9 |
| 11 | 2 | 2 | 1 | 2 | 2 | 1 | 10 |
| Всього | 21 | 15 | 11 | 26 | 22 | 15 | 110 |
Спеціальний алгоритм, звісно, класний, але досконалості немає меж. Тому поговоримо про модифіковану версію спеціального алгоритму (запропонував А. Пігіда: http://dspace.nbuv.gov.ua/bitstream/handle/123456789/90428/09-Pigida.pdf?sequence=1). Орієнтуватися будемо на суму квадратів відхилень. Вона розраховується так. Спочатку рахуємо суми за рядками і стовпчиками, окремо для округленого та неокругленого варіанту квот. В наведених вище таблицях є дві графи "Всього"; це і є суми за рядками та стовпчиками. Потім знаходимо різниці сум за рядками і стовпчиками (між округленим та неокругленим варіантами) і підносимо їх до квадрату. Після цього рахуємо загальну суму квадратів різниць за всіма рядками й стовчиками. Чим менша сума квадратів відхилень, тим краще вибіркова сукупність відповідає генеральній. Мета модифікованого алгоритму - знайти мінімальну суму квадратів відхилень.
В ідеалі, для того, щоб знайти мінімальну суму квадратів відхилень, потрібно: 1) виконати спеціальний алгоритм максимально можливу кількість разів (ітерацій). Вона становить n!, де n - кількість чисел в матриці з квотами. Вийде n! матриць, округлених за спеціальним алгоритмом; 2) порахувати n! сум квадратів відхилень; 3) обрати матрицю з найменшою сумою квадратів відхилень. Але це дуже трудомісткий процес. Так, наш вигаданий приклад неокруглених квот має 11 рядків та 6 стовпчиків, всього 66 чисел. Значить, кількість ітерацій становитиме 66!, а потім ще доведеться рахувати 66! сум квадратів відхилень. Вражає, так? Тому відійдемо від цього ідеалу подалі. Врешті-решт, для спеціального алгоритму можна виконати не 66!, а, наприклад, 1000 ітерацій. І з 1000 ітерацій знайти ту, де сума квадратів відхилень мінімальна. Подивимося, що з цього вийде.
# Модифікована версія спеціального алгоритму
# a_1000 - матриця, що містить неокруглені числа
# res - виводити початковий + округлений варіант матриці чи лише округлений?
# *res = T - виводимо початковий + округлений варіант
# *res = F - лише округлений варіант
r_round_1000 <- function(a_1000, res = F){
# Округлюємо 1000 разів, використовуючи спеціальний алгоритм
ar_1000 <- lapply(1:1000, function(x) r_round(a_1000))
# Квадрат різниці сум за рядками (до і після округлення)
sqbr <- lapply(1:1000, function(x) (apply(ar_1000[[x]], 1, sum) - apply(a_1000, 1, sum))^2)
# Квадрат різниці сум за стовпчиками (до і після округлення)
sqbc <- lapply(1:1000, function(x) (apply(ar_1000[[x]], 2, sum) - apply(a_1000, 2, sum))^2)
# Сума квадратів відхилень
total <- sapply(1:1000, function(x) sum(sqbr[[x]]) + sum(sqbc[[x]]))
# Шукаємо матрицю з найменшою сумою квадратів відхилень
nmbr <- which(is.na(match(total, min(total))) == F)
if (length(nmbr) == 1) {
final_1000 <- ar_1000[[nmbr]]
} else {
final_1000 <- ar_1000[[sample(nmbr, 1)]]
}
if (res == F) {
# Виводимо лише округлений варіант матриці
final_1000
} else {
# Виводимо початковий + округлений варіант матриці
both_1000 <- list(a_1000, final_1000)
names(both_1000) <- c("Початковий варіант", "Округлений варіант")
both_1000
}
}# Округлюємо за модифікованим алгоритмом
q_rrnd_1000 <- r_round_1000(q)
q_rrnd_total_1000 <- cbind(rbind(q_rrnd_1000, apply(q_rrnd_1000, 2, sum)), apply(rbind(q_rrnd_1000, apply(q_rrnd_1000, 2, sum)), 1, sum))
colnames(q_rrnd_total_1000)[7] <- "Всього"
rownames(q_rrnd_total_1000)[12] <- "Всього"
# Результати округлення за модифікованим алгоритмом
q_rrnd_total_1000 %>%
kbl(align = "c") %>%
kable_classic() %>%
column_spec(1, bold = T) %>%
row_spec(0, bold = T) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = F)| Чоловіки_18-30 | Чоловіки_31-54 | Чоловіки_55+ | Жінки_18-30 | Жінки_31-54 | Жінки_55+ | Всього | |
|---|---|---|---|---|---|---|---|
| 1 | 2 | 1 | 1 | 2 | 2 | 1 | 9 |
| 2 | 2 | 1 | 1 | 2 | 2 | 2 | 10 |
| 3 | 2 | 2 | 1 | 2 | 2 | 2 | 11 |
| 4 | 2 | 1 | 1 | 3 | 2 | 1 | 10 |
| 5 | 2 | 1 | 1 | 3 | 2 | 1 | 10 |
| 6 | 2 | 1 | 1 | 3 | 2 | 1 | 10 |
| 7 | 2 | 1 | 1 | 2 | 2 | 2 | 10 |
| 8 | 2 | 2 | 0 | 2 | 3 | 1 | 10 |
| 9 | 2 | 2 | 1 | 2 | 2 | 1 | 10 |
| 10 | 2 | 1 | 1 | 2 | 2 | 2 | 10 |
| 11 | 2 | 2 | 1 | 2 | 2 | 1 | 10 |
| Всього | 22 | 15 | 10 | 25 | 23 | 15 | 110 |
Порівняємо суму квадратів відхилень для спеціального та модифікованого алгоритму.
# Сума квадратів відхилень для спеціального алгоритму
brw <- (apply(q, 1, sum) - apply(q_rrnd, 1, sum))^2
bclmn <- (apply(q, 1, sum) - apply(q_rrnd, 1, sum))^2
total <- sum(brw + bclmn)
# Сума квадратів відхилень для модифікованого алгоритму
brw_1000 <- (apply(q, 1, sum) - apply(q_rrnd_1000, 1, sum))^2
bclmn_1000 <- (apply(q, 1, sum) - apply(q_rrnd_1000, 1, sum))^2
total_1000 <- sum(brw_1000 + bclmn_1000)
# Порівнюємо...
sq <- rbind(total, total_1000)
rownames(sq) <- c("Округлення за спеціальним алгоритмом", "Округлення за модифікованим алгоритмом")
colnames(sq) <- ""
sq %>%
kbl(align = "c") %>%
kable_classic_2 %>%
kable_styling(full_width = F)| Округлення за спеціальним алгоритмом | 24 |
| Округлення за модифікованим алгоритмом | 4 |
Ну і нарешті порівняємо розподіл за статтю і віком для неокругленого, спеціального та модифікованого алгоритму.
# Порівняння розподілу за статтю і віком
zr <- rbind(q_total[12, c(-7)], q_rrnd_total[12, c(-7)], q_rrnd_total_1000[12, c(-7)])
rownames(zr) <- c("Неокруглений варіант", "Округлення за спеціальним алгоритмом", "Округлення за модифікованим алгоритмом")
zr %>%
kbl(align = "c") %>%
kable_classic() %>%
column_spec(1, bold = T) %>%
row_spec(0, bold = T) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = F)| Чоловіки_18-30 | Чоловіки_31-54 | Чоловіки_55+ | Жінки_18-30 | Жінки_31-54 | Жінки_55+ | |
|---|---|---|---|---|---|---|
| Неокруглений варіант | 20.9 | 15.4 | 9.9 | 25.3 | 23.1 | 15.4 |
| Округлення за спеціальним алгоритмом | 21.0 | 15.0 | 11.0 | 26.0 | 22.0 | 15.0 |
| Округлення за модифікованим алгоритмом | 22.0 | 15.0 | 10.0 | 25.0 | 23.0 | 15.0 |
Так, щось забагато букв вийшло... Але ми побачили, що модифікована версія спеціального алгоритму дозволяє отримати ще кращу відповідність вибіркової сукупності генеральній. До речі, можна варіювати кількість ітерацій. Наприклад, їх може бути не 1000 (як у мене), а 10000. Чим більше ітерацій, тим більша імовірність отримати меншу суму квадратів відхилень. Разом з тим, збільшення кількості ітерацій сповільнюватиме роботу алгоритму.