<jupyter_start><jupyter_text>**BrainStation Capstone Project****Author: <NAME>**  **Crafting the Perfect Beer:** Predicting Beer Ratings with Sentiment Analysis and Clustering  Part III: K-Means Clustering ---  Introduction This notebook does not require a boosting environment and can be run on a standard Python kernel.Now that I have looked at supervised learning techniques to predict a beer rating, I want to use K-Means Clustering, an unsupervised learning technique, to see if we can uncover any similar groups or patterns in the text reviews that have yet to be seen. Unlike Logistic Regression, Decision Tree, and XGBoost where we try to optimize the models against specific metrics, we are looking for consistencies with clustering that tell us if the patterns found in the data are sound.  The technique that I will use is the K-Means Clustering algorithm. K-Means is one of the most commonly used clustering algorithm and works well on large data sets compared to other clustering techniques. Because our data is quite large and analysis on it takes a long time to run on my local machine, I will use the optimized version of K-Means to speed up the computation time.For the K-Means Clustering model, since I am only interested in the text reviews, I will not include the other columns that describe the ratings, review date, beer category, beer availability, abv, and brewery type.It should also be noted that scaling is usually required with K-Means Clustering since it is a distance based model, but because I used a TF-IDF Vectorizer on the text reviews, all of the values are between 0 and 1, which means they are already on the same scale. Therefore, I will not be using a scaler on the dataset. ---  Data Preprocessing Before clustering can begin, I need to first import the tools that will be used.<jupyter_code># trinity 
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# scaler
from sklearn.preprocessing import MinMaxScaler

# k-means clustering
from sklearn.cluster import KMeans, MiniBatchKMeans

# silhouette score & silhouette samples
from sklearn.metrics import silhouette_score, silhouette_samples<jupyter_output><empty_output><jupyter_text>Now I need to read in the data. This will be the same file that was used in the `beer_classification` notebook and can be found here.<jupyter_code># read in data

sample_df = pd.read_csv('data/beer_sample_cleaned.csv')
# show first few rows

pd.options.display.max_columns=200
sample_df.head()
# get shape

sample_df.shape<jupyter_output><empty_output><jupyter_text>The data was read in properly and has 44,506 rows and 3,682 columns. Because I will be modifying the dataset to run the K-Means model, I will make a copy of the original dataset first.<jupyter_code># create copy
sample_text = sample_df.copy()

# check
pd.options.display.max_columns=200
sample_df.head()<jupyter_output><empty_output><jupyter_text>Because we are only looking at the vectorized review text, I need to drop the first 35 columns.<jupyter_code># keep only review text columns
sample_text = sample_text.iloc[: , 35:]

# check
sample_text.head()<jupyter_output><empty_output><jupyter_text>The dataset is now ready. It has 44,506 rows and 3,647 features. Although scaling is not required for this dataset and overfitting is not an issue with K-Means, I will use Principle Component Analysis (PCA) to reduce the dimensionality in order to speed up modelling, since K-Means is a computationally expensive technique.  Reduce Dimensionality with PCA PCA takes linear combinations of existing features and returns new features known as principal components (PC). In order to apply PCA, we need to set the `n_components` parameter. This is the number of PCs that we want to keep. The default number is equal to the number of features in the dataset, so setting a smaller `n_components` would result in a reduction in dimensionality. To decide which features are most useful, we look at the ones with the highest variance, because if a feature has the same value for every data point, we can't use it to identify different classes. For my model, I will apply a pre-determined threshold for the proportion of variance we want to keep after dimensionality has been reduced. In this case, I will use `n_components = 0.8`, which ensures that we keep the number of PCs to explain 80% of variance.<jupyter_code># import PCA
from sklearn.decomposition import PCA

# instantiate PCA with `n_components=0.8` (keep 80% variance)
my_pca = PCA(n_components=0.8)

# fit and transform 
sample_text_pca = my_pca.fit_transform(sample_text)
# number of optimal PCs

sample_text_pca.shape<jupyter_output><empty_output><jupyter_text>Using PCA, I have reduced the number of columns from 3,647 to 981 columns, which means that we need 918 PCs to get the 80% variance. I can now fit the K-Means model to the PCA transformed data.  K-Means Clustering Due to the large size of our dataset, I will use `MiniBatchKMeans`, an optimized version of K-Means, to speed up computing time. `MiniBatchKMeans` works by creating mini-batches of randomly sampled small subsets of the data. The mini-batches are then assigned to the nearest cluster center (i.e. centroid), which is updated by taking the average of the sample and all previous samples that were assigned to that centroid. This process is repeated until convergence.The parameter that I will consider for this model is the number of clusters `n_clusters`, which is set at a default of 8. In order to select the ideal number, I will iterate through different a range of cluster sizes to get the indertia and silhouette score for each size. Inertia is the Within Sum of Squares (WSS), which is the squared distances between each point and its assigned cluster center. For K-Means, we want to minimize the WSS, so the lower the inertia the better. I will also look at the silhouette score to help decide on the ideal number of clusters. It measures how similar the points are to their assigned clusters in relation to other clusters that they are not assigned to. Therefore, unlike inertia, we are looking for higher silhouette scores.<jupyter_code># create empty list of metrics
metrics = list()

