マイブログ リスト

医療言語処理講座

2019年1月24日木曜日

Pyゼミ3.02　手書き画像からDICOM画像へ

手書き数字画像MNISTから医療画像DICOMへ発展させよう


Keyword：MNIST，DICOM，Neural Network，Chainer

前回のPyゼミ3.01にて，手書き数字画像を用いて，ニューラルネットワークの学習を行ないました。
前回のニューラルネットワークをDICOM画像に発展させることを考えてみます。


MNISTがどのようなデータ構造をもった訓練画像なのか調べてみる


対話モードで，Chainerのdatasetsをimportして，訓練データ，テストデータをそれぞれtrainとtestに読み込んでみます。

>>> from chainer import datasets
>>> from chainer.datasets import tuple_dataset
>>>
>>> train, test = datasets.get_mnist(ndim = 3)
次に，訓練データの数をlenメソッドを使って調べてみます。
>>> len(train)
60000

確かに６万件のデータ（手書き数字画像）が読み込まれていることがわかります。
それでは，最初の０番めのデータについて，同じくいくつ格納されているのかlenメソッドで調べてみます。
>>> len(train[0])
2

２つです。
つまり，train[0][0]とtrain[0][1]の２つです。

さらにtrain[0][0]について，調べてみます。
>>> len(train[0][0])
1

１つです。
train[0][0][0]は１つということです。
さらに調べてみると，
>>> len(train[0][0][0])
28

28と表示されました。
これはMNISTの手書き数字の画素サイズです。
さらに調べると。
>>> len(train[0][0][0][0])
28

また，28と表示されました。


少しまとめてみましょう。
最初の60000と表示された。これはデータ数（画像と教師ラベル）を表しています。
つまり，train[0]からtrain[59999]のデータがあることを示しています。

先頭の0番目のデータについて注目してみると
この中には２つの情報があるようです。
ひとつめの情報（tarin[0][0]）を表示してみましょう。

>>>  train[0][0]
array([[[ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
          0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
　　　       ・・・略・・・
          0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
          0.        ,  0.        ,  0.        ]]], dtype=float32)

かなり大きなリストのデータです。多分画像じゃないかって想像がつきます。

ふたつめの情報(train[0][1])を表示してみましょう。
>>> train[0][1]


5

５と表示された。これはラベルデータじゃないかと想像できます。
次に，画像の方をさらに見てみると，len(train[0][0])は１だったので，train[0][0][0]しかない。
これを表示してみると，
>>> train[0][0][0]
array([[ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
　　　　　         ・・・略・・・
         0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ]], dtype=float32)

