Gaming Behavior Analysis

Predicting Player Engagement using Decision Tree and SVM

import pandas as pd import numpy as np import matplotlib.pyplot as plt import seaborn as sns # Preprocessing from sklearn.preprocessing import StandardScaler, LabelEncoder from sklearn.model_selection import train_test_split, GridSearchCV from sklearn.model_selection import RandomizedSearchCV # Algorithms (add another clustering thing) from sklearn.tree import DecisionTreeClassifier from sklearn.svm import SVC from sklearn.cluster import KMeans from sklearn.mixture import GaussianMixture from sklearn.cluster import SpectralClustering # Metrics from sklearn.metrics import accuracy_score, precision_score, recall_score, f1_score, roc_auc_score from sklearn.metrics import confusion_matrix, classification_report, roc_curve, silhouette_score from sklearn.metrics import confusion_matrix import seaborn as sns # Training vs Validation from sklearn.model_selection import train_test_split, cross_val_score from sklearn.metrics import accuracy_score, classification_report import warnings warnings.filterwarnings("ignore")

Popular seaborn heatmap colors:

Diverging (good for correlation): Sequential: Misc:

1. Data Loading and Exploration

# Load gaming dataset (reload to get original data) df = pd.read_csv('online_gaming_behavior_dataset - Copy.csv') print(f"Dataset shape: {df.shape}") df.head()
# Basic info df.info() df.isnull().sum() # Shows count of missing values per column # Statistical summary df.describe()
# Filtered Statistical summary (numerical features only) df.describe included irrelevant columns num_features = ['Age', 'PlayTimeHours', 'InGamePurchases', 'SessionsPerWeek', 'AvgSessionDurationMinutes', 'PlayerLevel', 'AchievementsUnlocked'] # Run describe on the array of numerical fields to provide an actually useful output... df[num_features].describe().round(2)
# Check for missing values df.isnull().sum() # Unique values count df.nunique()

Filtered statistical summary analysis:

Demographics: Gaming Behaviour: Sessions: Session Duration: Engagement Indicators: Business Insights:

Data Quality Assessment

# Drop PlayerID (reload data set at top prior to running this cell) df = df.drop('PlayerID', axis=1) # Define categorical and numerical columns manually as select_dtypes() can be unreliable cat_cols = ['Gender', 'Location', 'GameGenre', 'GameDifficulty', 'EngagementLevel'] num_cols = ['Age', 'PlayTimeHours', 'InGamePurchases', 'SessionsPerWeek', 'AvgSessionDurationMinutes', 'PlayerLevel', 'AchievementsUnlocked'] print(f'Numerical columns: {num_cols}') print(f'Categorical columns: {cat_cols}') # Check unique values for categorical columns for col in cat_cols: print(f'{col}: {df[col].unique()}')
# Exploratory analysis and Data Preprocessing # Univariate analysis of data columns df[num_cols].describe().T # create histograms df[num_cols].hist(bins=20, figsize=(18, 15), alpha=0.7, color='skyblue', edgecolor='black') plt.tight_layout() plt.show()

Exploratory Data Analysis

