分析ノート - Part 51

2019-05-142019-05-14

pythonでひらがなとカタカナのリストを作成する

pythonでひらがなやカタカナの一覧を作成するスマートな方法を考えてみたので紹介です。
きっかけは、ある所でcsvファイルにあ、い、う、え、お、…と一行一文字書き並べたファイルを読み込んでいるコードを見かけたことです。

この記事でリストを作成する方法を書きますが、その前に本当にひらがな、カタカナのリストが必要なのかということは考えた方が良いと思います。
例えばある文字列がひらがなのみで成り立っているかどうかとか、カタカナが含まれているかと言った判定を行うのであれば、正規表現を使った方がよいです。
リストを準備して、for文で回して一文字ずつチェックするような効率が悪い処理を実装するべきではありません。
正直な所、このリストが必要になる場面は結構珍しいと思います。

それでもひらがなカタカナのリストが必要だという場合は、
ひらがな、カタカナにそれぞれ連続した文字コードを割り当てられていることを利用して生成するのが効率がいいのではないかと思います。
こちらも参照:pythonで文字と文字コードの相互変換

具体的には次のようなコードで作成できます。
生成結果は結合してprintしました。


hiragana = [chr(i) for i in range(ord("ぁ"), ord("ゖ")+1)]
print("".join(hiragana))
# ぁあぃいぅうぇえぉおかがきぎくぐけげこごさざしじすずせぜそぞただちぢっつづてでとどなにぬねのはばぱひびぴふぶぷへべぺほぼぽまみむめもゃやゅゆょよらりるれろゎわゐゑをんゔゕゖ
katakana = [chr(i) for i in range(ord("ァ"), ord("ヺ")+1)]
print("".join(katakana))
# ァアィイゥウェエォオカガキギクグケゲコゴサザシジスズセゼソゾタダチヂッツヅテデトドナニヌネノハバパヒビピフブプヘベペホボポマミムメモャヤュユョヨラリルレロヮワヰヱヲンヴヵヶヷヸヹヺ

“あ”より”ぁ”の方が先であること、最後に +1 しておかないと末尾の文字が抜けることなどに注意です。
“ゐ”や”ゑ”、”ヸ”、”ヹ”、”ヺ”などはいらないって場合もあると思うので必要に応じて微調整してください。
もしかしたら”゛”(濁点)や”゜”(半濁点)、その他”、”や”。”などの句読点なども必要だという場面もあるかもしれません。

ord(“ぁ”)などとせずに、数字で指定しても結果は同じになりますが、
意味が明確になるので、ord(“ぁ”)の方が良いと思います。

参考


hiragana = [chr(i) for i in range(12353, 12439)]
print("".join(hiragana))
# ぁあぃいぅうぇえぉおかがきぎくぐけげこごさざしじすずせぜそぞただちぢっつづてでとどなにぬねのはばぱひびぴふぶぷへべぺほぼぽまみむめもゃやゅゆょよらりるれろゎわゐゑをんゔゕゖ
katakana = [chr(i) for i in range(12449, 12539)]
print("".join(katakana))
# ァアィイゥウェエォオカガキギクグケゲコゴサザシジスズセゼソゾタダチヂッツヅテデトドナニヌネノハバパヒビピフブプヘベペホボポマミムメモャヤュユョヨラリルレロヮワヰヱヲンヴヵヶヷヸヹヺ

2019-05-132019-05-12

pandasで要素のユニークカウント

pandasのDataframeやSeriesの要素を重複を排除してカウントとする関数の紹介です。
nuniqueというのを使います。

pandas.DataFrame.nunique
pandas.Series.nunique

分析の中で、ユニークカウントする機会はよくあると思うのですが、この関数の知名度が低いようで、
少し回りくどい方法を取っている人をよく見かけます。

使い方は簡単です。試しに、4種類の値からなる要素10個の配列でやってみましょう。


import pandas as pd
data = pd.Series(['c', 'b', 'c', 'a', 'a', 'a', 'c', 'a', 'c', 'd'])
print(data.nunique())  # 4

簡単でした。

