Практика · Кластеризація k-means¶

Лекція: lecture.html · Тест: quiz.html · Домашнє: homework.md

Уся попередня половина курсу трималась на колонці з правильною відповіддю. Тут її немає: ми шукаємо сегменти ринку на дошці оголошень про вживані телефони — ту саму, що і в темі 08, — не маючи жодної мітки.

Що зробимо:

  1. зберемо ту саму дошку й лишимо дві ознаки — ціну й рік;
  2. напишемо k-means з нуля десятком рядків і подивимось, як падає інерція;
  3. звіримо свої числа з KMeans зі scikit-learn — вони мають зійтися;
  4. проженемо кілька випадкових стартів і побачимо розкид відповідей;
  5. побудуємо криву ліктя й силует для k від 2 до 10;
  6. подивимось, що буває без масштабування ознак;
  7. і наприкінці чесно зіставимо кластери з колонкою шахрайське, якої алгоритм не бачив.

Усі числа тут ті самі, що в лекції: генератор випадкових чисел зафіксовано зерном 42.

In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

from sklearn.preprocessing import StandardScaler
from sklearn.cluster import KMeans
from sklearn.metrics import silhouette_score, silhouette_samples, adjusted_rand_score
from sklearn.metrics import precision_score, recall_score, f1_score

# зерно фіксує всю випадковість: у тебе вийдуть точно ті самі числа, що в лекції
rng = np.random.default_rng(42)

pd.set_option("display.width", 120)
pd.set_option("display.max_columns", 12)
print("numpy", np.__version__, "· pandas", pd.__version__)
numpy 2.5.1 · pandas 2.3.3

1 · Збираємо ту саму дошку¶

Цей блок — код із практики теми 08 без змін: 1 200 оголошень і шість типових неприємностей, які там докладно розібрано. Тут ми лише відтворюємо таблицю, щоб числа збіглися до цифри.

In [2]:
кількість = 1200

моделі = ["Alfa A5", "Alfa A7", "Beta 12", "Beta 12 Pro", "Gamma X", "Gamma X Ultra"]
ціна_нового = {"Alfa A5": 5200, "Alfa A7": 7400, "Beta 12": 12000,
               "Beta 12 Pro": 17500, "Gamma X": 24000, "Gamma X Ultra": 34000}
частки_моделей = [0.24, 0.22, 0.18, 0.16, 0.12, 0.08]

модель = rng.choice(моделі, size=кількість, p=частки_моделей)
рік = rng.integers(2017, 2025, size=кількість)
стан = rng.choice(["нове", "дуже добре", "добре", "задовільне"],
                  size=кількість, p=[0.08, 0.32, 0.42, 0.18])
памʼять = rng.choice([64, 128, 256, 512], size=кількість, p=[0.30, 0.38, 0.24, 0.08])
вік_акаунта = np.round(rng.exponential(420, size=кількість) + 3).astype(int)

базова = np.array([ціна_нового[m] for m in модель])
знос = 0.82 ** (2024 - рік)
коефіцієнт_стану = np.array(
    [{"нове": 1.0, "дуже добре": 0.88, "добре": 0.75, "задовільне": 0.58}[s] for s in стан])
коефіцієнт_памʼяті = np.array(
    [{64: 0.85, 128: 1.0, 256: 1.15, 512: 1.32}[m] for m in памʼять])

типова_ціна_за_паспортом = базова * знос * коефіцієнт_стану * коефіцієнт_памʼяті
ціна = типова_ціна_за_паспортом * rng.lognormal(0, 0.13, size=кількість)

print("перші пʼять цін:", ціна[:5].round(0))
перші пʼять цін: [ 3085.  2810. 10822.  8346.  1852.]
In [3]:
# шахрай тим імовірніший, чим молодший акаунт; ціну він або занижує (приманка),
# або завищує під велику передоплату
шанс_шахрайства = 0.10 + 0.30 * np.exp(-вік_акаунта / 120)
шахрайське = rng.random(кількість) < шанс_шахрайства

ставить_дешево = rng.random(кількість) < 0.74
дешева_приманка = шахрайське & ставить_дешево
дорога_приманка = шахрайське & ~ставить_дешево

ціна[дешева_приманка] = (типова_ціна_за_паспортом[дешева_приманка]
                         * rng.uniform(0.20, 0.45, дешева_приманка.sum()))
