import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier
from sklearn.inspection import permutation_importance

# Ті самі клієнти банку, що в лекції про дерева, але тепер шість ознак замість двох.
# Дві справжні, дві майже-дублікати і дві чисто шумові — рівно як в інтерактиві 3 лекції.

def справжня_межа(стаж):
    """Сходинка, за якою насправді влаштовані дані. Модель її не знає."""
    return np.where(стаж < 9, 7.5, np.where(стаж < 17, 4.8, 2.2))


def згенерувати_клієнтів(скільки, генератор):
    стаж = генератор.uniform(0, 24, скільки)
    витрати = генератор.uniform(0, 10, скільки)

    лишиться = (витрати >= справжня_межа(стаж)).astype(int)
    перевернути = генератор.random(скільки) < 0.10   # 10% шуму в мітках
    лишиться[перевернути] = 1 - лишиться[перевернути]

    # майже-дублікати: несуть ту саму інформацію, але не збігаються з оригіналом
    транзакцій = витрати * 12 + генератор.normal(0, 3, скільки)
    вік_акаунта = стаж + генератор.normal(0, 2, скільки)

    # чистий шум: id сесії — випадкове число з мільйона різних значень,
    # канал звернення — теж шум, але всього двох значень. Різниця буде важлива.
    id_сесії = генератор.uniform(0, 1_000_000, скільки)
    канал = генератор.integers(0, 2, скільки).astype(float)

    ознаки = np.column_stack([стаж, витрати, транзакцій, вік_акаунта, id_сесії, канал])
    return ознаки, лишиться


назви_ознак = ["стаж", "витрати", "транзакцій", "вік акаунта", "id сесії", "канал"]

генератор = np.random.default_rng(61)
X_навч, y_навч = згенерувати_клієнтів(200, генератор)
X_тест, y_тест = згенерувати_клієнтів(1200, генератор)

print(f"навчальна вибірка: {X_навч.shape[0]} клієнтів × {X_навч.shape[1]} ознак")
print(f"тестова вибірка:   {X_тест.shape[0]} клієнтів")
print(f"різних значень «id сесії» серед 200 клієнтів: {len(np.unique(X_навч[:, 4]))}")
print(f"різних значень «канал»    серед 200 клієнтів: {len(np.unique(X_навч[:, 5]))}")

навчальна вибірка: 200 клієнтів × 6 ознак
тестова вибірка:   1200 клієнтів
різних значень «id сесії» серед 200 клієнтів: 200
різних значень «канал»    серед 200 клієнтів: 2

from math import comb


def ймовірність_правильної_більшості(кількість_моделей, точність_однієї):
    """Ймовірність, що більш ніж половина з M незалежних моделей дасть правильну відповідь."""
    ймовірність = 0.0
    # більшість — це строго більше половини голосів
    мінімум_голосів = кількість_моделей // 2 + 1
    for правильних in range(мінімум_голосів, кількість_моделей + 1):
        ймовірність += (comb(кількість_моделей, правильних)
                        * точність_однієї ** правильних
                        * (1 - точність_однієї) ** (кількість_моделей - правильних))
    return ймовірність


рядки = []
for точність in [0.45, 0.50, 0.60, 0.70]:
    рядок = {"точність одного дерева": точність}
    for моделей in [1, 11, 51, 201]:
        рядок[f"M={моделей}"] = ймовірність_правильної_більшості(моделей, точність)
    рядки.append(рядок)

print(pd.DataFrame(рядки).to_string(index=False, float_format=lambda v: f"{v:.4f}"))

 точність одного дерева    M=1   M=11   M=51  M=201
                 0.4500 0.4500 0.3669 0.2359 0.0774
                 0.5000 0.5000 0.5000 0.5000 0.5000
                 0.6000 0.6000 0.7535 0.9265 0.9979
                 0.7000 0.7000 0.9218 0.9986 1.0000

def частка_унікальних(розмір_вибірки, повторів, генератор):
    """Яка частка обʼєктів потрапляє в бутстреп-вибірку хоча б раз. Усереднено по повторах."""
    частки = []
    for _ in range(повторів):
        # витягуємо n індексів з n можливих — з поверненням, тому бувають дублікати
        індекси = генератор.integers(0, розмір_вибірки, розмір_вибірки)
        частки.append(len(np.unique(індекси)) / розмір_вибірки)
    return float(np.mean(частки))


