京橋のバイオインフォマティシャンの日常

データ分析、コマンドライン、プログラミングについての技術資料・自己アップデート・悩み事などをまとめています。最近、ディープラーニング関連のR言語の資料をまとめるべく注力してます。

grepコマンドで文字列処理をやってみた件【その2】検索語のヒット数カウントとか検索語の前後文字の抽出とか色々

f:id:skume:20220119233952p:plain

はじめに

grepコマンドによる文字列処理をやってみた」の続編です。

以前扱えていなかった内容をやっていきたいです。

grepコマンドの基本については過去の記事を参照のこと。

skume.hatenablog.com

まずは、サンプルデータをダウンロードする

$ svn export  https://github.com/kumeS/Blog/trunk/grep_practice_02

$ cd ./grep_practice_02

# test.txtを使っていく。
$ cat test.txt

#abcde
#ABCDE
#
#abcdefghijklmnopqrstuvwxyz
#ABCDEFGHIJKLMNOPQRSTUVWXYZ
#
#12345abcde67890fghij
#12345abcde67890fghij

内容

検索語に一致した文字列のみを抜き出して出力する

grep -oを使うと、検索後と同じものが、ヒットした文字列として出力されます。

#'abc'で検索すると、5個ヒット
$ grep -o 'abc' test.txt
#abc
#abc
#abc
#abc
#abc

# -o を入れずに検索すると、'abc'を含む行数の4個でヒット
$ grep 'abc' test.txt
#abcdeabcde
#abcdefghijklmnopqrstuvwxyz
#12345abcde67890fghij
#12345abcde67890fghij

#'fghij'で検索すると、3個ヒット
$ grep -o 'fghij' test.txt
#fghij
#fghij
#fghij

大文字・小文字を区別せず、検索語に一致した文字列のみを抜き出して出力する

grep -o -iを使うと、大文字・小文字が区別されずに、ヒットした文字列が出力されます。

#'abc' or 'ABC'で検索され、8個ヒット
$ grep -o -i 'abc' test.txt
#abc
#abc
#ABC
#ABC
#abc
#ABC
#abc
#abc

# -i のみで検索すると、6個ヒット
$ grep -i 'abc' test.txt
#abcdeabcde
#ABCDEABCDE
#abcdefghijklmnopqrstuvwxyz
#ABCDEFGHIJKLMNOPQRSTUVWXYZ
#12345abcde67890fghij
#12345abcde67890fghij

#'fghij' or 'FGHIJ' で検索され、4個ヒット
$ grep -o -i 'fghij' test.txt
#fghij
#FGHIJ
#fghij
#fghij

検索語に一致した文字列の数をカウントする【部分一致検索】

grep -ogrep -cを使うと、検索語に一致した文字列数が出力されます。

このとき注意することとしてコマンドを分けずに、grep -o -cにすると、行数カウントになります。

#'abc'で検索して、'abc'をカウントすると、5個ヒット
$ grep -o 'abc' test.txt | grep -c 'abc'
#5

# -o -c にすると、4個ヒットで、行数カウントになる
$ grep -o -c 'abc' test.txt
#4

#'fghij'で検索すると、3個ヒット
$ grep -o 'fghij' test.txt | grep -c 'fghij'
#3

検索語に一致した箇所の前後の行を出力する

grep -Cを使うと、前後の行(行数指定)を一緒に出力できます。

#前後1行を出力する
$ grep -C 1 'hijklmn' test.txt
#
#abcdefghijklmnopqrstuvwxyz
#ABCDEFGHIJKLMNOPQRSTUVWXYZ

#前1行、後3行を出力する
$ grep  -B 1  -A 3 'hijklmn' test.txt
#123456
#abcdefghijklmnopqrstuvwxyz
#ABCDEFGHIJKLMNOPQRSTUVWXYZ
#
#12345abcde67890fghij

grep -B -Aを使うと、-Bで検索語前の行数、-Aで検索語後の行数を指定できます。

検索語に一致した箇所の前後の文字列を出力する

上記で、行レベルで出力する方法を紹介したが、 正規表現. (ピリオド)を使うと、前後の文字列も出力できます。

#検索語を抜き出して表示
$ grep -o 'efghijk' test.txt
#efghijk

#検索語の前4文字も含めて出力
$ grep -o '....efghijk' test.txt
#abcdefghijk

#検索語の後4文字も含めて出力
$ grep -o 'efghijk....' test.txt
#efghijklmno

#検索語の前4文字、後4文字も含めて出力
$ grep -o '....efghijk....' test.txt
#abcdefghijklmno

また、検索語を「変数」で指定することもできる。

TERMS=$"efghijk" ; grep -o "....$TERMS...." test.txt
#OR
TERMS=$'efghijk' ; grep -o "....$TERMS...." test.txt

#abcdefghijklmno

やってみると、'"の組み合わせが重要らしいです。

;は、複数文を繋げて書く文法である。

