긍정적인 개발자의 꿈

coursera에서 진행하는 logistic regression을 구현해 보도록 하겠습니다.

logistic regression은 classification의 알고리즘이며 또한 supervisor learning입니다.

즉, 주어진 데이터를 여러형태로 분류하기 위해 사용합니다.

이전에 배운 linear regression 사용하지 못하고 새로운 접근법을 이용하여 알고리즘을 적용해야 합니다.

이번에 참고할 데이터 셋은 아래와 같이 3개로 구분되어 있습니다.

각각 x좌표, y좌표, 결과 라고 가정 할 수 있습니다.

파일이름 : ex2data1.txt

34.62365962451697,78.0246928153624,0

30.28671076822607,43.89499752400101,0

35.84740876993872,72.90219802708364,0

60.18259938620976,86.30855209546826,1

79.0327360507101,75.3443764369103,1

45.08327747668339,56.3163717815305,0

61.10666453684766,96.51142588489624,1

....

이 데이터의 1번째, 2번째를 data x로

3번째를 data y로 나눠 보겠습니다.

>> data = load('ex2data1.txt');

>> [m, n] = size(data);

>> X = data(:, 1:n-1);

>> y = data(:, n);

>> [m, n] = size(X);

나눠진 데이터를 각각 y 데이터를 보고 1은 positive, 0은 negative로 하여 

두 데이터를 나누어 그래프에 표시해 보겠습니다.

function plotData(X, y)

% create new chart

figure; hold on;

% find positive, negative position

pos = find(y == 1);

neg = find(y == 0);

% LineSpec is red / circle

plot(X(pos, 1), X(pos, 2), 'ro');

% LineSpec is blue / asterisk

plot(X(neg, 1), X(neg, 2), 'b*');

xlabel('x score')

ylabel('y score')

legend('positive', 'negative')

% end chart

hold off;

end

y가 1인 데이터의 인덱스를 pos에 담고, 0인 데이터의 인덱스를 neg에 담아

각각을 plot한 펑션입니다.

legend와 label를 붙여주고 실행하면 아래와 같이 그래프로 나타낼 수 있습니다.

>> plotData(X, y);

실행 하면 아래와 같은 데이터 셋으로 구성 된 것을 확인 해볼 수 있습니다.

이번에는 새터를 자동으로 구해주는 다른 알고리즘인 Stochastic Gradient descent를 알아보도록 하겠습니다.

앞서 배운 Gradient descent에 Stochastic의 접두사를 붙혀서 SGD라고 불립니다.

GD는 "batch" Gradient descent 라고도 불리며 이것은 매 이터레이션 마다 모든 트래이닝 데이터를 구해주기 때문에 계산의 코스트가 높다는 단점이 있습니다.

SGD는 위 단점을 해결하고자 Sum을 하지 않고 자동 새터를 구해주는 알고리즘 입니다.

위와 같이 SGD는 GD와 거의 동일한데 몇가지 다른 점이 있습니다.

1. 먼저 데이터를 셔플 해준다. (가장 큰 특징)

- 이유는 첫번째 부터 접근을 안하고 임의의 수 부터 접근하여 새터값을 빨리 줄이려는 이유인 것 같음. 

(적정한 수가 랜덤으로 선택되어 일정 새터값으로 수렴 되기 때문에)

2. 매 이터레이션 마다 모든 데이터를 더하는 sum of all 생략

3. m으로 나누지 않는다.

위와 같이 진행되며

매번 마다 sum을 안하기 때문에 코스트가 적게 듭니다.

코드는

function [theta, J_history] = stochasticGradientDescent(x, y, theta, alpha, iterations)

    m = length(y); % number of rows

    n = size(x, 2); % number of features

    h = zeros(n, 1);

    J_history = zeros(iterations, 1);

    % data shuffle

    myperm = randperm(m);

    Xshuffle = x(myperm , :);

    Yshuffle = y(myperm);

    for iter=1:iterations

        for i=1:m

            for j=1:n

                h(j) = (theta' * Xshuffle(i, :)' -  Yshuffle(i)) * Xshuffle(i, j);

            end

            for j=1:n

                theta(j) = theta(j) - alpha * h(j);

            end

        end

        % add history

        J_history(iter) = computeCostMulti(Xshuffle, Yshuffle, theta);

    end

