[Из песочницы] Как программист машину покупал

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

Как известно, для покупки авто на территории РФ существует несколько крупных авторитетных сайтов (auto.ru, drom.ru, avito.ru), поиску на которых я и отдал предпочтение. Моим требованиям отвечали сотни, а для некоторых моделей и тысячи, автомобилей, с перечисленных выше сайтов. Помимо того, что искать на нескольких ресурсах неудобно, так еще, прежде чем ехать смотреть авто «вживую», я хотел бы отобрать выгодные (цена которых относительно рынка занижена) предложения по априорной информации которую предоставляет каждый из ресурсов. Я, конечно, очень хотел решить несколько переопределенных систем алгебраических уравнений (возможно и нелинейных) высокой размерности вручную, но пересилил себя, и решил этот процесс автоматизировать.

Сбор данных


Данные я собирал со всех описанных выше ресурсов, а интересовали меня следующие параметры:

  • цена (price)
  • год выпуска (year)
  • пробег (mileage)
  • объем двигателя (engine.capacity)
  • мощность двигателя (engine.power)
  • тип двигателя (2 индикаторные взаимоисключающие переменные diesel и hybrid, принимающие значения 0 или 1, для дизельных и гибридных двигателей соответственно). Тип двигателя по умолчанию — бензиновый (не вынесен в третью переменную во избежании мультиколлинеарности).

    Таким образом:

    image

    Далее подобная логика для индикаторных переменных подразумевается по умолчанию.

  • тип кпп (индикаторная переменная mt (manual transmission), принимающая булево значение, для механической коробки передач). Тип кпп по умолчанию — автоматическая.

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

  • тип привода (2 индикаторные переменные front.drive и rear.drive, принимающие булевы значения). Тип привода по умолчанию — полный.
  • тип кузова (7 индикаторных переменных sedan, hatchback, wagon, coupe, cabriolet, minivan, pickup, принимающих булевы значения). Тип кузова по умолчанию — внедорожник/кроссовер


Незаполненные ленивыми продавцами данные я отмечал как NA (Not Available), чтобы можно было корректно обрабатывать эти значения с помощью R.

Визуализация полученных данных


Дабы не вдаваться в сухую теорию, давайте рассмотрим конкретный пример, будем искать выгодные Mercedes-Benz E-klasse не старше 2010 года выпуска, стоимостью до 1.5 млн. рублей в Москве. Для того, чтобы начать работать с данными, первым делом заполняем пропущенные значения (NA) на медианные, благо для этого в R есть функция median ().

dat <- read.csv("dataset.csv") # загружаем выборку в R

dat$mileage[is.na(dat$mileage)] <- median(na.omit(dat$mileage)) # например для пробега


Для остальных переменных процедура идентична, поэтому опущу этот момент.

Теперь посмотрим как цена зависит от регрессоров (визуализация индикаторных переменных на данном этапе нас не интересует).

image

Оказывается за эти деньги есть несколько машин 2013 года и даже одна 2014!

image

Очевидно, что чем меньше пробег, тем цена выше.

image

image

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

Хочу обратить ваше внимание на тот факт, что в большинстве статей по машинному обучению, которые попадались мне в последнее время, в том числе на «Хабре», очень мало внимания уделяется обоснованию правомерности использования выбранной модели, ее диагностики и ошибках.

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

Диагностика модели


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