ціна[дорога_приманка] = (типова_ціна_за_паспортом[дорога_приманка]
                         * rng.uniform(2.6, 3.8, дорога_приманка.sum()))
ціна = np.round(ціна, -1)

скарг = np.where(шахрайське, 1 + rng.poisson(3.0, кількість), rng.poisson(0.03, кількість))

дошка = pd.DataFrame({
    "модель": модель, "рік": рік, "стан": стан, "памʼять_гб": памʼять,
    "вік_акаунта": вік_акаунта, "скарг": скарг, "ціна": ціна,
    "шахрайське": шахрайське.astype(int),
})
print("шахрайських оголошень:", int(дошка["шахрайське"].sum()), "з", кількість)
шахрайських оголошень: 196 з 1200
In [4]:
# ті самі шість неприємностей із теми 08 — колекційні, одруки, памʼять текстом,
# пропуски в ціні та стані, дублікати
колекційні = дошка.index[дошка["модель"] == "Gamma X"][:4]
дошка.loc[колекційні, ["рік", "стан", "памʼять_гб"]] = [2017, "нове", 512]
дошка.loc[колекційні, "ціна"] = [82000.0, 88000.0, 91000.0, 95000.0]
дошка.loc[колекційні, ["шахрайське", "скарг"]] = 0

одруки = дошка.index[(дошка["ціна"] > 7000) & (дошка["ціна"] < 9600)
                     & (дошка["шахрайське"] == 0)][:2]
дошка.loc[одруки, "ціна"] = дошка.loc[одруки, "ціна"] * 10

памʼять_текстом = дошка["памʼять_гб"].astype(str)
із_одиницями = rng.random(len(дошка)) < 0.18
памʼять_текстом[із_одиницями] = памʼять_текстом[із_одиницями] + " ГБ"
дошка["памʼять_гб"] = памʼять_текстом

ймовірність_пропуску = np.where(дошка["шахрайське"] == 1, 0.25, 0.03)
дошка.loc[rng.random(len(дошка)) < ймовірність_пропуску, "ціна"] = np.nan
дошка.loc[rng.random(len(дошка)) < 0.04, "стан"] = np.nan

повтори = rng.choice(дошка.index, size=12, replace=False)
дошка = pd.concat([дошка, дошка.loc[повтори]], ignore_index=True)

print("таблиця як у темі 08:", дошка.shape)
assert дошка.shape == (1212, 8), "форма розійшлась із темою 22"
print("✅ форма збігається")
таблиця як у темі 08: (1212, 8)
✅ форма збігається

2 · Дві ознаки й чому саме такі¶

Кластеризація рахує відстані, тож їй потрібні числа й потрібні всі. Тому:

  • прибираємо повні дублікати рядків;
  • повертаємо памʼяті числовий вигляд (128 ГБ → 128);
  • викидаємо оголошення без ціни — для них головна ознака просто не існує;
  • беремо дві ознаки: логарифм ціни й рік випуску.

Чому логарифм? Різниця між 1 000 і 2 000 грн для покупця така сама, як між 10 000 і 20 000, а в гривнях це 1 000 проти 10 000. Логарифм робить ці дві відстані однаковими — і заразом приборкує хвіст із кількох оголошень по 80–95 тисяч (тема 10).

In [5]:
ринок = дошка.drop_duplicates().reset_index(drop=True)
ринок["памʼять_гб"] = ринок["памʼять_гб"].astype(str).str.replace(" ГБ", "", regex=False).astype(int)
ринок = ринок.dropna(subset=["ціна"]).reset_index(drop=True)

ціни = ринок["ціна"].to_numpy(float)
роки = ринок["рік"].to_numpy(int)
ознаки = np.c_[np.log(ціни), роки]

print("оголошень для кластеризації:", len(ринок))
print("шахрайських серед них:", int(ринок["шахрайське"].sum()),
      f"({100 * ринок['шахрайське'].mean():.1f} %)")
print("ціна: від", int(ціни.min()), "до", int(ціни.max()), "грн")
оголошень для кластеризації: 1100
шахрайських серед них: 139 (12.6 %)
ціна: від 130 до 95000 грн

3 · Масштаб: чому без нього не можна¶