# set range of Ks to search over
ks = np.arange(1, 21)

# iterating over ks
for k in ks:
    # instantiate
    kmm = MiniBatchKMeans(n_clusters=k, random_state=0, batch_size=1000)

    # fit on the PCA transformed data
    kmm.fit(sample_text_pca)
    
    # get cluster labels 
    kmm_labels = kmm.predict(sample_text_pca)
    
    # look at inertia
    inertia = kmm.inertia_
    
    # silhouette score
    if k > 1:
        silhouette = silhouette_score(sample_text_pca, kmm_labels)
    else:
        silhouette = np.NaN
        
    # List of dictionaries
    metrics.append({'inertia':inertia, 'silhouette':silhouette})
    
    print(f'{k} inertia:{kmm.inertia_}, silhouette:{silhouette}')<jupyter_output>1 inertia:9071.633619424741, silhouette:nan
2 inertia:8871.478275998386, silhouette:0.6703352137206224
3 inertia:8930.668567917508, silhouette:0.670678585397484
4 inertia:8870.637313537025, silhouette:0.6702349847195824
5 inertia:8930.528573914422, silhouette:0.6712018788536326
6 inertia:8928.16575030485, silhouette:0.6696227567405624
7 inertia:8929.378430754854, silhouette:0.6715893091650232
8 inertia:8922.010408904995, silhouette:0.6704734882408802
9 inertia:8930.428731876387, silhouette:0.6389885430727258
10 inertia:8858.675861373262, silhouette:0.6720936908223023
11 inertia:8929.160540913439, silhouette:0.46750427336068556
12 inertia:8928.428987358642, silhouette:0.6718459181327747
13 inertia:8856.282635451655, silhouette:0.6715836661848755
14 inertia:8856.037121693686, silhouette:0.6716194687391226
15 inertia:8858.04612818679, silhouette:-0.8038004388509319
16 inertia:8709.771782672866, silhouette:-0.11316024323576931
17 inertia:8916.464356947246, silhouette:0.6643064600534999
18 [...]<jupyter_text>The inertia and silhouette scores for each `n_cluster` have been identified and can be put into a dataframe and visualized for easier interpretation.<jupyter_code># put metrics in a dataframe
metrics_df = pd.DataFrame(metrics, index=ks)

# check
metrics_df
# visualize

plt.subplots(1, 2, figsize=(16, 8))
plt.subplot(1, 2, 1)
metrics_df['inertia'].plot(marker='o')
plt.title('Inertia Scree Plot for K Means', size = 16)
plt.xticks(ks)
plt.ylabel('Inertia', size = 14)
plt.xlabel('K - Number of Clusters', size = 14)
plt.subplot(1, 2, 2)
metrics_df['silhouette'].plot(marker='o')
plt.title('Silhouette Scree Plot for K Means', size = 16)
plt.ylabel('Mean Silhouette Score', size = 14)
plt.xlabel('K - Number of Clusters', size = 14)
plt.xticks(ks)
plt.show()<jupyter_output><empty_output><jupyter_text>Because we want to minimize inertia while maximizing the silhouette score, we should select the number of clusters at the "elbow" of the inertia score scree plot that also corresponds to a higher point in the silhouette scree plot. Based on the plots above, there isn't a very clear elbow in our inertia scree plot, but we see that the silhouette score is highest between 1 to 10 clusters, drops slightly at 12, then becomes erratic after `n_clusters=14`. Let's take a look at the cluster sizes to see if we can get a better idea of which `K` to select.<jupyter_code># range of Ks to search over
ks = np.arange(1, 15)

# iterating over ks
for k in ks:
    
    # instantiate
    kmm = KMeans(n_clusters=k)

    # fit on the whole dataset
    kmm.fit(sample_text_pca)
    
    # get cluster labels 
    kmm_labels = kmm.predict(sample_text_pca)
    
    # get cluster size
    size = np.unique(kmm_labels, return_counts=True)
    
    print(size)<jupyter_output>(array([0], dtype=int32), array([44506]))