end

위와 같이 데이터를 처음 한번 셔플해주고(row의 순서가 바뀜)

sum이 없고 m으로 나눠주는 곳도 없습니다.

위와 같이 구한것에 단항의 데이터를 이용하여 새터와 코스트 펑션을 구해보면

>>  data = load('ex1data1.txt');

>>  n = size(data, 2);

>>  x = data(:, 1:n-1);

>>  y = data(:, n);

>>  m = length(x);

>>  x = [ones(m, 1), x]; % add x0 1

>>  theta = zeros(n, 1);

>>  alpha = 0.001;

>>  iterations = 1000;

>>  [theta, J_history] = stochasticGradientDescent(x, y, theta, alpha, iterations);

>> theta

theta =

  -3.9598

   1.2511

>> computeCostMulti(x, y, theta)

ans =  4.5857

와 앞에서 구한 GD와 비슷한 것을 볼 수 있습니다.

여러번 진행해 보면 진행할 때 마다 새터값이 변경됩니다.

(랜덤으로 셔플하기 때문에 특정 구간으로 수렴하는 것 같은데 정확히 잘 모르겠네요)

일단 여기서 SGD 설명을 마치며 아직 의문이 드는 사항이 몇개 있는데 이것은 추후 공부하고 추가 하도록 하겠습니다.

(현재는 알파값은 GD보다 더 작게, iteration은 더 크게 해야 좋은 수가 나오는 것 같습니다.)

이전에 새터를 자동으로 구해주는 Gradient Descent 기법을 설명하였습니다.

하지만 적정한 알파(σ) 값을 지정해 주어야 하고

iteration 크기도 정해 주어야 값이 구해집니다.

(이말은 적정한 알파와 이터레이션 크기를 안해주면 잘 나오지 않는다는 의미도 됩니다.)

그리하여 알파와 이터레이션을 정해주지 않아도 자동으로 새터를 구해주는 Normal equation 기법에 대해 알아 봅시다.

새터를 구하는 공식은 위와 같습니다.

함수로 구현하면 

function theta = normalEquation(x, y)

