Implementing Logistic Regression from Scratch using Python

 Мы будем использовать набор данных make_classification от sklearn с четырьмя функциями.

import numpy as np 

from numpy import log,dot,exp,shape

import matplotlib.pyplot as plt

from sklearn.datasets import make_classification

X,y = make_classification(n_features=4)

from sklearn.model_selection import train_test_split  

X_tr,X_te,y_tr,y_te = train_test_split(X,y,test_size=0.1)

Стандартизация

Стандартизация — это процесс масштабирования данных вокруг среднего с единичным стандартным отклонением. Это означает, что мы фактически делаем среднее значение атрибута равным нулю, а результирующее распределение имеет стандартное отклонение, равное нулю. Некоторые алгоритмы уязвимы для функций с разными масштабами. Особенно, если мы используем градиентный спуск для оптимизации, тогда модели будет трудно дать точные результаты; например, если в наборе данных есть две характеристики: возраст и зарплата, то характеристика зарплаты с более высоким диапазоном, скорее всего, будет доминировать в результате. Таким образом, рекомендуется стандартизировать данные перед их подачей в алгоритм.

def standardize(X_tr):

    for i in range(shape(X_tr)[1]):

        X_tr[:,i] = (X_tr[:,i] - np.mean(X_tr[:,i]))/np.std(X_tr[:,i])

Инициализация параметров

Наборы данных всегда многомерны. Нам нужно будет использовать матрицы для любого вида вычислений. Итак, для ввода у нас есть две матрицы. Первый предназначен для векторов признаков, а второй - для параметров или весов. Наша первая матрица имеет размер m(x)n, где m - количество наблюдений, а n - размер наблюдений. А второй имеет размерность nx1. Здесь мы добавим столбец смещения из единиц в нашу матрицу векторов признаков и соответствующий термин параметра в вектор весов. Смещение важно, чтобы сделать модель более гибкой.        

def initialize(self,X):

        weights = np.zeros((shape(X)[1]+1,1))

        X = np.c_[np.ones((shape(X)[0],1)),X]

        return weights,X

 Функция стоимости

Функция стоимости или функция потерь — это та функция, которая описывает, насколько расчетное значение отклоняется от фактического значения. Линейная регрессия использует наименьший квадрат ошибки в качестве функции стоимости. Но функция наименьших квадратов ошибок для логистической регрессии невыпукла. При выполнении градиентного спуска шансов, что мы застрянем в локальном минимуме, больше. Поэтому вместо этого мы используем логарифмическую потерю в качестве функции стоимости.

h(x) = is the sigmoid function       

Градиентный спуск

Следующий шаг — градиентный спуск. Градиентный спуск — это алгоритм оптимизации, который отвечает за изучение наиболее подходящих параметров. 

Градиенты представляют собой вектор производной 1-го порядка функции стоимости. Это направление наискорейшего подъема или максимума функции. Для градиентного спуска мы движемся в направлении, противоположном градиентам. Мы будем обновлять веса на каждой итерации до сходимости. 

Что такое градиентный спуск?

Здесь альфа — это размер шага, отвечающий за скорость сходимости к глобальному минимуму. Если размер шага слишком мал, он будет сходиться медленно, но если он слишком велик, он может выйти за пределы минимума при спуске. Дифференцируя функцию стоимости, мы получаем выражение градиентного спуска

Это векторизованная форма выражения градиентного спуска, которую мы будем использовать в нашем коде.

def fit(self,X,y,alpha=0.001,iter=100):

        params,X = self.initialize(X)

        cost_list = np.zeros(iter,)

        for i in range(iter):

            params = params - alpha * dot(X.T, self.sigmoid(dot(X,params)) - np.reshape(y,(len(y),1)))

            cost_list[i] = cost(params)

        self.params = params

        return cost_list

Прогноз

Все, что мы сделали до сих пор, предназначено для этого шага. Мы обучили модель на обучающем наборе данных, и теперь мы будем использовать изученные параметры для прогнозирования невидимых данных.

def predict(self,X):

        z = dot(self.initialize(X)[1],self.weights)

        lis = []

        for i in self.sigmoid(z):

            if i>0.5:

                lis.append(1)

            else:

                lis.append(0)

        return lis

F1_score

Теперь, когда мы закончили прогнозирование, мы перейдем к разделу F1-показателей, где мы измерим, насколько хорошо наша модель прогнозирует невидимые данные. F1_score — это надежная метрика для оценки эффективности моделей классификации, а математически F1-score — это гармоническое среднее значение точности и полноты.

Напомним, точность, оценка F1 - простое объяснение метрики, машинное обучение

точность = точность - это количество истинных срабатываний по сравнению с суммой истинных срабатываний и ложных срабатываний.

precision = TP/(TP+FP)

отзыв = отзыв — это количество истинных положительных результатов в сумме истинных положительных и ложных отрицательных результатов.

recall = TP/(TP+FN)