検索語を行番号付きで出力する

grep -nを使うと、行番号付きで出力されます。

#'abc'行を行番号付きで出力
$ grep -o -n 'abc' test.txt
#1:abc
#abc
#4:abc
#7:abc
#8:abc

$ grep -n 'abc' test.txt
#1:abcdeabcde
#4:abcdefghijklmnopqrstuvwxyz
#7:12345abcde67890fghij
#8:12345abcde67890fghij

#'fghij'行を行番号付きで出力
$ grep -o -n 'fghij' test.txt
#4:fghij
#7:fghij
#8:fghij

検索語の行番号のみを出力する

行番号のみを出力するには、sedとの組み合わせで行います。

#'abc'行を行番号のみ出力
$ grep -n 'abc' test.txt | sed -e 's/:.*//g'
#1
#4
#7
#8

#'fghij'行を行番号のみ出力
$ grep -n 'fghij' test.txt | sed -e 's/:.*//g'
#4
#7
#8

grepのエラーメッセージを消す

そのままではあまり使わないけど、プログラムを組むときには役に立つかもです。

#フォルダに実行するとエラーが出る
$ grep 'abcd' test
grep: test: Is a directory

$ grep -s 'abcd' test
#OR
$ grep --no-messages 'abcd' test
#エラーメッセージの非表示

まとめ

grepコマンドは結構役に立つので、一通り覚えておくと、 いろいろなテキストの前処理とかに使えそうです。

補足

正規表現のまとめ

正規表現 (Regular Expression)は、文字列のパターン・マッチングなどの表記法として使用する。

また、メタキャラクタとは、文字列のパターンを表す特殊な記号のことを指す。

メタキャラクタ そのメタキャラクタが意味すること
. 任意の一文字
* 前にある文字の0回以上の繰り返し
+ 前にある文字の1回以上の繰り返し
? 前にある文字の0回あるいは1回の意味
^ 指定した文字が文頭
$ 指定した文字が文末
[ ] 指定した複数の文字の中のいずれか
{ } 前にある文字の指定した回数の繰り返し
\w 英数文字またはアンダーバー。[a-zA-Z0-9_]と同等
\W \w 以外の文字
\d 数値文字。[0-9]と同等
\D 数値以外の文字。[^0-9]と同等
\s 空白文字
\S 空白以外の文字
\b ワードの先頭、あるいは末尾にマッチ
\B ワード内の文字列にマッチ
\A 文字列の最初にマッチ
\Z 文字列の最後にマッチ
\t タブ
\n 改行
\r キャリッジリターン

参考

www.wakuwakubank.com

webplus8.com

www.gadgety.net

grepコマンドで文字列処理をやってみた件【その1】ファイル内のテキストに対する処理とか

f:id:skume:20220119224919p:plain

はじめに

grepコマンドは、ある特定の文字列を含むものを取り出す・検索するときに使用します。 大別して、ファイル内のテキストに対してgrepを行うか、ファイル名に対してgrepを行うかでやり方が変わってきます。

やってることは単純だが、いろいろと応用できるので、今回まとめてみました。

練習用のファイルのダウンロード

練習用のファイル・フォルダを以下のコマンドでダウンロードできるようにしています。

#ダウンロード
$ svn export  https://github.com/kumeS/Blog/trunk/grep_practice

#作業フォルダに移動
$ cd ./grep_practice

内容

ファイル内のテキストに対する処理

ある特定の文字列を含む行を表示する

#ファイル内の表示
$ cat test_01.txt

#'abc'を含む行を出力する
$ grep 'abc' test_01.txt
#OR
$ grep "abc" test_01.txt
#OR
$ cat test_01.txt | grep 'abc'

ある特定の文字列を含まない行を表示する

#'abc'を含まない行を出力する
$ grep -v 'abc' test_01.txt 
#OR
$ cat test_01.txt | grep -v "abc"

ある特定の文字列を含まない行を別ファイルに出力する

#'abc'を含まない行を「test_01_v.txt」として保存する
$ grep -v 'abc' test_01.txt > test_01_v.txt

#test_01_v.txtの表示
$ cat test_01_v.txt

オプション

-v: 文字列を含まない、あるいは、一致しないものを検索する

空白行を削除する

#空白行を含まない行の表示
$ cat test_02.txt | grep -v '^\s*$'

#別ファイルで保存する場合
$ cat test_02.txt | grep -v '^\s*$' > test_02_r.txt

オプション

^: 行頭

\s: 空白文字

^\s*: 行頭から0回以上繰り返しの空白

ファイル内の行数をカウントする

#ファイル内全体の行数のカウント
$ cat test_01.txt | grep '' -c

#grepを使わない方法
$ cat test_01.txt | wc -l

#'abc'を含む行数のみのカウント
$ cat test_01.txt | grep 'abc' -c

オプション

-c: 行数とかのカウント

複数条件で、文字列の検索をする場合

grepでのAND検索

$ cat test_02.txt

