Kaggle lung cancer detection: Tiền xử lý ảnh (Preprocessing)

kaggle-image-preprocessing

kaggle-image-preprocessing

Tiền xử lý ảnh (image preprocessing) là tiến trình cải thiện chất lượng ảnh sao cho ảnh đầu ra giữ lại được những đặc trưng (feature) quan trọng trong ảnh. Sau đó, làm đầu vào cho các quá trình xử lý tiếp theo. Ta có thể thay đổi độ sáng của ảnh (pixel brightness transformation), biến đổi hình học ảnh (geometric transformation), rút trích đặc trưng ảnh (feature extraction), nội suy ảnh (image interpolation), lọc ảnh (image filtering: blur, sharpen, noise removal), …

Trong bài viết này, ta sẽ học cách làm việc với tập dữ liệu ảnh số. Ta sẽ tìm hiểu về cách ảnh số được định dạng trong không gian như thế nào. Bắt đầu bằng việc đọc và hiển thị ảnh cho đến áp dụng tích chập (convolution) và lọc ảnh bằng Tensorflow. Ở đây, tôi sẽ sử dụng tập dữ liệu download từ Kaggle Data Science Bowl 2017. Do tập dữ liệu khá lớn, các bạn có thể download sample_images về làm việc.

Git source: basics.py

Đọc và hiển thị ảnh số

Đầu tiên, ta sẽ load các libraries cần thiết để load một bức ảnh bất kỳ.

import matplotlib.image as mpimg
import matplotlib.pyplot as plt

Ta dùng hàm imread() để đọc dữ liệu ảnh và imshow() để hiển thị ảnh.

img = mpimg.imread("imgs/dogs.jpg")

print "Image data", img.shape
print img

print "Show image"

plt.style.use('ggplot')
plt.colorbar()
plt.imshow(img)

plt.show()

# output
Image data (380, 499, 3)
[[[ 88  97 112]
  [ 82  91 106]
  [ 87  97 109]
  ...,
  [182 188 118]
  [177 183 113]
  [177 183 113]]

 [[194 201 209]
  [125 132 140]
  [ 25  29  38]
  ...,
  [134 137  66]
  [130 133  62]
  [131 134  63]]]
show-image

show-image

img là biến có kiểu dữ liệu ndarray, là mảng 3 chiều chứa các giá trị màu RGB của ảnh. Image data (380, 499, 3) nghĩa là bức ảnh này có chiều rộng (Height) 380 pixel, chiều dài (Width) 499 pixel và 3 Channel (kênh màu tương ứng Red, Green, Blue). Ta có thể thấy, mỗi điểm ảnh là tập hợp của ba màu cơ bản tạo thành. Ví dụ [88 97 112] là một điểm ảnh có Red:88, Green:97, Blue:112. Các điểm ảnh này nối tiếp nhau để hoàn thiện một bức ảnh.

Ta có thể quan sát từng kênh màu RGB như sau:

print "Show RGB channels"
plt.imshow(img[:, :, 0], cmap='gray')

plt.show()
plt.imshow(img[:, :, 1], cmap='gray')

plt.show()

plt.imshow(img[:, :, 2], cmap='gray')
plt.show()

Như vậy, với bất kỳ định dạng ảnh nào, ta chỉ cần thu được mảng dữ liệu ảnh là có thể làm việc được. Tương tự, đối với định dạng ảnh Y khoa DICOM, ta sử dụng thư viện pydicom để đọc và làm việc với định dạng ảnh này.

lung-example

lung-example

Giả sử bạn đã download thành công sample_images về máy. Tôi sẽ đọc folder ảnh của bệnh nhân có id 00cba091fa4ad62cc3200a657aeb957e. Sau đó, lưu tất cả ảnh đã đọc được vào mảng ảnh images như sau:

import glob
import os
import matplotlib.pyplot as plt
import numpy as np

from pydicom import dicomio

root_dir = "sample_images/00cba091fa4ad62cc3200a657aeb957e/"
os.chdir(root_dir)
images = []
for f in glob.glob("*.dcm"):
    ds = dicomio.read_file(f)
    img = ds.pixel_array
    images.append(img)

# convert to array
data = np.array(images)
print "Total images:", len(images)
print "Image dimensions:", images[0].shape
print "Combine dimensions:", data.shape

# output
Total images: 134
Image dimensions: (512, 512)
Combine dimensions: (134, 512, 512)

