雑記

Anaconda 環境を利用した Windows への Keras (TensorFlow backend) 導入

以下の記事の手順でできました。
qiita.com今回は Keras を利用するので Keras も導入しました。Anaconda インストール後、コマンドプロンプトを立ち上げて以下のように「仮想環境作成 → 仮想環境立ち上げ → TensorFlow導入 → Keras導入」しました。仮想環境を作成するとき Python のバージョンを指定できますが、3.5 にしないといけないようです。

C:\work> conda create -n tensorflow python=3.5 anaconda
C:\work> activate tensorflow
(tensorflow) C:\work> pip install --ignore-installed --upgrade https://storage.googleapis.com/tensorflow/windows/cpu/tensorflow-1.0.1-cp35-cp35m-win_amd64.whl
(tensorflow) C:\work> pip install keras

これで Keras (TensorFlow backend) がつかえるようになったことを確認するために以下のスクリプトを実行してみます。中間層が2層の全結合ニューラルネットワークで iris を分類します。

# -*- coding: utf-8 -*-
import numpy
import random
from keras.models import Sequential
from keras.layers.core import Dense, Activation
from keras.utils import np_utils
from sklearn import datasets

if __name__ == '__main__':
  # irisデータのロード
  iris = datasets.load_iris()
  x = iris.data
  y = np_utils.to_categorical(iris.target, num_classes=3)

  # シャッフルして訓練データとテストデータに分割
  n = iris.target.shape[0]
  indices = random.sample(range(n), n)
  x = x[indices]
  y = y[indices]
  n_train = int(n / 2)
  x_train = x[range(n_train)]
  y_train = y[range(n_train)]
  x_test = x[range(n_train, n)]
  y_test = y[range(n_train, n)]

  # モデルの構築と学習
  model = Sequential()
  model.add(Dense(10, input_shape=(x.shape[1],), activation='relu'))
  model.add(Dense(10, activation='relu'))
  model.add(Dense(3, activation='softmax'))
  model.compile(optimizer='adam', loss='categorical_crossentropy', metrics=['accuracy'])
  model.fit(x_train, y_train, epochs=100, batch_size=10, shuffle=True)

  # テスト結果表示
  scores = model.evaluate(x_train, y_train, verbose=0)
  print('Accuracy (train): %.2f%%' % (scores[1]*100))
  scores = model.evaluate(x_test,  y_test,  verbose=0)
  print('Accuracy (test) : %.2f%%' % (scores[1]*100))

以下のように学習できました。

(tensorflow) C:\work> python test.py
...
Epoch 99/100
75/75 [==============================] - 0s - loss: 0.3137 - acc: 0.9733
Epoch 100/100
75/75 [==============================] - 0s - loss: 0.3072 - acc: 0.9733
Accuracy (train): 97.33%
Accuracy (test) : 96.00%

さらに Anaconda 環境への OpenCV の導入

さらに画像処理ライブラリ OpenCV も導入します。
Anaconda で対象の仮想環境を立ち上げた状態で、以下のコマンドで導入できます。

(tensorflow) C:\work> conda install -c menpo opencv3

導入が上手くいったかどうか、以下のスクリプトで確認します。

# -*- coding: utf-8 -*-
import cv2
import numpy

if __name__ == '__main__':
  image = numpy.zeros((256, 256, 3), numpy.uint8)
  for i in range(256):
    for j in range(256):
      image[i, j, 0] = i
      image[i, j, 1] = j
  cv2.imwrite('hoge.png', image)

以下のような画像ファイル hoge.png が生成されます。

OpenCV で動画作成

OpenCV で動画を作成することができます。以下のスクリプトで movie.avi が生成されます。以下では各フレームを画素の行列から生成していますが、画像ファイルを読み込んでつなげて動画にすることもできます。

# -*- coding: utf-8 -*-
import cv2
import numpy

