盆暗の学習記録

データサイエンスを中心として,日々学んだことの備忘録としていく予定です。初心者であり独学なので内容には誤りが含まれる可能性が大いにあります。

cronでプログラムを定期的に自動実行させる

かなり詰まったのでメモ。

環境

cron

cronの設定方法は2つある

  • crontabコマンドを使う方法
  • /env/cron.dにcrontabファイルを置いておく方法

crontabコマンドを使う方法

概要

この方法のメリット

  • 簡単で、管理者でなくても設定できる

この方法のデメリット

  • ユーザーごとに単一のcrontabファイルを使うため、もし誤ってcrontab -rで削除してしまうと全部消える(人によってはこちらは非推奨とするらしい)

コマンド

# crontabを閲覧
crontab -l
# crontabを編集
crontab -e

設定方法

# m h  dom mon dow   command
分 時 日 月 曜日  コマンド

例えば次のような感じに書く。

# 全部* → 1分ごと
* *    * * *  [command]

# 毎時0分に実行
0 *    * * *  [command]

# 毎日0:00に実行
0 0    * * *  [command]

また、5つの*の部分の代わりに@monthlyなどの書き方をしても良い (Ubuntu Manpage: crontab - cron を駆動するための一覧表

/etc/cron.dに置く方法

概要

この方法のメリット

  • タスクのジャンルごとにcrontabファイルを分けられるので、本格的な使用に向く
  • あやまって削除してしまうなんてことは起こりにくい

この方法のデメリット

  • 管理者権限で操作する必要がある

crontabファイルの準備

/etc/crontabのフォーマットを使うのが楽なので、/etc/cron.dというディレクトリにコピーして編集する。

# test_cronという名前でコピーする場合
cp /etc/crontab /etc/cron.d/test_cron

この際、ファイル名に拡張子をつけない(.が含まれるファイルは実行されないらしい)

設定方法

crontabコマンドを使う場合と少し異なり、userのフィールドがある。

# m h dom mon dow user  command
分 時 日 月 曜日 ユーザー コマンド

例えば次のように書く。

# m h dom mon dow user  command
* *     * * *   user   echo "test" >> ~/test.log

私が詰まったポイント

私の場合、userをrootにしていたときはcronが動作しなかった。

自分がログイン中のユーザー名にしたところ、ちゃんと動作するようになった。

cronの再読み込み

cronデーモンは/etc/cron.d内のファイルや/etc/crontabの変更を監視しており、変更はすぐ再読み込みしてくれるはずだが、手動でやる場合は以下のようなコマンドが使える。

# cronの再起動
sudo service cron restart
# cronの状態を確認する
sudo service cron status

crontabの設定を楽にするサイト

f:id:nigimitama:20191110182341p:plain

f:id:nigimitama:20191110195920p:plain

参考になるサイト

cron力をつけよう!全てのcrontab入門者に贈る9個のテクニック · DQNEO起業日記

↑実践的な話でとてもよかった

スパースなデータをXGBoost.DMatrixに入れるときはpd.DataFrame/np.arrayを使ってはいけない

XGBoostにはDMatrixという独自のデータ保持用クラスがあります。Documentの説明では

optimized for both memory efficiency and training speed

と書いてあり、私は「自動で疎行列クラスとかにしてくれるんだろうなぁ」と思っていたのですが、そうでもない様子なのでメモしておきます。

以下では、「core.pyを読む」の節でxgboostのpythonで書かれた部分について調べた結果を述べ、それ以降の節ではC++で書かれた部分についてコードを動かすことで中身を推測していきます(私はC++が読めないため)。

(ここでは色々事情があってスパースなデータを使う場合のことを考えますが、決定木系のアルゴリズムならカテゴリカル変数はOne Hot Encoding(≒ダミー変数)ではなく、まずLabel Encodingを検討すべきです)

core.pyを読む

DMatrixクラスはpython-package/xgboost/core.pyの321行目あたりで定義されています。

f:id:nigimitama:20191106190007p:plain

DMatrixの__init__メソッドではデータのオブジェクトをdataで受け取るのですが、わりとすぐに_maybe_pandas_dataというメソッドに通しています。

# xgboost/core.py line 378
data, feature_names, feature_types = _maybe_pandas_data(data,
                                                        feature_names,
                                                        feature_types)

_maybe_pandas_data226行目あたりで定義されていて、pandasがインストール済みで、かつdataがDataFrameである場合にpd.DataFrameをnp.arrayに変換しています。

# xgboost/core.py line 226
def _maybe_pandas_data(data, feature_names, feature_types):
    """ Extract internal data from pd.DataFrame for DMatrix data """
    if not (PANDAS_INSTALLED and isinstance(data, DataFrame)):
        return data, feature_names, feature_types

    # (中略)
    
    data = data.values.astype('float')
    return data, feature_names, feature_types

DMatrix.__init__ではその下を見ていくと、dataのクラスに応じて別のinitメソッドを呼ぶ条件分岐があります。

scipy.sparse系の疎行列クラスであれば_init_from_csr_init_from_cscが呼ばれ、pd.DataFrameやnp.arrayであれば_init_from_npy2dが呼ばれることになります。

# xgboost/core.py line 393
if isinstance(data, (STRING_TYPES, os_PathLike)):
    handle = ctypes.c_void_p()
    _check_call(_LIB.XGDMatrixCreateFromFile(c_str(os_fspath(data)),
                                             ctypes.c_int(silent),
                                             ctypes.byref(handle)))
    self.handle = handle
elif isinstance(data, scipy.sparse.csr_matrix):
    self._init_from_csr(data)
elif isinstance(data, scipy.sparse.csc_matrix):
    self._init_from_csc(data)
elif isinstance(data, np.ndarray):
    self._init_from_npy2d(data, missing, nthread)
elif isinstance(data, DataTable):
    self._init_from_dt(data, nthread)
else:
    try:
        csr = scipy.sparse.csr_matrix(data)
        self._init_from_csr(csr)
    except:
        raise TypeError('can not initialize DMatrix from'
                        ' {}'.format(type(data).__name__))

小まとめ

_init_from_npy2dなどの中身はC++のメソッドを呼んでいて、そこから先はC++が読めない私ではわからない領域なのですが、ここまでのコードを読む限り、

(少なくともpython側では)疎行列クラスでDMatrixに投入しない限り、pd.DataFrameやnp.arrayが自動で疎行列クラスに変換されることはない

という事がわかります。

このことが分析上どう影響するかについては次節で調べてみます。

メモリ使用量の違い

私はc++で書かれている処理を読むことができないため、ここからはコードを動かしてみて 「pd.DataFrameからDMatrixオブジェクトを作った場合のメモリ使用量」と 「疎行列クラスからDMatrixオブジェクトを作った場合のメモリ使用量」 を比較して、間接的に両者の違いを見てみます。

実験方法

以下の方法で比較します。

  • 500,000行×1000列のスパースな特徴量を作る
  • DMatrixに投入して、元のデータオブジェクトを削除した時点のメモリ使用量で比較する

コード

import xgboost as xgb
import pandas as pd
import numpy as np
import psutil
import gc
import sys
from scipy.sparse import csr_matrix

# setting
nrows = 500_000
ncols = 1000
use_sparse = True if sys.argv[1] == "sparse" else False

# generate data
np.random.seed(0)
categories = pd.Series(np.random.randint(low=0, high=ncols, size=nrows))
X = pd.get_dummies(categories)

if use_sparse:
    X = csr_matrix(X)
else:
    # どうせ_init_from_npy2dでnp.float32のarrayに変換されるので、モニタリングしやすくするため予め変換しておく
    X = np.array(X, dtype=np.float32)

print(f"data size: {sys.getsizeof(X) / 1024**2:.2f} MB, shape: {X.shape}")
print(f"used memory before make DMatrix: {psutil.virtual_memory().used / 1024**3:.2f} GB")
dtrain = xgb.DMatrix(X)
del X
gc.collect()
print(f"used memory after delete X: {psutil.virtual_memory().used / 1024**3:.2f} GB")

結果

$ python3 DMatrix.py DataFrame
data size: 1907.35 MB, shape: (500000, 1000)
used memory before make DMatrix: 7.04 GB
used memory after delete X: 8.89 GB

$ python3 DMatrix.py sparse
data size: 0.00 MB, shape: (500000, 1000)
used memory before make DMatrix: 5.16 GB
used memory after delete X: 5.17 GB

考察

pd.DataFrame/np.arrayを使ったほうは、元のデータが約1.9GBあり、DMatrix変換後にメモリの使用量が約1.9GB増えています。DMatrixに投入したデータがそのままメモリに乗っかったような数値で推移しています。

つまり、python側からC++で書かれたメソッドにデータが渡った後も、特にnp.arrayを疎行列クラスに変換するような処理は行われていないと考えられます。

一方で疎行列クラスを使用した方は、DMatrixにデータを投じる前後でメモリ使用量がほとんど変化していません。

疎行列クラスのデータをDMatrixに投入すると、その先でも疎行列クラスとして受け取ってくれるのだと考えられます。

すなわち、スパースデータを使う場合は疎行列クラスに変換した後にDMatrixに投入すべきであると考えられます。

計算時間の違い

再びコードを動かしてみて 「pd.DataFrameから作ったDMatrixオブジェクトでの学習時間」と 「疎行列クラスから作ったDMatrixオブジェクトでの学習時間」 を比較して、間接的に両者の違いを見てみます。

実験方法

以下の条件で行います。

  • 回帰問題
  • 100,000行×100列の二値変数のみの特徴量
  • 予測を100回繰り返して、平均的な計算時間を測る
    • その際の時間を比較する
  • 疎行列クラスはscipy.sparse.csr_matrixを使用

なお、計測環境は以下のとおりです

コード

import xgboost as xgb
import pandas as pd
import numpy as np
import sys
from scipy.sparse import csr_matrix
import timeit

# setting
nrows = 100_000
ncols = 100
use_sparse = True if sys.argv[1] == "sparse" else False

# generate data
np.random.seed(0)
categories = pd.Series(np.random.randint(low=0, high=ncols, size=nrows))
X = pd.get_dummies(categories)
w = np.random.uniform(size=ncols)
y = np.dot(X, w) + np.random.normal(size=nrows)

if use_sparse:
    X = csr_matrix(X)
print(f"type of X: {type(X)}")

print(f"data size: {sys.getsizeof(X) / 1024**2:.2f} MB, shape: {X.shape}")
dtrain = xgb.DMatrix(X, label=y)

param = {"tree_method": "exact"}
t = timeit.Timer("xgb.train(param, dtrain)", globals=globals())
results = t.repeat(repeat=100, number=1)
print(f"train time: {np.mean(results):.1f} ± {np.std(results):.3f} sec. (n={len(results)})")

結果

$ python3 performance.py DataFrame
type of X: <class 'pandas.core.frame.DataFrame'>
data size: 9.54 MB, shape: (100000, 100)
train time: 1.5 ± 0.237 sec. (n=100)

$ python3 performance.py sparse
type of X: <class 'scipy.sparse.csr.csr_matrix'>
data size: 0.00 MB, shape: (100000, 100)
train time: 0.1 ± 0.015 sec. (n=100)
  • DataFrameからDMatrixを作って学習した場合は平均1.5秒程度
  • 疎行列クラスからDMatrixを作って学習した場合は平均0.1秒程度

の学習時間になりました。

考察

疎行列クラスはスパースデータの計算が高速化されることが知られています。

DataFrameを使ってDMatrixを作った場合は、xgboostの内部で疎行列クラスに変換されるようなことが無いため、疎行列クラスからDMatrixを作った場合に比べて圧倒的に計算が遅いという結果になったのだと思われます。

計算時間の点からも、スパースデータを使う場合は疎行列クラスに変換した後にDMatrixに投入すべきであると考えられます。

まとめ

  • pd.DataFrameやnp.arrayのデータをxgboost.DMatrixに入れると、DMatrix.__init__()ですべての要素がnp.float32のnp.arrayに変換される。

    • これの何が問題か:pd.get_dummies()の直後はuint8だがfloat32になるとメモリ使用量が約4倍増加することになる →DataFrame時点ではわからないがひっそりと大量のメモリを食う
    • np.arrayに変換されたデータはその後C++で書かれたAPIにデータが渡った後も疎行列クラスなどに変換されている様子はなく、float32のarrayのままであると考えられる。
      • (スパースなデータの場合、無駄にメモリ使用量と計算時間が多くなる)
  • scipy.sparseの疎行列クラスのデータをDMatrixに入れた場合、妙な変換はされず、疎行列クラスの特性(省メモリ・高速な計算)を持ったままXGBoostの学習に入ることができる

→スパースなデータを使う場合は、疎行列クラスのインスタンスにした後にDMatrixに入れるべき

XGBoostで自作の目的関数を使う

XGBoostの素敵なポイントの一つは、自分で定義した関数を目的関数に使うことができる点です。

でもどういう関数にしたらよいのかがわからなくて過去に戸惑ったことがあるのでメモしておきます。(詳しいやり方はXGBoostのdocumentationに書いてあります)

定義すべき関数

XGBoostなどの勾配ブースティング決定木アルゴリズムでは誤差関数・目的関数を2種類使います。

  1. early stopping1用の誤差関数(metric)
  2. 個々の木の学習で使う目的関数(objective)

誤差関数(metric function)

metricは誤差二乗和(Sum of Squared Error: SSE)

や平均二乗誤差(Mean of Squared Error: MSE)のことです。他の機械学習アルゴリズムでもよく見かけるやつです。

XGBoostでは、予測値のarrayと実測値が入ったDMatrixを受け取り、文字列と誤差を評価した実数を返すような関数を定義します。誤差二乗和なら以下のような感じです。

import numpy as np
import xgboost as xgb

def my_sse(y_hat: np.array, dtrain: xgb.DMatrix) -> [str, float]:
    """custom metric: sum of squared error"""
    y = dtrain.get_label()
    N = y_hat.shape[0]
    error = np.sum((y_hat - y)**2)
    return "sse", float(error)

目的関数(objective function)

XGBoostでは、objective functionは予測値のarrayと実測値が入ったDMatrixを受け取ってmetricの1次の勾配と2次の勾配を返すような関数として実装します。

metricがSSEなら、目的関数の出発点は総和を取る前の二乗誤差(squared error)であり、微分後の式を簡単にするために1/2を掛けて

とすれば、1次の勾配と2次の勾配は

となります。

def my_squared_error(y_hat: np.array, dtrain: xgb.DMatrix) -> [np.array, np.array]:
    """custom objective: squared error"""
    y = dtrain.get_label()
    N = y_hat.shape[0]
    gradient = y_hat - y
    hessian = np.ones(N)
    return gradient, hessian

使用例

上記の自作squared errorが正しく動いているか検証してみます。

from sklearn.datasets import load_boston
from sklearn.model_selection import train_test_split
import pandas as pd
import xgboost as xgb
import numpy as np


def my_sse(y_hat: np.array, dtrain: xgb.DMatrix) -> [str, float]:
    """custom metric: sum of squared error"""
    y = dtrain.get_label()
    N = y_hat.shape[0]
    error = np.sum((y_hat - y)**2)
    return "sse", float(error)

    
def my_squared_error(y_hat: np.array, dtrain: xgb.DMatrix) -> [np.array, np.array]:
    """custom objective: squared error"""
    y = dtrain.get_label()
    N = y_hat.shape[0]
    gradient = y_hat - y
    hessian = np.ones(N)
    return gradient, hessian


# load data
boston = load_boston()
X = boston["data"]
y = boston["target"]

# split data
X_train, X_test, y_train, y_test = train_test_split(X, y,
                                                    test_size=0.2, random_state=0)
# for early stopping
X_train, X_eval, y_train, y_eval = train_test_split(X_train, y_train,
                                                    test_size=0.2, random_state=0)

# make DMatrix data
dtrain = xgb.DMatrix(X_train, label=y_train)
deval = xgb.DMatrix(X_eval, label=y_eval)
dtest = xgb.DMatrix(X_test, label=y_test)

# 標準のmetricではmseは存在せずrmseしかないのでrmseにする
params = {'objective': 'reg:squarederror', 'metric': 'rmse'}
watchlist = [(deval, 'eval')]
n_trees = 100
# 標準のsquared errorで学習
model1 = xgb.train(params, dtrain, num_boost_round=n_trees,
                   evals=watchlist, early_stopping_rounds=2)
# カスタム目的関数で学習
model2 = xgb.train(params, dtrain, num_boost_round=n_trees,
                   evals=watchlist, early_stopping_rounds=2,
                   obj=my_squared_error, feval=my_sse)

# 予測結果が一致するかどうか確認
pred1 = model1.predict(dtest)
print(f"head of pred1: {pred1[0:5]}")
pred2 = model2.predict(dtest)
print(f"head of pred2: {pred2[0:5]}")
print(f"予測結果が一致するもの:{sum(pred1 == pred2)}件({sum(pred1 == pred2) / X_test.shape[0]:.0%})")
head of pred1: [23.791985 26.905416 23.237139 10.941677 21.959444]
head of pred2: [23.791985 26.905416 23.237139 10.941677 21.959444]
予測結果が一致するもの:102件(100%)

結果は一致しました。問題なくsquared error / mse(rmse)を実装できていると思われます。

自作の目的関数の使い所

例えば「(予測精度が多少下がってもいいから)予測値が実測値を上回らないようにしてほしい」とか、あるいはその逆とかいった非対称な評価関数を設定したいときには使えると思います。

また、Kaggleにおいてコンペの評価指標に向けて最適化するためにカスタム目的関数を利用する場合もあるみたいですね。


  1. early stoppingとは、決定木を作るごとに予測誤差を評価して、決定木を追加しても予測誤差が改善しなくなってきたときに決定木の追加を止める手法のことです。

データ分析時のメモリ使用量を減らす方法

最近少し覚えたことをまとめます。

(基本的にpythonのコードと共に述べていきますが、Rの場合についても少し触れていきます。)

不要なオブジェクトの削除

「使わなくなったオブジェクトを削除してメモリを開放する」という考えです。

del

例えば以下のコードを試してみます。

import gc
import psutil
import numpy as np
import sys
nrows = 500000
ncols = 100
print('データを生成します')
X = np.random.rand(nrows, ncols)
print(f'データのサイズ: {sys.getsizeof(X) / 1024**2: .1f} MB')
mem_before = psutil.virtual_memory().used
print(f'del前のメモリ使用量: {mem_before / 1024**3 :.1f} GB')
del X
mem_after = psutil.virtual_memory().used
mem_diff = mem_before - mem_after
print(f'del後のメモリ使用量: {mem_after / 1024**3 :.1f} GB ({mem_diff / 1024**2 :.1f} MB reduced)')

上のコードを実行すると、以下のような結果になります。

データを生成します
データのサイズ:  381.5 MB
del前のメモリ使用量: 11.5 GB
del後のメモリ使用量: 11.2 GB (362.4 MB reduced)

手動でのガベージコレクションは不要

delのあとに明示的にガベージコレクションgc.collect())を使う方法も考えられますが、基本的に自動でのガベージコレクションが働くのでやらなくて問題ないみたいです。