参考ですが、この関数を使わない方法も色々あります。


# value_counts　を使って、その長さを取る。
print(len(data.value_counts()))  # 4

# ユニークな値のリストを取得してその長さを取る。
print(len(data.unique()))  # 4

# 集合に変換してその要素数を取る
print(len(set(data)))  # 4

nunique の知名度を考えると、これらの書き方の方が適切な場面もあるかもしれません。

for文で実装してるような人も見たことがありますが流石にこういうのは避けた方がいいと思います。


# 悪い例
cnt = 0
tmp_set = set()
for d in data:
    if d not in tmp_set:
        cnt += 1
        tmp_set.add(d)
print(cnt)  # 4

2019-05-122021-09-21

標準化レーベンシュタイン距離は距離関数なのか

以前の記事で、標準化レーベンシュタイン距離(標準化編集距離)というのを紹介し、
自分も使っていたのですが挙動に少し違和感があったので確認しました。

参考:標準化レーベンシュタイン距離

レーベンシュタイン距離はその名の通り、距離関数なのですが、
これを標準化してしまうとどうも距離関数っぽくない動きをしてるように思えたのです。

念の為、距離関数というもの自体の定義をおさらいしておきましょう。

集合$X$に対して、$d:X\times X \rightarrow \mathbb{R}$ が距離関数であるとは、
$x,y,z \in X$に対して次の条件が成り立つ時に言います。
1. $d(x,y) \geq 0$ (非負性)
2. $d(x, y) = 0 \Leftrightarrow x = y$ （同一律）
3. $d(x, y) = d(y, x)$ (対象律)
4. $d(x, z) \leq d(x, y) + d(y, z)$ (三角不等式)

この条件のうち、 1. 2. 3. は特に問題ないのですが、
標準化レーベンシュタイン距離については、4. の三角不等式がちょっと怪しかったです。
で、反例を探してみたところ簡単に見つかりました。
$ld(*,*)$をレーベンシュタイン距離、$nld(*,*)$を標準化レーベンシュタイン距離とし、
x = ‘ab’, y= ‘aba’, z = ‘ba’ とおきます、
すると、
$ld(x, z) = 2$ なので、$nld(x, z) = 1$ ですが、
$ld(x, y) = ld(y, z) = 1$ なので、$nld(x, y) = nld(y, z) = \frac13$ です。

そのため、 $nld(x, z) > nld(x, y) + nld(y, z) = \frac23$ となり、
三角不等式を満たしません。

標準化レーベンシュタイン距離は標準化レーベンシュタインという名前の距離関数と考えるのは誤りで、
レーベンシュタイン距離という距離関数を標準化したもの(その結果距離関数ではなくなってしまったもの)と、
考える必要があります。

現状これで激しく困ったということはないのですが、
一部のライブラリにある、自分で作った距離関数を引数に渡せるようなものには、
標準化レーベンシュタイン距離は突っ込まない方が安全そうです。

2019-05-112021-09-21

pandasで指数平滑移動平均

昨日の記事が移動平均だったので、今日は指数平滑移動平均を扱います。
初めて知った日は衝撃だったのですが、pandasには指数平滑移動平均を計算する専用の関数が用意されています。
(pythonを使い始める前はExcel VBAでいちいち実装していたので非常にありがたいです。)

馴染みがない人もいると思いますので軽く紹介しておきます。
元のデータを${x_t}$とし、期間$n$に対して指数平滑移動平均${EWMA_t}$は次のように算出されます。
$$
\begin{align}\alpha &= \frac{2}{1+n}\\
EWMA_0 &= x_0\\
EWMA_t &= (1-\alpha)*EWMA_{t-1} + \alpha * x_t
\end{align}
$$

3番目の式を自分自身に逐次的に代入するとわかるのですが、
$EWMA_t$は、$x_t$から次のように算出されます。
$$
EWMA_t = \alpha\sum_{k=0}^{\infty}(1-\alpha)^k x_{t-k}
$$
$(1-\alpha)$の絶対値は1より小さいので、この無限級数の後ろの方の項は無視できるほど小さくなります。
結果的に、過程${x_t}$の最近の値に重みを置いた加重平均と見做せます。