# Correlation matrix for numerical features plt.figure(figsize=(10, 8)) correlation_matrix = df[num_cols].corr() sns.heatmap(correlation_matrix, annot=True, cmap='seismic', center=0) plt.title('Feature Correlation Matrix') plt.show()
# Distribution of categorical variables plt.Colorizer('viridis') fig, axes = plt.subplots(2, 3, figsize=(15, 10)) axes = axes.ravel() for i, col in enumerate(cat_cols): df[col].value_counts().plot(kind='bar', ax=axes[i]) axes[i].set_title(f'{col} Distribution') axes[i].tick_params(axis='x', rotation=45) # Hide unused subplot (find better way to do this) axes[5].set_visible(False) for i, col in enumerate(cat_cols): counts = df[col].value_counts() bar_colors = plt.cm.viridis(np.linspace(0, 1, len(counts))) counts.plot(kind='bar', ax=axes[i], color=bar_colors) plt.tight_layout() plt.show() # Target variable analysis print('Engagement Level Distribution:') print(df['EngagementLevel'].value_counts()) print(f'\nPercentages:') print(df['EngagementLevel'].value_counts(normalize=True) * 100) plt.figure(figsize=(8, 5)) sns.countplot(data=df, x='EngagementLevel', palette='viridis') plt.title('Player Engagement Distribution') plt.show() # Key insights by engagement level # Further break metrics down to mitigate poor formatting key_metrics = ['Age', 'PlayTimeHours', 'SessionsPerWeek', 'AvgSessionDurationMinutes', 'PlayerLevel'] engagement_stats = df.groupby('EngagementLevel')[key_metrics].mean().round(2) print(engagement_stats.to_string(float_format='%.2f'))
# Non-linear Relationship Exploration # Convert engagement to numeric from string df['Engagement_Score'] = df['EngagementLevel'].map({'Low': 0, 'Medium': 1, 'High': 2}) # Sample to reduce noise df_sample = df.sample(n=5000, random_state=42) # Scatter plots with engagement colors fig, axes = plt.subplots(2, 2, figsize=(12, 10)) features = ['Age', 'PlayTimeHours', 'SessionsPerWeek', 'PlayerLevel'] for i, feature in enumerate(features): ax = axes[i//2, i%2] sns.scatterplot(data=df_sample, x=feature, y='AvgSessionDurationMinutes', hue='EngagementLevel', alpha=0.7, ax=ax) ax.set_title(f'{feature} vs Session Duration by Engagement') plt.tight_layout() plt.show() # Feature interactions df['PlayTime_Sessions'] = df['PlayTimeHours'] * df['SessionsPerWeek'] df['Age_Level'] = df['Age'] * df['PlayerLevel'] # Age binned analysis df['Age_Bins'] = pd.cut(df['Age'], bins=5, labels=['Very Young', 'Young', 'Mid', 'Mature', 'Senior']) engagement_by_age = df.groupby(['Age_Bins', 'EngagementLevel']).size().unstack(fill_value=0) plt.figure(figsize=(10, 6)) engagement_by_age.plot(kind='bar', stacked=True) plt.title('Engagement Patterns Across Age Groups') plt.xticks(rotation=45) plt.show() # Correlation analysis print("Non-linear patterns found:") print(f"PlayTime alone: {df['PlayTimeHours'].corr(df['Engagement_Score']):.3f}") print(f"Sessions alone: {df['SessionsPerWeek'].corr(df['Engagement_Score']):.3f}") print(f"PlayTime*Sessions: {df['PlayTime_Sessions'].corr(df['Engagement_Score']):.3f}") print(f"Age*Level: {df['Age_Level'].corr(df['Engagement_Score']):.3f}")
# Multiple Non-linear Relationship Candidates df['Engagement_Score'] = df['EngagementLevel'].map({'Low': 0, 'Medium': 1, 'High': 2}) df_sample = df.sample(n=3000, random_state=42) fig, axes = plt.subplots(2, 2, figsize=(15, 10)) # 1. PlayTime (diminishing returns pattern) sns.boxplot(data=df, x='EngagementLevel', y='PlayTimeHours', ax=axes[0,0]) axes[0,0].set_title('PlayTime vs Engagement (Potential Threshold Effect)') # 2. PlayerLevel (S-curve pattern) sns.boxplot(data=df, x='EngagementLevel', y='PlayerLevel', ax=axes[0,1]) axes[0,1].set_title('Player Level vs Engagement (Progression Curve)') # 3. Session Duration (Sweet spot pattern) sns.boxplot(data=df, x='EngagementLevel', y='AvgSessionDurationMinutes', ax=axes[1,0]) axes[1,0].set_title('Session Duration vs Engagement (Optimal Range)') # 4. Age (U-shaped or inverted-U pattern) sns.boxplot(data=df, x='EngagementLevel', y='Age', ax=axes[1,1]) axes[1,1].set_title('Age vs Engagement (Life Stage Effect)') plt.tight_layout() plt.show() # Quantify non-linearity with polynomial features from sklearn.preprocessing import PolynomialFeatures from sklearn.linear_model import LinearRegression df['Engagement_Score'] = df['EngagementLevel'].map({'Low': 0, 'Medium': 1, 'High': 2}) features = ['PlayTimeHours', 'PlayerLevel', 'AvgSessionDurationMinutes', 'Age'] print("Non-linearity Evidence (R² improvement with polynomial features):") for feature in features: # Linear fit X_linear = df[[feature]].values linear_r2 = LinearRegression().fit(X_linear, df['Engagement_Score']).score(X_linear, df['Engagement_Score']) # Polynomial fit poly = PolynomialFeatures(degree=2) X_poly = poly.fit_transform(X_linear) poly_r2 = LinearRegression().fit(X_poly, df['Engagement_Score']).score(X_poly, df['Engagement_Score']) improvement = poly_r2 - linear_r2 print(f"{feature}: Linear R²={linear_r2:.3f}, Poly R²={poly_r2:.3f}, Improvement={improvement:.3f}")
df['Engagement_Score'] = df['EngagementLevel'].map({'Low': 0, 'Medium': 1, 'High': 2}) correlation_with_target = df[num_cols + ['Engagement_Score']].corr() plt.figure(figsize=(16, 6)) # Linear correlations plt.subplot(1, 2, 1) sns.heatmap(correlation_with_target, annot=True, cmap='seismic', center=0, fmt='.2f', square=True, cbar_kws={'shrink': 0.8}) plt.title('Linear Correlations\n(Decision Tree View)', fontsize=12, pad=20) plt.xticks(rotation=45, ha='right') plt.yticks(rotation=0) # Non-linear evidence plt.subplot(1, 2, 2) nonlinear_data = { 'Sessions': [0.606], 'PlayTime': [-0.002], 'PlayTime×Sessions': [0.406], 'Age×Level': [0.051] } nonlinear_df = pd.DataFrame(nonlinear_data, index=['Engagement']) sns.heatmap(nonlinear_df, annot=True, cmap='RdYlBu', center=0.3, fmt='.3f', square=True, cbar_kws={'shrink': 0.8, 'label': 'Correlation'}) plt.title('Feature Interactions\n(SVM View)', fontsize=12, pad=20) plt.xticks(rotation=45, ha='right') plt.yticks(rotation=0) plt.tight_layout(pad=3.0) plt.show()