となる。arrayの中で [ [ で囲まれているので二次元のリストだ，つまり画像が格納されていることがわかります。


画像かどうか表示してみましょう。OpenCVで簡単に見ることができます。
（OpneCVのインストールはPyゼミ1.02参照）

>>> import cv2
>>>
>>> cv2.imshow('image', (train[0][0][0]))
>>> cv2.waitKey(0)





教師ラベルが５であったように，確かに５の手書き数字です。

少し補足するとtrain[0][0][0]は28×28画素のモノクロの画像1枚（１プレーン，１チャンネル）なので，len( train[0][0] ）は1を出力しました。
もし，カラー画像の場合，len( train[0][0] )は3を出力します。
つまり，RGBそれぞれがtrain[0][0][0]，train[0][0][1]，train[0][0][2]に格納されています。

まとめると

訓練画像データは5次元のデータで，ラベルデータは2次元である。
1次元：データの順序番号を表す。(0～データ数-1）
2次元：画像またはラベルを表す。0のとき画像，1のときラベル
3次元：画像プレーン番号（モノクロは0の１つのみ，カラーは0,1または2の３つ）
4次元：画像の行数（縦の画素数，画像の高さ）
5次元：画像の列数（横の画素数，画像の幅）

画像の画素x，yにアクセスするには次のようになります。
　　　　train[データ番号][0][プレーン番号][画像の行（y）][画像の列（ｘ）]
たとえば，10番めの画像のx＝11，y＝12の画素を得るには次のように書きます。
　　　　pix = train[9][0][0][12][11]

ラベルについては下記のようにアクセスできます。
　　　　tarin[データ番号][1]

したがって，画像とそのラベルが与えられたら5次元の訓練データを作成するプログラムを作成すれば良いことになります。


DICOM画像から訓練データを作成する

DICOM画像（このサンプルは個人の画像なのでこの実験以外では使用しないでください）から訓練データを作成する仕様としては以下のように考えます。

  1. DICOM画像の入ったフォルダを指定して，その中の画像を訓練データとする。
  2. ラベル名はファイル名（Brain01, Chest02など）から抽出し，それにラベル番号を自動で（0から）付与する（サンプル画像のファイル名の命名がラベル名となっている）。
  3. DICOMファイル名，ラベル番号とラベル名のCSVファイルを作成する。
  4. 上記，CSVファイルから訓練データを作成する関数を作る


プログラムの概要


constractTrainData関数

次の２つの引数を受け取とります。

  • CSVファイル名（fname） ：DICOMファイル名，ラベル番号とラベル名のCSVファイル名を指定する。
  • リサイズ値（resz）：画像を任意のサイズに拡大または縮小する画素数を指定する。


プログラムははじめに画像とラベルを格納するリストを初期化します。
　　　　image = [ ]
　　　　label = [ ] 

次にfor文にて，CSVファイルから1行づつ読み込みlineに代入し，読み込んだ回数を0から変数n代入します。
変数lineの値はカンマ区切りのデータで，data（リスト）に分割して格納します。
　　　　data[0]にはDICOMファイル名が
　　　　data[1]にはラベル番号が
　　　　data[2]にはラベル名が格納されます。

readDicom2png関数（2.03参照）を用いてウィンドウニングした画像をimgに代入します。
このとき，imgの画素値は0から255の8ビットの値をとります。

もし，リサイズの値reszが0より大きければリサイズしますが，0の場合はリサイズしません。

画像imgは浮動小数点32ビット（np.float32）に変換し，さらに0から255の値を-1.0から1.0の値に変換します。
　　　　img = (img - 128)/128

そして，imageのリストにappendメソッドを用いてimgを追加します。
　　　　image.append( [ img]  )

次にラベル番号は整数32ビット（int32）に変換して，labelのリストに追加します。
　　　　t = np.array( int( data[1] ), dtype = np.int32 )
　　　　label.append( t )

最後に，画像imageとラベルlabelのリストをTupleDatasetメソッドを使って，統合してtrainデータを作成します。
　　　　train = tuple_dataset.TupleDataset( image, label )

TupleDatasetは複数の同じ長さのデータセットをひとつのtuple型のデータセットにまとめるメソッドです。


プログラムを入力して実行してみよう


makeDicomTrain.py


  1. import os, sys
  2. import numpy as np
  3. import cv2
  4. import chainer
  5. from chainer.datasets import tuple_dataset
  6. import pydicom
  7.  
  8. def constractTrainData( fname , resz ):
  9.     image = []
  10.     label = []
  11.     for n, line in enumerate(open(fname, 'r')):
  12.         data = line.split(",")          # 画像ファイル名と正解データに分割
  13.         print(n," File Name:", data[0] ," Label:", data[1])
  14.        
  15.         img = readDicom2png(data[0], "tmp.png")
  16.         if resz > 0:                    #画像のリサイズ       
  17.             img = cv2.resize(img,(resz, resz))       
  18.         img = img.astype(np.float32)   
  19.         img = (img - 128)/128
  20.         image.append([img])
  21.                
  22.         t = np.array(int(data[1]), dtype=np.int32)
  23.         label.append(t)
  24.    
  25.     train = tuple_dataset.TupleDataset(image, label)
  26.    
  27.     return train
  28.    
  29. def readDicom2png(dcmfnm, tmpfnm):
  30.     ds  = pydicom.read_file(dcmfnm)     
  31.     wc  = ds.WindowCenter               
  32.     ww  = ds.WindowWidth                        
  33.     img = ds.pixel_array                        
  34.     #表示画素値の最大と最小を計算する
  35.     max = wc + ww / 2                     
  36.     min = wc - ww / 2     
  37.     #ウインドニング処理
  38.     img = 255 * (img - min)/(max - min)  
  39.     img[img > 255] = 255                 
  40.     img[img < 0  ] = 0                   
  41.     img = img.astype(np.uint8)
  42.     cv2.imwrite(tmpfnm, img)
  43.     return img      
  44.  
  45. if __name__=='__main__':
  46.  
  47.     dirname = "../dcmdir1/"
  48.     fname   = "trainlist.csv"
  49.     resize  = 200
  50.     lblnum  = 0
  51.     lblname = {}
  52.     csvstr  = ""
  53.    
  54.     dlist = os.listdir(dirname)
  55.     print(dlist)
  56.     for nm in dlist:
  57.         bsnm = nm[:-2]              #ファイル名の最後の2文字（数字）を取り除く
  58.         if not bsnm in lblname:
  59.             lblname[bsnm] = lblnum
  60.             lblnum += 1
  61.         csvstr += dirname+nm+","+str(lblname[bsnm])+","+bsnm+"\n"
  62.     
  63.     f = open(fname, 'w')
  64.     f.write(csvstr)
  65.     f.close()
  66.       
  67.     #画像データ，教師データの構築
  68.     train = constractTrainData(fname, resize)
  69.    
  70.     print("\ntrain\n")
  71.     print("画像データ数    :",len(train))
  72.     print("画像/ラベル　　：",len(train[0]))
  73.     print("画像チャンネル数:",len(train[0][0]))
  74.     print("画像の高さ　　　:",len(train[0][0][0]))
  75.     print("画像の幅　　　　:",len(train[0][0][0][0]))
  76.     
  77.     cv2.imshow("image",(train[0][0][0]))
  78.     cv2.waitKey(0)
  79.     print("\nLabel Number:",train[0][1])


実行結果

makeDicomTrain.pyを実行すると，指定したフォルダ（../dcmdir1/）内のファイルを表示します。
8枚のDICOM画像のファイル名とラベル番号が表示されます。
最後に訓練データの構成が表示されます。
画像をリサイズしたので画像の高さと幅が200となっています。


$ python makeDicomTrain.py
['Abdomen01', 'Chest02', 'Chest01', 'Brain01', 'Abdomen03', 'Abdomen02', 'Brain02', 'Chest03']

0  File Name: ../dcmdir1/Abdomen01  Label: 0
1  File Name: ../dcmdir1/Chest02  Label: 1
2  File Name: ../dcmdir1/Chest01  Label: 1
3  File Name: ../dcmdir1/Brain01  Label: 2
4  File Name: ../dcmdir1/Abdomen03  Label: 0
5  File Name: ../dcmdir1/Abdomen02  Label: 0
6  File Name: ../dcmdir1/Brain02  Label: 2
7  File Name: ../dcmdir1/Chest03  Label: 1

train

画像データ数    : 8
画像/ラベル　　： 2
画像チャンネル数: 1
画像の高さ　　　: 200
画像の幅　　　　: 200

Label Number : 0


trainlist.csvの出力

CSVファイルにDICOMファイル名，ラベル番号とラベル名が記録されてます。









確認画像のPNG出力と表示

tmp.pngファイルがreadDicom2png関数から出力されています。





mnistNN.pyにconstractTrainData関数を実装してみる

実装方法は以下の要領で行います。
  1. mnistNN.pyをコピーしてdicomNN.pyに変更します。
  2. makeDicomTrain.pyのconstractTrainData関数，readDicom2png関数をdicomNN.pyにコピーします。
  3. mnistの入力部分は削除し，makeDicomTrain.pyのデータ入力部分を上記２つの関数のあとにコピーします。
  4. os, sys, cv2とpydicomのインポートを追加します。
  5. 今回，testデータはないので，trainデータの２つを使って実験をしてみましょう。

下記プログラムを参照してください。

dicomNN.py


  1. import os, sys
  2. import cv2
  3. import numpy as np
  4. import chainer
  5. from chainer import Function, report, training, utils, Variable
  6. from chainer import datasets, iterators, optimizers, serializers
  7. from chainer import Link, Chain, ChainList
  8. import chainer.functions as F
  9. import chainer.links as L
  10. from chainer.datasets import tuple_dataset
  11. from chainer.training import extensions
  12.  
  13. import pydicom
  14.  
  15. def constractTrainData( fname , resz ):
  16.     image = []
  17.     label = []
  18.     for n, line in enumerate(open(fname, 'r')):
  19.         data = line.split(",")          # 画像ファイル名と正解データに分割
  20.         print(n," File Name:", data[0] ," Label:", data[1])
  21.  
  22.         img = readDicom2png(data[0], "tmp.png")
  23.         if resz > 0:                    #画像のリサイズ       
  24.             img = cv2.resize(img,(resz, resz))       
  25.         img = img.astype(np.float32)   
  26.         img = (img - 128)/128
  27.         image.append([img])
  28.                
  29.         t = np.array(int(data[1]), dtype=np.int32)
  30.         label.append(t)
  31.        
  32.     train = tuple_dataset.TupleDataset(image, label)
  33.        
  34.     return train
  35.    
  36. def readDicom2png(dcmfnm, tmpfnm):
  37.     ds  = pydicom.read_file(dcmfnm)     
  38.     wc  = ds.WindowCenter               
  39.     ww  = ds.WindowWidth                        
  40.     img = ds.pixel_array                        
  41.     #表示画素値の最大と最小を計算する
  42.     max = wc + ww / 2                     
  43.     min = wc - ww / 2     
  44.     #ウインドニング処理
  45.     img = 255 * (img - min)/(max - min)  
  46.     img[img > 255] = 255                 
  47.     img[img < 0  ] = 0                   
  48.     img = img.astype(np.uint8)
  49.     cv2.imwrite(tmpfnm, img)
  50.     return img   
  51.  
  52. #訓練データとテストデータの取得
  53. dirname = "../dcmdir1/"
  54. fname   = "trainlist.csv"
  55. resize  = 28
  56. lblnum  = 0
  57. lblname = {}
  58. csvstr  = ""
  59.  
  60. dlist = os.listdir(dirname)
  61. print(dlist)
  62. for nm in dlist:
  63.     bsnm = nm[:-2]              #最後の2文字を取り除く
  64.     if not bsnm in lblname:
  65.         lblname[bsnm] = lblnum
  66.         lblnum += 1
  67.     csvstr += dirname+nm+","+str(lblname[bsnm])+","+bsnm+"\n"
  68. f = open(fname, 'w')
  69. f.write(csvstr)
  70. f.close()
  71.  
  72. #画像データ，教師データの構築
  73. train  = constractTrainData(fname, resize)
  74. print("\ntrain\n")
  75. print("画像データ数    :",len(train))
  76. print("画像/ラベル　　：",len(train[0]))
  77. print("画像チャンネル数:",len(train[0][0]))
  78. print("画像の高さ　　　:",len(train[0][0][0]))
  79. print("画像の幅　　　　:",len(train[0][0][0][0]))
  80.  
  81. class MyModel(Chain):                   #ネットワークモデルを定義する
  82.     def __init__(self, nLabel):                 #ネットワークの定義
  83.         super(MyModel, self).__init__(
  84.             l1=L.Linear(784,100),       #全結合ネットワーク
  85.             l2=L.Linear(100,100),           
  86.             l3=L.Linear(100,nLabel),
  87.         )
  88.        
  89.     def __call__(self, x):              #順伝播計算
  90.         h1 = F.relu(self.l1(x))         #活性化関数にReLuを定義
  91.         h2 = F.relu(self.l2(h1))
  92.         return self.l3(h2)
  93.        
  94. model = MyModel(lblnum)                 #ネットワークモデルを生成する
  95. model = L.Classifier(model, lossfun=F.softmax_cross_entropy)
  96. optimizer = optimizers.Adam()           #最適化アルゴリズムを設定する
  97. optimizer.setup(model)                  #モデルを最適化アルゴリズムにセットアップ
  98.  
  99. iterator = iterators.SerialIterator(train, 50, shuffle = True) #ミニバッチの設定
  100. updater = training.StandardUpdater(iterator, optimizer) #更新処理の設定
  101. trainer = training.Trainer(updater, (10, 'epoch'))      #訓練（エポック数）の設定
  102. trainer.extend(extensions.ProgressBar())                #進捗表示の設定
  103.  
  104. trainer.run()                           #訓練実行
  105.  
  106. test=train[:2]  #試しに2つのデータについて
  107.  
  108. #訓練結果からtestデータを用いて推定を行う
  109. ok = 0
  110. n  = 0
  111. for i in range(len(test)):              #テストデータ分繰り返す
  112.     x = Variable(np.array([ test[i][0] ], dtype=np.float32)) #入力データ
  113.     t = test[i][1]                      #正解データ
  114.     y = model.predictor(x)              #テストデータから推定結果を計算
  115.     out = F.softmax(y)
  116.     ans = np.argmax(out.data)           #推定結果からラベルを得る
  117.     n += 1
  118.    
  119.     np.set_printoptions(formatter={'float': '{: 0.2f}'.format})
  120.     print('\n#',i)
  121.     print('Infer :',out.data)           #推定結果を表示
  122.     print('Sprvs :','  ','----->'*(int(t)),t)#教師データの表示
  123.     if (ans == t):
  124.         ok += 1   
  125.         print('\tOK',ok/n)
  126.      
  127. print('平均正答率：', ok/len(test))


実行結果

mnistの手書きデータのサイズに合わせるために，DICOM画像を28×28にリサイズしています。
テストは未知のデータですべきですが，最初の２つの画像データでテストをしています。

今回はMNISTからDICOM画像への応用が成功しました。

$ python dicomNN.py
['Abdomen01', 'Chest02', 'Chest01', 'Brain01', 'Abdomen03', 'Abdomen02', 'Brain02', 'Chest03']
0  File Name: ../dcmdir1/Abdomen01  Label: 0
1  File Name: ../dcmdir1/Chest02  Label: 1
2  File Name: ../dcmdir1/Chest01  Label: 1
3  File Name: ../dcmdir1/Brain01  Label: 2
4  File Name: ../dcmdir1/Abdomen03  Label: 0
5  File Name: ../dcmdir1/Abdomen02  Label: 0
6  File Name: ../dcmdir1/Brain02  Label: 2
7  File Name: ../dcmdir1/Chest03  Label: 1

train

画像データ数    : 8
画像/ラベル　　： 2
画像チャンネル数: 1
画像の高さ　　　: 28
画像の幅　　　　: 28

# 0
Infer : [[ 0.67  0.22  0.11]]
Sprvs :     0
 OK 1.0

# 1
Infer : [[ 0.31  0.60  0.09]]
Sprvs :    -----> 1
 OK 1.0
平均正答率： 1.0



おわりに

MNISTはとても有名なデータセットなので，インターネットで調べると実はたくさん情報を得ることができます。
しかし，今回は解析的にデータ構造を追いかけることをしました。
こうした手法は，将来未知のデータセットに出会ったときに役に立つと思います。

今回訓練用データを作成するconstrauctTrain関数はこれからも使います。
また，画像ファイル名とラベル番号をCSVファイルに保存する考え方は，これからも訓練やテストのために使用します。




投稿者 時刻:
ラベル: ,

0 件のコメント:

コメントを投稿

登録: コメントの投稿 (Atom)