pexpectとclickで劇的ビフォーアフター

何故か最近GAMESSでFMOを計算することに全力を注いでおりまして。

GAMESSでFMO計算を行う場合にはfmoutilという便利ツールを使うことが多いのだが、これは対話的な環境の上にデフォルトというものが用意されていないので、毎回毎回手打ちしないといけない。

それから、きれいな複合体データでないと上手く動かない(途中に欠損があると駄目)とか実践で使っていく場合には色々と不都合があるので、今は自分用の 超便利有能空気読みすぎpython版のプリプロセッサ をガシガシ書いている。

で、ヘッダーチェックのためにfmoutilをすごい重宝してるんだけどなんせ毎回手打ちで1日何十回も対話していると流石に疲れる。

pexpectclickを使ってこんな感じに実行できるようにしてみた。

$ ./fuu --help
Usage: fuu [OPTIONS] IFILE

  Create FMO Input

Options:
  -o, --ofile TEXT         output option (default: fmo.inp)
  -n, --node INTEGER       node option (default: 4)
  -c, --core INTEGER       core option (default: 8)
  -m, --memory INTEGER     memory option (default: 48000)
  -w, --wavefunction TEXT  wavefunction option (default: RHF)
  -b, --basis-set TEXT     basis set option (default: STO-3G)
  --help                   Show this message and exit.

コードは至ってシンプルなのにモダンなfmoutilに生まれ変わって満足

#!/usr/bin/python
import pexpect
import click

bs = {"STO-3G": 1, "3-21G": 2, "6-31G": 3, "6-31G*":4, "6-311G*":5}
wf = {"RHF": 1, "DFT": 2, "MP2": 3, "CC": 4, "MCSCF": 5}

@click.command(help="Create FMO Input")
@click.argument("ifile")
@click.option("--ofile",  "-o", default="fmo.inp", help="output option (default: fmo.inp)")
@click.option("--node",  "-n", default=4, help="node option (default: 4)")
@click.option("--core",  "-c", default=8, help="core option (default: 8)")
@click.option("--memory",  "-m", default=48000, help="memory option (default: 48000)")
@click.option("--wavefunction",  "-w", default="RHF", help="wavefunction option (default: RHF)")
@click.option("--basis-set",  "-b", default="STO-3G", help="basis set option (default: STO-3G)")
def cli(ifile, ofile, node, core, memory, wavefunction, basis_set):
    c = pexpect.spawn("./fmoutil")
    c.expect(".*") #   Enter Job # :     1. Generate FMO input data for GAMESS
    c.sendline("1")
    c.expect(".*") #   Enter Input PDB File(s) :
    c.sendline(ifile)
    c.expect(".*") #   Enter Output File : 
    c.sendline(ofile)
    c.expect(".*") # 1> How many computer nodes do you want to use ?
    c.sendline("{}".format(node))
    c.expect(".*") # 1> How many CPU cores does each node have ?
    c.sendline("{}".format(core))
    c.expect(".*") # 1> How much memory does each node have (in megabyte) ?
    c.sendline("{}".format(memory))
    c.expect(".*") # 1> Choose runtyp (1:energy, 2:gradient, 3:geometry optimization)
    c.sendline("1") 
    c.expect(".*") # 1> How many layers do you want to define (1-5) ?
    c.sendline("1")
    c.expect(".*") # 1> Choose wavefunction type (1:RHF, 2:DFT, 3:MP2, 4:CC, 5:MCSCF)
    c.sendline("{}".format(wf[wavefunction]))
    c.expect(".*") # 1> Enter basis set (1:STO-3G, 2:3-21G, 3:6-31G, 4:6-31G*, 5:6-311G*)
    c.sendline("{}".format(bs[basis_set]))
    c.expect(".*") # 1> Do you add diffuse functions on COO- groups ? (1:yes, 2:no)
    c.sendline("2")
    c.expect(".*") # 1> Enter the desired n-body expansion (2 or 3) ?
    c.sendline("2")
    c.expect(".*") # 1> Would you like to perform the pair analysis (PIEDA) ? (1:yes, 2:no)
    c.sendline("1")
    c.expect(".*") # 1> Would you like to print Mulliken charges (1:yes, 2:no) ?
    c.sendline("1")
    c.expect(".*") # 1> Would you like to produce a cube file with the total electron density ? (1:no, 2:standard, 3:sparse)
    c.sendline("1")
    c.expect(".*") # 1> Whould you like to use PCM ? (1:yes, 2:no)
    c.sendline("2")
    c.expect(".*") # 1> Enter fragment size (1:1res-per-frg, 2:2res-per-frg)
    c.sendline("1")
    c.expect(".*") # 1> are S-S bonded CYSs combined to one ? (1:yes, 2:no)
    c.sendline("1")
    c.expect(".*") # 1> is GLY combined to the neighbor ? (1:yes, 2:no)
    c.sendline("2")
    c.expect(".*") # Enter Job # :
    c.sendline("")
    c.close()

