Nmatrix에 고통받는 개발자와 사라진 2주

February 23, 2017

얼마 전에 밑바닥부터 시작하는 딥러닝이라는 책을 읽고 심심풀이 삼아서 nmatrix를 사용해서 예제를 구현해보고 있었습니다. 그리고 2주 정도 씨름을 한 끝에 포기선언을 한 것이 2일 전. 여기에 대한 감상글이라고 쓰고 불평글이라고 읽습니다을 간단하게 정리해봅니다.

nmatrix 객체는 직렬화할 수 없습니다

첫번째 좌절 포인트. 데이터의 사이즈가 꽤나 묵직하기 때문에 초기 데이터 로딩이 시간이 꽤 많이 걸립니다. 구체적으로는 10~20초 내외. 책에서는 MNIST 데이터셋을 읽어온 뒤, numpy.array로 변환해 피클로 떨구고, 이후 예제에서는 이 파일을 읽어서 바로바로 사용합니다. 그래서 루비에서도 같은 방식으로 진행하려 했으나, gzip 푸는 데 한 세월 걸린건 둘째치고 nmatrix의 객체는 마셜 직렬화를 지원하지 않는다는 사실이 발각. 덕분에 매번 데이터를 루비의 일반 배열로 읽은 뒤, 이를 변환하게 됩니다. 컴파일 하는 것도 아니고 실행마다 매번 이를 기다리는 것은 상당한 고문이었습니다.

There is no Broadcast

소프트맥스 함수를 구현할 때에는 오버플로우를 피하기 위해서 각 데이터 셋의 가장 큰 값을 전체에서 빼는 연산을 넣고 있습니다. 다음과 같은 식인데요.

def softmax(x):
    if x.ndim == 2:
        x = x.T
        x = x - np.max(x, axis=0) # 여기
        y = np.exp(x) / np.sum(np.exp(x), axis=0)
        return y.T 

    x = x - np.max(x) # 여기
    return np.exp(x) / np.sum(np.exp(x))

이 예제에서 행렬의 차수에 따라서 구현이 갈리는 이유는 결과 데이터를 배치 처리하는 경우, 그렇지 않은 경우 계산 방식이 바뀌기 때문입니다. 문제가 되는 부분은 배치처리를 하는 부분(x = x - np.max(x, axis=0))인데, 각 데이터셋에서 최대값을 구한 뒤, 각각 매치되는 데이터셋에서 빼줍니다. 얼핏 보면 자연스럽게 느껴지는 이 코드의 이면에는 브로드캐스팅이라는 마법이 존재하는데요, 예시로 들어온 입력이 (2, 10)의 행렬이라고 가정합시다. 각 행이 데이터셋이라고 생각하면 10개의 값을 가지는 2개의 데이터셋이 되겠죠. 이 행렬의 크기를 잘 생각하면서 위의 코드를 살펴볼까요.

다른 크기의 행렬끼리 뺄셈을 해주고 있는 걸 확인하셨나요? 바로 여기가 브로드캐스팅이라는 마법이 발생하는 부분입니다. 좀 더 구체적으로 설명하자면, (1, 2)의 행렬을 행 방향으로 (10, 2) 크기로 확장하여 계산합니다. 구체적인 규칙이 궁금하다면 여기를 확인하세요.

자, 그럼 이제 이걸 브로드캐스팅과 유사한 기능이 없는 nmatrix로 구현하려면 어떻게 해야할까요? 별 수 있나요, 반복문을 돌려야지. XD

def softmax(x)
  x_row, x_col = x.shape
  x = x.reshape(x_col) if x_row == 1

  if x.dimensions == 2
    x = x.transpose
    # 1 => x - np.max(x, axis=0)
    x_max = x.max(0).to_a
    x.each_row do |row|
      (0...x_max.size).each { |i| row[i] -= x_max[i] }
    end

    y = np_exp(x)
    # 2 => np.exp(x) / np.sum(np.exp(x), axis=0)
    y_sum = y.sum(0).to_a
    y.each_row do |row|
      (0...y_sum.size).each { |i| row[i] /= y_sum[i] }
    end
    y.transpose 
  else
    x = x - x.max
    exp_x = np_exp(x)
    exp_x / np_sum(exp_x)
  end
end