генератор_бутстрепу = np.random.default_rng(7)

print("  n     симуляція   теорія 1-(1-1/n)^n   різниця")
for n in [10, 20, 50, 100, 500, 2000]:
    симуляція = частка_унікальних(n, 400, генератор_бутстрепу)
    теорія = 1 - (1 - 1 / n) ** n
    print(f"{n:5d}     {симуляція:.4f}       {теорія:.4f}          {abs(симуляція - теорія):.4f}")

межа = 1 - np.exp(-1)
print(f"\nграниця при n → ∞: 1 − e⁻¹ = {межа:.4f}")

симуляція_200 = частка_унікальних(200, 800, генератор_бутстрепу)
assert np.allclose(симуляція_200, 1 - (1 - 1 / 200) ** 200, atol=0.005), "симуляція розійшлася з теорією!"
print(f"✅ при n = 200 симуляція дала {симуляція_200:.4f}, теорія — {1 - (1 - 1 / 200) ** 200:.4f}")

  n     симуляція   теорія 1-(1-1/n)^n   різниця
   10     0.6500       0.6513          0.0013
   20     0.6406       0.6415          0.0009
   50     0.6347       0.6358          0.0011
  100     0.6329       0.6340          0.0011
  500     0.6310       0.6325          0.0015

 2000     0.6318       0.6322          0.0005

границя при n → ∞: 1 − e⁻¹ = 0.6321
✅ при n = 200 симуляція дала 0.6328, теорія — 0.6330

def навчити_беггінг(X, y, кількість_дерев, генератор, ознак_на_розріз=None):
    """Бутстреп + дерево на кожній вибірці. Повертає список навчених дерев і їхні OOB-індекси."""
    дерева = []
    невикористані = []
    n = len(y)

    for номер in range(кількість_дерев):
        індекси = генератор.integers(0, n, n)          # бутстреп: n з n з поверненням
        поза_вибіркою = np.setdiff1d(np.arange(n), індекси)   # ті, кого не витягли жодного разу

        # random_state різний для кожного дерева, щоб не було прихованої однаковості
        дерево = DecisionTreeClassifier(max_features=ознак_на_розріз, random_state=номер)
        дерево.fit(X[індекси], y[індекси])

        дерева.append(дерево)
        невикористані.append(поза_вибіркою)

    return дерева, невикористані


def проголосувати(дерева, X):
    """Середня ймовірність класу по всіх деревах — і мітка більшості з неї."""
    ймовірності = np.mean([дерево.predict_proba(X) for дерево in дерева], axis=0)
    return ймовірності.argmax(axis=1), ймовірності


генератор_беггінгу = np.random.default_rng(2024)
беггінг, _ = навчити_беггінг(X_навч, y_навч, 100, генератор_беггінгу)

одне_дерево = DecisionTreeClassifier(random_state=0).fit(X_навч, y_навч)
прогноз_беггінгу, _ = проголосувати(беггінг, X_тест)

print(f"одне дерево, тест     : {одне_дерево.score(X_тест, y_тест):.4f}")
print(f"беггінг зі 100 дерев  : {np.mean(прогноз_беггінгу == y_тест):.4f}")

одне дерево, тест     : 0.7925
беггінг зі 100 дерев  : 0.8692

беггінг_з_підмножиною, _ = навчити_беггінг(
    X_навч, y_навч, 100, np.random.default_rng(2024), ознак_на_розріз="sqrt")

прогноз_з_підмножиною, _ = проголосувати(беггінг_з_підмножиною, X_тест)

ліс = RandomForestClassifier(n_estimators=100, random_state=0).fit(X_навч, y_навч)

print(f"беггінг, усі 6 ознак у вузлі      : {np.mean(прогноз_беггінгу == y_тест):.4f}")
print(f"беггінг, √6 ≈ 2 ознаки у вузлі    : {np.mean(прогноз_з_підмножиною == y_тест):.4f}")
print(f"RandomForestClassifier (100 дерев): {ліс.score(X_тест, y_тест):.4f}")
print("\nНаш «беггінг з підмножиною ознак» — це і є випадковий ліс.")
print("Різниця з бібліотечним — тільки в тому, які саме випадкові числа випали.")