if __name__ == '__main__':
  vvw = cv2.VideoWriter('movie.avi', cv2.VideoWriter_fourcc('X','V','I','D'), 10, (200, 200))
  for i in range(10):
    for j in range(10):
      j_ = (9 - j) if (i % 2 == 1) else j
      image = numpy.zeros((10, 10, 3), numpy.uint8)
      image[i, j_, 0] = i  * 25
      image[i, j_, 1] = j_ * 25
      image = cv2.resize(image, (200, 200), interpolation=cv2.INTER_AREA)
      vvw.write(image)

以下のような動画ができます。.avi ははてなブログに貼り付けられないようなので、以下のページで .gif に変換しました。
AVI GIF 変換。オンライン フリー — Convertio

なお、cv2.VideoWriter の第2引数の4文字はコーデック指定です（ Video Codecs by FOURCC - fourcc.org ）。最適かどうかはわかりませんがとりあえず ('X','V','I','D') で生成できました。第3引数は1秒あたりのフレーム数、第4引数は画像サイズです。

雑記

データを分析するとき、取り込んだままのデータの列をそのままつかうのではなくて、元々あった列をもとに新しい列をつくることとかあると思います。でも Python とかよくわからないので Google 検索してどうやって加工するか調べると思います。

例えば、元々あった列のうち a, b, c 列の中央値をとり、新しい列 d に追加したいと思います。
以下の (1) (2) はどちらも望みを満たすようであることがわかります。でも (1) は遅いのでやってはいけないとわかります。
というか (2) のように一気にできるとわかれば (1) のように行ごとの処理にはしないと思いますが調べ方が難しいと思います。

# -*- coding: utf-8 -*-
import pandas
import numpy

if __name__ == "__main__":
  data = pandas.DataFrame({'name': numpy.repeat('Alice', 100000),
                           'a': range(0, 100000, 1),
                           'b': range(0, 200000, 2),
                           'c': range(0, 400000, 4)})

  # (1) 44.674 sec
  data['d'] = data.apply(lambda x: numpy.median(x[['a', 'b', 'c']]), axis=1) 
  # (2) 0.010 sec
  data['d'] = numpy.median(data[['a', 'b', 'c']].values, axis=1)

ともあれ、numpy.mean や numpy.median などで各行への処理を一気にできるのだったらそうすればいいことはわかります。
ただ、もっと自分で定義した処理をしたいときはベクトル演算が定義されているとは限らないので困ると思います。
そのような場合は numpy.vectorize で関数をベクトル化することにより処理速度を改善できると以下の記事にあります。
sinhrks.hatenablog.com

例えば、データに元々ある time 列をもとに、米国サマータイム期間かどうかを is_dst 列として付加したいとします。
apply をつかうと (1) のようにかけます。でも上の記事のアドバイス通りこの処理をベクトル化した (2) の方が速いです。

# -*- coding: utf-8 -*-
import pandas
import numpy

def is_dst(dt):
  return(pandas.Timestamp(dt).tz_localize('Asia/Tokyo').tz_convert('US/Eastern').dst().seconds > 0)

if __name__ == "__main__":
  data = pandas.DataFrame({
    'name': numpy.repeat('Alice', 100000),
    'time': numpy.repeat(['2017/03/12 15:59', '2017/03/12 16:00'], [50000, 50000])})
  data.time = pandas.to_datetime(data.time)

  # (1) 28.834 sec
  data['is_dst'] = data.apply(lambda x:
    x.time.tz_localize('Asia/Tokyo').tz_convert('US/Eastern').dst().seconds > 0, axis=1)
  # (2) 9.768 sec
  data['is_dst'] = numpy.vectorize(is_dst)(data.time)

ところで、上の実装で自作の is_dst 関数は受け取った日付時刻を pandas.Timestamp 型に変換していますが、time 列は元々 pandas.Timestamp 型です。なのになぜこの変換をしているのかというと、ベクトル化した関数を time 列に適用したときに time 列が勝手に numpy.datetime64 型になってしまったからです。なんでそうなるのかは Python とかよくわからないのでわかりません。

