Trinity を動かすのにハマる
なかなかアップされないこのブログですが、そもそもターミナルでコード書くとか解析とかしてハマらないとアップされない気もしますw。
さてさて、前段階として年末に本を出すのをお手伝いさせていただきました。
で、縁あって最近、またTrinityをかけていたのですが、コロナウイルスの影響でテレワークでございまして、セキュリティに気をつけようとか言っているわけです。で、ついうっかり
brew upgrade
とかやってしまったわけですね。
そうしたらば、Trinityのアップデートがかかってしまった（ように見えた）上に、Trinity動かすのに必要なJavaの環境が謎の挙動になってしまい、まぁとりあえず動かなくなってしまいました。。。
とりあえず、待てど暮らせどうんともすんとも。-vつけてもどうもこうもなのですが、--debugというオプションを見つけ
brew install --verbose --debug trinity
とやったところ、
java -version
のところで止まっていることが判明。これを手打ちしてもやっぱり動かない。
この間にHomebrewでJavaを入れたり抜いたりしていたので、なんかいろんなJavaが入っていたりw。
参考までに
$ brew config HOMEBREW_VERSION: 2.2.13 ... Java: 1.8.0_252, 1.8.0_73 ← ここにopenjdkの11とか14とかあった。 ...
みたいな感じ。で、動くやつと動かないやつがあるので、動くやつにしたい！でも変えられない。どこにこの設定が書かれているのか、隅々まで探したのがよくわからない。そんな中、
$ /usr/libexec/java_home -V
Matching Java Virtual Machines (2):
1.8.0_252, x86_64: "AdoptOpenJDK 8" /Library/Java/JavaVirtualMachines/adoptopenjdk-8.jdk/Contents/Home
1.8.0_73, x86_64: "Java SE 8" /Library/Java/JavaVirtualMachines/jdk1.8.0_73.jdk/Contents/Home
/Library/Java/JavaVirtualMachines/adoptopenjdk-8.jdk/Contents/Homeというのを発見。
うーん。よくわからんから、/Library/Java/JavaVirtualMachines/の下のadoptopenjdk-8.jdkみたいのを..._bakみたいにとりあえず退避したものの、この下にパスが通っているとここにリストされてしまうっぽいので、全然違うところに移したらば、無事にリストから消えてくれて、brew configしても消えてくれた。ｽﾊﾞﾗｼｽ。
これで、動かないJavaをとりあえず退避。
けっこう危険な方法に見えるのですが、
macOSで古いJDKをアンインストール - Qiita
を見ると、一応こうやってもOKということで敢行。
動くバージョンだけ残して、brew install trinity（ゴタゴタの中で一度 消してみたのである）。そうしたらば無事に入った！よっしゃ。
さてさて改めてTrinityを走らせたのですが、あれ、やっぱり動かん。これまたよく見るとやっぱりJavaのところで止まっている。（java -versionも動かない）。うーん。多分、何か別のが（依存関係で）入ったときに環境が壊れたのだろうと言うことでJavaを改めてインストール。今回は
$ brew cask install adoptopenjdk8
これでとりあえず動いたのでメデタシメデタシ。
（この後、FASTQのヘッダ行がとか、ホームディレクトリの~が展開されなかったりとかいろいろあるのだが、それはまた別の記事にしましょうか。）
Biopythonを使ってみる
俺は都会のマンションのベランダで、自分勝手に植物を育てる、ベランダーだ。
季節は春。俺は久しぶりに、ブログを更新することにした。
（植物男子ベランダー・一部改）
はてダがサービス終了ということで、はてブに移ってきました。
そもそもここにいろいろ書くのも最近はサボっていたので、まぁ、これを機会にまたちょこちょこ書くようにしていければとは思うんですけどねー（言うのは簡単）。
さてさて、ふとしたことでライフサイエンス系公共データをいじることとなり（ってまぁ、それが本業といえば本業なんですが）、いじるというか知識抽出なんですが、今までBioPerlでちゃらちゃらっとやっていたのを、せっかくだからBiopythonを使おうかと思い立ったわけです。
もともとpythonは10数年前から知ってはいてread-onlyだったです。今まで何回か、書きかけたこともあったですけどね。まぁ自分のころはperlが盛んにこの分野で使われていましたからねぇ。流行の最末端の自分ですが、世の中、だいぶpythonが隆盛してきたので（ってこともないですが）こういう機会でもないとpython使わねーしなというのもあり。
とまぁ、前置きが長いですが、
qiita.com
とか
biopython.org
を参考に書いてみたりして。
やったことはGenBankファイルのparse
元ファイル
LOCUS EU919294 1220 bp DNA linear INV 26-JUL-2016
DEFINITION Lampides boeticus voucher RMBR:102827 cytochrome c oxidase subunit
I (COI) gene, partial cds; mitochondrial.
ACCESSION EU919294
VERSION EU919294.1
KEYWORDS .
SOURCE mitochondrion Lampides boeticus (pea blue)
ORGANISM Lampides boeticus
Eukaryota; Metazoa; Ecdysozoa; Arthropoda; Hexapoda; Insecta;
Pterygota; Neoptera; Holometabola; Lepidoptera; Glossata; Ditrysia;
Papilionoidea; Lycaenidae; Polyommatinae; Lampides.
REFERENCE 1 (bases 1 to 1220)
AUTHORS Lohman,D.J., Peggie,D., Pierce,N.E. and Meier,R.
TITLE Phylogeography and genetic diversity of a widespread Old World
butterfly, Lampides boeticus (Lepidoptera: Lycaenidae)
JOURNAL BMC Evol. Biol. 8, 301 (2008)
PUBMED 18973689
REMARK Publication Status: Online-Only
REFERENCE 2 (bases 1 to 1220)
AUTHORS Lohman,D.J., Peggie,D., Pierce,N.E. and Meier,R.
TITLE Direct Submission
JOURNAL Submitted (23-JUL-2008) Department of Biological Sciences, National
University of Singapore, 14 Science Drive 4, Singapore 117543,
Republic of Singapore
FEATURES Location/Qualifiers
source 1..1220
/organism="Lampides boeticus"
/organelle="mitochondrion"
/mol_type="genomic DNA"
...書いたもの
#!/usr/env/python import sys import os import inspect from Bio import SeqIO gb_file = sys.argv[1] # input gb file path for record in SeqIO.parse(gb_file, 'genbank'): id = record.id name = record.name desc = record.description xrefs = record.dbxrefs # list annots = record.annotations # dic features = record.features # list seq = record.seq print('ID:\t' + id) print('name:\t' + name) print('desc:\t' + desc) print('seq:\t' + seq) # annotation ## source source = annots['source'] organism = annots['organism'] taxonomy = annots['taxonomy'] # list print('source:\t' + source) print('organism:\t' + organism) print('taxonomy:') for t in taxonomy: print('\t' + t) ## reference for ref in annots['references']: pmid = ref.pubmed_id medline = ref.medline_id journal = ref.journal print('pmid:\t' + pmid) print('medline:\t' + medline) print('journal:\t' + journal) # features for feat in features: type = feat.type print('type:\t' + type) for k in feat.qualifiers.keys(): print('subtype:\t' + k) vals = feat.qualifiers[k] for v in vals: print('\t' + v)
で、出てきたもの
ID: EU919294.1
name: EU919294
desc: Lampides boeticus voucher RMBR:102827 cytochrome c oxidase subunit I (COI) gene, partial cds; mitochondrial
seq: AACATTATATTTTATTTTTGGAATTTGAGCAGGAATATTGGGAACATCATTAAGAATTTTAATTCGAATAGAATTAGGAACTCCAGGATCTTTAATTGGTGATGATCAAATTTATAATACTATTGTAACAGCTCATGCTTTTATTATAATTTTTTTTATAGTTATACCTATTATAATTGGTGGATTTGGAAACTGATTAATCC
...
TTTCAGTTGATTAGCTACTATTTATGGAACTCAAATTAACTATAGACCTTCAATATTATGAAGTTTAGGTTTCATTTTTTTATTTACAGTAGGCGGATTAACAGGAGTAATTTTAGCTAACTCTTCAATTGATATTACTTTACATGATACTTATTATGTTGTAGCACATTTCCATTATGTTTTATCTATAGGAGCTGTATTTGCAATTTTTGGAGGATTTATTCATTGATATCCTTTATTTACTGGATTAATAATAAATCCTTATCTATTAAAAATTCAATTTATTATTATATTTATTGGAGTTAATT
source: mitochondrion Lampides boeticus (pea blue)
organism: Lampides boeticus
taxonomy:
Eukaryota
Metazoa
Ecdysozoa
Arthropoda
Hexapoda
Insecta
Pterygota
Neoptera
Holometabola
Lepidoptera
Glossata
Ditrysia
Papilionoidea
Lycaenidae
Polyommatinae
Lampides
pmid: 18973689
medline:
journal: BMC Evol. Biol. 8, 301 (2008)
...
type: source
subtype: organism
Lampides boeticus
subtype: organelle
mitochondrion
subtype: mol_type
genomic DNA
...
きったねーコードというのは承知しているんですけどね。
まぁ、ボチボチ手直ししつつ、さらにこれを使って先に進めれれば。
PlatanusをMacで使う
縁あって、ゲノム解析とか発現解析をしています。
de novoでゲノムデータをassembleするというので、Platanus ( http://platanus.bio.titech.ac.jp/ ) というのを使ってみることとしました。現在のバージョンは1.2.4のようで（去年に出たようだ）。サイトでは1.2.1ってのも配られていて、こちらはMacのバイナリもあるようですが、最新版の1.2.4はLinux 64 bit binary (precompiled)かソースで、という状況。
そこいらのマシンで自分のMacProが一番性能がよさげなので、これで動かすべくソースを落としてみました。
一応、trimmerも使おうかと思うので、まずはこちらから。ほどきまして、makeしますと、
nakazato@grouper:~/Downloads/Platanus_trim_v1.0.7$ make
g++ -o common.o -c common.cpp -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -DRUN_MODE=\"\"
clang: warning: argument unused during compilation: '-s'
In file included from common.cpp:2:
./common.h:13:10: fatal error: 'omp.h' file not found
#include
^
1 error generated.
make: *** [common.o] Error 1
はい、怒られました。もろもろ調べたところ、Macのg++はclangの方で（イミフ。でもよし）OpenMPは使わないようになっているとか。そこで（brewで brew install gcc49してからか?）Makefileを
CXX = g++
↓
CXX = g++-4.9
しますと
nakazato@grouper:~/Downloads/Platanus_trim_v1.0.7$ make
g++-4.9 -o common.o -c common.cpp -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -DRUN_MODE=\"\"
g++-4.9 -o trim.o -c trim.cpp -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -DRUN_MODE=\"\"
g++-4.9 -o main.o -c main.cpp -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -DRUN_MODE=\"\"
g++-4.9 -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -o platanus_trim common.o trim.o main.o
ld: warning: option -s is obsolete and being ignored
g++-4.9 -o main.oo -c main.cpp -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -DRUN_MODE=\"internal\"
g++-4.9 -s -std=c++0x -O3 -funroll-loops -fomit-frame-pointer -fopenmp -o platanus_internal_trim common.o trim.o main.oo
ld: warning: option -s is obsolete and being ignored
なにやら怒られましたが、platanus_internal_trimとplatanus_trimができたのでよしとしましょう。（ちゃんと動きました）。
はてさて、今度はPlatanus本体ですが
nakazato@grouper:~/Downloads/Platanus_v1.2.4$ make
g++ -o main.o -c main.cpp -std=c++0x -O3 -funroll-loops -Wall -fopenmp -finline-limit-50000 -lm -Dnullptr=0
clang: error: unknown argument: '-finline-limit-50000'
clang: warning: -lm: 'linker' input unused
make: *** [main.o] Error 1
はい、やっぱり怒られました。で、trimmerのと同様に直したわけですが
nakazato@grouper:~/Downloads/Platanus_v1.2.4$ make
g++ -o main.o -c main.cpp -std=c++0x -O3 -funroll-loops -Wall -fopenmp -finline-limit-50000 -lm -Dnullptr=0
clang: error: unknown argument: '-finline-limit-50000'
clang: warning: -lm: 'linker' input unused
make: *** [main.o] Error 1
nakazato@grouper:~/Downloads/Platanus_v1.2.4$ cp Makefile.org Makefile.arrange
nakazato@grouper:~/Downloads/Platanus_v1.2.4$ vi Makefile.arrange
nakazato@grouper:~/Downloads/Platanus_v1.2.4$ make clean
rm -f platanus main.o assemble.o scaffold.o scaffoldGraph.o gapClose.o common.o baseCommand.o seqlib.o mapper.o gapCloseOLC.o
nakazato@grouper:~/Downloads/Platanus_v1.2.4$ make -f Makefile.arrange
g++-4.9 -o main.o -c main.cpp -std=c++0x -O3 -funroll-loops -Wall -fopenmp -finline-limit-50000 -lm -Dnullptr=0
In file included from assemble.h:25:0,
from main.cpp:22:
counter.h:455:263: error: redeclaration of 'void Counter::countKmerOrWriteTemporary(bool&, const typename KMER::keyType&, DoubleHash *, FILE*, omp_lock_t*, const KMER&, unsigned int)' may not have default arguments [-fpermissive]
inline void Counter::countKmerOrWriteTemporary(bool &loopFlag, const typename KMER::keyType &key, DoubleHash tmpOccurrenceTable, FILE *unmappedFP, omp_lock_t lock, const KMER &kmer, const unsigned iterateTimes=32)
^
counter.h:684:190: error: redeclaration of 'void Counter::countKmerOrWriteTemporary(bool&, const typename KMER::keyType&, FILE*, omp_lock_t*, const KMER&, unsigned int)' may not have default arguments [-fpermissive]
inline void Counter::countKmerOrWriteTemporary(bool &loopFlag, const typename KMER::keyType &key, FILE *unmappedFP, omp_lock_t lock[], const KMER &kmer, const unsigned iterateTimes=32)
^
In file included from assemble.h:26:0,
from main.cpp:22:
graph.h:1764:80: error: redeclaration of 'void BruijnGraph::cutBranchIterative(long long unsigned int)' may not have default arguments [-fpermissive]
void BruijnGraph::cutBranchIterative(const unsigned long long numThread=1)
^
make: *** [main.o] Error 1
やっぱり怒られました。で、うねうね調べまして、むりくり先に進めるみたいなオプションをつけるとうまくいくとかで、上記のに加えて
CXXFLAGS = -std=c++0x -O3 -funroll-loops -Wall -fopenmp -finline-limit-50000 -lm -Dnullptr=0
↓
CXXFLAGS = -std=c++0x -O3 -funroll-loops -Wall -fopenmp -finline-limit-50000 -lm -Dnullptr=0 -fpermissive
としますと、まぁやっぱりくどくど怒られるんですが、無事にplatanusの実行ファイルができまして、ちゃんととりあえず、動き始めたんで一安心みたいなところですかね。
週末を越えて動かしているんですが、ちょっとファイルが大きかったのか、スレッド数とかメモリが小さかったのか、"WARNING: Sorry, memory exceeds specified value!!"と怒られつつ、まだ計算しております。。。（そのあたりはまた後日?
MacでのJavaのバージョンアップ
Javaって、まぁ、セキュリティのごにょごにょもあって最新版にせーとか何とか言われたりするわけでありますが、今般、8でないと動かないソフトを入れるにあたり、6から8にバージョンアップしようとしたわけであります。
そのソフト、ターミナルで入れるので、ならJavaのバージョンもターミナルで確認しますわな。
nakazato@gardeneel:~$ java -version
java version "1.6.X_XXX"
Java(TM) SE Runtime Environment (build 1.6.X_XXX-bXX)
Java HotSpot(TM) 64-Bit Server VM (build XX.XXX-bXX, mixed mode)
今はもう8になっているのでちとごまかしてますが。で、ググると、システム環境設定の中にJavaの項目があって、そこからアップデートできるというわけです。確かにそうなっていて、アップデートを促され、ポチッとやることで6から8にアップデートされたかのようになってめでたしめでたし、とは行かなくてですね。上のコマンドをたたいてもバージョンは6のままでございまして。なんかインストールするところが変わったのか、一部しかインストールされていないのか。
で、またごにょごにょ調べまして、インストールされているところにシンボリックリンクをはるかどうかの寸前まで行ったわけですが、結局、システム環境でインストールされるのはJREの方でございまして、コマンドラインの方で見てるのはJDKである、ということがわかりまして、普通にJavaのダウンロードページJava SE - Downloads | Oracle Technology Network | Oracleに行って、JDKのを落としてきてインストールしたらば
nakazato@gardeneel:~$ java -version
java version "1.8.0_101"
Java(TM) SE Runtime Environment (build 1.8.0_101-b13)
Java HotSpot(TM) 64-Bit Server VM (build 25.101-b13, mixed mode)
無事にアップデートできましたとさ。
サンタ
毎度おなじみ BioHackathon で山形は あつみ温泉に来ている。
外国からも交えて100人くらいか。今年は大所帯。
夜にNCBIとうちのと話していたらば、そのうちのが今日は何をやっていたかと訊かれていて、リストを作っていたと返事をしたところ、チェックは2回したか、お前はサンタか、ということを真顔で言われたわけである。実はこれは以下をふまえたジョークである。今度どこかで使いたい。
Santa Claus is Coming to Town（サンタが町にやってきた）
He's making a list
And checking it twice
Gonna find out Who's naughty and nice
Santa Claus is coming to town
昆虫学会＋応動昆＠大阪
来てます。
今回はポスター発表にしてます。内容は、まぁ、いつもどおり。本当はもっとモノに根ざした発表をしたいんですが。
せっかくなので、公共NGS中の昆虫データについて調べてみました。
↑ ポスターの図より。数字はプロジェクト数。1つのプロジェクトで複数の生物種が出てきた時は複数にカウントされているから、厳密にはこの表現はおかしいのかとも思う。
ハエが多いのはDrosophilaとかその近縁とかカが含まれているから、と思われ。
チョウ目はカイコとタバコスズメガか?、セミ目にはウンカやアブラムシが入っているから。甲虫が少ないですね、という話もあったが、生活環が長いから研究には不向きなのだ、という別の指摘が入った。
生物種は何種類か、という質問を受けたのでその場で調べてみたところ、1256種類だった。
上位の方。
770 7227 Drosophila melanogaster
83 7165 Anopheles gambiae
51 7460 Apis mellifera
40 7091 Bombyx mori
37 7159 Aedes aegypti
36 7240 Drosophila simulans
24 76194 Papilio polytes
24 7160 Aedes albopictus
23 7237 Drosophila pseudoobscura
21 7245 Drosophila yakuba
21 7029 Acyrthosiphon pisum
19 7173 Anopheles arabiensis
19 7038 Bemisia tabaci
15 7244 Drosophila virilis
15 7070 Tribolium castaneum
なんか、上の図と合わないんだが、なんで合わないかは、そのtaxonomy ID決め打ちか、下も見ているかの違いだったかと記憶している。（グラフのDrosophila melanogasterはDrosophila melanogaster なんちゃら/なんちゃらグループみたいのも含んでいる、ということ）
だいぶdryをやっていると忘れそうな感覚だが、この学会に来ている人はこういう学名を見てだいたいあれね、というのがわかる。（まぁ、自分も有名どころとか日ごろ図鑑を見つめている蝶関係とかならだいたいわかるが）。それは、wetの研究者が遺伝子発現で発現量が変動した遺伝子上位20とか見せられて、あー、あの遺伝子ね、とわかるのと同じことである。
oral聞いていると、なんちゃらアザミウマとかにわかに姿が浮かばない昆虫が乱発してニッチェ感が満載である。それに対して、私のやっているなんちゃらハマキガでは、とか細かいどうしで質問していてむしろこの人たちの研究対象をこうやってリスト化したいw。
BioProjectとかBioSampleを眺めてみた（続き）
普段の職場を飛び出してなので、（いまだに書類をこさえたり調整してメールしたりもありますが）コード書いたりするのがなんか久しぶりな気がします。ま、こうやってわいわいやった方が、ひとりでこもってやるよりシナジー効果が出ますわな。
昨日の続きとして、接頭辞のをちゃんと見たりもしてみた。（この後ろに数字が来る）
27780 SAMD
17 SAME
480372 SAMEA
3628607 SAMN
ついでにBioProjectでも同じことをしてみた。
288 PRJDA
3167 PRJDB
499 PRJEA
9964 PRJEB
143180 PRJNA
それから、仕事柄、文献とのリンクも見たいのだが、実際に文献で引用されている正例がとれずに、まだそこは野ざらしである。
本業の方だが、懸案の自動でAmazonのデータを更新の部分で、データベースをこちらでこさえて、転送するところまでは済。きちんとそれがcronでまわるのか（自分で打つ分には実行できるけど）、向こうできちんとロードできるかはまだなので、それはある程度 時を待ちたい。
あ、あとは15日をすぎたので、定例の月イチの統計情報の更新もdone
なかなか目に見えるものができないってのはアレですなぁ。
