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

Question

Прежде всего, прости меня, если название звучит немного запутанно.Английский не является моим родным языком, поэтому, хотя я бегло говорю на нем, некоторые технические термины могут быть неверными.Дайте мне знать, если и как я могу улучшить его.

Мне дали задание на курсе машинного обучения с Python.У меня есть таблица значений (я думаю, что я могу назвать ее матрицей), которую можно визуализировать с помощью следующего кода:

import numpy as np
from numpy.linalg import inv
import pandas as pd
import matplotlib.pyplot as plt


h = np.array([x for x in range(0,5500,500)])
v = np.array([x for x in range(0,65,5)])

print(h.shape,v.shape)

y = np.array([83.4, 78.7, 74.1, 69.7, 65.4, 61.3, 57.4, 53.6, 50.0, 46.5, 43.2,
 75.1, 70.8, 66.6, 63.3, 59.5, 55.7, 52.7, 49.2, 46.4, 43.2, 40.5,
 69.7, 65.7, 62.2, 58.8, 55.2, 52.0, 48.9, 45.9, 43.1, 40.3, 37.6,
 64.4, 61.0, 57.6, 54.2, 51.0, 48.0, 45.1, 42.4, 39.7, 37.1, 34.7,
 59.9, 56.6, 53.3, 50.3, 47.5, 44.6, 41.9, 39.3, 36.8, 34.4, 32.1,
 56.1, 53.0, 50.1, 47.2, 44.5, 41.9, 39.3, 36.9, 34.6, 32.3, 30.2,
 53.3, 50.4, 47.5, 44.8, 42.2, 39.8, 37.4, 35.1, 32.8, 30.7, 28.6,
 50.9, 48.1, 45.4, 42.8, 40.3, 38.0, 35.7, 33.4, 31.3, 29.3, 27.3,
 48.7, 46.0, 43.4, 40.9, 38.6, 36.3, 34.1, 31.9, 29.9, 27.9, 26.1,
 46.4, 43.9, 41.4, 39.0, 36.8, 34.6, 32.4, 30.4, 28.5, 26.6, 24.8,
 44.1, 41.7, 39.3, 37.1, 34.9, 32.8, 30.7, 28.8, 26.9, 25.2, 23.4,
 41.7, 39.4, 37.2, 34.9, 32.9, 30.9, 29.0, 27.2, 25.4, 23.7, 22.0,
 39.3, 37.2, 34.9, 32.9, 31.0, 29.1, 27.3, 25.6, 23.8, 22.3, 20.7])

# Confirm lengths all match
y.shape
len(y) == len(h)*len(v)

###################################
# Matrix visualization with Pandas
dataframe = pd.DataFrame(y.reshape(13,11), index=v, columns=h)

print(dataframe)

Для контекстуализации, строки - это скорость (м / с), а столбцы -высота (м) дирижабля.Значения - это тяга (N - Ньютон) самолета.

Задача, которую нужно решить:

Приблизительная матрица с функцией, построенной из элементов следующей базы:{1, v, h, v², h², vh, v³, h³, v²h, vh², v²h², v³h, vh³}

Прежде всего, я не совсем понимаю, чтолежит в основе вопроса.Это значения бета в многомерной линейной регрессии, верно?Каковы преимущества выполнения этого с таким количеством бета-версий?

С помощью друга я смог найти следующее решение:

base = []
i = 0
j = 0

for i in range(0,len(v)):
    for j in range(0,len(h)):
        base.append([1, v[i], h[j], v[i]**2 , h[j]**2, v[i] * h[j],
                    v[i]**3, h[j]**3, (v[i]**2) * h[j], v[i] * (h[j]**2),
                    (v[i]**2) * (h[j]**2), (v[i]**3) * h[j], v[i] * (h[j]**3)])

base = np.array(base)
base.shape
base_df = pd.DataFrame(base)
print(base_df)