#AND検索
#'AAA'と'ABB'を同行に含む行の表示
$ grep 'AAA' test_02.txt | grep 'ABB'

grepでのOR検索

$ cat test_02.txt

#OR検索(-e表記)
#'AAA'あるいは'BBB'を含む行の表示
$ grep -e 'AAA' -e 'BBB' test_02.txt

#OR検索(正規表現、\| を使う)
#'AAA'あるいは'BBB'を含む行の表示
$ grep 'AAA\|BBB' test_02.txt

オプション

-e: 検索パターンを指定する

-w: 単語全体でパターンと一致するものを検索する

-x: 行全体がパターンと一致するものを検索する

-i: 大文字と小文字を区別しない

-n: 検索結果に行番号を表示する

ファイル名に対する処理

ある特定の文字列をファイル名に含むファイルを検索する

#ディレクトリ内ファイルの表示
$ ls

#.txtファイルのみを表示する
$ ls | grep '.txt'
#.txtファイルのカウント
$ ls | grep '.txt' | wc -l

#ファイル内の行、ワード、バイトの各カウントの表示
$ ls | grep 'test_01.txt' | wc
#OR
$ wc test_01.txt

オプション

-l: ファイル名を検索対象とする

サブディレクトリも含めて検索する

#'.pdf'を含むファイル数のカウント
$ ls | grep '.pdf' | wc -l

#サブディレクトリを含めたファイル表示
$ du -a

#サブディレクトリを含めて、'.pdf'を含むファイルの表示
$ du -a | grep '.pdf'

#サブディレクトリを含めて、ファイル数のカウント
$ du -a | grep '.pdf' | wc -l

特定のファイルを検索して、消去する

#ディレクトリ内の'.pdf'を含むファイルを消去する
$ ls | grep '.pdf' | xargs rm -rf

#サブディレクトリも含めて削除する場合
$ du -a | grep '.pdf' | xargs rm -rf

ここで、 xargs は、前のコマンドの実行結果を引数にして、次のコマンド rm に渡すことを意味します。

補足

ターミナル・ショートカット
カーソルの移動
Ctrl + b 後退・一語後退
Ctrl + f 前進・一語前進
Ctrl + a 行頭へ移動
Ctrl + e 行末へ移動
削除
Ctrl + w 単語1個文削除
Ctrl + k 行末まで削除
Ctrl + u 行頭まで削除
Ctrl + d カーソル上の1文字削除
Ctrl + h カーソルの後方を1文字削除
履歴
一つ前のコマンド履歴
次のコマンド履歴
history (コマンド) 履歴
その他
Ctrl + c 実行中のコマンドを強制終了
Ctrl + z 実行中のコマンドを一時停止
Ctrl + d 終了、Logout
Ctrl + l 画面をクリアする
Ctrl + t 一つ前とカーソルの文字の入替
Ctrl + m, Ctrl + j, Ctrl + o Enter

cd : ディレクトリ移動

f:id:skume:20200607071106p:plain:w400

cat : ファイル内表示・結合・作成
#テキストファイル(test.txt)内の表示
$ cat test.txt

#行番号を付加して、ファイル表示
$ cat -n test.txt

#複数ファイルから、結合ファイルを作成(行方向に結合)
$ cat test1.txt test2.txt > test3.txt

#空ファイルの作成(Ctrl + d で終了すること)
$ cat > test4.txt

#すべての履歴表示、「~ (チルダ)」はHomeディレクトリ
$ cat ~/.bash_history
#OR
$ history

head / tail : ファイルの先頭や末尾を表示
#ファイルの先頭から、10行を表示
$ head -n 10 test.txt

#ファイルの末尾から、10行を表示
$ tail -n 10 test.txt

覚えておくと良いターミナルコマンド
#システムを終了する
$ shutdown -f now
              
#再起動する
$ reboot

#すべてのユーザーに実行権限を与える
$ chmod a+x test.command

#Rootユーザーでコマンドを実行する
$ sudo ...

その他のターミナルコマンド
#日時を表示するコマンド。
$ date

#2012年のカレンダーを表示する
$ cal 2020

#アクティブなジョブの表示
$ jobs

関連記事

skume.net

htmlwidgets for R のShowcaseにあるパッケージがCodePenでブログ表示できるかを調べた件

f:id:skume:20220119223900p:plain

はじめに

htmlwidgets for R パッケージは、Rでインタラクティブな図が作成できる王道的なパッケージであり、 それを使った色々な依存パッケージが開発されています。

www.htmlwidgets.org

今回、htmlwidgetsのshowcaseにある12パッケージをHatena Blog内の図表として表示するにあたって、(1)Html 出力の可否 、(2)CodePen対応できるか?、さらには(3)ブログ表示できるか?を調べてみました。

CodePenとは、ユーザーが作成したHTML、CSS、JavaScriptのコードスニペットのテスト、 および披露をするオンラインツールの1つです。

とりあえず、結果