k-means складає квадрати різниць по всіх ознаках. Подивись, який розкид має кожна колонка в своїх одиницях — і стане зрозуміло, чия різниця переважить.

In [6]:
# стандартне відхилення — це і є «типова відстань між двома об'єктами» по цій ознаці
print(f"розкид ціни у гривнях: {ціни.std():10.1f}")
print(f"розкид року в роках:   {роки.std():10.4f}")
print(f"відношення:            {ціни.std() / роки.std():10.0f} разів")

масштабувач = StandardScaler().fit(ознаки)
масштабовані = масштабувач.transform(ознаки)
print("\nпісля StandardScaler обидві ознаки мають розкид 1.0:",
      масштабовані.std(axis=0).round(6))
розкид ціни у гривнях:     8685.0
розкид року в роках:       2.3069
відношення:                  3765 разів

після StandardScaler обидві ознаки мають розкид 1.0: [1. 1.]

4 · k-means з нуля¶

Увесь алгоритм — це два кроки по черзі:

  1. розподілити — кожну точку до найближчого центра;
  2. пересунути — кожен центр у середнє своїх точок.

Повторюємо, поки центри перестануть рухатись. Нижче — рівно це, з журналом інерції на кожному кроці.

In [7]:
# Алгоритм Ллойда: повертає центри, мітки та інерцію
def кластеризувати(точки, центри, максимум_кроків=100, друкувати=False):
    центри = центри.copy()
    for крок in range(максимум_кроків):
        # квадрат відстані від кожної точки до кожного центра: (n, k)
        відстані = ((точки[:, None, :] - центри[None, :, :]) ** 2).sum(axis=2)
        мітки = відстані.argmin(axis=1)          # крок 1: до найближчого
        if друкувати:
            print(f"  крок {крок:2d}: інерція {відстані.min(axis=1).sum():9.2f}"
                  f"   розміри {np.bincount(мітки, minlength=len(центри))}")
        # крок 2: центр стає середнім своїх точок; порожній центр лишаємо на місці
        нові = np.array([точки[мітки == c].mean(axis=0) if (мітки == c).any() else центри[c]
                         for c in range(len(центри))])
        if np.allclose(нові, центри):
            break
        центри = нові
    відстані = ((точки[:, None, :] - центри[None, :, :]) ** 2).sum(axis=2)
    return центри, відстані.argmin(axis=1), відстані.min(axis=1).sum()


print("функція готова — 12 рядків, і це весь метод")
функція готова — 12 рядків, і це весь метод

Проженемо її на наших даних із трьох випадкових оголошень як стартових центрів. Дивись на стовпчик інерції: він жодного разу не зростає.

In [8]:
# ті самі три оголошення, що в лекції: старт A
старт_A = масштабовані[[960, 1003, 313]]
print("старт A, k = 3:")
центри_A, мітки_A, інерція_A = кластеризувати(масштабовані, старт_A, друкувати=True)
print(f"\nзупинилось на інерції {інерція_A:.4f}")
старт A, k = 3:
  крок  0: інерція   1486.79   розміри [285 169 646]
  крок  1: інерція    943.16   розміри [367 220 513]
  крок  2: інерція    852.77   розміри [419 243 438]
  крок  3: інерція    819.80   розміри [454 252 394]
  крок  4: інерція    806.63   розміри [463 254 383]
  крок  5: інерція    805.12   розміри [469 254 377]
  крок  6: інерція    804.35   розміри [473 256 371]
  крок  7: інерція    803.92   розміри [477 260 363]
  крок  8: інерція    803.55   розміри [480 261 359]
  крок  9: інерція    803.47   розміри [480 263 357]
  крок 10: інерція    803.44   розміри [481 264 355]
  крок 11: інерція    803.43   розміри [481 264 355]

зупинилось на інерції 803.4317

5 · Звірка з бібліотекою¶

Найцінніша клітинка практики: усередині KMeans немає магії. Дамо scikit-learn ті самі дані, ту саму кількість кластерів — і порівняємо два числа й два розбиття.

adjusted_rand_score міряє схожість двох розбиттів: 1.0 — це «ті самі групи», хай навіть номери в них переставлені.

In [9]:
бібліотечний = KMeans(n_clusters=3, n_init=10, random_state=42, tol=1e-10)
бібліотечний.fit(масштабовані)

