Машинное обучение с Python в розничной торговле : кластеризация магазинов по ассортиментной матрице продаж

Рассмотрим задачу кластеризации розничных магазинов одежды по ассортиментной матрице продаж. В качестве исходных данных возьмем ассортиментные матрицы продаж 72 магазинов одежды в осенне-зимний сезон  (массив данных сгенерирован с помощью sklearn.datasets.samples_generator). Доли в продажах разделены на 11 комбинаций по полу-товарной группе-сезон. В реальности различия в ассортиментных матрицах могут быть связаны с разными климатическими условиями (магазины могут находятся в разных городах), а также с разницей в расположении магазинов в городе и в торговом центре. В результате анализа необходимо разделить магазины на кластеры, чтобы в дальнейшем более правильно организовать их снабжение товарами.

Загружаем данные из файла Excel и проводим предварительный анализ :

df = pd.read_excel('am_shop.xlsx',index_col='shop')

# Размер данных

df.shape

(72, 11)

#Тип данных
df.dtypes

female pants universal       float64
female jacket demi-season    float64
female jacket winter         float64
female sweater winter        float64
female sweater universal     float64
male pants universal         float64
male jacket demi-season      float64
male jacket winter           float64
male sweater winter          float64
male sweater universal       float64
male hoody universal         float64
dtype: object

Наш датафрейм содержит 72 строки, выведем первые и последние пять строк

pd.concat([df.head().T, df.tail().T],axis=1)


sh1
sh2
sh3
sh4
sh5
sh68
sh69
sh70
sh71
sh72
female pants universal
0.1351
0.1609
0.1626
0.2029
0.1581
0.0554
0.0816
0.0769
0.0657
0.0567
female jacket demi-season
0.1499
0.0684
0.0534
0.0412
0.0562
0.0306
0.0441
0.0388
0.0316
0.0340
female jacket winter
0.0894
0.1188
0.1299
0.1622
0.1911
0.2684
0.2407
0.2631
0.2574
0.2225
female sweater winter
0.0329
0.0288
0.0430
0.0355
0.0425
0.0132
0.0196
0.0296
0.0263
0.0197
female sweater universal
0.0512
0.0453
0.0489
0.0564
0.0431
0.0095
0.0182
0.0227
0.0176
0.0134
male pants universal
0.1775
0.2136
0.2271
0.1915
0.1778
0.1261
0.1033
0.1049
0.1154
0.1381
male jacket demi-season
0.1441
0.0850
0.0831
0.0678
0.0654
0.0740
0.0588
0.0468
0.0427
0.0524
male jacket winter
0.0715
0.1038
0.1091
0.1116
0.1136
0.3182
0.3242
0.3237
0.3379
0.3572
male sweater winter
0.0309
0.0348
0.0421
0.0335
0.0388
0.0201
0.0260
0.0264
0.0329
0.0293
male sweater universal
0.0335
0.0505
0.0330
0.0287
0.0369
0.0306
0.0298
0.0290
0.0201
0.0304
male hoody universal
0.0841
0.0900
0.0677
0.0685
0.0765
0.0540
0.0535
0.0381
0.0524
0.0464

Сделаем копию введенных данных и наименование колонок

df_orig = df.copy()

colname=df_orig.columns.values.tolist()

С помощью корреляционной матрицы оценим взаимосвязь данных между собой, заменим единичные значения по диагонали нулями, чтобы в дальнейшем иметь возможность выбирать пары с наибольшим абсолютным значением корреляции

corr_mat = df.corr()

for x in range(corr_mat.shape[0]):

    corr_mat.iloc[x,x] = 0.0

    

corr_mat