雑記： sklearn.cluster.KMeans で k-means 法

sklearn.cluster.KMeans をつかってみたいときとかあると思います。iris データにつかうと以下のようになります。

# -*- coding: utf-8 -*-
import pandas
from sklearn import datasets
from sklearn.cluster import KMeans
from sklearn.metrics import confusion_matrix

if __name__ == "__main__":
  iris = datasets.load_iris()
  predict = KMeans(n_clusters=3).fit(iris.data)

  mat = pandas.DataFrame(confusion_matrix(iris.target, predict.labels_))
  mat.index = iris.target_names
  print(mat)

             0   1   2
setosa       0  50   0
versicolor  48   0   2
virginica   14   0  36

setosa はすべてクラスタ1に分類され（かつ、クラスタ1に他の品種は混ざっていない）、一方、virginica はクラスタ0とクラスタ2に分裂してしまったことがわかります。各品種が分類された最大クラスタが対角線上に並んでいないと見にくいという人は、以下のように列を入れ替えてください。

  mat = mat.loc[:,mat.idxmax(axis=1)]
  mat.columns = range(iris.target_names.shape[0])
  print(mat)

             0   1   2
setosa      50   0   0
versicolor   0  48   2
virginica    0  14  36

クラスタリングした結果、各クラスタの重心はどこになったかも確認できます（以下）。

  print(predict.cluster_centers_)

[[ 6.85        3.07368421  5.74210526  2.07105263]
 [ 5.006       3.418       1.464       0.244     ]
 [ 5.9016129   2.7483871   4.39354839  1.43387097]]

初期クラスタ重心を与えることもできます。
以下のようにすれば元データの 0 番目、50 番目、100 番目を初期クラスタ重心にできます。

# -*- coding: utf-8 -*-
import numpy
import pandas
from sklearn import datasets
from sklearn.cluster import KMeans
from sklearn.metrics import confusion_matrix

if __name__ == "__main__":
  iris = datasets.load_iris()
  initial_centers = numpy.array([iris.data[0], iris.data[50], iris.data[100]])
  predict = KMeans(n_clusters=3, init=initial_centers).fit(iris.data)

  mat = pandas.DataFrame(confusion_matrix(iris.target, predict.labels_))
  mat.index = iris.target_names
  print(mat)
  print(predict.n_iter_)

ただし、これを実行すると以下のように警告文が出てきます。これの意味は、KMeans はクラスタ初期化方法が k-means++（デフォルト）や random の場合はデフォルトで n_init=10 通りの初期値からのクラスタリングを試み、最もよいクラスタリング結果を出力します（モーメントがよく更新され続けた回をもっともよいクラスタリング結果とする）。しかし、初期クラスタ重心が明示的に与えられたときにはその他の初期化を試みるのはおかしいので、デフォルトの n_init=10 は無視して1回しかクラスタリングを試みませんよという意味です。その代わりにというか、クラスタリングのイテレーション回数 n_iter_ が返却されます。クラスタ重心は3回更新されるのみだったようです。

/usr/local/lib/python2.7/site-packages/sklearn/cluster/k_means_.py:889: RuntimeWarning:
 Explicit initial center position passed: performing only one init in k-means
 instead of n_init=10
  return_n_iter=True)
             0   1   2
setosa      50   0   0
versicolor   0  48   2
virginica    0  14  36
3

雑記： Keras で自作レイヤー（Dense）

Keras でオリジナルの自作レイヤーを追加したいときとかあると思います。
自作レイヤー自体は以下の記事でつかったことがありますが、これはウェイトをもつレイヤーではなく、最後にかぶせて損失関数のみをカスタマイズするためのレイヤーでした。
Keras で変分自己符号化器（VAE）を学習したい - クッキーの日記
ウェイトをもつレイヤーはどう追加するのか知りたいのですが、以下のドキュメントの通りにすればできるのか確認します。
Writing your own Keras layers - Keras Documentation
上のページの MyLayer は単純な Dense のようなので、実際に Dense に置き換えられるか確認します。