Key findings

The feature 'Sessions per week' is king and the strongest single predictor of engagement. Playtime alone is practically useless however, playtime combined with sessions per week provides a value of 0.406 correlation proving that non-linear relationships are present and possible in this dataset explaining why decision tree model is unable to perform as well as SVM

(Decision Tree View)

This is a correlation heatmap between various features:

Implication: Engagement is primarily driven by frequency of sessions and average session duration, not by age or total playtime.

(SVM View)

Shows interaction strength between features and engagement:

Implication: Engagement is highly influenced by session frequency, and there is a meaningful interaction between PlayTime and Sessions (players who play more frequently and longer in each session tend to be more engaged).

Extrapolation & Insights

2. Data Preprocessing

df_original = pd.read_csv('online_gaming_behavior_dataset - Copy.csv') # Allow analysis on the dataset without contaminating the original df_analysis = df_original.copy() # For preprocessing df_model = df_original.copy() df_model = df_model.drop('PlayerID', axis=1) # Handle categorical variables categorical_cols = df_model.select_dtypes(include=['object']).columns.tolist() categorical_cols.remove('EngagementLevel') for col in categorical_cols: le = LabelEncoder() df_model[col] = le.fit_transform(df_model[col]) le = LabelEncoder() df_model['EngagementLevel'] = le.fit_transform(df_model['EngagementLevel']) X = df_model.drop('EngagementLevel', axis=1) y = df_model['EngagementLevel'] scaler = StandardScaler() X_scaled = scaler.fit_transform(X) X_train, X_test, y_train, y_test = train_test_split(X_scaled, y, test_size=0.2, random_state=42, stratify=y)

Evaluate model training and validation