if __name__ == "__main__":
    cli()

以上、clickの提供でお送りしました☆

あと、clickで--helpしたときにデフォルトの値も表示させたかったんだけど、やりかたわからんかった。誰か知ってたら教えてください。尚、オプションにchoiceを入れてないのは僕の怠慢ですw

半自動スクレイピングに最適なファイル形式、それがmht

reCAPTCHAみたいな画像認証が組み込まれていてロボットも弾くようになっているサイトだけど、どうしてもスクレイピングしたいから、マニュアルアクセスでページを保存してスクレイピングのプログラムにまわしたい、さらにその際に表示されている画像も一緒にとりたいというなかなか難しそうな課題に取り組んでいた。

IEで保存すると(完全とかいうオプション選ぶんだっけ?忘れた)フォルダに画像ファイルとかまとめて保存してくれるんだけどスクレイピングにまわすの不便というか、ウェブサービスにしたいからuploadしやすいファイル形式がいいなぁと調べていたらIEだとmhtという形式で保存できるということに気がついた。

調べたら相当前からあったらしい。

そして単一のファイルでMIMEなのも便利なところ。

RFC 2557 は、MHTMLは電子メール以外のプロトコル (HTTPやFTP) でも一度に転送出来ることや、完全なHTML文書のアーカイブとして保存出来ることを示唆している。

unpackの実装見てみたら、emailモジュール使ってunpackしていた。まぁそうだろうなと。

そこまでわかればスクレイパーをサクッと書いてやりたい仕事は終了した。

もう、半自動スクレイピングは怖くないw

あとmhtをそのままメールで送ったらアクセスポリシーがとか画像の表示を許可しますか?とか言われなくていいんではないかと思ったけどどうなんだろう?

Python Testing読みなおそうっと

pygamessをopenbabelからRDKitに依存させるように書き換えてるんだけど、pyvowsのサイトにずっとアクセスできなくてnoseでいいじゃんってなってさっきpyvowsで書いてたテストを削除した。

pychembldbもunittest使ってるしpyvowsとはサヨナラだな。

そして、色々忘れていてあれなのでPython Testingを読み直すことにした。バッグが更に重くなるけど仕方ない。

ProductName Python Testing: Beginner's Guide
Daniel Arbuckle
Packt Publishing / ?円 ( 2010-01-19 )


pythonのformatでfloatの正負の記号を無視して揃える

ドットの前にスペース入れるらしい。解決するのに大分かかったのでメモ

>>> print "{:.1f}\n{:.1f}".format(1.0, -1.0)
1.0
-1.0
>>> print "{: .1f}\n{: .1f}".format(1.0, -1.0)
 1.0
-1.0

Shizuoka.pyお疲れ様でした

今回も盛り上がって良かったです。大いに楽しみました。

懇親会も美味しかったので良かった。席数が100を超えてる店で開店近くの飛び込みなのにギリギリだったのには驚いたけど人気店なんだなと納得。接客すごく良かったです。

1487591982

懇親会では色々ディープな話が聞けたけど、個人的には新しいmacbookproを思いとどまれたのが一番の収穫でした

早速帰ってから発注して、もう手元にあります(明日早起きして付け替え作業しようかなw)