さて、早速ですが計算してみましょう。

pandasのDataFrameおよび、Seriesに定義されているewm関数を使います。
pandas.DataFrame.ewm


import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
# データ作成
data = pd.Series(np.random.normal(0, 100, 200).cumsum() + 20000)
# 指数平滑移動平均の計算
data_ewm = data.ewm(span=10).mean()
# 可視化
plt.rcParams["font.size"] = 14
fig = plt.figure(figsize=(12, 7))
ax = fig.add_subplot(1, 1, 1)
ax.plot(data, label="元データ")
ax.plot(data_ewm, label="指数平滑移動平均")
plt.legend()
plt.show()

出力がこちら。

ここで一つ注意する点があります。
data_ewm = data.ewm(span=10).mean()
という風に、spanという変数名で期間$10$を渡しています。
ドキュメントを読んでいただくとわかるのですが、span=をつけないと、
comという別の変数に値が渡され、$\alpha$の計算が、
$\alpha=1/(1+com)$となり、結果が変わります。

また、spanやcomを使う以外にも、alpha=で$\alpha$のあたいを直接指定することも可能です。

2019-05-102019-05-08

pandasで移動平均や高値線、安値線を計算する

前回がローソク足だったので今回も市場データでよく使われるテクニックから移動平均を取り上げてみたいと思います。
ついでにHLバンド(ドンチャンチャンネル/高値線,安値線)も同様にもとまるので紹介します。

技術としては、window関数と呼ばれる種類の関数を使って算出します。

ドキュメントはこの辺り。
Window
pandas.DataFrame.rolling
pandas.Series.rolling

DataFrameとSeries両方に実装されていて、同じように使うことができます。
rolling() で指定期間ごとに区切ったデータを作り、それに対して、 meanやmax,minなどの関数を適用して
平均や最大値、最小値を算出して配列として返します。

実際に見た方が早いと思うので、乱数でランダムウォークデータを生成し、
計算して可視化してみましょう。


import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
# データ作成
data = pd.DataFrame(np.random.normal(0, 100, 200).cumsum() + 20000)
# 移動平均、期間の高値/安値線の計算
ma = data.rolling(10).mean()
high_band = data.rolling(20).max()
low_band = data.rolling(20).min()
# 可視化
fig = plt.figure(figsize=(12, 7))
ax = fig.add_subplot(1, 1, 1)
ax.plot(data, label="元データ")
ax.plot(ma, label="移動平均")
ax.plot(high_band, label="高値線")
ax.plot(low_band, label="安値線")
plt.legend()
plt.show()

結果がこちら。

期間の最初の方はデータ不足により線が途切れています。
この辺りの制御は min_periods などの引数で細かく調整できるので、
データ量やその時の目的に応じて調整して使いましょう。

2019-05-092019-05-08

pythonでローソク足を描く

以前(pythonを勉強し始めた頃)は、matplotlibでローソク足をかけたはずなのですが、最近は方法が変わってしまったようなのでそのメモです。
なお、ここでサンプルに使うデータはすでにcsvフィアルか何かで保存されているものとします。

以前は matplotlib.finance というのをimport でき、これを使ってかけたのですが、version 2.0 からなくなってしまったようです。
matplotlib.finance
This module is deprecated in 2.0 and has been moved to a module called mpl_finance.

そしてさらに良くないことに、移動先の mpl_finance ですが、あまりしっかり保守されてない様子。

githubのリポジトリに下記の文言があります。
The code is provided as is and is basically un-maintained.

ただ、一応動くようなので動かしてみましょう。
anacondaには含まれていないようなので、インストールから必要です。
pip install mpl_finance
これで、
mpl-finance==0.10.0
が入りました。

さて、使い方ですがun-maintainedの宣言通り、 mpl-finance の公式ドキュメントらしきものは見当たらず、
上の、matplotlib.finance時代のドキュメントを読んで使わないといけないようです。