パッケージ html output CodePen はてなブログ表示
Leaflet 473KB
Dygraphs 574KB
Plotly 3.5MB X X
rbokeh 955KB X X
Highcharter 995KB X X
visNetwork 1.1MB X X
networkD3 263KB
d3heatmap 229KB
DataTables 207KB
threejs 977KB X X
rglwidget X X X
DiagrammeR 2.5MB X X
MetricsGraphics 368KB

表示結果

Leaflet

Rスクリプト

if(!require("leaflet")){install.packages("leaflet")}; library(leaflet)
if(!require("magrittr")){install.packages("magrittr")}; library(magrittr)

m <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addMarkers(lng=174.768, lat=-36.852, popup="The birthplace of R") 
m

CodePenでの出力結果

See the Pen leaflet by skume (@kumes) on CodePen.

Dygraphs

Rスクリプト

if(!require("dygraphs")){install.packages("dygraphs")}; library(dygraphs)
if(!require("magrittr")){install.packages("magrittr")}; library(magrittr)

dygraph(nhtemp, main = "New Haven Temperatures") %>% 
  dyRangeSelector(dateWindow = c("1920-01-01", "1960-01-01"))

CodePenでの出力結果

See the Pen dygraph by skume (@kumes) on CodePen.

networkD3

Rスクリプト

if(!require("networkD3")){install.packages("networkD3")};
library(networkD3)
if(!require("magrittr")){install.packages("magrittr")}; library(magrittr)

data(MisLinks, MisNodes)
forceNetwork(Links = MisLinks, Nodes = MisNodes, Source = "source",
             Target = "target", Value = "value", NodeID = "name",
             Group = "group", opacity = 0.4)

CodePenでの出力結果(Run Pen クリック必要)

See the Pen networkD3 by skume (@kumes) on CodePen.

d3heatmap

Rスクリプト

if(!require("d3heatmap")){install.packages("d3heatmap")}; library(d3heatmap)

d3heatmap(mtcars, scale="column", colors="Blues")

CodePenでの出力結果

See the Pen d3heatmap by skume (@kumes) on CodePen.

DataTables

Rスクリプト

if(!require("DT")){install.packages("DT")}; library(DT)

datatable(iris, options = list(pageLength = 5))

CodePenでの出力結果

See the Pen datatable by skume (@kumes) on CodePen.

MetricsGraphics

Rスクリプト

if(!require("metricsgraphics")){install.packages("metricsgraphics")}; library(metricsgraphics)
if(!require("magrittr")){install.packages("magrittr")}; library(magrittr)

mjs_plot(mtcars, x=wt, y=mpg) %>%
  mjs_point(color_accessor=carb, size_accessor=carb) %>%
  mjs_labs(x="Weight of Car", y="Miles per Gallon")

CodePenでの出力結果

See the Pen metricsgraphics by skume (@kumes) on CodePen.

まとめ

CodePenの1MB制限で、Plotly、rbokeh、Highcharter、visNetwork、threejs、DiagrammeRがCodePenで保存できない。

うすうす気付いてたけど、Plotlyのhtml出力はやや重過ぎるので、Webコード共有には不向き。

visNetworkの結果は貼付できないけど、networkD3の結果は貼付できるという不思議な結果になった。。

CodePenでの保存ができると、だいたいブログでも表示できるっぽい。

MacOSXターミナルでのbz2形式の圧縮・解凍についてまとめてみた件

f:id:skume:20220119223039p:plain

.bz2について

bzip2では、 圧縮効率を良くするために、ブロックソート法などを用いています。 gzipやzipといったデータ圧縮に比べて、より高い圧縮率を示します。 また、bz2単独では、アーカイブ機能はありません。

Macでは、bzip2コマンドでbzip2圧縮、 bunzip2コマンドでbzip2解凍ができます。

それでは、Macのターミナルを起動して、実行してみます。

bzip2コマンドによるbzip2圧縮

bzip2コマンドは、 1つのファイルあるいはファイルパスに対して実行します。

bzip2 -z [File]

オプション

  • -z: 圧縮を行う(デフォルトなので省略化)

また、フォルダの圧縮はできない

bzip2圧縮を行い、任意のファイル名で保存する場合

bzip2 -z [File] > xxxx.bz2

bunzip2コマンドによるbzip2解凍

「XXXX.bz2」を解凍する場合

bunzip2 XXXX.bz2  

#OR

bzip2 -d XXXX.bz2 > [File]

オプション

  • -d: 伸張を行う

また、元ファイルを残して、解凍する場合には、 -kオプションを付けます。

bunzip2 -k XXXX.bz2  

.tar.bz2について

bz2圧縮ではアーカイブ機能がないので、 tarを組み合わせてアーカイブ化します。

.tar.bz2とは、tarでアーカイブ化して、 bzip2圧縮することを意味します。

tarコマンドによるtar.bz2圧縮

「XXXX.tar.bz2」として、tar.bz2圧縮アーカイブ化する場合

