Einführung

Viele Leute werden Erste Mustererkennung als Einführung in maschinelles Lernen lesen. Ich war einer von ihnen, und ich dachte, ich hätte es einmal gelesen und vollständig verstanden, aber wenn ich es noch einmal lese, verstehe ich immer noch nichts [^ 1], also habe ich beschlossen, es während der Implementierung erneut zu lesen. Dieser Artikel fasst den Code zusammen, der das im Buch in Python beschriebene Ausführungsbeispiel reproduziert (die Betriebsumgebung des beschriebenen Codes ist numpy 1.16.2, pandas 0.24.2). Diesmal ist Kapitel 4, "Wahrscheinlichkeitsmodell und Diskriminierungsfunktion".

Ausführungsbeispiel 4.1 Standardisierung

Standardisierung ist die Subtraktion des Durchschnittswertes und die Division durch die Standardabweichung. Unten ist der Code.

Laden der erforderlichen Bibliotheken

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style('whitegrid')

Irisdaten lesen

from sklearn.datasets import load_iris

iris = load_iris()
data = pd.DataFrame(iris['data'], columns=iris.feature_names)
target = pd.Series(iris['target']).map({0: 'setosa', 1:'versicolor', 2:'virginica'})

Durchschnittliche Länge und Breite der Irisblüten und der Kovarianzmatrix

print('mu = ', np.mean(data[['petal length (cm)', 'petal width (cm)']]))
print('Sigma = ', np.cov(data[['petal length (cm)', 'petal width (cm)']], rowvar=0))

#Ausgabe=>
# mu = 
#  petal length (cm)    3.758000
# petal width (cm)     1.199333
# dtype: float64
# Sigma = 
#  [[3.11627785 1.2956094 ]
#  [1.2956094  0.58100626]]

Streudiagramm vor der Standardisierung - Abb. 4.2 (a)

for t in target.unique():
    data_tmp = data[target==t]
    plt.scatter(data_tmp['petal length (cm)'], data_tmp['petal width (cm)'], label=t)
plt.legend()
plt.xlabel('petal length (cm)')
plt.ylabel('petal width (cm')
plt.xlim([1, 7])
plt.ylim([-1, 4])
plt.show()

Ausbreitungsdiagramm nach Standardisierung - Abb. 4.2 (b)

def standardize(X):
    return (X - np.mean(X)) / np.std(X)

std_data = standardize(data[['petal length (cm)', 'petal width (cm)']])
for t in target.unique():
    data_tmp = std_data[target==t]
    plt.scatter(data_tmp['petal length (cm)'], data_tmp['petal width (cm)'], label=t)
plt.legend()
plt.xlabel('petal length (cm)')
plt.ylabel('petal width (cm)')
plt.xlim([-2, 2])
plt.ylim([-2, 2])
plt.show()

Ausführungsbeispiel 4.2 Nicht korreliert

Rotieren über eine Matrix von Eigenvektoren, das heißt Unkorrelation. Unten ist der Code.

Ausbreitungsdiagramm nach Unkorrelation - Abb. 4.4 (b)

def no_correlation(X):
    cov_mat = np.cov(X, rowvar=0)
    eig_vals, eig_vecs = np.linalg.eig(cov_mat)
    return np.dot(X, eig_vecs)

no_corr_data = no_correlation(data[['petal length (cm)', 'petal width (cm)']])
for t in target.unique():
    data_tmp = no_corr_data[target==t]
    plt.scatter(data_tmp[:, 0], data_tmp[:, 1], label=t)
plt.legend()
plt.xlabel('petal length (cm)')
plt.ylabel('petal width (cm)')
plt.xlim([0, 8])
plt.ylim([-3, 3])
plt.show()

In diesem Ausführungsbeispiel ist das Vorzeichen der y-Komponente des Eigenvektors gegenüber dem im Buch beschriebenen umgekehrt, sodass auch die Richtung der y-Achse des Streudiagramms umgekehrt ist.