print(f"наша інерція:     {інерція_A:.4f}")
print(f"інерція KMeans:   {бібліотечний.inertia_:.4f}")
print(f"збіг розбиттів (ARI): {adjusted_rand_score(мітки_A, бібліотечний.labels_):.4f}")

assert np.isclose(інерція_A, бібліотечний.inertia_, rtol=1e-6), "інерція розійшлась!"
assert adjusted_rand_score(мітки_A, бібліотечний.labels_) == 1.0, "розбиття різні!"
print("✅ збігається")
наша інерція:     803.4317
інерція KMeans:   803.4317
збіг розбиттів (ARI): 1.0000
✅ збігається

6 · Той самий алгоритм, різні відповіді¶

Алгоритм гарантовано зупиняється, але не гарантовано в найкращій точці. Запустимо його з чотирьох різних стартів — тих самих, що в лекції, — і подивимось на інерцію.

In [10]:
старти = {"A": [960, 1003, 313], "B": [295, 925, 337],
          "C": [375, 333, 896], "D": [809, 340, 66]}

рядки = []
for назва, індекси in старти.items():
    центри, мітки, інерція = кластеризувати(масштабовані, масштабовані[індекси])
    у_гривнях = np.exp(масштабувач.inverse_transform(центри)[:, 0])
    рядки.append({"старт": назва, "інерція": round(інерція, 2),
                  "розміри": np.bincount(мітки).tolist(),
                  "центри, грн": np.sort(у_гривнях).round(0).tolist()})

таблиця_стартів = pd.DataFrame(рядки)
print(таблиця_стартів.to_string(index=False))

найкраща = таблиця_стартів["інерція"].min()
найгірша = таблиця_стартів["інерція"].max()
print(f"\nнайгірший старт гірший за найкращий на {найгірша - найкраща:.2f}"
      f" ({100 * (найгірша - найкраща) / найкраща:.1f} %)")
старт  інерція         розміри               центри, грн
    A   803.43 [481, 264, 355]  [1478.0, 5795.0, 7098.0]
    B   805.58 [445, 313, 342]  [1441.0, 5361.0, 7468.0]
    C   815.74 [520, 312, 268] [2380.0, 3201.0, 12283.0]
    D   880.78 [366, 467, 267]  [1725.0, 2981.0, 8509.0]

найгірший старт гірший за найкращий на 77.35 (9.6 %)

Ті самі дані, та сама функція — і чотири різні відповіді. Проти цього є два прийоми, і обидва вже стоять у scikit-learn за замовчуванням: розумний старт k-means++ і n_init — кілька запусків із вибором найменшої інерції.

Перевіримо, наскільки часто випадковий старт узагалі влучає в найкраще розбиття.

In [11]:
гсч_стартів = np.random.default_rng(42)
інерції = []
for спроба in range(30):
    індекси = гсч_стартів.choice(len(масштабовані), 3, replace=False)
    _, _, інерція = кластеризувати(масштабовані, масштабовані[індекси])
    інерції.append(інерція)

інерції = np.array(інерції)
найменша = інерції.min()
# «влучив» = зупинився в тому самому мінімумі, з точністю до похибки обчислень
влучань = int((інерції < найменша + 1e-6).sum())
print(f"30 випадкових стартів: від {інерції.min():.2f} до {інерції.max():.2f}")
print(f"у найкращий мінімум влучили {влучань} разів із 30")
print("унікальні зупинки:", np.unique(інерції.round(2)))
30 випадкових стартів: від 803.43 до 815.89
у найкращий мінімум влучили 7 разів із 30
унікальні зупинки: [803.43 803.46 805.59 815.74 815.89]

7 · Скільки брати кластерів¶

Тепер два інструменти вибору k — і подивимось, чи вони погодяться.

  • Лікоть: інерція для k від 1 до 10. Вона спадає завжди; шукаємо, де крива з обриву переходить у пологий схил.
  • Силует: для кожної точки рахують a — середню відстань до своїх, b — до найближчих чужих, і беруть (b − a) / max(a, b). Потім усереднюють по всіх точках. Ближче до +1 — краще.

Для кожного k робимо 50 випадкових стартів і лишаємо найкращий: інакше порівнюватимемо не значення k, а те, кому більше пощастило зі стартом.