tar -jcvf XXXX.tar.bz2 [File1] [File2] [File3]

#OR

tar -acvf XXXX.tar.bz2 [File1] [File2] [File3]

オプション

  • j: bzip2の意味

  • cvf: tarアーカイブ化

  • [File1] [File2] [File3]: 圧縮したいファイル(フォルダも可)

  • a: 拡張子による圧縮方式の自動判定

また、オプションの「-」は、tarの場合、省略できます。

tarコマンドで、フォルダごとの圧縮ができる

tarコマンドによるtar.bz2解凍

「XXXX.tar.bz2」をbz2解凍して、tarアーカイブを戻す場合

tar -jxvf xxxx.tar.bz2

#OR

tar -axvf xxxx.tar.bz2

オプション

  • j: bzip2の意味

  • xvf: tarアーカイブの解凍

  • a: 拡張子による解凍方式の自動判定

また、オプションの「-」は、tarの場合、省略できます。

R/Keras/TensorFlowでやる『ディープラーニング(Deep Learning)』のすゝめ【その2】教師なしニューラルネットワーク Autoencoder with 2D CNNの実装、そして色ムラ・ノイズ除去(Denoising)をやってみた件

f:id:skume:20220119221936p:plain

はじめに

「R/Keras/TensorFlowでやるディープラーニングのすゝめ」の連載2回目です(2022年1月19日アップデート版)。

【1】では、ベクトルデータに対する Autoencoderを取り上げましたが、 今回は、 2D Convolutional Neural Network (CNN: 畳み込みニューラルネットワーク) を使ったAutoencoderの実装について概説します。

skume.net

下記で登場する、CNNMaxプーリングUpサンプリング活性化関数 などについては、すでに多くの分かりやすい説明記事がありますので、そちらを参照のこと。

qiita.com

deepage.net

qiita.com

https://kansiho.hatenablog.com/entry/2018/05/07/%E5%9B%B3%E3%81%A7%E8%AA%AD%E3%82%80%E3%80%8C%E5%AE%9F%E8%A3%85%E3%83%87%E3%82%A3%E3%83%BC%E3%83%97%E3%83%A9%E3%83%BC%E3%83%8B%E3%83%B3%E3%82%B0%E3%80%8D%E3%81%9D%E3%81%AE%EF%BC%91kansiho.hatenablog.com

https://www.researchgate.net/publication/332057917_Transfer_learning_between_crop_types_for_semantic_segmentation_of_crops_versus_weeds_in_precision_agriculturewww.researchgate.net

今回の内容

R/Kerasのセットアップ

RStudioを起動して、Kerasパッケージをロードします。そして、Pythonを選択します。

library(keras)
#"/usr/local/bin/python"を選択する
reticulate::use_python("/usr/local/bin/python", required =T)

R/Kerasのセットアップの詳細については、前回の記事とかを参照のこと。

skume.hatenablog.com

www.slideshare.net

MNISTデータの準備

今回も、手書きアラビア文字であるMNISTデータセットを使用します。

MNISTデータセットをロードして、リスト型からアレイ型に変換します。

同データセットのトレーニングデータは6万枚の画像りますがが、 それを全て使用すると、残念ながら、CPU計算ではなかなか収束しない・・・

そこで、トレーニングデータの内で、1000枚の画像をランダムで ピックアップして、新たにデータセットを作ります。

また、2D CNNレイヤーを使用する場合、 inputデータは、4次元アレイ型(= 4次元テンソル)に変換しておくのがポイントです。

#MNISTのダウンロード
Data <- dataset_mnist()
str(Data)

#アレイ型に変換
xtrain <- Data$train$x
ytrain <- Data$train$y

#256階調 => 0〜1 に変換する
xtrain <- xtrain/255

#1000枚の画像をランダム・ピックアップ
Sam <- sample(1:dim(xtrain)[1], 1000, replace = F)
xtrain <- xtrain[Sam,,]
ytrain <- ytrain[Sam]

#4次元アレイ型に変換
x_train <- array_reshape(xtrain, dim=c(dim(xtrain)[1], 28, 28, 1))

dim(x_train)
#[1] 1000   28   28    1

このとき、4次元アレイ(= 4次元テンソル)は、 1次元目を画像番号、2と3次元目を2D画像データ、 4次元目をチャネル数( グレイスケールの場合、1 、RGBの場合、3 )で与えます。

Autoencoder with 2D CNN のモデル構築

2D CNNを使用する場合、実際の入力データは、4次元アレイの1次元目を除いた、2〜4次元目の3次元アレイとなります。 このアレイは、1と2次元目を2D画像データ、3次元目をチャネル数で与えます。 そのため、layer_inputshapeは、c(28, 28, 1) と3要素ベクトル(ピクセル数とチャネル数の組み合わせ)で記述します。

ちょっとした関数の説明

layer_conv_2d関数で、2次元畳み込み層(例えば、画像上の空間畳み込み)を与える。

