Files

1555 lines
69 KiB
Plaintext
Raw Permalink Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
{
"cells": [
{
"cell_type": "markdown",
"id": "7a19baa9",
"metadata": {},
"source": [
"Необходимо задать переменные, которые будут использоваться в ходе работы программы"
]
},
{
"cell_type": "code",
"execution_count": 38,
"id": "cf4829c9",
"metadata": {},
"outputs": [],
"source": [
"# 1) Названия файлов с исходными данными и пустого Excel-файла для вывода программы\n",
"# Если файлы лежат не в одной директории с программой, необходимо указать полный путь к файлам\n",
"# В целом, если конечный файл не создать, то программа его создаст сама, но лучше явно ей указать его.\n",
"INPUT_PATH = 'data.xlsx'\n",
"OUTPUT_PATH = 'output.xlsx'\n",
"# 2) Количество признаков в методе с включением\n",
"CNT_ENTER = 6\n",
"# 3) Количество признаков в методе с исключением\n",
"CNT_REMOVE = 6\n",
"\n",
"# 4) Обучающая выборка в формате словаря, где ключ - номер класса, значение - список номеров строк, принадлежащих этому классу\n",
"TRAIN_SAMPLE = {\n",
" 1: [4,23,28,53,59,69],\n",
" 2: [13,21,34,36,39,67,70,83],\n",
" 3: [14,20,24,40,45,61,65,68,74,78],\n",
" 4: [2,54],\n",
" 5: [18,31],\n",
" 6: [10,62],\n",
" 7: [47,79]\n",
"}\n"
]
},
{
"cell_type": "markdown",
"id": "ed63bb66",
"metadata": {},
"source": [
"## 1 Часть. Импорт и создание обучающей выборки"
]
},
{
"cell_type": "code",
"execution_count": 56,
"id": "74efa4f6",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"data:\n",
" Наименование X1 X2 X3 \\\n",
"0 Алтайский край -0.567424 -0.126713 -1.023246 \n",
"1 Амурская область 0.167298 -0.608339 -0.667051 \n",
"2 Архангельская область без автономного округа 0.246382 0.884700 -0.982147 \n",
"3 Астраханская область -0.560670 0.354912 -0.105358 \n",
"4 Белгородская область -0.116102 -0.030388 0.922129 \n",
"\n",
" X4 X5 X6 X7 X8 X9 \n",
"0 -0.917052 -0.432397 1.006771 -0.250080 -0.730154 0.163593 \n",
"1 -0.769444 0.622505 0.216527 -0.640826 -0.087726 1.593888 \n",
"2 -0.707543 1.135777 1.441228 -0.657606 -0.427835 -0.629321 \n",
"3 -0.221865 1.041252 -0.715829 -0.525759 0.894809 -0.162722 \n",
"4 0.778060 -0.533539 -0.390494 1.025235 -0.163306 -0.827550 \n"
]
}
],
"source": [
"import pandas as pd\n",
"import numpy as np\n",
"from scipy.stats import f\n",
"from sklearn.discriminant_analysis import LinearDiscriminantAnalysis\n",
"from scipy.spatial.distance import mahalanobis\n",
"\n",
"data = pd.read_excel(INPUT_PATH,\n",
" usecols=[i for i in range(0,10)])\n",
"data = data.iloc[:85]\n",
"# Стандартизация только числовых колонок (начиная со второй)\n",
"numeric_cols = data.columns[1:] # Все колонки кроме первой\n",
"data[numeric_cols] = (data[numeric_cols] - data[numeric_cols].mean(axis=0)) / data[numeric_cols].std(axis=0)\n",
"\n",
"print('data:')\n",
"print(data.head())\n",
"\n",
"\n",
"FEATURES = [f'X{i}' for i in range(1,10)]\n",
"# Сохраняем первую колонку (названия субъектов) и признаки X1-X9 для вывода в Excel\n",
"data_to_excel = data[[data.columns[0]] + FEATURES].copy()\n"
]
},
{
"cell_type": "code",
"execution_count": 40,
"id": "078c2afe",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" X1 X2 X3 X4 X5 X6 X7 \\\n",
"4 -0.116102 -0.030388 0.922129 0.778060 -0.533539 -0.390494 1.025235 \n",
"23 0.191144 -0.800989 0.470034 0.663782 -0.446576 -0.809216 0.423536 \n",
"28 -0.191100 -0.367526 1.018027 -0.064734 -0.955121 -0.997007 0.629696 \n",
"53 -0.808551 -0.560176 1.524921 0.468559 0.398480 -0.924936 0.001627 \n",
"59 -0.133486 -0.271201 0.771431 0.049543 -1.113924 -0.471701 -0.075084 \n",
"69 -0.450778 -0.271201 1.100226 0.240005 -0.400258 -0.893468 -0.012757 \n",
"13 -0.372445 -0.271201 -1.338342 -0.731351 -0.250908 -0.329589 0.080735 \n",
"21 -0.490466 -0.223038 -0.872548 -1.012282 0.369177 0.132275 -0.451446 \n",
"34 -0.464369 -0.656501 -1.379442 -0.331381 -0.420108 -0.242799 -0.247683 \n",
"36 -0.385119 -0.030388 -1.434241 -0.098065 0.190525 -0.204226 -0.489801 \n",
"39 -0.567507 -0.367526 -0.269756 -0.088542 -0.283047 -0.404197 -0.060701 \n",
"67 -0.304660 -0.319364 -0.584852 -0.807536 0.354053 -0.890931 0.078337 \n",
"70 -0.365525 -0.223038 -0.351955 -0.745636 -0.107230 0.055636 -0.111042 \n",
"83 -0.139989 0.065937 -0.762950 -0.669451 0.486389 0.402288 -0.041523 \n",
"14 -0.315123 0.017774 -0.242356 -0.431374 -0.080763 0.619009 -0.650414 \n",
"20 -0.443316 -0.126713 -0.310856 -0.550412 -0.622393 -0.317915 -0.290833 \n",
"24 0.014092 -0.800989 -0.091658 0.268574 -0.038226 0.004882 -0.700756 \n",
"40 -0.183554 -0.415689 -0.584852 0.397136 -0.064694 0.397720 -0.377132 \n",
"45 -0.417719 -0.656501 -0.119058 0.563790 -0.562842 0.910338 0.049571 \n",
"61 -0.306035 -0.800989 0.031640 0.782821 0.034558 0.468269 0.114296 \n",
"65 0.140534 -0.608339 -0.612252 -0.388520 -0.468316 0.399243 -0.417885 \n",
"68 -0.713376 0.065937 0.196038 -0.002834 -0.522196 0.016048 -0.055906 \n",
"74 -0.462868 -0.560176 -0.105358 0.025735 0.344601 0.278955 -0.082276 \n",
"78 -0.431018 -0.608339 0.168638 0.597121 -0.678162 0.488063 -0.132617 \n",
"2 0.246382 0.884700 -0.982147 -0.707543 1.135777 1.441228 -0.657606 \n",
"54 0.791628 0.981025 0.648132 -0.364712 0.940110 1.300639 -0.719694 \n",
"18 1.542984 -0.752826 0.374136 -1.755084 1.200054 0.153592 -0.718975 \n",
"31 1.073236 -0.126713 0.538534 -1.255121 0.591312 0.403811 -0.671989 \n",
"10 2.795688 -0.897314 2.127713 2.358894 2.356098 -0.758970 5.316241 \n",
"62 0.916945 -0.993639 0.291937 -0.240912 3.882681 1.425494 5.364185 \n",
"47 -0.431727 3.967103 -0.201257 1.549430 -0.964574 -0.246352 0.413947 \n",
"79 -0.588726 3.389152 0.127539 2.635063 -1.607346 -2.105480 0.833458 \n",
"\n",
" X8 X9 Class \n",
"4 -0.163306 -0.827550 1 \n",
"23 -0.201096 0.154444 1 \n",
"28 0.554701 0.194090 1 \n",
"53 -0.465625 -0.800103 1 \n",
"59 -0.314465 -0.998331 1 \n",
"69 -0.503415 0.401467 1 \n",
"13 0.781440 0.471610 2 \n",
"21 -0.314465 1.185232 2 \n",
"34 -0.238886 0.496007 2 \n",
"36 0.932599 1.532895 2 \n",
"39 0.592491 0.755229 2 \n",
"67 0.327962 0.230686 2 \n",
"70 -0.427835 0.599696 2 \n",
"83 0.365752 0.590547 2 \n",
"14 0.403541 -0.059033 3 \n",
"20 -0.352255 -0.080381 3 \n",
"24 0.025643 0.105649 3 \n",
"40 -0.465625 -0.559179 3 \n",
"45 -0.012147 -0.357901 3 \n",
"61 -0.654574 -0.116977 3 \n",
"65 -0.730154 -1.041027 3 \n",
"68 -1.070262 -0.424993 3 \n",
"74 -1.448161 -0.306056 3 \n",
"78 -0.352255 0.242884 3 \n",
"2 -0.427835 -0.629321 4 \n",
"54 -0.994683 -0.607974 4 \n",
"18 2.066294 0.904663 5 \n",
"31 0.781440 -0.138324 5 \n",
"10 -0.767944 -1.181311 6 \n",
"62 -0.049937 -0.873295 6 \n",
"47 -0.049937 -1.477129 7 \n",
"79 -0.654574 -2.733594 7 \n"
]
}
],
"source": [
"def get_train_data(data, features, train_samples=None): #объявление функции принимает 3 параметра последние не обязательно вводить (3 это обучающая выборка) фичерс это Х загаловки а дата это таблица с данными\n",
" train_data = pd.DataFrame() # создаётся переменная трейн дата( новый дата фрейм дата фрейм это таблица )\n",
" for cls, samples in train_samples.items(): # type: ignore # заходим в цикл переменная cls будет последовательно принемать значения номеров классов то есть каждый цикл это следующий номер класса, sempls бедет принимать список номеров строк которые принадлежат этому классу точка айтемс это пара значений (ключ и значение) для словаря train samples ключ это номер класса из обучающей выборки а значение это номера строк этого класса из обучающей выборки\n",
" train_samps = data[features].loc[samples] #train_samps это новая таблица в которую мы записываем значения(Х) которое входят в обучающаю выборку (значение номер класса и так до конца)\n",
" train_samps[\"Class\"] = cls # это доп колонка номер класса в которую вносится номер класса к которому принадлежит внесённый X\n",
" train_data = pd.concat([train_data, train_samps]) # сохраняем полоученные строчки в train_data а train_samps в котором они были до этого очищаем\n",
" train_data = train_data.astype({\"Class\": 'int32'}) # в train_data присваиваем итоговую таблицу с обучающей выборкой и показываем что в колонке класс данные типа интедгер\n",
" return train_data\n",
"\n",
"\n",
"train_data = get_train_data(data, FEATURES, TRAIN_SAMPLE) # type: ignore\n",
"print(train_data)\n",
"\n",
"# Создаём колонку с классами для всех строк (NaN для строк не из обучающей выборки)\n",
"train_sample_column = pd.Series(index=data.index, dtype='Int32')\n",
"train_sample_column[train_data.index] = train_data.Class\n",
"data_to_excel['Train sample'] = train_sample_column\n"
]
},
{
"cell_type": "code",
"execution_count": 41,
"id": "c05cdd57",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"<class 'pandas.core.frame.DataFrame'>\n",
"Index: 32 entries, 4 to 79\n",
"Data columns (total 10 columns):\n",
" # Column Non-Null Count Dtype \n",
"--- ------ -------------- ----- \n",
" 0 X1 32 non-null float64\n",
" 1 X2 32 non-null float64\n",
" 2 X3 32 non-null float64\n",
" 3 X4 32 non-null float64\n",
" 4 X5 32 non-null float64\n",
" 5 X6 32 non-null float64\n",
" 6 X7 32 non-null float64\n",
" 7 X8 32 non-null float64\n",
" 8 X9 32 non-null float64\n",
" 9 Class 32 non-null int32 \n",
"dtypes: float64(9), int32(1)\n",
"memory usage: 2.6 KB\n"
]
}
],
"source": [
"train_data.info()"
]
},
{
"cell_type": "markdown",
"id": "213ef4f8",
"metadata": {},
"source": [
"## 2 Часть. Преддискриминантный анализ"
]
},
{
"cell_type": "code",
"execution_count": 42,
"id": "95268435",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Ковариационная матрица\n",
" X1 X2 X3 X4 X5 X6 X7 \\\n",
"X1 0.131594 -0.011451 0.043813 0.088946 -0.062342 -0.066520 -0.002256 \n",
"X2 -0.011451 0.080149 0.010255 -0.041460 -0.008208 0.027059 -0.004493 \n",
"X3 0.043813 0.010255 0.236552 0.125222 -0.042629 -0.105374 0.015503 \n",
"X4 0.088946 -0.041460 0.125222 0.302019 -0.090391 -0.126444 0.039835 \n",
"X5 -0.062342 -0.008208 -0.042629 -0.090391 0.196549 0.090677 -0.022302 \n",
"X6 -0.066520 0.027059 -0.105374 -0.126444 0.090677 0.266100 -0.007099 \n",
"X7 -0.002256 -0.004493 0.015503 0.039835 -0.022302 -0.007099 0.084685 \n",
"X8 0.007213 0.003786 -0.085371 -0.054780 0.012072 0.045550 -0.008861 \n",
"X9 -0.004419 -0.001432 -0.024718 -0.030794 0.043951 0.043405 -0.034830 \n",
"\n",
" X8 X9 \n",
"X1 0.007213 -0.004419 \n",
"X2 0.003786 -0.001432 \n",
"X3 -0.085371 -0.024718 \n",
"X4 -0.054780 -0.030794 \n",
"X5 0.012072 0.043951 \n",
"X6 0.045550 0.043405 \n",
"X7 -0.008861 -0.034830 \n",
"X8 0.269534 0.101226 \n",
"X9 0.101226 0.232664 \n"
]
}
],
"source": [
"def scatter_matrix(samples):\n",
" # является ли подклассом?\n",
" if isinstance(samples, pd.Series): # проверка если по каким то причинам наши значения признаков имеют тип данных series то их конвертирует в data frame\n",
" samples = samples.to_frame()\n",
" d = samples - samples.mean() # вычитает из значений признаков средние значения признаков\n",
" res = np.zeros((d.shape[1], d.shape[1])) #создаёт матрицу нулей размерностью 9 на 9\n",
" # приводит к виду int: 32, 24, ...\n",
" for _, row in d.iterrows(): # проходимся циклом по каждой строке матрицы D датафрейм(там где от значений - срзнач)\n",
" col = row.to_frame() # берём строчку из таблицы D( которая сейчас имеет вид seria и мы приводим её к виду dataframe)\n",
" res += col @ col.T # матрица из нулей 9х9 берём строчку из датафрейма и умножаем её на неё же только транспонированную и так проходим по всем строчкам\n",
" return res\n",
"\n",
"\n",
"def classes_scatter_matrix(samples, labels): # передаём samples (значение признаков в обучающей выборке табличкой) и labels (номера классов) shape если 0 то число строк если 1 то число колонок\n",
" A = np.zeros((samples.shape[1], samples.shape[1])) # zeros создаёт матрицу shape на shape то есть создаёт матрицу число признаков на число призноков (9 на 9) заполненую нулями\n",
" for cls in labels.unique(): # переменная cls счётчик по классам принимает значения уникальных классов то есть у нас от 1 до 7\n",
" A += scatter_matrix(samples[labels == cls]) # В матрицу А прибавляем соответствующие значения из матрицы полученые в результате работы skater matrix для текущего класса\n",
" return A\n",
"\n",
"\n",
"\n",
"cov = pd.DataFrame(classes_scatter_matrix(train_data[FEATURES], train_data.Class) / (train_data.shape[0] - train_data.Class.unique().size), \\\n",
" index=FEATURES, columns=FEATURES)\n",
"\n",
"print('Ковариационная матрица')\n",
"print(cov)"
]
},
{
"cell_type": "code",
"execution_count": 43,
"id": "b2be0421",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Средние значения\n",
" X1 X2 X3 X4 X5 X6 X7 \\\n",
"1 -0.251479 -0.383580 0.967795 0.355869 -0.508490 -0.747804 0.332042 \n",
"2 -0.386260 -0.253140 -0.874261 -0.560531 0.042356 -0.185193 -0.155391 \n",
"3 -0.311838 -0.449402 -0.167007 0.126204 -0.265843 0.326461 -0.254395 \n",
"4 0.519005 0.932863 -0.167007 -0.536128 1.037944 1.370933 -0.688650 \n",
"5 1.308110 -0.439770 0.456335 -1.505103 0.895683 0.278701 -0.695482 \n",
"6 1.856316 -0.945477 1.209825 1.058991 3.119390 0.333262 5.340213 \n",
"7 -0.510227 3.678128 -0.036859 2.092247 -1.285960 -1.175916 0.623703 \n",
"\n",
" X8 X9 \n",
"1 -0.182201 -0.312664 \n",
"2 0.252382 0.732738 \n",
"3 -0.465625 -0.259701 \n",
"4 -0.711259 -0.618647 \n",
"5 1.423867 0.383169 \n",
"6 -0.408940 -1.027303 \n",
"7 -0.352255 -2.105362 \n",
" X1 X2 X3 X4 X5 X6 X7 \\\n",
"1 -0.251479 -0.383580 0.967795 0.355869 -0.508490 -0.747804 0.332042 \n",
"2 -0.386260 -0.253140 -0.874261 -0.560531 0.042356 -0.185193 -0.155391 \n",
"3 -0.311838 -0.449402 -0.167007 0.126204 -0.265843 0.326461 -0.254395 \n",
"4 0.519005 0.932863 -0.167007 -0.536128 1.037944 1.370933 -0.688650 \n",
"5 1.308110 -0.439770 0.456335 -1.505103 0.895683 0.278701 -0.695482 \n",
"6 1.856316 -0.945477 1.209825 1.058991 3.119390 0.333262 5.340213 \n",
"7 -0.510227 3.678128 -0.036859 2.092247 -1.285960 -1.175916 0.623703 \n",
"\n",
" X8 X9 \n",
"1 -0.182201 -0.312664 \n",
"2 0.252382 0.732738 \n",
"3 -0.465625 -0.259701 \n",
"4 -0.711259 -0.618647 \n",
"5 1.423867 0.383169 \n",
"6 -0.408940 -1.027303 \n",
"7 -0.352255 -2.105362 \n"
]
}
],
"source": [
"lda = LinearDiscriminantAnalysis()\n",
"lda.fit(train_data[FEATURES], train_data.Class)\n",
"means = pd.DataFrame(lda.means_, index=lda.classes_, columns=FEATURES) # type: ignore\n",
"print('Средние значения')\n",
"print(means)"
]
},
{
"cell_type": "code",
"execution_count": 44,
"id": "06d105f8",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Расстояние Махаланобиса (обучающая выборка)\n",
" 1 2 3 4 5 6 \\\n",
"1 0.0 22.343817 15.22444 87.857078 76.915784 580.080057 \n",
"2 22.343817 0.0 15.406909 72.484578 84.871051 620.59257 \n",
"3 15.22444 15.406909 0.0 60.045571 82.766139 660.598405 \n",
"4 87.857078 72.484578 60.045571 0.0 82.932649 618.932157 \n",
"5 76.915784 84.871051 82.766139 82.932649 0.0 562.34511 \n",
"6 580.080057 620.59257 660.598405 618.932157 562.34511 0.0 \n",
"7 356.469575 335.176705 355.369339 267.77654 526.670574 975.229758 \n",
"\n",
" 7 \n",
"1 356.469575 \n",
"2 335.176705 \n",
"3 355.369339 \n",
"4 267.77654 \n",
"5 526.670574 \n",
"6 975.229758 \n",
"7 0.0 \n"
]
}
],
"source": [
"def find_mahl_sqr_dist(centers, samples, covr): # функция принимает в себя дважды средние значения и один раз матрицу ковариций\n",
" res = pd.DataFrame(index=samples.index, columns=centers.index) # создаём новый датафрейм вверху заголовки это номера классов и слева заголовки тоже номера классов\n",
" for i in centers.index: # двойной цикл идёт по одной и той же таблице\n",
" for j in samples.index:\n",
" res[i][j] = mahalanobis(centers.loc[i], samples.loc[j], np.linalg.inv(covr)) ** 2 # вычисляется растояние махаланобиса в квадрате. np.linalg.inv(covr)-возвращает матрицу обратную матрице ковариации centers.loc[i] и samples.loc[j] возвращают i и j строки таблицы means(ср знач) и это значение записывается в ячейку ij\n",
" return res\n",
"\n",
"\n",
"cen_dis = find_mahl_sqr_dist(means, means, cov)\n",
"print('Расстояние Махаланобиса (обучающая выборка)')\n",
"\n",
"print(cen_dis)"
]
},
{
"cell_type": "markdown",
"id": "6c365f10",
"metadata": {},
"source": [
"## 3 Часть. Дискриминантный анализ"
]
},
{
"cell_type": "code",
"execution_count": 45,
"id": "a9ef0be3",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Коэффициенты дискриминантных функций:\n"
]
},
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" <th>3</th>\n",
" <th>4</th>\n",
" <th>5</th>\n",
" <th>6</th>\n",
" <th>7</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>X1</th>\n",
" <td>-5.234825</td>\n",
" <td>-3.148253</td>\n",
" <td>-4.217967</td>\n",
" <td>10.327559</td>\n",
" <td>16.656458</td>\n",
" <td>35.610185</td>\n",
" <td>-11.002064</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X2</th>\n",
" <td>-7.461316</td>\n",
" <td>-2.651336</td>\n",
" <td>-7.584118</td>\n",
" <td>13.105596</td>\n",
" <td>-11.178379</td>\n",
" <td>-4.201862</td>\n",
" <td>61.821871</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X3</th>\n",
" <td>5.873924</td>\n",
" <td>-4.182788</td>\n",
" <td>0.356458</td>\n",
" <td>-1.673889</td>\n",
" <td>9.525550</td>\n",
" <td>3.394379</td>\n",
" <td>-16.353702</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X4</th>\n",
" <td>-2.309080</td>\n",
" <td>-1.313906</td>\n",
" <td>0.911015</td>\n",
" <td>2.145491</td>\n",
" <td>-11.656789</td>\n",
" <td>-10.600878</td>\n",
" <td>22.281651</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X5</th>\n",
" <td>-3.567483</td>\n",
" <td>-1.664356</td>\n",
" <td>-4.297066</td>\n",
" <td>7.670232</td>\n",
" <td>4.942687</td>\n",
" <td>32.131610</td>\n",
" <td>5.652878</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X6</th>\n",
" <td>-1.022878</td>\n",
" <td>-3.339634</td>\n",
" <td>3.245656</td>\n",
" <td>5.188590</td>\n",
" <td>2.273415</td>\n",
" <td>-1.875044</td>\n",
" <td>-10.339378</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X7</th>\n",
" <td>2.264382</td>\n",
" <td>0.112277</td>\n",
" <td>-5.471875</td>\n",
" <td>-7.586796</td>\n",
" <td>-3.817930</td>\n",
" <td>77.090589</td>\n",
" <td>0.590885</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X8</th>\n",
" <td>1.619403</td>\n",
" <td>-1.400612</td>\n",
" <td>-1.324989</td>\n",
" <td>-2.957457</td>\n",
" <td>6.123916</td>\n",
" <td>-2.939071</td>\n",
" <td>1.797481</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X9</th>\n",
" <td>-0.671626</td>\n",
" <td>4.018564</td>\n",
" <td>-1.121021</td>\n",
" <td>-4.541955</td>\n",
" <td>-3.230166</td>\n",
" <td>2.291908</td>\n",
" <td>-7.498102</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Const</th>\n",
" <td>-7.817591</td>\n",
" <td>-6.087361</td>\n",
" <td>-5.803707</td>\n",
" <td>-23.736379</td>\n",
" <td>-34.669499</td>\n",
" <td>-289.317890</td>\n",
" <td>-153.089698</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" 1 2 3 4 5 6 \\\n",
"X1 -5.234825 -3.148253 -4.217967 10.327559 16.656458 35.610185 \n",
"X2 -7.461316 -2.651336 -7.584118 13.105596 -11.178379 -4.201862 \n",
"X3 5.873924 -4.182788 0.356458 -1.673889 9.525550 3.394379 \n",
"X4 -2.309080 -1.313906 0.911015 2.145491 -11.656789 -10.600878 \n",
"X5 -3.567483 -1.664356 -4.297066 7.670232 4.942687 32.131610 \n",
"X6 -1.022878 -3.339634 3.245656 5.188590 2.273415 -1.875044 \n",
"X7 2.264382 0.112277 -5.471875 -7.586796 -3.817930 77.090589 \n",
"X8 1.619403 -1.400612 -1.324989 -2.957457 6.123916 -2.939071 \n",
"X9 -0.671626 4.018564 -1.121021 -4.541955 -3.230166 2.291908 \n",
"Const -7.817591 -6.087361 -5.803707 -23.736379 -34.669499 -289.317890 \n",
"\n",
" 7 \n",
"X1 -11.002064 \n",
"X2 61.821871 \n",
"X3 -16.353702 \n",
"X4 22.281651 \n",
"X5 5.652878 \n",
"X6 -10.339378 \n",
"X7 0.590885 \n",
"X8 1.797481 \n",
"X9 -7.498102 \n",
"Const -153.089698 "
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"classes = np.unique(train_data.Class) # Получаем массив уникальных значений классов из обучающей выборки\n",
"\n",
"# Создаём список групп - для каждого класса выбираем все строки с признаками FEATURES, которые принадлежат этому классу\n",
"# Результат: список из DataFrame'ов, каждый содержит объекты одного класса\n",
"groups = [train_data[FEATURES][train_data.Class == cls] for cls in classes]\n",
"\n",
"\n",
"n = [len(g) for g in groups] # Создаём список n - количество объектов в каждом классе\n",
"\n",
"N = sum(n) # N - общее количество объектов в обучающей выборке (сумма всех элементов списка n)\n",
"\n",
"\n",
"p = train_data[FEATURES].shape[1] # количество признаков\n",
"\n",
"# Вычисляем объединённую (pooled) ковариационную матрицу:\n",
"# rowvar=False означает, что переменные в столбцах, ddof=1 - поправка на смещение\n",
"S_pooled = sum((ni - 1) * np.cov(g, rowvar=False, ddof=1) for g, ni in zip(groups, n)) / (N - len(classes))\n",
"\n",
"# Вычисляем обратную матрицу к объединённой ковариационной матрице\n",
"inv_S = np.linalg.inv(S_pooled)\n",
"\n",
"# Создаём список средних векторов для каждого класса\n",
"means = [g.mean(axis=0) for g in groups]\n",
"\n",
"# Вычисляем априорные вероятности классов (prior probabilities)\n",
"priors = np.array(n) / N\n",
"\n",
"# Создаём пустые словари для хранения коэффициентов дискриминантных функций\n",
"coef_stat = {} # Словарь для хранения векторов коэффициентов a для каждого класса\n",
"const_stat = {} # Словарь для хранения константных членов c для каждого класса\n",
"\n",
"# Проходим по всем классам и вычисляем коэффициенты дискриминантной функции для каждого\n",
"for cls, mu, p_j in zip(classes, means, priors):\n",
" # a - вектор коэффициентов при признаках X1...X9\n",
" a = inv_S @ mu # Матричное произведение\n",
"\n",
" # c - константный член дискриминантной функции\n",
" # p_j - априорная вероятность класса\n",
" c = -0.5 * mu.T @ inv_S @ mu + np.log(p_j)\n",
"\n",
" # Сохраняем вычисленные коэффициенты в словари по номеру класса\n",
" coef_stat[cls] = a\n",
" const_stat[cls] = c\n",
"\n",
"# Создаём DataFrame для удобного отображения коэффициентов\n",
"# Строки: X1, X2, ..., X9; Столбцы: номера классов\n",
"df_stat = pd.DataFrame(coef_stat, index=[f\"X{i+1}\" for i in range(p)])\n",
"\n",
"# Добавляем строку с константными членами в конец таблицы\n",
"df_stat.loc[\"Const\"] = const_stat\n",
"\n",
"print(\"Коэффициенты дискриминантных функций:\")\n",
"display(df_stat)\n"
]
},
{
"cell_type": "code",
"execution_count": 46,
"id": "7de23d86",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Распределение по классам\n",
" Class\n",
"0 3\n",
"1 2\n",
"2 4\n",
"3 2\n",
"4 1\n",
".. ...\n",
"79 7\n",
"80 1\n",
"81 5\n",
"82 5\n",
"83 2\n",
"\n",
"[84 rows x 1 columns]\n"
]
}
],
"source": [
"lda = LinearDiscriminantAnalysis()\n",
"lda.fit(train_data[FEATURES], train_data.Class)\n",
"means = pd.DataFrame(lda.means_, index=lda.classes_, columns=FEATURES) #type: ignore\n",
"\n",
"def LDA_predict(lda, x):# принимает в себя результаты линейного дискр анализа и значения признаков всех объектов( исходная таблица только с значениями X1...X9)\n",
" return pd.DataFrame( # функция считает распределение по классам (классификация масива тестовых векторов Х )\n",
" lda.predict(x),\n",
" columns=[\"Class\"],\n",
" index=x.index\n",
" )\n",
"\n",
"\n",
"lda_predict = LDA_predict(lda, data[FEATURES])\n",
"print('Распределение по классам')\n",
"print(lda_predict)\n",
"data_to_excel['Result Lda'] = lda_predict"
]
},
{
"cell_type": "code",
"execution_count": 47,
"id": "f37c3126",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Расстояние Махланобиса\n",
" 1 2 3 4 5 6 \\\n",
"0 35.548618 16.784286 11.032019 67.652466 106.380096 685.635422 \n",
"1 44.397398 15.17559 28.556746 72.867432 63.671602 621.00928 \n",
"2 92.049277 64.122689 58.571719 3.93205 97.092489 629.959981 \n",
"3 45.929629 37.680543 49.035861 60.870665 92.950853 658.258753 \n",
"4 10.384607 36.685522 29.484569 81.018346 95.923694 485.761693 \n",
".. ... ... ... ... ... ... \n",
"79 361.914041 347.840907 367.588069 300.470231 556.010044 989.726755 \n",
"80 11.323179 32.848072 20.11476 85.940301 75.565423 551.195134 \n",
"81 529.682935 494.161709 479.158855 331.787176 276.596841 630.305397 \n",
"82 369.683371 396.241564 330.687362 230.499407 183.163732 621.077741 \n",
"83 26.559382 6.11464 19.063954 43.463375 66.032665 555.489917 \n",
"\n",
" 7 \n",
"0 362.911582 \n",
"1 439.360623 \n",
"2 247.054034 \n",
"3 248.019142 \n",
"4 299.084459 \n",
".. ... \n",
"79 6.031622 \n",
"80 372.395174 \n",
"81 977.516974 \n",
"82 891.232175 \n",
"83 309.482441 \n",
"\n",
"[84 rows x 7 columns]\n"
]
}
],
"source": [
"samp_dist = find_mahl_sqr_dist(means, data[FEATURES], cov)\n",
"print('Расстояние Махланобиса')\n",
"print(samp_dist)"
]
},
{
"cell_type": "code",
"execution_count": 48,
"id": "8ba2b145",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Вероятности\n",
" 1 2 3 4 5 \\\n",
"0 2.724508e-06 4.313697e-02 9.568603e-01 9.702972e-14 3.778389e-22 \n",
"1 3.380276e-07 9.984489e-01 1.550797e-03 7.406988e-14 7.353441e-12 \n",
"2 2.201448e-19 3.402740e-13 6.825218e-12 1.000000e+00 5.894780e-21 \n",
"3 1.193238e-02 9.838574e-01 4.207907e-03 2.265694e-06 2.449497e-13 \n",
"4 9.998788e-01 2.592478e-06 1.186566e-04 1.530751e-16 8.876653e-20 \n",
".. ... ... ... ... ... \n",
"79 1.578468e-77 2.393959e-74 1.541647e-78 1.157359e-64 3.747166e-120 \n",
"80 9.798381e-01 2.767082e-05 2.013425e-02 2.047062e-17 3.664435e-15 \n",
"81 3.312640e-55 2.282701e-47 5.166393e-44 1.036504e-12 1.000000e+00 \n",
"82 9.438354e-41 2.151774e-46 4.619775e-32 5.262459e-11 1.000000e+00 \n",
"83 2.720782e-05 9.980490e-01 1.923774e-03 1.936053e-09 2.432529e-14 \n",
"\n",
" 6 7 \n",
"0 6.217431e-148 7.450791e-78 \n",
"1 6.952668e-133 1.934890e-93 \n",
"2 1.147524e-136 1.609652e-53 \n",
"3 4.305525e-136 5.205061e-47 \n",
"4 1.976985e-104 6.799073e-64 \n",
".. ... ... \n",
"79 2.473522e-214 1.000000e+00 \n",
"80 1.915701e-118 1.283056e-79 \n",
"81 1.560137e-77 6.267969e-153 \n",
"82 8.094335e-96 1.757483e-154 \n",
"83 1.264055e-120 3.323616e-67 \n",
"\n",
"[84 rows x 7 columns]\n"
]
}
],
"source": [
"def LDA_predict_probab(lda, x): # принимает в себя результаты линейного дискр анализа и значения признаков всех объектов( исходная таблица только с значениями X1...X9)\n",
" return pd.DataFrame( # возвращает апостериорные вероятности классификации в соответствии с каждым классом в массиве тестовых векторов\n",
" lda.predict_proba(x),\n",
" columns=lda.classes_,\n",
" index=x.index\n",
" )\n",
"\n",
"lda_post_prob = LDA_predict_probab(lda, data[FEATURES])\n",
"print('Вероятности')\n",
"print(lda_post_prob)"
]
},
{
"cell_type": "markdown",
"id": "df1e683e",
"metadata": {},
"source": [
"## 4 Часть. Пошаговый ДА с включением"
]
},
{
"cell_type": "code",
"execution_count": 49,
"id": "243acaa1",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Forward stepwise\n",
"Step: 0\n",
"Empty DataFrame\n",
"Columns: [Wilk's lmbd, Partial lmbd, F to enter, P value]\n",
"Index: []\n",
"\n",
"Step: 1\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 1.0 0.034339 117.172492 4.570418e-17\n",
"\n",
"Step: 2\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.055327 0.034633 111.495293 2.539633e-16\n",
"X2 0.034339 0.055802 67.682511 7.462368e-14\n",
"\n",
"Step: 3\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.009143 0.053446 67.889761 1.792580e-13\n",
"X2 0.007714 0.063345 56.682141 1.240902e-12\n",
"X5 0.001916 0.255016 11.198399 7.409067e-06\n",
"\n",
"Step: 4\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.002715 0.055551 62.338586 1.092269e-12\n",
"X2 0.002441 0.061807 55.658064 3.490100e-12\n",
"X5 0.000531 0.283947 9.246521 4.132537e-05\n",
"X1 0.000489 0.308693 8.211369 9.737300e-05\n",
"\n",
"Step: 5\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.000903 0.066279 49.307354 2.671740e-11\n",
"X2 0.000960 0.062367 52.619515 1.421293e-11\n",
"X5 0.000181 0.331033 7.072970 3.186012e-04\n",
"X1 0.000168 0.356066 6.329628 6.410529e-04\n",
"X3 0.000151 0.396708 5.322606 1.782148e-03\n",
"\n",
"Step: 6\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.000346 0.079099 38.807979 5.450728e-10\n",
"X2 0.000470 0.058251 53.890403 2.662696e-11\n",
"X5 0.000064 0.427580 4.462485 5.054419e-03\n",
"X1 0.000078 0.349448 6.205518 8.358543e-04\n",
"X3 0.000077 0.354847 6.060375 9.604737e-04\n",
"X6 0.000060 0.457851 3.947056 9.143395e-03\n",
"\n",
"Step: 7\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.000094 0.135646 20.178341 2.661233e-07\n",
"X2 0.000196 0.064997 45.553320 2.823152e-10\n",
"X5 0.000023 0.542947 2.665705 4.761796e-02\n",
"X1 0.000039 0.322062 6.665822 6.423386e-04\n",
"X3 0.000038 0.333369 6.332330 8.660969e-04\n",
"X6 0.000028 0.453227 3.820269 1.146526e-02\n",
"X4 0.000027 0.463805 3.660921 1.382036e-02\n",
"\n",
"Step: 8\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.000061 0.135927 19.070677 6.712073e-07\n",
"X2 0.000090 0.091516 29.780979 2.085605e-08\n",
"X5 0.000015 0.548214 2.472314 6.373804e-02\n",
"X1 0.000026 0.321306 6.336900 1.015968e-03\n",
"X3 0.000022 0.379230 4.910770 3.881752e-03\n",
"X6 0.000016 0.501614 2.980693 3.350271e-02\n",
"X4 0.000016 0.501861 2.977754 3.362441e-02\n",
"X9 0.000013 0.648374 1.626957 1.970044e-01\n",
"\n",
"Step: 9\n",
" Wilk's lmbd Partial lmbd F to enter P value\n",
"X7 0.000039 0.137183 17.820289 1.783586e-06\n",
"X2 0.000053 0.100352 25.400542 1.346467e-07\n",
"X5 0.000010 0.560292 2.223553 9.119062e-02\n",
"X1 0.000015 0.369418 4.836395 4.725286e-03\n",
"X3 0.000015 0.349955 5.262941 3.138798e-03\n",
"X6 0.000010 0.529227 2.520391 6.244848e-02\n",
"X4 0.000009 0.565830 2.174061 9.722558e-02\n",
"X9 0.000008 0.643040 1.572819 2.151830e-01\n",
"X8 0.000008 0.650606 1.521583 2.304240e-01\n",
"\n",
"['X7', 'X2', 'X5', 'X1', 'X3', 'X6']\n"
]
}
],
"source": [
"def wilks_lambda(samples, labels): # samples - на каждой итерации мы передаём в функцию wilks lmbd dataframe сначала с одной колонкой Х1 постепенно увеличивая число колонок пока не дойдём до конца labels - колонка с номерами классов\n",
" if isinstance(samples, pd.Series):\n",
" samples = samples.to_frame()\n",
" # определитель матрицы рассеивания\n",
" dT = np.linalg.det(scatter_matrix(samples))\n",
" # определитель классовой матрицы рассеивания\n",
" dE = np.linalg.det(classes_scatter_matrix(samples, labels))\n",
" return dE / dT\n",
"\n",
"\n",
"def f_p_value(lmbd, n_obj, n_sign, n_cls): # sign это число признаков вошедших в модель lmbd - значение лямбды( число) n_obj - число объектов в обучающей выборке n_cls - число классов в обучающей выборке\n",
" num = (1-lmbd)*(n_obj - n_cls - n_sign)\n",
" den = lmbd * (n_cls - 1)\n",
" f_value = num / den\n",
" p = f.sf(f_value, n_cls-1, n_obj-n_cls-n_sign)\n",
" return f_value, p\n",
"\n",
"def forward(samples, labels): # samples - значение признаков в обучаюзей выборке labels - колонка с номерами классов f_in точность( это F to inter в статистике)\n",
" st_columns = [\"Wilk's lmbd\", \"Partial lmbd\", \"F to enter\", \"P value\"] # создаётся список названий колонок таблицы\n",
" n_cls = labels.unique().size #число уникальных классов нашей обучающей выборки\n",
" n_obj = samples.shape[0] # число объектов в обучающей выборке\n",
" # хранение пременных вне и в модели(е)\n",
" out = {0: pd.DataFrame(columns=st_columns, index=samples.columns, dtype=float)} # создаётся словарик ключу(ключ показывает какой у нас шаг метода) 0 ставится в соответствие дата фрейм(пустой) с колонками из переменной st_colums и индексами(строки таблицы) Х1,,,Х9\n",
" into = {0: pd.DataFrame(columns=st_columns, dtype=float)} # создаётся словарик ключу(ключ показывает какой у нас шаг метода) 0 ставится в соответствие дата фрейм(пустой) с колонками из переменной st_colums\n",
" step = 0 # шаг нашего метода\n",
"\n",
" while True:\n",
" model_lmbd = wilks_lambda(samples[into[step].index], labels)\n",
" # расчёт характеристик элементов вне модели\n",
" for el in out[step].index: # el переменная счётчик по списку индексов датафрейма внутри out (Х1 Х2,,, Х9)\n",
" lmbda = wilks_lambda(samples[into[step].index.tolist() + [el]], labels) # мы с датафрейма samples берём значение из колонок в квадратных скобках список колонок который есть в inta для текущего шага + текущая колонку из счётчика el и эти значения помещаются в функцию вилкс лямбда\n",
" partial_lmbd = lmbda / model_lmbd #\n",
" f_lmbd, p_value = f_p_value(partial_lmbd, n_obj, into[step].index.size, n_cls)\n",
" out[step].loc[el] = lmbda, partial_lmbd, f_lmbd, p_value # type: ignore\n",
" # расчёт характеристик элементов в моделе\n",
" for el in into[step].index:\n",
" lmbda = wilks_lambda(samples[into[step].index.drop(el)], labels)\n",
" partial_lmbd = model_lmbd / lmbda\n",
" f_lmbd, p_value = f_p_value(partial_lmbd, n_obj, into[step].index.size-1, n_cls)\n",
" into[step].loc[el] = lmbda, partial_lmbd, f_lmbd, p_value # type: ignore\n",
"\n",
" if out[step].index.size == 0:\n",
" break\n",
"\n",
" # добавление нового элемента\n",
" el_to_enter = out[step][\"F to enter\"].idxmax()\n",
" into[step+1] = pd.concat([into[step], out[step].loc[[el_to_enter]]])\n",
" out[step+1] = out[step].drop(index=el_to_enter)\n",
"\n",
" step += 1\n",
" return into, out\n",
"\n",
"\n",
"\n",
"into, out = forward(train_data[FEATURES], train_data.Class)\n",
"print(\"Forward stepwise\")\n",
"for i, tab in into.items():\n",
" print(\"Step: \", i)\n",
" print(tab, end=\"\\n\\n\")\n",
"\n",
"forw_stepwise = into[CNT_ENTER].index.tolist() # смотрим каждый шаг выбираем последний шаг где p-value не превышает 0,05. смотрим сколько признаков на этом шаге вошло в модель\n",
"print(forw_stepwise)"
]
},
{
"cell_type": "code",
"execution_count": 50,
"id": "6bb39a9c",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Pi: [0.1875 0.25 0.3125 0.0625 0.0625 0.0625 0.0625]\n",
"Распределение\n",
" Class\n",
"0 3\n",
"1 2\n",
"2 4\n",
"3 2\n",
"4 1\n"
]
}
],
"source": [
"forw_stepwise_lda = LinearDiscriminantAnalysis().fit(train_data[forw_stepwise], train_data.Class)\n",
"\n",
"print(\"Pi: \", forw_stepwise_lda.priors_)\n",
"forw_stepwise_pred = LDA_predict(forw_stepwise_lda, data[forw_stepwise])\n",
"print(\"Распределение\")\n",
"print(forw_stepwise_pred.head())\n",
"data_to_excel[\"Result forward\"] = forw_stepwise_pred"
]
},
{
"cell_type": "code",
"execution_count": 51,
"id": "c8e6ad7e",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Функции Фишера ПДАсВ:\n"
]
},
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" <th>3</th>\n",
" <th>4</th>\n",
" <th>5</th>\n",
" <th>6</th>\n",
" <th>7</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>X1</th>\n",
" <td>1.736280</td>\n",
" <td>-1.743108</td>\n",
" <td>-4.758759</td>\n",
" <td>-5.226456</td>\n",
" <td>-6.519691</td>\n",
" <td>72.846915</td>\n",
" <td>10.519265</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X2</th>\n",
" <td>-6.145526</td>\n",
" <td>-2.381514</td>\n",
" <td>-8.110648</td>\n",
" <td>11.978715</td>\n",
" <td>-4.772944</td>\n",
" <td>0.412523</td>\n",
" <td>51.767905</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X3</th>\n",
" <td>-3.172120</td>\n",
" <td>-0.849988</td>\n",
" <td>-4.681666</td>\n",
" <td>6.484264</td>\n",
" <td>6.900413</td>\n",
" <td>34.438056</td>\n",
" <td>0.270260</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X4</th>\n",
" <td>-5.767018</td>\n",
" <td>-3.575317</td>\n",
" <td>-4.300118</td>\n",
" <td>9.991338</td>\n",
" <td>13.505700</td>\n",
" <td>31.253687</td>\n",
" <td>-3.123795</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X5</th>\n",
" <td>4.403649</td>\n",
" <td>-4.202633</td>\n",
" <td>1.213440</td>\n",
" <td>0.308135</td>\n",
" <td>2.874115</td>\n",
" <td>0.482478</td>\n",
" <td>-8.547150</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X6</th>\n",
" <td>-0.755884</td>\n",
" <td>-2.768609</td>\n",
" <td>2.925547</td>\n",
" <td>4.204502</td>\n",
" <td>3.521690</td>\n",
" <td>-0.577544</td>\n",
" <td>-13.660256</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Const</th>\n",
" <td>-7.086063</td>\n",
" <td>-4.589115</td>\n",
" <td>-5.259904</td>\n",
" <td>-18.973678</td>\n",
" <td>-19.159540</td>\n",
" <td>-280.003432</td>\n",
" <td>-110.069841</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" 1 2 3 4 5 6 \\\n",
"X1 1.736280 -1.743108 -4.758759 -5.226456 -6.519691 72.846915 \n",
"X2 -6.145526 -2.381514 -8.110648 11.978715 -4.772944 0.412523 \n",
"X3 -3.172120 -0.849988 -4.681666 6.484264 6.900413 34.438056 \n",
"X4 -5.767018 -3.575317 -4.300118 9.991338 13.505700 31.253687 \n",
"X5 4.403649 -4.202633 1.213440 0.308135 2.874115 0.482478 \n",
"X6 -0.755884 -2.768609 2.925547 4.204502 3.521690 -0.577544 \n",
"Const -7.086063 -4.589115 -5.259904 -18.973678 -19.159540 -280.003432 \n",
"\n",
" 7 \n",
"X1 10.519265 \n",
"X2 51.767905 \n",
"X3 0.270260 \n",
"X4 -3.123795 \n",
"X5 -8.547150 \n",
"X6 -13.660256 \n",
"Const -110.069841 "
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"forw_stepwise_lda = LinearDiscriminantAnalysis().fit(train_data[forw_stepwise], train_data.Class)\n",
"\n",
"classes = np.unique(train_data.Class)\n",
"groups = [train_data[forw_stepwise][train_data.Class == cls] for cls in classes]\n",
"n = [len(g) for g in groups]\n",
"N = sum(n)\n",
"p = train_data[forw_stepwise].shape[1]\n",
"\n",
"# Общая (pooled) ковариация как в Statistica\n",
"S_pooled = sum((ni - 1) * np.cov(g, rowvar=False, ddof=1) for g, ni in zip(groups, n)) / (N - len(classes))\n",
"inv_S = np.linalg.inv(S_pooled)\n",
"\n",
"# Средние по классам и априорные вероятности\n",
"means = [g.mean(axis=0) for g in groups]\n",
"priors = np.array(n) / N\n",
"\n",
"# Классификационные функции (Statistica)\n",
"coef_stat = {}\n",
"const_stat = {}\n",
"\n",
"for cls, mu, p_j in zip(classes, means, priors):\n",
" a = inv_S @ mu\n",
" c = -0.5 * mu.T @ inv_S @ mu + np.log(p_j)\n",
" coef_stat[cls] = a\n",
" const_stat[cls] = c\n",
"\n",
"df_stat = pd.DataFrame(coef_stat, index=[f\"X{i+1}\" for i in range(p)])\n",
"df_stat.loc[\"Const\"] = const_stat\n",
"print(\"Функции Фишера ПДАсВ:\")\n",
"display(df_stat)"
]
},
{
"cell_type": "markdown",
"id": "45a0da72",
"metadata": {},
"source": [
"## 5 Часть. Пошаговый ДА с исключением"
]
},
{
"cell_type": "code",
"execution_count": 52,
"id": "74093d18",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Backward stepwise\n",
"Step: 0\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.000015 0.369418 4.836395 4.725286e-03\n",
"X2 0.000053 0.100352 25.400542 1.346467e-07\n",
"X3 0.000015 0.349955 5.262941 3.138798e-03\n",
"X4 0.000009 0.565830 2.174061 9.722558e-02\n",
"X5 0.000010 0.560292 2.223553 9.119062e-02\n",
"X6 0.000010 0.529227 2.520391 6.244848e-02\n",
"X7 0.000039 0.137183 17.820289 1.783586e-06\n",
"X8 0.000008 0.650606 1.521583 2.304240e-01\n",
"X9 0.000008 0.643040 1.572819 2.151830e-01\n",
"\n",
"Step: 1\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.000026 0.321306 6.336900 1.015968e-03\n",
"X2 0.000090 0.091516 29.780979 2.085605e-08\n",
"X3 0.000022 0.379230 4.910770 3.881752e-03\n",
"X4 0.000016 0.501861 2.977754 3.362441e-02\n",
"X5 0.000015 0.548214 2.472314 6.373804e-02\n",
"X6 0.000016 0.501614 2.980693 3.350271e-02\n",
"X7 0.000061 0.135927 19.070677 6.712073e-07\n",
"X9 0.000013 0.648374 1.626957 1.970044e-01\n",
"\n",
"Step: 2\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.000039 0.322062 6.665822 6.423386e-04\n",
"X2 0.000196 0.064997 45.553320 2.823152e-10\n",
"X3 0.000038 0.333369 6.332330 8.660969e-04\n",
"X4 0.000027 0.463805 3.660921 1.382036e-02\n",
"X5 0.000023 0.542947 2.665705 4.761796e-02\n",
"X6 0.000028 0.453227 3.820269 1.146526e-02\n",
"X7 0.000094 0.135646 20.178341 2.661233e-07\n",
"\n",
"Step: 3\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.000082 0.284888 8.367167 1.277550e-04\n",
"X2 0.000340 0.068884 45.057411 1.394937e-10\n",
"X3 0.000069 0.337744 6.536070 6.131245e-04\n",
"X4 0.000064 0.365254 5.792738 1.247026e-03\n",
"X6 0.000064 0.366742 5.755716 1.293539e-03\n",
"X7 0.000347 0.067374 46.141445 1.121008e-10\n",
"\n",
"Step: 4\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.000232 0.275052 9.224878 5.243354e-05\n",
"X2 0.000908 0.070269 46.308408 4.898174e-11\n",
"X3 0.000187 0.341989 6.734261 4.357539e-04\n",
"X4 0.000181 0.353031 6.414144 5.907362e-04\n",
"X7 0.000947 0.067399 48.429295 3.179135e-11\n",
"\n",
"Step: 5\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.000655 0.275907 9.622842 3.073058e-05\n",
"X2 0.003223 0.056083 61.712393 1.211773e-12\n",
"X3 0.000531 0.340281 7.108733 2.619443e-04\n",
"X7 0.003487 0.051840 67.063885 5.141883e-13\n",
"\n",
"Step: 6\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X1 0.001916 0.277240 9.993420 1.834980e-05\n",
"X2 0.009351 0.056813 63.638909 3.595709e-13\n",
"X7 0.010934 0.048586 75.064819 6.044530e-14\n",
"\n",
"Step: 7\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X2 0.034339 0.055802 67.682511 7.462368e-14\n",
"X7 0.055327 0.034633 111.495293 2.539633e-16\n",
"\n",
"Step: 8\n",
" Wilk's lmbd Partial lmbd F to remove P value\n",
"X7 1.0 0.034339 117.172492 4.570418e-17\n",
"\n",
"Step: 9\n",
"Empty DataFrame\n",
"Columns: [Wilk's lmbd, Partial lmbd, F to remove, P value]\n",
"Index: []\n",
"\n",
"['X1', 'X2', 'X3', 'X4', 'X6', 'X7']\n"
]
}
],
"source": [
"def backward(samples, labels):\n",
" st_columns = [\"Wilk's lmbd\", \"Partial lmbd\", \"F to remove\", \"P value\"]\n",
" n_cls = labels.unique().size\n",
" n_obj = samples.shape[0]\n",
" # хранение пременных вне и в модели(е)\n",
" into = {0: pd.DataFrame(columns=st_columns, index=samples.columns, dtype=float)}\n",
" out = {0: pd.DataFrame(columns=st_columns, dtype=float)}\n",
" step = 0\n",
"\n",
" while True:\n",
" # print(step)\n",
" model_lmbd = wilks_lambda(samples[into[step].index], labels)\n",
" # расчёт характеристик элементов вне модели\n",
" for el in out[step].index:\n",
" lmbda = wilks_lambda(samples[into[step].index.tolist() + [el]], labels)\n",
" partial_lmbd = lmbda / model_lmbd\n",
" f_lmbd, p_value = f_p_value(partial_lmbd, n_obj, into[step].index.size, n_cls)\n",
" out[step].loc[el] = lmbda, partial_lmbd, f_lmbd, p_value # type: ignore\n",
" # расчёт характеристик элементов в моделе\n",
" for el in into[step].index:\n",
" lmbda = wilks_lambda(samples[into[step].index.drop(el)], labels)\n",
" partial_lmbd = model_lmbd / lmbda\n",
" f_lmbd, p_value = f_p_value(partial_lmbd, n_obj, into[step].index.size-1, n_cls)\n",
" into[step].loc[el] = lmbda, partial_lmbd, f_lmbd, p_value # type: ignore\n",
"\n",
" if into[step].index.size == 0:\n",
" break\n",
"\n",
" # удаление элемента\n",
" el_to_remove = into[step][\"F to remove\"].idxmin()\n",
" out[step+1] = pd.concat([out[step], into[step].loc[[el_to_remove]]])\n",
" into[step+1] = into[step].drop(index=el_to_remove)\n",
"\n",
" step += 1\n",
" return into, out\n",
"\n",
"\n",
"into, out = backward(train_data[FEATURES], train_data.Class)\n",
"print(\"Backward stepwise\")\n",
"for i, tab in into.items():\n",
" print(\"Step: \", i)\n",
" print(tab, end=\"\\n\\n\")\n",
"\n",
"back_stepwise = into[len(into) - 1 - CNT_REMOVE].index.tolist()\n",
"print(back_stepwise)"
]
},
{
"cell_type": "code",
"execution_count": 53,
"id": "f804180f",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Pi: [0.1875 0.25 0.3125 0.0625 0.0625 0.0625 0.0625]\n",
"Распределение\n",
" Class\n",
"0 3\n",
"1 3\n",
"2 4\n",
"3 2\n",
"4 1\n"
]
}
],
"source": [
"back_stepwise_lda = LinearDiscriminantAnalysis().fit(train_data[back_stepwise], train_data.Class)\n",
"\n",
"print(\"Pi: \", back_stepwise_lda.priors_)\n",
"back_stepwise_pred = LDA_predict(back_stepwise_lda, data[back_stepwise])\n",
"print(\"Распределение\")\n",
"print(back_stepwise_pred.head())\n",
"data_to_excel[\"Result backward\"] = back_stepwise_pred"
]
},
{
"cell_type": "code",
"execution_count": 54,
"id": "75ae71e1",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Функции Фишера ПДАсИ:\n"
]
},
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" <th>3</th>\n",
" <th>4</th>\n",
" <th>5</th>\n",
" <th>6</th>\n",
" <th>7</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>X1</th>\n",
" <td>-3.771640</td>\n",
" <td>-2.857112</td>\n",
" <td>-3.319760</td>\n",
" <td>7.084028</td>\n",
" <td>16.643087</td>\n",
" <td>25.314092</td>\n",
" <td>-12.520287</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X2</th>\n",
" <td>-6.063287</td>\n",
" <td>-2.536821</td>\n",
" <td>-6.090944</td>\n",
" <td>10.678588</td>\n",
" <td>-12.177667</td>\n",
" <td>-15.673274</td>\n",
" <td>60.676258</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X3</th>\n",
" <td>4.937764</td>\n",
" <td>-3.905883</td>\n",
" <td>0.357235</td>\n",
" <td>0.196979</td>\n",
" <td>7.971863</td>\n",
" <td>7.845312</td>\n",
" <td>-16.315901</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X4</th>\n",
" <td>-1.851979</td>\n",
" <td>-0.914431</td>\n",
" <td>1.743032</td>\n",
" <td>1.116366</td>\n",
" <td>-13.010087</td>\n",
" <td>-15.719599</td>\n",
" <td>21.180522</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X5</th>\n",
" <td>-1.973156</td>\n",
" <td>-3.164523</td>\n",
" <td>1.875504</td>\n",
" <td>6.247621</td>\n",
" <td>3.328162</td>\n",
" <td>6.663230</td>\n",
" <td>-10.050651</td>\n",
" </tr>\n",
" <tr>\n",
" <th>X6</th>\n",
" <td>3.300529</td>\n",
" <td>-1.165738</td>\n",
" <td>-4.143743</td>\n",
" <td>-7.414054</td>\n",
" <td>-3.475812</td>\n",
" <td>69.419338</td>\n",
" <td>2.432008</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Const</th>\n",
" <td>-6.656663</td>\n",
" <td>-4.606433</td>\n",
" <td>-3.962780</td>\n",
" <td>-16.111413</td>\n",
" <td>-29.617907</td>\n",
" <td>-216.567018</td>\n",
" <td>-146.680111</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" 1 2 3 4 5 6 \\\n",
"X1 -3.771640 -2.857112 -3.319760 7.084028 16.643087 25.314092 \n",
"X2 -6.063287 -2.536821 -6.090944 10.678588 -12.177667 -15.673274 \n",
"X3 4.937764 -3.905883 0.357235 0.196979 7.971863 7.845312 \n",
"X4 -1.851979 -0.914431 1.743032 1.116366 -13.010087 -15.719599 \n",
"X5 -1.973156 -3.164523 1.875504 6.247621 3.328162 6.663230 \n",
"X6 3.300529 -1.165738 -4.143743 -7.414054 -3.475812 69.419338 \n",
"Const -6.656663 -4.606433 -3.962780 -16.111413 -29.617907 -216.567018 \n",
"\n",
" 7 \n",
"X1 -12.520287 \n",
"X2 60.676258 \n",
"X3 -16.315901 \n",
"X4 21.180522 \n",
"X5 -10.050651 \n",
"X6 2.432008 \n",
"Const -146.680111 "
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"back_stepwise_lda = LinearDiscriminantAnalysis().fit(train_data[back_stepwise], train_data.Class)\n",
"\n",
"classes = np.unique(train_data.Class)\n",
"groups = [train_data[back_stepwise][train_data.Class == cls] for cls in classes]\n",
"n = [len(g) for g in groups]\n",
"N = sum(n)\n",
"p = train_data[back_stepwise].shape[1]\n",
"\n",
"# Общая (pooled) ковариация как в Statistica\n",
"S_pooled = sum((ni - 1) * np.cov(g, rowvar=False, ddof=1) for g, ni in zip(groups, n)) / (N - len(classes))\n",
"inv_S = np.linalg.inv(S_pooled)\n",
"\n",
"# Средние по классам и априорные вероятности\n",
"means = [g.mean(axis=0) for g in groups]\n",
"priors = np.array(n) / N # можно заменить на np.ones(len(classes))/len(classes), если Statistica = равные априоры\n",
"\n",
"# Классификационные функции (Statistica)\n",
"coef_stat = {}\n",
"const_stat = {}\n",
"\n",
"for cls, mu, p_j in zip(classes, means, priors):\n",
" a = inv_S @ mu\n",
" c = -0.5 * mu.T @ inv_S @ mu + np.log(p_j)\n",
" coef_stat[cls] = a\n",
" const_stat[cls] = c\n",
"\n",
"df_stat = pd.DataFrame(coef_stat, index=[f\"X{i+1}\" for i in range(p)])\n",
"df_stat.loc[\"Const\"] = const_stat\n",
"print(\"Функции Фишера ПДАсИ:\")\n",
"display(df_stat)"
]
},
{
"cell_type": "code",
"execution_count": 55,
"id": "aee7756f",
"metadata": {},
"outputs": [],
"source": [
"# Сохраняем в Excel без номеров строк (индексов)\n",
"data_to_excel.to_excel(OUTPUT_PATH, index=False)\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": ".venv",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.14.0"
}
},
"nbformat": 4,
"nbformat_minor": 5
}