In [12]:
гсч_кривої = np.random.default_rng(42)
результати = {}
for k in range(2, 11):
    найкраще = None
    for спроба in range(50):
        індекси = гсч_кривої.choice(len(масштабовані), k, replace=False)
        центри, мітки, інерція = кластеризувати(масштабовані, масштабовані[індекси])
        if найкраще is None or інерція < найкраще[2]:
            найкраще = (центри, мітки, інерція)
    результати[k] = найкраще

# при k = 1 центр один — це середнє всієї вибірки
інерція_при_одному = ((масштабовані - масштабовані.mean(axis=0)) ** 2).sum()

крива = pd.DataFrame({
    "k": [1] + list(range(2, 11)),
    "інерція": [інерція_при_одному] + [результати[k][2] for k in range(2, 11)],
})
крива["падіння"] = -крива["інерція"].diff()
крива["силует"] = [np.nan] + [silhouette_score(масштабовані, результати[k][1])
                              for k in range(2, 11)]
print(крива.round(4).to_string(index=False))
 k   інерція   падіння  силует
 1 2200.0000       NaN     NaN
 2 1069.6716 1130.3284  0.4438
 3  803.4317  266.2399  0.3841
 4  573.8607  229.5710  0.3867
 5  480.6250   93.2358  0.3642
 6  392.2461   88.3789  0.3836
 7  347.8271   44.4190  0.4012
 8  309.4528   38.3743  0.3771
 9  275.0980   34.3548  0.3768
10  248.0538   27.0442  0.3788
In [13]:
фігура, (ліворуч, праворуч) = plt.subplots(1, 2, figsize=(11, 3.6))

ліворуч.plot(крива["k"], крива["інерція"], marker="o", color="#c2185b")
ліворуч.set_title("Лікоть: інерція")
ліворуч.set_xlabel("кількість кластерів k")
ліворуч.set_ylabel("інерція (менше = щільніше)")
ліворуч.grid(alpha=0.25)

праворуч.plot(крива["k"][1:], крива["силует"][1:], marker="o", color="#0f766e")
праворуч.set_title("Силует: середній по всіх точках")
праворуч.set_xlabel("кількість кластерів k")
праворуч.set_ylabel("силует (більше = чіткіше)")
праворуч.grid(alpha=0.25)

plt.tight_layout()
plt.show()

найкраще_за_силуетом = int(крива.loc[крива["силует"].idxmax(), "k"])
print("найбільший силует — при k =", найкраще_за_силуетом,
      f"({крива['силует'].max():.4f})")
print("падіння інерції з k=3 на k=4:", round(крива.loc[3, "падіння"], 1),
      "· з k=4 на k=5:", round(крива.loc[4, "падіння"], 1))
No description has been provided for this image
найбільший силует — при k = 2 (0.4438)
падіння інерції з k=3 на k=4: 229.6 · з k=4 на k=5: 93.2

Два чесно пораховані критерії показують на різні числа: силует найбільший при k = 2, а крива інерції помітно вирівнюється після k = 4. Це нормально, і третьої метрики як арбітра шукати не треба.

Останнє слово має задача. Нам потрібні сегменти, які можна назвати словами й показати менеджерові дошки, тому далі беремо k = 3 — і записуємо це як рішення, а не як відкриття.

In [14]:
центри3, мітки3, інерція3 = результати[3]

# упорядкуємо кластери від найдешевшого до найдорожчого, щоб номери не стрибали
порядок = np.argsort(центри3[:, 0])
перепис = {старий: новий for новий, старий in enumerate(порядок)}
мітки3 = np.array([перепис[м] for м in мітки3])
центри3 = центри3[порядок]

у_натуральних = масштабувач.inverse_transform(центри3)

сегменти = pd.DataFrame({
    "сегмент": ["дешеві старі", "міцна середина", "свіжі дорогі"],
    "оголошень": [int((мітки3 == c).sum()) for c in range(3)],
    "центр, грн": np.exp(у_натуральних[:, 0]).round(0),
    "центр, рік": у_натуральних[:, 1].round(1),
    "медіанна ціна": [np.median(ціни[мітки3 == c]).round(0) for c in range(3)],
    "медіанний рік": [int(np.median(роки[мітки3 == c])) for c in range(3)],
})
print(сегменти.to_string(index=False))
       сегмент  оголошень  центр, грн  центр, рік  медіанна ціна  медіанний рік
  дешеві старі        355      1478.0      2018.9         1570.0           2019
