[수치해석] 최소 제곱 방법(Least square approximation)

회귀 분석에서 많이 쓰이는 최소 제곱 오차에 대한 수학적 설명입니다. 

 

회귀 분석은 분포된 데이터 값을 한줄로 표현할 수 있는 적절한 함수를 찾아내는 것인데 

최소 제곱 오차를 쓰면 괜찮은 함수를 찾아낼 수 있습니다. 

가장 직관적인 접근을 적어놓겠습니다.

 

Linear case

어떤 곡선을 넣을 것인가에 따라 방향이 달라질 수 있는데 가장 쉬운 방법인 최적의 직선을 만들어 보겠습니다. 

즉, 직선은 $ y= a+bx $ 로 표현되는데 적절한 a와 b를 잘 찾아내어 최적의 함수를 구해보는 것입니다.

최적의 함수를 찾아내기 위해 최소 제곱 방법을 사용합니다. 

 

최적의 직선을 $ y = a+bx $ 라고 합시다.

데이터 $ (x_i,f_i), \ \ i = 0,...,n $ 에 대해서, 각 데이터와 함수간의 에러를 $ e_i $ 라 하면 에러함수는 다음과 같은 식이 됩니다.

$$ e_i = f_i - (a +bx_i) $$ 

모든 데이터에 대해서 최소의 오차를 내야하는데 1) 각각의 오차를 최소를 할것인지 2) 전체에러의 최소화를 할지에 따라 진행방향이 달라질겁니다.

그런데 각 오차를 최소화한다는 건 사실상 불가능합니다.

왜냐하면 하나의 오차를 최소화하기 위해 직선을 옮기면 다른 쪽의 오차는 분명히 커지기 때문에 결국 모두 최소화하기 위해서 하나씩 건드리면 균형이 무너져 어려운 작업이 되기 때문입니다.

그래서 전체 에러의 합이 가장 작은 경우를 최적화 되었다고 여기고 진행하게 됩니다.

 

목표는 에러의 합인 $ S = \sum_{i=0}^{n} e_i $ 의 최소화할 수 있는 a,b를 구하는 것입니다.

하지만, $ e_i $는 + 와 - 모두가 가능합니다. 위 식대로 S를 구하면 오차가 상쇄되어 문제가 생깁니다. 

모든 $ e_i $ 가 양수가 되게 하기 위해 제곱을 합니다.

즉, S를 구하는 문제를 다시 쓰면 다음과 같습니다. 

$$ S = \min_{a,b} {\left[E(a,b) = \sum_{i=0}^{n} (f_i-a-bx_i)^2 \right]} \tag{1} $$

 

a,b 를 구해야하니 E(a,b) 로써 표현하였습니다.

 

(1) 식은 두 변수 a,b의 2차식으로 되어있고 2차식의 최소값은 극소값 중에 가장 작은 값으로 구할 수 있습니다.

따라서, $$ \frac{ \partial E}{ \partial a} = \frac{ \partial E}{ \partial b} = 0 $$ 을 만족해야 합니다.

 

실제로 (1)식을 미분하면 다음과 같이 나옵니다.

$$ \begin{align} 0=\frac{ \partial E}{ \partial a} & = -2 \sum_{i=0}^{n} (f_i-a-bx_i) \\ 0= \frac{ \partial E}{ \partial b} & = -2 \sum_{i=0}^{n} x_i(f_i-a-bx_i) \end{align} $$

 

정리하면 다음과 같습니다. 

$$ \begin{align} (n+1)a + \sum_{i=0}^{n} x_ib & = \sum_{i=0}^{n} f_i \tag{2} \\ \sum_{i=0}^{n} x_ia + \sum_{i=0}^{n} x_i^2b & = \sum_{i=0}^{n} x_if_i \tag{3} \end{align} $$

 

식 (2)에서 (n+1)으로 나눈다면 시그마로 이루어진 식은 평균이 되어버립니다. 따라서 평균 표시를 $ \overline { } $ 이라 하면, $$ a = \overline{f} - \overline{x}b $$ 가 됩니다.

구한 a 를 식 (3) 에 적용하면 b를 구할 수 있습니다.

$$ b = \frac{ \overline{x} \overline{f} - \overline{xf}} { \overline{x}^2 - \overline{x^2}} $$

 

이런식으로 a,b를 구해 데이터의 적합한 오차가 최소화된 직선을 구하게 됩니다.

 

파이썬으로 구현

a,b 구하는 방법을 알았으니 파이썬으로 한번 구현해보겠습니다.

import numpy as np
import matplotlib.pyplot as plt

# 데이터
x = np.array([0.0, 0.1, 0.3, 0.4, 0.5])
f = np.array([0.5, 1.4, 2.0, 2.5, 3.1])

# 최소제곱으로 a,b 구하기
def least_square_line(x,f):
    b = (np.mean(x)*np.mean(f) - np.mean(x*f))/(np.mean(x)*np.mean(x)-np.mean(x*x))
    a = np.mean(f)-np.mean(x)*b
    return a,b
    
a,b = least_square_line(x,f)
print(f'y = {a:.2f}+{b:.2f}x')

 

그림으로 보면 다음과 같습니다.

plt.plot(x, a+b*x)
plt.scatter(x,f)

 

이런 식으로 최소제곱 방법으로 데이터의 적합한 직선을 만들 수 있습니다.

 

Non-linear case

이번엔 직선으로 표현이 안 되는 경우입니다. 

곡선이니 보통 증가율 모델로 접근해야하므로 지수함수 모양으로 나타나야 합니다.

따라서 최적의 함수는 $ y = ae^{bx} $ 가 됩니다.

이걸 로그를 이용하고 각 로그를 대문자로 치환하면 다음과 같이 됩니다.

$$ Y = A+bx \ \ , A = \ln{a} , Y = \ln{y} $$

그럼 구하는 방식이 Linear case 와 같아지고 데이터셋은 $ \{(x_i,Y_i)\} $ 로 바뀌게 됩니다.

Linear case 같은 방식으로 A와 b를 구하면 원하는 함수를 구하게 됩니다.

$ a= e^{A} $ 로 놓아서 원래 식으로 출력을 하면 되겠습니다.

 

파이썬으로 구현

x = np.array([0,1,2,4,6,7,9])
t = np.array([100, 56.8, 38.6, 13.5, 4.78, 2.9, 1.16])
Y = np.log(t)

def non_linear_least_square(x,Y):
    b = (np.mean(x)*np.mean(Y) - np.mean(x*Y))/(np.mean(x)*np.mean(x)-np.mean(x*x))
    A = np.mean(Y)-np.mean(x)*b
    a = np.exp(A)
    return a,b
    
a,b = non_linear_least_square(x,Y)
plt.plot(x, a*np.exp(b*x))
plt.scatter(x,t)

관련 포스팅

[데이터 사이언스/머신러닝 딮러닝] - [Python] 회귀(Regression)

 

 

참고문헌 

Introduction to numerical analysis/alastair wood