読書記録(2020年8,9月)

行動を変えるデザイン

元々ファスト&フローとか予測どおりに不合理とか行動経済学の話が好きなので、それをプロダクトにどう適用していくのかが気になって読んだ。行動変容が必要なアプリをつかってもらうために、どういった障壁があり、人の行動特性(ここが行動経済学の部分)を踏まえてその障壁を超える人を増やすにはどういうことをしたら良いかということが書いてある。そこそこページ数が多いが、後半は前半でいったことの繰り返しが多かったように感じた。

コンテナ物語

流通コストを大きく下げたコンテナ技術がどのような経緯で使われるようになったかという話。今まで個別に詰めていたものをコンテナに詰め替えるだけでは当然済むわけではなく、コンテナを効率的に使用するための船や港の設計、この技術により職がなくなる港の労働者や労働組合との交渉やコンテナの規格についての取り決めなど、有効であるとわかっていてもアイデアが実際に広く使用されるまでには多くの障壁があり労力と時間が必要。はんこがなかなか撤廃できないのと同じような事象かな思う。

OKR

OKRは目標管理方法の1つでObjectives(目標)とKey Results(指標)のこと。前半部分は小説で後半部分が解説。基本的には実施することに優先順をつけることと、4半期に1回という短いサイクルで行うこと、また人事考査と切り離して達成が半々くらい目標を設定することで、各個人の能力を発揮しかつストレートに会社業績に直結させるような仕組み。人事考査とOKRが切り離されている、かつ組織に楽をしようとしている人がいない等の性善説が前提としてあるので、自走力が高い人が多い組織では有効そうかなという印象を受けた。

Pythonクローリング&スクレイピング[増補改訂版] -データ収集・解析のための実践開発ガイド

今までスクレイピングはweb上の資料をかいつまんでやっていたが、本著をもっと前に買っておけばよかった。スクレイピング自体の説明に留まらず、unixコマンドを用いた基本的なスクレイピングの考えたかや、python自体に関する説明、サーバにおけるスクレイピングコードの運用まで網羅されており一連の流れが理解しやすかった。

たった一人の分析から事業は成長する 実践 顧客起点マーケティング

主題として適切なn=1を選定しその人に合うように解像度を上げて戦略を考えるという話がある。以前読んだ「平均思考は捨てなさい」にも書かれていたが、セグメントで考えてしまうと特徴は平均値にまとまるが、実はすべての要素が平均的な人物というのは存在しないために、ありえないペルソナをつくる可能性もあるので、理にかなっていると思った。ただ同時にn=1の選び方を誤ると筋違いの仮説を立ててしまう可能性があるので、そこは丁寧に行う必要がありそうだなと思った。 またマーケティングという言葉の意味を今まであまり理解していなかったが、売り上げを最大化するためにプロダクトをいかにして顧客に届けるかという活動であり、必要であればプロダクト側への働きかけも行うものだということを理解できた。

われわれはなぜ嘘つきで自信過剰でお人好しなのか

進化心理学という分野の本。要は人間が遺伝子を残す上で、嘘つきで自信過剰でお人好しであることが重要な要素だったからであるという内容。 その性質のためにリーダーに適任だと思えた人が合理的に考えると不適であったり、立場がその人の性質を変えてしまったりもするので、 この人類のもつ特性を踏まえた上で、それが現在の世の中において必要なものなのかどうかは考える必要がある。

UNIXという考え方―その設計思想と哲学

UNIXという考え方―その設計思想と哲学

UNIXという考え方―その設計思想と哲学

  • 作者:Mike Gancarz
  • 発売日: 2001/02/01
  • メディア: 単行本

大きな1つの作品を作るのではなく、小さな部品の組み合わせを考える、試作をできるだけ早く作る、その時の最大効率を考えるのではなく、周囲の環境が変わっても使用できるような移植性を重視する、等のunixの設計思想はビジネスの展開の仕方と似ているなあと思った。

