投稿

4月, 2018の投稿を表示しています

macで外付け英字キーボード

iMac純正のキーボードが打ちにくくて腱鞘炎になりそうだったのでlogicoolの英字キーボードを買った。キーマップを変えるのにkalabinarを導入した。 https://github.com/tekezo/Karabiner-Elements の You can download the latest Karabiner-Elements from https://pqrs.org/osx/karabiner/ をダウンロード。Launchpadからkalabinerを開いて、 add items で以下のように変更していく。 元のキー(印字) 編集後のキー caps_lock left_control left_command(Windows) fn left_control left_command left_option(Alt) left_option これを機に ・Ctrl-V,C,Xを親指使ってやるのは諦めた ・日本語の切り替えは覚え直し とした。 とした。 追記 karabinerはいらなかった。OS標準のkey mapperでいけた。

Inverse-normal transformation

Inverse normal transformationとは正規分布しないデータをランク変換してそのランクを正規分布の累積分布関数で変換すること。変換後の値は正規分布する。Rでやると qnorm((rank(x)-0.5)/length(x)) となる。正規分布しない量的形質に正規性を持たせるためにこのような処理を行う。 BMIを目的とし、年齢・性別・バッチでで補正してQTLをしたい場合、まずlmで残差を計算して、この残差をInverse normal transformationして形質として回帰分析を行う。 res = lm(BMI〜Age+Sex+Batch)$residuals std_res = qnorm((rank(res)-0.5)/length(res)) のような感じ。 バラ色の解決策でもない様子 https://www.ncbi.nlm.nih.gov/pubmed/19526352

optparseの使い方

pythonでコマンドラインツールを作っている。普通にpositionalに指定してもよいのだが、できれば機能を大幅に変えたいと思っている。そうするとpositionalな引数の指定では不都合が多い。簡単にコマンドラインオプションの設定ができるoptparseというパッケージがある。 # option.py from optparse import OptionParser parser = OptionParser() parser.add_option('-i', '--infile', action="store', dest='inputfile', help="name of input file') options, args = parser.parse_args() こうしておくと辞書形式で options にそれぞれのoptionを保存してくれる。また自動でヘルプドキュメントを作成してくれる。 $ python option.py -h Usage: test2.py [options] Options: -h, --help show this help message and exit -i INPUTFILE, --infile=INPUTFILE name of input file その他に設定可能なオプションたち。 add_option( dest = varname, action = 'store/store_const/store_true/store_false', type = 'int/float/long/complex', default = 'string/True/False', choices = ['a', 'b', 'c'] )

joinに関してshellとRとpythonと

結論 import pandas as pd leftfile = pd.read_csv("leftfile.txt", sep="\t") rightfile = pd.read_csv("rightfile.txt", sep="\t") mergedfile = pd.merge(leftfile, rightfile, how='outer', sort=True) 欠損値を含むいくつかのテーブルをある値をキーにしてjoinしたい。もちろんヘッダ込みで、欠損はNA埋めて。こういう場合いくつか方法があるように思いますが、環境的にpythonしか使えなかったという。 Rで Rではdplyrが便利です dplyr::full_join(leftfile, rightfile, by="V1") しかし権限がなくてライブラリのインストールができなかったという・・・。 Shellで shellは以下のように実現できます join -a 1 -a 2 -e "NA" -o auto leftfile.txt rightfile.txt ヘッダは先にのぞいておく必要がある -o auto が効かないバージョンだと-o 1.1,1.2,1.3…2.2,2.3…と延々オプションを指定することになる などという問題が。 pythonで pandasは使えた() import pandas as pd leftfile = pd.read_csv("leftfile.txt", sep="\t") rightfile = pd.read_csv("rightfile.txt", sep="\t") mergedfile = pd.merge(leftfile, rightfile, how='outer', sort=True) なんとなくshellでjoinを使うのが早いのかなと思っていましたが、dplyrは早いですし、pandasも割と軽快でした。

GRS

Mega et al Lancet 2015 によると、I個のマーカーによるjさんのGRSは以下の式で計算される。 G R S j = ∑ i I [ l o g ( O R i ) × G e n o t y p e i j ] GRS_j = \sum_i^I \big[ log(OR_i) \times Genotype_{ij} \big] G R S j ​ = i ∑ I ​ [ l o g ( O R i ​ ) × G e n o t y p e i j ​ ] 複数のGRSを比較するときは、pupulation全体のmean GRSを引いて、sd GRSで割るという方法が用いられるようです(Tada et al. EHJ 2017)。

変数展開

いつも忘れるのでメモ。#は前置パターン、%は後置パターンと覚えよう。そして二つ置くと最長、一つだと最短と覚える。つまり var=/tmp/jorijori/txt.long.sep echo var /tmp/jorijori/txt.long.sep echo {var#*/} tmp/jorijori/txt.long.sep echo {var##*/} txt.long.sep echo {var$.*} /tmp/jorijori/txt.long echo {var$.*} /tmp/jorijori/txt のとなる。