# Split data into training and validation sets def split_data(X, y, test_size=0.2, random_state=42): return train_test_split(X, y, test_size=test_size, random_state=random_state, stratify=y) # Train and validate model def train_validate_model(model, X_train, X_val, y_train, y_val): # Train model model.fit(X_train, y_train) # Validate model y_pred = model.predict(X_val) accuracy = accuracy_score(y_val, y_pred) print(f"Validation Accuracy: {accuracy:.4f}") print("\nClassification Report:") print(classification_report(y_val, y_pred)) return model, accuracy # Cross-validation for robust evaluation def cross_validate_model(model, X, y, cv=5): scores = cross_val_score(model, X, y, cv=cv, scoring='accuracy') print(f"Cross-validation scores: {scores}") print(f"Mean CV accuracy: {scores.mean():.4f} (+/- {scores.std() * 2:.4f})") return scores
# Apply training/validation print("=== Training/Validation Results ===") # Decision Tree with validation print("\n--- Decision Tree ---") dt_model, dt_val_acc = train_validate_model(DecisionTreeClassifier(random_state=42), X_train, X_test, y_train, y_test) # SVM with validation print("\n--- SVM ---") svm_model, svm_val_acc = train_validate_model(SVC(random_state=42), X_train, X_test, y_train, y_test) # Cross-validation comparison print("\n--- Cross-Validation Comparison ---") print("Decision Tree CV:") cross_validate_model(DecisionTreeClassifier(random_state=42), X_scaled, y) print("\nSVM CV:") cross_validate_model(SVC(random_state=42), X_scaled, y)

How the process works:

The data gets 'split'

40,034 total players

Training phase

Validation phase

Cross fold validation applied:

Split data into 5 equal parts:

This is important as:

3. Decision Tree Model

# Basic Decision Tree dt = DecisionTreeClassifier(random_state=42) dt.fit(X_train, y_train) dt_pred = dt.predict(X_test) # Evaluate print("Decision Tree Results:") print(f"Accuracy: {accuracy_score(y_test, dt_pred):.4f}") print(f"Precision: {precision_score(y_test, dt_pred, average='weighted'):.4f}") print(f"Recall: {recall_score(y_test, dt_pred, average='weighted'):.4f}") print(f"F1-Score: {f1_score(y_test, dt_pred, average='weighted'):.4f}")
# Hyperparameter tuning for Decision Tree dt_params = { 'max_depth': [3, 5, 10, None], 'min_samples_split': [2, 5, 10], 'min_samples_leaf': [1, 2, 4] } dt_grid = GridSearchCV(DecisionTreeClassifier(random_state=42), dt_params, cv=5, scoring='accuracy', n_jobs=-1) dt_grid.fit(X_train, y_train) print(f"Best DT params: {dt_grid.best_params_}") print(f"Best DT score: {dt_grid.best_score_:.4f}")
# Hypertuned tree # Pass in params max_depth=10, {min_samples_split=4, min_samples_leaf=10} dt = DecisionTreeClassifier(max_depth=10, min_samples_split=4, min_samples_leaf=10, random_state=42) dt.fit(X_train, y_train) dt_pred = dt.predict(X_test) # Evaluate print("Decision Tree Results:") print(f"Accuracy: {accuracy_score(y_test, dt_pred):.4f}") print(f"Precision: {precision_score(y_test, dt_pred, average='weighted'):.4f}") print(f"Recall: {recall_score(y_test, dt_pred, average='weighted'):.4f}") print(f"F1-Score: {f1_score(y_test, dt_pred, average='weighted'):.4f}") cm = confusion_matrix(y_test, dt_pred) # Visualize plt.figure(figsize=(8, 6)) sns.heatmap(cm, annot=True, fmt='d', cmap='Blues', xticklabels=['Low', 'Medium', 'High'], yticklabels=['Low', 'Medium', 'High']) plt.xlabel('Predicted') plt.ylabel('Actual') plt.title('Confusion Matrix') plt.show()

4. SVM Model