「意思決定」の科学 なぜ、それを選ぶのか

表が出続ける限り続けられるコイン投げのように、極めて少ない確率で極めて大きな利益が得られるような事例では期待値が発散するというサンクトペテルブルクパラドックスを例としてリスク回避性などの人間の持つバイアスを効用関数として定式化し意思決定の説明がされている。 この本では内容の検証を実験ベース行なっているが、個人的には大規模データに適用することに興味がある。 ちなみにサンクトペテルブルクはロシアの西側、北欧に近い場所に位置している。

確率思考 不確かな未来から利益を生みだす

長年ポーカーで生計を立てていた著者が、ポーカーで学んだ確率的な意思決定についての思考を、ポーカーだけでなく政治やスポーツの事例について記している。意思決定を行う場合、選択によって確率を考えることに加え、過去より現在を重視する、結論ありきで理由を考える、という人間自身のもつバイアスを考慮する必要があるという話。可能な限り公正な視点であらゆる確率とその結果を考えた方が良いという話。

DTW(Dynamic Time Warping)で台風軌道をクラスタリングする

はじめに

多次元時系列データのクラスタリングがしたいと思って探していたところ、 ちょうどこちらのブログの題材が台風軌道のクラスタリングという、多次元時系列かつ系列長の異なるデータをクラスタリングするというものだったので、理解を兼ねて同じ内容をpythonで実施してみたのが今回の内容になります。

参考資料

題材と内容を参考にさせていただいたブログ

DTWについてのわかりやすい資料

気象庁の台風データ

tsleanのドキュメント

実装

気象庁から2019年に発生した台風のデータを取得します。 下記のように、各台風に対して6時間毎の緯度経度の時系列データが存在します。

f:id:rmizutaa:20200926085018p:plain

foliumで可視化するとこのような感じ。2019年は全部で29の台風がありました。

import folium

map_osm = folium.Map(location=[35.6,139.7],zoom_start=3)
for num in df["台風番号"].unique():
    point=np.array(df[df["台風番号"]==num][["緯度","経度"]])
    folium.PolyLine(
            locations=point,                   
            ).add_to(map_osm)
map_osm

f:id:rmizutaa:20200926085113p:plain

この台風軌道を3つのグループにクラスタリングします。 tslearnのドキュメントがわかりやすいですが、 系列長が同じかつある特定の時間での比較が重要である場合はユークリッド距離を用いたkmeansでもクラスタリングができますが、 例えば1期分ずれているが、全体的に動きが類似しているような時系列をクラスタリングする場合はDTW(Dynamic Time Warping)が有効になります。

DTWの考え方としては、2つの時系列データに対しすべての点の組み合わせの距離を計算し、系列同士の距離が最短となる経路を探索するものになります。 サンプルとして出てくるものは1次元のデータが多いですが、組み合わせ計算のときはユークリッド距離を使用することが多いので、多次元データにも対応できます。 今回のケースだと緯度と経度という2次元の特徴に対してDTWを用いたクラスタリングを実施する形となります。

調べるとpythonでもいくつかパッケージが見つかるのですが、tsleanというパッケージの信頼性が高そうだったので、今回はこのパッケージを用いてクラスタリング及び結果の可視化を行います。

from tslearn.clustering import TimeSeriesKMeans
from tslearn.utils import to_time_series_dataset
from tslearn.barycenters import dtw_barycenter_averaging

points=[]
for num in df["台風番号"].unique():
    points.append(np.array(df[df["台風番号"]==num][["緯度","経度"]]))
points = to_time_series_dataset(points)

dba_km = TimeSeriesKMeans(n_clusters=3,
                          n_init=5,
                          metric="dtw",
                          verbose=True,
                          max_iter_barycenter=10,
                          random_state=22)

pred = dba_km.fit_predict(points)


