kuroの覚え書き

96の個人的覚え書き

カテゴリープロット、ボックスプロット、バイオリンプロット

最近Nature系の論文は投稿時に平均値+エラーバーだけでなく、全データ点をプロットしたグラフを要求するようになった。
投稿サイトでもそういうグラフをオンラインで作ってくれるサイトが紹介されているのでそれを使えばいいのだが、みんなそのサイトを使ってグラフを書き直しているためか、どれもこれもみんな似たようなグラフになっていてキモい。
ということで、自前で描画できるアプリを作ろうということでいつものようにPythonでプログラムを書く。
と言ってもmatplotlibとseabornで描かせるだけなので自作というのもおこがましいのだが。

f:id:k-kuro:20200804233920p:plain

サクサクっと完成。これをexcelで作ったグラフに重ねてしまえばOK。

ついでに、最近良く見るようになったバイオリンプロットも描けるようにしてみた。
f:id:k-kuro:20200804234038p:plain

Flaskのプログラムにいつものように継ぎ足しているので、ライブラリのインポートやformとかの定義は省略しているが、ソースは以下のような感じ。

@app.route('/cat_graph', methods=['GET', 'POST'])
def cat_graph_index():
    user = g.user.name
    form=StatForm()
    dir = "./user/" + user + "/graph"
    if os.path.exists(dir):
        shutil.rmtree(dir)
    if not os.path.exists(dir):
        os.makedirs(dir)
    input = dir + "/input.txt"
    output = dir + "/output.txt"


    if request.method == 'POST':
        if 'file' not in request.files:
            flash('No file part')
            return redirect(request.url)
        file = request.files['file']
        if file.filename == '':
            if form.samp1.data:
                samp = form.samp1.data.strip()
                f = open(input, 'w')
                f.write(samp)
                f.close()
            else:
                flash('No data input')
                return redirect(request.url)
        if file and allowed_file(file.filename):
            file.save(input)

    if os.path.exists(input):
        if form.fig_x.data:
            fig_x=form.fig_x.data
        if form.fig_y.data:
            fig_y=form.fig_y.data

        if form.data_f.data=='csv':
            df = pd.read_csv(input)
        if form.data_f.data=='tsv':
            df = pd.read_table(input)

        graph_figa = dir + "/graph_" + str(time.time()) + ".eps"
        graph_fig = dir + "/graph_" + str(time.time()) + ".png"

        # X = df.iloc[:, 1:].apply(lambda x: (x-x.mean())/x.std(), axis=0)

        if form.graph_type.data == "cat":
            if form.jitter.data:
                sns.catplot(x=df.columns.values[0],y=df.columns.values[1], data=df)
            else:
                sns.catplot(x=df.columns.values[0],y=df.columns.values[1], data=df, jitter=False)

            plt.savefig(graph_figa)
            plt.savefig(graph_fig)

            img_url = "../../static/" + graph_fig

            return render_template('/tools/cat_graph.html', form=form, img_url=img_url)

        if form.graph_type.data == "box":
            sns.catplot(x=df.columns.values[0],y=df.columns.values[1], data=df, kind='box')
            plt.savefig(graph_fig)
            plt.savefig(graph_figa)

            img_url = "../../static/" + graph_fig

            return render_template('/tools/cat_graph.html', form=form, img_url=img_url)

        if form.graph_type.data == "violin":
            sns.catplot(x=df.columns.values[0],y=df.columns.values[1], data=df, kind='violin')
            plt.savefig(graph_fig)
            plt.savefig(graph_figa)

            img_url = "../../static/" + graph_fig

            return render_template('/tools/cat_graph.html', form=form, img_url=img_url)

    return render_template('/tools/cat_graph.html', form=form)

@app.route('/graph/dl/')
def dl_graph_index():
    user = g.user.name
    dir = "./user/" + user + "/graph"
    file = "./user/" + user + "/graph_figs"
    if os.path.exists(dir):
        shutil.make_archive(file, 'zip', dir)

    downloadFileName = 'figs' + str(time.time()) + '.zip'
    downloadFile = "./static/user/" + user + '/graph_figs.zip'

    return send_file(downloadFile, as_attachment = True, attachment_filename = downloadFileName)