міцна середина        264      5795.0      2018.6         5160.0           2019
  свіжі дорогі        481      7098.0      2022.7         6890.0           2023

Назви в першій колонці придумали ми. Алгоритм видав лише номери — і це важлива різниця: описати групу словами можна тільки після того, як подивишся, що в ній лежить.

In [15]:
гсч_малюнка = np.random.default_rng(1)   # тільки щоб «розтрусити» точки по вертикалі
кольори = ["#c2185b", "#0f766e", "#c2620f"]
plt.figure(figsize=(8, 4.2))
for c in range(3):
    у_кластері = мітки3 == c
    plt.scatter(ціни[у_кластері], роки[у_кластері] + гсч_малюнка.uniform(-0.25, 0.25, у_кластері.sum()),
                s=7, alpha=0.55, color=кольори[c], label=сегменти["сегмент"][c])
plt.scatter(np.exp(у_натуральних[:, 0]), у_натуральних[:, 1],
            s=180, marker="D", color="white", edgecolor="black", zorder=5, label="центри")
plt.xscale("log")
plt.xlabel("ціна, грн (логарифмічна шкала)")
plt.ylabel("рік випуску")
plt.title("Три сегменти дошки")
plt.legend(fontsize=8)
plt.grid(alpha=0.2)
plt.show()
print("рік трохи «розтрушено» по вертикалі — інакше 1 100 точок лягли б вісьмома лініями")
No description has been provided for this image
рік трохи «розтрушено» по вертикалі — інакше 1 100 точок лягли б вісьмома лініями

8 · Що буває без масштабування¶

Той самий алгоритм, ті самі дані — але ознаки в їхніх власних одиницях: ціна в гривнях, рік у роках. Розкид ціни більший за розкид року в 3 765 разів, тож у формулі відстані рік просто зникне.

In [16]:
сирі_ознаки = np.c_[ціни, роки]

гсч_сирих = np.random.default_rng(7)
найкраще_сире = None
for спроба in range(30):
    індекси = гсч_сирих.choice(len(сирі_ознаки), 3, replace=False)
    центри, мітки, інерція = кластеризувати(сирі_ознаки, сирі_ознаки[індекси])
    if найкраще_сире is None or інерція < найкраще_сире[2]:
        найкраще_сире = (центри, мітки, інерція)

сирі_мітки = найкраще_сире[1]
порядок_сирих = np.argsort([np.median(ціни[сирі_мітки == c]) for c in range(3)])
перепис_сирих = {int(старий): новий for новий, старий in enumerate(порядок_сирих)}
сирі_мітки = np.array([перепис_сирих[int(мітка)] for мітка in сирі_мітки])

без_масштабу = pd.DataFrame({
    "кластер": range(3),
    "оголошень": [int((сирі_мітки == c).sum()) for c in range(3)],
    "ціна від": [int(ціни[сирі_мітки == c].min()) for c in range(3)],
    "ціна до": [int(ціни[сирі_мітки == c].max()) for c in range(3)],
    "рік від": [int(роки[сирі_мітки == c].min()) for c in range(3)],
    "рік до": [int(роки[сирі_мітки == c].max()) for c in range(3)],
})
print(без_масштабу.to_string(index=False))
print("\nзбіг із масштабованим розбиттям (ARI):",
      round(adjusted_rand_score(сирі_мітки, мітки3), 4))
 кластер  оголошень  ціна від  ціна до  рік від  рік до
       0        956       130    11390     2017    2024
       1        136     11660    45320     2017    2024
       2          8     52620    95000     2017    2024

збіг із масштабованим розбиттям (ARI): -0.0038

Кожен із трьох кластерів містить оголошення всіх років з 2017-го по 2024-й: рік не вплинув ні на що. А сама ціна поділена безглуздо — 956 оголошень в одній купі й вісім дорогих викидів в окремому кластері. Збіг двох розбиттів між собою — на рівні випадкового.

Правило без винятків: будь-який метод, який рахує відстані, вимагає приведення ознак до спільної шкали (тема 09).

9 · Силует по точках: кому в кластері незручно¶

