読者です 読者をやめる 読者になる 読者になる

EmacsでPythonとTeXとMarkdownを書いてる動画と使ってるパッケージ

はじめに

今年は頑張ってアウトプットを増やすことを目標にした。 とはいえ、特に書くこともなかったので、Emacsについて書く。

正月やることがなかったので、またEmacsの設定を見なおした。 Emacs Rocksという動画に触発されて、動画で説明することにした。 普段pythonTeXMarkdownしか書かないので、それらについての設定しかできていない。

環境

github.com

python

youtu.be

python modeはこのころ とあまり変わっていない。

  • yasnippet
    予め登録されているnpとかmainとかの後にtabを打つと、スニペットを挿入してくれる。 自分でスニペットを登録できるので、独自のライブラリに関するスニペットとかも作れて便利。

  • jedi
    jediはオムニ補完、つまり文法的な部分をある程度汲んだ上で補完をしてくれる。 動画ではやらなかったが関数定義にジャンプしたりもできる。 importした直後とかはたまに補完候補をだすのに数秒かかることもある。

  • flymake-python-pyflakes
    文法エラーやwarningの表示。 flymake-cursorを入れないとカーソル位置のエラーが表示されなかったようなそうでもなかったような気がする。

  • py-autopep8
    pep8準拠のコードに自動的に修正。 演算子の両側にスペースを入れたり、改行を入れたり。 after-save-hookに登録して、保存時に修正するといい感じ。

  • iedit
    複数箇所の同時編集。変数名を変えたいときに便利。 範囲を関数内に限定したり、手動で範囲を決めたりできる。 multiple-cursorsとどっちがいいんでしょうか。

  • smartparens
    括弧などの補完。 補完もそうだけど、範囲選択して括弧を入力すると、その範囲を括弧で囲ってくれるのがいい。

  • highlight-parentheses
    括弧の色を深さで変える。 趣味。

TeX

youtu.be

  • AUCTeX
    以前までYaTeXを使っていたけど、どちらがいいかわからなかったので、とりあえず使ってみている。 数式に関するショートカットがなくなったけど、yasnippet使うしどっちでもいい気がする。 C-c C-c でlatexmkを呼ぶようにできるらしいので、その点はこっちのほうがいいかも。

  • yasnippet
    機能はpythonと同じだが、複雑な構文が多いのでかなり役に立つ。 figureの構文とか覚えてられないし。 TeX固有の機能として、cite + Tabとかref+Tabでreftexを呼ぶようになっている。 デフォルトのままなので、自分の環境に合わせてもうちょっといじったほうがいいかなと思っている。

  • reftex
    bibファイルは使い回しなので、どういう名前で登録したか覚えていないし、 長い文章を書いていると、\labelで設定した名前なんか忘れる。 reftexは参考文献入れるときや\refを使うときに、名前を補完してくれる。 参考文献は特に著者名前でもキーワードでもヒットするので便利。
    同じbibファイルを使いまわすときは、default-bibliographyに登録しておくと便利。

  • magic latex buffer
    \beginとか\endを▽△に変えたり、数式を表示したり、見やすいように色々やってくれる。コマンドのタイポに気づきやすい。 数式以外の場所で"_"を打っても添字扱いされるところが少し困る。

Markdown

youtu.be

  • OrgTbl
    markdownで表を書いてると、自分で合わせない限り縦棒がずれているので見にくい。 orgの表を書くモードをMarkdownモードでも使えるので、それを使う。 orgとMarkdownで一部書式が違うので、保存時に変換している。

  • flymd
    Markdownのリアルタイムプレビューはいくつかパッケージがあるけど、今はflymdを使っている。 画像や数式もちゃんと表示されるのでいい。

  • misc
    研究ノート的なものをmarkdownで書いているので、C-x C-nで常に同じファイルが開くようになっている。 org-modeでいいじゃないかという気もしている。

その他

git-gutter-fringeの見た目をちょっといじって使っている。 いい感じ。

f:id:ksknw:20170108204246p:plain

今後の課題

未だに改行をEnterでしか打てないので、C-mとかC-jとか使っていきたい。

参考

関係データ学習の実装 ツイッターデータのスペクトルクラスタリングとSBM

概要

関係データ学習の学習のために,自分で実装して理解する. ツイッターのフォローフォロワー関係を使って,グラフラプラシアンを求めスペクトルクラスタリングを行った. その結果,なんとなくクラスタリングできた. また,確率的ブロックモデルによる非対称データクラスタリングをStanによって実装しようとした. これはうまくいっていない.

はじめに

関係データ学習という本を買って読んでいる.

www.kspub.co.jp

本の内容は前半と後半に分かれていて,前半は関係データをスペクトルクラスタリングしたり,確率的ブロックモデルでクラスタリングしたりする話.後半は行列分解やテンソル分解の話になっている. まだ前半の途中までしか読めていないが,予想していたよりも数式が簡単だったこともあり,実際のデータに適用してみたくなった. 数年前に書いたツイッターのフォローフォロワー関係をダウンロードするスクリプトを見つけたので,これを使ってデータを作り,自分のフォロー中の人たちをクラスタリングしてみることにした. (カット最小化とスペクトルクラスタリングの関係とか,まだよく理解できてないので,理論的な説明はしない.)

ツイッターデータのダウンロード

tweepyを使って関係データをダウンロードする. スクリプトは何年も前に書いたものであり,なんとなく恥ずかしいのでここでは省略するが,

  • 自分がフォローしている人のフォローしている人を全てダウンロード.
  • ダウンロードした中から自分がフォローしている人たちの関係を抜き出す.

ようになっている. ツイッターapi制限があるので,100人ちょっとだけど,それなりに時間がかかる.

ダウンロードしたのは以下のようなデータ。 つながっている部分が青、ない部分が白。鍵垢のデータはとれてない部分もある。 ツイッターではフォローしているがされていないなどの関係もあるので,非対称データになる.

import pandas as pd

data = pd.read_csv("./combinationTable.csv")
uname = data[:1].get_values()[0]

%matplotlib inline
import pylab as plt
import seaborn
import numpy as np

X = (data.get_values()[1:]=="True")

def plot_matrix(matrix, clusters=None):
    f, ax = plt.subplots(figsize=(10, 10))
    plt.pcolor(matrix, cmap=plt.cm.Blues )

    if not clusters is None:
        clusters_diff = np.r_[[0], np.diff(clusters)]
        print clusters_diff
        for i,c in enumerate(clusters_diff):
            if c==1:
                ax.axhline(i, c="grey") #, linewidth=1)
                ax.axvline(i, c="grey")#, linewidth=1)
    plt.savefig("../../temp.svg")
    plt.show()

plot_matrix(X)

f:id:ksknw:20161225112150p:plain

すでにある程度のクラスタができているように見える. これは自分が今まで所属した(ツイッター的意味での)クラスタを表している. 自分がフォローした順にアカウントが並んでいるので,右上からだいたいポケモン,大学,大学研究室,研究者みたいな感じ. 研究者クラスタ以外は知り合いをフォローしていることが多いので,数が少ない.

はじめはこの関係データを対称データに変換し,2章を参考にしつつスペクトルクラスタリングを行う. 以下のように対称データに変換した.

symmetry_X = (X.T + X)
plot_matrix(symmetry_X)

f:id:ksknw:20161225112224p:plain

非正規化グラフラプラシアンに基づくスペクトルクラスタリング

対称データ{ \displaystyle
X}クラスタ{ \displaystyle
K}が与えられたとき、非正規化グラフラプラシアン{ \displaystyle
L}は以下のように定義される。

{ \displaystyle
L=D-X}

ここで{ \displaystyle
D}は次数行列で、オブジェクトから伸びるリンクの総数。 { \displaystyle
L}をプロットすると以下。

D = np.eye(symmetry_X.shape[0]) * symmetry_X.sum(axis=0)
K = 10
L = D - X
plot_matrix(L)

f:id:ksknw:20161225125154p:plain

先ほどまでとはグラフの色が異なるので注意. 対角成分だけ1以上の値をとるのでこんな感じの図になる。

スペクトルクラスタリングでは、以下の手順でクラスタを抽出する。

クラスタ数は適当に10とした。本来なら色々チューニングしないといけないが、面倒なので固定する。

import scipy.linalg