werkzeugの使い方が変わった件

$ python3 manage.py runserver
Traceback (most recent call last):
  File "manage.py", line 2, in <module>
    from app import app, manager
  File "/Users/kuro/database/app/__init__.py", line 11, in <module>
    from app.models import *
  File "/Users/kuro/database/app/models.py", line 4, in <module>
    from werkzeug import check_password_hash, generate_password_hash
ImportError: cannot import name 'check_password_hash' from 'werkzeug' (/usr/local/lib/python3.7/site-packages/werkzeug/__init__.py)

となってFlaskが動かない。

以前はmodels.pyの中で

from werkzeug import check_password_hash, generate_password_hash

としておけばよかったのだが、werkzeugの構成に変更があったらしく

from werkzeug.security import check_password_hash, generate_password_hash

このように書くようになった。

CentOS7のKVMにUbuntuを乗せる

$ sudo virt-install --name kvm3_ubuntu --memory 8192 --disk /home/kvm/images/kvm3,size=100  --location /tmp/ubuntu-20.04-desktop-amd64.iso --network default --vcpu 8 
WARNING  オペレーティングシステムを検出できません。仮想マシンのパフォーマンスが低下する可能性があります。最適なパフォーマンスを得るには、--os-variant で OS を指定する必要があります。

インストールの開始中...
ファイル .treeinfo を読出中...                      |    0 B  00:00     
ファイル content を読出中...                        |    0 B  00:00     
ファイル info を読出中...                           |   57 B  00:00     
ファイル info を読出中...                           |   57 B  00:00     
ファイル info を読出中...                           |   57 B  00:00     
ERROR    Ubuntu ツリーに対する hvm カーネルが見つかりませんでした。
仮想マシンのインストールが成功したように見えません。
成功していれば、次を実行すると、仮想マシンを再起動できます:
  virsh --connect qemu:///system start kvm3_ubuntu
そうでなければ、インストールをやり直してください。

およよ。

どういうわけかダウンロードしたisoイメージからインストールができない。

$ sudo virt-install --name kvm3_ubuntu --memory 8192 --disk /home/kvm/images/kvm3,size=100  --location http://ftp.jaist.ac.jp/pub/Linux/ubuntu/dists/bionic/main/installer-amd64/ --network default

これでいけた。

追記
ローカルに置いたisoイメージでインストールするときは--locationでなく--cdromとする。

KVM環境をホストの電源オンオフと連動

ホストの電源を入れたら自動で起動させるには

$ sudo virsh list --all
 Id    Name                           State
----------------------------------------------------
 1     centos7                        running
 -     centos8                        shut off

$ sudo virsh autostart centos7

ホストをシャットダウンするとき自動的にゲストもシャットダウンする。
Ubuntuの場合

$ sudo nano /etc/default/libvirt-guests 
# URIs to check for running guests
# example: URIS='default xen:/// vbox+tcp://host/system lxc:///'
#URIS=default

# action taken on host boot
# - start   all guests which were running on shutdown are started on boot
#           regardless on their autostart settings
# - ignore  libvirt-guests init script won't start any guest on boot, however,
#           guests marked as autostart will still be automatically started by
#           libvirtd
ON_BOOT=ignore

# Number of seconds to wait between each guest start. Set to 0 to allow
# parallel startup.
#START_DELAY=0

# action taken on host shutdown
# - suspend   all running guests are suspended using virsh managedsave
# - shutdown  all running guests are asked to shutdown. Please be careful with
#             this settings since there is no way to distinguish between a
#             guest which is stuck or ignores shutdown requests and a guest
#             which just needs a long time to shutdown. When setting
#             ON_SHUTDOWN=shutdown, you must also set SHUTDOWN_TIMEOUT to a
#             value suitable for your guests.
ON_SHUTDOWN=shutdown

# Number of guests will be shutdown concurrently, taking effect when
# "ON_SHUTDOWN" is set to "shutdown". If Set to 0, guests will be shutdown one
# after another. Number of guests on shutdown at any time will not exceed number
# set in this variable.
PARALLEL_SHUTDOWN=10