keras.layers.Dense をつかった場合

まず、 普通に Keras の Dense をつかって MNIST の分類器をつくると例えば以下になると思います。モデルは適当です。
→ 下記スクリプトの設定で 95% 程度のテスト精度になりました。

# -*- coding: utf-8 -*-
import numpy
from keras.layers import Input, Dense, Dropout
from keras.models import Model
from keras.datasets import mnist
from keras.utils import np_utils

if __name__ == "__main__":
  (x_train, y_train), (x_test, y_test) = mnist.load_data()
  x_train = x_train.astype('float32') / 255.
  x_test = x_test.astype('float32') / 255.
  x_train = x_train.reshape((len(x_train), numpy.prod(x_train.shape[1:])))
  x_test = x_test.reshape((len(x_test), numpy.prod(x_test.shape[1:])))
  y_train = np_utils.to_categorical(y_train, num_classes=10)
  y_test = np_utils.to_categorical(y_test, num_classes=10)

  input = Input(shape=(x_train.shape[1],))
  converted = Dense(20, activation='relu')(input)
  converted = Dropout(0.2)(converted)
  converted = Dense(20, activation='relu')(converted)
  converted = Dropout(0.2)(converted)
  converted = Dense(10, activation='softmax')(converted)
  classifier = Model(input, converted)
  classifier.compile(optimizer='adam', loss='categorical_crossentropy',
                     metrics=['accuracy'])
  print(classifier.summary())

  classifier.fit(x_train, y_train, epochs=20, batch_size=100, shuffle=True)

  scores = classifier.evaluate(x_train, y_train, verbose=0)
  print('Accuracy (train): %.2f%%' % (scores[1]*100))
  scores = classifier.evaluate(x_test, y_test, verbose=0)
  print('Accuracy (test) : %.2f%%' % (scores[1]*100))

_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
input_1 (InputLayer)         (None, 784)               0         
_________________________________________________________________
dense_1 (Dense)              (None, 20)                15700     
_________________________________________________________________
dropout_1 (Dropout)          (None, 20)                0         
_________________________________________________________________
dense_2 (Dense)              (None, 20)                420       
_________________________________________________________________
dropout_2 (Dropout)          (None, 20)                0         
_________________________________________________________________
dense_3 (Dense)              (None, 10)                210       
=================================================================
Total params: 16,330
Trainable params: 16,330
Non-trainable params: 0
_________________________________________________________________

Accuracy (train): 95.69%
Accuracy (test) : 94.92%

自作 Dense をつかった場合

次に、import Dense をつかわずに自作 Dense ( class MyDense ) に置き換えたのが以下です。
→ これも 95% 程度のテスト精度になったので、ちゃんと Dense になっていると思います。
　 パラメータ数も上と同じ 16,330 です。

Keras ドキュメント（Writing your own Keras layers - Keras Documentation）にあるとおり、自作レイヤーのクラスでは build メソッドでウェイトを定義し、call メソッドにウェイトを用いたこのレイヤーの計算処理を実装すればよいはずです。compute_output_shape メソッドは次のレイヤーのために出力 shape を返すようにしておきます。

# -*- coding: utf-8 -*-
import numpy
from keras import backend as K
from keras.engine.topology import Layer
from keras.layers import Input, Dropout
from keras.models import Model
from keras.datasets import mnist
from keras.utils import np_utils

class MyDense(Layer):
  def __init__(self, output_dim, activation, **kwargs):
    self.output_dim = output_dim
    self.activation = activation
    super(MyDense, self).__init__(**kwargs)

  def build(self, input_shape):
    self.kernel = self.add_weight(name='kernel',
                                  shape=(input_shape[1], self.output_dim),
                                  initializer='glorot_uniform')
    self.bias   = self.add_weight(name='bias',
                                  shape=(1, self.output_dim),
                                  initializer='zeros')
    super(MyDense, self).build(input_shape)

  def call(self, x):
    if self.activation == 'relu':
      return(K.relu(K.dot(x, self.kernel) + self.bias))
    elif self.activation == 'softmax':
      return(K.softmax(K.dot(x, self.kernel) + self.bias))
    else:
      return(K.dot(x, self.kernel) + self.bias)

  def compute_output_shape(self, input_shape):
    return(input_shape[0], self.output_dim)