Như vậy, folder này có tất cả 134 ảnh lát cắt của phổi và mỗi ảnh như vậy có kích thước 512×512. Bây giờ, ta đã có tất cả dữ liệu ảnh của bệnh nhân. Từ đây, ta có thể áp dụng nhiều phép biến đổi trên ảnh. Ví dụ như dưới đây, ta sẽ lấy trung bình và phương sai của các điểm ảnh để tạo thành một ảnh mới.

print "Calculating mean images"
mean_img = np.mean(data, axis=0)
plt.imshow(mean_img.astype(np.uint8))
plt.show()

print "Calculating deviation images"
std_img = np.std(data, axis=0)
plt.imshow(std_img.astype(np.uint8))
plt.show()
mean-of-images

mean-of-images

deviation-of-images

deviation-of-images

Với cách biểu diễn này, ta đã thu gọn tập dữ liệu gồm 134 bức ảnh xuống còn 1 bức ảnh như trên nhờ vậy ta đã giảm tải cho mô hình học máy sau này.

Histograms

Ta có thể sử dụng histograms để quan sát bức ảnh theo một hướng khác. Đầu tiên, ta cần biến đổi các bức ảnh sang mảng 1 chiều thay vì mảng 4 chiều như lúc đầu (N x H x W x C) và dùng hàm hist() để hiển thị histogram của ảnh.

Ghi chú:

  • N: số lượng ảnh trong mảng (không đề cập nếu đang xét 1 ảnh)
  • H: chiều rộng ảnh (height)
  • W: chiều dài ảnh (width)
  • C: kênh màu (channel)
# convert to flattened array
flattened = data.ravel()
print "First image:", data[:1]
print "First 10 values:", flattened[:10]

print "Histogram"
plt.hist(flattened, 255)
plt.show()

# output
First image: [[[-2000 -2000 -2000 ..., -2000 -2000 -2000]
  [-2000 -2000 -2000 ..., -2000 -2000 -2000]
  [-2000 -2000 -2000 ..., -2000 -2000 -2000]
  ...,
  [-2000 -2000 -2000 ..., -2000 -2000 -2000]
  [-2000 -2000 -2000 ..., -2000 -2000 -2000]
  [-2000 -2000 -2000 ..., -2000 -2000 -2000]]]
First 10 values: [-2000 -2000 -2000 -2000 -2000 -2000 -2000 -2000 -2000 -2000]
histogram

histogram

Ta có thể thấy giá trị điểm ảnh phân bố nhiều xung quanh các giá trị -2000, 0, 1000. Tiếp theo, ta quan sát histogram của trung bình ảnh sẽ trông như thế nào

print "Histogram Equalization"
plt.hist(mean_img.ravel(), 255)
plt.show()
histogram-equalization

histogram-equalization

Histogram cho thấy các giá trị phân bố đều vào khoảng giá trị từ 0 đến 1000. Ngoài ra, khi chúng ta lấy một bức ảnh bất kỳ trừ cho trung bình ảnh, ta đã loại bỏ mọi thông tin không cần thiết ra khỏi ảnh, còn lại là các thông tin quan trọng miêu tả cho bức ảnh hiện tại.

print "Normalizing our data"
bins = 20
fig, axs = plt.subplots(1, 3, figsize=(12, 6), sharey=True, sharex=True)
axs[0].hist(data[0].ravel(), bins)
axs[0].set_title("img distribution")
axs[1].hist(mean_img.ravel(), bins)
axs[1].set_title("mean distribution")
axs[2].hist((data[0] - mean_img).ravel(), bins)
axs[2].set_title("(img - mean) distribution")
plt.show()
histogram-normalizing

histogram-normalizing

Khởi động TensorFlow

Trong phần này, ta sẽ làm quen với cách làm việc của TensorFlow. Về cơ bản mọi hoạt động đều được diễn ra trong Graph và các phép tính toán được thực hiện thông qua Session. Nhưng trước hết, ta hãy import library tensorflow vào chương trình của chúng ta và tạo mảng 10 phần tử nằm trong khoảng -2.0 đến 2.0 với kiểu dữ liệu số thực.

x = tf.linspace(-2.0, 2.0, 10)

print x

# output
Tensor("LinSpace:0", shape=(10,), dtype=float32)