# Number of seconds we're willing to wait for a guest to shut down. If parallel
# shutdown is enabled, this timeout applies as a timeout for shutting down all
# guests on a single URI defined in the variable URIS. If this is 0, then there
# is no time out (use with caution, as guests might not respond to a shutdown
# request). The default value is 300 seconds (5 minutes).
SHUTDOWN_TIMEOUT=120

# If non-zero, try to bypass the file system cache when saving and
# restoring guests, even though this may give slower operation for
# some file systems.
#BYPASS_CACHE=0

# If non-zero, try to sync guest time on domain resume. Be aware, that
# this requires guest agent with support for time synchronization
# running in the guest. For instance, qemu-ga doesn't support guest time
# synchronization on Windows guests, but Linux ones. By default, this
# functionality is turned off.
#SYNC_TIME=1

CentOS7の場合

$ sudo nano /etc/sysconfig/libvirt-guests

Linuxbrew + python3でglibcのバージョンが足りない件(解決)

新規構築したCentOS7環境にLinuxbrewでpython3.7をインストールするとmatplotlibがインポートできない現象が出て、どうにかpython3.6を入れてお茶を濁そうとしていたわけだ。
Mac OSXのHomebrewではgit checkoutで古いバージョンのインストールができるようなことがあちこちにかかれていたのだけれど、Linuxbrewでは結局成功しなかった。

途方に暮れかけ、ネットを検索していて、これは!という記事があった、とおもったらなんと自分の書き込み(つまりこのサイト)だったw

$ brew install glibc

たったこれだけのことで解決してしまった。

CentOS7はGLIBC1.17までしか対応せず、1.18はインストールのしようがないということだったけど、brewでは2.23がインストールできちゃうのね。

SRAToolkit

公開されているSRAファイルをダウンロードして再解析がしたいとき、データのダウンロードが結構面倒なので(サイトの構造が複雑すぎてなかなかファイル本体にたどり着けない)専用のツールを利用したい。
とおもってインストールをbrewからやってみた。

$ brew install sratoolkit

ところがだ。

$ fastq-dump --gzip --split-files SRRxxxxxxxx

localeがおかしいとかのwarningはまあいいとして、どうもperlのライブラリがぶつかって動かないらしい。


というわけで、バイナリを手動でインストールしてみた。