hi = K-1
lo = 0
eigen_value,eigen_vector = scipy.linalg.eigh(L,eigvals=(lo,hi))
def plot_eigen_vec(eigen_vector):
    plt.figure(figsize=(30,10))
    plt.plot(np.abs(eigen_vector))
    plt.show()

from sklearn.cluster import KMeans
def kmeans(eigen_vector, unames, n_clusters):
    kmean = KMeans(n_clusters=n_clusters)
    clusters = kmean.fit_predict(eigen_vector)
    cluster_uname = zip(clusters, uname)
    cluster_uname.sort(key=lambda x:x[0])
    print cluster_uname
    return clusters

def sort_by_cluster(matrix, clusters):
    sorted_mat = zip(clusters, matrix)
    sorted_mat.sort(key=lambda x:x[0])
    _,sorted_mat = zip(*sorted_mat)
    
    sorted_mat = zip(clusters, np.array(sorted_mat).T)
    sorted_mat.sort(key=lambda x:x[0])
    sorted_clusters,sorted_mat = zip(*sorted_mat)
    return np.array(sorted_mat).T, sorted_clusters
plot_eigen_vec(eigen_vector)
clusters = kmeans(eigen_vector, uname, K)
plot_matrix(*sort_by_cluster(symmetry_X, clusters))
sorted_mat, sorted_clusters = sort_by_cluster(symmetry_X, clusters)

クラスタリング結果は以下。 クラスタの番号順に並び替えて、境界に線を引いた。

f:id:ksknw:20161225125552p:plain

でかいクラスタが1つと1アカウントしかないクラスタがいくつかできた。 正規化しない単純なグラフラプラシアンだとこういう結果になるらしい。

対称正規化グラフラプラシアンに基づくスペクトルクラスタリング

もう少し望ましいクラスタリングを得るために、グラフラプラシアンの正規化を行う。 これは(緩和した)グラフカットの正規化と等価になる(?)らしいが、細かい説明は本(とその参考文献)に譲る。 読まねば。

対象正規化グラフラプラシアン{ \displaystyle
L_\mathrm{sym}}は、上で定義したグラフラプラシアン{ \displaystyle
L}から以下のように求まる。

{ \displaystyle
L_\mathrm{sym} =D^{-\frac{1}{2}} L D^{-\frac{1}{2}}}

L = D - X
L_sym = np.dot(np.dot(scipy.linalg.sqrtm(np.linalg.inv(D)) , L) , scipy.linalg.sqrtm(np.linalg.inv(D))) # D**(-0.5) * L * D**(-0.5)
plot_matrix(L_sym)
eigen_value,eigen_vector = scipy.linalg.eigh(L_sym,eigvals=(lo,hi))
eigen_vector = eigen_vector/np.linalg.norm(eigen_vector)

対象正規化グラフラプラシアンのプロットは以下。 非正規化グラフラプラシアンに比べて対角項が小さくなっている。 また、クラスタリングする前に固有ベクトルを正規化している。

f:id:ksknw:20161225144212p:plain

plot_eigen_vec(eigen_vector)
clusters = kmeans(eigen_vector, uname, K)
plot_matrix(*sort_by_cluster(symmetry_X, clusters))

クラスタリング結果は以下。 先程に比べればだいぶ良さそうだが、真ん中の巨大なクラスタにほとんど周りとつながっていないものがいくつか含まれている。

f:id:ksknw:20161225144619p:plain

酔歩正規化グラフラプラシアン

別の正規化から導かれるグラフラプラシアンを使ってクラスタリングをする。

{ \displaystyle
L_\mathrm{rw} =D^{-1} L }

L = D - X
L_rw = np.dot(np.linalg.inv(D) , L)  # D**(-1) * L 

plot_matrix(L_rw)
eigen_value,eigen_vector = scipy.linalg.eigh(L_rw,eigvals=(lo,hi))

plot_eigen_vec(eigen_vector)
clusters = kmeans(eigen_vector, uname, K)
plot_matrix(*sort_by_cluster(symmetry_X, clusters))

f:id:ksknw:20161225145741p:plain

でかいクラスタが少し分割されて、まあこっちのほうがいいかなぁという感じ。

固有ベクトルをPCAで2次元に落としてプロットする。 左でポケモンクラスタが孤立している。 右に研究室のクラスタがきた。ドクターの人と先生のアカウントは研究者のクラスタとエッジをたくさんもっている。 研究室が同じかつ大学の同級生は大学同級生クラスタと研究室クラスタの間に配置された。 このような感じで、まあまあ妥当なクラスタリングができた。

from sklearn.decomposition import PCA
# from sklearn.manifold import TSNE
# decomp = TSNE()
decomp = PCA()
decomped = decomp.fit_transform(eigen_vector)

%matplotlib inline
plt.figure(figsize=(20,20))
for i,row in enumerate(X):
    for j,is_connected in enumerate(row):
        if is_connected:
            plt.plot(decomped[[i,j], 0], decomped[[i,j], 1], alpha=0.2,linewidth=0.3, c="gray")
plt.scatter(decomped[:,0], decomped[:,1], c=clusters, linewidths=0, cmap=plt.cm.jet)
#for d,u in zip(decomped, uname):
#    plt.text(d[0], d[1], u)
plt.show()

f:id:ksknw:20161229121731p:plain

確率的ブロックモデルによる非対称データクラスタリング

ここまでは、本来非対称データであるツイッターの関係データを対称データとしてクラスタリングした。 ここからは、確率的ブロックモデル(SBM)を使って非対称データのままクラスタリングする。

{ \displaystyle
i} 番目のアカウントが { \displaystyle
j} 番目のアカウントをフォローしているかどうかを { \displaystyle
x_{i,j}
} と書く。 また、{ \displaystyle
z_{1,i}}{ \displaystyle
i}番目のアカウントが所属するフォローする側のクラスタ{ \displaystyle
z_{2,i}}をフォローされる側のクラスタとする。 フォローする側のクラスタとか意味がわからないかもしれないが、例えば機械学習関連のアカウントをフォローしまくっているクラスタとかそういう感じだと思う。

このときSBMでは以下の生成モデルを仮定する。

{ \displaystyle
\pi_1 \mid \alpha_1 \sim \mathrm{Dirichlet}(\alpha_1) }
{ \displaystyle
\pi_2 \mid \alpha_2 \sim \mathrm{Dirichlet}(\alpha_2)}
{ \displaystyle
z_{1,i} = k\mid \pi_1 \sim \mathrm{Discrete}(\pi_1)}
{ \displaystyle
z_{2,i} = l\mid \pi_2 \sim \mathrm{Discrete}(\pi_2)}
{ \displaystyle
\theta_{k,l} \mid a_0,b_0 \sim \mathrm{Beta}(a_0,b_0)}
{ \displaystyle
x_{i,j} \mid {\theta_{k,l}} , z_{1,i}, z_{2,j} \sim \mathrm{Bernoulli}(\theta_{z_{1,i}, z_{2,j}})
}

SBMを使うと、「フォローする(1)」と「フォローされる(2)」にわけてクラスタを考えることになる。 例えば「フォロワーするクラスタ」は有名な研究者をフォローするクラスタだが、「フォローされるクラスタ」は研究室の人とかそういうことになる。 うっ

SBMでは、クラスタ間の関係の強さ  \theta_{k,l} とそれぞれのアカウントがどのクラスタに所属するか  z_{1,i} z_{2,i}を同時に推定する。 本では周辺化ギブスサンプラーによる推論方法の導出を行っている。 ここでは、Stanを用いて、SBMを解くことを試みる。 まずは、上の数式をそのままStanに変換すると以下のようになる。

data{
  int N;
  int K;
  int x[N, N];
}

parameters{
  vector<lower=0>[K] alpha1;
  vector<lower=0>[K] alpha2;
  int z1[N];
  int z2[N];
  vector[K] pi1;
  vector[K] pi2;
  real<lower=0> a0;
  real<lower=0> b0;
  matrix[K,K] theta;
}

model{
  pi1 ~ dirichlet(alpha1);
  pi2 ~ dirichlet(alpha2);

  for (i in 1:N){
    z1[i] ~ categorical(pi1);
    z2[i] ~ categorical(pi2);
  }

  for (k in 1:K){
    for(l in 1:K){
      theta[k,l] ~ beta(a0, b0);
    }
  }

  for (i in 1:N){
    for (j in 1:N){
      x[i,j] ~ bernoulli(theta[z1[i], z2[j]]);
    }
  }
}

