Julia を使って Fasta ファイルの中に、空のシーケンスのエントリーがあるか調べる

# 概要

あるソフトウェアを使うときに fasta ファイルが入力になるのだが

```
>dummy_1.1 dummy data 1
MSS

>empty_1.1 empty data 1

>dummy_2.1 dummy data 2
AAA
```

こんなファイルがあったとしたら、

```
>empty_1.1 empty data 1

```

このような行がある場合に処理を進めたくない。

このように、空のシーケンスのエントリーがあるかを調べたい。
できたら、高速に。

# 手法

ワンライナーでかいていたのだが、たくさんあるしなぁ、高速化とか、考えるとなぁとか思っていたら、以下の記事がちょうどよくでてきました。

* [正規表現からJuliaのコードを生成する - Qiita](https://qiita.com/bicycle1885/items/b9ae737bb1655a3b7252) 

Julia のコードから書くことはできないが、
正規表現ならなんとかかけるかなぁというかんじがするので、試してみることに。

今後のことを考えると、できるだけ速く動くようにできる余地を残しておきたかったので、Julia を使ってみることに。

## 以下作業記録

[BioJulia/Automa.jl: A julia code generator for regular expressions](https://github.com/BioJulia/Automa.jl) 

の中のサンプルにfastaについのものがあったので、

[Automa.jl/fasta.jl at master · BioJulia/Automa.jl](https://github.com/BioJulia/Automa.jl/blob/master/example/fasta.jl) 

### まずは、現状だとどうなるのか確認

まずは上記の `fasta.jl` でダミーデータを使うとどうなるのかを調べた

```
$ bin/julia -qL check-fasta-entry-has-no-sequence.jl 
julia> records
3-element Array{FASTARecord,1}:
 FASTARecord("dummy_1.1", "dummy data 1", UInt8[0x4d, 0x53, 0x53])
 FASTARecord("empty_1.1", "empty data 1", UInt8[])                
 FASTARecord("dummy_2.1", "dummy data 2", UInt8[0x41, 0x41, 0x41])
```

空の部分、`empty_1.1` も、parseされている。

目標はこの `empty_1.1` があれば、まずはエラーとしたい。
本当は、このようなエントリーがあるところを全部だしたいが、まずは１つでもあったらエラーが目標

### 簡単になおせるとおもった

ソースをながめていたら、シーケンスが以下のように組み立てられていたので

```julia
    letters     = re"[A-Za-z*-]*"
    sequence    = re.cat(letters, re.rep(newline * letters))
```

letters のところの正規表現の最後の `*` を `+` にすれば一気に解決じゃないか！
とおもって試してみました。

```diff
--- a/check-fasta-entry-has-no-sequence.jl
+++ b/check-fasta-entry-has-no-sequence.jl
@@ -12,7 +12,7 @@ fasta_machine = (function ()
     newline     = re"\r?" * lf
     identifier  = re"[!-~]*"
     description = re"[!-~][ -~]*"
-    letters     = re"[A-Za-z*-]*"
+    letters     = re"[A-Za-z*-]+"
     sequence    = re.cat(letters, re.rep(newline * letters))
     record      = re.cat('>', identifier, re.opt(re" " * description), newline, sequence)
     fasta       = re.rep(record)
```

#### 結果

期待通り、エラーになっている！


```console
$ bin/julia -qL check-fasta-entry-has-no-sequence.jl 
ERROR: LoadError: failed to parse on line 3
parse_fasta(::String) at /home/vagrant/julia-0d7248e2ff/check-fasta-entry-has-no-sequence.jl:73
include_from_node1(::String) at ./loading.jl:576
include(::String) at ./sysimg.jl:14
(::Base.Distributed.##135#136{Base.#include,Tuple{String},Array{Any,1}})() at ./distributed/remotecall.jl:314
run_work_thunk(::Base.Distributed.##135#136{Base.#include,Tuple{String},Array{Any,1}}, ::Bool) at ./distributed/process_messages.jl:56
#remotecall_fetch#140(::Array{Any,1}, ::Function, ::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
remotecall_fetch(::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
#remotecall_fetch#144(::Array{Any,1}, ::Function, ::Function, ::Int64, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:367
remotecall_fetch(::Function, ::Int64, ::String, ::Vararg{Any,N} where N) at ./distributed/remotecall.jl:367
(::Base.##507#509{Base.JLOptions})() at ./task.jl:335
while loading /home/vagrant/julia-0d7248e2ff/check-fasta-entry-has-no-sequence.jl, in expression starting on line 81
#remotecall_fetch#140(::Array{Any,1}, ::Function, ::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:340
remotecall_fetch(::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
#remotecall_fetch#144(::Array{Any,1}, ::Function, ::Function, ::Int64, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:367
remotecall_fetch(::Function, ::Int64, ::String, ::Vararg{Any,N} where N) at ./distributed/remotecall.jl:367
(::Base.##507#509{Base.JLOptions})() at ./task.jl:335
Stacktrace:
 [1] sync_end() at ./task.jl:287
 [2] macro expansion at ./task.jl:303 [inlined]
 [3] process_options(::Base.JLOptions) at ./client.jl:279
 [4] _start() at ./client.jl:371
```

## あとは、正常系で試すだけ

`empty_1.1` を外して実行すれば、問題がないことを確認して、終了とおもい、実行してみた。

```
$ bin/julia -qL check-fasta-entry-has-no-sequence.jl 
ERROR: LoadError: failed to parse on line 9
parse_fasta(::String) at /home/vagrant/julia-0d7248e2ff/check-fasta-entry-has-no-sequence.jl:73
include_from_node1(::String) at ./loading.jl:576
include(::String) at ./sysimg.jl:14
(::Base.Distributed.##135#136{Base.#include,Tuple{String},Array{Any,1}})() at ./distributed/remotecall.jl:314
run_work_thunk(::Base.Distributed.##135#136{Base.#include,Tuple{String},Array{Any,1}}, ::Bool) at ./distributed/process_messages.jl:56
#remotecall_fetch#140(::Array{Any,1}, ::Function, ::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
remotecall_fetch(::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
#remotecall_fetch#144(::Array{Any,1}, ::Function, ::Function, ::Int64, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:367
remotecall_fetch(::Function, ::Int64, ::String, ::Vararg{Any,N} where N) at ./distributed/remotecall.jl:367
(::Base.##507#509{Base.JLOptions})() at ./task.jl:335
while loading /home/vagrant/julia-0d7248e2ff/check-fasta-entry-has-no-sequence.jl, in expression starting on line 81
#remotecall_fetch#140(::Array{Any,1}, ::Function, ::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:340
remotecall_fetch(::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
#remotecall_fetch#144(::Array{Any,1}, ::Function, ::Function, ::Int64, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:367
remotecall_fetch(::Function, ::Int64, ::String, ::Vararg{Any,N} where N) at ./distributed/remotecall.jl:367
(::Base.##507#509{Base.JLOptions})() at ./task.jl:335
Stacktrace:
 [1] sync_end() at ./task.jl:287
 [2] macro expansion at ./task.jl:303 [inlined]
 [3] process_options(::Base.JLOptions) at ./client.jl:279
 [4] _start() at ./client.jl:371
```

エラーがでている。。。



## 中の人にアドバイスをもらう。

幸運にも、中の人にアドバイスをもらえる機会があったので、きいてみたところ。

そもそも、 `empty_1.1` があったときのエラーをよくみると、

```
ERROR: LoadError: failed to parse on line 3
```

このようになっている。
実際空としてエラーにしてほしいのは、５行目。

```
1 >dummy_1.1 dummy data 1
2 MSS
3 
4 >empty_1.1 empty data 1
5 
6 >dummy_2.1 dummy data 2
7 AAA
```

3行目は最初の空行部分なので、ここでは引っかからないようにすればよいというアドバイスをいただき、最初の空行にあたらないで、エントリーの部分だけがあたるように、
さらにコードの改修案を教えていただいた。

### 修正方法

```diff
--- a/check-fasta-entry-has-no-sequence.jl
+++ b/check-fasta-entry-has-no-sequence.jl
@@ -12,8 +12,9 @@ fasta_machine = (function ()
     newline     = re"\r?" * lf
     identifier  = re"[!-~]*"
     description = re"[!-~][ -~]*"
-    letters     = re"[A-Za-z*-]*"
-    sequence    = re.cat(letters, re.rep(newline * letters))
+    letters1     = re"[A-Za-z*-]+"
+    letters2     = re"[A-Za-z*-]*"
+    sequence    = re.cat(letters1, re.rep(newline * letters2))
     record      = re.cat('>', identifier, re.opt(re" " * description), newline, sequence)
     fasta       = re.rep(record)
 
@@ -23,8 +24,10 @@ fasta_machine = (function ()
     identifier.actions[:exit]   = [:identifier]
     description.actions[:enter] = [:mark]
     description.actions[:exit]  = [:description]
-    letters.actions[:enter]     = [:mark]
-    letters.actions[:exit]      = [:letters]
+    letters1.actions[:enter]     = [:mark]
+    letters1.actions[:exit]      = [:letters]
+    letters2.actions[:enter]     = [:mark]
+    letters2.actions[:exit]      = [:letters]
     record.actions[:exit]       = [:record]
 
     # Finally, compile the final FASTA pattern into a state machine.
```

### 実行結果

```console
$ bin/julia -qL check-fasta-entry-has-no-sequence.jl 
ERROR: LoadError: failed to parse on line 5
parse_fasta(::String) at /home/vagrant/julia-0d7248e2ff/check-fasta-entry-has-no-sequence.jl:76
include_from_node1(::String) at ./loading.jl:576
include(::String) at ./sysimg.jl:14
(::Base.Distributed.##135#136{Base.#include,Tuple{String},Array{Any,1}})() at ./distributed/remotecall.jl:314
run_work_thunk(::Base.Distributed.##135#136{Base.#include,Tuple{String},Array{Any,1}}, ::Bool) at ./distributed/process_messages.jl:56
#remotecall_fetch#140(::Array{Any,1}, ::Function, ::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
remotecall_fetch(::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
#remotecall_fetch#144(::Array{Any,1}, ::Function, ::Function, ::Int64, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:367
remotecall_fetch(::Function, ::Int64, ::String, ::Vararg{Any,N} where N) at ./distributed/remotecall.jl:367
(::Base.##507#509{Base.JLOptions})() at ./task.jl:335
while loading /home/vagrant/julia-0d7248e2ff/check-fasta-entry-has-no-sequence.jl, in expression starting on line 84
#remotecall_fetch#140(::Array{Any,1}, ::Function, ::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:340
remotecall_fetch(::Function, ::Base.Distributed.LocalProcess, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:339
#remotecall_fetch#144(::Array{Any,1}, ::Function, ::Function, ::Int64, ::String, ::Vararg{String,N} where N) at ./distributed/remotecall.jl:367
remotecall_fetch(::Function, ::Int64, ::String, ::Vararg{Any,N} where N) at ./distributed/remotecall.jl:367
(::Base.##507#509{Base.JLOptions})() at ./task.jl:335
Stacktrace:
 [1] sync_end() at ./task.jl:287
 [2] macro expansion at ./task.jl:303 [inlined]
 [3] process_options(::Base.JLOptions) at ./client.jl:279
 [4] _start() at ./client.jl:371
```

重要な部分は２行目のところにでていて

```console
$ bin/julia -qL check-fasta-entry-has-no-sequence.jl 
ERROR: LoadError: failed to parse on line 5
```

これをみると、エラーの発生行が３行目から５行目に変化していることがわかる。

## 問題のないケースでは、エラーにならないことの確認

最後に、`empty_1.1` を取り除いて、どこにもエラーがないときは、
正しく終了することの確認

与えるデータは、以下のように、`empty_1.1` をはぶいたものとする。

```
>dummy_1.1 dummy data 1
MSS

>dummy_2.1 dummy data 2
AAA
```
### 実行結果

問題がないことが確認できる

```
$ bin/julia -qL check-fasta-entry-has-no-sequence.jl 
julia> records
2-element Array{FASTARecord,1}:
 FASTARecord("dummy_1.1", "dummy data 1", UInt8[0x4d, 0x53, 0x53])
 FASTARecord("dummy_2.1", "dummy data 2", UInt8[0x41, 0x41, 0x41])
julia>
```

## さらなるアドバイス

* いらないアクションをはずしてしまう

など、いくつかアドバイスをいただいたので、それを反映させてみたい。


## ワンライナー版

もともとのワンライナーはこんなかんじ。
`a.fasta` についてチェックしている。

```
pcregrep -M "^>.*\n\n" a.fasta
```