colors=["blue","yellow","green","red"]
c0=[]
c1=[]
c2=[]

map_osm = folium.Map(location=[35.6,139.7],zoom_start=3)
for num,c in zip(df["台風番号"].unique(),pred):
    if c==0:
        c0.append(np.array(df[df["台風番号"]==num][["緯度","経度"]]))
    elif c==1:
        c1.append(np.array(df[df["台風番号"]==num][["緯度","経度"]]))
    else:
        c2.append(np.array(df[df["台風番号"]==num][["緯度","経度"]]))
 
    point=np.array(df[df["台風番号"]==num][["緯度","経度"]])
    folium.PolyLine(
            locations=point,     
            color=colors[c],
            ).add_to(map_osm)

map_osm

f:id:rmizutaa:20200926101217p:plain

右下の太平洋上で発生して、そこから東南アジアへ向かうもの、日本の左側を通るもの、日本の東側を通るもの、という別れ方をしているようです。 それほど別れ方に違和感はないように感じます。 元記事ではDTWを使った距離計算をした後にk-medoidsを用いてクラスタリングを行なっているのですが、 今回の手法はDTWで距離計算をした後にDBA(dtw barycenter averaging)で重心計算を行い、その重心を用いてkmeansでクラスタリングを行なっています。

クラスタの重心点をDBAで計算したものを赤線で引いてみます。

dba_c0=dtw_barycenter_averaging(c0)
dba_c1=dtw_barycenter_averaging(c1)
dba_c2=dtw_barycenter_averaging(c2)

folium.PolyLine(
            locations=dba_c0,     
            color=colors[3],
            ).add_to(map_osm)

folium.PolyLine(
            locations=dba_c1,     
            color=colors[3],
            ).add_to(map_osm)

folium.PolyLine(
            locations=dba_c2,     
            color=colors[3],
            ).add_to(map_osm)

map_osm

f:id:rmizutaa:20200926101235p:plain

かなりガタガタしていますが、大まかな特徴を捉えるくらいには使えるでしょうか。

まとめ

DTWを用いて台風軌道のクラスタリングを行いました。 系列長が不定で多次元なデータもクラスタリングを行うことができるので、例えばgpsデータなど探せば応用できる対象が結構あるのではないかなぁと思っています。

使用したコードは下記になります。

blog/dtw_typhoon at master · rmizuta3/blog · GitHub

Dashで更新可能な地図を表示する

はじめに

選択した項目に対して、インタラクティブに地図を更新するwebアプリを作成したいと思いました。 最初はstreamlitのpydeckで試していたのですが、地図を複数回レンダリングすると発生するバグが解消できなかったため、今回はDashで実装しました。

実施事項とデータについて

住民基本台帳の都道府県間の住所移動のデータを用いて、 移動後の都道府県を選択した場合に、各都道府県の移動前の人数を地図上に可視化します。

地図上にレンダリングするのに必要な都道府県のgeojsonファイルはこちらのサイトのものを使用させていただきました。

地図描写に使用する都道府県間の移動データはこのような形になります。 全国・都道府県・大都市が移動先の都道府県、valueが人数になります。

f:id:rmizutaa:20200815212925p:plain

実装

import plotly.express as px
import dash
import dash_html_components as html
import dash_core_components as dcc
from dash.dependencies import Input, Output
import pandas as pd
import json
import geopandas as gpd
import numpy as np


GEOJONFILE = "japan.geojson"
DATAFILE = "FEH.csv"

df = pd.read_csv(DATAFILE)
df = df[(df["性別"] == "総数") & (df["時間軸(年次)"] == "2019年") & (df["国籍"] == "移動者")]
df = df[df["value"] != "-"]
df["value"] = df["value"].astype(int)  # colorに使用するカラムは数値型に

jsonfile = gpd.read_file(GEOJONFILE)
dfjson = json.loads(jsonfile.to_json())