また、ガベージコレクションを行ったとしても、メモリーがOSに返されるとは限らないようです1

Rの場合

一方でRはオブジェクトを削除しただけではメモリ使用量が変わらないため、手動でのガベージコレクションが有効なようです2

rm(X)
gc(reset = TRUE)

データ型の最適化

numpyにはuint8float32int64などの細かいデータ型が用意されています。

例えば以下のようにnp.iinfo()np.finfo()のメソッドを使うと、そのデータ型でもてる最小値・最大値などの情報を見ることができます。

import numpy as np
print(np.iinfo(np.int8))
print(np.iinfo(np.int16))
print(np.iinfo(np.int32))
print(np.iinfo(np.int64))
print(np.finfo(np.float16))
print(np.finfo(np.float32))
print(np.finfo(np.float64))
Machine parameters for int8
---------------------------------------------------------------
min = -128
max = 127
---------------------------------------------------------------

Machine parameters for int16
---------------------------------------------------------------
min = -32768
max = 32767
---------------------------------------------------------------
(以下略)

8bitのint8は-128~127の範囲の値しか表現できませんが、例えば二値変数(0, 1)や月(1~12)や日(1~31)のデータなどであればこのくらいの表現力で十分ですよね。なので最適化しましょう、というのがこのアプローチの考え方です。