layer_conv_2dfiltersは、出力空間の次元(すなわち、畳み込みの出力フィルタの数)を整数値で与える。

layer_conv_2dkernel_size(sizeで与える)は、2次元畳み込みウィンドウの幅と高さを、整数値で与える。

layer_conv_2dpaddingでは、"valid" か "same" を指定する。"same"の場合、同じアレイを返す。

layer_max_pooling_2d関数は、空間データの最大プーリング処理を行う層を与える。

layer_max_pooling_2dpool_sizeは、ダウンスケールする係数 (垂直, 水平)を整数値で与える。例えば、c(2, 2) を指定すると、垂直・水平の空間次元で入力が半分になる。

layer_upsampling_2d関数は、2D入力用のアップサンプリング層を与える。

layer_upsampling_2dsizeでは、行と列のアップサンプリング係数を整数値で与える。

これらのKeras内の関数を使って、autoencoderモデルを構築してみる。

#インプット層の設定
input <- layer_input(shape = c(28, 28, 1))
kernel_size <- c(3,3)
filters <- 32

#アウトプット層の設定
output = input %>%
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_max_pooling_2d(pool_size=c(2,2), padding="same") %>% 
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_max_pooling_2d(pool_size=c(4,4), padding="same") %>% 
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_upsampling_2d(size=c(4,4))  %>% 
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="valid")  %>% 
  layer_upsampling_2d(size=c(2,2))  %>% 
  layer_conv_2d(filters=1, kernel_size=kernel_size, activation="sigmoid", padding="same")
  
Autoencoder2DCNN <- keras_model(input, output)
summary(Autoencoder2DCNN)

#Model: "model"
#_____________________________________________________
#Layer (type)           Output Shape          Param # 
#=====================================================
#input_1 (InputLayer)   [(None, 28, 28, 1)]   0       
#_____________________________________________________
#conv2d (Conv2D)        (None, 28, 28, 32)    320     
#_____________________________________________________
#max_pooling2d (MaxPool (None, 14, 14, 32)    0       
#_____________________________________________________
#conv2d_1 (Conv2D)      (None, 14, 14, 32)    9248    
#_____________________________________________________
#max_pooling2d_1 (MaxPo (None, 4, 4, 32)      0       
#_____________________________________________________
#conv2d_2 (Conv2D)      (None, 4, 4, 32)      9248    
#_____________________________________________________
#up_sampling2d (UpSampl (None, 16, 16, 32)    0       
#_____________________________________________________
#conv2d_3 (Conv2D)      (None, 14, 14, 32)    9248    
#_____________________________________________________
#up_sampling2d_1 (UpSam (None, 28, 28, 32)    0       
#_____________________________________________________
#conv2d_4 (Conv2D)      (None, 28, 28, 1)     289     
#=====================================================
#Total params: 28,353
#Trainable params: 28,353
#Non-trainable params: 0
#_____________________________________________________

このとき、input_1 (InputLayer) [(None, 28, 28, 1)]conv2d_4 (Conv2D) (None, 28, 28, 1)のサイズを一致させる必要があります。

このモデルでは、28x28 => 14x14 => 4x4 => 16x16 => 14x14 => 28x28 とピクセルサイズがダウンスケール・アップスケールしています。

このとき、4x4以下までピクセル数を下げると、精度が悪くなる。また、16x16 => 14x14は、layer_conv_2dpaddingvalidとすることでピクセル数を微調整しています。

あと、MNISTデータは、28x28ピクセルであるが、モデル構築の観点からは、使用する2D画像のピクセル数は、2N x 2N ピクセル(例えば、25 = 32、27 = 128)が望ましいです。

DLモデルの出力

source("https://gist.githubusercontent.com/kumeS/41fed511efb45bd55d468d4968b0f157/raw/0f64b83700ac578d0c39abd420da5373d4317083/DL_plot_modi_v1.1.R")

Autoencoder2DCNN %>% plot_model_modi(width=1, height=1.25)

plot_model_modiの出力結果)

https://kumes.github.io/Blog/ConvolutionalAutoencoder/1_model.html

TensorFlowでもやってみると

tf <- reticulate::import(module = "tensorflow")
py_plot_model <- tf$keras$utils$plot_model
py_plot_model(Autoencoder2DCNN, to_file='Autoencoder2DCNN_tf.png', 
              show_shapes=T, show_layer_names=T, 
              expand_nested=T, dpi=100)

py_plot_modelの出力結果)

https://kumes.github.io/Blog/ConvolutionalAutoencoder/1_Autoencoder2DCNN_tf.png

モデルが確認できれば、compile & fit を実行する。

# compile
Autoencoder2DCNN %>% 
  compile(optimizer="rmsprop", loss="mean_squared_error")

# Fit
Autoencoder2DCNN %>% 
  fit(x_train, x_train, epochs=200, batch_size=32, 
      shuffle = T, verbose=1)