Середній силует ховає найцікавіше. Подивимось на розподіл силуету всередині кожного з трьох сегментів — і на те, скільки точок сидять не в тій групі (силует нижчий за нуль).

In [17]:
силует_точок = silhouette_samples(масштабовані, мітки3)

розбір = pd.DataFrame({
    "сегмент": сегменти["сегмент"],
    "оголошень": [int((мітки3 == c).sum()) for c in range(3)],
    "середній силует": [силует_точок[мітки3 == c].mean().round(3) for c in range(3)],
    "силует < 0": [int((силует_точок[мітки3 == c] < 0).sum()) for c in range(3)],
})
print(розбір.to_string(index=False))
print(f"\nсередній силует усього розбиття: {силует_точок.mean():.4f}")
print(f"оголошень із відʼємним силуетом:  {int((силует_точок < 0).sum())} із {len(силует_точок)}")
       сегмент  оголошень  середній силует  силует < 0
  дешеві старі        355            0.358           8
міцна середина        264            0.385           0
  свіжі дорогі        481            0.403           2

середній силует усього розбиття: 0.3841
оголошень із відʼємним силуетом:  10 із 1100

10 · Найчесніша перевірка: а чи це шахраї?¶

У таблиці все-таки є колонка шахрайське — розмітка, з якою ми розбирали дошку ще в темі 08. Алгоритм її не бачив: він працював тільки з ціною й роком. Подивімось, чи має знайдене розбиття хоч якийсь стосунок до неї.

In [18]:
шахрайське_справді = ринок["шахрайське"].to_numpy()

порівняння = pd.DataFrame({
    "сегмент": сегменти["сегмент"],
    "усього": [int((мітки3 == c).sum()) for c in range(3)],
    "чесних": [int(((мітки3 == c) & (шахрайське_справді == 0)).sum()) for c in range(3)],
    "шахрайських": [int(((мітки3 == c) & (шахрайське_справді == 1)).sum()) for c in range(3)],
})
порівняння["частка шахрайських"] = (100 * порівняння["шахрайських"] / порівняння["усього"]).round(1)
print(порівняння.to_string(index=False))
print(f"\nпо дошці загалом: {100 * шахрайське_справді.mean():.1f} %")
print(f"збіг кластерів із розміткою (ARI): {adjusted_rand_score(шахрайське_справді, мітки3):.4f}")
       сегмент  усього  чесних  шахрайських  частка шахрайських
  дешеві старі     355     282           73                20.6
міцна середина     264     246           18                 6.8
  свіжі дорогі     481     433           48                10.0

по дошці загалом: 12.6 %
збіг кластерів із розміткою (ARI): 0.0072

Сигнал є — і він слабкий: у найдешевшому сегменті шахрайських утричі більше, ніж у «міцній середині». Спокуса очевидна: оголосити цей сегмент детектором шахрайства. Порахуймо, що з цього вийде, тими самими метриками, що й у темі 05.

In [19]:
# «детектор»: усе, що потрапило в найдешевший сегмент, вважаємо підозрілим
прогноз_за_кластером = (мітки3 == 0).astype(int)

print(f"precision: {precision_score(шахрайське_справді, прогноз_за_кластером):.3f}")
print(f"recall:    {recall_score(шахрайське_справді, прогноз_за_кластером):.3f}")
print(f"F1:        {f1_score(шахрайське_справді, прогноз_за_кластером):.3f}")

# орієнтир: ткнути пальцем у стільки ж рядків навмання, нічого не рахуючи.
# точність такого тику дорівнює частці шахрайських на дошці, а повнота —
# частці позначених рядків. рахуємо ці два числа прямо, без жеребкування.
позначено = int(прогноз_за_кластером.sum())
точність_навмання = шахрайське_справді.mean()
повнота_навмання = позначено / len(шахрайське_справді)
f1_навмання = (2 * точність_навмання * повнота_навмання
               / (точність_навмання + повнота_навмання))

print(f"\nдля порівняння — {позначено} рядків навмання:")
print(f"precision: {точність_навмання:.3f}")
print(f"recall:    {повнота_навмання:.3f}")
print(f"F1:        {f1_навмання:.3f}")
precision: 0.206
recall:    0.525
F1:        0.296

для порівняння — 355 рядків навмання:
precision: 0.126
recall:    0.323
F1:        0.182