беггінг, усі 6 ознак у вузлі      : 0.8692
беггінг, √6 ≈ 2 ознаки у вузлі    : 0.8733
RandomForestClassifier (100 дерев): 0.8708

Наш «беггінг з підмножиною ознак» — це і є випадковий ліс.
Різниця з бібліотечним — тільки в тому, які саме випадкові числа випали.

# rf.estimators_ — це звичайні DecisionTreeClassifier, їх можна викликати окремо
наші_ймовірності = np.mean([дерево.predict_proba(X_тест) for дерево in ліс.estimators_], axis=0)
ймовірності_sklearn = ліс.predict_proba(X_тест)

print(f"дерев у лісі: {len(ліс.estimators_)}")
print(f"наші ймовірності, перші 3 клієнти:\n{наші_ймовірності[:3]}")
print(f"sklearn,        перші 3 клієнти:\n{ймовірності_sklearn[:3]}")
print(f"максимальна різниця: {np.max(np.abs(наші_ймовірності - ймовірності_sklearn)):.2e}")

assert np.allclose(наші_ймовірності, ймовірності_sklearn), "усереднення розійшлося!"
print("\n✅ прогноз лісу = середнє прогнозів його дерев, до останнього знаку")

дерев у лісі: 100
наші ймовірності, перші 3 клієнти:
[[0.96 0.04]
 [0.66 0.34]
 [0.79 0.21]]
sklearn,        перші 3 клієнти:
[[0.96 0.04]
 [0.66 0.34]
 [0.79 0.21]]
максимальна різниця: 0.00e+00

✅ прогноз лісу = середнє прогнозів його дерев, до останнього знаку

def середня_розбіжність(прогнози):
    """Середня по всіх парах запусків частка обʼєктів, де прогнози розійшлися."""
    розбіжності = []
    for i in range(len(прогнози)):
        for j in range(i + 1, len(прогнози)):
            розбіжності.append(np.mean(прогнози[i] != прогнози[j]))
    return float(np.mean(розбіжності))


генератор_підвибірок = np.random.default_rng(11)
прогнози_дерева = []
прогнози_лісу = []
рядки_стабільності = []

for спроба in range(8):
    # прибираємо чверть клієнтів навмання — «трохи інші дані»
    індекси = генератор_підвибірок.choice(len(y_навч), size=150, replace=False)

    дерево = DecisionTreeClassifier(random_state=0).fit(X_навч[індекси], y_навч[індекси])
    ліс_спроби = RandomForestClassifier(n_estimators=25, random_state=0)
    ліс_спроби.fit(X_навч[індекси], y_навч[індекси])

    прогнози_дерева.append(дерево.predict(X_тест))
    прогнози_лісу.append(ліс_спроби.predict(X_тест))

    рядки_стабільності.append({
        "спроба": спроба + 1,
        "test · дерево": дерево.score(X_тест, y_тест),
        "test · ліс": ліс_спроби.score(X_тест, y_тест),
    })

таблиця_стабільності = pd.DataFrame(рядки_стабільності)
print(таблиця_стабільності.to_string(index=False, float_format=lambda v: f"{v:.3f}"))

 спроба  test · дерево  test · ліс
      1          0.748       0.843
      2          0.782       0.839
      3          0.785       0.846
      4          0.783       0.852
      5          0.795       0.856
      6          0.785       0.857
      7          0.805       0.865
      8          0.779       0.849

розбіжність_дерева = середня_розбіжність(прогнози_дерева)
розбіжність_лісу = середня_розбіжність(прогнози_лісу)

print("СТАБІЛЬНІСТЬ (менше — краще)")
print(f"  дерево: два запуски розходяться на {розбіжність_дерева:.1%} тестових клієнтів")
print(f"  ліс   : два запуски розходяться на {розбіжність_лісу:.1%} тестових клієнтів")
print(f"  ліс стабільніший у {розбіжність_дерева / розбіжність_лісу:.1f} раза\n")

print("РОЗКИД ТОЧНОСТІ по восьми запусках")
print(f"  дерево: від {таблиця_стабільності['test · дерево'].min():.3f} "
      f"до {таблиця_стабільності['test · дерево'].max():.3f}, "
      f"std = {таблиця_стабільності['test · дерево'].std():.4f}")
