正解データがあるときの教師なし文書クラスタリングの評価

# 概要
異なる方法で教師なしクラスタリングしたとき、どれが一番優れているかを検証するため、正解データを使って教師なしクラスタリングの評価をする方法を検討、実装しました。
具体的にはランド指数に近い指標から算出したクラスタリング類似度の評価することにして、TF-IDF、fastText、SentenceBERTを用いた場合の予測結果と正解データの類似度の傾向から、指標の妥当性を判断しました。
# 背景
アンケートの自由記述からどんな意見があったかを抽出するようなケースを想定します。このとき、教師なしクラスタリング（以下、単に「クラスタリング」と記述）をして、群ごとにざっくりと目を通す方法が考えられます。しかし、クラスタリングは一発でうまく行かないことも多く、いろんな手法のうちどれが適しているかを試行錯誤したくなることがあります。

しかし、クラスタリングの評価は難しいです。群内の凝集性（同じ群の要素同士がなるべく近くにあること）と群間の離散性（違う群の要素同士がなるべく遠くにあること）を使うことが多いですが、埋め込み空間のみでの評価が、人間から見て妥当な評価になるとは限らず、また、異なるクラスタリング手法の比較にも必ずしも適していません。

そこで正解データがある場合に、正解データとの類似度をみることで、クラスタリング手法の妥当性を評価することを考えます。正解データは例えば、アンケートデータの一部に人がラベリングしたようなものを想定します。正解データがあれば教師ありの学習をすればよいのではと思うかもしれませんが、自由記述の探索的な分析など、付与すべきラベルの種類が事前に定まっていない場合には、正解データにすべての種類のラベルが含まれるとは限らないため、教師ありよりも教師なしのほうが分析に適している可能性が高いです。

# ２つのクラスタリングの類似度評価
正解クラスタリングと予測クラスタリングの類似度を測ることで、予測クラスタリングの評価とすることを考えます。２つのクラスタリングの類似度を測る方法はいくつかありますが、今回は以下を試してみます。