Không giống như cơ chế của numpy, ở đây ta không in ra được các giá trị của biến x, thay vào đó là các thông tin liên quan đến biến này gồm tên của biến x “LinSpace:0”, số chiều của x là 10, và mỗi phần tử có kiểu dữ liệu số thực. Lý do là bởi Tensor chưa thực hiện tính toán cho khai báo này mà phải thông qua graph để thực thi lệnh.

Tiếp theo, ta sẽ truy xuất Graph mặc định và quan sát các operations bên trong

g = tf.get_default_graph()
print [op.name for op in g.get_operations()]

# output
[u'LinSpace/start', u'LinSpace/stop', u'LinSpace/num', u'LinSpace']

Ta có thể thấy Tensorflow đặt tên cho mỗi phép tính (operations) của chúng ta để phản ánh lại những gì sắp được tính toán. Để có thể truy xuất đến một operation bất kỳ (tensor) ta thực hiện như sau

print g.get_tensor_by_name('LinSpace' + ':0')

# output
Tensor("LinSpace:0", shape=(10,), dtype=float32)

Sau cùng, để có thể thực hiện tính toán trong Tensorflow, ta cần phải tạo một tf.Session. Sau đó, ta có thể sử dụng hàm run() để tiến hành tính toán hoặc sử dụng chính biến tensor x vừa rồi để thực thi bằng cách truyền vào tham số là session vừa tạo.

# Create Session
sess = tf.Session()

# Tell session to compute
print "Session computes"
computed_x = sess.run(x)
print(computed_x)

# Evaluate itself using this session
print "Variable evaluates"
computed_x = x.eval(session=sess)
print(computed_x)

print "Tensor shapes"
print(x.get_shape())
# convert to list format
print(x.get_shape().as_list())

# Close the session
sess.close()

# output
Session computes
[-2.         -1.55555558 -1.11111116 -0.66666663 -0.22222221  0.22222233
  0.66666675  1.11111116  1.55555558  2.        ]
Variable evaluates
[-2.         -1.55555558 -1.11111116 -0.66666663 -0.22222221  0.22222233
  0.66666675  1.11111116  1.55555558  2.        ]

Ngoài ra, ta có thể chỉ định session làm việc với graph nào

sess = tf.Session(graph=g)

sess.close()

Thông thường, ta sẽ làm việc với Graph mặc định. Tuy nhiên, ta vẫn có thể tạo ra Graph mới như sau

g2 = tf.Graph()

Để đơn giản hóa quá trình prototype, ta có thể sử dụng tf.InteractiveSession

sess = tf.InteractiveSession()
print x.eval()

# output
[-2.         -1.55555558 -1.11111116 -0.66666663 -0.22222221  0.22222233
  0.66666675  1.11111116  1.55555558  2.        ]

Convolution với 2-D Gaussian Kernel

Convolution là cách chúng ta filter ảnh để đạt được một vài hiệu ứng như làm mờ ảnh, làm sắc nét ảnh, phát hiện biên cạnh trên ảnh, … Trong phần này, ta sẽ tạo một kernel từ hàm Gaussian. Đầu tiên, ta sẽ tạo đường cong Gaussian dựa vào dãy số x vừa tạo cùng với công thức phát sinh đường cong

f(x, \mu, \sigma^2) = \frac{1}{\sqrt{2\sigma^2\pi}} e^{-\frac{(x - \mu)^2}{2\sigma^2}}

mean = 0.0
sigma = 1.0

z = (tf.exp(tf.neg(tf.pow(x - mean, 2.0) /
                   (2.0 * tf.pow(sigma, 2.0)))) *
     (1.0 / (sigma * tf.sqrt(2.0 * 3.1415))))

res = z.eval()
plt.style.use("ggplot")
plt.plot(res)
plt.show()
gaussian-curve

gaussian-curve

Ta tạo 2D Gaussian kernel bằng cách nhân hai vector lại với nhau. Nếu ta có 2 vector N x 1 thì khi (N, 1) x (1, N) ta sẽ thu được ma trận (N, N) tương đương với một 2D kernel

# store the number of values in our Gaussian curve.
ksize = z.get_shape().as_list()[0]

# multiply the two to get a 2d gaussian
z_2d = tf.matmul(tf.reshape(z, [ksize, 1]), tf.reshape(z, [1, ksize]))

# Execute the graph
plt.imshow(z_2d.eval())
plt.colorbar()
plt.show()
gaussian-kernel

gaussian-kernel