ローソク足を書く関数は次の4つあり、それぞれデータの渡し方が違います。
.candlestick2_ochl(ax, opens, closes, highs, lows, width=4, colorup=’k’, colordown=’r’, alpha=0.75)
.candlestick2_ohlc(ax, opens, highs, lows, closes, width=4, colorup=’k’, colordown=’r’, alpha=0.75)
.candlestick_ochl(ax, quotes, width=0.2, colorup=’k’, colordown=’r’, alpha=1.0)
.candlestick_ohlc(ax, quotes, width=0.2, colorup=’k’, colordown=’r’, alpha=1.0)

今回は手元のデータと相性が良いので .candlestick_ohlc を使います。
quotes に日付、始値、高値、安値、終値、の5列のデータがデータ件数行だけ並んだ配列を渡してあげる必要があります。
ここで面倒なのは日付の渡し方で、float型で渡す必要があります。
ドキュメントに time must be in float days format – see date2numとある通り、専用の関数があるのでそれを使います。

matplotlib.dates.date2num(d)

また、この関数は引数がdatetime型なので、元々が2019-05-07 のような文字列になっているならば、
datetime型に変換しておく必要があります。
それにはpandasの to_datetimeを使います。
pandas.to_datetime
(いつもならそれぞれ1記事使ってるようなテクニックですね。to_datetimeの方は便利なのでそのうち専用記事書くかも。)

前置きが長くなりましたが、ここまでの情報でできるので日経平均のcsvデータからローソク足を書いてみましょう。


import pandas as pd
import mpl_finance
import matplotlib.pyplot as plt
from matplotlib.dates import date2num

# データの読み込み
df = pd.read_csv("./日経平均データ.csv")
print(df.head())

'''
        date      open      high       low     close
0  2019-3-11  21062.75  21145.94  20938.00  21125.09
1  2019-3-12  21361.61  21568.48  21348.81  21503.69
2  2019-3-13  21425.77  21474.17  21198.99  21290.24
3  2019-3-14  21474.58  21522.75  21287.02  21287.02
4  2019-3-15  21376.73  21521.68  21374.85  21450.85
'''

# dateの型変換
# まずdatetime型にする
print(df["date"].dtypes)  # object
df["date"] = pd.to_datetime(df["date"])
print(df["date"].dtypes)  # datetime64[ns]
# 続いて float型へ
df["date"] = matplotlib.dates.date2num(df["date"])
print(df["date"].dtypes)  # float64

print(df.head())
'''
       date      open      high       low     close
0  737129.0  21062.75  21145.94  20938.00  21125.09
1  737130.0  21361.61  21568.48  21348.81  21503.69
2  737131.0  21425.77  21474.17  21198.99  21290.24
3  737132.0  21474.58  21522.75  21287.02  21287.02
4  737133.0  21376.73  21521.68  21374.85  21450.85
'''

# 可視化
fig = plt.figure(figsize=(13, 7))
ax = fig.add_subplot(1, 1, 1)
mpl_finance.candlestick_ohlc(ax, df.values)
plt.show()

こうして出来上がるチャートが次です。

正直これ単体では手間の割に可視化するメリットがないなーというのが正直なところです。
ただ、matplotlibの仕組みに乗っかっているので、
自分のオリジナルの指標などを追加していくことができます。

2019-05-082019-05-06

matplotlibのpcolorとpcolormesh

先日のトピックモデルの記事中で、試しにヒートマップでの可視化を試みた時、
matplotlibのpcolorって関数を使用しました。

参考:pythonでトピックモデル(LDA)

matplotlibでヒートマップを描こうと思うと少々無理やりな実装になるものも含めて、
imshowや、contourf、pcolorなど複数の方法が考えられ、結構迷いますが一番自然に書ける気がして採用しました。
しかしどうやらこのpcolor、あまり評判がよろしくないようです。

公式ドキュメントを見ても次ように記載があります。

matplotlib.axes.Axes.pcolor

Hint

pcolor() can be very slow for large arrays. In most cases you should use the similar but much faster pcolormesh instead. See there for a discussion of the differences.

要するに pcolormesh を使う方が良いようです。