Sigma = np.cov(data[['petal length (cm)', 'petal width (cm)']], rowvar=0)
Lambda, S = np.linalg.eig(Sigma)
print('Lambda = \n', Lambda)
print('S = \n', S)
#Ausgabe=>
# Lambda = 
#  [3.66123805 0.03604607]
# S = 
#  [[ 0.92177769 -0.38771882]
#  [ 0.38771882  0.92177769]]

Kovarianzmatrix nach Unkorrelation

print('Sigma = \n', np.cov(no_corr_data, rowvar=0))
#Ausgabe=>
# Sigma = 
#  [[3.66123805e+00 6.01365790e-16]
#  [6.01365790e-16 3.60460707e-02]]

Ausführungsbeispiel 4.3 Bleaching

Bleaching ist Zentralisierung, Unkorrelation und Normalisierung der Standardabweichung auf 1. Unten ist der Code.

Streudiagramm nach dem Aufhellen - Abb. 4.6 (b)

def whitening(X):
    cov_mat = np.cov(X, rowvar=0)
    eig_vals, eig_vecs = np.linalg.eig(cov_mat)
    diag_sqrt_eig_vals = np.diag(np.sqrt(eig_vals))

    return np.dot(np.dot(X-np.mean(X, axis=0), eig_vecs), np.linalg.inv(diag_sqrt_eig_vals.T))

whitened_data = whitening(data[['petal length (cm)', 'petal width (cm)']])
for t in target.unique():
    data_tmp = whitened_data[target==t]
    plt.scatter(data_tmp[:, 0], data_tmp[:, 1], label=t)
plt.legend()
plt.xlabel('petal length (cm)')
plt.ylabel('petal width (cm)')
plt.xlim([-3, 3])
plt.ylim([-3, 4])
plt.show()

Kovarianzmatrix nach dem Aufhellen

print('Sigma = \n', np.cov(whitened_data, rowvar=0))
#Ausgabe=>
# Sigma = 
#  [[1.00000000e+00 1.64411249e-15]
#  [1.64411249e-15 1.00000000e+00]]

Ausführungsbeispiel 4.4 Sekundäre Diskriminanzfunktion / lineare Diskriminanzfunktion

Eine Funktion, die eine Klasse identifiziert, vorausgesetzt, dass die klassenbedingten Wahrscheinlichkeiten einer Normalverteilung folgen. Unten ist der Code.

Lesen von Pima Indian Daten

Die Daten stammen aus R-Datensätzen.

pima_train = pd.read_csv('https://raw.githubusercontent.com/vincentarelbundock/Rdatasets/master/csv/MASS/Pima.tr.csv')
pima_test = pd.read_csv('https://raw.githubusercontent.com/vincentarelbundock/Rdatasets/master/csv/MASS/Pima.te.csv')

Streudiagramm der Trainingsdaten - Abb. 4.8

for t in ['Yes', 'No']:
    pima_tmp = pima_train[pima_train['type']==t]
    plt.scatter(pima_tmp['glu'], pima_tmp['bmi'], label=t)
plt.legend()
plt.xlabel('glu')
plt.ylabel('bmi')
plt.xlim([50, 200])
plt.ylim([15, 50])
plt.show()

Implementierung der sekundären Identifikationsfunktion