(トレーニングの結果)

https://kumes.github.io/Blog/ConvolutionalAutoencoder/2_Autoencoder.html

結果の評価

次に、このモデルによるTrainデータセットの変換結果を実際の画像と比較してみます。

まずは、アレイ型を4次元から3次元に変換します。

library(EBImage)
pred_imgs <- Autoencoder2DCNN %>% predict(x_train)
pred_imgsR <- array_reshape(pred_imgs, dim=c(dim(pred_imgs)[1], 28, 28))
dim(pred_imgsR)

par(mfrow=c(3,2))
for (i in 1:6) {
  m <- sample(1:dim(xtrain)[1], 1, replace = F)
  display(combine(t(xtrain[m,,]), t(pred_imgsR[m,,])), 
          method="raster", nx=2, all=TRUE, spacing = 0.01, margin = 2)
}

左側が入力画像(オリジナル画像)、右側が予測出力画像の結果を示します。

https://kumes.github.io/Blog/ConvolutionalAutoencoder/3_Pred.png

さらに、Autoencoderの圧縮特徴量層を取り出して、t分布型確率的近傍埋め込み法 (t-SNE; T-distributed Stochastic Neighbor Embedding) で次元圧縮を行い、2D 表示をしてみます。

#圧縮特徴量層までの変換モデルを取り出して、予測する
summary(Autoencoder2DCNN)
intermediate_layer <- keras_model(inputs = Autoencoder2DCNN$input,
                                  outputs = get_layer(Autoencoder2DCNN, "conv2d_2")$output)
summary(intermediate_layer)
intermediate_output <- predict(intermediate_layer, x_train)

#圧縮特徴量層をベクトルに変換
str(intermediate_output)
intermediate_output_re <- array_reshape(intermediate_output, dim=c(dim(intermediate_output)[1], dim(intermediate_output)[2]*dim(intermediate_output)[3]*dim(intermediate_output)[4]))
xy <- data.frame(ytrain, intermediate_output_re)

#t-SNEによる圧縮特徴量ベクトルの次元圧縮
#install.packages("Rtsne")
library(Rtsne)
xy.tSNE <- Rtsne(as.matrix(xy[,-1]), check_duplicates = FALSE, verbose=TRUE, 
                 dims = 2, perplexity = 25, theta = 0.5, max_iter = 5000)
plot(xy.tSNE$Y, cex=0.5, pch=19, col=rainbow(10)[xy[,1]+1],
     xlab="t-SNE PC1", ylab="t-SNE PC2")

#結果表示
par(mfrow=c(1,1), mai=c(0.75,0.75,0.2,0.2), mgp = c(2,1,0))
xy1 <- data.frame(ytrain, xy.tSNE$Y)
plot(xy1[,2:3], cex=0.5, pch=19, col="white",
     xlab="t-SNE PC1", ylab="t-SNE PC2")
a <- range(xy1[,2][is.finite(xy1[,2])])
b <- range(xy1[,3][is.finite(xy1[,3])])
a1 <- diff(a)*0.0125
b1 <- diff(b)*0.0125

for(n in 1:nrow(xy1)){
  #n <-1
  v <- col2rgb(rainbow(10)[xy1[n,1] + 1]) / 255
  img = channel(xtrain[n,,], 'rgb')
  img[,,1] <- img[,,1]*v[1]
  img[,,2] <- img[,,2]*v[2]
  img[,,3] <- img[,,3]*v[3]
  ff <- t(as.raster(img))
  ff[ff == "#000000"] <- "#00000000"
  rasterImage(ff, xy1[n,2]-a1, xy1[n,3]-b1, 
              xy1[n,2]+a1, xy1[n,3]+b1)
}

https://kumes.github.io/Blog/ConvolutionalAutoencoder/4_tSNE.png

CNNモデルの方*1が、良い感じに分離境界がでています。

色ムラに対するDenoising Autoencoder

このセクションでは、Autoencoderを用いた、色ムラのノイズ除去を行ってみます。

ゴマシオノイズを消す事例は多くあったので、色ムラ戻す事例をやってみます。

まずは、続けてやるとよくないので、RStudioを初期化してみます。

.rs.restartR()
library(keras)
reticulate::use_python("/usr/local/bin/python", required =T)

色ムラがある手書き文字の生成

ランダムな方向に線形の色ムラを施して、4次元アレイ型に変換します。

str(xtrain)
xtrain_noise <- xtrain
abc <- sample(1:4, 1000, replace=T)

for(n in 1:dim(xtrain_noise)[1]){
  m <- abc[n]
  ab <- matrix(1, 28, 28)*seq(from = 0, to = 1, by = 1/27)
  if(m == 1){
    xtrain_noise[n,,] <- xtrain[n,,]*ab
  }
  if(m == 2){
    xtrain_noise[n,,] <- xtrain[n,,]*t(ab)
  }
  if(m == 3){
    xtrain_noise[n,,] <- xtrain[n,,]*ab[28:1,]
  }
  if(m == 4){
    xtrain_noise[n,,] <- xtrain[n,,]*t(ab[28:1,])
  }
}

