Pandasでログ統計処理スクリプト（その１）
特にブログに書くことがなくネタ切れであります。今回はコンテストでもなく移動運用でもなく自作ネタでもない、ソフトウエアの話題にしましょうか。
長年無線運用をやっていて交信ログがある程度の量になってくると、交信統計というのをやってみたくなるものです。例えば、バンド・モードごとやDX/国内エリア別の交信数の統計は運用スタイルやアクティビティを反映するものなので分析するとなかなか興味深いものです。
Turbo HAMLOGやCTestWinなどのログソフトでも統計分析が可能ですが、細かいところには手が届いていません。例えば、移動局と固定局の割合を調べたり、記念局の交信数を調べること（簡単には）できないようになっています。ログデータを解析するのはさほど難しくないので、PerlやRuby, Pythonなどの文字列解析が得意なスクリプト言語に親和性のある方は、自分で解析スクリプトを作るのが良いかと思いますし、当局もRubyやPythonの勉強がてら分析スクリプトを自分で作っています。
実際にプログラミングしたり開発環境を調べているとデータ分析用のライブラリとして有名なPandasを使うとなかなか便利だということに気が付きました。PandasはPythonで動作するデータ解析支援ライブラリで、科学技術数値データを扱うときに便利です。今回は交信ログの解析に使うので、かなりもったいない使い方ではありますが、後に述べるようにメソッドが充実しているので、こちらのプログラミングが非常に楽になります。せっかくですから最新プログラミングツールを無線活動に活用してみましょう。では例を少しお見せしたいと思います。
データはTurbo HAMLOGのCSV形式とします。CTestWinが出力するHAMLOG形式でもOKです。ご自分の交信ログからCSV形式でログを変換してみてください。他のログソフトでも少しカスタマイズすれば使えます。
次にPythonの開発環境を手に入れてください。WindowsやMacでも良いですが、Unix(Linux)のほうがいろいろと便利です。Pythonはバージョン2と3があるのですが、2でもOKですが、できれば3を選んでください。
次に、IPythonを用意してください。Pythonのインタラクティブツールで、これがあると割と便利です。IPythonも（できれば）バージョン3を選んでください。次にライブラリのPandasも必要となります。これらツールとライブラリはPIPで入手できます。すでにインストールされているライブラリはシェル上でこのようにして確認できます。
~$ pip3 list
pandasがこのリストに含まれていなかったら、
~$ pip3 install pandas
としてお使いの環境にPandasをPython3用にインストールしてください。
次に、以下のスクリプトをlogstat.pyという名前のファイルにセーブします。そのスクリプトと同ディレクトリにHAMLOGのcsvファイルを置いておきます。
# logstat.py import sys import re import pandas as pd from datetime import datetime, timedelta columns_list = ['call', 'datetime', 'time', 'urrst', 'myrst', 'freq', 'mode', 'jcc', 'gl', 'qsl', 'name', 'qth', 'urnr', 'mynr', 'rem3'] log = pd.read_csv(sys.argv[1], names=columns_list)
シェル上で、ipython3を実行します。
~$ ipython3 Python 3.5.2 (default, Sep 10 2016, 08:21:44) Type "copyright", "credits" or "license" for more information. IPython 5.1.0 -- An enhanced Interactive Python. ? -> Introduction and overview of IPython's features. %quickref -> Quick reference. help -> Python's own help system. object? -> Details about 'object', use 'object??' for extra details. In [1]:
となるので、先のスクリプトとCSVログファイル(2012WPXCW.csv)を指定して実行します。
In [1]: run logstat.py 2012WPXCW.csv In [2]:
となると成功。次にインタラクティブモードでログデータをいろいろ見てみましょう。log.head()はログデータ（2次元データフレーム）の先頭5行を表示するメソッド。
In [2]: log.head()
Out[2]:
call datetime time urrst myrst freq mode jcc gl qsl name qth \
0 HS7AT 12/05/26 04:04U 599 599 21.0 CW NaN NaN J NaN NaN
1 R9UT 12/05/26 04:13U 599 599 21.0 CW NaN NaN J NaN NaN
2 E21EIC 12/05/26 04:14U 599 599 21.0 CW NaN NaN J NaN NaN
3 NH2T 12/05/26 04:15U 599 599 21.0 CW NaN NaN J NaN NaN
4 PT5T 12/05/26 04:16U 599 599 21.0 CW NaN NaN J NaN NaN
urnr mynr rem3
0 1 137 8
1 2 96 8
2 3 236 8
3 4 559 8
4 5 471 8
In [3]:
このように正しくCSVファイルを読み込んでpandasのデータフレームにおさめられていることがわかります。なお、これはCTestWinから出力させたCSVファイルの例で、HAMLOGから直接出力したものだとちょっと異なるかもしれませんので、カラム名などを適当にアレンジしてみてください。
・・・つづく・・・
Pandasでログ統計処理スクリプト（その2）
pandasを使うとたった一行でCSV形式のログファイルをデータフレームに収めることができます。データフレームは行と列のインデックスで管理するExcelファイルのようなものです。まずはQSO数を確かめてみましょう。これはデータフレームの行数を調べれば良いので、
In [3]: len(log) Out[3]: 500
などとします。他にもやり方がありますので調べてみてください。
さて、HAMLOG形式の文字列を取り込んでいるので不都合なことがあります。それは時刻です。'datetime'と'time'という列インデックスにとりあえず日付と時刻を当てはめてみましたが、このままでは文字列データに過ぎないので、大小の比較などができません。またJSTとUTCが混在している場合もあります。そこで、以下のようなスクリプトでUTCへ変換しかつdatetime形式にして'datetime'列に上書きしてしまいましょう。
date_list = log.datetime time_list = log.time datetime_list = [] for i, s in enumerate(date_list): date_str = s + '-' + time_list[i] logdate = datetime.strptime(date_str.rstrip('JZU'), '%y/%m/%d-%H:%M') if time_list[i].find('J') == 5: logdate = logdate - timedelta(hours = 9) datetime_list.append(logdate) log.datetime = datetime_list del log['time']; del log['jcc']; del log['gl']; del log['qsl'] del log['name']; del log['qth']; del log['rem3']
元のスクリプトlogstat.pyの後ろに上記スクリプトを書き加えてもう一度ipythonでrunさせます。log.head()でログデータの一部を見てみると、datetime列が若干変わっているはずです。GLやQTHなどの列もこの場合不要なのでdelコマンドで削除してしまいました。
In [5]: log.head()
Out[5]:
call datetime urrst myrst freq mode urnr mynr
0 HS7AT 2012-05-26 04:04:00 599 599 21.0 CW 1 137
1 R9UT 2012-05-26 04:13:00 599 599 21.0 CW 2 96
2 E21EIC 2012-05-26 04:14:00 599 599 21.0 CW 3 236
3 NH2T 2012-05-26 04:15:00 599 599 21.0 CW 4 559
4 PT5T 2012-05-26 04:16:00 599 599 21.0 CW 5 471
datetime列がpythonで時刻として扱える型になっています。確認したい場合は、列だけをリスト表示させます。長いのでやはり先頭の5つだけ表示させたいので、.head()をつけます。
In [7]: log.datetime.head() Out[7]: 0 2012-05-26 04:04:00 1 2012-05-26 04:13:00 2 2012-05-26 04:14:00 3 2012-05-26 04:15:00 4 2012-05-26 04:16:00 Name: datetime, dtype: datetime64[ns]
のようにリストの要素データの型がdatetime64となっているのを確認できます。
次に、コールサインからエリアを判別する関数call_area()を作っておいて（後述）、ログにあるコールを一つ一つその関数を呼び出してエリアの文字列をリストを作って、データフレームに'area'という列を新しく作って、加えるということをします。
callarea_list = [] for c in log.call: callarea_list.append(call_area(c)) log['area'] = callarea_list
In [11]: log.tail(30).head()
Out[11]:
call datetime urrst myrst freq mode urnr mynr area
470 P40H 2012-05-27 09:33:00 599 599 14.0 CW 473 1742 DX
471 JA2VHO 2012-05-27 11:19:00 599 599 7.0 CW 474 706 JA2
472 K7RF 2012-05-27 11:20:00 599 599 7.0 CW 475 849 DX
473 N6WIN 2012-05-27 11:23:00 599 591 7.0 CW 476 440 DX
474 KT5J 2012-05-27 11:27:00 599 599 7.0 CW 477 1353 DX
のように、area列に日本のコールにはエリア、日本以外のコールにはDXと識別されたデータを付与しておきます。以上で、データの加工は終わりました。次からログデータの統計処理です。
・・・つづく・・・
Pandasでログ統計処理スクリプト（その3）
コールサイン文字列から国内エリアと日本以外の局を判別する関数です。あまり厳密な方法ではないです（例えば/MMやゲストオペレータのコールなど）。pythonは正規表現のマッチングがシンプルに書けないし、後方参照が仕方がまどろっこしいので困ったものです。
def call_area(callsign): jacall_re = re.compile('(^J[AD-S]\d)|(^[78][J-N]\d)') area_ja_re = re.compile('^J[AD-S](\d)') area_7_re = re.compile('^7[K-N][1-4]') area_7j_re = re.compile('^7J(\d)') area_8_re = re.compile('^8[J-N](\d)') area_jd1_re = re.compile('^JD1') call_sp = callsign.split('/') if len(call_sp) >= 2: if jacall_re.search(call_sp[0]): m_areanum = re.match('(\d)', call_sp[1]) if m_areanum: area = 'JA' + m_areanum.group(1) else: area = 'DX' else: area = 'DX' else: if jacall_re.search(callsign): m_ja = area_ja_re.search(callsign) m_7 = area_7_re.search(callsign) m_7j = area_7j_re.search(callsign) m_8 = area_8_re.search(callsign) if m_ja: area = 'JA' + m_ja.group(1) elif m_7: area = 'JA1' elif m_7j: area = 'JA' + m_7j.group(1) elif m_8: area = 'JA' + m_8.group(1) else: area = 'DX' m_jd1 = area_jd1_re.search(callsign) if m_jd1: area = 'DX' return area
この関数を使って例えば、
In [12]: call_area('JF1DIR')
Out[12]: 'JA1'
In [13]: call_area('JF1DIR/9')
Out[13]: 'JA9'
In [14]: call_area('KH0/JF1DIR')
Out[14]: 'DX'
とエリアに対応する文字列を返します。前回の（その２）でこの関数を使って新たに'area'列にエリア文字列の列を作っておきました。
さて、ここからが本番です。
まずはこのログからどのバンドの運用があったのかを調べるのは、
In [15]: set(log.freq)
Out[15]: {3.5, 7.0, 14.0, 21.0, 28.0}
とします。log.freqでfreq列のリストが得られ、そのset（重複を除いたリスト）を表示させればよいのです。
次に、ログのすべてのデータから、エリア別の統計を取りたい場合は、
In [16]: log.area.value_counts() Out[16]: DX 419 JA1 29 JA3 12 JA2 11 JA0 9 JA7 7 JA6 6 JA4 4 JA8 2 JA5 1 Name: area, dtype: int64
で得られます。python/pandasの偉いところはIPythonのインタラクティブ表示やスクリプト内のprint()文で表示させるとちゃんと整形されて表示されるということです。ちなみに、このメソッドの返り値はSeries型になっているので、
In [17]: a = log.area.value_counts() In [18]: a[0] Out[18]: 419 In [19]: a.index[0] Out[19]: 'DX'
で個々の値に参照することができます。
例えば7MHzの運用分だけでエリア別統計を取りたい場合は、freq列にフィルターをかけてカウントします。
In [20]: log[log.freq==7].area.value_counts() Out[20]: DX 59 JA1 16 JA2 7 JA3 6 JA0 5 JA7 5 JA4 1 JA6 1 JA5 1 Name: area, dtype: int64
QSO数だけならば、
In [21]: log[log.freq==7].area.count() Out[21]: 101
でよいです。
QSO時間でログを切り取る、例えば、2012年5月26日18時〜24時までのデータだけを切り出す場合は、
In [22]: log[(log.datetime>=datetime(2012,5,26,18,0,0,0))&(log.datetime
のようにします。
call列をキーにログ全体をソートするには、
In [23]: log.sort('call').head()
Out[23]:
call datetime urrst myrst freq mode urnr mynr area
66 3V8BB 2012-05-26 06:25:00 599 599 21.0 CW 67 503 DX
33 3Z2X 2012-05-26 05:20:00 599 599 21.0 CW 34 68 DX
41 6M0NR 2012-05-26 05:32:00 599 599 21.0 CW 42 172 DX
165 6M0NR 2012-05-26 11:38:00 599 599 7.0 CW 167 192 DX
141 7J1YAJ 2012-05-26 10:55:00 599 599 7.0 CW 143 3039 JA1
ただし元のlogデータフレームには影響しない非破壊のメソッドです。
以上のようなデータフレームに対するメソッドを駆使すればやりたいことは大体出来てしまいます。位置から自分でプログラミングすると結構面倒ですが、上記メソッドをスクリプト内に組み込んであげれば各種統計解析が可能となります。いかがでしょう？pythonの概念がわかってしまうと意外と簡単ですよね。
