Изучая и внедряя курс Эндрю Нг в Python, я обнаружил, что двойная проверка моей работы и ответов с людьми в Интернете была чрезвычайно полезной. Я нашел эту реализацию особенно полезной, но она закончилась после Упражнения 5. Итак, вот мой минимальный рабочий пример Упражнения 5, Регуляризованная линейная регрессия и смещение против дисперсии. Обязательно попробуйте исходные инструкции, а затем проверяйте, только если вы застряли!

1. Регулярная линейная регрессия

1.1 Визуализация набора данных

Набор данных разделен на 3 части: обучение (n = 12), перекрестная проверка (n = 21) и тест (n = 21). Начните с построения графика данных: вода, вытекающая из плотины (y), в зависимости от изменения уровня воды (x).

import numpy as np
import matplotlib.pyplot as plt
from scipy.io import loadmat
import scipy.optimize as opt
df = loadmat("C:/Users/name/Desktop/ex5/ex5data1.mat")
X, y = df["X"], df["y"]
Xval, yval = df["Xval"], df["yval"]
Xtest, ytest = df["Xtest"], df["ytest"]
def plotData():
    plt.figure(figsize=(8, 5))
    plt.ylabel("Water flowing out of the dam")
    plt.xlabel("Change in water level")
    plt.plot(X, y, "rx")
    plt.grid(True)

plotData()

Давайте также определим наши входные данные со столбцом единиц, включенным для термина смещения тета (0).

m = len(y)
X_ones = np.insert(X, 0, 1, axis=1)
Xval_ones = np.insert(Xval, 0, 1, axis=1)
Xtest_ones = np.insert(Xtest, 0, 1, axis=1)

1.2 Функция стоимости регуляризованной линейной регрессии

Теперь напишите регуляризованную функцию стоимости, следя за тем, чтобы не регуляризовать тета(0).

def computeCost(theta, X, y, lmbda):
    unreg_cost = (1 / (2*m)) * np.sum(np.square((X @ theta) - y))
    reg = (lmbda / (2*m)) * np.sum(np.square(theta[1:len(theta)]))
    cost = unreg_cost + reg
    return cost

Установив нашу начальную тету, мы получим начальную стоимость 303,951.

theta = np.ones((2, 1))
computeCost(theta, X_ones, y, 0)

1.3 Градиент регуляризованной линейной регрессии

Давайте сделаем то же самое для начального градиента, получив [-15,303, 598,167]

def computeGradient(theta, X, y, lmbda):
    grad = (1 / m) * (X.T @ ((X @ theta) - y))  # 2x12 @ 12x1 = 2x1
    reg_grad = (lmbda / m) * theta[1:len(theta)]
    grad[1:] = grad[1:] + reg_grad
    return grad

computeGradient(theta, X_ones, y, 0)

1.4 Подгонка линейной регрессии

Теперь давайте вычислим оптимальные значения тета, используя opt.fmin_cg.

def optimizeTheta(theta, X, y, lmbda, print_output=True):
    opt_theta = opt.fmin_cg(f=computeCost, x0=theta.flatten(),
                            fprime=computeGradient,
                            args=(X, y.flatten(), lmbda),
                            disp=print_output, epsilon=1.49e-12,
                            maxiter=1000)
    opt_theta = opt_theta.reshape((theta.shape[0], 1))
    return opt_theta

Помните, что нам нужно сгладить определенные входные данные, чтобы это работало правильно. Расчет дает оптимизированную тета [13 088, 0,368] и оптимизированную стоимость 22,384. График данных показывает, что это не очень подходит. Запомните формулу функции гипотезы для построения линии.

fit_theta = optimizeTheta(theta, X_ones, y, 0)
computeCost(fit_theta, X_ones, y, 0)
plotData()
plt.plot(X, (X_ones @ fit_theta))

2. Смещение-дисперсия

Здесь мы построим обучение и тестовые ошибки на кривой обучения, чтобы исследовать систематическую ошибку и дисперсию.

2.1 Кривые обучения

Кривая обучения отображает ошибку как функцию размера обучающей выборки. Взяв подмножество обучающего набора, мы вычисляем тета-значения, вычисляем соответствующую ошибку обучения и ошибку перекрестной проверки, затем повторяем увеличение размера обучающего набора, чтобы увидеть форму кривой обучения. Внимательно рассмотрите каждую строку, чтобы разобрать ее.

def learningCurve(theta, X_ones, y, lmbda):
    xvalue, train_error_set, cv_error_set = [], [], []
    for i in range(1, m+1):
        X_train = X_ones[:i]
        y_train = y[:i]
        xvalue.append(y_train.shape[0])
        theta_train = optimizeTheta(theta, X_train, y_train, lmbda,
                                    print_output=False)
        train_error_val = computeCost(theta_train.reshape((-1, 1)),
                                      X_train, y_train, lmbda)
        train_error_set.append(train_error_val)
        cv_error_val = computeCost(theta_train.reshape((-1, 1)), Xval_ones,
                                   yval, lmbda)
        cv_error_set.append(cv_error_val)
    plt.figure(figsize=(8, 5))
    plt.plot(xvalue, train_error_set, label="Training")
    plt.plot(xvalue, cv_error_set, label="Cross validation")
    plt.legend(loc="upper right")
    plt.xlabel("No. of training examples")
    plt.ylabel("Error")
    plt.title("Learning curve for linear regression")
    plt.grid(True)