app = dash.Dash(__name__)
app.layout = html.Div([
    html.Div(
        html.H1('都道府県間の住所移動者(2019)',
                style={'textAlign': 'center'})
    ),
    dcc.Dropdown(
        id='selectplace',
        options=[{'label': i, 'value': i}
                 for i in df["全国・都道府県・大都市"].unique()],
        placeholder="移動先都道府県を選択",
        value="東京都"
    ),
    dcc.Loading(
        dcc.Graph(id='japanmap')
    )
])

@ app.callback(
    dash.dependencies.Output('japanmap', 'figure'),
    [dash.dependencies.Input('selectplace', 'value')])
def update_output(selectplace_value):
    selectdf = df[df["全国・都道府県・大都市"] == selectplace_value]
    fig = px.choropleth_mapbox(selectdf,
                               geojson=dfjson,
                               locations=selectdf['移動前の住所地'],
                               color=selectdf["value"],
                               featureidkey='properties.nam_ja',
                               color_continuous_scale="Viridis",
                               mapbox_style="carto-positron",
                               zoom=5,
                               center={"lat": 36, "lon": 138},
                               opacity=0.5,
                               labels={"value": "人数"}
                               )

    fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0})
    return fig

if __name__ == '__main__':
    app.run_server(debug=True)

https://github.com/rmizuta3/blog/tree/master/dash_chropleth

実装メモ

  • px.choropleth_mapboxについて
    • selectdfのlocationsとdfjsonのfeatureidkeyのカラムを紐づけて地図表示をしている。
    • colorに入る値は数値型にしておかないとレンダリングが終わらなくなる。
    • 地図のレンダリング速度は渡すjsonfileの大きさに依存する。今回の都道府県単位のデータだと2~3秒くらいだが、市区町村単位などより細かい地図を使用したい場合はjsonfileの方もデータを絞る必要が出てくる。

結果

このように、selectboxに入力した都道府県に対する各都道府県からの移動人数で色分けされた図ができました。 マウスオーバで各県からの移動人数の数値も確認できます。

f:id:rmizutaa:20200816131422g:plain

Whooshで検索機能の向上を図る

はじめに

前に類似本検索システムを作成したのですが、その中で数万ある本の候補の中から探したい本の検索する部分があります。 そのときは入力された単語に対し検索を全書籍に対して行う、という最も単純な手法を実装したのですが、 もう少しいいやり方がないかなーといくつか資料を読んで改善を実施したのでその過程を記述します。

参考資料

現状の問題点

  • 複数のワードを入力することができず、OR検索やNOT検索もできない
  • 登録されている書籍を全検索しているので、件数が増えた場合に検索時間が線形に増える。
  • 検索が一致した後のリストの返し方に何も優先順位をつけていない

解決方法

pythonで利用できる全文検索パッケージのwhooshを使います。商用ではサーバ機能を併せ持つElasticsearch等が使われることが多そうだったのですが、今回はheroku上に実装したいという意図もありこちらを使用しました。ちなみにドキュメントの類はElasticsearchの方が豊富そうでした。

検索を実施する最も簡単な方法は、検索したい文字列を検索対象すべてに対し順に検索をかけることですが、これだと検索対象の数だけ時間がかかってしまいます。それを防ぐために単語単位でインデックスを持たせるという手法があります。 例えば「犬も棒に当たる」という検索文字列に対し、「犬」、「棒」、「当たる」というそれぞれの単語に対して事前に対応する検索対象が何かを別々に記録しておくことで、走査が必要な数が全体での単語のユニーク数に収まるようになりますし、複数の単語の検索結果を比較することでAND検索やOR検索も行えるようになります。これを行うためには形態素解析を事前に行う必要があります。

実は今回実験の対象とするデータだと、検索先の文書件数が2万弱とデータ量が多くなく全検索でも速度はあまり気にならないので、Whooseを使うメリットはANDやORを使うことができるようになることと、ランキングを考慮できることになります。