$ wget https://ftp-trace.ncbi.nlm.nih.gov/sra/sdk/current/sratoolkit.current-centos_linux64.tar.gz
$ tar xvzf sratoolkit.current-centos_linux64.tar.gz
$ mv sratoolkit.2.10.7-centos_linux64 ~/src/
$ ln -s ~/src/sratoolkit.2.10.7-centos_linux64/bin/* ~/bin

homebrewでpython3.6(前のバージョン)をインストールする

$ cd /home/linuxbrew/.linuxbrew/Homebrew/Library/Taps/homebrew/homebrew-core/Formula
$ git log python.rb
commit b02eee23d9ac8e1027a212f269da6edd4c1f0f5f
Author: BrewTestBot <homebrew-test-bot@lists.sfconservancy.org>
Date:   Mon Apr 20 06:02:17 2020 +0000

    python: update 3.7.7_1 bottle.

commit 8bcd19011df9b5d42f5286bd5f3ec65e2aafc8c7
Author: Michka Popoff <michkapopoff@gmail.com>
Date:   Mon Apr 20 08:01:09 2020 +0200

    python: bump revision for libffi (#20152)

commit 4f455d5f8c2d45682448416d53d0c47355445177
Author: Michka Popoff <michkapopoff@gmail.com>
Date:   Sat Apr 18 18:17:50 2020 +0200

    python: fix style for Linuxbrew

・・・・・

こんな感じにgit logを開く
とにかく3.7を入れちゃうとモジュール関係も最新を入れてしまうのでGLIBC_2.18が必要となってしまうので、3.6に戻ってもらう。

commit ac64d11bca23f0d221967fef3b89e1d30053e067
Author: LinuxbrewTestBot <testbot@linuxbrew.sh>
Date:   Sun Jun 17 15:53:29 2018 +0000

    python: update 3.6.5_1 bottle for Linuxbrew.

commit 3ec730ba2ea8fcbd8fb853778b07236c3edb02c1
Merge: d07628b798 597da37dcc
Date:   Sun Jun 17 17:23:04 2018 +0200

    Merge branch homebrew/master into linuxbrew/master
    
     Conflicts:
            Formula/gdbm.rb
            Formula/pipenv.rb
            Formula/pre-commit.rb
            Formula/python.rb
            Formula/python@2.rb
            Formula/sip.rb
            Formula/vim@7.4.rb

commit f2a764ef944b1080be64bd88dca9a1d80130c558
Author: BrewTestBot <brew-test-bot@googlegroups.com>
Date:   Sun Jun 17 13:24:46 2018 +0000

    python: update 3.6.5_1 bottle.

commit b60c0325575ed4c3c90ba061880a79bc0c208cf0
Author: commitay <commitay@users.noreply.github.com>
Date:   Sun Jun 17 18:01:44 2018 +1000

    python: revision for gdbm
    
    Also, update `setuptools` and `wheel` resources

commit 5f5331172328b5d21592518fece5c332c42bc770
Merge: dbc1429eab c1f561a772
Author: Michka Popoff <michkapopoff@gmail.com>
Date:   Wed Jun 13 23:07:55 2018 +0200

    Merge branch homebrew/master into linuxbrew/master

commit c1f561a772a747d711657afa2644287aa0cc9a40
Author: commitay <commitay@users.noreply.github.com>
Date:   Wed Jun 13 19:26:30 2018 +1000

    python 3.7.0rc1 (devel) (#28975)

ここでやっと3.6登場

commit 0036460b7cf3f16e5e02880d893717ee86404ef3
Author: Chuancong Gao <chuanconggao@users.noreply.github.com>
Date:   Thu Mar 29 19:40:41 2018 -0700

    python 3.6.5

このへんで行くか。

$ git checkout 0036460b7cf3f16e5e02880d893717ee86404ef3 python.rb
$ HOMEBREW_NO_AUTO_UPDATE=1 brew install python
Warning: Calling 'devel' blocks in formulae is deprecated! Use 'head' blocks or @-versioned formulae instead.
Please report this issue to the homebrew/core tap (not Homebrew/brew or Homebrew/core), or even better, submit a PR to fix it:
  /home/linuxbrew/.linuxbrew/Homebrew/Library/Taps/homebrew/homebrew-core/Formula/python.rb:15

==> Downloading https://files.pythonhosted.org/packages/72/c2/c09362ab29338413ab
######################################################################## 100.0%
==> Downloading https://files.pythonhosted.org/packages/c4/44/e6b8056b6c8f2bfd14
######################################################################## 100.0%
==> Downloading https://files.pythonhosted.org/packages/fa/b4/f9886517624a4dcb81
######################################################################## 100.0%
==> Downloading https://www.python.org/ftp/python/3.6.5/Python-3.6.5.tar.xz
##############################                                   [  938.835230] serial8250: too much work for irq4
######################################################################## 100.0%
Warning: Calling 'devel' blocks in formulae is deprecated! Use 'head' blocks or @-versioned formulae instead.
Please report this issue to the homebrew/core tap (not Homebrew/brew or Homebrew/core), or even better, submit a PR to fix it:
  /home/linuxbrew/.linuxbrew/Homebrew/Library/Taps/homebrew/homebrew-core/Formula/python.rb:15

Error: An exception occurred within a child process:
  NoMethodError: undefined method `installed?' for OS::Mac::CLT:Module
Did you mean?  instance_of?

うーん、やっぱりだめか?

次はanacondaでトライ

Python3.7.7/CentOS7でエラー

Traceback (most recent call last):
  File "manage.py", line 2, in <module>
    from app import app, manager
  File "/home/kuro/database/app/__init__.py", line 23, in <module>
    from app.views import samples, genes, anno, fasta, blast, blast_nb, crispr, clustalw, seeds, plasmid, stat
  File "/home/kuro/database/app/views/crispr.py", line 8, in <module>
    import matplotlib
  File "/home/linuxbrew/.linuxbrew/Cellar/python/3.7.7_1/lib/python3.7/site-packages/matplotlib/__init__.py", line 205, in <module>
    _check_versions()
  File "/home/linuxbrew/.linuxbrew/Cellar/python/3.7.7_1/lib/python3.7/site-packages/matplotlib/__init__.py", line 190, in _check_versions
    from . import ft2font
ImportError: /lib64/libc.so.6: version `GLIBC_2.18' not found (required by /home/linuxbrew/.linuxbrew/lib/libstdc++.so.6)

matplotlibをインポートしようとするとGLIBC_2.18が無いと言って落ちる。
CentOS7ではGLIBC_2.18はサポートしないらしいので、そろそろCentOS7も8に上げるべき時期が来たのかもしれない。

KVM仮想環境でテストを開始するか。
それともbrewで入れずanacondaにしてみるか?

brewで新規にインストールすると3.7.7がインストールされてしまうので、brewの環境を一旦アンインストールし、yumで標準リポジトリから3.6系を入れてみることにした。そうすると問題はなくなることは確認できたが、yumで入れるといちいちsudoでパッケージを追加したりと若干使い勝手が悪い。

次はbrewで明示的に3.6系をインストールする方法を試してみる。
参照
macOS Mojave に Python 3.6 環境を構築する - Satoshi Oikawa - Medium

cDNA FASTAファイルから最長ORFを抽出し、5'UTR/CDS/3'UTRに分割してそれぞれのFASTAファイルを作成する(改訂版)

#DNA FASTAファイルから最長のORFを抽出し5UTR,CDS,3UTRに分割して保存する。

import sys, os, re
from Bio import SeqIO
from Bio.Alphabet import IUPAC
from Bio.Seq import Seq

fasta_file = sys.argv[1]
faname = os.path.basename(fasta_file)
fdir = os.path.dirname(fasta_file)
fname = os.path.splitext(faname)[0]
fext = os.path.splitext(fasta_file)[1]
utr_5_fasta = os.path.splitext(fasta_file)[0] + '_5utr' + fext
cds_fasta = os.path.splitext(fasta_file)[0] + '_cds' + fext
utr_3_fasta = os.path.splitext(fasta_file)[0] + '_3utr' + fext
for record in SeqIO.parse(fasta_file, 'fasta'):
    cds = re.findall('(ATG(?:\S{3})*?T(?:AG|AA|GA))', str(record.seq).upper())

    if cds:
        match = max(cds, key = len, , default=0)
        seq = Seq(match, IUPAC.ambiguous_dna)
        utr_5 = re.sub(str(seq.strip()+'[ACGT]*'), '', str(record.seq).upper())
        utr_3 = re.sub(str(utr_5 + seq.strip()), '', str(record.seq).upper())
        with open(utr_5_fasta, 'a') as f:
            f.write(">" + record.id +'\n')
            f.write(str(utr_5) + '\n')
        with open(cds_fasta, 'a') as f:
            f.write(">" + record.id +'\n')
            f.write(str(seq) + '\n')
        with open(utr_3_fasta, 'a') as f:
            f.write(">" + record.id +'\n')
            f.write(str(utr_3) + '\n')
    else:
        t_cds = re.findall('(ATG[ACGT]*)', str(record.seq).upper())
        if t_cds:
            match = max(t_cds, key = len)
            seq = Seq(match, IUPAC.ambiguous_dna)
            utr_5 = re.sub(str(seq.strip()+'[ACGT]*'), '', str(record.seq).upper())
            utr_3 = re.sub(str(utr_5 + seq.strip()), '', str(record.seq).upper())
            with open(utr_5_fasta, 'a') as f:
                f.write(">" + record.id +'\n')
                f.write(str(utr_5) + '\n')
            with open(cds_fasta, 'a') as f:
                f.write(">" + record.id +'\n')
                f.write(str(seq) + '\n')
            with open(utr_3_fasta, 'a') as f:
                f.write(">" + record.id +'\n')
                f.write(str(utr_3) + '\n')

CDSがstopまで入っていないとき5UTRも切り出せていなかったところを改善してみた。