ドキュメントはこちら。
matplotlib.axes.Axes.pcolormesh

Differences between pcolor() and pcolormesh()
の部分を読んでも、あまり pcolorにメリットを感じないので、
いっそのこと pcolor 自体を pcolormesh のエイリアスか何かに変えてしまっても良さそうなのですが、
戻り値の型が違うこともありますし、何か同じように使えない事情もあるようですね。

とりあえず今後はこれまでpcolorを使っていた場面ではpcolormeshを使うようにしようと思います。
ちなみに、x軸y軸が等間隔の場合はimshowの方がさらに速いそうです。

imshow
If X and Y are each equidistant, imshow can be a faster alternative.

2019-05-072019-05-05

標準化レーベンシュタイン距離

以前の記事で、レーベンシュタイン距離を計算できるライブラリを紹介しました。
参考:pythonで編集距離(レーベンシュタイン距離)を求める

どちらかというと、ライブラリよりもアルゴリズム側の説明の続きなのですが、
標準化されたレーベンシュタイン距離(normalized Levenshtein distance)というものも提案されています。
これは、二つの文字列のレーベンシュタイン距離を、文字数が多い方の文字数で割った値として定義されます。
固有名詞の名寄せなどでレーベンシュタイン距離を使う場合、
こちらを使った方がうまく行くことが多いようです(個人的な経験から。)

以前紹介した、 python-Levenshtein に
実装されているんじゃないかと期待してしばらく調べていたのですが、どうやらこれは実装されてないようです。
特にLevenshtein.ratio という関数に期待したのですがこれは全然違いました。

ということで自分で実装しましょう。


import Levenshtein


def normalized_distance(text0, text1):
    dist = Levenshtein.distance(text0, text1)
    max_len = max(len(text0), len(text1))
    return dist / max_len

関数名は normalized_levenshtein_distance にしたかったが流石に長すぎるので少し短縮。
ただ、これでも長いですね。(pep8的につらい)

これによって、例えば通常のレーベンシュタイン距離では、
“バニラ” と “アイス”　の組み合わせも “チョコレート” と “チョコレートアイス” もどちらも距離は3でしたが、
標準化した距離を使うことで、
前者は距離1、後者は 1/3(=0.333…)と算出されるようになり、前者の方が離れてると見なせます。

使ってみた結果がこちら。


print(Levenshtein.distance("アイス", "バニラ"))  # 3
print(Levenshtein.distance("チョコレートアイス", "チョコレート"))  # 3
print(normalized_distance("アイス", "バニラ"))  # 1.0
print(normalized_distance("チョコレートアイス", "チョコレート"))  # 0.3333333333333333

2019-05-062019-05-04

pythonでトピックモデル(LDA)

この前の記事で、scikit-learnのニュース記事のジャンルをロジスティック回帰で予測するというモデルを作ってみました。

参考:scikit-learnのニュースデータをロジスティック回帰で分類

今回はアプローチを変えて、トピックモデルを試してみようと思います。
どちらかというと、20newsのデータセットでもう少し何かやりたいというのが主目的で、
トピックモデルの理論的な説明については今回は省略します。
興味のあるかたへは、講談社から出ている岩田具治先生の、トピックモデル (機械学習プロフェッショナルシリーズ)
が非常にわかりやすかったのでおすすめです。ページ数も少なめでありがたい。
(数式が多くて書くのが大変なのですがゆくゆくは時系列分析みたいにこのブログでも説明したい。)

さて、pythonでトピックモデルを実装するには gensim を使うのが一般的のようです。
gensim topic modelling for humans
ただ、今回はいつも使っているscikit-learnでやってみました。
(gensimはword2vec等で使ってるのですがscikit-learnに比べると少し苦手。)

scikit-learnでトピックモデルを実装するために読むドキュメントはこちら。
sklearn.decomposition.LatentDirichletAllocation
Topic extraction with Non-negative Matrix Factorization and Latent Dirichlet Allocation

サンプルコードと同じことをしてもしょうがないので、少し工夫をしています。
・サンプルデータのカテゴリーを前回の記事同様に5個に絞る(その代わりそのカテゴリの全データを使用)
・カテゴリーごとに各文章のトピッックを可視化