Чтобы модель была корректной, необходимо выполнение условий теоремы Гаусса-Маркова:

  1. Модель данных правильно специфицирована, т.е.:
    • отсутствуют пропущенные значения.
      Подробнее

      Это условие выполнено (см. раздел Визуализация полученных данных), пропущенные значения заменены на медианные.


    • отсутствует мультиколлинеарность между регрессорами.
      Подробнее
      Проверим, выполняется ли это условие:
      dat_cor <- as.matrix(cor(dat)) # рассчет корреляции между переменными
      
      dat_cor[is.na(dat_cor)] <- 0 # заменяем пропущенные значения на 0 (т.к. например в кузове пикап или минивен искомого авто не бывает)
      
      library(corrplot) # подключаем библиотеку corrplot, для красивой визуализации
      
      palette <-colorRampPalette(c("#7F0000","red","#FF7F00","yellow","#7FFF7F", "cyan", "#007FFF", "blue","#00007F"))
      
      corrplot(dat_cor, method="color", col=palette(20), cl.length=21,order = "AOE", addCoef.col="green") # рисуем таблицу зависимостей между переменными
      
      

      image

      Видно что присутствует большая корреляция (>0.7) между объемом и мощностью двигателя, поэтому при дальнейшем анализе мы не будем учитывать переменную engine.capacity, т.к. именно мощность двигателя позволит более точно построить регрессионную модель по сравнению с объемом двигателя (атмосферный бензиновый, бензиновый с турбонаддувом, дизельный моторы — при одном и том же объеме могут имеют разную мощность).


    • отсутствуют выбросы.
      Подробнее
      Выбросы — показатели выделяющиеся из общей выборки (см. раздел Визуализация полученных данных), оказывают существенное воздействие на оценки коэффициентов регрессионной модели. Статистический метод, способный действовать в условиях выбросов, называется робастным — линейная регрессия к ним не относится, в отличие, например, от робастной регрессии Хубера или метода урезанных квадратов.

      Меры влияния выбросов на оценки модели можно подразделить на общие и специфические. Общие меры, такие как расстояние Кука, dffits, ковариацинное отношение (covratio), расстояние Махаланобиса, показывают как i-тое наблюдение влияет на положение всей регрессионной зависимости, их мы и будем использовать для идентификации выбросов. Специфические меры влияния, такие как dfbetas, показывают влияние i-того наблюдения на отдельные параметры регрессионной модели.

      Ковариацинное отношение (covratio) — общая мера влияния наблюдения. Представляет собой отношение детерминанта ковариационной матрицы с удаленным i-ым наблюдением к детерминанту ковариацинной матрицы для всего набора данных.

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

      Расстояние Махаланобиса — мера удаления наблюдения от центра системы, но т.к. она не учитывает зависимый или независимый характер переменной и рассматривает их в облаке рассеяния равнозначно, данная мера не предназначена для регрессионного анализа.
      Таким образом, для обнаружения выбросов будем использовать меры dffits и covratio.

      model <- lm(price ~ year + mileage  + diesel + hybrid + mt + front.drive + rear.drive + engine.power + sedan + hatchback + wagon + coupe + cabriolet + minivan + pickup, data = dat) # линейная модель
      
      model.dffits <- dffits(model) # рассчитаем меру dffits
      
      

      Значимыми для меры dffits являются показатели превосходящие величину 2*sqrt (k/n) = 0.42, поэтому следует их отбросить (k — количество переменных, n — число строк выборки).
      model.dffits.we <- model.dffits[model.dffits < 0.42]
      
      model.covratio <- covratio(model)  # рассчитаем ковариационное отношения для модели
      
      

      Значимые для меры covratio показатели можно найти из неравенства | model.covratio[i] -1 | > (3*k)/n.
      model.covratio.we <- model.covratio[abs(model.covratio -1) < 0.13]
      
      dat.we <- dat[intersect(c(rownames(as.matrix(model.dffits.we))), c(rownames(as.matrix(model.covratio.we)))),] # наблюдения без выбросов
      
      model.we <- lm(price ~ year + mileage  + diesel + hybrid + mt + front.drive + rear.drive + engine.power + sedan + hatchback + wagon + coupe + cabriolet + minivan + pickup, data = dat.we) # линейная модель построенная по выборке без выбросов
      
      

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

      image

      image

      image

      Итого, выбранный нами путь идентификации выбросов позволил исключить из общий выборки 18 наблюдений, что несомненно положительно скажется на точности определения коэффициентов линейной модели с помощью МНК.



  2. Все регрессоры детерминированы и не равны.
    Подробнее

    Данное условие выполнено.


  3. Ошибки не носят систематического характера, дисперсия ошибок одинакова (гомоскедастичность)
    Подробнее
    Посмотрим как распределены ошибки модели (для расчета ошибок модели в R есть функция resid ()).
    plot(dat.we$year, resid(model.we))
    plot(dat.we$mileage, resid(model.we))
    plot(dat.we$engine.power, resid(model.we))
    
    

    image

    image

    image

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


  4. Ошибки распределены нормально.
    Подробнее
    Для проверки данного условия построим график квантилей остатков против квантилей, которые можно было бы ожидать при условии, что остатки нормально распределены.
    qqnorm(resid(model.we)) 
    qqline(resid(model.we))
    
    

    image