実装メモ

検索部分の実装はこのような形になりました。 各機能の詳しい説明については公式を参照するのが確実です。

## indexの型の作成
schema = Schema(title=TEXT(stored=True), content=TEXT(stored=True, analyzer=StandardAnalyzer(stoplist=None)),count=NUMERIC(stored=True, sortable=True))
if not os.path.exists("heroku_index"):
    os.mkdir("heroku_index")
ix = create_in("heroku_index", schema)

# index作成
writer = ix.writer()
for num in range(len(bookdict_all_sort_upper10)):
    #表層形
    titlewords = set([m.surface() for m in tokenizer_obj.tokenize(
        bookdict_all_sort_upper10[num][0], mode)])
    #正規化
    titlewords = titlewords.union(set([m.normalized_form(
    ) for m in tokenizer_obj.tokenize(bookdict_all_sort_upper10[num][0], mode)]))
    #辞書表記
    titlewords = titlewords.union(set([m.dictionary_form(
    ) for m in tokenizer_obj.tokenize(bookdict_all_sort_upper10[num][0], mode)]))
    writer.add_document(title=bookdict_all_sort_upper10[num][0], content=" ".join(list(titlewords-remove_words)),
                        count=bookdict_all_sort_upper10[num][1])
writer.commit()

# 検索
outarr=[]
q = qp.parse(searchword_parsed)
with ix.searcher() as s:  # with内で処理する必要あり
    results = s.search(q, limit=999, sortedby="count",
                        reverse=True)  # 本の登場数でソート
    for i in results:
        outarr.append(i.values()[2]) #検索結果を格納

実装メモ

  • 形態素解析について

    • sudachiを使用
    • sudachiでは形態素解析の方法に3つの選択肢があるが、固有名詞を抽出できそうなSplitMode.Cを使用している。
    • 表記揺れや類似語に対応できることを期待して表層型と辞書型と正規化型の3つの形態素解析の結果を用いている。
  • 全文検索エンジンについて

    • whooshを使用
    • 後で格納した値を取り出したい場合はスキーマ作成時にstored_Trueとしておく必要がある。
    • 元が英語ベースでデフォルトだと形態素解析後に1文字のものはstopwordとして除外されてしまうため、1文字の検索を有効化したいときはスキーマ作成時にcontext=analyzer=StandardAnalyzer(stoplist=None)としておく必要がある。
    • 検索後に値を取り出す際はwith内で実施する必要がある。
    • whooshはデフォルトだとBM25というTFIDFに似たアルゴリズムでランキングを返しているが、今回は検索ワードの類似性よりも本の登場回数を用いたかったため、検索のsortedbyの部分を本の登場回数に変更している。

結果

OR検索の例 f:id:rmizutaa:20200803085943p:plain

NOT検索の例 f:id:rmizutaa:20200803090138p:plain

今回の更新によりスペース区切りの複数ワード入力や、OR検索やNOT検索ができるようになりました。 また、本の登録冊数でソートしているので、目的の本が上位に来やすくなったかと思います。

まとめ

whooshを使用してサイトの検索機能の向上を図りました。

いろいろなサイトの検索サジェストを見ていると、 入力された単語とは直接関係はないが検索したい対象に近い結果が出てくるような仕組みもあったので、 また時間があればそのあたりも調べてみたいです。

使用したコードはこちらになります。(コード部のみ)

github.com

読書記録(2020年6,7月)

今まで読んだ本の感想をたまにtwitterで流していましたが、140文字制限がきついのと後で見返しにくいので、こちらにまとめることにしてみます。

図書館情報学オタクと学ぶ検索エンジニア入門

booth.pm

情報検索の分野は初見だったが、具体的な例を元に情報検索とは何か、どのようなツールを使うのか、評価はどうするのか、どうビジネス的な価値を出すのかが記述されており、この分野でやっていることの片鱗を知ることができた。