if __name__ == "__main__":
  (x_train, y_train), (x_test, y_test) = mnist.load_data()
  x_train = x_train.astype('float32') / 255.
  x_test = x_test.astype('float32') / 255.
  x_train = x_train.reshape((len(x_train), numpy.prod(x_train.shape[1:])))
  x_test = x_test.reshape((len(x_test), numpy.prod(x_test.shape[1:])))
  y_train = np_utils.to_categorical(y_train, num_classes=10)
  y_test = np_utils.to_categorical(y_test, num_classes=10)

  input = Input(shape=(x_train.shape[1],))
  converted = MyDense(20, activation='relu')(input)
  converted = Dropout(0.2)(converted)
  converted = MyDense(20, activation='relu')(converted)
  converted = Dropout(0.2)(converted)
  converted = MyDense(10, activation='softmax')(converted)
  classifier = Model(input, converted)
  classifier.compile(optimizer='adam', loss='categorical_crossentropy',
                     metrics=['accuracy'])
  print(classifier.summary())

  classifier.fit(x_train, y_train, epochs=20, batch_size=100, shuffle=True)

  scores = classifier.evaluate(x_train, y_train, verbose=0)
  print('Accuracy (train): %.2f%%' % (scores[1]*100))
  scores = classifier.evaluate(x_test, y_test, verbose=0)
  print('Accuracy (test) : %.2f%%' % (scores[1]*100))

_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
input_1 (InputLayer)         (None, 784)               0         
_________________________________________________________________
my_dense_1 (MyDense)         (None, 20)                15700     
_________________________________________________________________
dropout_1 (Dropout)          (None, 20)                0         
_________________________________________________________________
my_dense_2 (MyDense)         (None, 20)                420       
_________________________________________________________________
dropout_2 (Dropout)          (None, 20)                0         
_________________________________________________________________
my_dense_3 (MyDense)         (None, 10)                210       
=================================================================
Total params: 16,330
Trainable params: 16,330
Non-trainable params: 0
_________________________________________________________________

Accuracy (train): 95.83%
Accuracy (test) : 95.14%

class MyDense の実装は Keras ドキュメント（Writing your own Keras layers - Keras Documentation）の MyLayer を参考にしましたが、必要に応じて keras.engine.topology.Layer のソースを参照して修正しました。
keras/topology.py at master · fchollet/keras · GitHub

• 当初 Keras ドキュメントの MyLayer のとおりに実装したところ、パラメータ数からバイアス項がないことに気付いたので、バイアス項を追加しました。
• add_weight の引数に name が必要だったので追加しました。name をつけて何につかうのかわかりません。→ 後から自作 MyDense にはバイアス項がないことに気付いたとき、ユニークな名前でいくつでも weight 行列を追加できることに気付きました。
• kernel と bias の初期化方法は本家 Dense に合わせました（初期化が uniform だといまいち到達精度が 94% 程度にしかなりませんでした）（誤差もあるかもしれませんが）。
• 今回のモデルでは Dense を ReLU または softmax で活性化するので、Dense よろしく MyDense インスタンス生成時に活性化関数も指定するようにしました（ReLU と softmax 以外未対応）（活性化関数自体を引数に取ればいいんだろうとは思います）。

---

これで Dense が自作できることは確認できたのですが、自分がつくりたいのは Dense ではないのでこれからもっと調べていきます。

雑記： 対象日付が米国サマータイム期間か判定（R, Python）