[２つのクラスタリング結果がどのくらい似ているかの指標 | takemikami's note](https://takemikami.com/2019/02/25/clustdiff.html)

クラスタリング予測結果において各要素のペアを作り、同じクラスタに属するかどうかを判定します。この判定が正解データの同じペアに対する同様の判定とどの程度一致するかを見積もる方法です。一致度の指標にはF値を用いています。なお、accuracyの場合は、[ランド指数（rand index）](https://scikit-learn.org/stable/modules/generated/sklearn.metrics.adjusted_rand_score.html)という名前が付いているようですが、F値の場合の名前はよくわかりませんでした。しかし、accuracyよりもF値のほうが指標として良さそうな気がするので、今回はF値でいきます。

# 文書クラスタリング
正解が存在する文書データに対して、複数の方法でクラスタリングを実施し、上記類似度を計算します。文書クラスタリングは、文書の埋め込みとクラスタリングの２ステップでされることが多いです。そこで今回は文書の埋め込み手法の違いが上記類似度に与える影響を見てみます。クラスタリングはkmeansに固定し、埋め込みの方法について以下を比較します。

* トリグラムのTF-IDF
* 単語のTF-IDF
* 単語の分散表現（fastText）の平均
* SentenceBERT

また参考として教師あり学習による分類も評価します。
仮説としては教師あり学習が最も類似度が高く、SentenceBERT、fastText、単語TF-IDF、トリグラムTF-IDFと続くのではないかと思いますので、実際に類似度がそのように変化したら、ある程度妥当な指標と言えるのではないかと思います。

# 実装（関数の用意）
## 検証用データの用意
[livedoorニュースコーパス](https://www.rondhuit.com/download.html)をダウンロードし、タイトルとカテゴリの情報だけ抽出し、適当なパス（今回は`text/titles.csv`）に保存します。

```Python
# extract titles
DIR_PATH = "text" # 解凍したディレクトリへのパス
SUBDIR_NAMES = ["dokujo-tsushin", "it-life-hack", "kaden-channel", "livedoor-homme", "movie-enter", "peachy", "smax", "sports-watch", "topic-news"]

rows = []
for subdir_name in SUBDIR_NAMES:
  dir_path = Path(DIR_PATH, subdir_name)
  for file in dir_path.iterdir():
    # 条件に合わないファイル名のとき処理をスキップ
    if not file.stem.startswith(subdir_name): continue

    with open(file) as f:
      lines = f.read().splitlines()
    # ３行目がタイトル
    title = lines[2]
    rows.append({
      "title": title
      , "category": subdir_name
    })

# データフレームに変換してテキストに出力
pd.DataFrame(rows).to_csv(str(Path(DIR_PATH).joinpath("titles.csv")), index=False)
```

以下のようなファイルができます。

```:text/titles.csv
title,category
タニタに続き、第二弾！ 話題の社員食堂は家庭薬膳,dokujo-tsushin
マスクで顔を隠す独女たち,dokujo-tsushin
薄っぺらい親友関係が増殖中？ 女の友情はきょわい〜！？,dokujo-tsushin
【オトナ女子映画部】“隙だらけ”のウザい女に学ぶ男心のつかみ方『乱暴と待機』,dokujo-tsushin
女医が教えるオンナの体のウソホント vol.10「体の変化とホルモン分泌の周期」presented by ゆるっとcafe,dokujo-tsushin
結婚相手に妥協できるもの、できないもの　（男性編）,dokujo-tsushin
...
```

## 埋め込みの取得
埋め込みを取得するための関数を定義します。
まずライブラリを読み込みます。環境に存在しなければpipでインストールしておきます。

```Python
import pandas as pd
from pathlib import Path
from sklearn.model_selection import train_test_split
from sklearn.cluster import KMeans
from sklearn.feature_extraction.text import TfidfVectorizer
from sklearn import metrics
from sklearn.naive_bayes import MultinomialNB
from sklearn.pipeline import make_pipeline
from transformers import BertJapaneseTokenizer, BertModel
import torch
from sudachipy import dictionary
from sudachipy import tokenizer
import pickle
import fasttext
import fasttext.util
```

### トリグラムTF-IDF
単語分かち書きをせず、トリグラムでTF-IDFを計算します。scikit-learnのTfidfVectorizerを使います。

```Python
# tfidf
def ngram_tfidf(texts, *, ngram_range = (3,3)):
  vectorizer = TfidfVectorizer(
                    analyzer="char"
                    , ngram_range=ngram_range
                    , max_df=0.9
                    , min_df = 5)
  return vectorizer.fit_transform(texts)
```

### 単語TF-IDF
`Sudachi`で単語分かち書きしてから、TF-IDFを計算します。

```Python
def word_tfidf(texts, *, ngram_range = (1,1)):
  tokenizer_obj = dictionary.Dictionary(dict="full").create()
  mode = tokenizer.Tokenizer.SplitMode.A
  wakachi_texts = [" ".join([m.surface() for m in tokenizer_obj.tokenize(text, mode)]) for text in texts]
  vectorizer = TfidfVectorizer(
    analyzer = "word"
    , ngram_range = ngram_range
    , max_df = 0.9
    , min_df = 5
  )
  return vectorizer.fit_transform(wakachi_texts)
```

### fastText

* Sudachiで分かち書き
* 単語ごとの分散表現の平均値を計算

というステップです。
fastTextのモデルは事前にダウンロードしておいてください。
[Word vectors for 157 languages | fastText](https://fasttext.cc/docs/en/crawl-vectors.html)
ダウンロートおよびfasttext.load_modelの実行にまあまあ時間がかかります。

```Python
def fasttext_vector(texts, *, model=None, model_path = "fasttext/cc.ja.300.bin"):
  ft = model or fasttext.load_model(model_path)
  tokenizer_obj = dictionary.Dictionary(dict="full").create()
  mode = tokenizer.Tokenizer.SplitMode.A
  vectors = []
  for text in texts:
    tokens = tokenizer_obj.tokenize(text)
    words = [token.surface() for token in tokens]
    vec = ft.get_word_vector(words[0])
    for w in words[1:]:
      vec += ft.get_word_vector(w)
    mean_vec = vec / len(words)
    vectors.append(mean_vec)
  return vectors
```

### SentenceBERT
モデルは以下からリンクされているsentence-bert-base-ja-mean-tokens-v2を使わせていただきました。コードも同記事およびhuggingfaceのページを参考にさせていただきました。
[【日本語モデル付き】2020年に自然言語処理をする人にお勧めしたい文ベクトルモデル](https://qiita.com/sonoisa/items/1df94d0a98cd4f209051)

```Python
# 参考 https://qiita.com/sonoisa/items/1df94d0a98cd4f209051
class SentenceBertJapanese:
    def __init__(self, model_name_or_path, device=None):
        self.tokenizer = BertJapaneseTokenizer.from_pretrained(model_name_or_path)
        self.model = BertModel.from_pretrained(model_name_or_path)
        self.model.eval()

        if device is None:
            device = "cuda" if torch.cuda.is_available() else "cpu"
        self.device = torch.device(device)
        self.model.to(device)

    def _mean_pooling(self, model_output, attention_mask):
        token_embeddings = model_output[0] #First element of model_output contains all token embeddings
        input_mask_expanded = attention_mask.unsqueeze(-1).expand(token_embeddings.size()).float()
        return torch.sum(token_embeddings * input_mask_expanded, 1) / torch.clamp(input_mask_expanded.sum(1), min=1e-9)


    def encode(self, sentences, batch_size=8):
        all_embeddings = []
        iterator = range(0, len(sentences), batch_size)
        for batch_idx in iterator:
            batch = sentences[batch_idx:batch_idx + batch_size]

            encoded_input = self.tokenizer.batch_encode_plus(batch, padding="longest", 
                                           truncation=True, return_tensors="pt").to(self.device)
            model_output = self.model(**encoded_input)
            sentence_embeddings = self._mean_pooling(model_output, encoded_input["attention_mask"]).to('cpu')

            all_embeddings.extend(sentence_embeddings)

        # return torch.stack(all_embeddings).numpy()
        return torch.stack(all_embeddings)

def sentencebert(texts, *, model=None):
    MODEL_NAME = "sonoisa/sentence-bert-base-ja-mean-tokens-v2"  # <- v2です。
    model = model or SentenceBertJapanese(MODEL_NAME)
    sentence_embeddings = model.encode(texts, batch_size=8)
    return sentence_embeddings.detach().numpy() # 他の関数が扱いやすいようにnumpyに変換して返す
```

## クラスタリング用の関数

埋め込みベクトルのリストを入力として、kmeansクラスタリングの結果（ラベル）を返す関数です。クラスタ数はとりあえずライブドアコーパスの9をデフォルトにしています。

```Python
# k-meansでクラスタ分析。とりあえず9つのグループに分けてみる
def kmeans_clustering(vectors, *, n_clusters=9):
  km_model = KMeans(n_clusters=n_clusters, random_state = 0)
  km_model.fit(vectors)
  return km_model.labels_
```

## 教師あり学習（トリグラム、ナイーブベイズ）
比較用に教師あり学習したモデルで予測する関数も用意しておきます。
管理を楽にするため、関数内でモデルの学習も行っています。
（トリグラムのナイーブベイズなので、そこまで時間はかからなかったので、そうしました）

```Python
# 参考 
def supervised_naivebayse_vector(texts, *, model=None, train_text = train_df["title"], train_category=train_df["category"]):
  if model is None:
    model = make_pipeline(
      TfidfVectorizer(
                      analyzer="char"
                      , ngram_range=(3,3)
                      , max_df=0.9
                      , min_df = 5)
      , MultinomialNB()
    )
    model.fit(train_text, train_category)
  return model.predict(texts)
```
## 類似度評価用の関数
正解ラベルと予測ラベルを入力として類似度を返す関数です。
`get_pair_label`はラベルリストを入力として、要素ペア（組み合わせ）が同じクラスタに属しているかどうかのリストを返します。
`cluster_similarity`では`get_pair_label`で取得したリストのF値値やその他参考情報を出力します。

```Python
def get_pair_label(cluster):
  labels = []
  for i0, v0 in enumerate(cluster):
    for i1, v1 in enumerate(cluster):
      if i1<=i0: continue
      labels.append(v0==v1)
  return labels

def cluster_similarity(correct_cluster, test_cluster):
  correct_pairs = get_pair_label(correct_cluster)
  test_pairs = get_pair_label(test_cluster)
  combined_pairs = [(v0,v1) for v0, v1 in zip(correct_pairs, test_pairs)]

  correct_true, correct_false = correct_pairs.count(True), correct_pairs.count(False)
  test_true, test_false = test_pairs.count(True), test_pairs.count(False)
  true_positive = combined_pairs.count((True, True))
  false_positive = combined_pairs.count((False, True))
  true_negative = combined_pairs.count((False, False))
  false_negative = combined_pairs.count((True, False))

  scores = {
    "ct_cf_tt_tf": (correct_true, correct_false, test_true, test_false)
    , "tp_fp_tn_fn": (true_positive, false_positive, true_negative, false_negative)
    , "precision": metrics.precision_score(correct_pairs, test_pairs)
    , "recall": metrics.recall_score(correct_pairs, test_pairs)
    , "f1": metrics.f1_score(correct_pairs, test_pairs)
    , "accuracy": metrics.accuracy_score(correct_pairs, test_pairs)
  }
  return scores    
```

# 検証
## データの読み込み
事前に作った検証用データを読み込みます。
全部で7000件以上あり多いので、１割を検証に使います。
残りの9割は教師あり学習に使います。

```Python
FILE_PATH = "text/titles.csv"

# テストデータの読み込み
df = pd.read_csv(FILE_PATH)
# 検証を素早くできるようにテストデータ数を制限
train_df, test_df = train_test_split(df, train_size=0.9, random_state = 0, shuffle=True, stratify=df["category"])
# indexをリセット
train_df, test_df = train_df.reset_index(drop=True), test_df.reset_index(drop=True)

print("all")
print(df["category"].value_counts())
print("")
print("train")
print(train_df["category"].value_counts())
print("")
print("test")
print(test_df["category"].value_counts())
```

```
all
sports-watch      900
dokujo-tsushin    870
it-life-hack      870
movie-enter       870
smax              870
kaden-channel     864
peachy            842
topic-news        770
livedoor-homme    511
Name: category, dtype: int64

train
sports-watch      810
smax              783
it-life-hack      783
dokujo-tsushin    783
movie-enter       783
kaden-channel     777
peachy            758
topic-news        693
livedoor-homme    460
Name: category, dtype: int64

test
sports-watch      90
kaden-channel     87
smax              87
dokujo-tsushin    87
movie-enter       87
it-life-hack      87
peachy            84
topic-news        77
livedoor-homme    51
Name: category, dtype: int64
```

## トリグラムTF-IDF

```Python
# ngram tfidf, kmeans
X = ngram_tfidf(test_df["title"])
test_labels = kmeans_clustering(X)
cluster_similarity(test_df["category"], test_labels)
```

```
{'ct_cf_tt_tf': (30397, 240819, 172013, 99203),
 'tp_fp_tn_fn': (21785, 150228, 90591, 8612),
 'precision': 0.12664740455663234,
 'recall': 0.7166825673586209,
 'f1': 0.21525616323304184,
 'accuracy': 0.4143413367942894}
```

## 単語TF-IDF

```Python
# word tfidf, kmeans
X = word_tfidf(test_df["title"])
test_labels = kmeans_clustering(X)
cluster_similarity(test_df["category"], test_labels)
```

```
{'ct_cf_tt_tf': (30397, 240819, 85461, 185755),
 'tp_fp_tn_fn': (12802, 72659, 168160, 17595),
 'precision': 0.14979932366810592,
 'recall': 0.4211599828930487,
 'f1': 0.22099466588409952,
 'accuracy': 0.667224647513421}
```

## fastText

```Python
ft = fasttext.load_model('fasttext/cc.ja.300.bin')
vectors = fasttext_vector(test_df["title"], model = ft)
test_labels = kmeans_clustering(vectors)
cluster_similarity(test_df["category"], test_labels)
```

```
{'ct_cf_tt_tf': (30397, 240819, 33096, 238120),
 'tp_fp_tn_fn': (6382, 26714, 214105, 24015),
 'precision': 0.19283297075175249,
 'recall': 0.20995492976280555,
 'f1': 0.2010300348069866,
 'accuracy': 0.812957200165182}
```

## SentenceBERT
CPUの場合、推論に時間がかかるので、結果をpickleで保存しておき、２回目以降はそちらを読み込むだけですむようにしています。

```Python
# sentence bert, kmeans
embedding_binary_path = "embedding/sentencebert_embedding.pickle"
if Path(embedding_binary_path).exists():
  with open(embedding_binary_path, "rb") as f:
    sentence_embeddings = pickle.load(f)
else:
  sentence_embeddings = sentencebert(test_df["title"])
  with open("embedding/sentencebert_embedding.pickle", "wb") as f:
    pickle.dump(sentence_embeddings.detach().numpy(), f)

test_labels = kmeans_clustering(sentence_embeddings)
cluster_similarity(test_df["category"], test_labels)
```

```
{'ct_cf_tt_tf': (30397, 240819, 32532, 238684),
 'tp_fp_tn_fn': (8420, 24112, 216707, 21977),
 'precision': 0.25882208287224884,
 'recall': 0.277001019837484,
 'f1': 0.2676031718285687,
 'accuracy': 0.8300653353784437}
```

## 教師あり学習
教師ありでいきなりラベルを予測しているので、kmeansは不要です。

```Python
test_labels = supervised_naivebayse_vector(test_df["title"])
cluster_similarity(test_df["category"], test_labels)
```

```
{'ct_cf_tt_tf': (30397, 240819, 31432, 239784),
 'tp_fp_tn_fn': (19530, 11902, 228917, 10867),
 'precision': 0.6213413082209214,
 'recall': 0.6424976148962068,
 'f1': 0.6317423862588752,
 'accuracy': 0.9160484632175093}
```

# 結果と考察
F値だけ並べると以下の結果でした。

|手法|F値|
|---|---|
|トリグラムTF-IDF|0.215|
|単語TF-IDF|0.220|
|fastText|0.209|
|SentenceBERT|0.267|
|教師あり学習|0.631|

fastTextが予想外にTF-IDFよりも低かったですが、ほかは仮説通りの順位になりました。
よって、類似度の指標として悪くはなさそうな気がします。もちろん、一つのドメインのデータだけで結論を出すのは早計で、今後、いろんな手法や正解データで試して、頑健性をチェックする必要はありますが。

また、SentenceBERTでも教師あり学習のF値とは大きく開きがあり、正解データを再現できているかという観点ではあまり質が良くないことがわかります。今回は正解データに関するヒューリスティクスを全く入れていないので、実応用上無理のない範囲でそうしたものを取り込んでいくことで、精度が高まるかもしれません。

# おわりに
正解データが存在する場合の教師なしクラスタリングの評価方法について検討しました。ランド指数のaccuracyをF値で置き換えた指標で試してみたところ、埋め込み手法の差によるクラスタリングの評価をある程度仮説通りに行えたように思える結果となりました。

今後の課題として以下を考えています。

* 埋め込みのバリエーション：LDAなど
* クラスタリングのバリエーション：ward法など。合わせて最適なクラスタ数の決定方法についても要実装
* 評価指標のバリエーション：purity inversed-purityのF値、BCubedなど
* 純粋な教師なしだと限界がありそうなので、正解データに関するヒューリスティクスを導入することで結果が向上するか検討

以上です。
同じ興味を持つ人に向けて、なにかの参考になれば幸いです。

# 参考
今回の実装のGitHub
* https://github.com/JiroShimaya/ClusteringSimilarity

文章のベクトル化
* [SudachiPy](https://github.com/WorksApplications/SudachiPy/blob/develop/docs/tutorial.md)
* [tf-idfでベクトル化したラジオ感想ツイートをクラスタリングして可視化する](https://note.com/himaratsu/n/necefee6e5454)
* [機械学習 〜 テキスト分類（ナイーブベイズ分類器） 〜](https://qiita.com/fujin/items/39d450b910bf2be866b5)
* [fastTextとDoc2Vecのモデルを作成してニュース記事の多クラス分類の精度を比較する](https://qiita.com/kazuki_hayakawa/items/ca5d4735b9514895e197)
* [【日本語モデル付き】2020年に自然言語処理をする人にお勧めしたい文ベクトルモデル](https://qiita.com/sonoisa/items/1df94d0a98cd4f209051)
* [https://huggingface.co/sonoisa/sentence-bert-base-ja-mean-tokens-v2](https://huggingface.co/sonoisa/sentence-bert-base-ja-mean-tokens-v2)

エラー・デバッグ関係
* [Pytorch: Can't call numpy() on Variable that requires grad. Use var.detach().numpy() instead](https://stackoverflow.com/questions/55466298/pytorch-cant-call-numpy-on-variable-that-requires-grad-use-var-detach-num)
  * sentence-bertの出力をkmeansに入力したらエラーが出たときの解消方法


クラスタリング結果の比較方法
* [２つのクラスタリング結果がどのくらい似ているかの指標 | takemikami's note](https://takemikami.com/2019/02/25/clustdiff.html)
* [Precision and recall for clustering?](https://stats.stackexchange.com/questions/15158/precision-and-recall-for-clustering/80194)
* [Evaluation of clustering](https://nlp.stanford.edu/IR-book/html/htmledition/evaluation-of-clustering-1.html)
* [クラスタリング＋二次元マッピングによる可視化を定量評価する](https://qiita.com/muraokaz/items/e02c2e440d285062b6d2)
* https://scikit-learn.org/stable/modules/clustering.html
* （注意：いきなりダウンロードされます）[特徴抽出を目的とした文書クラスタからの 一貫性阻害要素除去](https://ipsj.ixsq.nii.ac.jp/ej/index.php?action=pages_view_main&active_action=repository_action_common_download&item_id=92744&item_no=1&attribute_id=1&file_no=1&page_id=13&block_id=8)