データ分析のための数理モデル入門-本質をとらえた分析のために

統計や機械学習だけでなく、強化学習やエージェントベースモデルまで、モデルを作成する方法について、かなり網羅的に触れられている。1項目あたり数ページなので、各手法については触りを紹介するような感じ。既に知識がある場合は、説明資料を作成する際の参考になるし、自分の知らない分野もあるはずなので、そこに知識のindexをはる助けにもなるかなと思った。

測りすぎ――なぜパフォーマンス評価は失敗するのか

測定した指標を政治やビジネスなどに適用した場合に結果がどうなったかという例示集。測定したものをどう判断するかということが大事であり、測定すること自体が目的となったり、指標をハックするような行為がまかり通るようになるのは望ましくないとのことだった。確かに単純に売上を指標としてもそれは長期的な投資を軽視することになるように、本来複雑なことを単一の数字に落とすのは見逃す要素があることを意味する。どこまで測定が必要でその測定によりわかることはどこまでか、ということを認識しておくことが重要。

Team-Geek-―Googleギークたちはいかにしてチームを作るのか

googleエンジニアの技術ではなく人間に焦点を当てた本。 技術者といえど、会社員として働く上での生産性を高めるためにはチームやプロジェクトの他メンバとの協業が必須であり、 その対応が重要であることはわかるんだが人への働きかけは難しいんだよな…。

ザ・ゴール

ザ・ゴール

ザ・ゴール

適切な目標を設定し、その目標を達成するためのボトルネックに着目する。それ以外の部分については工数に遊びが出てもかまわない、という内容。工場における生産効率を例にしているが、他分野でも同じことが言えると思う。作業をするのは楽しいが、他に優先度の高い項目はないか、この作業により一体いくらの金額を得ることができるんだろうか?という意識は頭の片隅に常に置いておきたい。あとがきの著者はある自社商品の宣伝目的でこの本を書いたが、この本の内容を正しく理解すると自社製品は適していない、という結果になったという話が一番面白かった。

データの見えざる手:ウエアラブルセンサが明かす人間・組織・社会の法則

ウェアラブルバイスで取得したセンサデータから、仕事効率や人の幸せを測るという取り組みや、その結果について記述されており、内容は非常に興味深かった。ただ個人的にはウェアラブルから取得できるデータからわかるのは人の性質のほんの1側面にすぎないと感じており、このセンサデータだけから汎用的な結論を出す部分については若干疑問が残った。

創作人物と実在人物の誕生日を比較する

はじめに

TLで下記のマンガやアニメ等のキャラクターの誕生日には偏りがあることがわかるヒートマップが流れてきました。

これを見て結構偏りがあるんだなーと思うと同時に、疑り深い私は実在の人物についてもある程度なら 出生日はコントロールできるはずで、本当にこの特性はマンガやアニメ等の創作人物に限ったものなのか? 実在人物と本当に差があるといえるのか?という疑問をもちました。 そこで実在人物と創作人物の誕生日を比較し、本当に違いがあるかの検証を行いました。

データ収集

実在人物の出生日については政府の人口動態調査の結果を利用します。ここから1999~2018年の20年分の出生日の情報を使用しました。 現状2歳〜22歳の誕生日の人のデータになりますが、他の年齢の人との差はないと仮定します。 取得できるcsvファイルに少し癖があるので、集計のためには少し整形を行う必要があります。 創作人物の出生日は元ツイと同じくこちらのサイトから取得しました。

取得データの確認

実在人物については21,373,093人、創作人物については39,267人のデータが集まりました。 このデータを集計し、実在人物、創作人物それぞれで誕生者数のTOP5,WORST5の日付のデータを確認します。

TOP5

  • 実在人物
日付 誕生者数 誕生者割合(%) 誕生者数/誕生者数平均
9月25日 66897 0.31 1.15
1月5日 66896 0.31 1.15
12月25日 66848 0.31 1.14
9月26日 66827 0.31 1.14
5月2日 66363 0.31 1.14
  • 創作人物