예쁘지는 않지만 적당히 동작하는 코드입니다. 이쪽이 설명하기 용이하니 그대로 방치했다는 변명과 함께 잠깐 살펴보면, 브로드캐스팅을 사용할 수 없어서 직접 연산한 부분이 두 곳 있음을 확인할 수 있네요. 까짓거 손으로 계산좀 할 수 있지, 라고 생각한 독자분에게 질문을 하나. 만약 CNN을 배치처리하는 경우, 입력이 3차원이고 배치 처리라면 한 차원을 추가해야하므로 레이어에서는 4차원 행렬을 처리해야 합니다. 매번 브로드캐스팅이 필요할 때마다 직접 구현할 생각이 있으신가요?

연산자 오버로딩의 불완전 지원

시그모이드 함수는 numpy의 지원을 받으면 아래와 같이 표현할 수 있습니다.

def sigmoid(x):
    return 1 / (1 + np.exp(-x))

잊으면 안되는 부분은 x가 다차원 배열이라는 점입니다. 스칼라값과 행렬간의 연산이라고 하면 머릿속에 물음표 부호가 떠오를만 하지만, element-wise한 연산이라는 점을 이해한다면 그렇게 복잡하지도 않습니다. 요는 스칼라값을 행렬의 각 요소에 대해서 전부 연산하므로 위의 코드는 다음과 같이 정리할 수 있겠네요.

그럼 이제 이걸 nmatrix로 구현해보죠. 일단 np.exp에 해당하는 유틸리티 함수가 없으므로 우리는 각 요소에 대해서 필요한 연산을 해주는 함수가 필요합니다. 여기 설명에서는 그렇게 중요한 부분이 아니므로 이 연산을 처리하는 np_exp라는 함수가 있다고 가정하죠.

def sigmoid(x)
  1 / (1 + np_exp(-x))
end

return이 없어서 한결 깔끔하고 좋아보이는 것은 기분탓입니다. 왜냐하면 실행이 되지 않으니까요.

TypeError: NMatrix can't be coerced into Integer
  from (irb):4:in `+'
  from (irb):4
  from /.rbenv/versions/2.4.0/bin/irb:11:in `<main>'

nmatrix에서는 행렬 (연산자) 스칼라값은 지원하지만 스칼라값 (연산자) 행렬을 지원하지 않기 때문에 발생하는 문제입니다. 그러므로 우선 1 + np_exp(-x)np_exp(-x) + 1로 바꿔보죠.

def sigmoid(x)
  1 / (np_exp(-x) + 1)

당연하지만 여전히 에러가 발생합니다. 왜냐하면 아직 스칼라값을 행렬로 나누는 연산이 있기 때문이죠. 나누는 행렬과 같은 사이즈의 행렬을 만드는 수밖에 없지만, 위에서 말했듯 브로드캐스팅이 지원되지 않으므로 선택지는 ‘손으로 구현한다’ 뿐입니다. 다행히도 넘겨받은 행렬과 같은 크기의 행렬을 만들고 1로 초기화해주는 ones_like라는 함수를 가져다 쓸 수 있습니다.

def sigmoid(x)
  b = np_exp(-x) + 1
  t = NMatrix.ones_like(b)
  t / b
end

이예! 축하합니다. 시그모이드 함수 구현에 성공했습니다. 어떠신가요?

결론

처음 numpy를 사용하면서 느낀 점은 혼란스럽기 그지없고 동작의 일관성을 파악하기 어려운 라이브러리네, 였지만 가면 갈수록 이게 없으면 신경망을 만드는건 사람이 할 일이 아니라는 생각만 들게 됩니다.

덕분에 이 라이브러리를 보기 전이라면 모를까, 본 이후인 지금으로서는 nmatrix는 신경망과는 관계없는 라이브러리로 보입니다.

나름대로 머리를 굴려가며 짠 코드인데, 더 나은 방법은 분명 있을 겁니다. nmatrix를 그대로 활용하든 보완하는 코드를 덧붙이든 말이죠. 문제는 관련 정보가 너무 부족해서 더 좋은 설명이나 코드를 볼 수 없다는 점이라고도 할 수 있겠네요. 하다못해 코드를 읽고 싶어도 C라서 읽기도 어렵고…

narray의 새 버전은 numpy를 많이 의식했는지 브로드캐스팅 등등 여러 기능을 지원한다는 모양이라 좀 기대되지만 아직 개발 도중이므로 뭐라고 평가하기가 어렵네요.

결론으로 현재 시점으로 그럭저럭 쓸만하다는 다차원 행렬 처리 라이브러리가 이 정도라면 루비로 신경망과 관련된 뭔가를 해보겠다는 생각은 버리는게 좋겠다는 답에 도달했습니다.

Reference