これには色々問題があるが、特にint型の変数がparametersの中に含まれていることが問題である。 Stanでは離散的なパラメータを扱うことができないので、zで周辺化する必要がある。 zで周辺化したり、遅すぎないようにベクトルにしたりで以下のようになる。(たぶんあっていると思う。)間違っている。

data{
  int N;
  int K;
  int x[N*N];
}

parameters{
  vector<lower=0>[K] alpha1;
  vector<lower=0>[K] alpha2;
  simplex[K] pi1;
  simplex[K] pi2;
  real<lower=0> a0;
  real<lower=0> b0;
  vector<lower=0,upper=1>[K*K] theta;
}

model{
  real gamma[K*K];

  pi1 ~ dirichlet(alpha1);
  pi2 ~ dirichlet(alpha2);
  theta ~ beta(a0, b0);

  for (k in 1:K){
    for (l in 1:K){
      gamma[k+ K*(l-1)] = log(pi1[k]) + log(pi2[l]) + bernoulli_lpmf(x | theta[k+ K*(l-1)]);
    }
  }
  target += log_sum_exp(gamma);
}

サンプリングする。

import pystan
data_dict= {"N":len(X),
            "K":10,
            "x":X.astype(int).reshape(-1)}
fit = pystan.stan("sbm_vector.stan", data=data_dict)
fit
            mean se_mean     sd   2.5%    25%    50%    75%  97.5%  n_eff   Rhat
alpha1[0]  2.6e7   3.1e7  5.3e7 3429.6 4934.1  3.9e4  1.2e7  1.8e8    3.0    3.2
alpha1[1]  3.2e8   3.5e8  5.0e8  3.0e5  1.4e6  3.7e7  5.4e8  1.4e9    2.0   6.44
alpha1[2]  1.0e7   6.9e6  9.7e6 460.13  1.1e6  6.9e6  1.5e7  3.4e7    2.0    2.7
alpha1[3]  2.9e7   3.6e7  5.0e7 9752.8  2.5e4  7.6e5  4.8e7  1.4e8    2.0   9.99
alpha1[4]  9.3e6   1.2e7  1.6e7   1.29  5.7e4  3.8e5  9.9e6  5.1e7    2.0   5.75
...
theta[97]   0.58    0.14   0.25   0.16    0.3    0.6   0.74   0.98    3.0   3.63
theta[98]   0.29    0.04   0.16   0.02   0.23   0.23   0.35   0.78   17.0   1.12
theta[99]   0.35    0.09   0.16   0.11   0.19   0.34   0.54   0.57    3.0   2.23
lp__       -4015   40.35  57.06  -4105  -4053  -4019  -3975  -3929    2.0  12.14
fit.plot()

f:id:ksknw:20161229165851p:plain

うまく収束しなかった。 特に対策をしていないために、ラベルスイッチングが起こっているのが原因かなと思ったので、  \thetaの片側をorderに変えた。 0~1の制約も残さないといけないので、すこしごちゃごちゃした。

data{
  int N;
  int K;
  int x[N*N];
}

parameters{
  vector<lower=0>[K] alpha1;
  vector<lower=0>[K] alpha2;
  simplex[K] pi1;
  simplex[K] pi2;
  real mu;
  real<lower=0> sigma;
  ordered[K] theta_l[K];
}

transformed parameters{
  vector<lower=0, upper=1>[K] theta[K];
  for (k in 1:K){
    for (l in 1:K){
      theta[k][l] = inv_logit(theta_l[k][l]);
    }
  }
}

model{
  real gamma[K*K];
  pi1 ~ dirichlet(alpha1);
  pi2 ~ dirichlet(alpha2);
  for (k in 1:K)
    theta_l[k] ~ normal(mu, sigma);

  for (k in 1:K){
    for (l in 1:K){
      gamma[k+ K*(l-1)] = log(pi1[k]) + log(pi2[l]) + bernoulli_lpmf(x | theta[k][l]);
    }
  }
  target += log_sum_exp(gamma);
}

サンプリングする。

fit2 = pystan.stan("sbm2.stan", data=data_dict)
fit2
               mean se_mean     sd   2.5%    25%    50%    75%  97.5%  n_eff   Rhat
alpha1[0]     2.3e7   1.3e7  2.3e7 6866.1  1.7e6  1.7e7  3.5e7  7.8e7    3.0   2.35
alpha1[1]     4.4e7   2.0e7  3.4e7  2.1e4  7.3e6  4.3e7  6.1e7  1.1e8    3.0    2.5
alpha1[2]     3.8e8   5.0e8  7.0e8  2.4e5  1.4e6  4.0e6  4.0e8  2.2e9    2.0   2.96
alpha1[3]     3.6e7   3.6e7  5.1e7  5.1e6  7.1e6  9.3e6  3.9e7  1.7e8    2.0   3.64
alpha1[4]     6.1e7   6.1e7  1.1e8  5.6e5  1.1e6  8.5e6  5.7e7  3.5e8    3.0   4.63
alpha1[5]     3.9e7   3.6e7  5.1e7 2659.3 9569.2  8.9e6  6.2e7  1.5e8    2.0   4.23
...
theta[8,9]     0.95    0.07   0.13   0.51    1.0    1.0    1.0    1.0    4.0   1.66
theta[9,9]     0.93    0.08   0.16   0.38    1.0    1.0    1.0    1.0    4.0    1.8
lp__          -4060    6.59  13.17  -4087  -4070  -4060  -4049  -4039    4.0   1.98

が、これもだめ。

諦めて一端ここで終わり。

おわりに

関係データ学習の学習のために,対称データのスペクトルクラスタリングと非対称データのSBMを実装した。 スペクトルクラスタリングの結果、なんとなくクラスタを抽出できた。 SBMをStanで実装しようとしたが、いまいち収束しなかった。

とりあえず、最後まで読んで、頑張れたら頑張ります。

参考