それから、ドキドキするので発表時間は過少申告しないように気をつけましょうw

次回は夏くらいにまた西の方でやりたいですね。浜松の方でうなぎを食べるとかの美味しいものをセットにしたい。

週末はShizuoka.pyです

週末はShizuoka.pyですので、SLに乗るなり、でラーメンを堪能するなり、楽しんでください。

尚僕はまだ資料に取り掛かってない、ヤバイ…

Shizuoka.pyを2/18にやります

一度年末に調整したのですが、タイミングが合わずにこの日にやることになりました。Shizdevのほうでも人工知能まわりをやっているので今回は共催となります。

python使いが静岡東部に多い+コミュニティFの利用料金が無料なのに綺麗という理由で最近のShizuoka.pyは富士でばかりやっていましたが、流石に他の地域で美味しいものを食べながらの懇親会をしなきゃなという理由で今回島田開催です。

尚、懇親会場は決まってないのでリクエストがあれば@fmkz___@yajuまでよろしくお願いします。

島田のとなりの金谷には大井川鉄道があってSLが走っているので、鉄成分高めの方は満足するんじゃないかなぁと思います。僕は鉄成分がホメオパシー並に薄いのですが五和で撮ったとき(connpassの写真、今回は料理ではない)には「SLウォー!!!!」ってなりました。五和の中屋酒店に行った時にたまたま遭遇したんだけどねw

PyConJP2016に行ってきた

まず予習しておいた戸山公園にはピカチュウはいなかった。東側も行ってみたけどブーバーくらいしかいなかったし、高低差があるので歩くのがめんどくさかった。西側はオニスズメばかりw

PyConJPは去年は不参加だったので二年ぶりの参加でした。初中級向けの演題が多かった気がして、間口が広げたかったのかなぁと思ったが一方でディープな話題が少なくて初期の頃とは色々変わっているんだなぁと感じた。これは選考プロセスのせいもあるのかな?

あとディープラーニング関連の演題が多かった。一昨年は機械学習の話題も多かった気がするけど。僕の最近の興味もディープラーニングなので、Pythonというよりはディープラーニングの演題を聴きに行った(というかほとんどそれしか聴いてないw)

尚、狙った演題を椅子に座って聞くためには一つ前の演題をスキップしないとしんどいことになるのでなかなか厳しかった。もう少し自由に出入りできると良いのかなと思うが、大学では仕方ないかな(逆にカジュアルに出入りできる教室の設計だと困るw)

深層学習の次の一手

次の一手ということで、なんとなくディープな話が聞けるのかなと思ったけど、クラウドサービスの話が多かった印象。 自社サービスの宣伝するならもっとグイグイやってくれたほうがよかったかな。

GPUを使える仮想環境のニーズって個人レベルだとかなりあると思うんだよね。気軽に試したいし。

尚、CNTKはちょっと興味を持った。

ニューラルネットワークのフレームワークであるChainerで始める対話Botの作成

自然言語処理系のDLはあまり追いかけてなくて、 word2vecくらいしか理解してないんだけど、あれはディープじゃないから。RNNとかLSTMをつかってやろうとしたことの話で、具体的なコードが出てこなかったけど、やろうとしていることは理解できて面白かった。

LSTMと化合物のフラグメントってどう対応とったらうまいくのかなぁとか思った。ゲノムのほうがやりやすい気はするけどそういう論文は見たことないなぁ。

確率的ニューラルネットの学習と Chainer による実装

とても勉強になった。ツイートでメモとってたので、帰りの新幹線で読み返してた。

画像生成と化合物合成って結構似ていると思っているので、こういう方法で化合物生成するのって悪くはないんじゃないかなぁーと。

そのうち試す。

Pythonでpyftpdlibを使ってFTPサーバーを作る際に使ったテクニックの紹介

要件定義と言うか仕様を固めるの難しいよなぁと思いながら聞いていた。パッケージを作りたいヒトは聞いておいて良い内容だったと思います。僕も、あのときああいう風にすればよかったのかとか思いながら聞いていた。

Deep Learning with Python & TensorFlow