自動で型変換するコード

こんな感じのやつを使っています。

import pandas as pd
import numpy as np

def reduce_mem_usage(df: pd.DataFrame) -> pd.DataFrame:
    """データフレームのすべてのカラムのデータ型を最適化してメモリ使用量を抑えるやつ"""
    start_mem = df.memory_usage().sum() / 1024**2
    print(f'Memory usage of dataframe is {start_mem:.2f} MB')
    df = df.apply(_optimize_dtype, axis=0)
    end_mem = df.memory_usage().sum() / 1024**2
    print(f'Memory usage after optimization is: {end_mem:.2f} MB')
    print(f'Decreased by {100 * (start_mem - end_mem) / start_mem:.1f}%')
    return df


def _optimize_dtype(col):
    """あるカラム/変数のデータ型を最適化するやつ"""
    col_type = col.dtype
    if col_type != 'object' and col_type != 'datetime64[ns]':
        c_min = col.min()
        c_max = col.max()
        # intにできるなら変換するが、NAがあればfloatのままにする(intはNAを保持できない)
        is_int = False
        if _is_int(col):
            if col.isnull().sum() == 0:
                is_int = True
            else:
                print(f'{col.name} was not convert to int because it has NA')
        # bit数の最適化
        if is_int:
            if c_min > np.iinfo(np.int8).min and c_max < np.iinfo(np.int8).max:
                col = col.astype(np.int8)
            elif c_min > np.iinfo(np.int16).min and c_max < np.iinfo(np.int16).max:
                col = col.astype(np.int16)
            elif c_min > np.iinfo(np.int32).min and c_max < np.iinfo(np.int32).max:
                col = col.astype(np.int32)
            elif c_min > np.iinfo(np.int64).min and c_max < np.iinfo(np.int64).max:
                col = col.astype(np.int64)
        else:
            if c_min > np.finfo(np.float16).min and c_max < np.finfo(np.float16).max:
                # feather-formatがfloat16を受け入れないためfloat32にする
                col = col.astype(np.float32)
            elif c_min > np.finfo(np.float32).min and c_max < np.finfo(np.float32).max:
                col = col.astype(np.float32)
            else:
                col = col.astype(np.float64)
    return col