時刻カラムがあるようなデータのプレ処理で、各レコードが米国サマータイム期間かどうか自動判定したいときとかあると思います。サマータイムかそうでないかで東京とニューヨークの時差が違うからね。
とりわけ、時刻カラムのタイムゾーンは日本時間なのですが、その時刻がニューヨークでサマータイムかどうかがほしいことがあると思います。以下の IS_DST 列を追加したいということです（サマータイム突入時には、現地では深夜2時を3時にジャンプさせますが、このジャンプ時刻は日本時間だと16時になります）。

DATETIME	IS_DST
2017/03/12 15:59:59	0
2017/03/12 16:00:00	1

時刻文字列を入れると、日本時間のその時刻にニューヨークではサマータイムかどうか返す関数があればよく、R だと例えば以下のようにかけます（参考： R の日付時刻オブジェクト - クッキーの日記）。

is_dst <- function(datetime_string='2017/03/12 15:59:59') {
  time_lt <- strptime(datetime_string, format='%Y/%m/%d %H:%M:%S')
  time_ct <- as.POSIXct(time_lt, tz='Japan')
  return(as.POSIXlt(time_ct, tz="EST5EDT")$isdst == 1)
}

print(is_dst('2017/03/12 15:59:59')) # FALSE
print(is_dst('2017/03/12 16:00:00')) # TRUE

ミソは、isdst のようなプロパティを取り出せるのは POSIXlt 型ですが、POSIXlt 型のままタイムゾーン間を行き来することはできないので、一旦 POSIXct 型を経由するところです（もっとスマートな方法があるかもしれないですがわからない）。

一方、Python では例えば以下のようにかけます（だいたい pandas をつかうので pandas.Timestamp 型をつかいます）。

# -*- coding: utf-8 -*-
import pandas

def is_dst(datetime_string='2017/03/12 15:59:59'):
  datetime = pandas.Timestamp(datetime_string)
  datetime = datetime.tz_localize('Asia/Tokyo')
  datetime = datetime.tz_convert('US/Eastern')
  return(datetime.dst().seconds > 0)

if __name__ == "__main__":
  print(is_dst('2017/03/12 15:59:59')) # False
  print(is_dst('2017/03/12 16:00:00')) # True

時刻文字列をタイムスタンプに変換して、タイムゾーンを日本時間にして、アメリカ東部標準時に変換するだけです。dst() というメソッドは Standard Time との時差（つまり、サマータイムのときは1時間で、そうでなければ0時間）が取得できるようなので、0 より大きければサマータイムです。

---

上の Python スクリプトでいま家では上手くいったんだけど、会社ではある1つの日付文字列から上の方法でサマータイムかどうかを上手く取得できなかった気がするんだよね（何かエラーが出た）。Pandas DataFrameで日時データのタイムゾーン変換 - pLog という記事によると、データフレーム型ではタイムスタンプをインデックスにしないと tz_convert ができないとあります。となると、ある1つの文字列に対して判定を行いたい場合、以下のようにすべきなのでしょうか…？

def is_dst(datetime_string='2017/03/12 15:59:59'):
  df = pandas.DataFrame({'datetime': [datetime_string]})
  df.index = pandas.DatetimeIndex(df.datetime, name='datetime')
  df.index = df.index.tz_localize('Asia/Tokyo')
  df.index = df.index.tz_convert('US/Eastern')
  return(df.index[0].dst().seconds > 0)

というかある1つの文字列に対して判定を行いたいだけなら pandas をつかわなくてもいいですね。
Python: How do you convert datetime/timestamp from one timezone to another timezone? - Stack Overflow

# -*- coding: utf-8 -*-
import datetime
import pytz

def is_dst(datetime_string='2017/03/12 15:59:59'):
  dt0 = datetime.datetime.strptime(datetime_string, '%Y/%m/%d %H:%M:%S')
  tz0 = pytz.timezone('Asia/Tokyo')
  tz1 = pytz.timezone('US/Eastern')
  dt1 = tz0.localize(dt0).astimezone(tz1)
  return(dt1.dst().seconds > 0)