theta = pinv(x' * x) * x' * y;

end

이럴경우 최소 오차가 GD와 얼마나 차이가 있을지를 확인해 보겠습니다.

다항일 경우

>> data = load('ex1data2.txt');

>> n = size(data, 2);

>> x = data(:, 1:n-1);

>> y = data(:, n);

>> m = length(x);

>> [x, mu, sigma] = meanNormalization(x);

>> x = [ones(m, 1), x]; % add x0 1

>> theta = normalEquation(x, y)

ans =

  3.4041e+005

  1.1063e+005

  -6.6495e+003

cost function을 구해보면

>> computeCostMulti(x, y, theta)

ans =   2.0433e+009

이며 이전에 구한 값 

theta =

  3.4041e+005

  1.1063e+005

  -6.6493e+003

코스트 펑션은

ans =   2.0433e+009

와 비슷합니다.

단항일 경우 어떨지 다시 비교해보면

>> data = load('ex1data1.txt');

>> y = data(:, 2);

>> m = length(y);

>> x = [ones(m, 1), data(:,1)];

>> theta = zeros(2, 1);

>> theta = normalEquation(x, y)

theta =

  -3.8958

   1.1930

>>computeCost(x, y, theta)

ans =  4.4770

이전에 구한 값은

theta =

  -3.6303

   1.1664

ans =  4.4834

와 같이 비슷합니다.

NE 기법이 GD와 결과도 비슷하고 계산도 빠르다는 것을 볼 수 있습니다.

각각의 장단점은 아래 그림과 같습니다.

GD나 NE를 적절히 사용하라고 하네요.

이전에 Gradient Descent라는 기법을 이용하여 단항에서 자동으로 새터를 구해 주었습니다.

이번에는 다항일 때의 Gradient Descent 기법을 알아 보겠습니다.

다항일 때는 별로 달라지는게 없고 x0이 원래 존재 했었지만 1이여서 생략 했던 것이고

이전에 구한 multiple variable cost function에 h를 구한 식을 가져다가 

n(number of features) 만큼 iteration을 돌리면 됩니다.

식을 그림과 같으며 코드를 구현해 보면

* Gradient Descent Multi

function [theta, J_history] = gradientDescentMulti(x, y, theta, alpha, iterations)

m = length(y); % number of rows

n = size(x, 2); % number of features

    h = zeros(n, 1);

    J_history = zeros(iterations, 1);

    for i=1:iterations

    % sum of all hypotesis * x

    for j=1:n

H = (theta' * x')';

h(j) = sum((H - y)' * x(:, j));

    end

    % update theta

    for j=1:n

    theta(j) = theta(j) - alpha * (h(j) / m);

    end

    % add history

    J_history(i) = computeCostMulti(x, y, theta);

    end

end

다항 코스트 펑션 구하는 코드 + 단항 GD 코드로 구해지는 것을 알 수 있습니다.

>> data = load('ex1data2.txt');

>> n = size(data, 2);

>> x = data(:, 1:n-1);

>> y = data(:, n);

>> m = length(x);

>> [x, mu, sigma] = meanNormalization(x);

>> x = [ones(m, 1), x]; % add x0 1

>> theta = zeros(n, 1)

>> alpha = 0.01;

>> iterations = 1500;

>> theta = gradientDescentMulti(x, y, theta, alpha, iterations)

theta =

  3.4041e+005

  1.1063e+005

  -6.6493e+003

새터가 새로 변경된 것을 볼 수 있습니다.

그럼 코스트 펑션을 다시 구해보면

>> computeCostMulti(x, y, theta)

ans =   2.0433e+009

위와 같이 코스트 펑션이 나옵니다.

이 값이 reasonable 한지는 차차 알아보도록 하겠습니다.

다항일경우 각 feature에 대한 스케일 조정이 필요합니다.

예를들어 size는 4자리의 수를 표현되고 방의개수는 1자리로 표현됩니다.

이렇게 스케일이 안맞을 경우 Gradient Descent에서 찾는데 많이 걸리기 때문에 스케일을 먼저 조정해 줄 필요가 있습니다.

* 그림 설명 - 왼쪽그림은 스케일 조정이 되지 않아 GD가 오래 걸리니 오른쪽 그림처럼 스케일을 조정해야 속도가 빠릅니다.

* 그림 설명 - Mean normalization이라는 기법을 이용하여 x의 범위를 조정해주는 기법

그럼 Mean normalization을 이용하여 x의 변화에 대해 알아보도록 합시다.

새로운 데이터 셋(ex1data2.txt) 파일이다. 

x1, x2, y 순으로 값이 되어 있으며 x1은 size, x2는 방의개수, y는 가격입니다.

2104,3,399900

1600,3,329900

2400,3,369000

1416,2,232000

3000,4,539900

먼저 현재 데이터를 스케일 조정 없이 데이터를 그려보도록 합니다.

>> data = load('ex1data2.txt');

>> n = size(data, 2);

>> x = data(:, 1:n-1);

>> y = data(:, n);

...

>> plot(x(:, 1), y, 'rx')

>> hold on

>> plot(x(:, 2), y, 'bx')

그래프와 같이 0을 중심으로 데이타가 모여 있는 것을 확인 할 수 있습니다.

이후에 앞에서 설명한 cost function 을 구해보면

이전까지는 단항(One Variable) 일 때 선형회귀분석에 대해서 알아 보았습니다.

하지만 현실에서는 단항보다는 여러개의 feature 들(Multiple Variable)로 되어 있는 문제를 많이 접하게 됩니다.

위와 같이 가격을 결정하는데 사이즈, 방의개수, 층수, 년식등 많은 정보들을 필요로 합니다.

이러한 문제 처럼 차례대로 각각을 x1, x2, x3, x4 라고 price를 y라고 설정한 후 단항에서와 같이 똑같이 새터를 구합니다.

다항일 때도 이전에 했던 것과 같이 이 값들을 대표하는 새터를 구하고, 임의의 값에 대한 y를 예측해 보는 것을 설명 하겠습니다.

거의 단항과 비슷하며 몇개 추가된 개념과, feature 개수에 따른 기존 로직 수정이 조금 필요 합니다.

위와 같이 다항일 때는 hypothesis가 구해지는 식입니다.

사실 x0라는 값이 존재하는데 이 값은 1이기 때문에 생략되서 새터 0만 나오는 것입니다.

다항일 경우 위와 같이 간편하게 게산할 수 있습니다. (theta transpose * x)

기존에서 x가 하나 더 추가된 새로운 데이터 셋(ex1data2.txt) 파일 입니다. 

x1, x2, y 순으로 값이 되어 있으며 x1은 size, x2는 방의개수, y는 가격입니다.

2104,3,399900

1600,3,329900

2400,3,369000

1416,2,232000

3000,4,539900

그럼 앞서 설명한 코스트 펑션을 Multiple variable로 확장하여 가장 작은 오차를 찾아 보겠습니다.

* cost function

J = 1/2m * 1 to m 까지 (h(x)- y)제곱 의 더하기

h = theta transpose * x

function J = computeCostMulti(x, y, theta)

m = length(y);

H = (theta' * x')'; % theta transpose * x transpose

S = sum((H - y) .^ 2);

J = S/(2*m);

end

위의 transpose 식을 이용하여 멀티 코스트 펑션을 구해 보았습니다.

참고로 transpose 는 ' 기호를 사용 합니다.

행렬 곱, 행렬 뺄셈 등을 이용하여 여러번 행렬의 트랜스포즈 해주었습니다.

>> data = load('ex1data2.txt');

>> n = size(data, 2);

>> x = data(:, 1:n-1);

>> y = data(:, n);

>> m = length(x);

>> x = [ones(m, 1), x]; % add x zero

>> theta = zeros(n, 1)

>> computeCostMulti(x, y, theta)

ans =   6.5592e+010

위와 같이 새터 [0;0;0] 일때 코스트 펑션을 구해 보았습니다.

이전 블로그에서 데이터에 맞는 theta를 구하였습니다.

theta =

  -3.6303

   1.1664

이 새터를 가지고 그럼 임의의 x값으로 y를 예측 해봅시다.

* 20일 때 y값은 얼마일지?

>> newX = [1, 20]

newX * theta를 통해 y값을 구할 수 있다.

>> predict = newX * theta

predict =  19.697

이와 같이 선형 그래프의 새터를 구하고 x값을 대입하여 y를 구할 수 있습니다.

이것을 그래프로 나타내면

>> plotManualTheta(x, y, theta)

>> hold on

>> plot(20 , predict, 'rx', 'MarkerSize', 25);

와 같이 라인에 정확한 값을 구할 수 있습니다.

앞서 새터로 인하여 cost function와 그래프 모양이 어떻게 변하는지 알아 보았습니다.

여기서는 새터를 자동으로 구해주는 알고리즘인 Gradient Descent 알고리즘에 대해 설명합니다.

* Gradient Descent

gradient descent의 알고리즘은 위와 같습니다.

간단하게 설명하면 

: 입력한 iteration 개수 만큼

: 각 새터를 구하고 그값을 이용하여 iteration 을 돌립니다.

: 알파 값도 입력해줍니다.

자세한 이론은 심화 강의를 통해 학습하기 바랍니다.

alpha 값과 iteration 값을 입력해줍니다.

코드는 아래와 같습니다.

function theta = gradientDescent(x, y, theta, alpha, iterations)

    m = length(y); % number of rows

    for i=1:iterations

        h1 = 0;

        h2 = 0;

        for j=1:m

            h1 = h1 + (theta(1) + theta(2) * x(j, 2) - y(j));

            h2 = h2 + (theta(1) + theta(2) * x(j, 2) - y(j)) * x(j, 2);

        end

        theta(1) = theta(1) - alpha * (h1 / m);

        theta(2) = theta(2) - alpha * (h2 / m);

        % computeCost(x, y, theta); % print cost function

    end

end

: 앞서 구한 cost function과 유사하며 iterations, alpha 값을 더 입력해준 것이 특징입니다.

그럼 알파, 이터레이션 값까지 하여 알고리즘을 통해 새터를 구한 값을 대입해서 계산합니다.

>> data = load('ex1data1.txt');

>> y = data(:, 2);

>> m = length(y);

>> x = [ones(m, 1), data(:,1)];

>> theta = zeros(2, 1);

>> alpha = 0.01;

>> iterations = 1500;

>> theta = gradientDescent(x, y, theta, alpha, iterations)

theta =

  -3.6303

   1.1664

구한 새터값은 위와 같다. cost function을 구해봅니다.

>> computeCost(x, y, theta)

ans =  4.4834

값이 4.4로 좋아 보입니다.

그래프로 나타나면

>> plotManualTheta(x, y, theta)

위와 같이 꽤 괜찮은 그래프가 그려지는 것을 볼 수 있습니다.

우리의 목표는 cost function값이 최소가 되는 새터를 찾는다고 하였습다.

먼저 새터를 찾는다는 것이 무엇인지 확인해 볼 필요가 있을 것 같습니다.

앞서 plot을 통해 x, y를 표시하던 그래프에 cost function으로 오차율을 나타내고

실제 새터를 그래프로 그려보면서 수동으로 새터를 찾아 보겠습니다.

먼저 새터를 그리는 함수부터 정의해봅니다.

* 새터를 표시하는 plot 함수

function plotManualTheta(x, y, theta)

plot(x(:,2) , y, 'rx', 'MarkerSize', 5); % first plot

hold on; % keep previous plot visible

plot(x(:,2), x*theta, '-') % second plot

legend('Training data', 'Linear regression') % legend

hold off % don't overlay any more plots on this figure

end

: hold on - plot을 고정 시키는 명령어

: legend - 범례

: hold off - plot 고정 해제

: plot(x(:,2), x*theta, '-') 

 - 두번째 plot을 그리는 x, y 데이터 입니다.

 - 여기서 두번째 파라미터에는 x와 theta에서 들어온 값을 곱해 줍니다.

 - 행렬 곱을 이용하여 x = 1 * 2, y = 2 * 1 행렬을 곱하여 1*1행렬로 계산하여 그래프를 표시 합니다.

예제를 통해 알아 보겠습니다.

>> data = load('ex1data1.txt');

>> y = data(:, 2);

>> m = length(y);

>> x = [ones(m, 1), data(:,1)];

>> theta = zeros(2, 1);

* 새터가 0,0 일때

>> computeCost(x, y, theta)

h =  6222.1

ans =  32.073

>> plotManualTheta(x, y, theta)

* 새터가 0, 1일때

>> theta = [0; 1]

>> computeCost(x, y, theta)

h =  1444.5

ans =  7.4459

>> plotManualTheta(x, y, theta)

위와 같이 계속 하다 보면 최적의 새터를 찾을 수 있습니다.

문서 초반에 설명을 못했지만 이 글을 쓰고 있는 목적은 선형회귀분석을 구하는 것을 octave를 이용하여 구하는게 목적이였습니다.

추후 스칼라를 이용하여 머신러닝을 구현하는것이 최종 Goal.

최종 골은 가장작은 Cost Function을 구하는거고 (오차율 최소화)

거기에 맞는 새터를 구하는것입니다.

(선형 회귀 분석에서 가장 잘 맞는 라인을 찾는 것 이라고 생각하면 됩니다.)

자세한 이론은 추후 설명하고

이번에는 코스트 펑션을 구하도록 구현 해보겠습니다.

식은 아래와 같습니다.

1/2m * 1 to m 까지 ((새터0 + 새터1 * x) - y)제곱 의 더하기

cost function을 구현해보면

* 데이터 로딩

>> data = load('ex1data1.txt')

* y 변수 선언

>>  y = data(:, 2);

* m 변수 선언

>> m = length(y)

m = 97

: 행의 개수

* x 변수 선언

>> x = [ones(m, 1), data(:,1)];

: x변수에 첫번째 컬럼에 1값을 대입시킴

: 이유는 cost function의 h에 새터0 + 새터1 * x의 식인데 사실 새터 0에도 x를 곱해야 하는데 1이라 식에서 생략해서 그렇습니다.

x =

    1.0000    6.1101

    1.0000    5.5277

    1.0000    8.5186

    1.0000    7.0032

* 새터 선언

>> theta = zeros(2, 1);

: 0으로 초기화

theta =

   0

   0

* cost function 구현

function J = computeCost(x, y, theta)

m = length(y); % number of rows

J = 0;

h = 0;

for i=1:m

h = h + (theta(1) + theta(2) * x(i, 2) - y(i))^2;

end

h % print hypothesis

J = h / (2*m);

end

* cost function 호출

>> computeCost(x, y, theta)

h =  6222.1

ans =  32.073

32.073이라는 값이 리턴됩니다.

이 값을 최대한으로 줄이는 새터를 찾는게 이 프로그램의 최종 목표입니다.