Tiếp theo, ta sẽ áp dụng phương pháp nhân tích chập (Convolution) bằng cách sử dụng Gaussian kernel vừa mới tạo lên một ảnh bất kỳ. Đầu tiên, ta cần chuyển ảnh ban đầu về định dạng 1 x H x W x 1 (4D input) và kernel về định dạng H x W x I x O (filters)

Ghi chú:

  • 4D input: [batch, in_height, in_width, in_channels]
  • Filters: [filter_height, filter_width, in_channels, out_channels]
# use tensorflow to reshape matrix
img = mean_img.astype(np.float32)
img_4d = tf.reshape(img, [1, img.shape[0], img.shape[1], 1])
print "Tensorflow image shape:", img_4d.get_shape().as_list()

# Reshape with 4d format: H x W x I x O
z_4d = tf.reshape(z_2d, [ksize, ksize, 1, 1])
print "Tensorflow kernel shape:", z_4d.get_shape().as_list()

# output
Tensorflow image shape: [1, 512, 512, 1]
Tensorflow kernel shape:', [100, 100, 1, 1]

Sau đó, ta có thể thực hiện nhân tích chập cho ảnh với kernel đã tạo

convolved = tf.nn.conv2d(img_4d, z_4d, strides=[1, 1, 1, 1], padding='SAME')
res = convolved.eval()

plt.imshow(np.squeeze(res), cmap='gray')
plt.imshow(res[0, :, :, 0], cmap='gray')
plt.show()
convolution-gaussian

convolution-gaussian

Ta có thể thấy hiệu ứng của Gaussian filter là làm mờ ảnh. Tham số strides cho biết các kernel di chuyển dọc theo ảnh. [1,1,1,1] nghĩa là ta sẽ nhân tích chập cho từng ảnh, từng pixel, từng kênh màu, và cho bất kỳ kernel nào. Padding cho biết cách chúng ta xử lý biên ảnh như thế nào. “SAME” nghĩa là ta muốn kernel phải được apply vào biên ảnh, lúc này biên ảnh sẽ được thêm vào các giá trị zeros để kernel có thể apply. “VALID” nghĩa là không sử dụng padding, lúc này kích thước ảnh sẽ bị thay đổi.

Do matplotlib không xử lý được ảnh 4D nên ta cần chuyển đổi sang số chiều ảnh ban đầu. Ta có thể dùng hàm squeeze() hoặc chỉ định chiều nào của ma trận được sử dụng để hiển thị ảnh.

Ngoài ra, các bạn có thể tự tạo các filter đơn giản để quan sát sự thay đổi cho mỗi kernel.

def convolution_with_filter(img_4d, filter):
    convolved = tf.nn.conv2d(img_4d, filter, strides=[1, 1, 1, 1], padding='SAME')
    res = convolved.eval()

    plt.imshow(np.squeeze(res), cmap='gray')
    plt.imshow(res[0, :, :, 0], cmap='gray')
    plt.show()

# apply sharpen filter
sharpen_filter = np.zeros([3, 3, 1, 1])
sharpen_filter[1, 1, :, :] = 5
sharpen_filter[0, 1, :, :] = -1
sharpen_filter[1, 0, :, :] = -1
sharpen_filter[2, 1, :, :] = -1
sharpen_filter[1, 2, :, :] = -1

convolution_with_filter(img_4d, sharpen_filter)

# apply top sobel filter
top_sobel_filter = np.zeros([3, 3, 1, 1])
top_sobel_filter[0, 0, :, :] = 1
top_sobel_filter[0, 1, :, :] = 2
top_sobel_filter[0, 2, :, :] = 1
top_sobel_filter[2, 0, :, :] = -1
top_sobel_filter[2, 1, :, :] = -2
top_sobel_filter[2, 2, :, :] = -1

convolution_with_filter(img_4d, top_sobel_filter)
sharpen-filter

sharpen-filter

top-sobel-filter

top-sobel-filter

Tham khảo thêm:

Advertisements

Trả lời

Mời bạn điền thông tin vào ô dưới đây hoặc kích vào một biểu tượng để đăng nhập:

WordPress.com Logo

Bạn đang bình luận bằng tài khoản WordPress.com Đăng xuất / Thay đổi )

Twitter picture

Bạn đang bình luận bằng tài khoản Twitter Đăng xuất / Thay đổi )

Facebook photo

Bạn đang bình luận bằng tài khoản Facebook Đăng xuất / Thay đổi )

Google+ photo

Bạn đang bình luận bằng tài khoản Google+ Đăng xuất / Thay đổi )

Connecting to %s