日付 誕生者数 誕生者割合(%) 誕生者数/誕生者数平均
7月7日 904 2.30 8.43
3月3日 553 1.40 5.15
5月5日 440 1.12 4.10
1月1日 393 1.00 3.66
4月1日 361 0.92 3.36

WORST5

  • 実在人物
日付 誕生者数 誕生者割合(%) 誕生者数/誕生者数平均
2月29日 12704 0.06 0.22
1月1日 40602 0.19 0.70
1月2日 40884 0.19 0.70
12月31日 42460 0.20 0.73
1月3日 43607 0.20 0.75
  • 創作人物
日付 誕生者数 誕生者割合(%) 誕生者数/誕生者数平均
5月26日 42 0.11 0.39
12月11日 43 0.11 0.40
1月26日 49 0.12 0.46
12月17日 50 0.13 0.47
5月29日 54 0.14 0.50

日にち自体が少ない2/29を除いて実在人物は平均値からの乖離が0.7~1.15倍なのに対し、創作人物は0.39~8.43倍とかなり幅があることがわかります。 また、TOPとWORSTで現れる日付も全く違ったものになります。

割合をプロットしてみます。

f:id:rmizutaa:20200729230604p:plain

x軸の実在人物は左端の2/29の一点をのぞいてほぼ同じ箇所に固まっていますが、 y軸の創作人物はかなり幅があり、特に最も上の7/7の点はかなり外れ値的な位置付けとなっています。 また、本来日数が少ないはずの2/29にその影響が伺えません。

最後に、もう結果は明らかな気もしていますが念の為、 実在人物の誕生日と創作人物の誕生日には関連があるかという仮説に対しての独立性の検定も行ってみました。 結果としてはp値は0.0で帰無仮説は棄却され、有意水準を1%とすると実在人物の誕生日と創作人物の誕生日には関連があるとはいえませんでした。 ただ自由度が366、有意水準を1%とした場合はカイ二乗の値が431を越えると棄却域に入るのですが、 最も値が大きい7月7日の創作人物のセルのカイ二乗の値が5188なので、これだけで棄却域に入っているんですよね。 366*2のセルの内1つでも異常な値があったら仮説が棄却されてしまうので、 やってはみましたが、自由度が高すぎる場合にこの検定を使うのはあまり妥当ではないかもしれません…。

結論

創作人物の誕生日には、現実的ではない偏りがありました。

使用したコードは以下になります。 github.com

スペクトラルクラスタリングを用いてグラフ構造のデータをクラスタリングする

はじめに

グラフ構造のデータをクラスタリングする方法について調べていて、 スペクトラルクラスタリングという手法が使えそうだったので、その実験結果を記述します。

参考資料

使用データ

20年4月にtwitter apiで取得したホロライブ(vtuber)のtwitterの共通フォロワーを元に作成した相関行列を利用します。 (少し古いデータですがこれ結構データ取得に時間がかかるので…)

f:id:rmizutaa:20200720214824p:plain

データはこのような感じで、ライバー間で共通のフォロワーが多い組み合わせのときに値が大きくなるようなデータとなっています。(フォロワーが全くかぶっていない場合に0、完全に一致している場合に1になります。)

共通フォロワーが0の組み合わせはなく全結合のグラフなので、グラフネットワークとして図示するとこんな感じです。なんの情報量もない…。

f:id:rmizutaa:20200719203623p:plain

実験

このグラフをいくつかの類似した特性をもつグループに分割するためにスペクトラルクラスタリングを行います。

参考資料[1][2]の記述を参考に計算を行います。 手順はこのようになります。

  1. 隣接行列を元に正規化ラプラシアン行列を作成

  2. ラプラシアン行列を固有値分解し、固有値の小さい方からクラスタリング分割数と同じ次元数分の固有ベクトルを抽出

  3. 抽出した固有ベクトルに対しkmeansを実施