learningCurve(theta, X_ones, y, 1)

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

3. Полиномиальная регрессия

Здесь мы добавляем больше функций, используя большие мощности существующей функции x (уровень воды). Давайте создадим функцию, которая берет матрицу с уже включенными единицами для члена тета (0) и возвращает нашу матрицу более высокой размерности.

def polyFeature(X_ones, p):
    X_poly = X_ones.copy()
    for i in range(1, p):
        X_poly = np.insert(X_poly, X_poly.shape[1],
                           np.power(X_poly[:, 1], i+1), axis=1)
    return X_poly

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

def featureNormalize(X_ones):
    X_norm = X_ones.copy()
    means = np.mean(X_norm, axis=0)
    stds = np.std(X_norm, axis=0)
    X_norm[:, 1:] = (X_norm[:, 1:] - means[1:]) / stds[1:]
    return X_norm, means, stds

Теперь посчитаем нашу начальную стоимость. Мы выбрали степень 5 вместо 8, потому что она лучше иллюстрирует проблему переобучения.

test_power = 5
Xp_ones = polyFeature(X_ones, test_power)
Xp_ones_norm, Xp_mean, Xp_std = featureNormalize(Xp_ones)
thetap = np.ones((Xp_ones_norm.shape[1], 1))
L_poly = 0
computeCost(thetap, Xp_ones_norm, y, L_poly)
# Initial cost = 104.835

Мы также можем увидеть, как модель работает с наборами данных проверки и тестирования.

Xval_poly_norm, d1, d2 = featureNormalize(polyFeature(Xval_ones, test_power))
Xtest_poly_norm, d1, d2 = featureNormalize(polyFeature(Xtest_ones, test_power))
fit_thetap = optimizeTheta(thetap, Xp_ones_norm, y, L_poly)
computeCost(fit_thetap, Xp_ones_norm, y, L_poly)  # J = 0.208
computeCost(fit_thetap, Xval_poly_norm, yval, L_poly)  # J = 14.313
computeCost(fit_thetap, Xtest_poly_norm, ytest, L_poly)  # J = 11.226

Давайте попробуем построить нашу полиномиальную гипотезу и кривую обучения.

def plotFit(theta, X_mean, X_std):
    xvalue = np.linspace(-55, 55, 50)
    xvalue_ones = np.insert(xvalue.reshape((-1, 1)), 0, 1, axis=1)
    xpoly = polyFeature(xvalue_ones, test_power)
    xpoly[:, 1:] = (xpoly[:, 1:] - X_mean[1:]) / X_std[1:]
    plotData()
    plt.plot(xvalue, (xpoly @ theta), "b--")
    plt.title("Predicted water flowing out for polynomial of order
              5")

plotFit(fit_thetap, Xp_mean, Xp_std)

И для кривой обучения

plotPolyLearningCurve(thetap, X_ones, y, L_poly)

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

3.2 Настройка параметра регуляризации

Затем мы можем настроить параметр регуляризации, чтобы увидеть, как это влияет на дисперсию смещения. Установка lambda = 1 дает следующую подгонку.

L_poly = 1

При регуляризации ошибка обучения и перекрестной проверки сходится к аналогичному относительно низкому значению.

3.4 Выбор лямбды с использованием набора перекрестной проверки

Кроме того, можно оценить, какое значение лямбда дает минимальную ошибку для обучения и перекрестной проверки.

def bestLambda():
    lmbdas = [0, 0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1, 3, 10]
    train_cost_set, cv_cost_set = [], []
    for i in range(len(lmbdas)):
        fit_thetap = optimizeTheta(thetap, Xp_ones_norm, y,
                                   lmbdas[i])
        train_cost = computeCost(fit_thetap, Xp_ones_norm, y,
                                 lmbdas[i])
        train_cost_set.append(train_cost)
        cv_cost = computeCost(fit_thetap, Xval_poly_norm, yval,
                              lmbdas[i])
        cv_cost_set.append(cv_cost)
    plt.figure(figsize=(8, 5))
    plt.plot(lmbdas, train_cost_set, label="Training")
    plt.plot(lmbdas, cv_cost_set, label="Cross validation")
    plt.legend(loc="upper left")
    plt.xlabel("Lambda")
    plt.ylabel("Error")
    plt.title("Comparison of cost and value of lambda")
    plt.grid(True)

bestLambda()

В нашем случае мы видим, что значение лямбда, которое соответствует нашему набору для обучения и перекрестной проверки, составляет около 1.

На этом пока все, дайте мне знать, если у вас есть предложения по улучшению!