print(f"  ліс   : від {таблиця_стабільності['test · ліс'].min():.3f} "
      f"до {таблиця_стабільності['test · ліс'].max():.3f}, "
      f"std = {таблиця_стабільності['test · ліс'].std():.4f}")

СТАБІЛЬНІСТЬ (менше — краще)
  дерево: два запуски розходяться на 19.0% тестових клієнтів
  ліс   : два запуски розходяться на 6.7% тестових клієнтів
  ліс стабільніший у 2.8 раза

РОЗКИД ТОЧНОСТІ по восьми запусках
  дерево: від 0.748 до 0.805, std = 0.0165
  ліс   : від 0.839 до 0.865, std = 0.0084

ліс_oob = RandomForestClassifier(n_estimators=200, oob_score=True, random_state=0)
ліс_oob.fit(X_навч, y_навч)

n = len(y_навч)
сума_голосів = np.zeros((n, 2))   # накопичуємо ймовірності класів
скільки_дерев_голосувало = np.zeros(n)

# estimators_samples_ каже, які індекси потрапили в бутстреп кожного дерева
for дерево, індекси_бутстрепу in zip(ліс_oob.estimators_, ліс_oob.estimators_samples_):
    поза_вибіркою = np.setdiff1d(np.arange(n), індекси_бутстрепу)
    # це дерево цих клієнтів не бачило — його голос про них чесний
    сума_голосів[поза_вибіркою] += дерево.predict_proba(X_навч[поза_вибіркою])
    скільки_дерев_голосувало[поза_вибіркою] += 1

наш_oob_прогноз = сума_голосів.argmax(axis=1)
наша_oob_точність = float(np.mean(наш_oob_прогноз == y_навч))

print(f"кожен клієнт отримав від {int(скільки_дерев_голосувало.min())} "
      f"до {int(скільки_дерев_голосувало.max())} чесних голосів "
      f"(у середньому {скільки_дерев_голосувало.mean():.1f} з 200)")
print(f"\nнаша OOB-точність   : {наша_oob_точність:.6f}")
print(f"sklearn oob_score_  : {ліс_oob.oob_score_:.6f}")

assert np.allclose(наша_oob_точність, ліс_oob.oob_score_), "OOB розійшлася!"
print("\n✅ збігається — жодної магії всередині oob_score_ немає")

кожен клієнт отримав від 56 до 100 чесних голосів (у середньому 73.9 з 200)

наша OOB-точність   : 0.800000
sklearn oob_score_  : 0.800000

✅ збігається — жодної магії всередині oob_score_ немає

import warnings

кількості_дерев = [5, 10, 20, 40, 80, 150, 300]
oob_криві = []

for скільки in кількості_дерев:
    ліс_кроку = RandomForestClassifier(n_estimators=скільки, oob_score=True, random_state=0)
    # при 5-10 деревах частина клієнтів не отримує жодного чесного голосу,
    # і sklearn чесно про це попереджає. Нам ця межа й потрібна — глушимо попередження.
    with warnings.catch_warnings():
        warnings.simplefilter("ignore")
        ліс_кроку.fit(X_навч, y_навч)
    oob_криві.append({
        "дерев": скільки,
        "OOB": ліс_кроку.oob_score_,
        "test": ліс_кроку.score(X_тест, y_тест),
        "розбіжність": abs(ліс_кроку.oob_score_ - ліс_кроку.score(X_тест, y_тест)),
    })

таблиця_oob = pd.DataFrame(oob_криві)
print(таблиця_oob.to_string(index=False, float_format=lambda v: f"{v:.4f}"))

# OOB рахується всього на 200 обʼєктах — у неї є власна похибка, і чимала
похибка_oob = np.sqrt(0.8 * 0.2 / len(y_навч))
print(f"\nстандартна похибка OOB-оцінки на {len(y_навч)} обʼєктах: ±{похибка_oob:.3f}")

 дерев    OOB   test  розбіжність
     5 0.6250 0.8083       0.1833
    10 0.6900 0.8300       0.1400
    20 0.7750 0.8642       0.0892
    40 0.8150 0.8600       0.0450
    80 0.8000 0.8675       0.0675
   150 0.8050 0.8717       0.0667
   300 0.8050 0.8750       0.0700