(array([0, 1], dtype=int32), array([34655,  9851]))
(array([0, 1, 2], dtype=int32), array([34702,  2560,  7244]))
(array([0, 1, 2, 3], dtype=int32), array([34702,  7242,  2561,     1]))
(array([0, 1, 2, 3, 4], dtype=int32), array([ 4364,  2063, 33917,  3543,   619]))
(array([0, 1, 2, 3, 4, 5], dtype=int32), array([33922,  4500,  1948,   190,  3825,   121]))
(array([0, 1, 2, 3, 4, 5, 6], dtype=int32), array([ 3810, 33902,     2,  2023,   906,  3859,     4]))
(array([0, 1, 2, 3, 4, 5, 6, 7], dtype=int32), array([33906,  3190,  2119,  3217,     1,    60,    62,  1951]))
(array([0, 1, 2, 3, 4, 5, 6, 7, 8], dtype=int32), array([33902,  2317,     1,    14,    94,  1926,  3087,    52,  3113]))
(array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9], dtype=int32), array([33937,    52,    95,  4031,     1,    20,     1,     1,  2102,
        4266]))
(array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10], dtype=int32), array([33877,  2724,     1,  1705,     1,  1789,  2300, [...]<jupyter_text>We can see that the cluster sizes are quite uneven, with many including clusters of only 1 or 2. With a 3 clusters, the sizes are slightly better, and with `n_clusters=3`, inertia is relatively low and the silhouette score is high. Therefore, I will use `n_clusters=3` for my K-Means model, which will be fit to the PCA-transformed data. ---  Model v.01<jupyter_code># fit K Means (3 clusters)
# instantiate
kmm_batch = KMeans(n_clusters = 3)

# fit 
kmm_batch.fit(sample_text_pca)

# get cluster labels
kmm_labels_batch = kmm_batch.predict(sample_text_pca)<jupyter_output><empty_output><jupyter_text>I can check the cluster size for our K-Means model.<jupyter_code># check cluster size
np.unique(kmm_labels_batch, return_counts=True)<jupyter_output><empty_output><jupyter_text>We can see that the the three cluster sizes are 34,701 for cluster 1, 7,244 for cluster 2, and 2,561 for cluster 3. It should be noted that there isn't a hierarchy with the cluster labels and it doesn't make a difference which cluster the data points belong to, just that they are in the different clusters.  Clusters  I can now add the cluster labels back to our original dataset to see which cluster each feature belongs to. Because I am adding a column to the data, I will make a copy of the original dataframe first and then creating a new dataframe with the clusters as the categories and the original columns as the rows.<jupyter_code># create copy of original data
sample_text_copy = sample_text.copy()

# add `clusters` column
sample_text_copy['clusters'] = kmm_labels_batch

# create dataframe by grouping clusters and getting means of each column
text_clusters_df = sample_text_copy.groupby('clusters').mean()

# check
text_clusters_df.head()<jupyter_output><empty_output><jupyter_text>We can visualize the top 20 features for each cluster for easier interpretation.<jupyter_code># get top 20 features for cluster 1
cluster_1 = text_clusters_df.T.sort_values(0, ascending = False).head(20)

# visualize
plt.figure(figsize=(10, 8))
cluster_1[0].plot(kind='barh', color='orange')
plt.xlabel('mean', size=14)
plt.ylabel('word', size=14)
plt.xticks(size=12)
plt.yticks(size=12)
plt.title('Top 20 Words in Cluster 1', size=16)
plt.show()<jupyter_output><empty_output><jupyter_text>We can see that "beer" is the top word for cluster 1, which should probably be removed to get better insights. I will repeat the process for the remaining 2 clusters.<jupyter_code># get top 20 features for cluster 2
cluster_2 = text_clusters_df.T.sort_values(1, ascending = False).head(20)

# visualize
plt.figure(figsize=(10, 8))
cluster_2[1].plot(kind='barh', color='orange')
plt.xlabel('mean', size=14)
plt.ylabel('word', size=14)
plt.xticks(size=12)
plt.yticks(size=12)
plt.title('Top 20 Words in Cluster 2', size=16)
plt.show()<jupyter_output><empty_output><jupyter_text>"Beer" is also the top word for cluster 2. Let's look at the final cluster.<jupyter_code># get top 20 features for cluster 3
cluster_3 = text_clusters_df.T.sort_values(2, ascending = False).head(20)

# visualize
plt.figure(figsize=(10, 8))
cluster_3[2].plot(kind='barh', color='orange')
plt.xlabel('mean', size=14)
plt.ylabel('word', size=14)
plt.xticks(size=12)
plt.yticks(size=12)
plt.title('Top 20 Words in Cluster 3', size=16)
plt.show()<jupyter_output><empty_output><jupyter_text>Although "beer" isn't the top word for cluster 3, it still appears on the list. Let's remove "beer" from our dataset and repeat the whole process again. I'll make a copy of our original dataset. ---  Model v.02<jupyter_code># create copy of data
sample_text_2 = sample_text.copy()

# check
sample_text_2.head()<jupyter_output><empty_output><jupyter_text>Recall that our original dataset had 3,647 columns. I will now remove the `beer` column, so we should end up with 3,646 columns.<jupyter_code># remove `beer` column
sample_text_2 = sample_text_2.drop(columns = 'beer')

# check
sample_text_2.shape<jupyter_output><empty_output><jupyter_text>The `beer` column was succesfully dropped. I will do PCA with `n_components=0.8` on this dataset like before.<jupyter_code># instantiate PCA with `n_components=0.8` (keep 80% variance)
my_pca = PCA(n_components=0.8)

# fit and transform 
sample_text_2_pca = my_pca.fit_transform(sample_text_2)

# get number of optimal PCs
sample_text_2_pca.shape<jupyter_output><empty_output><jupyter_text>The number of PCs for the updated dataset is 925. I can fit a K-Means model with `n_clusters=3` again like in model v.01.<jupyter_code># fit K Means (3 clusters)
# instantiate
kmm_2 = KMeans(n_clusters=3)

# fit 
kmm_2.fit(sample_text_2_pca)

# get cluster labels
kmm_2_labels = kmm_2.predict(sample_text_2_pca)
# look at cluster sizing

np.unique(kmm_2_labels, return_counts=True)<jupyter_output><empty_output><jupyter_text>Our new cluster sizes are 35,955 for cluster 1, 6,985 for cluster 2, and 2,566 for cluster 3. Let's make a copy of `sample_text_2` and add the `cluster` label column.<jupyter_code># create copy of dataframe
sample_text_2_copy = sample_text_2.copy()

# add clusters column
sample_text_2_copy['clusters'] = kmm_2_labels<jupyter_output><empty_output><jupyter_text>I can now create a new dataframe with just the clusters and the original columns.<jupyter_code># create dataframe by grouping clusters and getting means of each column
text_clusters_2 = sample_text_2_copy.groupby('clusters').mean()

# check
text_clusters_2.head()<jupyter_output><empty_output><jupyter_text>Clusters  Like in v.01, we can extract the top 20 features for each cluster using some bar graphs.<jupyter_code># get top 20 features for cluster 1
v2_cluster_1 = text_clusters_2.T.sort_values(0, ascending = False).head(20)

# visualize
plt.figure(figsize=(10, 8))
v2_cluster_1[0].plot(kind='barh', color='orange')
plt.xlabel('mean', size=14)
plt.ylabel('word', size=14)
plt.xticks(size=12)
plt.yticks(size=12)
plt.title('Top 20 Words in Cluster 1', size=16)
plt.show()<jupyter_output><empty_output><jupyter_text>Looking at cluster 1, we see words that suggest a positive preference, such as "great", "good", "best", "better", "favorite", as well as words that describe how the beer is served, such as "bottle" and "tap", and the beer style, such as "ipa", and "ale." This suggests that IPAs are generally associated with positive sentiment, which supports our earlier EDA where we found IPAs ranking third in average `score` rating.Let's take a look at cluster 2 now.<jupyter_code># get top 20 features for cluster 2
v2_cluster_2 = text_clusters_2.T.sort_values(1, ascending = False).head(20)

# visualize
plt.figure(figsize=(10, 8))
v2_cluster_2[1].plot(kind='barh', color='orange')
plt.xlabel('mean', size=14)
plt.ylabel('word', size=14)
plt.xticks(size=12)
plt.yticks(size=12)
plt.title('Top 20 Words in Cluster 2', size=16)
plt.show()<jupyter_output><empty_output><jupyter_text>Looking at cluster 2, we see words that mostly describe the beer flavour, such as "hop", "taste", "light", "citrus", "orange", "malt", and "fruit." Looking at the beer categories, these words align most with Wild/Sour beers, that generally have the characteristics listed here. There are also some sentiment words, such as "nice", "good", and "like." Again, we saw earlier in EDA that the wild/sour beer category had the second highest average `score` rating.We can look at the final cluster now.<jupyter_code># get top 20 features for cluster 3
v2_cluster_3 = text_clusters_2.T.sort_values(2, ascending=False).head(20)

# visualize
plt.figure(figsize=(10, 8))
v2_cluster_3[2].plot(kind='barh', color='orange')
plt.xlabel('mean', size=14)
plt.ylabel('word', size=14)
plt.xticks(size=12)
plt.yticks(size=12)
plt.title('Top 20 Words in Cluster 3', size=16)
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09105708704145865, 'fraction_numerical': 0.0664749294551769, 'mean_word_length': 5.1800804828973845, 'pattern_counts': {'":': 0, '<': 92, '<?xml version=': 0, '>': 93, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '23', 'dir': 'jupyter-structured-clean-dedup', 'id': '582910', 'license': 'MIT', 'n_tokens_mistral': 6852, 'n_tokens_neox': 5876, 'n_words': 2244, 'path': '3 - beer_clustering.ipynb', 'repo_name': 'lynntlho/BeerCraft', 'size': '18429'}