1

У меня есть продажи товаров с января 2013 года по октябрь 2015 года. Я просто хочу спрогнозировать общий объем продаж на следующий месяц. Просто для того, чтобы научиться, я хотел бы превратить его в многократную регрессионную модель, закодированную с нуля, без каких-либо библиотек. Пока что мне удалось получить бета-версии, но я не знаю, как получить прогноз на следующий месяц.

Вот исторические данные ежемесячных продаж с января 2013 года по октябрь 2015 года,

Ниже приведены исторические данные ежемесячных продаж с января 2013 года по октябрь 2015 года, ts:

date_block_num
0     131479.0
1     128090.0
2     147142.0
3     107190.0
4     106970.0
5     125381.0
6     116966.0
7     125291.0
8     133332.0
9     127541.0
10    130009.0
11    183342.0
12    116899.0
13    109687.0
14    115297.0
15     96556.0
16     97790.0
17     97429.0
18     91280.0
19    102721.0
20     99208.0
21    107422.0
22    117845.0
23    168755.0
24    110971.0
25     84198.0
26     82014.0
27     77827.0
28     72295.0
29     64114.0
30     63187.0
31     66079.0
32     72843.0
33     71056.0

Я пытался сделать простую линейную регрессию:

$$y_t = \alpha + \beta x_t +\varepsilon$$

Сначала я попытался оценить alpha и beta, а затем использовал predict(alpha,beta,34). Так я и сделал:

import random

def predict(alpha, beta, x_i):
  return alpha+ beta * x_i

def error(alpha, beta, x_i, y_i):
  """the error from predicting beta * x_i + alpha
  when the actual value is y_i"""
  return y_i - predict(alpha, beta, x_i)

def sum_of_squarred_errors(alpha, beta, x, y):
  return sum(errors(alpha, beta, x_i, y_i)**2
             for x_i, y_i in zip(x,y))
  
def correlation(x,y):
  stdev_x = standard_deviation(x)
  stdev_y = standard_deviation(y)
  if stdev_x > 0 and stdev_y >0:
    return covariance(x,y)/ stdev_x/ stdev_y
  else:
    return 0

def squared_error(x_i, y_i, theta):
  alpha, beta = theta
  return error(alpha, beta, x_i, y_i) ** 2

def squared_error_gradient(x_i, y_i, theta):
  alpha, beta = theta
  return [-2 * error(alpha, beta, x_i, y_i),
          -2 * error(alpha, beta, x_i, y_i) * x_i]

def in_random_order(data):
  """generator that returns the elements if data in random order"""
  indexes = [i for i, _ in enumerate(data)] # create a list of indexes
  random.shuffle(indexes) # suffle them
  for i in indexes:
    yield data[i]

def scalar_multiply(c, v):
    """c is a number, v is a vector"""
    return [c*v_i for v_i in v]

def minimize_stochastic(target_fn, gradient_fn, x,y, theta_0, alpha_0=0.01):
  print("x: ", x, "\ny: ",y.tolist())
  data = zip(x,y)
  theta = theta_0  #initial guess
  alpha = alpha_0  # initial step size
  min_theta, min_value = None, float('inf') # the minimum so far
  iterations_with_no_improvment = 0

  # if we ever go 100 iterations with no improvment, stop
  while iterations_with_no_improvment < 100:
    value = sum(target_fn(x_i, y_i, theta) for x_i, y_i in data)
    # print("value: ", value)

    if value < min_value:
      # if we've found a new minimum, remember it
      # and go back to the original step size
      min_theta, min_value = theta, value
      iterations_with_no_improvment = 0
      alpha = alpha_0
    else:
      # otherwise we're not improving, so try shrinking the step size
      iterations_with_no_improvment +=1
      alpha *=0.9

    # and take a gradient step for each of the data points
    for x_i, y_i in in_random_order(data):
      gradient_i = gradient_fn(x_i, y_i, theta)
      theta = vector_substract(theta, scalar_multiply(alpha, gradient_i))
  return min_theta

# choose random value to start
random.seed(0)
theta = [random.random(), random.random()]