#4次元アレイ型に変換
x_train <- array_reshape(xtrain, dim=c(dim(xtrain)[1], 28, 28, 1))
x_train_noise <- array_reshape(xtrain_noise, dim=c(dim(xtrain_noise)[1], 28, 28, 1))

そして、色ムラ有無の2D画像を並べて確認してみます。

library(EBImage)

par(mfrow=c(3,2))
for (i in 1:6) {
  m <- sample(1:dim(x_train)[1], 1, replace = F)
  EBImage::display(combine(t(x_train_noise[m,,,]), t(x_train[m,,,])), 
          method="raster", nx=2, all=TRUE, spacing = 0.01, margin = 2)
}

https://kumes.github.io/Blog/ConvolutionalAutoencoder/5_Denoise.png

左側が色ムラ画像(入力画像)、右側がオリジナル画像(出力画像)を示します。

Autoencoder for denoising モデルの構築

上記モデルとほぼ同じで、filtersは64に変更して構築しました。

input <- layer_input(shape = c(28, 28, 1))
filters <- 64
kernel_size <- c(3,3)

output = input %>%
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_max_pooling_2d(pool_size=c(2,2), padding="same") %>% 
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_max_pooling_2d(pool_size=c(2,2), padding="same") %>% 
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_upsampling_2d(size=c(2,2))  %>% 
  layer_conv_2d(filters=filters, kernel_size=kernel_size, activation="relu", padding="same") %>% 
  layer_upsampling_2d(size=c(2,2))  %>% 
  layer_conv_2d(filters=1, kernel_size=kernel_size, activation="sigmoid", padding="same")

AutoencoderDenoising <- keras_model(input, output)
summary(AutoencoderDenoising)

モデル表示を行ってみると

source("https://gist.githubusercontent.com/kumeS/41fed511efb45bd55d468d4968b0f157/raw/0f64b83700ac578d0c39abd420da5373d4317083/DL_plot_modi_v1.1.R")
AutoencoderDenoising %>% plot_model_modi(width=1, height=1.25)

plot_model_modiの出力結果)

https://kumes.github.io/Blog/ConvolutionalAutoencoder/6_AutoencoderDenoising_model.html

tf <- reticulate::import(module = "tensorflow")
py_plot_model <- tf$keras$utils$plot_model
py_plot_model(AutoencoderDenoising, to_file='AutoencoderDenoising_tf.png', 
              show_shapes=T, show_layer_names=T, 
              expand_nested=T, dpi=100)

py_plot_modelの出力結果)

https://kumes.github.io/Blog/ConvolutionalAutoencoder/6_AutoencoderDenoising_model_tf.png

Compileでは、最適化アルゴリズムをadam、損失関数をbinary_crossentropyにした*2

# compile
AutoencoderDenoising %>% 
  compile(optimizer="adam", loss="binary_crossentropy")

# Fit
AutoencoderDenoising %>% 
  fit(x_train_noise, x_train, epochs=200, batch_size=32, 
      shuffle = T, verbose=1)

https://kumes.github.io/Blog/ConvolutionalAutoencoder/7_DenoiseRes.html

このモデルによる色補正の結果を見てみます。

library(EBImage)
pred_imgs <- AutoencoderDenoising %>% predict(x_train_noise)
pred_imgsR <- array_reshape(pred_imgs, dim=c(dim(pred_imgs)[1], 28, 28))
dim(pred_imgsR)

par(mfrow=c(3,2))
for (i in 1:6) {
  m <- sample(1:dim(xtrain_noise)[1], 1, replace = F)
  EBImage::display(combine(t(xtrain_noise[m,,]), t(pred_imgsR[m,,]), t(xtrain[m,,])), 
          method="raster", nx=3, all=TRUE, spacing = 0.01, margin = 2)
}

https://kumes.github.io/Blog/ConvolutionalAutoencoder/8_Denoise_results.png

左側が色ムラ画像(入力画像)、真ん中がAutoencoderによるノイズ除去変換後(予測画像)、右側がオリジナル画像(出力画像)です。

トレーニング画像内の内挿であるが、ほぼ色ムラが除去されています。

まとめ

Autoencoderでノイズ除去ができることが分かりました。

また、元画像を入力にして、ノイズ画像を出力にすると、ノイズ付加するAutoencoderモデルとなります。

ただ、実際には、タスクとなるちょうど良いノイズがはいった対応画像を手に入れるのがやや大変そうに思いました。

R/Kerasを用いたDeep Learningの推薦図書

参考文献

www.datatechnotes.com

Autoencoders with Keras, TensorFlow, and Deep Learning - PyImageSearch

elix-tech.github.io

*1:t-SNEが良いのかもだけど

*2:損失関数が、mean_squared_error ではうまく収束しなかった