# Basic 'sampled' svm for benchmark print("Training SVM on sample for faster execution...") X_sample, _, y_sample, _ = train_test_split(X_train, y_train, train_size=0.3, random_state=42, stratify=y_train) svm = SVC(random_state=42, probability=True) svm.fit(X_sample, y_sample) # Train on 30% sample faster with minimal performance loss svm_pred = svm.predict(X_test) # Test on full test set print("SVM Results:") print(f"Accuracy: {accuracy_score(y_test, svm_pred):.4f}") print(f"Precision: {precision_score(y_test, svm_pred, average='weighted'):.4f}") print(f"Recall: {recall_score(y_test, svm_pred, average='weighted'):.4f}") print(f"F1-Score: {f1_score(y_test, svm_pred, average='weighted'):.4f}") print(f"ROC AUC: {roc_auc_score(y_test, svm.predict_proba(X_test), multi_class='ovr'):.4f}") print(f"Training sample size: {len(X_sample)} (30% of training data)")
# SVM produced from GridSearch: # Params: {'C': 1, 'gamma': 'scale', 'kernel': 'rbf'} # Try adjusting to C = 10? svm = SVC(C=10, gamma='scale', kernel='rbf', random_state=42, probability=True) svm.fit(X_train, y_train) svm_pred = svm.predict(X_test) print("SVM Results:") print(f"Accuracy: {accuracy_score(y_test, svm_pred):.4f}") print(f"Precision: {precision_score(y_test, svm_pred, average='weighted'):.4f}") print(f"Recall: {recall_score(y_test, svm_pred, average='weighted'):.4f}") print(f"F1-Score: {f1_score(y_test, svm_pred, average='weighted'):.4f}") print(f"ROC AUC: {roc_auc_score(y_test, svm.predict_proba(X_test), multi_class='ovr'):.4f}") # First run with out-of-box params (GridSearch BEST) # Params for these results = {'C': 1, 'gamma': 'scale', 'kernel': 'rbf'}
# SVM produced from RandomSearch: # Params: {'C': 1, 'gamma': 'scale', 'kernel': 'rbf'} svm = SVC(C=10, gamma='auto', kernel='rbf', random_state=42, probability=True) svm.fit(X_train, y_train) svm_pred = svm.predict(X_test) print("SVM Results:") print(f"Accuracy: {accuracy_score(y_test, svm_pred):.4f}") print(f"Precision: {precision_score(y_test, svm_pred, average='weighted'):.4f}") print(f"Recall: {recall_score(y_test, svm_pred, average='weighted'):.4f}") print(f"F1-Score: {f1_score(y_test, svm_pred, average='weighted'):.4f}") print(f"ROC AUC: {roc_auc_score(y_test, svm.predict_proba(X_test), multi_class='ovr'):.4f}") # First run with out-of-box params RandomSearch BEST # Params {'C': 1, 'gamma': 'scale', 'kernel': 'rbf'}
# Hyperparameter tuning for SVM (calibrate runtime as it's running GridSearch on 40k plus samples) # SVM scales poorly with large datasets # An RBF Kernel is quite computationally expensive # added parallel processing via 'n_jobs=-1' (uses all available cpu cores) # This now runs 96.1% quicker than previous iterations! Down from 36 minutes to 1 minute 24.5 seconds # Use subset for hyperparameter tuning (0.2) which means sample size is now only 20% # Try with x% of dataset find best compromise from that X_sample, _, y_sample, _ = train_test_split(X_train, y_train, train_size=0.2, random_state=42, stratify=y_train) svm_params = { 'C': [0.1, 1, 10], 'kernel': ['rbf', 'linear'], 'gamma': ['scale', 'auto'] } # CV folds: Take X amount of training data and splits it into equal parts, i.e. for 5 folds it will train on 4 folds and test on 1 and repeats this 5 times # Each fold gets to be used as test data exactly once # Averages the X amount of scores in this case 5 to produce the final evaluation metric # Adjust levels of cv folds to find acceptable level (3-5 produce similar results, 3 is a good compromise) svm_grid = GridSearchCV(SVC(random_state=42, probability=True), svm_params, cv=3, scoring='accuracy', n_jobs=-1) svm_grid.fit(X_sample, y_sample) print(f"Best SVM params: {svm_grid.best_params_}") print(f"Best SVM score: {svm_grid.best_score_:.4f}")
# Alternative hyperparam tuning method via random search: # Runs 50% quicker on average X_sample, _, y_sample, _ = train_test_split(X_train, y_train, train_size=0.2, random_state=42, stratify=y_train) svm_params = { 'C': [0.1, 1, 10], 'kernel': ['rbf', 'linear'], 'gamma': ['scale', 'auto'] } # Test only 6 random combinations instead of all 12 svm_random = RandomizedSearchCV(SVC(random_state=42, probability=True), svm_params, n_iter=6, cv=3, scoring='accuracy', n_jobs=-1, random_state=42) svm_random.fit(X_sample, y_sample) print(f"Best SVM params: {svm_random.best_params_}") print(f"Best SVM score: {svm_random.best_score_:.4f}") # Comparing model comparison methods (GridSearch - Tests all combinations) and RandomizedSearch(tests only a random subset of the dataset)

