Одномерная линейная регрессия в R —  Различные подходы

Цель этой статьи — продемонстрировать различные подходы, которые мы можем использовать для построения одномерной линейной регрессии в R. Мы увидим формулы, используемые для каждого подхода, и сравним их предсказания.

Введение

Линейная регрессия — это метод прогнозирования непрерывной переменной путем нахождения линейной связи между независимыми переменными и зависимыми переменными.

Линейное уравнение:

Линейная регрессия с одной независимой переменной называется одномерной линейной регрессией. Для одномерной линейной регрессии m всегда равно 1 (т. е. m = 1)

Есть два подхода к прогнозированию в LR. Они есть

  1. Нормальное уравнение
  2. Градиентный спуск

Нормальное уравнение

Есть два способа добиться предсказания с помощью нормального уравнения.

  1. Матричная форма
  2. Ковариация/дисперсия

Матричная форма

Ковариация/дисперсия

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

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

Функция гипотезы:

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

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

Импорт библиотек и загрузка данных

Установка библиотек

# norm_pkgs <- c('moments','nortest')
# install.packages(norm_pkgs)
# install.packages('matlib')
library(MASS) # for matrix operations

Загрузка данных

data <- read.csv('input/ex1data1.txt',header=FALSE, col.names = c('population','profit'))
head(data)
##   population  profit
## 1     6.1101 17.5920
## 2     5.5277  9.1302
## 3     8.5186 13.6620
## 4     7.0032 11.8540
## 5     5.8598  6.8233
## 6     8.3829 11.8860
X = as.matrix(data$population)
y = as.matrix(data$profit)

Преобразование переменных X и Y в массивы numpy (матрицы/векторы), а также добавление столбца 1 к X (поскольку значение X0 равно 1 в w0. x0, т.е. предполагается, что коэффициент перехвата равен 1)

Нормальное уравнение — ковариация/дисперсия

getCoeffs <- function(X,y)
{
  covariance_xy = cov(X,y)
  variance_x = var(X)
  coeffs = covariance_xy/variance_x
  y_intercept = mean(y) - coeffs * mean(X)
  return(list('intercept' = y_intercept[1], 'coefficient'=coeffs[1]))
  
}
getPred <- function(X,y)
{
  res <- getCoeffs(X,y)
  print(res)
  preds = res$intercept +(res$coefficient * X)
  return(preds)
}
preds = getPred(X,y)
## $intercept
## [1] -3.895781
## 
## $coefficient
## [1] 1.193034
## Model Evaluation
SSE = sum((preds-y)^2) # Sum of squared error
SST = sum((y-mean(y))^2) # sum of squared total
n = nrow(X) # number of observations
q = 1 # number of coefficients
MSE =  SSE/(n-q) # Mean squared error
MST = SST/(n-1) # Mean squared total
MAPE = sum(abs(preds-y)^2)/n # Mean Absolute Percentage Error
R_squared = 1-(SSE/SST) # R Square
A_rsquare = 1-(MSE/MST) # Adjusted R Square
std_err = MSE^1/2 # Standard Error/Root Mean Squared error
MSR = (SST-SSE)/(q-1)
f_static = MSR/MSE # F static
# paste('Mean Squared Percentage Error :',MSE)
# paste('Mean Absolute Percentage Error :',MAPE)
# print(paste('R Square :',R_squared))
# print(paste('Adjusted R Square :',A_rsquare))
# print(paste('Standard Error :',std_err))
print(paste('F Statistic :',f_static))
## [1] "F Statistic : Inf"

Среднеквадратическая ошибка в процентах: 9,047213
Средняя абсолютная ошибка в процентах: 8,9539428
Квадрат R: 0,7020316
Скорректированный квадрат R: 0,7020316
Стандартная ошибка: 4,5236065
F Статистика :

См. также:  Как использовать ИИ для улучшения управления качеством продукта и превзойти ожидания ваших клиентов

lm() в пакете статистики

lrModel <- lm(profit~population, data=data)
print(lrModel)
## 
## Call:
## lm(formula = profit ~ population, data = data)
## 
## Coefficients:
## (Intercept)   population  
##      -3.896        1.193
summary(lrModel)
## 
## Call:
## lm(formula = profit ~ population, data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.8540 -1.9686 -0.5407  1.5360 14.1982 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -3.89578    0.71948  -5.415 4.61e-07 ***
## population   1.19303    0.07974  14.961  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.024 on 95 degrees of freedom
## Multiple R-squared:  0.702,  Adjusted R-squared:  0.6989 
## F-statistic: 223.8 on 1 and 95 DF,  p-value: < 2.2e-16

Наблюдение:

P-значение модели является значительным.

Обе оценки, R-квадрат и скорректированный R-квадрат, выше и лучше Независимая переменная «популяция» значима

modsummary <- summary(lrModel)
modcoeffs <- modsummary$coefficients
std.err <- modcoeffs['population','Std. Error']
beta.estim <- modcoeffs['population','Estimate']
t_valve <- modcoeffs['population','t value']
# print(paste('t-statistic :',t_valve))
# print(paste('Std. Error :',std.err))
# print(paste('Beta Estimate/Coefficient :',beta.estim))

t-статистика: 14,9608056
стандарт. Ошибка: 0,0797439
Бета-оценка/коэффициент: 1,1930336

Наблюдение:

t-статистика хорошая (должно быть больше 1,96, чтобы p-значение было меньше 0,05), а стандартная ошибка меньше (близка к 0)

print(paste('AIC :',AIC(lrModel)))
## [1] "AIC : 493.907190165596"
print(paste('BIC :',BIC(lrModel)))
## [1] "BIC : 501.631323101106"

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

X = cbind(1,matrix(data$population))
y = as.matrix(data$profit)

gradDescent <- function(x,y,alpha=0.005,maxIter=100,threshold=0.0001)
{
  thetas <- list(maxIter)
  costs <- double(maxIter)
  converged <- FALSE
  theta <- matrix(c(0,0),nrow = 2)
  m <- nrow(y)
  i <- 1
  while(converged==FALSE)
  {
    preds <- x %*% theta
    theta <- theta - alpha * (t(x) %*% (preds-y))/m
    cost = mean((preds-y)^2)/2
    
    if (i > 1)
    {
      if ((costs[i-1]-cost) == threshold)
      {   converged <- TRUE      }
    }
    if (i == maxIter)
    {      converged <- TRUE    }
    costs[i] <- cost
    thetas[[i]] <- theta
    i <- i+1
  }
  print(paste('theta :',thetas[which.min(costs)]))
  print(paste('Cost :',min(costs)))
  print(paste('Iteration :',which.min(costs)))
  #return(thetas[which.min(costs)])
}
gradDescent(X,y,maxIter=15000)
## [1] "theta : c(-3.89577556561726, 1.1930331104718)"
## [1] "Cost : 4.47697137597775"
## [1] "Iteration : 15000"

Вывод

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

См. также:  Расширенные хранимые процедуры в MySQL

Ссылка:

https://sebastianraschka.com/faq/docs/closed-form-vs-gd.html

Понравилась статья? Поделиться с друзьями:
IT Шеф
Добавить комментарий

;-) :| :x :twisted: :smile: :shock: :sad: :roll: :razz: :oops: :o :mrgreen: :lol: :idea: :grin: :evil: :cry: :cool: :arrow: :???: :?: :!: