Category: Python

CellProfiler의 UnMix 모듈로 색을 뽑아내는게 가장 대조가 좋아서, 어떻게 해서든 구현해 보고 싶었다. 다운로드 화면에서 코드가 GitHub에 있다는 내용을 봐서 코드를 보고 참고해서 구현해 보기로 했다. 따라서, 이 내용은 다음의 라이센스를 따른다.

The BSD 3-Clause License
Copyright ⓒ 2003 – 2020 Broad Institute, Inc. All rights reserved.

Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:

1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
3. Neither the name of the Broad Institute, Inc. nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED “AS IS.” BROAD MAKES NO EXPRESS OR IMPLIED REPRESENTATIONS OR WARRANTIES OF ANY KIND REGARDING THE SOFTWARE AND COPYRIGHT, INCLUDING, BUT NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, CONFORMITY WITH ANY DOCUMENTATION, NON-INFRINGEMENT, OR THE ABSENCE OF LATENT OR OTHER DEFECTS, WHETHER OR NOT DISCOVERABLE. IN NO EVENT SHALL BROAD, THE COPYRIGHT HOLDERS, OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF, HAVE REASON TO KNOW, OR IN FACT SHALL KNOW OF THE POSSIBILITY OF SUCH DAMAGE.

If, by operation of law or otherwise, any of the aforementioned warranty disclaimers are determined inapplicable, your sole remedy, regardless of the form of action, including, but not limited to, negligence and strict liability, shall be replacement of the software with an updated version if one exists.

Development of CellProfiler has been funded in whole or in part with federal funds from the National Institutes of Health, the National Science Foundation, and the Human Frontier Science Program.

변환에 앞서 알아야 하는 것은 다음과 같다.

• 이미지에 포함된 모든 염색 종류를 포함해야 한다. 그러니까 H&E 염색이라면 H와 E를 모두 필요하다. CellProfiler에서 기본 제공하는 것이 아니라면 비슷한 것을 고르거나, Custom으로 입력할 수 있다. Custom에는 분석에 사용된 RGB 흡광도를 추정할 수 있는 메뉴가 있다.
• 2가지 색상이 있다면 색상 순서의 선후 관계는 없다.
• 일련의 변환 과정을 거치면 1채널 gray 값으로 변환된다.
• 0~255 범위가 아닌 0~1 범위를 이용하는 것으로 추정된다.
• 변환 과정 중에서 지수와 로그를 이용하기 때문에 0인 경우를 처리하기 위하여 작은 수 epsilon을 더해주며, 기본값은 1/256/2 이다. 0을 제외한 최소값 1/256의 1/2 값이다.
• 나도 이 결과에 책임져 줄 수는 없으나, 꽤 비슷했다는 경험이 있다.

코드 설명

col1 = (.. ,.. ,.. )
col2 = (.. ,.. ,.. )
stain = np.array([col1, col2])
stain = stain / np.sqrt(np.sum(stain ** 2))
ext = np.matrix(stain).I[:,0].flatten().tolist()

im = crop.copy() / 255
re_im = np.exp(np.sum(np.log(im + eps) * ext, 2))
re_im -= eps
re_im[re_im < 0] = 0; re_im[re_im > 1] = 1
re_im = (1 - re_im) * 255
im_th = cv2.threshold(re_im, 120, 255, cv2.THRESH_BINARY)[1]
im_th = np.array(im_th, dtype=np.uint8)

• 소스 코드에서 선택한 염색의 흡광도를 이용한다. 여러개를 하나의 array로 묶은 다음 일련의 변환 과정을 거친다.
• 나의 경우는 col1 값의 이미지만 필요하기 때문에 0번째 값을 선택하여 ext로 저장했다.
• 수치를 0~1로 변환한 다음 복잡한 과정을 거쳐 변환하여 re_im에 저장한다.
• 초기에 더한 epsilon을 제거 뺀다.
• 혹시 변환값이 0보다 작거나 1보다 크면 0과 1로 되도록 한다.
• 1에서 뺀 값을 이용한다(inverse).
• findContours에서 이상하게 오류가 나서 255를 다시 곱했다. 그리고 uint8으로 변환하여 주면 동작한다.

cvtColor로 BGR을 RGB로 바꾸었다면 순서대로 RGB로 배열되어 있다. 하나씩 분리한다.

import cv2
import numpy as np
import matplotlib.pyplot as plt

src = cv2.cvtColor(cv2.imread("seq2337.tif"), cv2.COLOR_BGR2RGB)
im_r = im[:,:,0]
im_g = im[:,:,1]
im_b = im[:,:,2]

자료 타입이 uint8으로 되어 있다. 일정한 정도를 빼거나 더할 경우 아마 255 범위를 넘어가면 다시 반복한다. 그래서 일반 정수형으로 변환이 필요하다. 내가 시도한 것은 일정한 정도를 뺐을 때 0보다 작으면 0으로 치환하는 것이다.

im_r1 = np.array(im_r, dtype=np.int64) - 255
im_r1[im_r1[:] < 0] = 0

이렇게 변환한 값을 다시 RGB로 보고 싶으면 uint8으로 변환한 다음에 쌓아준다. dstack을 이용한다.

im_r1 = np.array(im_r1, dtype=np.uint8)
re_im = np.dstack((im_r1, im_g, im_b))

이건 굉장히 쉽다. 원본 이미지 변수에서 행과 열의 순서대로 범위를 지정하면 된다. 어려운 점은 0부터 시작한다는 것이다. 내가 잘 모르는 것이라면.. 100픽셀로 하려면 하나 더 더해주어야 한다는 것이다. 다음의 식으로 쭉 훑는 방법을 적용해 보았다.

src = cv2.cvtColor(cv2.imread("seq2337.tif"), cv2.COLOR_BGR2RGB)
src_height = src.shape[0]
src_weight = src.shape[1]
for i in range(0, np.int(src_height/100)+1):
    for ii in range(0, np.int(src_weight/100)+1):
        crop = src[(i*100):(i*100+99+1), (ii*100):(ii*100+99+1)]

OpenCV에서 이미지를 복제할 때 R 처럼 src2 = src 했더니 상호 연동이 되었다. .copy()를 이용해서 해야 독립적이다.

src_color = cv2.cvtColor(cv2.imread("seq2337_1.tif"), cv2.COLOR_BGR2RGB)
src_color2 = src_color.copy()