стандартна похибка OOB-оцінки на 200 обʼєктах: ±0.028

fig, ax = plt.subplots(figsize=(9, 4.5))

ax.plot(таблиця_oob["дерев"], таблиця_oob["OOB"], marker="o", lw=2,
        color="crimson", label="OOB (усередині навчальних даних)")
ax.plot(таблиця_oob["дерев"], таблиця_oob["test"], marker="o", lw=2,
        color="teal", label="test (1200 нових клієнтів)")

ax.set_xscale("log")
ax.set_xlabel("кількість дерев (лог. шкала)")
ax.set_ylabel("частка правильних відповідей")
ax.set_title("OOB проти тестової вибірки")
ax.legend()
ax.grid(alpha=.25)
plt.tight_layout()
plt.show()

ліс_важливостей = RandomForestClassifier(n_estimators=300, random_state=0)
ліс_важливостей.fit(X_навч, y_навч)

mdi = pd.Series(ліс_важливостей.feature_importances_, index=назви_ознак)
print("MDI (feature_importances_), сума = 1:")
print(mdi.sort_values(ascending=False).to_string(float_format=lambda v: f"{v:.4f}"))

print(f"\n⚠️  «id сесії» — випадкове число, яке не несе жодної інформації, —")
print(f"    отримало {mdi['id сесії']:.1%} важливості.")
print(f"    «канал» — теж чистий шум, але лише двох значень — отримав {mdi['канал']:.1%}.")
print(f"    Різниця між ними {mdi['id сесії'] / mdi['канал']:.0f}× — і вона повністю")
print(f"    пояснюється кількістю різних значень, а не корисністю ознаки.")

MDI (feature_importances_), сума = 1:
витрати       0.2780
транзакцій    0.2419
стаж          0.1937
вік акаунта   0.1669
id сесії      0.1052
канал         0.0143

⚠️  «id сесії» — випадкове число, яке не несе жодної інформації, —
    отримало 10.5% важливості.
    «канал» — теж чистий шум, але лише двох значень — отримав 1.4%.
    Різниця між ними 7× — і вона повністю
    пояснюється кількістю різних значень, а не корисністю ознаки.

# спершу зробимо це руками для однієї ознаки, щоб було видно, що всередині
базова_точність = ліс_важливостей.score(X_тест, y_тест)
генератор_перемішування = np.random.default_rng(5)

X_зіпсований = X_тест.copy()
номер_id_сесії = назви_ознак.index("id сесії")
# перемішуємо колонку між обʼєктами: значення ті самі, звʼязок з міткою зруйновано
X_зіпсований[:, номер_id_сесії] = генератор_перемішування.permutation(X_зіпсований[:, номер_id_сесії])

print(f"точність лісу на цілих даних          : {базова_точність:.4f}")
print(f"точність із перемішаним «id сесії»    : {ліс_важливостей.score(X_зіпсований, y_тест):.4f}")
падіння = базова_точність - ліс_важливостей.score(X_зіпсований, y_тест)
print(f"падіння (це і є permutation importance): {падіння:+.4f}")
print("\nЧисло вийшло майже нульове й навіть трохи відʼємне: зіпсувавши цю колонку,")
print("ми нічому не нашкодили. Відʼємне значення — звичайний шум вимірювання,")
print("а не «ознака заважає моделі».")

точність лісу на цілих даних          : 0.8750
точність із перемішаним «id сесії»    : 0.8758

падіння (це і є permutation importance): -0.0008

Число вийшло майже нульове й навіть трохи відʼємне: зіпсувавши цю колонку,
ми нічому не нашкодили. Відʼємне значення — звичайний шум вимірювання,
а не «ознака заважає моделі».

# те саме для всіх ознак одразу, з десятьма повторами — бо перемішування випадкове
результат_pi = permutation_importance(ліс_важливостей, X_тест, y_тест,
                                      n_repeats=10, random_state=0)

порівняння = pd.DataFrame({
    "MDI (train)": ліс_важливостей.feature_importances_,
    "PI (test)": результат_pi.importances_mean,
    "PI · розкид": результат_pi.importances_std,
}, index=назви_ознак).sort_values("PI (test)", ascending=False)