alpha, beta = minimize_stochastic(squared_error,
                                  squared_error_gradient, ts.index.values,
                                  ts.values,
                                  theta,
                                  0.001)

print("alpha: ", alpha, "beta: ", beta)

Но у него супернизкие альфа и бета:

alpha:  0.8444218515250481 beta:  0.7579544029403025

Таким образом, общий объем продаж за 34 (ноябрь 2015 года): 26,614871551495334 что выглядит невозможным по сравнению с 33 (октябрь 2015 г.): 71056,0

Так я испортил алгоритм линейной регрессии? Я думаю, что мои случайные значения для начала, возможно, слишком низкие:

theta = [random.random(), random.random()]

Тем не менее, они должны увеличиваться в любом случае до тех пор, пока больше нет входа, не так ли?

Так как же выбрать исходные теории для простой линейной регрессии?

5
  • А не слишком ли широко вы шагаете по градиенту? -2 не многовато? Не факт, что здесь проблема, просто бросилось в глаза.
    – CrazyElf
    27 окт 2020 в 13:17
  • А вообще для time series простая регрессия точно не взлетит. У вас явный спад в январе, нужно выделять месяц года как фичу, это как минимум. Просто из сквозного номера месяца модель ничего не поймёт. Ну и третий год у вас весь год спад, значит год году рознь, надо и год выделять как фичу отдельно.
    – CrazyElf
    27 окт 2020 в 13:20
  • NameError: name 'in_random_order' is not defined
    – CrazyElf
    27 окт 2020 в 13:27
  • NameError: name 'scalar_multiply' is not defined
    – CrazyElf
    27 окт 2020 в 13:29
  • @CrazyElf Да, это так. Прости меня за то, что я забыл, дорогой друг. Я только что узнал о том, чего не хватает. 27 окт 2020 в 15:53

1 ответ 1

1

Давайте по порядку.

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

Во-вторых, если вы хотите использовать именно одномерную линейную регрессию, и хотите ее написать вручную, без привлечения сторонних модулей, то надо использовать "классику жанра" и искать уравнения в явном виде

numerator=np.sum((x-np.mean(x))*(y-np.mean(y)))
denominator=np.sum((x-np.mean(x))**2)
b1=numerator/denominator
b0=np.mean(y)-b1*np.mean(x)
plt.plot(x,y,c='b')
x=np.arange(1,40)
y_pred=b0+b1*x
plt.plot(x,y_pred,c='r')

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

Идея запускать градиентный поиск для данной задачи - немного сомнительна. Если вам поможет, то значения коэффициентов b0 и b1 для вашего примера равны 141311.93048128343 и -1943.5288006111537 соответственно.

P.S. Что такое "многократная регрессионная модель" я, увы, не знаю. Расшифруйте пожалуйста сей термин.

4
  • Какая-нибудь авторегрессия может имелась в виду? )
    – CrazyElf
    27 окт 2020 в 14:01
  • @CrazyElf не знаю, сначала я подумал про экспоненциальное сглаживание, но формула явно не та. А так, конечно авторегрессионная модель должна быть, SARIMA, если повезет - то с нулевыми i и m параметрами.
    – passant
    27 окт 2020 в 14:05
  • @passant Большое спасибо за ваше предложение. Я считаю, что использование «классики жанра» и поиск уравнений в явном виде - лучшее решение для одномерной линейной регрессии. Как вы думаете, можно ли отремонтировать мою модель в любом случае или она действительно испорчена? 28 окт 2020 в 18:23
  • Что значит "испорчена". Я ваш код честно говоря не разбирал подробно, но фрагменты полностью (с точностью до имени переменной) совпадающие с моим примером нашел. Вот и отбросьте все лишнее, что "портит". А хотите что-то добавить "от себя" - добавляйте, но потом, когда убедитесь, что "ядро" работает правильно.
    – passant
    28 окт 2020 в 20:36

Ваш ответ

By clicking “Отправить ответ”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Всё ещё ищете ответ? Посмотрите другие вопросы с метками или задайте свой вопрос.