female pants universal	female jacket demi-season	female jacket winter	female sweater winter	female sweater universal	male pants universal	male jacket demi-season	male jacket winter	male sweater winter	male sweater universal	male hoody universal
female pants universal	0.0000	0.2054	-0.6937	0.4556	0.6439	0.6669	0.1764	-0.8025	0.1905	0.3398	0.4886
female jacket demi-season	0.2054	0.0000	-0.2618	0.0197	0.4025	-0.0166	0.7095	-0.5104	-0.1300	0.1757	0.2676
female jacket winter	-0.6937	-0.2618	0.0000	-0.3881	-0.5272	-0.8632	-0.4227	0.7438	-0.3969	-0.5925	-0.7436
female sweater winter	0.4556	0.0197	-0.3881	0.0000	0.5526	0.2676	-0.1644	-0.4614	0.6918	0.4062	0.4126
female sweater universal	0.6439	0.4025	-0.5272	0.5526	0.0000	0.2919	0.3717	-0.7324	0.2527	0.4531	0.5513
male pants universal	0.6669	-0.0166	-0.8632	0.2676	0.2919	0.0000	0.1946	-0.6835	0.3387	0.4929	0.5926
male jacket demi-season	0.1764	0.7095	-0.4227	-0.1644	0.3717	0.1946	0.0000	-0.5417	-0.1481	0.2164	0.4420
male jacket winter	-0.8025	-0.5104	0.7438	-0.4614	-0.7324	-0.6835	-0.5417	0.0000	-0.3110	-0.5596	-0.7395
male sweater winter	0.1905	-0.1300	-0.3969	0.6918	0.2527	0.3387	-0.1481	-0.3110	0.0000	0.3619	0.4562
male sweater universal	0.3398	0.1757	-0.5925	0.4062	0.4531	0.4929	0.2164	-0.5596	0.3619	0.0000	0.5978
male hoody universal	0.4886	0.2676	-0.7436	0.4126	0.5513	0.5926	0.4420	-0.7395	0.4562	0.5978	0.0000

#Выведем пары с максимальными корреляциями c их значениями
pc=corr_mat.abs().idxmax()
cormax=corr_mat.loc[pc.index,pc.values]
df_cor=pd.DataFrame({'colname ': pc,'cor': np.diagonal(cormax)})
df_cor.sort_values(by=['cor'])

colname
cor
female jacket winter
male pants universal
-0.8632
male pants universal
female jacket winter
-0.8632
female pants universal
male jacket winter
-0.8025
male jacket winter
female pants universal
-0.8025
male hoody universal
female jacket winter
-0.7436
female sweater universal
male jacket winter
-0.7324
male sweater universal
male hoody universal
0.5978
female sweater winter
male sweater winter
0.6918
male sweater winter
female sweater winter
0.6918
female jacket demi-season
male jacket demi-season
0.7095
male jacket demi-season
female jacket demi-season
0.7095

Проведем некоторую предобработку данных: приведем данные к одному масштабу 

from sklearn.preprocessing import MinMaxScaler

scaler = MinMaxScaler()

df = scaler.fit_transform(df)

df=pd.DataFrame(df,columns=colname)

Визуализируем зависимости в данных

sns.set_context('notebook')

sns.set_style('white')

sns.pairplot(df)

Далее используем метод главных компонент для уменьшения размерности наших данных.

from sklearn.decomposition import PCA

pca_list = list()

feature_weight_list = list()

# Цикл по количеству компонент

for n in range(1, 12):

    

    # Создаем модель

    PCAmod = PCA(n_components=n)

    PCAmod.fit(df)

    

    # Сохраняем модель и дисперсию

    pca_list.append(pd.Series({'n':n, 'model':PCAmod,

                               'var': PCAmod.explained_variance_ratio_.sum()}))

    

    # Расчет и сохранение важности переменных

    abs_feature_values = np.abs(PCAmod.components_).sum(axis=0)

    feature_weight_list.append(pd.DataFrame({'n':n, 

                                             'features': df.columns,

                                             'values':abs_feature_values/abs_feature_values.sum()}))

    

#Выводим суммарную долю объясненной дисперсии по компонентам

pca_df = pd.concat(pca_list, axis=1).T.set_index('n')

pca_df

model	var
n
1	PCA(n_components=1)	0.5530
2	PCA(n_components=2)	0.6870
3	PCA(n_components=3)	0.7863
4	PCA(n_components=4)	0.8653
5	PCA(n_components=5)	0.9171
6	PCA(n_components=6)	0.9454
7	PCA(n_components=7)	0.9656
8	PCA(n_components=8)	0.9808
9	PCA(n_components=9)	0.9920
10	PCA(n_components=10)	0.9999
11	PCA(n_components=11)	1.0000

