{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Clustering text documents using k-means\n\nThis is an example showing how the scikit-learn API can be used to cluster\ndocuments by topics using a [Bag of Words approach](https://en.wikipedia.org/wiki/Bag-of-words_model).\n\nTwo algorithms are demonstrated, namely :class:`~sklearn.cluster.KMeans` and its more\nscalable variant, :class:`~sklearn.cluster.MiniBatchKMeans`. Additionally,\nlatent semantic analysis is used to reduce dimensionality and discover latent\npatterns in the data.\n\nThis example uses two different text vectorizers: a\n:class:`~sklearn.feature_extraction.text.TfidfVectorizer` and a\n:class:`~sklearn.feature_extraction.text.HashingVectorizer`. See the example\nnotebook `sphx_glr_auto_examples_text_plot_hashing_vs_dict_vectorizer.py`\nfor more information on vectorizers and a comparison of their processing times.\n\nFor document analysis via a supervised learning approach, see the example script\n`sphx_glr_auto_examples_text_plot_document_classification_20newsgroups.py`.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Authors: The scikit-learn developers\n# SPDX-License-Identifier: BSD-3-Clause"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Loading text data\n\nWe load data from `20newsgroups_dataset`, which comprises around 18,000\nnewsgroups posts on 20 topics. For illustrative purposes and to reduce the\ncomputational cost, we select a subset of 4 topics only accounting for around\n3,400 documents. See the example\n`sphx_glr_auto_examples_text_plot_document_classification_20newsgroups.py`\nto gain intuition on the overlap of such topics.\n\nNotice that, by default, the text samples contain some message metadata such\nas `\"headers\"`, `\"footers\"` (signatures) and `\"quotes\"` to other posts. We use\nthe `remove` parameter from :func:`~sklearn.datasets.fetch_20newsgroups` to\nstrip those features and have a more sensible clustering problem.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import numpy as np\n\nfrom sklearn.datasets import fetch_20newsgroups\n\ncategories = [\n    \"alt.atheism\",\n    \"talk.religion.misc\",\n    \"comp.graphics\",\n    \"sci.space\",\n]\n\ndataset = fetch_20newsgroups(\n    remove=(\"headers\", \"footers\", \"quotes\"),\n    subset=\"all\",\n    categories=categories,\n    shuffle=True,\n    random_state=42,\n)\n\nlabels = dataset.target\nunique_labels, category_sizes = np.unique(labels, return_counts=True)\ntrue_k = unique_labels.shape[0]\n\nprint(f\"{len(dataset.data)} documents - {true_k} categories\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Quantifying the quality of clustering results\n\nIn this section we define a function to score different clustering pipelines\nusing several metrics.\n\nClustering algorithms are fundamentally unsupervised learning methods.\nHowever, since we happen to have class labels for this specific dataset, it is\npossible to use evaluation metrics that leverage this \"supervised\" ground\ntruth information to quantify the quality of the resulting clusters. Examples\nof such metrics are the following:\n\n- homogeneity, which quantifies how much clusters contain only members of a\n  single class;\n\n- completeness, which quantifies how much members of a given class are\n  assigned to the same clusters;\n\n- V-measure, the harmonic mean of completeness and homogeneity;\n\n- Rand-Index, which measures how frequently pairs of data points are grouped\n  consistently according to the result of the clustering algorithm and the\n  ground truth class assignment;\n\n- Adjusted Rand-Index, a chance-adjusted Rand-Index such that random cluster\n  assignment have an ARI of 0.0 in expectation.\n\nIf the ground truth labels are not known, evaluation can only be performed\nusing the model results itself. In that case, the Silhouette Coefficient comes in\nhandy. See `sphx_glr_auto_examples_cluster_plot_kmeans_silhouette_analysis.py`\nfor an example on how to do it.\n\nFor more reference, see `clustering_evaluation`.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from collections import defaultdict\nfrom time import time\n\nfrom sklearn import metrics\n\nevaluations = []\nevaluations_std = []\n\n\ndef fit_and_evaluate(km, X, name=None, n_runs=5):\n    name = km.__class__.__name__ if name is None else name\n\n    train_times = []\n    scores = defaultdict(list)\n    for seed in range(n_runs):\n        km.set_params(random_state=seed)\n        t0 = time()\n        km.fit(X)\n        train_times.append(time() - t0)\n        scores[\"Homogeneity\"].append(metrics.homogeneity_score(labels, km.labels_))\n        scores[\"Completeness\"].append(metrics.completeness_score(labels, km.labels_))\n        scores[\"V-measure\"].append(metrics.v_measure_score(labels, km.labels_))\n        scores[\"Adjusted Rand-Index\"].append(\n            metrics.adjusted_rand_score(labels, km.labels_)\n        )\n        scores[\"Silhouette Coefficient\"].append(\n            metrics.silhouette_score(X, km.labels_, sample_size=2000)\n        )\n    train_times = np.asarray(train_times)\n\n    print(f\"clustering done in {train_times.mean():.2f} \u00b1 {train_times.std():.2f} s \")\n    evaluation = {\n        \"estimator\": name,\n        \"train_time\": train_times.mean(),\n    }\n    evaluation_std = {\n        \"estimator\": name,\n        \"train_time\": train_times.std(),\n    }\n    for score_name, score_values in scores.items():\n        mean_score, std_score = np.mean(score_values), np.std(score_values)\n        print(f\"{score_name}: {mean_score:.3f} \u00b1 {std_score:.3f}\")\n        evaluation[score_name] = mean_score\n        evaluation_std[score_name] = std_score\n    evaluations.append(evaluation)\n    evaluations_std.append(evaluation_std)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## K-means clustering on text features\n\nTwo feature extraction methods are used in this example:\n\n- :class:`~sklearn.feature_extraction.text.TfidfVectorizer` uses an in-memory\n  vocabulary (a Python dict) to map the most frequent words to features\n  indices and hence compute a word occurrence frequency (sparse) matrix. The\n  word frequencies are then reweighted using the Inverse Document Frequency\n  (IDF) vector collected feature-wise over the corpus.\n\n- :class:`~sklearn.feature_extraction.text.HashingVectorizer` hashes word\n  occurrences to a fixed dimensional space, possibly with collisions. The word\n  count vectors are then normalized to each have l2-norm equal to one\n  (projected to the euclidean unit-sphere) which seems to be important for\n  k-means to work in high dimensional space.\n\nFurthermore it is possible to post-process those extracted features using\ndimensionality reduction. We will explore the impact of those choices on the\nclustering quality in the following.\n\n### Feature Extraction using TfidfVectorizer\n\nWe first benchmark the estimators using a dictionary vectorizer along with an\nIDF normalization as provided by\n:class:`~sklearn.feature_extraction.text.TfidfVectorizer`.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from sklearn.feature_extraction.text import TfidfVectorizer\n\nvectorizer = TfidfVectorizer(\n    max_df=0.5,\n    min_df=5,\n    stop_words=\"english\",\n)\nt0 = time()\nX_tfidf = vectorizer.fit_transform(dataset.data)\n\nprint(f\"vectorization done in {time() - t0:.3f} s\")\nprint(f\"n_samples: {X_tfidf.shape[0]}, n_features: {X_tfidf.shape[1]}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "After ignoring terms that appear in more than 50% of the documents (as set by\n`max_df=0.5`) and terms that are not present in at least 5 documents (set by\n`min_df=5`), the resulting number of unique terms `n_features` is around\n8,000. We can additionally quantify the sparsity of the `X_tfidf` matrix as\nthe fraction of non-zero entries divided by the total number of elements.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "print(f\"{X_tfidf.nnz / np.prod(X_tfidf.shape):.3f}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We find that around 0.7% of the entries of the `X_tfidf` matrix are non-zero.\n\n\n### Clustering sparse data with k-means\n\nAs both :class:`~sklearn.cluster.KMeans` and\n:class:`~sklearn.cluster.MiniBatchKMeans` optimize a non-convex objective\nfunction, their clustering is not guaranteed to be optimal for a given random\ninit. Even further, on sparse high-dimensional data such as text vectorized\nusing the Bag of Words approach, k-means can initialize centroids on extremely\nisolated data points. Those data points can stay their own centroids all\nalong.\n\nThe following code illustrates how the previous phenomenon can sometimes lead\nto highly imbalanced clusters, depending on the random initialization:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from sklearn.cluster import KMeans\n\nfor seed in range(5):\n    kmeans = KMeans(\n        n_clusters=true_k,\n        max_iter=100,\n        n_init=1,\n        random_state=seed,\n    ).fit(X_tfidf)\n    cluster_ids, cluster_sizes = np.unique(kmeans.labels_, return_counts=True)\n    print(f\"Number of elements assigned to each cluster: {cluster_sizes}\")\nprint()\nprint(\n    \"True number of documents in each category according to the class labels: \"\n    f\"{category_sizes}\"\n)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "To avoid this problem, one possibility is to increase the number of runs with\nindependent random initiations `n_init`. In such case the clustering with the\nbest inertia (objective function of k-means) is chosen.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "kmeans = KMeans(\n    n_clusters=true_k,\n    max_iter=100,\n    n_init=5,\n)\n\nfit_and_evaluate(kmeans, X_tfidf, name=\"KMeans\\non tf-idf vectors\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "All those clustering evaluation metrics have a maximum value of 1.0 (for a\nperfect clustering result). Higher values are better. Values of the Adjusted\nRand-Index close to 0.0 correspond to a random labeling. Notice from the\nscores above that the cluster assignment is indeed well above chance level,\nbut the overall quality can certainly improve.\n\nKeep in mind that the class labels may not reflect accurately the document\ntopics and therefore metrics that use labels are not necessarily the best to\nevaluate the quality of our clustering pipeline.\n\n### Performing dimensionality reduction using LSA\n\nA `n_init=1` can still be used as long as the dimension of the vectorized\nspace is reduced first to make k-means more stable. For such purpose we use\n:class:`~sklearn.decomposition.TruncatedSVD`, which works on term count/tf-idf\nmatrices. Since SVD results are not normalized, we redo the normalization to\nimprove the :class:`~sklearn.cluster.KMeans` result. Using SVD to reduce the\ndimensionality of TF-IDF document vectors is often known as [latent semantic\nanalysis](https://en.wikipedia.org/wiki/Latent_semantic_analysis) (LSA) in\nthe information retrieval and text mining literature.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from sklearn.decomposition import TruncatedSVD\nfrom sklearn.pipeline import make_pipeline\nfrom sklearn.preprocessing import Normalizer\n\nlsa = make_pipeline(TruncatedSVD(n_components=100), Normalizer(copy=False))\nt0 = time()\nX_lsa = lsa.fit_transform(X_tfidf)\nexplained_variance = lsa[0].explained_variance_ratio_.sum()\n\nprint(f\"LSA done in {time() - t0:.3f} s\")\nprint(f\"Explained variance of the SVD step: {explained_variance * 100:.1f}%\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Using a single initialization means the processing time will be reduced for\nboth :class:`~sklearn.cluster.KMeans` and\n:class:`~sklearn.cluster.MiniBatchKMeans`.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "kmeans = KMeans(\n    n_clusters=true_k,\n    max_iter=100,\n    n_init=1,\n)\n\nfit_and_evaluate(kmeans, X_lsa, name=\"KMeans\\nwith LSA on tf-idf vectors\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We can observe that clustering on the LSA representation of the document is\nsignificantly faster (both because of `n_init=1` and because the\ndimensionality of the LSA feature space is much smaller). Furthermore, all the\nclustering evaluation metrics have improved. We repeat the experiment with\n:class:`~sklearn.cluster.MiniBatchKMeans`.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from sklearn.cluster import MiniBatchKMeans\n\nminibatch_kmeans = MiniBatchKMeans(\n    n_clusters=true_k,\n    n_init=1,\n    init_size=1000,\n    batch_size=1000,\n)\n\nfit_and_evaluate(\n    minibatch_kmeans,\n    X_lsa,\n    name=\"MiniBatchKMeans\\nwith LSA on tf-idf vectors\",\n)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Top terms per cluster\n\nSince :class:`~sklearn.feature_extraction.text.TfidfVectorizer` can be\ninverted we can identify the cluster centers, which provide an intuition of\nthe most influential words **for each cluster**. See the example script\n`sphx_glr_auto_examples_text_plot_document_classification_20newsgroups.py`\nfor a comparison with the most predictive words **for each target class**.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "original_space_centroids = lsa[0].inverse_transform(kmeans.cluster_centers_)\norder_centroids = original_space_centroids.argsort()[:, ::-1]\nterms = vectorizer.get_feature_names_out()\n\nfor i in range(true_k):\n    print(f\"Cluster {i}: \", end=\"\")\n    for ind in order_centroids[i, :10]:\n        print(f\"{terms[ind]} \", end=\"\")\n    print()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### HashingVectorizer\nAn alternative vectorization can be done using a\n:class:`~sklearn.feature_extraction.text.HashingVectorizer` instance, which\ndoes not provide IDF weighting as this is a stateless model (the fit method\ndoes nothing). When IDF weighting is needed it can be added by pipelining the\n:class:`~sklearn.feature_extraction.text.HashingVectorizer` output to a\n:class:`~sklearn.feature_extraction.text.TfidfTransformer` instance. In this\ncase we also add LSA to the pipeline to reduce the dimension and sparcity of\nthe hashed vector space.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from sklearn.feature_extraction.text import HashingVectorizer, TfidfTransformer\n\nlsa_vectorizer = make_pipeline(\n    HashingVectorizer(stop_words=\"english\", n_features=50_000),\n    TfidfTransformer(),\n    TruncatedSVD(n_components=100, random_state=0),\n    Normalizer(copy=False),\n)\n\nt0 = time()\nX_hashed_lsa = lsa_vectorizer.fit_transform(dataset.data)\nprint(f\"vectorization done in {time() - t0:.3f} s\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "One can observe that the LSA step takes a relatively long time to fit,\nespecially with hashed vectors. The reason is that a hashed space is typically\nlarge (set to `n_features=50_000` in this example). One can try lowering the\nnumber of features at the expense of having a larger fraction of features with\nhash collisions as shown in the example notebook\n`sphx_glr_auto_examples_text_plot_hashing_vs_dict_vectorizer.py`.\n\nWe now fit and evaluate the `kmeans` and `minibatch_kmeans` instances on this\nhashed-lsa-reduced data:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fit_and_evaluate(kmeans, X_hashed_lsa, name=\"KMeans\\nwith LSA on hashed vectors\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fit_and_evaluate(\n    minibatch_kmeans,\n    X_hashed_lsa,\n    name=\"MiniBatchKMeans\\nwith LSA on hashed vectors\",\n)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Both methods lead to good results that are similar to running the same models\non the traditional LSA vectors (without hashing).\n\n## Clustering evaluation summary\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\nimport pandas as pd\n\nfig, (ax0, ax1) = plt.subplots(ncols=2, figsize=(16, 6), sharey=True)\n\ndf = pd.DataFrame(evaluations[::-1]).set_index(\"estimator\")\ndf_std = pd.DataFrame(evaluations_std[::-1]).set_index(\"estimator\")\n\ndf.drop(\n    [\"train_time\"],\n    axis=\"columns\",\n).plot.barh(ax=ax0, xerr=df_std)\nax0.set_xlabel(\"Clustering scores\")\nax0.set_ylabel(\"\")\n\ndf[\"train_time\"].plot.barh(ax=ax1, xerr=df_std[\"train_time\"])\nax1.set_xlabel(\"Clustering time (s)\")\nplt.tight_layout()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        ":class:`~sklearn.cluster.KMeans` and :class:`~sklearn.cluster.MiniBatchKMeans`\nsuffer from the phenomenon called the [Curse of Dimensionality](https://en.wikipedia.org/wiki/Curse_of_dimensionality) for high dimensional\ndatasets such as text data. That is the reason why the overall scores improve\nwhen using LSA. Using LSA reduced data also improves the stability and\nrequires lower clustering time, though keep in mind that the LSA step itself\ntakes a long time, especially with hashed vectors.\n\nThe Silhouette Coefficient is defined between 0 and 1. In all cases we obtain\nvalues close to 0 (even if they improve a bit after using LSA) because its\ndefinition requires measuring distances, in contrast with other evaluation\nmetrics such as the V-measure and the Adjusted Rand Index which are only based\non cluster assignments rather than distances. Notice that strictly speaking,\none should not compare the Silhouette Coefficient between spaces of different\ndimension, due to the different notions of distance they imply.\n\nThe homogeneity, completeness and hence v-measure metrics do not yield a\nbaseline with regards to random labeling: this means that depending on the\nnumber of samples, clusters and ground truth classes, a completely random\nlabeling will not always yield the same values. In particular random labeling\nwon't yield zero scores, especially when the number of clusters is large. This\nproblem can safely be ignored when the number of samples is more than a\nthousand and the number of clusters is less than 10, which is the case of the\npresent example. For smaller sample sizes or larger number of clusters it is\nsafer to use an adjusted index such as the Adjusted Rand Index (ARI). See the\nexample\n`sphx_glr_auto_examples_cluster_plot_adjusted_for_chance_measures.py` for\na demo on the effect of random labeling.\n\nThe size of the error bars show that :class:`~sklearn.cluster.MiniBatchKMeans`\nis less stable than :class:`~sklearn.cluster.KMeans` for this relatively small\ndataset. It is more interesting to use when the number of samples is much\nbigger, but it can come at the expense of a small degradation in clustering\nquality compared to the traditional k-means algorithm.\n\n"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.9.21"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}