SVM Results:

SVM scores and what they mean

Accuracy

Precision

Recall

F1-Score

ROC AUC

SVM Model (89.8% accuracy):

5. Model Comparison

# Compare best models best_dt = dt_grid.best_estimator_ best_svm = svm_grid.best_estimator_ best_svmr = svm_random.best_estimator_ dt_final_pred = best_dt.predict(X_test) svm_final_pred = best_svm.predict(X_test) svmr_final_pred = best_svmr.predict(X_test)
# Results comparison results = pd.DataFrame({ 'Model': ['Decision Tree', 'SVM'], 'Accuracy': [accuracy_score(y_test, dt_final_pred), accuracy_score(y_test, svm_final_pred)], 'Precision': [precision_score(y_test, dt_final_pred, average='weighted'), precision_score(y_test, svm_final_pred, average='weighted')], 'Recall': [recall_score(y_test, dt_final_pred, average='weighted'), recall_score(y_test, svm_final_pred, average='weighted')], 'F1-Score': [f1_score(y_test, dt_final_pred, average='weighted'), f1_score(y_test, svm_final_pred, average='weighted')] }) print(results)
# Compare all model variants models = { 'Decision Tree (Default)': dt, 'Decision Tree (Tuned)': best_dt, 'SVM (Default)': svm, 'SVM (GridSearch)': best_svm, 'SVM (RandomSearch)': best_svmr } # Get predictions for all models predictions = {} for name, model in models.items(): predictions[name] = model.predict(X_test) # Create comprehensive results DataFrame results_data = [] for name, pred in predictions.items(): results_data.append({ 'Model': name, 'Accuracy': accuracy_score(y_test, pred), 'Precision': precision_score(y_test, pred, average='weighted'), 'Recall': recall_score(y_test, pred, average='weighted'), 'F1-Score': f1_score(y_test, pred, average='weighted'), 'ROC_AUC': roc_auc_score(y_test, models[name].predict_proba(X_test), multi_class='ovr') if hasattr(models[name], 'predict_proba') else 'N/A' }) results_df = pd.DataFrame(results_data) print("=== Complete Model Performance Comparison ===") print(results_df.round(4)) # Hyperparameter comparison print("\n=== Hyperparameter Comparison ===") hyperparams = pd.DataFrame({ 'Model': ['DT Default', 'DT Tuned', 'SVM Default', 'SVM GridSearch', 'SVM RandomSearch'], 'Key Parameters': [ 'max_depth=None, min_samples_split=2', f'max_depth={best_dt.max_depth}, min_samples_split={best_dt.min_samples_split}', 'C=1.0, kernel=rbf, gamma=scale', f'C={best_svm.C}, kernel={best_svm.kernel}, gamma={best_svm.gamma}', f'C={best_svmr.C}, kernel={best_svmr.kernel}, gamma={best_svmr.gamma}' ], 'Tuning Method': ['None', 'GridSearchCV', 'None', 'GridSearchCV', 'RandomizedSearchCV'] }) print(hyperparams) # Performance improvement analysis print("\n=== Performance Improvements ===") dt_improvement = results_df[results_df['Model'] == 'Decision Tree (Tuned)']['Accuracy'].iloc[0] - \ results_df[results_df['Model'] == 'Decision Tree (Default)']['Accuracy'].iloc[0] svm_grid_improvement = results_df[results_df['Model'] == 'SVM (GridSearch)']['Accuracy'].iloc[0] - \ results_df[results_df['Model'] == 'SVM (Default)']['Accuracy'].iloc[0] svm_random_improvement = results_df[results_df['Model'] == 'SVM (RandomSearch)']['Accuracy'].iloc[0] - \ results_df[results_df['Model'] == 'SVM (Default)']['Accuracy'].iloc[0] print(f"Decision Tree tuning improvement: {dt_improvement:.4f}") print(f"SVM GridSearch improvement: {svm_grid_improvement:.4f}") print(f"SVM RandomSearch improvement: {svm_random_improvement:.4f}")

Analysis on results

Overall based on this dataset SVM is a better model

Algorithm differences

Decision tree: SVM: SVM boasts superior pattern recognition