tensorflowの触りだけだった。初めてのヒト向けに具体的な事例を紹介する内容だったと思う。知らない人向けには面白かったと思うが、ある程度使っているヒトにとってはちょっと物足りなかったかも。

二年ぶりの参加だったけど、色々と刺激を受けて楽しかった二日間だった。来年も参加します。

DNN Regressionの自分の実装がダメな件

TensorflowのDNN Regressionを使っていて便利なんだけどハイレベルなAPIばっかり使っていると自分でもっと別なの実装したくなった時に困るだろうなぁと、自分で実装してみたら…

$python boston.py
 RF: R2: 0.736624, RMSE:4.631010
SVM: R2: 0.515467, RMSE:6.281300
TFL: R2: 0.647992, RMSE:5.353827
 MY R2: -1.007581, RMSE:12.785702

精度悪すぎ… sofでも同じような質問を見つけたが同業だろう…w

ログを取ってみたところ、どうも収束が悪いっぽい。

import numpy as np
import tensorflow as tf
from sklearn import cross_validation
from sklearn import datasets
from sklearn.svm import SVR
from sklearn.ensemble import RandomForestRegressor
from sklearn.metrics import mean_squared_error
from sklearn.metrics import r2_score
from math import sqrt
import random

tf.logging.set_verbosity(tf.logging.ERROR)

def test_svm(x_train, x_test, y_train, y_test):
    clf = SVR(kernel='linear').fit(x_train, y_train)
    y_pred = clf.predict(x_test)
    r2 = r2_score(y_test, y_pred)
    rmse = sqrt(mean_squared_error(y_test, y_pred))
    print('SVM: R2: {0:f}, RMSE:{1:f}'.format(r2, rmse))

def test_rf(x_train, x_test, y_train, y_test):
    rlf = RandomForestRegressor().fit(x_train, y_train)
    y_pred = rlf.predict(x_test)
    r2 = r2_score(y_test, y_pred)
    rmse = sqrt(mean_squared_error(y_test, y_pred))
    print(' RF: R2: {0:f}, RMSE:{1:f}'.format(r2, rmse))

def test_tf_learn_dnn(x_train, x_test, y_train, y_test, steps=10000, hidden=[20, 20]):
    feature_columns = [tf.contrib.layers.real_valued_column("", dimension=13)]
    tfl = tf.contrib.learn.DNNRegressor(hidden_units=hidden, feature_columns=feature_columns)
    tfl.fit(x=x_train, y=y_train, steps=steps)
    y_pred = tfl.predict(x_test)
    r2 = r2_score(y_test, y_pred)
    rmse = sqrt(mean_squared_error(y_test, y_pred))
    print('TFL: R2: {0:f}, RMSE:{1:f}'.format(r2, rmse))

def inference(x):
    hidden1 = 20
    hidden2 = 20

    with tf.name_scope("l1") as scope:
        w1 = tf.Variable(tf.truncated_normal([13, hidden1]), name="w1")
        b1 = tf.Variable(tf.zeros([hidden1]), name="b1")
        h1 = tf.nn.relu(tf.matmul(x, w1) + b1)
        #h1 = tf.nn.dropout(h1, 0.9)

    with tf.name_scope("l2") as scope:
        w2 = tf.Variable(tf.truncated_normal([hidden1, hidden2]), name="w1")
        b2 = tf.Variable(tf.zeros([hidden2]), name="b2")
        h2 = tf.nn.relu(tf.matmul(h1, w2) + b2)
        #h2 = tf.nn.dropout(h2, 0.9)

    with tf.name_scope("l3") as scope:
        w3 = tf.Variable(tf.truncated_normal([hidden2, 1]), name="w2")
        b3 = tf.Variable(tf.zeros([1]), name="b3")
        y = tf.matmul(h2, w3) + b3
    return y

def loss(model, y_):
    return tf.reduce_mean(tf.square(tf.sub(model, y_)))

def training(loss, rate):
    return tf.train.AdagradOptimizer(rate).minimize(loss)