Ось і вся ціна нашого «детектора»: половину шахраїв упіймано, але з кожних пʼяти піднятих тривог чотири хибні. Проти випадкового тику він виграє (0.206 проти 0.126), і виграє мало. Скоригований індекс Ренда 0.007 каже те саме одним числом: збіг із розміткою на рівні випадкового.

Скажемо прямо: кластеризація не знайшла шахраїв. Вона й не мала. Шахрайство на цій дошці — не окрема область простору, а співвідношення ціни оголошення з типовою ціною такого самого телефона. В ознаках «ціна» й «рік» цього відношення немає — його треба спершу побудувати руками (тема 10). Далі в курсі, коли дійдемо до моделей із учителем, ми саме так і зробимо.

Кластеризація знаходить те, що є в геометрії даних, а не те, що тобі цікаво. Вона знайшла цінові сегменти — вони там справді є. Шахраїв не знайшла — їх у цій геометрії немає.

11 · Перевірка на стійкість¶

Останній звичай, який варто завести: перш ніж комусь показувати кластери, перевір, чи переживуть вони втрату частини даних. Приберемо 10 % рядків навмання, перезапустимо — і подивимось, чи лишилось розбиття тим самим на решті точок.

In [20]:
гсч_стійкості = np.random.default_rng(42)
збіги = []
for спроба in range(5):
    залишені = гсч_стійкості.choice(len(масштабовані),
                                    size=int(0.9 * len(масштабовані)), replace=False)
    підвибірка = масштабовані[залишені]
    найкраще = None
    for старт in range(10):
        індекси = гсч_стійкості.choice(len(підвибірка), 3, replace=False)
        центри, мітки, інерція = кластеризувати(підвибірка, підвибірка[індекси])
        if найкраще is None or інерція < найкраще[2]:
            найкраще = (центри, мітки, інерція)
    # порівнюємо з повним розбиттям на тих самих рядках
    збіги.append(adjusted_rand_score(мітки3[залишені], найкраще[1]))

print("ARI повного розбиття з розбиттям на 90 % рядків:")
for i, v in enumerate(збіги, 1):
    print(f"  спроба {i}: {v:.4f}")
print(f"\nсередній збіг: {np.mean(збіги):.4f}")
ARI повного розбиття з розбиттям на 90 % рядків:
  спроба 1: 0.9818
  спроба 2: 0.8229
  спроба 3: 0.9942
  спроба 4: 0.9941
  спроба 5: 0.9809

середній збіг: 0.9548

Чотири спроби з пʼяти дали збіг понад 0.98, одна — 0.82: у ній випадково зникло досить оголошень із межі між «дешевими старими» й «міцною серединою», і межа зсунулась. Загалом це добрий результат: сегменти описують дані, а не конкретний запуск. Якби середній збіг вийшов 0.6, показувати це розбиття менеджерові було б рано.


Завдання¶

🟢 Рівень 1 — База¶

Побудуй розбиття на k = 4 тією самою функцією кластеризувати (не забудь про кілька випадкових стартів) і опиши словами, на що саме розпався один із трьох наших сегментів.

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

🟡 Рівень 2 — Плюс¶

Додай третю ознаку — памʼять_гб (не забудь масштабувати всі три разом!) і повтори криву ліктя та силует для k від 2 до 10.

Питання, на які треба відповісти числами: чи змінилось k із найбільшим силуетом? Чи змінились розміри груп при k = 3? Чи стали сегменти зрозумілішими — тобто чи можеш ти й далі назвати кожен словами?

🔴 Рівень 3 — Виклик¶

Реалізуй k-means++ — розумний вибір стартових центрів:

  1. перший центр — випадкове оголошення;
  2. для кожного оголошення порахуй d² — квадрат відстані до найближчого вже обраного центра;
  3. наступний центр витягни випадково з імовірністю, пропорційною цьому d² (стане в пригоді rng.choice(..., p=...));
  4. повторюй, поки центрів не стане k.

Потім порівняй: 30 запусків із випадковим стартом проти 30 запусків із k-means++. Скільки разів кожен спосіб влучив у найкраще розбиття? Наскільки менший розкид інерції?

Зроблено, якщо є два числа «влучань із 30» і два значення розмаху інерції — і одне речення про те, чи варта ця складність результату.