Decision tree's weaknesses

6. Clustering Analysis

inertias = [] kmeans_silhouette = [] gmm_silhouette = [] spc_silhouette = [] k_range = range(2, 11) # KMeans evaluation (Reference point) for k in k_range: kmeans = KMeans(n_clusters=k, random_state=42) kmeans.fit(X_scaled) inertias.append(kmeans.inertia_) kmeans_silhouette.append(silhouette_score(X_scaled, kmeans.labels_)) # GaussianMixture evaluation for k in k_range: gmm = GaussianMixture(n_components=k, random_state=42) gmm.fit(X_scaled) gmm_silhouette.append(silhouette_score(X_scaled, gmm.predict(X_scaled))) # Plot plt.figure(figsize=(15, 5)) plt.subplot(1, 3, 1) plt.plot(k_range, inertias, 'bo-') plt.title('KMeans Elbow Method') plt.xlabel('Number of Clusters') plt.ylabel('Inertia') # Use silhouette score instead more clear plt.subplot(1, 3, 2) plt.plot(k_range, kmeans_silhouette, 'ro-') plt.title('KMeans Silhouette Score') plt.xlabel('Number of Clusters') plt.ylabel('Silhouette Score') plt.subplot(1, 3, 3) plt.plot(k_range, gmm_silhouette, 'go-') plt.title('GMM Silhouette Score') plt.xlabel('Number of Clusters') plt.ylabel('Silhouette Score') plt.tight_layout() plt.show() # Kmeans optima optimal_k = k_range[np.argmax(kmeans_silhouette)] print(f"Optimal number of clusters (kmeans): {optimal_k}") # GMM optima optimal_k = k_range[np.argmax(gmm_silhouette)] print(f"Optimal number of clusters (GMM): {optimal_k}") # Use optimal K value as found via silhouette score as its clearer than elbow final_kmeans = KMeans(n_clusters=optimal_k, random_state=42) kmeansclusters = final_kmeans.fit_predict(X_scaled) final_Gmm = GaussianMixture(n_components=optimal_k, random_state=42) gmmclusters = final_Gmm.fit_predict(X_scaled) # Kmeans final print(f"Silhouette Score (Kmeans): {silhouette_score(X_scaled, kmeansclusters):.4f}") print(f"Inertia: {final_kmeans.inertia_:.2f}") # GMM final print(f"Silhouette Score (gmm): {silhouette_score(X_scaled, gmmclusters):.4f}") # Add clusters to original dataframe df['Cluster'] = gmmclusters print("\nCluster distribution:") print(df['Cluster'].value_counts().sort_index())

7. Business Insights and Recommendations

# Analyze distribution of engagement levels across clusters cluster_engagement_analysis = df.groupby(['Cluster', 'EngagementLevel']).size().unstack(fill_value=0) print("Distribution of Engagement Levels by Cluster:") print(cluster_engagement_analysis) # Visualize cluster engagement distribution with heatmap plt.figure(figsize=(12, 6)) sns.heatmap(cluster_engagement_analysis, annot=True, fmt='d', cmap='YlOrRd', cbar_kws={'label': 'Number of Players'}) plt.title('Distribution of Player Engagement Levels Across Clusters') plt.xlabel('Engagement Level') plt.ylabel('Cluster') plt.tight_layout() plt.show()

Clusters:

How to utilise these clusters to target niche groups:

It is now possible to split these groups and design experiences tailor made for specific audiences, i.e. cluster 3 are far more inclined to purchase additional content or premium content such as cosmetics and all efforts should be made to market these items heavily toward them

Cluster one being more so for weekend players due to intense play style but moderate session count with long durations would likely be more inclined to enjoy time limited events

Cluster zero would be more so favoured toward accessible smaller purchases due to less commitment

High risk clusters

Clusters zero and one are the highest risk clusters for churn.

Cluster zero is the highest risk for players disengaging with the game:

Risk factors:

Cluster one: Shows far more inconsistent play patterns

Based on this analysis these clusters can be used as early warning for player retention

Based on which cluster a player falls into appropriate business action could be taken.

Based on the results from the correlation and subsequent correlation analysis it can be said that session frequency (not total playtime) is a better predictor of engagement, as such clusters showing a drastic decrease in this metric will be the highest risk of churn.