Least Squares Problem
Kong, Seokkyu
2018 10 13
Summary
Least Squares Problem 에 대한 기초 지식을 복습한다.
참고: MIT 선형대수학 과정 Lecture 15, PRML(Pattern Recognition and machine learning) 1장
Gilbert Strang 교수님은 15 강의 내용중 (1, 1), (2, 2), (3, 2) 에 대한 최소 제곱문제를 설명하고 있다.
여기서 주목할 것은 어떤 데이터를 입력으로 보고 어떤 데이터를 출력으로 볼것이냐이다. 이와 같은 형태의 데이터들은 주변에서 많이 볼 수 있다. 강의에서는 앞에 있는 1, 2, 3을 입력값으로 보고 뒤의 1, 2, 2 를 출력값으로 본다.
선형회귀에서 가설에 대한 모델(hypothesis)은 y = a0 + a1 * x로 주어진다.
위의 데이터를 방정식으로 다시 쓰면,
1 = a0 + a1 * 1 2 = a0 + a1 * 2 2 = a0 + a1 * 3
이것을 행렬형태로 쓰면
library(ggplot2)
x <- c(1, 2, 3)
y <- c(1, 2, 2)
A <- matrix(c(1, 1, 1, x), nrow=3)
A## [,1] [,2]
## [1,] 1 1
## [2,] 1 2
## [3,] 1 3우리가 찾는 것은 a0, a1 이다. 하지만, Ax = b 형태에서 이 방정식을 만족하는 해 x = (a0, a1)는 없다. 그리고 b도 A의 컬럼스페이스 공간에 있지 못해서 projection을 해야 한다.
이럴때 최소제곱해를 찾아보자.
# At A x_hat = At b => x_hat = (At A)-1 At b
At <- t(A)
x_hat <- solve(At %*% A) %*% At %*% y
intercept <- x_hat[1]
slope <- x_hat[2]
ggplot(data = data.frame(x, y), aes(x = x, y = y)) + geom_point() +
scale_x_continuous(limits = c(0, 5)) + scale_y_continuous(limits = c(0, 5)) +
geom_abline(intercept = intercept, slope = slope, color = 'red')
까만 점은 측정값이고 빨간선은 선형회귀선이 된다.
그런데, 만약 측정된 데이터가 직선형태가 아닌 sin 주기 데이터를 가지면 어떻게 되나? 노이즈를 가진 데이터를 생성해보자.
x <- seq(from = 0, to = 1, by = 0.1) # 11개의 입력값을 생성한다.
y <- sin(2*pi*x) + rnorm(11, 0, 0.02) # 11개의 출력값을 생성한다. 노이즈를 약간 섞는다.
y## [1] -0.0005595046 0.6004035159 0.9551742214 0.9448357758 0.6045452376
## [6] 0.0158118639 -0.6141137936 -0.9951025018 -0.9444906971 -0.5986758223
## [11] 0.0268994705plot(x, y)
이 데이터를 y = a0 + a1x 라는 직선으로 추정을 하면
A <- matrix(c(rep(1, 11), x), nrow=11)
A## [,1] [,2]
## [1,] 1 0.0
## [2,] 1 0.1
## [3,] 1 0.2
## [4,] 1 0.3
## [5,] 1 0.4
## [6,] 1 0.5
## [7,] 1 0.6
## [8,] 1 0.7
## [9,] 1 0.8
## [10,] 1 0.9
## [11,] 1 1.0At <- t(A)
x_hat <- solve(At %*% A) %*% At %*% y # 계수를 구한다.
b_hat <- A %*% x_hat # 예측치를 구한다. Normal Equation
intercept <- x_hat[1]
slope <- x_hat[2]
# 추정된 기울기와 절편을 확인한다.
intercept## [1] 0.7020913slope## [1] -1.405141ggplot(data = data.frame(x, y), aes(x = x, y = y)) + geom_point() +
scale_x_continuous(limits = c(0, 2)) + scale_y_continuous(limits = c(-2, 2)) +
geom_abline(intercept = intercept, slope = slope, color = 'red')
# rmse 값을 구한다.
rmse <- sqrt(mean((b_hat - y)^2))
rmse## [1] 0.5203128그래프와 같이 전혀 맞지 않는다. 그래서 모델을 고차다항식으로 재설정한다.
y = a0 + a1 * x + a2 * x^2 + a3 * x^3
가설의 모델의 위와 같다면, 이때 우리가 찾는 것은 계수 a0, a1, a2, a3가 된다. 하지만, 이 역시도 Ax = b 에 대한 정확한 솔루션은 없다. x = [a0 a1 a2 a3], b = [y1 ~ y11]
그래서 솔루션(계수)을 찾기 위해 Least Squares를 이용할 수 있고, 계수를 구하면 예측치도 구할 수 있다.
A <- matrix(c(rep(1, 11), x, x^2, x^3), nrow=11)
A## [,1] [,2] [,3] [,4]
## [1,] 1 0.0 0.00 0.000
## [2,] 1 0.1 0.01 0.001
## [3,] 1 0.2 0.04 0.008
## [4,] 1 0.3 0.09 0.027
## [5,] 1 0.4 0.16 0.064
## [6,] 1 0.5 0.25 0.125
## [7,] 1 0.6 0.36 0.216
## [8,] 1 0.7 0.49 0.343
## [9,] 1 0.8 0.64 0.512
## [10,] 1 0.9 0.81 0.729
## [11,] 1 1.0 1.00 1.000At <- t(A)
x_hat <- solve(At %*% A) %*% At %*% y # 계수를 구한다.
b_hat <- A %*% x_hat # 예측치를 구한다.
# 추정된 계수항을 확인한다.
x_hat## [,1]
## [1,] -0.06664807
## [2,] 10.90994054
## [3,] -32.41285780
## [4,] 21.65708631# 측정데이터와 예측데이터를 가시화한다.
data_set <- data.frame(x = x, y = y, type='real')
predict_set <- data.frame(x = x, y = b_hat, type='predict')
mydata <- rbind(data_set, predict_set)
ggplot(data = mydata, aes(x = factor(x), y = y, group = type, colour = type)) + geom_point() 
RMSE (Root mean square error) 를 구해보면…
rmse <- sqrt(mean((b_hat - y)^2))
rmse## [1] 0.08530647자, 이제 우리는 짐작할 수 있다. 데이터가 (입력값, 출력값) 과 같은 형태로 주어질때 Least Squares 문제로 선형회귀식을 구할 수 있었다. 하지만, 데이터의 패턴에 따라, 1차 직선으로 fitting이 되기도 하고, 주기성을 갖는 데이터라면 고차 다항식으로 모델을 설정해야 했다. 그럼, 3차 이상의 고차 다항식으로 모델을 설정했다면, 행렬 A를 구성하는 방법을 알 수 있겠지?
그런데, 데이터가 (입력값1, 입력값2, ,,, 출력값)과 같은 형태라면, 고차 다항식을 어떻게 설정해야 할까? 음, 위와 비슷한 방식을 취한다. 각 feature 별 고차 다항식 계수를 만들고 각 feature를 서로 묶어서 고차 다항식 계수를 만든다. 이럴 경우 계수 항목이 무척 많아지겠지?
이런 문제들로 여러 가설에 대한 모델식이 있다면 어떤 것이 효율적인지 구하고, 각 모델에 있어서도 계수의 중요도에 따라 계수항이 있기도 하고 없어지기도 한다.
모델이 bias라면, 좀더 많은 고차 다항식의 모델을 고려할 수 있고, variance라면, 모델의 항을 줄이는 방안을 고려해야 한다.