ラプラシアン行列は正規化する場合としない場合があるようですが、正規化しなかった場合にはクラスタリングをうまく行えなかったため、今回は正規化したものを使用しています。

import pandas as pd
import numpy as np
from scipy import linalg
from sklearn.cluster import KMeans
#相関行列の読み込み
cor=pd.read_csv("./cor_matrix_hololive.csv")

#隣接行列
A=cor-np.diag([1 for _ in range(len(cor))])

#対角行列
D= np.diag(A.sum(axis=1))

#ラプラシアン行列
L = D-A

#正規化ラプラシアン行列
L_norm=((np.linalg.inv(D)**(1/2))) @ L @ ((np.linalg.inv(D)**(1/2)))

#固有値分解(昇順で結果が出る)
vals, vecs = linalg.eig(L_norm) 

#kmeans
kmeans = KMeans(n_clusters=6,random_state=0)
kmeans.fit(vecs[:,1:6])

PCAだと大きい固有値をもつ固有ベクトルだけを抽出することで次元削減をしているので、 この手法で小さい方の固有値から使用するのはなぜかと思っていたのですが、 固有値が小さいということは、その成分の分散が小さい、つまり全体に対する影響量が小さいということになります。今回は全結合しているグラフのうち、切断しても全体への影響が少ない部分で分割したいという意図があるので、小さい方から選択していると理解しました。また、固有値が0である個数がグラフが完全に分離している数に相当するために、使用する固有ベクトルの次元数はクラスタリングの分割数(=グラフカットを行いたい数)と同じにしているようです。

sklearnでもスペクトラルクラスタリングのモジュールがあります。 隣接行列を入力とする場合はaffinity='precomputed'とすればクラスタリングが行えます。

from sklearn.cluster import SpectralClustering
sc = SpectralClustering(6, affinity='precomputed', n_init=100,assign_labels='kmeans',random_state=0)
sc.fit(A) 

分割数6でクラスタリングを行った場合、クラスタリング結果は以下のようになりました。 今回の条件だと、手計算したものと、sklearnで計算した結果は一致しました。

クラスタ1 クラスタ2 クラスタ3 クラスタ4 クラスタ5 クラスタ6
不知火フレア 常闇トワ 大空スバル 友人A AZKi 紫咲シオン
さくらみこ 姫森ルーナ 大神ミオ ロボ子さん アキ・ローゼンタール 赤井はあと
星街すいせい 角巻わため 猫又おかゆ ホロライブ公式 夜空メル 夏色まつり
潤羽るしあ 天音かなた 百鬼あやめ ときのそら 癒月ちょこ 湊あくあ
兎田ぺこら 桐生ココ 戌神ころね 白上フブキ
白銀ノエル
宝鐘マリン

この結果をざっくり解釈するとこのような形になります。

ホロライブのメンバーはそれぞれの加入時期によって1~4期生等の名称がつけられているのですが、 今回の結果だと、この加入時期が大きな影響を与えているようです。 加入時期が同じということは、同時に認知される可能性が高いですし、コラボ配信などで絡む機会も多い傾向がありますので、 今回の結果がtwitterの共通フォロワー数を用いたことから考えると、比較的妥当性の高いものではないかと思われます。

ちなみに細かいパラメータなどで違いはあるようで、クラスタ数によっては手計算とsklearnの結果にある程度差異は出ます。

まとめ

今回は全結合かつ、各パスにそれほど顕著な差があるわけではないデータだったのですが、 スペクトラルクラスタリングを用いることで、思ったよりも解釈可能なクラスタリング結果を出すことができました。 グラフネットワークの分析は今までやってこなかったですが、結構面白そうです。 あと計算を追っていて線形代数を再履修したくなってきました…。

使用したコードはこちらになります。 github.com