前置きが長くなりましたが、やってみましょう。
必要ライブラリーのインポートとデータの読み込み


from sklearn.feature_extraction.text import CountVectorizer
from sklearn.decomposition import LatentDirichletAllocation
from sklearn.datasets import fetch_20newsgroups
import matplotlib.pyplot as plt
import numpy as np

remove = ('headers', 'footers', 'quotes')
categorys = [
        "rec.sport.hockey",
        "soc.religion.christian",
        "sci.med",
        "comp.windows.x",
        "talk.politics.mideast",
    ]
twenty_news = fetch_20newsgroups(
                                subset='all',
                                remove=remove,
                                categories=categorys
                            )
X = twenty_news.data

続いて、単語の出現頻度を数え、LDAのモデルを構築して学習します。
トピック数はカテゴリー数と同じ5でも試したのですが、
どうやら6か7にして、あまり重要でない単語を引き受けるトピックを作った方が納得性の高いものになりました。
サンプルコードは20カテゴリーを10トピックでうまく処理できているのに何故だろう？


# 単語の出現頻度データを作成
tf_vectorizer = CountVectorizer(max_df=0.90, min_df=5, stop_words='english')
tf = tf_vectorizer.fit_transform(X)
len(tf_vectorizer.get_feature_names())

# LDAのモデル作成と学習
lda = LatentDirichletAllocation(
                            n_components=7,
                            learning_method='online',
                            max_iter=20
                        )
lda.fit(tf)

それでは、学習した7個のトピックについて、それぞれの頻出語をみてみます。


features = tf_vectorizer.get_feature_names()

for tn in range(7):
    print("topic #"+str(tn))
    row = lda.components_[tn]
    words = ', '.join([features[i] for i in row.argsort()[:-20-1:-1]])
    print(words, "\n")

出力は下記の通りです。(乱数の影響で、モデルの学習をやり直すと結果は変わります。)

topic #0
god, people, think, don, know, just, like, does, say, believe, jesus, church, time, way, did, christ, things, good, christian, question

topic #1
25, 10, 11, 12, 14, 16, 15, 17, 20, 13, 18, 19, 55, 30, la, period, 24, 21, pit, 92

topic #2
armenian, armenians, turkish, people, turkey, armenia, turks, greek, genocide, russian, azerbaijan, government, history, muslim, university, soviet, war, 000, ottoman, killed

topic #3
game, don, said, team, just, didn, hockey, like, know, went, year, time, games, think, got, people, going, did, ll, came

topic #4
israel, jews, jewish, israeli, arab, state, people, world, right, public, arabs, rights, human, war, anti, peace, adl, states, country, palestinian

topic #5
medical, health, disease, cancer, patients, use, new, hiv, doctor, season, good, treatment, years, aids, high, drug, number, time, information, vitamin

topic #6
edu, use, file, window, com, server, program, dos, windows, available, motif, using, version, widget, sun, set, display, mit, x11, information

#1があまり意味のない数値を引き受けてくれていますが、
それ以外は、トピックごとに、宗教や国際的な話題、スポーツに医療に、コンピューターなどの単語が分類されています。

最後に、元の各テキストが、カテゴリーごとに妥当なトピック(話題)を持つと判定さているのか可視化してみてみましょう。
どんな可視化方法が一番わかりやすいか色々試したのですが、カラーマップが比較的良さそうでしたので紹介します。
(このほか箱ひげ図などもそこそこ綺麗に特徴が出ましたが。)


topic_data = lda.transform(tf)
fig = plt.figure(figsize=(6, 25))
for i in range(5):
    ax = fig.add_subplot(6, 1, 1+i)
    im = ax.pcolor(topic_data[twenty_news.target == i], vmax=1, vmin=0)
    fig.colorbar(im)
    # 軸の設定
    ax.set_xticks(np.arange(7) + 0.5, minor=False)
    ax.set_xticklabels(np.arange(7))
    ax.set_title(twenty_news.target_names[i])
plt.show()

出力がこちら。