class QuadraticDiscriminantFunction:
    def __init__(self, threshold=0):
        self._S = None
        self._c = None
        self._F = None
        self._threshold = threshold
    
    def fit(self, X, Y):
        X_pos = X[Y==1]
        X_neg = X[Y==0]

        mu_pos = np.mean(X_pos, axis=0)
        mu_neg = np.mean(X_neg, axis=0)
        Sigma_pos = np.cov(X_pos, rowvar=0)
        Sigma_neg = np.cov(X_neg, rowvar=0)
        p_pos = len(Y[Y==1])/len(Y)
        p_neg = len(Y[Y==0])/len(Y)
        
        self._S = np.linalg.inv(Sigma_pos) - np.linalg.inv(Sigma_neg)
        self._c = np.dot(mu_neg, np.linalg.inv(Sigma_neg)) - np.dot(mu_pos, np.linalg.inv(Sigma_pos))
        self._F = np.dot(np.dot(mu_pos, np.linalg.inv(Sigma_pos)), mu_pos) - np.dot(np.dot(mu_neg, np.linalg.inv(Sigma_neg)), mu_neg)\
                    + np.log(np.linalg.det(Sigma_pos)/np.linalg.det(Sigma_neg))\
                    - 2*np.log(p_pos/p_neg)
    
    def predict(self, X):
        xSx = np.diag(np.dot(np.dot(X, self._S), X.T))
        cx = np.dot(self._c, X.T)
        Y_hat = xSx + 2*cx + self._F
        return np.where(Y_hat < self._threshold, 1, 0)

Implementierung der linearen Diskriminanzfunktion

class LinearDiscriminantFunction:
    def __init__(self, threshold=0):
        self._c = None
        self._F = None
        self._threshold = threshold
    
    def fit(self, X, Y):
        X_pos = X[Y==1]
        X_neg = X[Y==0]

        mu_pos = np.mean(X_pos, axis=0)
        mu_neg = np.mean(X_neg, axis=0)
        Sigma_pos = np.cov(X_pos, rowvar=0)
        Sigma_neg = np.cov(X_neg, rowvar=0)
        p_pos = len(Y[Y==1])/len(Y)
        p_neg = len(Y[Y==0])/len(Y)
        Sigma_pool = p_pos*Sigma_pos + p_neg*Sigma_neg
        
        self._c = np.dot(mu_neg, np.linalg.inv(Sigma_pool)) - np.dot(mu_pos, np.linalg.inv(Sigma_pool))
        self._F = np.dot(np.dot(mu_pos, np.linalg.inv(Sigma_pool)), mu_pos) - np.dot(np.dot(mu_neg, np.linalg.inv(Sigma_pool)), mu_neg)\
                    - 2*np.log(p_pos/p_neg)
    
    def predict(self, X):
        cx = np.dot(self._c, X.T)
        Y_hat = 2*cx + self._F
        return np.where(Y_hat < self._threshold, 1, 0)

Lernen und Vorhersagen zum Zeichnen von Bayes'schen Identifikationsgrenzen

X_train = pima_train[['glu', 'bmi']].values
Y_train = pima_train['type'].map({'Yes':1, 'No':0})

qdf = QuadraticDiscriminantFunction(threshold=0)
qdf.fit(X_train, Y_train)
ldf = LinearDiscriminantFunction(threshold=0)
ldf.fit(X_train, Y_train)

X1_min, X1_max = np.min(X_train[:, 0])-1, np.max(X_train[:, 0])+1
X2_min, X2_max = np.min(X_train[:, 1])-1, np.max(X_train[:, 1])+1
X1, X2 = np.meshgrid(np.arange(X1_min, X1_max, 1), np.arange(X2_min, X2_max, 0.5))
X = np.c_[X1.ravel(), X2.ravel()]

Y_qdf_pred = qdf.predict(X)
Y_ldf_pred = ldf.predict(X)

Durch die sekundäre Diskriminanzfunktion erhaltene Bayes-Diskriminierungsgrenze - Abb. 4.9 (a)

plt.contour(X1, X2, Y_qdf_pred.reshape(X1.shape), colors='limegreen')
for t in ['Yes', 'No']:
    pima_tmp = pima_train[pima_train['type']==t]
    plt.scatter(pima_tmp['glu'], pima_tmp['bmi'], label=t)
plt.legend()
plt.xlabel('glu')
plt.ylabel('bmi')
plt.xlim([50, 200])
plt.ylim([15, 50])
plt.show()