def f1_score(y,y_hat):

    tp,tn,fp,fn = 0,0,0,0

    for i in range(len(y)):

        if y[i] == 1 and y_hat[i] == 1:

            tp += 1

        elif y[i] == 1 and y_hat[i] == 0:

            fn += 1

        elif y[i] == 0 and y_hat[i] == 1:

            fp += 1

        elif y[i] == 0 and y_hat[i] == 0:

            tn += 1

    precision = tp/(tp+fp)

    recall = tp/(tp+fn)

    f1_score = 2*precision*recall/(precision+recall)

    return f1_score

Теперь, когда мы закончили с каждой частью, мы собираем все вместе

(.env) boris@boris-All-Series:~/GRADDECENTLGREG$ cat gradDescLogReg.py

import numpy as np 

from numpy import log,dot,exp,shape

import matplotlib.pyplot as plt

from sklearn.datasets import make_classification

X,y = make_classification(n_features=4)

from sklearn.model_selection import train_test_split  

X_tr,X_te,y_tr,y_te = train_test_split(X,y,test_size=0.1)

def standardize(X_tr):

    for i in range(shape(X_tr)[1]):

        X_tr[:,i] = (X_tr[:,i] - np.mean(X_tr[:,i]))/np.std(X_tr[:,i])

def F1_score(y,y_hat):

    tp,tn,fp,fn = 0,0,0,0

    for i in range(len(y)):

        if y[i] == 1 and y_hat[i] == 1:

            tp += 1

        elif y[i] == 1 and y_hat[i] == 0:

            fn += 1

        elif y[i] == 0 and y_hat[i] == 1:

            fp += 1

        elif y[i] == 0 and y_hat[i] == 0:

            tn += 1

    precision = tp/(tp+fp)

    recall = tp/(tp+fn)

    f1_score = 2*precision*recall/(precision+recall)

    return f1_score

class LogidticRegression:

    def sigmoid(self,z):

        sig = 1/(1+exp(-z))

        return sig

    def initialize(self,X):

        weights = np.zeros((shape(X)[1]+1,1))

        X = np.c_[np.ones((shape(X)[0],1)),X]

        return weights,X

    def fit(self,X,y,alpha=0.001,iter=400):

        weights,X = self.initialize(X)

        def cost(theta):

            z = dot(X,theta)

            cost0 = y.T.dot(log(self.sigmoid(z)))

            cost1 = (1-y).T.dot(log(1-self.sigmoid(z)))

            cost = -((cost1 + cost0))/len(y)

            return cost

        cost_list = np.zeros(iter,)

        for i in range(iter):

            weights = weights - alpha*dot(X.T,self.sigmoid(dot(X,weights))-np.reshape(y,(len(y),1)))

            cost_list[i] = cost(weights)

        self.weights = weights

        return cost_list

    def predict(self,X):

        z = dot(self.initialize(X)[1],self.weights)

        lis = []

        for i in self.sigmoid(z):

            if i>0.5:

                lis.append(1)

            else:

                lis.append(0)

        return lis

standardize(X_tr)

standardize(X_te)

obj1 = LogidticRegression()

model= obj1.fit(X_tr,y_tr)

y_pred = obj1.predict(X_te)

y_train = obj1.predict(X_tr)

#Let's see the f1-score for training and testing data

f1_score_tr = F1_score(y_tr,y_train)

f1_score_te = F1_score(y_te,y_pred)

print("Logistic regression fares")

print(f1_score_tr)

print(f1_score_te)

print("Logistic regression fares in comparison to sklearn’s logistic regression")

from sklearn.linear_model import LogisticRegression

from sklearn.metrics import f1_score

model = LogisticRegression().fit(X_tr,y_tr)

y_pred = model.predict(X_te)

print(f1_score(y_te,y_pred))

Update fit(self,X,y,alpha=0.001,iter=400) as follows :

    def fit(self,X,y,alpha=0.001,iter=400):

        weights,X = self.initialize(X)

        def cost(theta):

            z = dot(X,theta)

            cost0 = y.T.dot(log(self.sigmoid(z)))

            cost1 = (1-y).T.dot(log(1-self.sigmoid(z)))

            cost = -((cost1 + cost0))/len(y)

            return cost

        cost_list = np.zeros(iter,)

        for i in range(iter):

            weights = weights - alpha*dot(X.T,self.sigmoid(dot(X,weights))-np.reshape(y,(len(y),1)))

            cost_list[i] = cost(weights)

        self.weights = weights

        # Plotting Cost curve 

        xs = [i for i in range(iter)]

        ys = [cost_list[i] for i in range(iter)]

        plt.plot(xs,ys)

        plt.xlabel('xs - axis')

        plt.ylabel('ys - axis')

        plt.title('Cost function plot')

        plt.draw()

        return cost_list

.  .  .  . .  .  .  .  .  .  .  . 

standardize(X_tr)

standardize(X_te)

obj1 = LogidticRegression()

model= obj1.fit(X_tr,y_tr)

y_pred = obj1.predict(X_te)

y_train = obj1.predict(X_tr)

#Let's see the f1-score for training and testing data

f1_score_tr = F1_score(y_tr,y_train)

f1_score_te = F1_score(y_te,y_pred)

print("Logistic regression fares")

print(f1_score_tr)

print(f1_score_te)

plt.show()

Теперь выполним тот же код выполняя плот убывания функции потерь