def _is_int(col):
    """あるカラム/変数が整数かどうかを判定する"""
    fractional, _ = np.modf(col)
    if np.sum(fractional) == 0: # 小数部がすべて0
        return True
    else:
        return False

使用例は次のような感じです。

import numpy as np
import pandas as pd
import sys

# データを生成
nrows = 500_000
ncols = 100
X = [np.random.randint(low=0, high=50, size=nrows) for x in range(ncols)]
X = pd.DataFrame(np.array(X).T)
print('データを生成 ---------------------------')
print('dtypes: ')
print(X.dtypes.value_counts().to_frame(name='rows'))
print(f'データ型変換前のオブジェクトのサイズ: {sys.getsizeof(X) / 1024**2 : .1f} MB')

print('\nデータ型を変換 -------------------------')
X = reduce_mem_usage(X)
print('dtypes: ')
print(X.dtypes.value_counts().to_frame(name='rows'))
print(f'データ型変換後のオブジェクトのサイズ: {sys.getsizeof(X) / 1024**2 : .1f} MB')
データを生成 ---------------------------
dtypes:
       rows
int64   100
データ型変換前のオブジェクトのサイズ:  381.5 MB

データ型を変換 -------------------------
Memory usage of dataframe is 381.47 MB
Memory usage after optimization is: 47.68 MB
Decreased by 87.5%
dtypes:
      rows