Fehlerrate der sekundären Diskriminanzfunktion

X_test = pima_test[['glu', 'bmi']].values
Y_test = pima_test['type'].map({'Yes':1, 'No':0})

Y_train_qdf_pred = qdf.predict(X_train)
print('error rate (train) = {:.1f}%'.format(np.mean(Y_train_qdf_pred != Y_train)*100))
Y_test_qdf_pred = qdf.predict(X_test)
print('error rate (test) = {:.1f}%'.format(np.mean(Y_test_qdf_pred != Y_test)*100))

#Ausgabe=>
# error rate (train) = 24.5%
# error rate (test) = 23.8%

Durch die lineare Diskriminanzfunktion erhaltene Bayes-Diskriminierungsgrenze - Abb. 4.9 (b)

plt.contour(X1, X2, Y_ldf_pred.reshape(X1.shape), colors='limegreen')
for t in ['Yes', 'No']:
    pima_tmp = pima_train[pima_train['type']==t]
    plt.scatter(pima_tmp['glu'], pima_tmp['bmi'], label=t)
plt.legend()
plt.xlabel('glu')
plt.ylabel('bmi')
plt.xlim([50, 200])
plt.ylim([15, 50])
plt.show()

Fehlerrate der linearen Diskriminanzfunktion

Y_train_ldf_pred = ldf.predict(X_train)
print('error rate (train) = {:.1f}%'.format(np.mean(Y_train_ldf_pred != Y_train)*100))
Y_test_ldf_pred = ldf.predict(X_test)
print('error rate (test) = {:.1f}%'.format(np.mean(Y_test_ldf_pred != Y_test)*100))
#Ausgabe=>
# error rate (train) = 24.0%
# error rate (test) = 22.0%

ROC-Kurve der quadratischen Diskriminanzfunktion / linearen Diskriminanzfunktion - Abb. 4.10

qdf_tpr_list = []
qdf_fpr_list = []
ldf_tpr_list = []
ldf_fpr_list = []

for th in np.arange(-9, 7, 0.1):
    qdf = QuadraticDiscriminantFunction(threshold=th)
    qdf.fit(X_train, Y_train)
    Y_qdf_pred = qdf.predict(X_test)

    tpr_qdf = len(Y_test[(Y_test==0)&(Y_qdf_pred==0)]) / len(Y_test[Y_test==0])
    qdf_tpr_list.append(tpr_qdf)
    fpr_qdf = len(Y_test[(Y_test==1)&(Y_qdf_pred==0)]) / len(Y_test[Y_test==1])
    qdf_fpr_list.append(fpr_qdf)
                        
    ldf = LinearDiscriminantFunction(threshold=th)
    ldf.fit(X_train, Y_train)
    Y_ldf_pred = ldf.predict(X_test)
    
    tpr_ldf = len(Y_test[(Y_test==0)&(Y_ldf_pred==0)]) / len(Y_test[Y_test==0])
    ldf_tpr_list.append(tpr_ldf)
    fpr_ldf = len(Y_test[(Y_test==1)&(Y_ldf_pred==0)]) / len(Y_test[Y_test==1])
    ldf_fpr_list.append(fpr_ldf)

plt.plot(qdf_fpr_list, qdf_tpr_list, label='Quadratic')
plt.plot(ldf_fpr_list, ldf_tpr_list, label='Linear')
plt.legend()
plt.xlabel('false positive ratio')
plt.ylabel('true positive ratio')
plt.show()

abschließend

Ich habe versucht, das Ausführungsbeispiel von Kapitel 4 von Hajipata in Python zu reproduzieren. Was die Richtigkeit des Codes betrifft, stimmt die gezeichnete Figur ungefähr mit der Figur im Buch überein, also Yoshi! Wir überprüfen nur den Abschluss. Wenn Sie also Fehler finden, würden wir uns freuen, wenn Sie eine Bearbeitungsanfrage stellen könnten.