MacのExcelで改行コードが変

情弱なのでexcelを使う時があります。そしてこれを再度csv/tsvに書き出してスクリプトを組みたいことがあります。しかし、Excel for Macでは書き出した時に改行コードが ^M とかになります。 wc -l しても行数0になります。UNIXの改行コードは通常 \n ですが、excelは \r で返してきます、trを使って変換しましょう。 tr unreadable.tsv | tr "\r" "\n" > readable.tsv

tailで先頭行を読み飛ばす

先頭行を簡単に読み飛ばすのに tail -n +2 hoge.txt が使える。逆に最後の一行を読み飛ばすのは tail -n -1 hoge.txt

awkの正規表現にシェル変数を

awkで行を抜き出す時などにシェルから変数を渡して、それを正規表現内で用いて処理したいと常から思っていたがうまくいっていなかった。下記のようにすれば良いらしい。変数宣言のところでダブルクオーテーションで囲んで正規表現を作っておいてこれを~で呼び出せば良い。 tmp=22 awk -v pattern="chr{tmp}_" '1~pattern {print $0}' これでchr22_という文字列が1番目のカラムに含まれる行のみ出力される、

SnpEff

Javaで書かれている。何も書かれていないVCFファイル(test.chr22.vcf)を処理する。 #CHROM POS ID REF ALT QUAL FILTER INFO 22 17071756 . T C . . . 22 17072035 . C T . . . 22 17072258 . C A . . . 22 17072674 . G A . . . 22 17072747 . T C . . . 22 17072781 . C T . . . 22 17073043 . C T . . . 22 17073066 . A G . . . 22 17073119 . C T . . . 以下のスクリプトで処理すると java -Xmx4g -jar snpEff.jar GRCh37.75 examples/test.chr22.vcf > test.chr22.ann.vcf こうなる ##SnpEffVersion="4.3t (build 2017-11-24 10:18), by Pablo Cingolani" ##SnpEffCmd="SnpEff GRCh37.75 examples/test.chr22.vcf " ##INFO=<ID=ANN,Number=.,Type=String,Description="Functional annotations: 'Allele | Annotation | Annotation_Impact | Gene_Name | Gene_ID | Feature_Type | Feature_ID | Transcript_BioType | Rank | HGVS.c | HGVS.p | cDNA.pos / cDNA.length | CDS.pos / CDS.length | AA.pos / AA.length | Distance | ERRORS / WARNINGS / INFO' "> ##INFO=<ID=LOF,Numbe...

localでnvimをbuild

管理者権限なしでのnvimのbuild(RHLE7) 前にもよく似た記事を書いたが、今回はRHLE7でやったためうまくいかないところがあったので修正している。 git clone https://github.com/neovim/neovim.git cd neovim make CMAKE_EXTRA_FLAGS="-DCMAKE_INSTALL_PREFIX=HOME/neovim" make install export PATH="HOME/neovim/bin:$PATH" だけでいけた。 ~/.config/nvim/init.vim は自分で作成すること。以下、deinを用いたプラグインのインストール。init.vimと同じディレクトリにdein.tomlというファイルを置いておく init.vim let s:dein_dir = expand('~/.cache/dein') if &runtimepath !~# '/dein.vim' let s:dein_repo_dir = s:dein_dir '/repos/github.com/Shougo/dein.vim' if !isdirectory(s:dein_repo_dir) call system('git clone https://github.com/Shougo/dein.vim ' . shellescape(s:dein_repo_dir)) endif execute 'set runtimepath^=' . s:dein_repo_dir endif if dein#load_state(s:dein_dir) call dein#begin(s:dein_dir) let s:toml_dir = expand('~/.config/nvim') call dein#load_toml(s:toml_dir . '/dein.toml', {'lazy': 0}) call dein#end() call dein#sav...

rsync

rsyncを使ったことがなかった。使い方のメモ。 rsync -av -e ssh USER@HOST:ORIGIN DIST ORIGINはHOSTにおける絶対パスもしくはUSERのホームディレクトリからの相対パス。DISTはローカルのディレクトリ。これでうまくいくはずが、なぜかリモートサーバに接続できない・・・。調べるとログイン時にtmuxを立ち上げる設定にしているとうまくいかないようだ。 -t オプションで解決できる。 rsync -av -e ssh -t USER@HOST:DIST オプション 効果 -delete 削除も同期する –existing 更新ファイルのみ同期 –exclude 除外対象の指定 –progress progress -v 経過を表示 -r 再帰的にコピー -l シンボリックリンクもコピー -p パーミッションもコピー -g グループもコピー -o ファイル所有者もコピー -a -rlpgoD -D –devices --specialsと同じ –devices ブロックデバイスをコピー –specials 特殊ファイルをコピー -z データ転送時に圧縮