int8   100
データ型変換後のオブジェクトのサイズ:  47.7 MB

このコードはkaggleのkernelを元にえじさんが改造したものを自分用に改造したものです。

えじさんのものからの変更点は

  • loggerに関する部分をprintに置き換え(私がまだloggerに慣れていないため)
  • floatであっても小数点以下が全部0でNAの無いカラムはintに変換する
  • forループでなくapply()を使う(applyの方が何倍か速い3

です。

疎行列クラスを使う

ダミー変数(One Hot Encodingした列)が多いデータの場合、非ゼロ要素の情報しか保持しない疎行列用のクラスを使うことでメモリ使用量を大幅に削減できます。

pythonの場合はscipy.sparseを、Rの場合は{Matrix}パッケージを使うことで効率よく疎行列を保持できます。

import pandas as pd
import numpy as np
import sys
import scipy.sparse
nrows = 500000
ncols = 500
categories = pd.Series(np.random.randint(low=0, high=ncols, size=nrows))
X = pd.get_dummies(categories)
print(f'dtypes: {X.dtypes.value_counts().to_dict()}')
print(f'データのサイズ: {sys.getsizeof(X) / 1024**2: .1f} MB')
X_csr = scipy.sparse.csr_matrix(X, shape=X.shape)
print(f'csr_matrixのサイズ: {sys.getsizeof(X_csr): .1f} Byte')
dtypes: {dtype('uint8'): 500}
データのサイズ:  238.4 MB
csr_matrixのサイズ:  56.0 Byte

  1. アプリケーションのメモリー使用量の最適化

  2. Rでメモリを解放したい - まずは蝋の翼から。

  3. 使用例で使った50万行x100列のデータの変換を10回行って実行時間の平均をとると、applyを使う場合1.62±0.01秒、forループの場合19.668±0.292秒でした。

pythonで塗り分け地図を描く

f:id:nigimitama:20190928225127p:plain

基本的なplotの仕方を覚えたのでメモ。

geopandasによるplot

データの用意

市区町村レベルの行政区域データは国土数値情報の行政区域データから、 町丁レベルの行政区域データはe-Statの境界データダウンロードの小地域データから取得できます。

e-Statの小地域データを例にやってみます

input_path = "shapefile" # 展開先フォルダ名

# e-stat 国勢調査 小地域(町丁・字等別) 東京都全域
shapefile_url = "https://www.e-stat.go.jp/gis/statmap-search/data?dlserveyId=A002005212015&code=13&coordSys=1&format=shape&downloadType=5"
shapefile_name = "tokyo.zip"

# ダウンロード
from urllib.request import urlretrieve
urlretrieve(url=shapefile_url, filename=shapefile_name)

# 解凍
import zipfile
with zipfile.ZipFile(shapefile_name) as existing_zip:
    existing_zip.extractall(input_path)

# ファイル名を取得
import os
files = os.listdir(input_path)
shapefile = [file for file in files if ".shp" in file][0]
print(f"downloaded shapefile: {shapefile}")

# 読み込み
import geopandas as gpd
shapefile_path = os.path.join(input_path, shapefile)
df = gpd.read_file(shapefile_path, encoding='cp932')
print(f"{shapefile_path} is loaded")

# 東京都の島嶼部を除く
import pandas as pd
islands = ['大島町', '利島村', '新島村', '神津島村', '三宅村', '御蔵島村', '八丈町', '青ヶ島村', '小笠原村']
is_not_islands = [df["CITY_NAME"] != island for island in islands]
is_not_islands = pd.concat(is_not_islands, axis=1).all(axis=1)
df = df.loc[is_not_islands, :]

# 陸地だけにする
df = df.loc[df["HCODE"] == 8101, :]

こんな感じのデータです(詳細はe-Statから定義書をご参照ください)

f:id:nigimitama:20190928224557p:plain

plot

全体図

島嶼部を除いた、いつもの東京都の図。

# 全体図
df.plot(figsize=[10,10])

f:id:nigimitama:20190928224655p:plain

塗り分け図

# 人口で塗り分け
df.plot(column="SETAI", legend=True, figsize=[30,10], cmap="Oranges")

f:id:nigimitama:20190928224751p:plain

# 人口密度で塗り分け
df["pop_density"] = df["JINKO"] / df["AREA"]
df.plot(column="pop_density", legend=True, figsize=[30,10], cmap="Blues")

f:id:nigimitama:20190928224833p:plain

一部の地域ですごく高い値がある関係で、ほとんどの地域で薄い色に…

値のラベルをつけた図

# 値ラベル用にgeometryから当該ポリゴン内のある地点を取得
df["coords"] = df["geometry"].apply(lambda x: x.representative_point().coords[:])
df["coords"] = [coords[0] for coords in df["coords"]]

# 文京区の人口密度
temp = df.query("CITY_NAME == '文京区'")
temp.plot(column="pop_density", legend=True, figsize=[30,10], cmap="Blues")

# 値ラベル
import matplotlib.pyplot as plt
for i, row in temp.iterrows():
    plt.annotate(s=round(row["pop_density"], 2), xy=row["coords"], horizontalalignment="center")

f:id:nigimitama:20190928225127p:plain

[WSL Ubuntu]pythonで塗り分け地図を描くための環境構築

f:id:nigimitama:20190928230707p:plain

Ubuntuへのgeopandasというライブラリのインストールまでの流れと、塗り分け図のプロットの例を書いていきます。

私のPCの環境について

WindowsですがWSLを使います。

geopandasのインストール

windowsへのインストールはめちゃくちゃ面倒くさい上に何度やっても失敗したのでUbuntuに入れます。

sudo apt update
sudo apt upgrade

# pipがない場合、インストール
sudo apt install python3-pip
# Geopandasの依存ライブラリをインストール
sudo pip3 install numpy pandas shapely fiona pyproj six descartes
# Install Geopandas
sudo pip3 install git+https://github.com/geopandas/geopandas.git

jupyterが入っていない場合は

sudo pip3 install jupyterlab

でjupyterを入れ、$ jupyter labで起動します。

起動して

import matplotlib.pyplot as plt
%matplotlib inline

import geopandas as gpd
world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))
world.plot()

