Практика · PCA — метод головних компонент¶
Лекція: lecture.html · Тест: quiz.html · Домашнє: homework.md
Той самий приклад, що й у лекції: дошка оголошень про вживані телефони з теми 08. Сім числових ознак, які неможливо намалювати разом.
Що ми зробимо:
- зберемо ту саму таблицю й дістанемо з неї сім числових ознак;
- знайдемо пару колонок, які майже дублюють одна одну;
- порахуємо PCA вручну на NumPy — центрування, коваріаційна матриця,
np.linalg.eigh; - звіримо результат із
sklearn.decomposition.PCA— числа мають зійтися; - повторимо ручний приклад із шести точок із лекції й перевіримо власні числа 28 і 4;
- подивимось на пояснену дисперсію: скільки компонент треба для 90 %;
- намалюємо дошку у двох компонентах і пофарбуємо за
шахрайське, якого PCA не бачив; - побачимо, у що перетворюється перша компонента без
StandardScaler; - поміряємо витік, який дає PCA, порахований до поділу даних.
Зерно генератора зафіксовано (np.random.default_rng(42)), тож числа збігатимуться
з лекцією до цифри.
In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from sklearn.model_selection import train_test_split
# зерно фіксує всю випадковість: у тебе вийдуть точно ті самі числа, що в лекції
rng = np.random.default_rng(42)
pd.set_option("display.width", 130)
pd.set_option("display.max_columns", 14)
np.set_printoptions(suppress=True, linewidth=130)
print("numpy", np.__version__, "· pandas", pd.__version__)
numpy 2.5.1 · pandas 2.3.3
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)
assert дошка.shape == (1212, 8), "форма розійшлась із темою 22"
print("таблиця як у темі 08:", дошка.shape)
таблиця як у темі 08: (1212, 8)
2 · Сім числових ознак¶
PCA працює тільки з числами й тільки без пропусків. Тому робимо чотири речі:
- прибираємо суфікс «ГБ» і повні дублікати рядків — це та сама безумовна чистка, що й у практиці теми 08;
- перетворюємо стан на бал від 1 до 4 (порядкова шкала, тема 09);
- додаємо дві оцінки «скільки цей телефон приблизно коштує»: медіану реальних цін для пари «модель + рік» і оцінку з каталогу (ціна нового мінус 18 % за кожен рік);
- залишаємо тільки оголошення з відомою ціною — їх 1 100 із 1 200.
Колонку скарг не беремо: це витік, знайдений ще в
темі 08, — скарги зʼявляються вже після публікації.
In [5]:
дошка["памʼять_гб"] = pd.to_numeric(дошка["памʼять_гб"].str.replace(" ГБ", "", regex=False))
дошка = дошка.drop_duplicates().reset_index(drop=True)
# стан — порядкова шкала: чим більше, тим кращий апарат
бали_стану = {"задовільне": 1, "добре": 2, "дуже добре": 3, "нове": 4}
# пропуск у стані заповнюємо найчастішим балом, щоб не викидати рядок цілком
дошка["стан_бал"] = дошка["стан"].map(бали_стану).fillna(2).astype(int)
# без ціни неможливо порахувати жодну з цінових ознак, тож ці рядки відкладаємо
оголошення = дошка.dropna(subset=["ціна"]).reset_index(drop=True)
# перша оцінка вартості: медіана реальних цін усередині групи «модель + рік»
оголошення["типова_ціна"] = (оголошення.groupby(["модель", "рік"])["ціна"]
.transform("median"))
# друга оцінка вартості: ціна нового апарата мінус 18 % зношення за кожен рік
оголошення["оцінка_каталогу"] = (оголошення["модель"].map(ціна_нового)
* 0.82 ** (2024 - оголошення["рік"])).round(0)
назви_ознак = ["ціна", "типова_ціна", "оцінка_каталогу", "рік",
"памʼять_гб", "стан_бал", "вік_акаунта"]
X = оголошення[назви_ознак].to_numpy(float)
таргет = оголошення["шахрайське"].to_numpy()
print("матриця ознак:", X.shape)
print("шахрайських оголошень:", int(таргет.sum()),
f"({таргет.mean() * 100:.1f} %)")
print()
print(оголошення[назви_ознак].agg(["mean", "std", "min", "max"]).T.round(1))
матриця ознак: (1100, 7) шахрайських оголошень: 139 (12.6 %)
mean std min max ціна 6498.2 8689.0 130.0 95000.0 типова_ціна 5771.6 5077.0 1020.0 31915.0 оцінка_каталогу 7394.5 6149.9 1296.0 34000.0 рік 2020.5 2.3 2017.0 2024.0 памʼять_гб 169.5 123.8 64.0 512.0 стан_бал 2.3 0.9 1.0 4.0 вік_акаунта 404.7 406.7 4.0 3519.0
3 · Дві колонки, які кажуть те саме¶
типова_ціна й оцінка_каталогу відповідають на одне питання різними способами.
Подивимось, наскільки вони насправді різні.
In [6]:
матриця_кореляцій = оголошення[назви_ознак].corr()
# шукаємо найсильнішу пару, щоб не вибирати її на око
найсильніша_пара, найсильніша_кореляція = None, 0.0
for i in range(len(назви_ознак)):
for j in range(i + 1, len(назви_ознак)):
значення = матриця_кореляцій.iloc[i, j]
if abs(значення) > abs(найсильніша_кореляція):
найсильніша_пара = (назви_ознак[i], назви_ознак[j])
найсильніша_кореляція = значення
print(матриця_кореляцій.round(3))
print()
print("найсильніша пара:", найсильніша_пара, "· кореляція", round(найсильніша_кореляція, 3))
ціна типова_ціна оцінка_каталогу рік памʼять_гб стан_бал вік_акаунта
ціна 1.000 0.596 0.586 0.247 0.210 0.155 -0.030
типова_ціна 0.596 1.000 0.989 0.501 0.020 -0.008 -0.030
оцінка_каталогу 0.586 0.989 1.000 0.527 0.015 -0.028 -0.031
рік 0.247 0.501 0.527 1.000 -0.030 -0.027 -0.032
памʼять_гб 0.210 0.020 0.015 -0.030 1.000 -0.027 -0.025
стан_бал 0.155 -0.008 -0.028 -0.027 -0.027 1.000 0.063
вік_акаунта -0.030 -0.030 -0.031 -0.032 -0.025 0.063 1.000
найсильніша пара: ('типова_ціна', 'оцінка_каталогу') · кореляція 0.989
Кореляція 0.989 означає, що друга колонка майже не додає нового виміру. PCA має це помітити сам — перевіримо наприкінці розділу 5.
4 · PCA руками на NumPy¶
Три дії, і жодної бібліотечної магії:
- стандартизація — кожну колонку до середнього 0 і відхилення 1 (без цього перемагає ознака з найбільшими числами, розділ 8);
- коваріаційна матриця 7×7 —
np.cov; - власні числа й вектори —
np.linalg.eigh.
eigh (а не eig) — тому що коваріаційна матриця симетрична: для таких матриць
є швидший і точніший алгоритм, і він завжди повертає дійсні числа.
Повертає він їх за зростанням, тому список доведеться перевернути.
In [7]:
масштабувальник = StandardScaler().fit(X)
Z = масштабувальник.transform(X) # центровано й поділено на стандартне відхилення
# rowvar=False означає «колонки — це ознаки, рядки — обʼєкти»
коваріація = np.cov(Z, rowvar=False)
print("коваріаційна матриця 7×7 (вона ж матриця кореляцій, бо дані стандартизовані):")
print(np.round(коваріація, 3))
коваріаційна матриця 7×7 (вона ж матриця кореляцій, бо дані стандартизовані): [[ 1.001 0.596 0.586 0.248 0.21 0.155 -0.03 ] [ 0.596 1.001 0.99 0.501 0.02 -0.008 -0.03 ] [ 0.586 0.99 1.001 0.528 0.015 -0.028 -0.031] [ 0.248 0.501 0.528 1.001 -0.03 -0.027 -0.032] [ 0.21 0.02 0.015 -0.03 1.001 -0.027 -0.025] [ 0.155 -0.008 -0.028 -0.027 -0.027 1.001 0.063] [-0.03 -0.03 -0.031 -0.032 -0.025 0.063 1.001]]
In [8]:
власні_числа, власні_вектори = np.linalg.eigh(коваріація)
# eigh віддає власні числа за зростанням — нам потрібен зворотний порядок
порядок = np.argsort(власні_числа)[::-1]
власні_числа = власні_числа[порядок]
# власні вектори стоять у СТОВПЦЯХ, тому переставляємо стовпці, а потім транспонуємо:
# так кожен рядок стане однією компонентою — як у sklearn
компоненти_вручну = власні_вектори[:, порядок].T
частка = власні_числа / власні_числа.sum()
print("власні числа :", власні_числа.round(4))
print("частка розкиду, % :", (частка * 100).round(2))
print("накопичено, % :", (np.cumsum(частка) * 100).round(2))
власні числа : [2.7963 1.1127 1.0755 0.9463 0.6691 0.3962 0.0103] частка розкиду, % : [39.91 15.88 15.35 13.51 9.55 5.66 0.15] накопичено, % : [ 39.91 55.79 71.14 84.65 94.2 99.85 100. ]
5 · Звіряємо з бібліотекою¶
Тепер те саме через sklearn.decomposition.PCA. Числа мають збігтися — з однією
обмовкою про знак.
Про знак. Якщо вектор u задає напрямок найбільшого розкиду, то −u задає той
самий напрямок, просто дивиться в інший бік. Обидва однаково правильні, і різні
бібліотеки (і навіть різні версії однієї) можуть повернути різні. Тому порівнювати
компоненти треба за модулем, а картинку в компонентах не дивує, коли вона
виявляється дзеркальною.
In [9]:
pca = PCA().fit(Z)
# власні числа — це і є explained_variance_ бібліотеки
assert np.allclose(власні_числа, pca.explained_variance_), "власні числа розійшлись!"
# компоненти збігаються з точністю до знака, тому порівнюємо модулі
assert np.allclose(np.abs(компоненти_вручну), np.abs(pca.components_)), "напрямки розійшлись!"
print("✅ ручний розрахунок збігається з sklearn")
print()
print("наші власні числа :", власні_числа.round(4))
print("sklearn :", pca.explained_variance_.round(4))
print()
print("знак першої компоненти: наш", np.sign(компоненти_вручну[0, 0]),
"· sklearn", np.sign(pca.components_[0, 0]))
✅ ручний розрахунок збігається з sklearn наші власні числа : [2.7963 1.1127 1.0755 0.9463 0.6691 0.3962 0.0103] sklearn : [2.7963 1.1127 1.0755 0.9463 0.6691 0.3962 0.0103] знак першої компоненти: наш -1.0 · sklearn 1.0
In [10]:
навантаження = pd.DataFrame(pca.components_.round(3), columns=назви_ознак,
index=[f"PC{i + 1}" for i in range(len(назви_ознак))])
print("з якими вагами вихідні ознаки входять у кожну компоненту:")
print(навантаження)
print()
print("остання компонента цілком складається з двох колонок,")
print("які ми в розділі 3 знайшли як найсхожішу пару:")
print(навантаження.loc["PC7"].sort_values(key=abs, ascending=False).head(3))
з якими вагами вихідні ознаки входять у кожну компоненту:
ціна типова_ціна оцінка_каталогу рік памʼять_гб стан_бал вік_акаунта
PC1 0.438 0.571 0.572 0.387 0.056 0.018 -0.034
PC2 0.379 -0.069 -0.093 -0.281 0.557 0.627 0.247
PC3 0.093 -0.042 -0.040 -0.107 0.653 -0.468 -0.575
PC4 -0.060 0.028 0.037 0.019 0.318 -0.535 0.778
PC5 -0.345 -0.178 -0.150 0.813 0.330 0.239 -0.004
PC6 0.730 -0.378 -0.363 0.312 -0.224 -0.208 0.028
PC7 0.005 -0.702 0.712 -0.025 0.002 0.013 -0.001
остання компонента цілком складається з двох колонок,
які ми в розділі 3 знайшли як найсхожішу пару:
оцінка_каталогу 0.712
типова_ціна -0.702
рік -0.025
Name: PC7, dtype: float64
PC7 — це різниця між типова_ціна й оцінка_каталогу, і на неї припадає 0.15 %
розкиду. PCA знайшов дублікат сам, без жодної підказки.
6 · Шість точок із лекції¶
У лекції ми рахували все руками на шести оголошеннях і двох ознаках. Перевіримо, що арифметика на папері дає ті самі числа, що NumPy.
In [11]:
# ціна й оцінка каталогу шести оголошень, тисячі гривень
шість_точок = np.array([[4.0, 8.0], [8.0, 5.0], [9.0, 8.0],
[13.0, 10.0], [16.0, 9.0], [16.0, 14.0]])
середні = шість_точок.mean(axis=0)
центровані = шість_точок - середні
C = np.cov(центровані, rowvar=False)
# власні числа матриці 2×2 через дискримінант — рівно як у лекції
слід = C[0, 0] + C[1, 1]
визначник = C[0, 0] * C[1, 1] - C[0, 1] ** 2
дискримінант = слід ** 2 - 4 * визначник
лямбда_1 = (слід + np.sqrt(дискримінант)) / 2
лямбда_2 = (слід - np.sqrt(дискримінант)) / 2
print("середні :", середні)
print("матриця C :", C.round(2).tolist())
print(f"слід {слід:.1f} · визначник {визначник:.1f} · дискримінант {дискримінант:.0f}")
print(f"власні числа: {лямбда_1:.1f} і {лямбда_2:.1f} "
f"(перша компонента тримає {лямбда_1 / слід * 100:.1f} %)")
середні : [11. 9.] матриця C : [[23.2, 9.6], [9.6, 8.8]] слід 32.0 · визначник 112.0 · дискримінант 576 власні числа: 28.0 і 4.0 (перша компонента тримає 87.5 %)
In [12]:
# напрямок першої компоненти: розвʼязок (C − λ₁·I)·u = 0
напрямок = np.array([C[0, 1], лямбда_1 - C[0, 0]])
напрямок = напрямок / np.linalg.norm(напрямок)
# проєкція кожної точки на цей напрямок — скалярний добуток
проєкції = центровані @ напрямок
print("одиничний вектор:", напрямок.round(3),
f"· кут {np.degrees(np.arctan2(напрямок[1], напрямок[0])):.1f}°")
print("проєкції шести точок:", проєкції.round(2))
print(f"дисперсія проєкцій: {проєкції.var(ddof=1):.1f} (мала вийти {лямбда_1:.1f})")
# та сама перевірка, що в лекції: власне число дорівнює розкиду проєкцій
assert np.isclose(проєкції.var(ddof=1), лямбда_1), "розкид проєкцій не дорівнює власному числу!"
assert np.isclose(лямбда_1, 28.0) and np.isclose(лямбда_2, 4.0), "числа розійшлись із лекцією!"
print("✅ 28 і 4 — рівно ті числа, що в лекції")
одиничний вектор: [0.894 0.447] · кут 26.6° проєкції шести точок: [-6.71 -4.47 -2.24 2.24 4.47 6.71] дисперсія проєкцій: 28.0 (мала вийти 28.0) ✅ 28 і 4 — рівно ті числа, що в лекції
7 · Скільки компонент лишити¶
Повертаємось до семи ознак. Питання практичне: якщо ми хочемо зберегти 90 % розкиду, скільки компонент доведеться взяти?
In [13]:
накопичено = np.cumsum(pca.explained_variance_ratio_) * 100
# searchsorted знаходить першу позицію, де накопичена частка дотягує до порога
скільки_треба = int(np.searchsorted(накопичено, 90)) + 1
for k, значення in enumerate(накопичено, start=1):
позначка = "← тут проходимо 90 %" if k == скільки_треба else ""
print(f"компонент {k}: накопичено {значення:6.2f} % {позначка}")
print()
print(f"для 90 % потрібно {скільки_треба} компонент із {len(назви_ознак)}")
print(f"дві компоненти, які можна намалювати, тримають {накопичено[1]:.1f} %")
компонент 1: накопичено 39.91 % компонент 2: накопичено 55.79 % компонент 3: накопичено 71.14 % компонент 4: накопичено 84.65 % компонент 5: накопичено 94.20 % ← тут проходимо 90 % компонент 6: накопичено 99.85 % компонент 7: накопичено 100.00 % для 90 % потрібно 5 компонент із 7 дві компоненти, які можна намалювати, тримають 55.8 %
Пʼять із семи — це майже не стиснення. Наша таблиця стискається погано, і це чесна відповідь, а не помилка: сім ознак справді розповідають про різні речі.
8 · Дошка в двох компонентах¶
Тепер картинка, заради якої все й починалось. Фарбуємо точки за колонкою
шахрайське, якої PCA не бачив жодного разу.
In [14]:
координати = pca.transform(Z) # 1 100 × 7, нам потрібні перші два стовпці
# цінова сімʼя моделі: Alfa дешеві, Beta середні, Gamma дорогі
сімʼя = оголошення["модель"].str.split().str[0]
фігура, (зліва, справа) = plt.subplots(1, 2, figsize=(13, 5.5))
for назва, колір in [("Alfa", "#0f766e"), ("Beta", "#c2620f"), ("Gamma", "#c2185b")]:
маска = (сімʼя == назва).to_numpy()
зліва.scatter(координати[маска, 0], координати[маска, 1],
s=9, alpha=0.6, c=колір, label=назва)
зліва.set_title("колір — цінова сімʼя моделі")
зліва.legend()
справа.scatter(координати[таргет == 0, 0], координати[таргет == 0, 1],
s=9, alpha=0.35, c="#8899aa", label="чесні")
справа.scatter(координати[таргет == 1, 0], координати[таргет == 1, 1],
s=14, alpha=0.85, c="#c2185b", label="шахрайські")
справа.set_title("колір — колонка «шахрайське»")
справа.legend()
for вісь in (зліва, справа):
вісь.set_xlabel("PC1")
вісь.set_ylabel("PC2")
вісь.axhline(0, lw=0.8, c="#cccccc")
вісь.axvline(0, lw=0.8, c="#cccccc")
plt.tight_layout()
plt.show()
print("медіана PC1 за ціновою сімʼєю:")
for назва in ["Alfa", "Beta", "Gamma"]:
маска = (сімʼя == назва).to_numpy()
print(f" {назва:6s} {np.median(координати[маска, 0]):+.2f}")
print()
print(f"медіана PC1 у шахрайських: {np.median(координати[таргет == 1, 0]):+.2f}")
print(f"медіана PC1 у чесних : {np.median(координати[таргет == 0, 0]):+.2f}")
медіана PC1 за ціновою сімʼєю: Alfa -1.01 Beta -0.32 Gamma +1.31 медіана PC1 у шахрайських: -0.38 медіана PC1 у чесних : -0.43
Ліворуч структура є: перша компонента впорядкувала оголошення за ціновим класом моделі, хоча колонки «модель» у неї не було. Праворуч структури немає: шахрайські оголошення розсипані так само, як чесні, і медіани в них практично однакові.
Це не поломка. PCA максимізує розкид, а не користь для задачі. Найрозкиданіший напрямок у цій таблиці — «дорогий телефон проти дешевого», його метод і повернув. А шахрая видає відношення ціни до типової — дія нелінійна, і жодна зважена сума ознак її не відтворює.
9 · Без масштабування перша компонента вироджується¶
Порахуємо PCA ще раз, але на сирій матриці X, де ціна виміряна в гривнях,
а стан — у балах від 1 до 4.
In [15]:
pca_без_масштабу = PCA().fit(X)
порівняння = pd.DataFrame({
"стандартне відхилення": X.std(axis=0, ddof=1).round(2),
"вага в PC1 без масштабу": pca_без_масштабу.components_[0].round(6),
"вага в PC1 після StandardScaler": pca.components_[0].round(3),
}, index=назви_ознак)
print(порівняння)
print()
print(f"частка PC1 без масштабу : {pca_без_масштабу.explained_variance_ratio_[0] * 100:.1f} %")
print(f"частка PC1 після StandardScaler: {pca.explained_variance_ratio_[0] * 100:.1f} %")
# сума квадратів ваг усього вектора дорівнює одиниці — подивимось, як вона розподілена
цінові = pca_без_масштабу.components_[0][:3]
решта = pca_без_масштабу.components_[0][3:]
print()
print(f"на три цінові колонки припадає {(цінові ** 2).sum() * 100:.4f} % довжини вектора,")
print(f"на рік, памʼять, стан і вік акаунта разом — {(решта ** 2).sum():.7f}")
стандартне відхилення вага в PC1 без масштабу вага в PC1 після StandardScaler ціна 8688.98 0.757664 0.438 типова_ціна 5076.98 0.417038 0.571 оцінка_каталогу 6149.89 0.502015 0.572 рік 2.31 0.000090 0.387 памʼять_гб 123.78 0.001648 0.056 стан_бал 0.85 0.000007 0.018 вік_акаунта 406.74 -0.001321 -0.034 частка PC1 без масштабу : 79.5 % частка PC1 після StandardScaler: 39.9 % на три цінові колонки припадає 99.9996 % довжини вектора, на рік, памʼять, стан і вік акаунта разом — 0.0000045
Формально метод відпрацював. Змістовно він побудував складну конструкцію, яка означає «ціна»: решти ознак у першій компоненті просто немає.
10 · Витік: PCA, порахований до поділу¶
PCA не бачить таргета — і саме тому здається, що рахувати його на всіх даних безпечно. Це не так: середні, відхилення й напрямки компонент запамʼятовують тестові рядки, і тест перестає бути незнайомим.
Поміряємо ефект прямо. Метрика проста: яку частку розкиду тестових рядків утримують дві компоненти. Робимо це двічі — чесно й з витоком.
In [16]:
# Частка розкиду ТЕСТОВИХ рядків, яку утримують перші компоненти.
# Масштабувальник і PCA вчаться тільки на навчальних рядках — так, як має бути.
def частка_утриманого(навчальні, тестові, скільки_компонент=2):
свій_масштаб = StandardScaler().fit(навчальні)
свій_pca = PCA(n_components=скільки_компонент).fit(свій_масштаб.transform(навчальні))
тест_Z = свій_масштаб.transform(тестові)
відновлено = свій_pca.transform(тест_Z) @ свій_pca.components_
return 1 - ((тест_Z - відновлено) ** 2).sum() / (тест_Z ** 2).sum()
# Те саме, але масштабувальник і PCA бачили всі рядки, зокрема тестові.
def частка_утриманого_з_витоком(усі, тестові, скільки_компонент=2):
спільний_масштаб = StandardScaler().fit(усі)
спільний_pca = PCA(n_components=скільки_компонент).fit(спільний_масштаб.transform(усі))
тест_Z = спільний_масштаб.transform(тестові)
відновлено = спільний_pca.transform(тест_Z) @ спільний_pca.components_
return 1 - ((тест_Z - відновлено) ** 2).sum() / (тест_Z ** 2).sum()
# один поділ дав би шумне число, тому усереднюємо по тридцяти різних поділах
чесні, з_витоком = [], []
for зерно in range(30):
навч_індекси, тест_індекси = train_test_split(
np.arange(len(X)), test_size=0.25, random_state=зерно)
чесні.append(частка_утриманого(X[навч_індекси], X[тест_індекси]))
з_витоком.append(частка_утриманого_з_витоком(X, X[тест_індекси]))
print("усі 1 100 оголошень, дві компоненти, середнє по 30 поділах:")
print(f" чесно : {np.mean(чесні) * 100:.1f} %")
print(f" з витоком : {np.mean(з_витоком) * 100:.1f} %")
усі 1 100 оголошень, дві компоненти, середнє по 30 поділах: чесно : 54.6 % з витоком : 56.1 %
Півтора відсотка — не катастрофа, і саме тому помилку рідко помічають. Але різниця завжди в один бік, і чим менше даних, тим вона більша. Перевіримо це на маленьких вибірках: 40 оголошень, половина в тест, 200 повторів.
In [17]:
генератор = np.random.default_rng(42)
чесні_виміри, виміри_з_витоком = [], []
for повтор in range(200):
# беремо випадкові 40 рядків — маленька вибірка, як буває в реальних задачах
вибірка = X[генератор.choice(len(X), size=40, replace=False)]
навч, тест = train_test_split(np.arange(40), test_size=0.5, random_state=повтор)
чесні_виміри.append(частка_утриманого(вибірка[навч], вибірка[тест]))
виміри_з_витоком.append(частка_утриманого_з_витоком(вибірка, вибірка[тест]))
print("40 оголошень, 200 повторів, дві компоненти:")
print(f" чесно : {np.mean(чесні_виміри) * 100:.1f} %")
print(f" з витоком : {np.mean(виміри_з_витоком) * 100:.1f} %")
print()
розрив = (np.mean(виміри_з_витоком) - np.mean(чесні_виміри)) * 100
print(f"Розрив — {розрив:.1f} відсоткового пункту, і він вигаданий:")
print("PCA припасувався саме до тих рядків, на яких його потім перевіряють.")
40 оголошень, 200 повторів, дві компоненти: чесно : 53.4 % з витоком : 62.8 % Розрив — 9.4 відсоткового пункту, і він вигаданий: PCA припасувався саме до тих рядків, на яких його потім перевіряють.
Правильно так: StandardScaler і PCA кладуться в Pipeline разом із моделлю
(тема 09). Тоді на кожному згині крос-валідації
вони перераховуються заново й тільки на навчальній частині згину.
Завдання¶
🟢 Рівень 1¶
Побудуй графік накопиченої поясненої дисперсії: по горизонталі кількість компонент від 1 до 7, по вертикалі накопичений відсоток, горизонтальна пунктирна лінія на рівні 90 %. Підпиши точку, у якій крива перетинає поріг.
🟡 Рівень 2¶
Викинь із таблиці колонку оцінка_каталогу — ту саму, що майже дублює
типова_ціна. Порахуй PCA на шести ознаках, що лишились, і порівняй із семи:
скільки компонент тепер треба для 90 % і що сталося з останнім власним числом.
Поясни результат одним реченням.
🔴 Рівень 3¶
Заміни три цінові колонки їхніми логарифмами (np.log) і повтори весь розділ 8.
Подивись на PC1–PC2 знову: чи стало видно шахрайські оголошення? Знайди компоненту,
у якій ваги лог_ціна та лог_типова_ціна мають різні знаки — і поясни, чому саме
вона й є «відхилення ціни від типової» з теми 10,
та який у неї внесок у загальний розкид.