print(порівняння.to_string(float_format=lambda v: f"{v:.4f}"))
print(f"\n«id сесії»: MDI = {ліс_важливостей.feature_importances_[4]:.4f}, "
      f"PI = {результат_pi.importances_mean[4]:.4f}")
print("Permutation importance майже нульова — і вона права.")

             MDI (train)  PI (test)  PI · розкид
витрати           0.2780     0.1774       0.0109
стаж              0.1937     0.1195       0.0093
транзакцій        0.2419     0.1078       0.0086
вік акаунта       0.1669     0.0406       0.0033
id сесії          0.1052     0.0015       0.0026
канал             0.0143     0.0006       0.0016

«id сесії»: MDI = 0.1052, PI = 0.0015
Permutation importance майже нульова — і вона права.

fig, ax = plt.subplots(figsize=(10, 4.5))

позиції = np.arange(len(назви_ознак))
ширина = 0.38

# MDI нормована на суму 1, PI — у пунктах точності. Щоб порівняти форму,
# нормуємо обидві на власний максимум: нас цікавить порядок, а не абсолют.
mdi_норм = ліс_важливостей.feature_importances_ / ліс_важливостей.feature_importances_.max()
pi_норм = результат_pi.importances_mean / результат_pi.importances_mean.max()

ax.bar(позиції - ширина / 2, mdi_норм, ширина, color="crimson", label="MDI (на навчальних даних)")
ax.bar(позиції + ширина / 2, pi_норм, ширина, color="teal", label="Permutation (на тестових)")

ax.set_xticks(позиції)
ax.set_xticklabels(назви_ознак, rotation=15)
ax.set_ylabel("важливість, нормована на максимум")
ax.set_title("Дві важливості на одному лісі: дві останні ознаки — чистий шум")
ax.axhline(0, color="gray", lw=1)
ax.legend()
ax.grid(alpha=.25, axis="y")
plt.tight_layout()
plt.show()

# лишаємо тільки стаж, витрати й обидва шуми — без майже-дублікатів
без_дублікатів = [0, 1, 4, 5]
назви_без_дублікатів = [назви_ознак[i] for i in без_дублікатів]

ліс_без_дублікатів = RandomForestClassifier(n_estimators=300, random_state=0)
ліс_без_дублікатів.fit(X_навч[:, без_дублікатів], y_навч)

порівняння_дублікатів = pd.DataFrame({
    "MDI з дублікатами": [mdi["стаж"], mdi["витрати"]],
    "MDI без дублікатів": [ліс_без_дублікатів.feature_importances_[0],
                           ліс_без_дублікатів.feature_importances_[1]],
}, index=["стаж", "витрати"])

print(порівняння_дублікатів.to_string(float_format=lambda v: f"{v:.4f}"))
print(f"\nточність з 6 ознаками : {ліс_важливостей.score(X_тест, y_тест):.4f}")
print(f"точність з 4 ознаками : {ліс_без_дублікатів.score(X_тест[:, без_дублікатів], y_тест):.4f}")
print("\nВажливість оригіналів підскочила, хоча якість моделі майже не змінилась:")
print("дублікати нічого не додавали, зате «крали» частку важливості.")
print("Тому важливості ознак — інструмент для гіпотез, а не для остаточних висновків.")

         MDI з дублікатами  MDI без дублікатів
стаж                0.1937              0.3464
витрати             0.2780              0.4649

точність з 6 ознаками : 0.8750
точність з 4 ознаками : 0.8625

Важливість оригіналів підскочила, хоча якість моделі майже не змінилась:
дублікати нічого не додавали, зате «крали» частку важливості.
Тому важливості ознак — інструмент для гіпотез, а не для остаточних висновків.

Практика 12 · Випадковий ліс¶

1. Теорема Кондорсе: коли голосування допомагає¶

2. Бутстреп: звідки береться 63.2%¶

3. Беггінг вручну¶

4. Ліс — це рівно середнє своїх дерев¶

5. Наскільки ліс стабільніший — у числах¶

6. OOB-оцінка: валідація задарма¶

7. Пастка MDI: як шум отримує 10% важливості¶

Друга пастка: корельовані ознаки¶

💻 Завдання¶

🟢 Рівень 1 — разом¶

🟡 Рівень 2 — самостійно¶

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

🧪 Самоперевірка¶