と書いて

f:id:nigimitama:20190928222903p:plain

のように図が出ればインストールは成功です。

記事冒頭の図は

world.plot(column="pop_est", figsize=[10,10], cmap="Oranges")

というコードで描くことができます。色が濃い国は人口が多い国です。

f:id:nigimitama:20190928230724p:plain

Chromeの標準フォントを変更してギザギザ文字を無くす

f:id:nigimitama:20190825120233p:plain

はじめに

WindowsMacに比べるとデフォルトのフォントが汚いことが多いですが,CSSが適切に設定されていないWebサイトのフォントについてはChrome側のデフォルトのフォント設定を変更することで解決できるみたいです。

やり方をメモしておきます。

設定方法

Chromeの「設定」→「フォントをカスタマイズ」から設定できます。

f:id:nigimitama:20190825124607p:plain

お好みのフォント・フォントサイズを設定してください。

f:id:nigimitama:20190825124654p:plain

おすすめフォント

游ゴシックにすると文字が細すぎて見づらいので,「源ノ角ゴシック(source han sans)」というオープンソースのフォントにしたらいい感じでした

f:id:nigimitama:20190825122852p:plain

源ノ角ゴシックCode (Source Han Code JP)というプログラミング用のもあるので固定幅フォントはRictyの代わりにこちらにしてもいいかも。

游ゴシックの場合、Medium(やや太字)に設定するのがおすすめです。