Опробирование выбранной модели


Итак, настал момент, когда с бюрократией покончено, корректность использования модели линейной регрессии не вызывает сомнений. Нам остается только посчитать коэффициенты линейного уравнения, по которым мы сможем получить свои, предсказанные цены на автомобили и сравнить их с ценами из объявлений.

model.we <- lm(price ~ year + mileage  + diesel + hybrid + mt + front.drive + rear.drive + engine.power + sedan + hatchback + wagon + coupe + cabriolet + minivan + pickup, data = dat.we)

coef(model.we) # коэффициенты линейной модели
  (Intercept)      year           mileage    diesel    rear.drive  engine.power         sedan 
  -1.76e+08   8.79e+04  -1.4e+00  2.5e+04  4.14e+04     2.11e+03        -2.866407e+04        
    
predicted.price <- predict(model.we, dat) # предскажем цену по полученным коэффициентам

real.price <- dat$price # вектор цен на автомобили полученный из объявлений

profit <- predicted.price - real.price # выгода между предсказанной нами ценой и ценой из объявлений


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

plot(real.price,profit)
abline(0,0)

image

И на какую выгоду в процентном соотношении можно рассчитывать?


sorted <- sort(predicted.price /real.price, decreasing = TRUE)
sorted[1:10]
      69       42      122      248      168       15      244      271      109      219 
1.590489 1.507614 1.386353 1.279716 1.279380 1.248001 1.227829 1.209341 1.209232 1.204062 


Да, экономия 59% это очень здорово, но весьма сомнительно, нужно смотреть авто «вживую», т.к. бесплатный сыр обычно в мышеловке, ну или продавцу срочно нужны деньги. А вот начиная с 4-го места (экономия 28%) и далее, результат кажется вполне реалистичным.

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

Напоследок


Дорогие друзья, я сам, прочитав бы подобную статью задал бы в первую очередь 3 вопроса:

  1. Какая у модели точность?
  2. Почему выбрана простейшая линейная регрессия и где сравнение с другим моделями?
  3. А почему бы не сделать сервис для поиска выгодных автомобилей?


Поэтому отвечаю:

  1. Протестируем модель по схеме 80/20.
    split <- runif(dim(dat.we)[1]) > 0.2 # разделяем нашу выборку
    
    train <- dat.we[split,] # выборка для обучения
    
    test <- dat.we[!split,] # тестовая выборка
    
    train.model <- lm(price ~ year + mileage  + diesel + hybrid + mt + front.drive + rear.drive + engine.power + sedan + hatchback + wagon + coupe + cabriolet + minivan + pickup, data = train) # линейная модель построенная по выборке для обучения
    
    predictions <- predict(train.model, test) # проверим точность на тестовой выборке
    
    print(sqrt(sum((as.vector(predictions - test$price))^2))/length(predictions)) # точность предсказания цены (в рублях)
    [1] 11866.34 
    
    

    Много это или мало — можно узнать только в сравнении.
  2. Линейная регрессия, как один из азов машинного обучения позволяет не отвлекаться на изучение алгоритма, а полностью погрузиться в саму задачу.

    Поэтому я решил разделить своё повествование на 2 (а там как пойдет) статьи.

    В следующей статье речь пойдет о сравнении моделей для решения нашей задачи с целью выявления победителя.

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

    И вот, что из этого получилось — Mercedes-Benz E-klasse не старше 2010 года выпуска, стоимостью до 1.5 млн. рублей в Москве.

    image

    Ещё
    image


Ссылки


1. Исходные данные в формате csv — dataset.txt

© Habrahabr.ru