Диаграмма долей объясненной дисперсии по компонентам

plt.bar(range(1,12),PCAmod.explained_variance_ratio_)

plt.title('Percentage of Variance Explained')

plt.xlabel('Number of Components')

plt.ylabel('Percentage of Variance Explained'

Диаграмма суммарной доли объясненной дисперсии по компонентам

sns.set_context('talk')

ax = pca_df['var'].plot(kind='bar')

ax.set(xlabel='Number of Components',

       ylabel='Percent explained variance',

       title='Explained Variance vs Components')

Выведем веса показателей  по компонентам

features_df = (pd.concat(feature_weight_list)

               .pivot(index='n', columns='features', values='values'))

features_df.columns=colname

features_df

female pants universal	female jacket demi-season	female jacket winter	female sweater winter	female sweater universal	male pants universal	male jacket demi-season	male jacket winter	male sweater winter	male sweater universal	male hoody universal
n
1	0.0978	0.0308	0.1314	0.0696	0.0897	0.1196	0.0466	0.1361	0.0586	0.1024	0.1174
2	0.0623	0.0837	0.0806	0.1217	0.0651	0.0770	0.1195	0.1131	0.1285	0.0803	0.0683
3	0.0501	0.0928	0.0901	0.1305	0.1017	0.1295	0.0954	0.0927	0.1036	0.0619	0.0516
4	0.0899	0.0755	0.0707	0.1161	0.0981	0.1075	0.0928	0.0858	0.0916	0.1038	0.0680
5	0.0896	0.0694	0.0669	0.0960	0.0940	0.0878	0.0975	0.0709	0.1109	0.1349	0.0822
6	0.0783	0.0873	0.0651	0.0914	0.0965	0.0813	0.0926	0.0656	0.1101	0.1198	0.1120
7	0.0786	0.0841	0.0871	0.0794	0.1121	0.0711	0.0895	0.0840	0.0962	0.1079	0.1099
8	0.0708	0.0883	0.0920	0.0898	0.1136	0.0691	0.0859	0.0849	0.1032	0.0951	0.1072
9	0.0772	0.0928	0.0902	0.1017	0.1016	0.0733	0.0923	0.0813	0.1052	0.0867	0.0977
10	0.0891	0.0986	0.0843	0.0935	0.1029	0.0801	0.1016	0.0782	0.0989	0.0830	0.0898
11	0.0905	0.0998	0.0916	0.0874	0.0976	0.0865	0.0996	0.0910	0.0924	0.0773	0.0862

Визуализируем состав и доли данных по компонентам

ax = features_df.plot(kind='bar', figsize=(13,8))

ax.legend(loc='upper right')

ax.set(xlabel='Number of dimensions',

       ylabel='Relative importance',

       title='Feature importance vs Dimensions')

Обычно предпочтительнее выбирать такое количество измерений, которое

соответствует достаточно большой порции дисперсии (например 95%).

Но так как мы хотим визуализировать расположение кластеров, уменьшим размер данных до двух компонент, как было видно ранее, они объясняют ~69% дисперсии

pca = PCA(n_components=2, random_state=11)

pca.fit(df)
df_pca = pca.transform(df)
df_pca = pd.DataFrame(df_pca,columns=['Component1', 'Component2'])

df_pca.shape

(72, 2)


Выведем первые и последние пять строк полученного датафрейма

pd.concat([df_pca.head(), df_pca.tail()])


Component1
Component2
0
0.9431
0.7945
1
0.9818
0.2189
2
0.8057
-0.0369
3
0.6521
0.0999
4
0.6104
-0.0842
67
-0.7478
0.2571
68
-0.6098
0.1086
69
-0.6888
-0.0883
70
-0.7305
-0.1522
71
-0.6697
-0.0407




Перейдем к кластеризации методом k-средних. В качестве оценки качества 
кластеризации используем коэффициент силуэта, который показывает, насколько
среднее расстояние до объектов в нашем кластере отличается от среднего расстояния до объектов в других кластерах. Коэффициент силуэта ограничен диапазоном 
от -1 до 1. Значения, близкие к -1, соответствуют плохой (разбросанной) кластеризации,
значения, близкие к нулю, указывают на то, что кластеры пересекаются и
перекрываются, значения, близкие к 1, соответствуют «плотные» четко очерченные
скопления. Таким образом, чем крупнее силуэт, тем отчетливее выделяются
скопления, и они представляют собой компактные, плотно сгруппированные облака
точек. Рассмотрим варианты с количеством кластеров от 2 до 8 и
рассчитаем коэффициент силуэта для каждого варианта.

plt.figure(figsize=(20, 20))
plt.subplot(4, 2, 1)
plt.title('Instances')
plt.scatter(df_pca.loc[:,'Component1'],df_pca.loc[:,'Component2'])

colors = ['b', 'g', 'r', 'c', 'm', 'y', 'k', 'b']
markers = ['>', 's', 'D', 'v', '^', 'p', '*', '+']
tests = [2, 3, 4, 5, 6,7,8]
subplot_counter = 1
for t in range(2,9):
    subplot_counter += 1
    plt.subplot(4, 2, subplot_counter)
    kmeans_model = KMeans(n_clusters=t).fit(df_pca)
    for i, l in enumerate(kmeans_model.labels_):
        plt.plot(df_pca.loc[i,'Component1'],df_pca.loc[i,'Component2'], 
        color=colors[l], marker=markers[l],ls='None')
        plt.title('K = %s, Silhouette Coefficient = %.03f' % 
        (t,metrics.silhouette_score(df_pca, kmeans_model.labels_,
        metric='euclidean')))
        plt.scatter(kmeans_model.cluster_centers_[:,0],
        kmeans_model.cluster_centers_[:,1],s=200,c='black',marker='*')

По этим графикам можно определить, что оптимальным с точки зрения показателя 
силуэта будет разделение на три кластера. Выведем значения силуэта по кластерам

s = []
x=[]
for n_clusters in range(2,9):
    kmeans = KMeans(n_clusters=n_clusters)
    kmeans.fit(df_pca)
    labels = kmeans.labels_
    centroids = kmeans.cluster_centers_
    s.append(metrics.silhouette_score(df_pca, labels, metric='euclidean'))
    x.append(n_clusters)
plt.plot(x,s)
plt.ylabel("Silouette")
plt.xlabel("n_clusters")
plt.title("Silouette for K-means")
plt.grid(True)
sns.despine()




Для количественной оценки качества также кластеризации используется
 внутрикластерное значение SSE (искажение). В случае применения библиотеки 
scikit-leam внутрикластерное значение SSE после подгонки модели КМеаns доступно
через атрибут inertia_.На основе внутрикластерного значения SSE построен метод локтя - графический инструмент для оценки оптимального количества кластеров k.  Идея 
метода заключается  в том, чтобы идентифицировать значение k, при котором
 искажение начинает  увеличиваться быстрее всего, что прояснится, если построить
 график искажения  для разных значений k.
from sklearn.cluster import KMeans
distortions = []
for i in range(1, 11):
    km = KMeans(n_clusters=i, 
                init='k-means++', 
                n_init=10, 
                max_iter=300, 
                random_state=0)
    km.fit(df_pca)
    distortions.append(km.inertia_)
plt.plot(range(1, 11), distortions, marker='o')
plt.xlabel('Number of clusters')
plt.ylabel('Distortion')
plt.tight_layout()
plt.grid(True)



Как видно на графике  локоть находится в точке k = 3, свидетельствуя о том, 
что k = 3 - хороший выбор для этого набора данных.

При кластеризацию важную роль играет возможность ее визуализации и с этой
 точки зрения наилучшую возможность дают методы иерархической кластеризации и
ее инструмент - дендрограмма, позволяющая обрабатывать и представлять
многомерные данные.

Построим  дендограмму  с помощью библиотеки SciPy.
В SciPy используется функция, которая принимает массив данных X в качестве
 аргумента и вычисляет массив связей (linkage array) с записанными сходствами между
 кластерами. Затем  этот массив передается функции SciPy dendrogram, чтобы
построить дендрограмму  

# импортируем функцию dendrogram и функцию кластеризации ward из SciPy
from scipy.cluster.hierarchy import dendrogram, ward
# применяем кластеризацию ward к массиву данных X
# функция SciPy ward возвращает массив с расстояниями
# вычисленными в ходе выполнения агломеративной кластеризации
linkage_array = ward(df)
#строим дендрограмму для массива связей, содержащего расстояния
# между кластерами
dendrogram(linkage_array)

# делаем отметки на дереве, соответствующие двум или трем кластерам
ax = plt.gca()

bounds = ax.get_xbound()
ax.plot(bounds, [0.7, 0.7], '--', c='k')
ax.plot(bounds, [0.45, 0.45], '--', c='k')
ax.plot(bounds, [0.3, 0.3], '--', c='k')
ax.text(bounds[1], 0.7, ' two clusters', va='center', fontdict={'size': 15})
ax.text(bounds[1], 0.45, ' three clusters', va='center', fontdict={'size': 15})
ax.text(bounds[1], 0.3, ' four clusters', va='center', fontdict={'size': 15})
plt.xlabel("observation index")
plt.ylabel("cluster distance")




В библиотеке scikit-learn имеется также реализация AgglomerativeClustering, 
позволяющая выбирать количество кластеров, которые желательно возвратить.
Класс AgglomerativeClustering полезен, если мы хотим подрезать иерархическое
дерево кластеров. Устанавливая число кластеров n_cluster=3 и еще два параметра :
метрика, используемая для вычисления связи affinity ='euclidean' и linkage='ward' 
 критерий связи, определяет, какое расстояние использовать между наборами
 наблюдений 
 Алгоритм объединит пары кластеров, которые минимизируют этот критерий,
'ward' минимизирует дисперсию объединяемых кластеров. Добавим, что они также 
установлены по умолчанию.

from sklearn.cluster import AgglomerativeClustering
h_clusters= AgglomerativeClustering(n_clusters=3,affinity='euclidean',linkage='ward')
h_clusters.fit(df_pca)

Выведем полученный коэффициент силуэта :

metrics.silhouette_score(df_pca, h_clusters.labels_, metric='euclidean')

0.5120138550006529

Он получился выше, чем у KMeans

Также выведем номера кластеров

np.unique(h_clusters.labels_)

array([0, 1, 2], dtype=int64)

Визуализируем кластеры

plt.figure(figsize=(10, 5))
plt.scatter(df_pca.loc[h_clusters.labels_==0,'Component1'],
df_pca.loc[h_clusters.labels_==0,'Component2'],
            c='lightblue',edgecolor='black',marker='o',s=40,label='cluster 1')
plt.scatter(df_pca.loc[h_clusters.labels_==1,'Component1'],
df_pca.loc[h_clusters.labels_==1,'Component2'],
            c='red',edgecolor='black',marker='o',s=40,label='cluster 2')
plt.scatter(df_pca.loc[h_clusters.labels_==2,'Component1'],
df_pca.loc[h_clusters.labels_==2,'Component2'],
            c='green',edgecolor='black',marker='o',s=40,label='cluster 3')
plt.legend()




Добавим в наши данные номера кластеров

df_orig["cluster"] = h_clusters.labels_
df_orig.head().T

shop
sh1
sh2
sh3
sh4
sh5
female pants universal
0.1351
0.1609
0.1626
0.2029
0.1581
female jacket demi-season
0.1499
0.0684
0.0534
0.0412
0.0562
female jacket winter
0.0894
0.1188
0.1299
0.1622
0.1911
female sweater winter
0.0329
0.0288
0.0430
0.0355
0.0425
female sweater universal
0.0512
0.0453
0.0489
0.0564
0.0431
male pants universal
0.1775
0.2136
0.2271
0.1915
0.1778
male jacket demi-season
0.1441
0.0850
0.0831
0.0678
0.0654
male jacket winter
0.0715
0.1038
0.1091
0.1116
0.1136
male sweater winter
0.0309
0.0348
0.0421
0.0335
0.0388
male sweater universal
0.0335
0.0505
0.0330
0.0287
0.0369
male hoody universal
0.0841
0.0900
0.0677
0.0685
0.0765
cluster
2.0000
2.0000
2.0000
2.0000
2.0000

Для того, чтобы понять, как отличаются ассортиментные матрицы кластеров
сгруппируем данные и выведем средние показатели по кластерам

df_by_cl = df_orig.groupby('cluster')
df_cl=pd.DataFrame(df_by_cl.mean())
df_cl['cluster']=df_cl.index
df_cl.T

cluster
0
1
2
female pants universal
0.1101
0.0804
0.1477
female jacket demi-season
0.0473
0.0384
0.0667
female jacket winter
0.2079
0.2552
0.1334
female sweater winter
0.0303
0.0247
0.0342
female sweater universal
0.0283
0.0189
0.0425
male pants universal
0.1760
0.1263
0.2230
male jacket demi-season
0.0629
0.0504
0.0846
male jacket winter
0.2077
0.3041
0.1141
male sweater winter
0.0352
0.0287
0.0364
male sweater universal
0.0347
0.0270
0.0403
male hoody universal
0.0598
0.0459
0.0771
cluster
0.0000
1.0000
2.0000


Перед тем, как выводить показатели кластеров в разных разрезах переведем данные
из широкого формата в длинный

df_cl_new=pd.melt(df_cl,id_vars=['cluster'],
        value_vars=['female pants universal','female jacket demi-season',
        'female jacket winter','female sweater winter','female sweater universal',
        'male pants universal','male jacket demi-season','male jacket winter',
        'male sweater winter','male sweater universal','male hoody universal'])

df_cl_new.head()

cluster
variable
value
0
0
female pants universal
0.1101
1
1
female pants universal
0.0804
2
2
female pants universal
0.1477
3
0
female jacket demi-season
0.0473
4
1
female jacket demi-season
0.0384

Выделим в отдельные колонки показатели пола, товарной группы и сезона

df_cl_new[['sex','tg','season']] =df_cl_new['variable'].
apply(lambda x: pd.Series(str(x).split()))

df_cl_new.head()

cluster
variable
value
sex
tg
season
0
0
female pants universal
0.1101
female
pants
universal
1
1
female pants universal
0.0804
female
pants
universal
2
2
female pants universal
0.1477
female
pants
universal
3
0
female jacket demi-season
0.0473
female
jacket
demi-season
4
1
female jacket demi-season
0.0384
female
jacket
demi-season

Первая группировка по сезону

df_gr_season=df_cl_new.groupby(['cluster','season'])
df_gr_season.sum()

value
cluster
season
0
demi-season
0.1101
universal
0.4088
winter
0.4810
1
demi-season
0.0887
universal
0.2985
winter
0.6127
2
demi-season
0.1513
universal
0.5305
winter
0.3181

По сезону кластеры условно можно разделить : 0 - центральный регион, 

1 - северный (большая зимнего товара) и 2 - южный регион

Вторая группировка по полу :

df_gr_sex=df_cl_new.groupby(['cluster','sex'])
df_gr_sex.sum()

value
cluster
sex
0
female
0.4239
male
0.5762
1
female
0.4176
male
0.5824
2
female
0.4246
male
0.5754

Принципиальной разницы в долях мужского и женского товара нет

Третья группировка по товарным группам

df_gr_sex=df_cl_new.groupby(['cluster','tg'])
df_gr_sex.sum()

value
cluster
tg
0
hoody
0.0598
jacket
0.5257
pants
0.2861
sweater
0.1284
1
hoody
0.0459
jacket
0.6480
pants
0.2067
sweater
0.0994
2
hoody
0.0771
jacket
0.3988
pants
0.3707
sweater
0.1534

По товарным группам можно условно разделить кластеры также, как и по сезону. 

В северном максимальная доля курток, в южном - брюки, толстовки, свитера.

Таким образом, проведенная кластеризация и ее анализ показали целесообразность

разделения магазинов на кластеры, что должно привести к более корректной поставке

товаров в магазины.