def test_my_dnn(x_train, x_test, y_train, y_test, batch_size=32, epoch=10000, shuffle=True):
    max_size = x_train.shape[0]
    n = max_size - batch_size
    idx = list(range(x_train.shape[0]))
    x = tf.placeholder(tf.float32, shape=[None, 13])
    y_ = tf.placeholder(tf.float32, shape=[None])

    model = inference(x)
    loss_value = loss(model, y_)
    train_op = training(loss_value, 0.1)

    init = tf.initialize_all_variables()
    sess = tf.Session()
    sess.run(init)

    for e in range(epoch):
        if shuffle:
            random.shuffle(idx)
            x_train = x_train[idx]
            y_train = y_train[idx]
        for i in range(n / batch_size):
            batch = batch_size * i
            x_train_b = x_train[batch:batch + batch_size]
            y_train_b = y_train[batch:batch + batch_size]
            _, l = sess.run([train_op, loss_value], feed_dict={x: x_train_b, y_: y_train_b})
        #if e % 100 == 0:
        #    print e, l

    y_pred = sess.run(model, feed_dict={x: x_test})
    y_pred = y_pred.T[0]
    r2 = r2_score(y_test, y_pred)
    rmse = sqrt(mean_squared_error(y_test, y_pred))
    print(' MY R2: {0:f}, RMSE:{1:f}'.format(r2, rmse))

if __name__ == "__main__":
    boston = datasets.load_boston()
    x_train, x_test, y_train, y_test = cross_validation.train_test_split(
        boston.data, boston.target, test_size=0.2, random_state=0)

    test_rf(x_train, x_test, y_train, y_test)
    test_svm(x_train, x_test, y_train, y_test)
    test_tf_learn_dnn(x_train, x_test, y_train, y_test)
    test_my_dnn(x_train, x_test, y_train, y_test, shuffle=True, epoch=30000)

というわけでソースコードを読んでみたんだけど、わかりにくかったので途中までしか追えてなく、結局原因がつかめていない状態。

以下メモ

  • 実際に使われているクラスはDNNLinearCombinedRegressorで、名前の由来がわからなかったけど、TensorFlow Wide & Deep Learning Tutorialを読めばわかる

  • fitで訓練するんだけど、_get_input_fnで入力、出力のplaceholderを返すinput_fnとバッチ用のデータセットを返すfeed_fnという関数が返される。

  • バッチは_get_input_fnでランダムシャッフルされているし、デフォルトのサイズは32

  • デフォルトのoptimizerはAdaGrad

  • clip-gradientsというオプションがあるのだがどこで使われているのかわからなかった

DNNをRandom Forest (RF)やSupport Vector Machine (SVM)と比較したい

TensorFlowのDNNチュートリアルだとトレーニングセットとテストセットをファイルから読みだすので、実用的にはちょっと面倒くさい。scikit-learnのよろしく分割してくれるメソッド使ったほうが楽でしょう。

またこScikit-learnとTensorFlowを組み合わせることでそれぞれのアルゴリズムの精度を比較することが簡単にできるので便利。

import tensorflow as tf
import numpy as np
from sklearn import datasets
from sklearn import svm
from sklearn.ensemble import RandomForestClassifier
from sklearn import cross_validation

iris = datasets.load_iris()
x_train, x_test, y_train, y_test = cross_validation.train_test_split(iris.data, iris.target, test_size=0.4, random_state=0)

classifier = tf.contrib.learn.DNNClassifier(hidden_units=[10, 20, 10], n_classes=3)
classifier.fit(x=x_train, y=y_train, steps=200)
dnn_accuracy_score = classifier.evaluate(x=x_test, y=y_test)["accuracy"]
print('DNN Accuracy: {0:f}'.format(dnn_accuracy_score))

clf = svm.SVC(kernel='linear').fit(x_train, y_train)
svm_accuracy_score = clf.score(x_test, y_test)
print('SVM Accuracy: {0:f}'.format(svm_accuracy_score))

rlf = RandomForestClassifier().fit(x_train, y_train)
rf_accuracy_score = rlf.score(x_test, y_test)
print('RF Accuracy: {0:f}'.format(rf_accuracy_score))