[http://www.kspub.co.jp/book/detail/1529212.html:title]

追記(2017年1月5日)

上記はたぶん間違っているので、以下のようなモデルを作った。 しかし、これもうまくいっていない。

data{
  int N;
  int K;
  int x[N,N];
}

parameters{
  simplex[K] pi1;
  simplex[K] pi2;
  real<lower=0> a0;
  real<lower=0> b0;
  /* vector<lower=0>[K] alpha1; */
  /* vector<lower=0>[K] alpha2; */
  matrix<lower=0, upper=1>[K,K] theta;
}

transformed parameters{
  matrix<upper=0>[K,K] soft_z[N, N];
  for (i in 1:N){
    for (j in 1:N){
      for (k in 1:K){
        for (l in 1:K){
          soft_z[i,j, k,l] = log(pi1[k]) + log(pi2[l]) + bernoulli_lpmf(x[i,j] | theta[k,l]);
        }
      }
    }
  }
}

model{
  /* pi1 ~ dirichlet(alpha1); */
  /* pi2 ~ dirichlet(alpha2); */
  to_vector(theta) ~ beta(a0, b0);

  for (i in 1:N){
    for (j in 1:N){
      target += log_sum_exp(soft_z[i,j]);
    }
  }
}

Stanで気温変動のフィッティング

はじめに

この記事は Stan Advent Calendar 2016n日目の記事ではありません。 21日目の記事です。

話題のアヒル本(言われるまで何がアヒルなのかわからなかった)を読んだ。

statmodeling.hatenablog.com

Stanの文法を知るために買ったけれど、全体的に統計モデリング、解析について書いてあって、勉強になった。 本を読んで勉強になったことはたくさんあるけど、特にグラフを死ぬほど書かないといけないという認識を持った。 自分ではそれなりにグラフを書いているつもりだったけれど、pariplotとかviolinplotとかちゃんと書いたことがなかった。

せっかく本を読んだので、適当にデータをとってきて解析してAdvent Calendarなるものに挑戦しようと思っていたけど、もう埋まっていた。悲しいので一人で勝手にやることにする。 よく見たら21日が空いていたので、フライングということにさせてください。

今回は特に12章を見ながら時系列データをいじくることにした。

知識がないので、間違ったことをしていたら教えていただけると嬉しいです。

データ

時系列データとして天気に関するデータを使うことにした。 気象庁のサイト から天気に関するいろんなデータがダウンロードできるので利用する。

今回は愛知県の気温とか降水量とかを過去5年分ぐらいダウンロードした。 日本語のファイルがやや面倒だった。

import pandas as pd
aichi = pd.read_csv("./weather_aichi.csv", encoding='SHIFT-JIS')
aichi["年"] = aichi[u"年月日"].apply(lambda x: x.split("/")[0])
aichi["月_s"] = aichi[u"年月日"].apply(lambda x: "%02d"%int(x.split("/")[1]))
aichi["日_s"] = aichi[u"年月日"].apply(lambda x: "%02d"%int(x.split("/")[2]))
aichi["m"] = aichi[u"年月日"].apply(lambda x: int(x.split("/")[1]))
aichi["d"] = aichi[u"年月日"].apply(lambda x: int(x.split("/")[2]))
aichi["is_rain"]  = aichi[u"降水量の合計(mm)"]>0
aichi = aichi.set_index(u"年月日")
month_groups = aichi.groupby("月_s")
aichi
Unnamed: 0 平均気温(℃) 降水量の合計(mm) 月_s 日_s m d is_rain
年月日
2011/12/3 1 12.2 8.5 2011 12 03 12 3 True
2011/12/4 2 11.9 0.0 2011 12 04 12 4 False
2011/12/5 3 10.4 0.0 2011 12 05 12 5 False
2011/12/6 4 8.0 0.5 2011 12 06 12 6 True
... ... ... ... ... ... ... ... ... ...
2016/11/30 1825 8.2 0.0 2016 11 30 11 30 False
2016/12/1 1826 12.0 0.0 2016 12 01 12 1 False
2016/12/2 1827 11.0 0.0 2016 12 02 12 2 False
2016/12/3 1828 9.8 0.0 2016 12 03 12 3 False

1828 rows × 9 columns

今回は、平均気温をモデリングしてみることにした。 とりあえずグラフを書く。

%matplotlib inline 
import pylab as plt
import seaborn as sns
import matplotlib.font_manager as fm
font = {'family' : 'TakaoGothic'}
plt.rc('font', **font)

ax = aichi[[ u"平均気温(℃)"]].plot(figsize=(20,10))

f:id:ksknw:20161211134350p:plain

当たり前だけど春夏秋冬がある。 雨のデータも後で使うのでプロットする。

ax = aichi[ u"平均気温(℃)"]["2011/12/3": "2012/12/3"].plot(figsize=(20,10))
aichi[ u"降水量の合計(mm)"]["2011/12/3": "2012/12/3"].plot( kind="bar", ax=ax)
plt.show()

f:id:ksknw:20161211134429p:plain

その他色々プロットを書く。

sns.pairplot(aichi[[u"平均気温(℃)",   u"降水量の合計(mm)", "m"]], hue="m", size=4)

f:id:ksknw:20161211134443p:plain

アヒル本では対角線の上と下で別のプロットを書いていたけど、正直手でやるのはちょっと面倒。 変数の型とか見て自動でやるやつをそのうち作りたい。

fig, ax = plt.subplots(figsize=(20,20))
sns.violinplot(data=aichi, x="m", y =u"平均気温(℃)")
fig, ax = plt.subplots(figsize=(20,20))
sns.violinplot(data=aichi, x="m", y =u"平均気温(℃)", hue="is_rain", split=True)
plt.show()

f:id:ksknw:20161211134455p:plain

f:id:ksknw:20161211134516p:plain

雨が降ると夏は気温が下がるけど、冬は逆に気温が高い。そりゃそうだな。

temperatures = aichi[u"平均気温(℃)"].get_values()
data = {"Y":temperatures,
       "N":len(aichi),
       "is_rain":aichi["is_rain"].astype(int).get_values()+1}

状態空間モデルでのフィッティング

アヒル本12章を見ながら、いくつかのモデルを使って気温の変動をモデル化する。

状態空間モデル(1階差分)

観測された気温({ \displaystyle
Y
})は 真の気温({ \displaystyle
\mu
})を平均とする正規分布に従い、 真の気温は、前日の真の気温を平均とする正規分布に従うというモデルを仮定する。 式的には以下。

{ \displaystyle
 \mu_n \sim \mathcal{N}(\mu_{n-1}, s_{\mu})
}

{ \displaystyle
 Y_n \sim \mathcal{N}(\mu_{n}, s_{Y})
}

以下のようにStanでモデルを書ける。

data{
  int N;
  vector[N] Y;
}

parameters{
  vector[N] mu;
  real<lower=0> s_mu;
  real<lower=0> s_Y;
}

model{
  mu[2:N] ~ normal(mu[1:(N-1)], s_mu);
  Y ~ normal(mu, s_Y);
}

python からStanをよんでサンプリングして、色々とプロットする。

import pystan
import numpy as np
def fit(model_filename):
    fit = pystan.stan(model_filename, data=data)
    ms = fit.extract()
    fit.plot()
    return fit, ms

def fit_plot(fit, ms, list_max_x=[-1,500,100,50], ylim=None):
    mu_mean = ms["mu"].mean(axis=0)
    mu_std  = ms["mu"].std(axis=0)

    def plot_line(max_x):
        fig = plt.figure(figsize=(10,10))
        plt.plot(mu_mean[:max_x], linewidth=0.5)
        
        plt.scatter(range(len(temperatures[:max_x])), temperatures[:max_x], marker=".",# linestyle="None", 
                    c=aichi["is_rain"].get_values()[:max_x].astype(float))
    
        plt.fill_between(range(len(mu_mean[:max_x])), (mu_mean-mu_std)[:max_x], (mu_mean+mu_std)[:max_x], alpha=0.3)
        plt.fill_between(range(len(mu_mean[:max_x])), (mu_mean-mu_std*2)[:max_x], (mu_mean+mu_std*2)[:max_x], alpha=0.2)

        plt.xlim(0, len(temperatures[:max_x]))
        if not ylim==None:
            plt.ylim(*ylim)
        plt.show()
   
    for max_x in list_max_x:
        plot_line(max_x)
fit1 = fit("./model.stan") 
fit_plot(*fit1, list_max_x=[-1,100])

f:id:ksknw:20161211140401p:plain

1000のあたりとか、見るからにサンプリングが安定していない感じがある。 サンプリング結果を見てみる。

fit1[0]
           mean se_mean     sd   2.5%    25%    50%    75%  97.5%  n_eff   Rhat
mu[0]     12.18  6.1e-3   0.39  11.34  11.99  12.19  12.36  12.98 4000.0    1.0
mu[1]     11.84  6.0e-3   0.38  10.98  11.66  11.87  12.04  12.63 4000.0    1.0
mu[2]     10.36  6.1e-3   0.38   9.49  10.16  10.39  10.54  11.15 4000.0    1.0
mu[3]      8.21  6.5e-3   0.41   7.47   7.96   8.12   8.46   9.13 4000.0   1.05
mu[4]      9.99  6.2e-3   0.39    9.1   9.77  10.04   10.2  10.74 4000.0   1.02

mu[1825]  11.77  6.5e-3   0.41  10.85  11.53  11.87  12.04  12.48 4000.0   1.07
mu[1826]  10.97  6.1e-3   0.39  10.14  10.78  10.99  11.18  11.76 4000.0    1.0
mu[1827]   9.86  6.3e-3    0.4   9.08   9.67   9.83  10.06  10.75 4000.0    1.0
s_mu       1.72    0.04   0.07   1.59   1.67   1.72   1.78   1.85    3.0   1.71
s_Y        0.37    0.12   0.18   0.07    0.2   0.42   0.49   0.63    2.0   3.15
lp__     -622.4  903.47 1277.7  -1830  -1483  -1218 121.56 1970.3    2.0   5.28

アヒル本によると、

n_effはStanが自己相関等から判断した実効的なMCMCサンプル数である。

とあるので、今回のサンプリングではs_muやs_Yはまともにサンプリングされていないことがわかる。

また、Rhatはパラメータが収束したかどうかを表す値であり、

「chain数が3以上ですべてのパラメータでRhat < 1.1 となること」を「収束した」と見なすことにする。

とある。 ということで、このモデルではうまく収束しなかった。 サンプリング数を増やしたりしてもいいのかもしれないが、今回は別のモデルを作ることにする。

状態空間モデル(2階差分)

真の気温の変化はその前の日の真の気温の変化を平均とする正規分布に従うというモデル。 式的には以下。

{ \displaystyle
 \mu_n-\mu_{n-1} \sim \mathcal{N}(\mu_{n-1}-\mu_{n-2}, s_{\mu})
}

{ \displaystyle
 Y_n \sim \mathcal{N}(\mu_{n}, s_{Y})
}

Stanのコードは以下。

data{
  int N;
  vector[N] Y;
}

parameters{
  vector[N] mu;
  real<lower=0> s_mu;
  real<lower=0> s_Y;
}

model{
  mu[3:N] ~ normal(2 * mu[2:(N-1)] - mu[1:(N-2)], s_mu);
  Y ~ normal(mu, s_Y);
}

サンプリングする。

fit2 = fit("./model2.stan") 
fit_plot(*fit2, list_max_x=[-1,100])

f:id:ksknw:20161211140436p:plain

サンプリング結果

fit2[0]
           mean se_mean     sd   2.5%    25%    50%    75%  97.5%  n_eff   Rhat
mu[0]     12.12    0.02   1.06  10.07  11.42  12.13  12.83  14.17 4000.0    1.0
mu[1]     11.26    0.01   0.76   9.76  10.75  11.25  11.78  12.75 4000.0    1.0
mu[2]     10.41    0.01   0.67   9.12   9.95  10.41  10.86  11.73 4000.0    1.0
mu[3]      9.64    0.01   0.68   8.28   9.19   9.65  10.09  10.97 4000.0    1.0
mu[4]      9.04    0.01   0.68    7.7   8.59   9.04   9.49  10.38 4000.0    1.0

mu[1825]  10.28    0.01   0.67   8.99   9.82  10.28  10.73  11.59 4000.0    1.0
mu[1826]   10.4    0.01   0.75   8.94    9.9  10.38   10.9  11.89 4000.0    1.0
mu[1827]  10.42    0.02   1.04   8.42   9.72  10.42  11.14  12.39 4000.0    1.0
s_mu       0.49  9.5e-3   0.06   0.35   0.45   0.49   0.53   0.61   42.0    1.1
s_Y        1.43  7.5e-3   0.06   1.33   1.39   1.43   1.46   1.55   55.0   1.07
lp__      -1158   28.55 180.58  -1459  -1281  -1171  -1060 -711.0   40.0   1.09

Rhatの値を見てみる1.1ギリギリという感じ。 n_effも50とかしかないので、厳しそう。

一応{ \displaystyle \mu}ベイズ信頼区間のグラフを見る。 白丸が晴れの日、黒丸は雨の日を表す。 濃い青は68%ぐらい、薄い青は95%ぐらい。 青線は平均値を表す。

f:id:ksknw:20161211181516p:plain

なんとなく雨の日が上手くあっていないような気がする。 雨が降ると前日との温度差が大きくなることが予想されるので、そのへんのモデルをちゃんと作ったほうが良さそう。

状態空間モデル(2階差分 降水量を考慮するモデル)

はじめの方に描いたグラフなどから天気を考慮する必要がありそう。 そこで、この題材に対するモチベーションとモデルの当てはまりから以下のようなモデルを仮定する。

{ \displaystyle
 \mu_n-\mu_{n-1} \sim \mathcal{N}(\mu_{n-1}-\mu_{n-2}, s_{\mu})
}

雨のとき

{ \displaystyle
      Y_n \sim \mathcal{N}(\mu_{n}, s_{Y|\mathrm{rain}})
}

雨がふらなかったとき

{ \displaystyle
      Y_n \sim \mathcal{N}(\mu_{n}, s_{Y|\mathrm{no\ rain}})
}

雨が振った日と振らなかった日では、真の気温から観測値へ乱数の分散が変わるというモデル 本当は分散だけでなく平均も変わったほうがいい気がするが、その変わり方が季節ごとに違うので面倒でやめた。

fit3 = fit("./model3.stan") 
fit_plot(*fit3, list_max_x=[-1,100])

f:id:ksknw:20161211181548p:plain

fit3[0]
       mean se_mean     sd   2.5%    25%    50%    75%  97.5%  n_eff   Rhat
mu[0]     11.41    0.04   2.55   6.51   9.72  11.38  13.13   16.4 4000.0    1.0
mu[1]     10.83    0.03   1.59   7.72   9.74  10.84   11.9  13.97 4000.0    1.0
mu[2]     10.26    0.01   0.93   8.46   9.64  10.25  10.89  12.09 4000.0    1.0
mu[3]      9.71    0.01   0.72   8.32   9.22    9.7  10.18  11.11 4000.0    1.0

mu[1826]  10.52  9.9e-3   0.63   9.31  10.08  10.52  10.96  11.77 4000.0    1.0
mu[1827]  10.38    0.01   0.88   8.65   9.78  10.39  10.99  12.07 4000.0    1.0
s_mu       0.69  7.5e-3   0.07   0.55   0.64   0.69   0.74   0.82   88.0   1.03
s_Y[0]     1.06  5.5e-3   0.06   0.95   1.01   1.05    1.1   1.18  123.0   1.02
s_Y[1]      1.7  2.9e-3   0.07   1.57   1.65    1.7   1.75   1.84  583.0   1.01
lp__      -1508    13.2 122.38  -1712  -1595  -1520  -1433  -1247   86.0   1.03

前のモデルよりちゃんと収束してる感はある。 s_Y[1]が雨が降った日なので、雨の日は晴れの日よりも、真の気温からの分散が大きくなる。

前のモデルとの違いがぱっとわからないので、重ねてグラフをかく。

fig = plt.figure(figsize=(10,10))
max_x = 100
colors = ["b", "g"]
for i,ms in enumerate([fit2[1], fit3[1]]): 
    mu_mean = ms["mu"].mean(axis=0)
    mu_std  = ms["mu"].std(axis=0)

    plt.plot(mu_mean[:max_x], linewidth=0.5)    
    plt.scatter(range(len(temperatures[:max_x])), temperatures[:max_x], marker=".",# linestyle="None", 
                c=aichi["is_rain"].get_values()[:max_x].astype(float))
    
    #plt.fill_between(range(len(mu_mean[:max_x])), (mu_mean-mu_std)[:max_x], (mu_mean+mu_std)[:max_x], alpha=0.3)
    plt.fill_between(range(len(mu_mean[:max_x])), (mu_mean-mu_std*2)[:max_x], (mu_mean+mu_std*2)[:max_x], alpha=0.2
                    ,color=colors[i])

plt.xlim(0, len(temperatures[:max_x]))
plt.ylim([-5, 20])
plt.show()

f:id:ksknw:20161211141729p:plain

青が単純な2階差分、緑が雨で条件わけしたモデル。 正直良くなっているのかどうかよくわからない。

まとめと今後

アヒル本を参考にして、気温変動のモデリングを行った。 Stanの使い方はなんとなくわかってきたけど、 統計的に正しいのかとか、いまいちよくわかっていない。

今回は年ごとの関係を考慮しないでモデルを作った。 sinとかでモデル化できそうな気がする。 また、雨の日と晴れの日でsinの強さが変わるようにすると、もっといい感じにモデル化できるように思う。

Stanでまだわかってないことの1つにディリクレ混合過程がある。 アヒル本には載っていなかったけど、このへんを見ながらやってみたい。

fit.plot()のfig_sizeを大きくしたいんだけど、どうやらpython側でなくstan側からプロットされている(?)ようで、いまいちやり方がわからなかった。

(この記事には、グラフを死ぬほどかかないといけないという認識があまり反映されていないのではないか。)

参考

jupyter notebookで触れるプロットを描く

昔、Emacsから使えないからjupyter notebookは使わないという旨の記事を書いた。 実際まだ、Emacsからは使えていないけど、最近はだいたいのコードをEmacs+pythonで書いて、実行やグラフを描く部分をjupyter notebookでやるということをよくやっている。

jupyter notebook内にグラフを描く方法として

%matplotlib inline

という文を実行しておくという方法がある

大体の場合にはこれで間に合うが、3Dグラフを回したり、アニメーションを動かしたりすることができなかった。 調べてみると他にも色々できることがわかったのでまとめる。 この記事では以下のようなプロットを描く。

  • 回転できる3D
  • アニメーション
  • スライドバー付き
  • もっと綺麗なグラフを描く

markdownで出力すると全く触れなくなったので動画を作った。

www.youtube.com

回転できる3Dプロット

%matplotlib inline

でプロットを作成すると、ノートブックの中に画像として出力されてしまうので回せない。

%matplotlib notebook

を使うとノートブックの外に出力したときと同じように、インタラクティブなグラフが描ける。 こちらの3次元プロットを例として描く。

%matplotlib notebook
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import numpy as np

x = np.arange(-3, 3, 0.25)
y = np.arange(-3, 3, 0.25)
X, Y = np.meshgrid(x, y)
Z = np.sin(X)+ np.cos(Y)

fig = plt.figure()
ax = Axes3D(fig)
ax.plot_wireframe(X,Y,Z)

plt.show()

f:id:ksknw:20160804225927p:plain

スライドバーつき

ipywidgetsというやつを使えばできる。 jupyter notebookをインストールしただけで使おうとすると以下のコマンドを実行するように表示されたので、おとなしく従うと動くようになった。

sudo jupyter nbextension enable --py --sys-prefix widgetsnbextension

以下を実行すると、スライドバーが表示される。 このバーをいじると変数の値が動的に変わるので、コールバックにプロットを登録しておくと、スライドバーをいじるだけで勝手にプロットが変わる。

自分の環境では

%matplotlib notebook

とするとカクカクでいまいちだったが、

%matplotlib inline

にするといい感じ。

%matplotlib inline
from ipywidgets import interact
import numpy as np

def scatter(num_data):
    x = range(num_data)
    y = [np.sin(t/5.0) for t in x]
    plt.plot(x, y)
    plt.show()

interact(scatter, num_data=(1,200, 1), value=2)

f:id:ksknw:20160804225709p:plain

ipywidgetsには他にもボタンとかプログレスバーとか色々あって、ノートを簡単なアプリにできる。

%matplotlib inline

import pylab as plt
import time
from ipywidgets import FloatProgress
from IPython import display

prg = FloatProgress(min=0, max=99, value=1)
display.display(prg)

for i in range(100):
    prg.value = i
    time.sleep(0.01)

もっと綺麗なグラフを描く

matplotlibではなくて、plotlyというライブラリを使えばできる。 matplotlibと書き方が結構違うので、改めて覚えてまでやるほどかと言われると微妙だが、グラフは綺麗でかっこいい。 コードやグラフ、データをオンラインで共有したりもできるらしい。

デモのサイトを見ていると可能性を感じる。

from plotly.offline import iplot, init_notebook_mode
import plotly.plotly as py
from plotly.graph_objs import Scatter, Data


init_notebook_mode(connected=True) 

trace0 = Scatter(
    x=[1, 2, 3, 4],
    y=[10, 15, 13, 17]
)

trace1 = Scatter(
    x=[1, 2, 3, 4],
    y=[16, 5, 11, 9]
)

data = Data([trace0, trace1])
unique_url = iplot(data, filename = 'basic-line')

f:id:ksknw:20160804230515p:plain

参考

DCGANで湖岸道路ぽい画像を生成する

f:id:ksknw:20160625232128g:plain:w300

はじめに

先日、小倉であった某学会に参加した。 面白そうな話はいくつかあったけど、中でもDNNを使った画像の生成について興味を持った。 深層学習ウェイ系の某先生もマルチモーダルとか生成とか言ってた気がするし、 判別するより生成するほうが見た目に楽しそうなので、 こちらの発表でも使われていたDCGANを使って、画像を生成してみることにした。

Deep Convolutional Generative Adversarial Networks

論文はこちらとてもわかりやすいこちらのブログがわかりやすい。

GANはAdversarialという名の通り、2つのネットワークを競合させて学習を行うアルゴリズム。 GANでは普通のAutoEncoderなんかで画像を生成する場合と違って、Discriminatorというやつを作る。 Discriminatorは入力された画像が、Generatorの生成した画像か元画像かを判別する。 GeneratorはDiscriminatorを騙せるように画像の生成を学習し、Discriminatorはそれを判別できるように学習を行う。 これによって、最終的に単に入力画像だけからGeneratorを学習するよりも良くなるということらしい。

具体的には以下のような式を最適化する。

f:id:ksknw:20160626130749p:plain

右辺第1項は入力画像に対してDが出力を正しく出したときに大きくなり、第2項はDがGの画像に反応しないときに大きくなる。 Gはこの評価関数を小さくする、つまりDが騙される方に学習し、Dは逆に正しく見分けられるように学習する。

面白いのは、画像をGeneratorに直接学習させて、例えば入力画像との誤差を最小化するのではなく、 Discriminatorを騙せさえすればいいというところ(と思う)。 実際には入力画像にはなくても、Dが騙される「ぽい」画像を作れば評価関数の値を下げることができる。

DCGANはDiscriminatorにConvNet、GaneratorにDeconvNetを使う

mattyaさんによるchainerの実装

とてもおもしろい研究だし、arxivを見る感じ、最初にDCGANの論文がでたのが2015年の11月と少し古いのもあって、 既に色々なライブラリを使った実装が公開されていて、自分でコードを書く必要はなさそう。 chainerがお気に入りなので、その中から、mattyaさんによる実装を使わせてもらう。

git clone して動かそうとすると、残念ながら2GしかないGPUではメモリ足りなくて動かなかった。 コードを見ていると順次メモリを解放してやればギリギリいけそうだったので、210行目あたりに

del z
del x
del x2
del yl
del yl2
del L_dis
del L_gen

と書いてやると何とか動くようになった。

入力データ

入力データとして琵琶湖(マザーレイク)周辺の道路画像を使う。 動画を撮影しに行くのは大変なので、これを使ってGoogleMapのストリートビューから画像を取ってくる。

とってきた画像は以下のようなもの。 f:id:ksknw:20160612114759g:plain

道路に対して前だけでなく、横を向いている画像もある。 画像データは合計で600枚ぐらい。 この画像を237x119にリサイズしたあと、96x96の画像をランダムに切り出して4560枚の画像を作った。 入力画像からランダムに10枚選んだものを以下に示す。 思ったよりちゃんと道が写っているものが少ない。

f:id:ksknw:20160626005852p:plainf:id:ksknw:20160626005853p:plainf:id:ksknw:20160626005854p:plainf:id:ksknw:20160626005855p:plainf:id:ksknw:20160626005916p:plainf:id:ksknw:20160626005917p:plainf:id:ksknw:20160626005918p:plainf:id:ksknw:20160626005919p:plainf:id:ksknw:20160626010125p:plainf:id:ksknw:20160626010126p:plain

画像の中にはかなり類似したものもあるだろうし、そもそもDNNの学習に数千枚の画像というのは少ないように思うけど、JSAIの発表でもそれぐらいの枚数だったように思うので、まずはやってみる。

結果

GTX750Tiで学習。うるさいパソコンとともに数日過ごした。 10万ぐらいのグラボがほしい。

学習途中で生成された画像は以下。 上半分(50枚)の画像はエポック間で同じzベクトルで生成している画像。下半分は毎回ランダムにサンプリングしなおしたzベクトルで生成している。

最初 f:id:ksknw:20160626000823p:plain

1エポック後 f:id:ksknw:20160626000902p:plain

5エポック後 f:id:ksknw:20160626000936p:plain

10エポック後 f:id:ksknw:20160626001032p:plain

61エポック後 f:id:ksknw:20160625230705p:plain

62エポック後 f:id:ksknw:20160626003455p:plain

63エポック後 f:id:ksknw:20160626003511p:plain

大体10エポックぐらいから遠目にはそんなに変わってないように思う。 あと、上半分は同じzベクトルから生成しているはずなのに、61から63エポックで毎回結構大きく変わっているのも気になる。学習が収束している感じがあまりしない。画像が少ないとかそういうことが関係しているのだろうか。

なんか赤いなにかを作ったりしているし、まだもうちょっと学習してほしい感もある。あと1週間ぐらい回したほうがいいのかもしれない。

zベクトルをいじって色々画像を作る

DCGANでは学習時に乱数のベクトル(zベクトル)をDeconvNetworkに入力して画像を生成する。 生成した全ての画像において、Discriminatorを騙せるように学習するので、 理想的には学習後はどのzベクトルを選択しても、それっぽい画像を生成できるようになるはずである。 しかも、論文を読んでるとzベクトルを動かしていくと連続的に画像が変化していくらしい。楽しそう。

ということで適当に画像を生成した後、そこから3枚選んでそれらの間にある画像を生成してみる。 以下のように適当にプログラムを書く。

def interpolate(index1, index2, img_i):
    import pandas as pd

    epoch = 61

    all_z = pd.read_csv("z.csv").get_values()
    print all_z.shape
    sub = all_z[index2] - all_z[index1]
    gen = Generator().to_gpu()
    dis = Discriminator().to_gpu()

    serializers.load_hdf5("%s/dcgan_model_dis_%d.h5" % (out_model_dir, epoch), dis)
    serializers.load_hdf5("%s/dcgan_model_gen_%d.h5" % (out_model_dir, epoch), gen)
    gen.to_gpu()
    dis.to_gpu()
    z = xp.random.uniform(-1, 1, (batchsize, nz), dtype=np.float32)
    for i in range(100):
        z[i, :] = xp.array(all_z[index1] + sub / 100.0 * i)
    pylab.rcParams['figure.figsize'] = (0.96, 0.96)

    print z
    z = Variable(z)
    x = gen(z, test=True)
    x = x.data.get()

    for i_ in range(100):
        tmp = ((np.vectorize(clip_img)(x[i_, :, :, :]) + 1) / 2).transpose(1, 2, 0)
        pylab.imshow(tmp)
        pylab.axis('off')
        pylab.savefig('interpolate_imgs/%09d.png' % (img_i))
        img_i += 1

        pylab.clf()
    f.close()

interpolate(4, 10, 0)  # index+1がtest_imgsの画像の番号と対応する
interpolate(10, 23, 100)
interpolate(23, 4, 200)

深い理由はあまりないけど、なんとなく色々な画像を作っていて良さそうなので、61エポックのときの学習モデルを使うことにした。 以下のような3枚を選んで、それらの間のzベクトルを線形につないで、それぞれの間ごとに100枚ずつ画像を生成する。 f:id:ksknw:20160625235309p:plain:w300

生成された画像をつなげて動かしてみると以下のようになった。 f:id:ksknw:20160625232128g:plain:w300

なんか、ぽいものが生成されている。

最初の何枚かを静止画で見る。

f:id:ksknw:20160626012429p:plainf:id:ksknw:20160626012430p:plainf:id:ksknw:20160626012431p:plainf:id:ksknw:20160626012432p:plainf:id:ksknw:20160626012433p:plainf:id:ksknw:20160626012434p:plainf:id:ksknw:20160626012435p:plainf:id:ksknw:20160626012436p:plainf:id:ksknw:20160626012437p:plainf:id:ksknw:20160626012438p:plainf:id:ksknw:20160626012439p:plainf:id:ksknw:20160626012440p:plainf:id:ksknw:20160626012441p:plainf:id:ksknw:20160626012442p:plainf:id:ksknw:20160626012443p:plainf:id:ksknw:20160626012444p:plainf:id:ksknw:20160626012445p:plainf:id:ksknw:20160626012446p:plainf:id:ksknw:20160626012447p:plainf:id:ksknw:20160626012448p:plainf:id:ksknw:20160626012449p:plainf:id:ksknw:20160626012450p:plainf:id:ksknw:20160626012451p:plainf:id:ksknw:20160626012452p:plain f:id:ksknw:20160626012453p:plainf:id:ksknw:20160626012454p:plainf:id:ksknw:20160626012455p:plainf:id:ksknw:20160626012456p:plainf:id:ksknw:20160626012457p:plainf:id:ksknw:20160626012458p:plainf:id:ksknw:20160626012459p:plainf:id:ksknw:20160626012500p:plainf:id:ksknw:20160626012501p:plainf:id:ksknw:20160626012502p:plainf:id:ksknw:20160626012503p:plainf:id:ksknw:20160626012504p:plainf:id:ksknw:20160626012505p:plainf:id:ksknw:20160626012506p:plainf:id:ksknw:20160626012507p:plainf:id:ksknw:20160626012508p:plainf:id:ksknw:20160626012509p:plainf:id:ksknw:20160626012510p:plainf:id:ksknw:20160626012511p:plainf:id:ksknw:20160626012512p:plainf:id:ksknw:20160626012513p:plainf:id:ksknw:20160626012514p:plainf:id:ksknw:20160626012515p:plainf:id:ksknw:20160626012516p:plainf:id:ksknw:20160626012517p:plainf:id:ksknw:20160626012518p:plainf:id:ksknw:20160626012519p:plainf:id:ksknw:20160626012520p:plainf:id:ksknw:20160626012521p:plainf:id:ksknw:20160626012522p:plainf:id:ksknw:20160626012523p:plainf:id:ksknw:20160626012524p:plainf:id:ksknw:20160626012525p:plainf:id:ksknw:20160626012526p:plainf:id:ksknw:20160626012527p:plainf:id:ksknw:20160626012528p:plainf:id:ksknw:20160626012529p:plainf:id:ksknw:20160626012530p:plainf:id:ksknw:20160626012531p:plainf:id:ksknw:20160626012532p:plainf:id:ksknw:20160626012533p:plainf:id:ksknw:20160626012534p:plainf:id:ksknw:20160626012535p:plainf:id:ksknw:20160626012536p:plainf:id:ksknw:20160626012537p:plainf:id:ksknw:20160626012538p:plain

木が建物っぽい何かになりつつ、道の向きも気づいたら変わっている(気がする)。

まとめ

DCGANを使って、琵琶湖周辺の道路っぽい画像を生成した。 zを連続的に変化させると画像も連続的に変化していった。 ぱっと見はそれなりにできているように思うけれど、1枚1枚ちゃんと見ると、全体的にまだ線がぼやっとしていたり、自然な画像はまだ少し遠い。

誰か早くLSTM-GAN発表してコードをくれ。どうせ半年後までには誰かやるんだろ、わかってんだぞ

参考

キーボード上のいらないキーをよく使うキーに変える

はじめに

日本語配列のキーボードには変換、無変換とかカナとかほぼ使わないキーがいくつかある。 これらのキーは親指のホームポジションあたりにあるので、Altなどのよく使うキーに置き換えて有効活用することにした。

入れ替えたキーは以下

  • 無変換->全角半角
  • 変換->Alt
  • かな->Alt
  • Insert -> Delete

Xmodmap

XmodmapはXorgでキーマップを変更するためのやつ。 こちらが詳しい。

~/.Xmodmap に以下を書き込む

! 無変換->全角半角
keycode 102 = Zenkaku_Hankaku
! 変換->Alt
keycode 100 = Alt_R
! かな->alt
keycode 101 = Alt_R


! 以下はHHKBではいらない
! Insert -> Delete
! keycode 118 = Delete

!はコメント。HHKBを使っているのでInsertを入れ替える必要はなかったけど、会社では普通のキーボードを使っているので書いておく。

あとは~/.zshrcなどに

xmodmap ~/.Xmodmap

と書いておけば、ログインするたびに上の設定を読み込んでくれる。

C-nでDownしたい、が、できない。

Emacsユーザーなので、Chromeなんかを使ってるとC-nで大量のウィンドウを生成したり、C-pで印刷プレビューを開いてしまうことがよくある。

Xmodmapを使えばできるかと思ったけど、できなかった。 試したのは以下。

keycode 57 = n N Down Down

この右側2つがてっきりCtrl押した時の挙動とかそういうのかと思ったけど、どうやら違うらしい。 こちらによると、これはMode_switchというキーを押した時の挙動らしい。

さらにこちらには

"Usually the Mode_switch key is used on non-US keyboards for a few selected keys only"

とか書いてあって、なんじゃそりゃ。

左CtrlをこのMode_switchキーに割り当てることはできるけど、それをやるともはやCtrlキーでなくなってしまうので、本末転倒であった。

とりあえず、Chromeのアドオンで対応することにする。

参考

PythonをEmacsで書く+α

はじめに

この記事を書いてから早1年。 暇なときにちょこちょこといじっているうちに、helmを導入したり、tabbarを入れたりと色々変わっていた。 全部書いていると多すぎるので、中でも一番変わったpythonを書くための設定について書く。 加えて気に入っているパッケージについても書く。

全体の設定のgithubのレポジトリはこちら。 github.com elispのコードを多少書いたけど、その部分だけでなく、common.elに書いてある部分がないと動かないものもあると思う。

環境

python mode

1年前はelpyとか使ってなんとかかんとかpythonを書いていた。 今思うとよく頑張って書いてたなって思う。 今はpythonを書くために主に以下のパッケージに頼っている。 これらは補完、コード規約準拠、文法チェック、テンプレート展開の機能。

以下を使うためには、

    $ sudo apt-get install pyflakes
    $ sudo pip install jedi epc autopep8

をする必要がある。

jedi

jediはオムニ補完、つまり文法的な部分をある程度汲んだ上で補完をしてくれるパッケージ。 補完だけでなく関数の定義にジャンプする機能もある。 公式はこちら

こちらの設定を参考にして、設定した。

  (jedi:setup)
  (define-key jedi-mode-map (kbd "<C-tab>") nil) ;;C-tabはウィンドウの移動に用いる
  (setq jedi:complete-on-dot t)
  (setq ac-sources
    (delete 'ac-source-words-in-same-mode-buffers ac-sources)) ;;jediの補完候補だけでいい
  (add-to-list 'ac-sources 'ac-source-filename)
  (add-to-list 'ac-sources 'ac-source-jedi-direct)
  (define-key python-mode-map "\C-ct" 'jedi:goto-definition)
  (define-key python-mode-map "\C-cb" 'jedi:goto-definition-pop-marker)
  (define-key python-mode-map "\C-cr" 'helm-jedi-related-names)

何も考えずに設定すると、元々のauto-completeの補完も表示されてしまうので、それらは消しておくといいと思う。

こんな感じで補完される。 f:id:ksknw:20160501185304p:plain

autopep8

コードを書くとき、演算子の両側にスペースを入れるとか、無駄な改行を入れ過ぎないとか、色々気にしてないと見た目に汚いコードになってしまう。 手で直すぐらいなら、「動くからまあいいか」となるかもしれないけど、自動で直してくれるならそれに越したことはない。 autopep8はpep8というコーディング規約と比較してダメなところを勝手に修正してくれる機能。

例えばこんな感じのコードを書いて保存すると

# -*- coding: utf-8 -*-
def test(a):



    return a+1




if __name__ == '__main__':



    print test(   10)

勝手にこのようなコードに変換されて保存される。

# -*- coding: utf-8 -*-
def test(a):

    return a + 1


if __name__ == '__main__':

    print test(10)

ちなみに、","の後ろにスペースを置くことは規約違反ではないので、 これは修正されるけど、

    test(10, 20, 30)
    test(1  ,2  ,3)

こういうのは修正されない。

    test(10, 20, 30)
    test(1,  2,  3)

ので、縦方向にインデントを揃えたいときはスペースをどこに入れるかをちょっとだけ気にしないといけない。

設定は以下。pep8の規定に1行79字以内というのがあるが、さすがに厳しすぎないかと思うので、200字に変えてある。

(require 'py-autopep8)
(setq py-autopep8-options '("--max-line-length=200"))
(setq flycheck-flake8-maximum-line-length 200)
(py-autopep8-enable-on-save)

pyflakes (flymake-cursor)

pythonコンパイルが不要なので、逆にエラーがあったとしても、その箇所が実行されるまでわからない。 pyflakesはとても便利で、エラーやwarningを表示してくれる。 これをEmacsから呼ぶことで、画面内にエラーなどを表示できる。

    (flymake-mode t)
    ;;errorやwarningを表示する
    (require 'flymake-python-pyflakes)
    (flymake-python-pyflakes-load)

f:id:ksknw:20160501185309p:plain

こんな感じでpyflakesを勝手に実行して結果を表示してくれる。

yasnippet

他のエディタになかなか移れない原因の一端がこの機能。 短いフレーズを入力した後、tabを押すと予め登録されたフレーズを挿入してくれる。 githubのレポジトリはこちら。 単純なフレーズの挿入だけでなく、

class $0:
    def __init__(self, $1)
        $2

とか書いておくと、$のところを順番にカーソル移動させて穴埋めみたいな感じでコードを書くこともできる。

以下のようなものがたくさん登録してある。

  • . + Tab

    self.

  • np + Tab

    import numpy as np

  • ifmain + Tab

    if __name__ == '__main__':

気に入っているパッケージ

helm

1年前はanythingを使っていた。 anything-filelist+には大変お世話になったんだけど、helmにはhelm-locateというものがあって(anythingにもあるのかもしれないが)、そっちのほうが便利そうだったので、helmに乗り換えた。 helm-locateはデフォルトではand検索ができなくて不便だったけれど、こちらの設定を使うと、and検索できるようになって、all.filelistとかも作らなくていいしfilelist+いらないなってなった。

全般的な設定はこちらを参考にした 。

locateはファイルを作った直後なんかは更新されていないので、 $ sudo updatedb とかやる必要がたまにある。

見た目

本質的ではないけど、なんやかんやこれをいじるのが一番楽しいかもしれない。

最近追加したもののなかでTabbarがいい感じだった。 ファイル開くときはだいたいhelm使ってるから別にいらないっちゃいらないけど。

;;tabバーを追加する。
; tabbar.el http://cloverrose.hateblo.jp/entry/2013/04/15/183839
(require 'tabbar)
(tabbar-mode 1)
;; グループ化しない
(setq tabbar-buffer-groups-function nil)
;;画像はいらない
(setq tabbar-use-images nil)
;; 左に表示されるボタンを無効化
(dolist (btn '(tabbar-buffer-home-button
               tabbar-scroll-left-button
               tabbar-scroll-right-button))
  (set btn (cons (cons "" nil)
                 (cons "" nil))))
;; タブ同士の間隔
;; http://ser1zw.hatenablog.com/entry/2012/12/31/022359

(setq tabbar-separator '(0.8))

(defun my-tabbar-buffer-list ()
  (delq nil
        (mapcar #'(lambda (b)
                    (cond
                     ;; Always include the current buffer.
                     ((eq (current-buffer) b) b)
                     ((buffer-file-name b) b)
                     ((char-equal ?\  (aref (buffer-name b) 0)) nil)
;;          ((equal "*scratch*" (buffer-name b)) b) ; *scratch*バッファは表示する
             ((equal "*eww*" (buffer-name b)) b) ; *eww*バッファは表示する
             ((char-equal ?* (aref (buffer-name b) 0)) nil) ; それ以外の * で始まるバッファは表示しない
                     ((buffer-live-p b) b)))
                (buffer-list))))
(setq tabbar-buffer-list-function 'my-tabbar-buffer-list)

;; tabbar外観変更
(set-face-attribute
 'tabbar-default nil
 :family (face-attribute 'default :family)
 :background (face-attribute 'mode-line-inactive :background)
 :height 1.0)
(set-face-attribute
 'tabbar-unselected nil
 :background (face-attribute 'mode-line-inactive :background)
 :foreground (face-attribute 'mode-line-inactive :foreground)
 :box nil)
(set-face-attribute
 'tabbar-selected nil
 :background (face-attribute 'mode-line :background)
 :foreground (face-attribute 'mode-line :foreground)
 :box nil)

その他細々と設定して、こんな感じの見た目になる。 f:id:ksknw:20160501185314p:plain

Emacsを再起動したときに、終了時の状態に戻す

Emacsを起動して、画面を分割して、昨日開いていたファイルを開いてってやるのめんどくさい。 こちらを参考にして、desktop-save-modeを設定した。 以下は設定後に自動的にinit.elに追加されていた。

;;再起動時に色々復元
(custom-set-variables
 ;; custom-set-variables was added by Custom.
 ;; If you edit it by hand, you could mess it up, so be careful.
 ;; Your init file should contain only one such instance.
 ;; If there is more than one, they won't work right.
 '(desktop-save-mode t))
(custom-set-faces
 ;; custom-set-faces was added by Custom.
 ;; If you edit it by hand, you could mess it up, so be careful.
 ;; Your init file should contain only one such instance.
 ;; If there is more than one, they won't work right.
 )

今後やろうと思ってる設定

  • quick run + pdb

    Emacs上でpythonを起動してEmacs上でデバッグしたい。 パスの設定とか面倒でまだやってない。

  • ein mode

    ipython(jupyter) notebookをEmacsから使うやつ。jupyterのバージョンが上がったせいか上手く使えない。

参考