plotly (JS) お試し (Scatterplot Matrix もどき)

直接サポートされてなくても

Subplots で模擬出来るはずだよなぁ、て話。

散布図行列ってネタがなかなかいいのがないのよね

なんで散布図行列なんて描きたいか、てのは、統計解析に慣れてるなら自明ではあるんだけれど、簡単に説明すると、「なんかわからんけど関連しそうかどうか、まずは眺めてみたい」というような、「本格的な分析を始める前に真っ先にやってみること」として行うもの。なので相関係数を計算してみたりだとか、その他もろもろ「高等な」分析を始める前に、手っ取り早く視覚化するのに使う。(「風が吹いたら桶屋が儲かる、のような思わぬ相関は、ありそうか、なさそうか」てこと。)

作為的過ぎる「説明のためだけのデータ」をご用意して「ほーらね」してみたところで、何言いたいのかわからんと思うので、少しは現実的な例で見せたいが、現実にそんなに都合のいいデータなんか、そうそうあるはずもなく。だけれどもそれでも「現実の例」にしたい。

というわけで、「比較的」簡単に個人が入手出来るデータと、簡単に計算出来るほかの物理量の中で、こんなデータで遊んでみようと思う:

要するに「月齢と交通事故に関連はありそうか、気象条件との関連はありそうか」ということを「知りたい」というストーリーね。明確な関連が出るとは思えないけれど、とにかく「そう思ったとしてやってみたとしたらどうか」という話。

「人身交通事故日別発生状況」がとにかく扱いづらい形でなぁ、PDF で入手したにせよ上のリンクの Google がキャッシュした html であるにせよ、データを簡単に取り出せるようになってない。仕方ないのでキャッシュの html をお取り寄せ後に以下の Python プログラムでデータを取り出すようにしてみた:

 1 # -*- coding: utf-8 -*-
 2 from __future__ import unicode_literals
 3 
 4 import io
 5 import re
 6 from datetime import datetime
 7 
 8 
 9 _PAT = re.compile(r'<div style="position:absolute;top:(\d+);left:(\d+)"><nobr>(.*)</nobr></div>')
10 
11 
12 def to_structured(fn):
13     hy = int(re.match(r"H(\d+)\.txt", fn).group(1))
14     year = (2016 - 28) + hy
15 
16     starting = False
17     by_top = {}  # top: {left: col}
18     for line in io.open(fn, encoding="utf-8").readlines():
19         m = _PAT.match(line)
20         if m:
21             top, left, col = m.group(1, 2, 3)
22             top, left = int(top), int(left)
23             #
24             if not starting and col == u"1日":
25                 starting = True
26             elif starting and col == u"合計":
27                 break
28             if starting:
29                 #print(repr((top, left, col)))
30                 if top not in by_top:
31                     by_top[top] = {}
32                 by_top[top][left / 30 * 30] = col
33     #
34     result = {}
35     collefts = None  # left, left, ...
36     for d, top in enumerate(sorted(by_top.keys())):
37         if not collefts:
38             collefts = list(sorted(by_top[top].keys()))
39         tmp = [
40             by_top[top].get(left)
41             for left in collefts[1:]]
42         result.update(dict(
43             [
44                 ("{:04d}-{:02d}-{:02d}".format(year, (i + 3) / 3, d + 1),
45                  dict(zip(("件数", "死者数", "負傷者数"), tuple(map(int, tmp[i:i + 3])))))
46                 for i in range(0, len(tmp), 3)
47                 if tmp[i] is not None
48              ]
49             ))
50     return result
51 
52 def to_structured_all():
53     from glob import glob
54     result = {}
55     for fn in glob("H*.txt"):
56         result.update(to_structured(fn))
57     return result
58 
59 #import json
60 #import sys, codecs ; sys.stdout = codecs.getwriter("utf-8")(sys.stdout)
61 #print(json.dumps(to_structured_all(), indent=2, sort_keys=True, ensure_ascii=False))

気象庁のデータの方もそんなに扱いやすい形ではないけれど、一応ちゃんと csv なので、こっちは説明はしない。

月齢、というか「phase as percent of the surface illuminated」は pyEphem を使うと例えばこんなんで求まる:

1 date = "2017-01-03"
2 
3 from datetime import datetime
4 import ephem
5 
6 moon = ephem.Moon()
7 moon.compute(ephem.Date(datetime.strptime(date, "%Y-%m-%d")))
8 moon.phase

この phase は、まぁ女性なら「月経」に関係するんだろうし、そうでなくても「潮汐」に直結するので、「これが何かに関係しているという想像することは、そう大外ししてない」ものではあるんだけれど、まぁ交通事故と関連付けるのは難しいだろうな、というネタなのね。まさにこういう「関係ありそうな気がするもの」が本当にそうなのかをざっくり知りたい、という目的に、散布図行列を使う、というわけだ。