base_tp = np.transpose(base)
d = np.dot(base_tp, base)
inv_d = inv(d)
# print(inv_d.shape, base_tp.shape)
x = np.dot(inv_d, base_tp)
beta = np.dot(x, y)

print(beta.shape,beta)

И это то, что ядо сих пор.Я верю, что это правильно, но я не буду высовываться из-за этого.Должен ли я изменить исходный массив Y с y.reshape(143,1), прежде чем продолжить?Если так, почему?

Answer 1

Прежде всего, я не очень понимаю, что лежит в основе вопроса.Это значения бета в многомерной линейной регрессии, верно?Каковы преимущества выполнения этого с таким количеством бета-версий?

Задача, стоящая перед нами, состоит в том, чтобы предположить, что значения y могут быть смоделированы в терминах значений v и hв качестве линейной комбинации указанных базовых элементов.Упомянутые вами бета - это коэффициенты, с которыми участвуют отдельные базисные элементы.

Формально ищется решение для A . betas = y, где y входит как одномерный массив (вектор) и матричные элементы A_ij содержат j-й базовый элемент, оцененный для значений v и h, соответствующих y_i (т. е. значениям v и h, для которых y_iбыл «измерен»).

Поскольку проблема переопределена (матрица A имеет больше строк, чем столбцов), нужно как-то ослабить определение «решения» - одним из возможных кандидатов является решение наименьших квадратовкоторый может быть формально получен (как вы уже реализовали в своем вопросе) как решение системы уравнений (A^T.A) . betas = A^T . y.

Однако матрица (A^T.A) имеет тенденцию быть довольно плохо обусловленной, поэтому вместоЯвно вычисляя обратное, я бы предложил либо использовать решатель наименьших квадратов, предоставленный numpy:

beta, res, rank, s = np.linalg.lstsq(base, y, rcond = -1)

, либо использовать псевдообратное instead:

beta = np.dot(np.linalg.pinv(base), y)

Теперь, если вы проверите максимальную относительную погрешность полученной «подгонки» как:

np.max(np.abs((np.dot(base, beta) - y) / y))

, то подход основан на вычислении прямой обратной матрицы A^T.A дает приблизительно 0.05, тогда как решатель наименьших квадратов предоставляет значение 0.01.

Должен ли я изменить исходный массив y с помощью y.reshape (143,1), прежде чем продолжить?

Если вы не измените y, решение beta будет иметь форму (13,), а если вы сделаете y.reshape(143, 1), то форма beta будет(13, 1), тем не менее оба случая одинаково действительны ...

Answer 2

Насколько я понимаю, это проблема множественной регрессии, когда скорость самолета (v) и высота (h) используются для моделирования тяги (t), то есть "тяга = модель (скорость, высота)".Это означает, что некоторая математическая формула должна быть смоделирована с использованием заданных данных и предоставленных математических функций, таких как v в квадрате, высота в кубах и т. Д. Самый простой способ - использовать «линейную регрессию», где различные комбинации данных функций используются какв «t = a + b h + c v» и т. д. Заданное значение «1» будет означать «a * 1», которое является просто «a», это обычно называют смещениемпараметр.

Я сделал трехмерную диаграмму рассеяния данных, и у нее, кажется, есть некоторый резкий разрыв, вот изображение, показывающее это:

однако с помощью этого трехмерного анимированного GIF-файла (6,6 МБ) данных, вращающихся в трехмерном пространстве, легче увидеть: http://zunzun.com/temp/static_images/rotation.gif

Чтобы помочь вам в создании трехмерной диаграммы рассеяния, трехмерного графика поверхностифункции, и здесь трехмерный контурный график - это некоторый код Python для не -линейного подбора, который создает эти графики с помощью matplotlib.

import numpy, scipy, scipy.optimize
import matplotlib
from mpl_toolkits.mplot3d import  Axes3D
from matplotlib import cm # to colormap 3D surfaces from blue to red
import matplotlib.pyplot as plt

