[Scipy+Numpy] 3. 2D Histogram + [Matplotlib] 11. Plot Main: Pcolorm...

글의 내용은 [SciPy] 2. Density Estimation (Application to Scatter Plot)에서 이어집니다.

전 편에서 산점도 Scatter Plot은 단점으로 너무 많은 데이타를 그릴 경우 겹쳐지는 부분 때문에 밀도를 제대로 가늠할 수 없다는 얘기를 했습니다. 그래서 데이타의 밀도를 계산하여 Contour 를 위에 올려놨었죠. 이렇게요:

위 그림에는 X, Y, 그리고 C (Color) 이렇게 3가지 변수가 얽혀있는 모습입니다.

이후 데이타의 특성을 더 잘 파악하기 위하여 C에 따라 4 그룹으로 나눴습니다. 흔히 말하는 Quartile로 나눴어요. Quartile은 값을 한줄로 줄세운 뒤 (sort) 균등한 머리수(?)로 그룹을 나누는 방식이죠.
1Q: 0-25%, 2Q: 25-75%, 3Q: 75-100%

C값에 따라 나눴으므로 색깔이 확연히 구분되는 건 당연한건데, 밀도 계산 만으로는 뭐랄까, 좀 명확하지가 않아서요. 보이는 모양도 Contour를 몇 개 그리느냐에 따라 조금 달라 보이기도 하고. 그래서 원 상태의 분포는 과감히 포기하고 밀도에 집중하기 위해 2D 히스토그램을 이용하기로 합니다.

Numpy.histogram2d

numpy.histogram2d(x, y, bins=10, range=None, normed=None, weights=None, density=None)
"Compute the bi-dimensional histogram of two data samples."

Returns:
H : ndarray, shape(nx, ny)
xedges : ndarray, shape(nx+1,)
yedges : ndarray, shape(ny+1,)

Numpy.histogram2d의 사용법은 Numpy.histogram (1d)와 유사합니다. 처음의 인풋만 데이타 1줄이 아닌 2가지 데이타 (위에서 x와 y)를 넣어주면 됩니다. 히스토그램을 나누는 경계값도 x와 y에 따라 각각 입력하면 더 정확한 결과를 얻을 수 있겠죠.

이렇게 해서 나온 결과 H를 그림으로 표현해야 하는데, Bin 경계값이 일률적이지 않은 상황까지 고려하면 pcolormesh가 적절하다고 합니다. 간단한 사용법은 다음과 같습니다.

>>> H, xedges, yedges = np.histogram2d(x, y, bins=(xedges, yedges))
>>> H = H.T  # Let each row list bins with common y range.
>>> X, Y = np.meshgrid(xedges, yedges)
>>> ax.pcolormesh(X, Y, H)

여기서 주의해야 할 부분이 3가지가 있는데,

첫째는 histogram2d의 결과물인 H를 Transpose 해줘야 한다는 겁니다. 왜냐하면, 위에 적혀있듯이 히스토그램 함수를 통해 나오는 H는 (nx,ny) 크기의 행렬이거든요. 그런데 보통 그림그릴 때 들어가는 행렬은 (ny,nx)의 형태여야 합니다. 이것은 Python이 C와 같이 Row-major order 여서 그런건데, 그냥 간단히 외우자면, Python에서 모든 행렬은 가장 바깥쪽(오른쪽)이 먼저입니다. 보통 그림은 x 먼저, 그리고 y의 순서이므로 x가 가장 오른쪽에 위치하여야 합니다.

둘째는 pcolormesh는 먼저 각 데이타를 그림으로 표시할 좌표계 그물망을 계산해서 제공해야 합니다. 그러면 그 그물망에 맞춰서 데이타를 그림으로 표시해요. 그물망 만드는 방법은 Numpy.meshgrid 함수를 이용하는데, 이건 또 (X,Y) 형태로 합니다. 헷갈리죠. ^^

셋째는, pcolormesh 라는 함수 자체가 격자점이 아니라 격자의 4 꼭지점과 모서리로 둘러싸인 사각형 공간을 적당한 색으로 채우는 함수입니다. 따라서 격자 그물망은 데이타 행렬보다 크기가 1씩 커야합니다. 3x3 의 격자점 안에는 사각형이 4개 들어있으니까요.

그렇게해서 나온 그림입니다.

"Colormap.set_under()" 함수를 이용하여 값이 0일때는 회색 ("0.8")으로 처리하였습니다. 그런데 아래에 Colorbar를 추가한다고 해도 값을 구별하기가 어려워서 가장 큰 값 근처에는 직접 숫자를 적어봅니다.

숫자를 표시하는 방법은 다음과 같습니다.

    for n in range(ny):
        for m in range(nx):
            if H[n,m]> some good condition:
                xloc=(xedges[m]+xedges[m+1])/2.
                yloc=(yedges[n]+yedges[n+1])/2.
                ax.annotate("{:4.1f}".format(H[n,m]),xy=(xloc,yloc),ha='center',va='center',stretch='semi-condensed',fontsize=12)

자 이제 다 되었을까요?

그림 그리는 것만 따지자면 다 되었습니다. 하지만 어떤 그림이 가장 "예쁘냐"는 측면에서는 아직 멀었죠. 이 히스토그램의 경우 Bin을 몇개로 정하느냐에 따라 느낌이 확연히 다르니까요.

Bin 개수가 너무 작으면 선명도가 떨어지고, 너무 많으면 분산되어서 집중이 안됩니다. 이중에서 제가 전달하고자 하는 메세지를 가장 잘 표현할 만한걸 정해야죠. 하나하나 보면서...

[Scipy+Numpy] 3. 2D Histogram + [Matplotlib] 11. Plot Main: Pcolormesh

Matplotlib List

F2PY List

Scipy+Numpy List