ともあれ、ここまでの作業で、こんなデータを手にしてみた:

 1 {
 2   "2012-01-01": {
 3     "illuminated": 47.6, 
 4     "件数": 24, 
 5     "平均気温(℃)": 6.3, 
 6     "平均風速(m/s)": 2.4, 
 7     "日照時間(時間)": 6.5, 
 8     "死者数": 0, 
 9     "負傷者数": 47, 
10     "降水量の合計(mm)": 0.0, 
11     "降雪量合計(cm)": 0.0
12   }, 
13   "2012-01-02": {
14     "illuminated": 57.1, 
15     "件数": 35, 
16     "平均気温(℃)": 6.4, 
17     "平均風速(m/s)": 5.0, 
18     "日照時間(時間)": 5.3, 
19     "死者数": 0, 
20     "負傷者数": 48, 
21     "降水量の合計(mm)": 0.0, 
22     "降雪量合計(cm)": 0.0
23   }, 
24   /* ... snip ... */
25   "2016-12-31": {
26     "illuminated": 2.9, 
27     "件数": 37, 
28     "平均気温(℃)": 7.2, 
29     "平均風速(m/s)": 2.2, 
30     "日照時間(時間)": 8.7, 
31     "死者数": 5, 
32     "負傷者数": 48, 
33     "降水量の合計(mm)": 0.0, 
34     "降雪量合計(cm)": 0.0
35   }
36 }

これをどうにか散布図行列にしてみよう。(とはいえ全部やると図が大き過ぎてわけわからなくなるので、実例では全部ではやらない。)実際散布図行列に使いたいデータは、以下のように、日付情報ではなく子のキーでまとめあげたこれである:

 1 {
 2   "illuminated": [
 3     6.6, 
 4     86.9, 
 5     18.1, 
 6     25.9, 
 7     54.8, 
 8     45.1, 
 9     34.5, 
10     /* ... */
11   ], 
12   "件数": [
13     101, 
14     92, 
15     84, 
16     45, 
17     104, 
18     126, 
19     80, 
20     53, 
21     87, 
22     /* ... */
23   ], 
24   /* ... */
25 }

出来てみた

貼り付けるとバカでか過ぎて開くのが大変なので、データは平成28年分だけ:

「フレームのソースを表示」で見てくれ、というのは正直キビシイし、はっきりいってアタシ自身が超絶に迷走したのでちゃんと説明しないわけにはいかない。

とにかく「散布図マトリクスの直接のサポート」がないからには「自力でサブプロットを駆使する」ことになるわけなんだけれど、その「サブプロット」のやり方が独特というか、「プア」というか。結構長い時間試行錯誤しててやっとこ出来た。「わかってしまえば簡単」かというと、そうでもない気がするが、こんな感じ:

「traces」内の xaxis, yaxis にはとにかく連番
「layout」内の domain が位置決め(%)
「layout」内、anchor は 1. で連番にしてるのでシンプルに「xaxis10 なら anchor は y10」

ただし「うまく出来たパターン」というだけであって、多分「こうせねばならん」てことではないと思う。(どうせ単に名前で対応付けているんであろうし。) にしても「anchor」の意味がさっぱりわからん。なんなんだろうかこれは。

散布図をちょっとだけ解読してみる

「plotly.js で散布図マトリクス」という技術ネタなので別に分析はどうでもいいちゃぁどうでもいいんだけれど、ちょっとは見とこうかと。

illuminated (月齢に近いもの) と事故件数なんだけど、よーくみると 0 と 100 の付近の「母集団」が多いのよね。これは「事故件数が多い/少ない」とは関係なくて、たぶん「交通安全週間」のようなものが月の特定の日付(1日、15日など)にあることに関係してるんじゃないかと思う。

風速と件数は、そもそも風速が大きい場合母数が小さくなるのは、「ヒドい天気だから外に出ない」ことに関係してんだろね。事故件数そのものも右肩下がりにみえなくもない。

とにかく詳細にみれば何か発見出来るかもしれないけれど、「見慣れてないと何がなんだか」だとも思うし、あと注意して欲しいのは、散布図マトリクスは「相関してて当たり前」のものも出てしまいがちてことね。よくあるのが「あるパラメータ z は x と y に基づく計算結果」の場合に z と x、z と y に図から相関を読み取ってしまって「うぁ、大発見」と思ってしまったりすること。非常に基礎的な話だけど、「わかっている従属関係がないもの」どうしだけを選ばないとダメですわよ(直交・独立)。