graphWidth = 800 # units are pixels
graphHeight = 600 # units are pixels

# 3D contour plot lines
numberOfContourLines = 16


def SurfacePlot(func, data, fittedParameters):
    f = plt.figure(figsize=(graphWidth/100.0, graphHeight/100.0), dpi=100)

    matplotlib.pyplot.grid(True)
    axes = Axes3D(f)

    x_data = data[0]
    y_data = data[1]
    z_data = data[2]

    xModel = numpy.linspace(min(x_data), max(x_data), 20)
    yModel = numpy.linspace(min(y_data), max(y_data), 20)
    X, Y = numpy.meshgrid(xModel, yModel)

    Z = func(numpy.array([X, Y]), *fittedParameters)

    axes.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=cm.coolwarm, linewidth=1, antialiased=True)

    axes.scatter(x_data, y_data, z_data) # show data along with plotted surface

    axes.set_title('Surface Plot (click-drag with mouse)') # add a title for surface plot
    axes.set_xlabel('X Data') # X axis data label
    axes.set_ylabel('Y Data') # Y axis data label
    axes.set_zlabel('Z Data') # Z axis data label

    plt.show()
    plt.close('all') # clean up after using pyplot or else thaere can be memory and process problems


def ContourPlot(func, data, fittedParameters):
    f = plt.figure(figsize=(graphWidth/100.0, graphHeight/100.0), dpi=100)
    axes = f.add_subplot(111)

    x_data = data[0]
    y_data = data[1]
    z_data = data[2]

    xModel = numpy.linspace(min(x_data), max(x_data), 20)
    yModel = numpy.linspace(min(y_data), max(y_data), 20)
    X, Y = numpy.meshgrid(xModel, yModel)

    Z = func(numpy.array([X, Y]), *fittedParameters)

    axes.plot(x_data, y_data, 'o')

    axes.set_title('Contour Plot') # add a title for contour plot
    axes.set_xlabel('X Data') # X axis data label
    axes.set_ylabel('Y Data') # Y axis data label

    CS = matplotlib.pyplot.contour(X, Y, Z, numberOfContourLines, colors='k')
    matplotlib.pyplot.clabel(CS, inline=1, fontsize=10) # labels for contours

    plt.show()
    plt.close('all') # clean up after using pyplot or else thaere can be memory and process problems


def ScatterPlot(data):
    f = plt.figure(figsize=(graphWidth/100.0, graphHeight/100.0), dpi=100)

    matplotlib.pyplot.grid(True)
    axes = Axes3D(f)
    x_data = data[0]
    y_data = data[1]
    z_data = data[2]

    axes.scatter(x_data, y_data, z_data)

    axes.set_title('Scatter Plot (click-drag with mouse)')
    axes.set_xlabel('X Data')
    axes.set_ylabel('Y Data')
    axes.set_zlabel('Z Data')

    plt.show()
    plt.close('all') # clean up after using pyplot or else thaere can be memory and process problems


def func(data, a, alpha, beta):
    t = data[0]
    p_p = data[1]
    return a * (t**alpha) * (p_p**beta)


if __name__ == "__main__":
    xData = numpy.array([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0])
    yData = numpy.array([11.0, 12.1, 13.0, 14.1, 15.0, 16.1, 17.0, 18.1, 90.0])
    zData = numpy.array([1.1, 2.2, 3.3, 4.4, 5.5, 6.6, 7.7, 8.0, 9.9])

    data = [xData, yData, zData]

    # this example uses curve_fit()'s default initial paramter values
    fittedParameters, pcov = scipy.optimize.curve_fit(func, [xData, yData], zData)

    ScatterPlot(data)
    SurfacePlot(func, data, fittedParameters)
    ContourPlot(func, data, fittedParameters)

    print('fitted prameters', fittedParameters)