Спектральное вложение для нелинейного понижения размерности.

Спектральное вложение (Spectral Embedding) - это подход к вычислению нелинейного вложения. Scikit-learn реализует собственные лапласианские карты, которые находят низкоразмерное представление данных с помощью спектрального разложения лапласиана графа. Сгенерированный граф можно рассматривать как дискретную аппроксимацию многообразия низкой размерности в пространстве высокой размерности. Минимизация функции стоимости на основе графика гарантирует, что точки, близкие друг к другу на многообразии, отображаются близко друг к другу в низкоразмерном пространстве с сохранением локальных расстояний.

Алгоритм Spectral Embedding (Laplacian Eigenmaps)

Входными данными алгоритма являются:

• обучающая выборка $D = \{ x_1 ,x_2 ,...,x_p \} \subset \Re ^n$

• функция близости $K (·, ·)$: $\Re ^n \times \Re ^n \to [0,\infty )$. Например, $\;K (x, x^*)=0$ для $\left| {x - x^* } \right| < \varepsilon$ и $0$ – иначе, где $\varepsilon > 0$.Однако, функция близости может зависеть и от объектов обучающей выборки.

Выходом алгоритма являются низкоразмерные описания точек обучающей выборки: $\;y_1 ,...,y_p \in \Re ^m ,\;\;m < n$

1.По выборке $D = \{ x_1 ,x_2 ,...,x_p \} \subset \Re ^n$ строим матрицу близости (affinity matrix) $M = (M{}_{ij})_{i,j = 1}^p $ размера $p \times p$ $M{}_{ij} = K(x_i ,x_j )$

2.Определим диагональную матрицу $Q$ : $Q_{ii} = \sum\nolimits_{j = 1}^p {M_{ij} }$ , и вычисляем матрицу:

$$\tilde M = (\tilde M_{ij} )_{i,j}^p = Q^{ - 0.5} \times M \times Q^{ - 0.5} {\rm \;\;\;\;(1)}$$

3.Находим $m + 1$ наименьшее собственное значение $\lambda _0 = 0 < \lambda _1 \le ... \le \lambda _m$ матрицы $\tilde M$ и соответствующие им собственные векторы $V_0 ,...,V_m \in \Re ^p$

4.Используя полученную матрицу собственных векторов, получаем вложение: $Y = (y_1 \left| {y{}_2} \right|...\left| {y_p } \right.)^T = (V_1 \left| {V{}_2} \right|...\left| {V_m } \right.)$ , где $y_1 ,...,y_p \in \Re ^m$

Этим алгоритмом решается задача минимизации:

$$\sum\limits_{i,j = 1}^p {K(x_i } ,x_j ) \cdot \left| {y{}_i - y_j } \right|^2 {\rm \;\;\;\;(2)}$$

При ограничениях $\sum\nolimits_{i = 1}^p {y_i } = 0\;\;$ и $\;\;\sum\nolimits_{i = 1}^p {y_i^2 } = 1$

Пример реализации алгоритма $LE$ на $python$ и $R$

Для рачетов использовались данный о Российских регионах. В $python$ реализации исходным было $3D$ пространство, регионы размечены на три класса.Для R реализаци исходное пространство $4D$, а регионы размечены на четыре класса.

В sklearn.manifold.SpectralEmbedding мы задаем следующие параметры n.components = 2 количество координат для коллектора (размерность выходного редуцируемого пространства) n.neighbors=5 число ближайших соседей для построения графа ближайших_соседей,random.state=0 -начальное значение генератор случайных чисел для воспроизводимости результатов остальные параметры принимаются по умолчанию

Код $python$ :

import warnings
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
import scipy.stats as st
import matplotlib.pyplot as plt
from pandas import Series, DataFrame
import pandas as pd
#import umap
#import umap.plot
#from sklearn.manifold import Isomap
#from sklearn.manifold import LocallyLinearEmbedding
from sklearn.manifold import SpectralEmbedding

warnings.simplefilter('ignore')
result =pd.read_table('dan04.txt',sep='\s+')
print(result)

X =DataFrame(result, columns=['P10','P11','P14'])
X=np.array (X)

SUMX=X[:,0] +X[:,1] + X[:,2]
X[:,0] =X[:,0] /SUMX
X[:,1] =X[:,1] /SUMX
X[:,2] =X[:,2] /SUMX



YT=result['P18']
y=np.array(YT)

XX1=result.query("P18 in [0]")
XX2=result.query("P18 in [1]")
XX3=result.query("P18 in [2]")

XX1 =DataFrame(XX1, columns=['P10','P11','P14'])
XX1=np.array (XX1)
XX2 =DataFrame(XX2, columns=['P10','P11','P14'])
XX2=np.array (XX2)
XX3 =DataFrame(XX3, columns=['P10','P11','P14'])
XX3=np.array (XX3)

SUMX1=XX1[:,0] +XX1[:,1] + XX1[:,2]
XX1[:,0] =XX1[:,0] /SUMX1
XX1[:,1] =XX1[:,1] /SUMX1
XX1[:,2] =XX1[:,2] /SUMX1

SUMX2=XX2[:,0] +XX2[:,1] + XX2[:,2]
XX2[:,0] =XX2[:,0] /SUMX2
XX2[:,1] =XX2[:,1] /SUMX2
XX2[:,2] =XX2[:,2] /SUMX2

SUMX3=XX3[:,0] +XX3[:,1] + XX3[:,2]
XX3[:,0] =XX3[:,0] /SUMX3
XX3[:,1] =XX3[:,1] /SUMX3
XX3[:,2] =XX3[:,2] /SUMX3

# Преобразование данных, уменьшающее размерность до 2
#результат имеет два измерения
mod_LE  = SpectralEmbedding(n_components=2, n_neighbors = 5,
 random_state=0 )
X_LE = mod_LE.fit_transform(X)
X_LE.shape



XL =DataFrame(X_LE, columns=['1','2'])
XL = XL.join(YT)

XL1=XL.query("P18 in [0]")
XL2=XL.query("P18 in [1]")
XL3=XL.query("P18 in [2]")

XL1 =DataFrame(XL1, columns=['1','2'])
XL1=np.array (XL1)
XL2 =DataFrame(XL2, columns=['1','2'])
XL2=np.array (XL2)
XL3 =DataFrame(XL3, columns=['1','2'])
XL3=np.array (XL3)

#построение исходного графика 3 D
fig = plt.figure(1, figsize=(10, 8))
ax = fig.add_subplot(111, projection='3d')

xs1 = XX1[:,0] 
ys1 = XX1[:,1] 
zs1 = XX1[:,2] 

xs2 = XX2[:,0] 
ys2 = XX2[:,1] 
zs2 = XX2[:,2] 

xs3 = XX3[:,0] 
ys3 = XX3[:,1] 
zs3 = XX3[:,2] 

ax.scatter(xs1, ys1, zs1,c='g', cmap=plt.cm.nipy_spectral, edgecolor='k',s=80,label='A')
ax.scatter(xs2, ys2, zs2,c='r', cmap=plt.cm.nipy_spectral, edgecolor='k',s=80,label='P')
ax.scatter(xs3, ys3, zs3,c='b', cmap=plt.cm.nipy_spectral, edgecolor='k',s=80,label='D')
ax.set_title("Исходное пространство предикторов",fontsize=16, fontweight='bold')
ax.set_xlabel('P10',fontsize=14, fontweight='bold')
ax.set_ylabel('P11',fontsize=14, fontweight='bold')
ax.set_zlabel('P14',fontsize=14, fontweight='bold')
ax.legend(loc="best")

y = np.choose(y, [1, 2, 0]).astype(np.float)
fig = plt.figure(2,figsize=(10, 8))
plt.clf()
plt.cla()
plt.scatter(XL1[:, 0], XL1[:, 1], c='g', cmap=plt.cm.nipy_spectral, edgecolor='k',s=100,label='A')
plt.scatter(XL2[:, 0], XL2[:, 1], c='r', cmap=plt.cm.nipy_spectral, edgecolor='k',s=100,label='P')
plt.scatter(XL3[:, 0], XL3[:, 1], c='b', cmap=plt.cm.nipy_spectral, edgecolor='k',s=100,label='D')
plt.title("Снижение размерности методом LE",fontsize=16, fontweight='bold')
plt.xlabel("LE1",fontsize=14, fontweight='bold')
plt.ylabel("LE2",fontsize=14, fontweight='bold')
plt.legend(loc='best')

plt.show()

          P10      P11     P14        P16  P18
0    148863.0   710829  257038   336148.5    0
1       262.0   218544   85146   253157.2    0
2      5005.0   448428   29651   225731.2    1
3      7726.0   448223  219151   552288.4    0
4       972.0   152792   16085   164827.4    0
..        ...      ...     ...        ...  ...
76  1016799.0    60300   11147   148497.4    2
77    10330.0     6632    5772    24075.7    2
78    67502.0     1055    1334     9573.8    2
79  1808823.0  6413663    7308  4798454.0    2
80    24035.0  2615910       0  1412406.0    1

[81 rows x 5 columns]

В в функции llle (LE) R пакета 'Rdimtools' мы задаем следующие параметры ndim =2 целочисленное целевое измерение., K= 6 размер ближайшего окружения для каждой точки данных,P=3 размер искусственного соседства,bandwidth = (0.1,0.5,0.9) параметр масштаба для гауссова ядра. Он должен быть в (0, 1), остальные параметры по умолчанию.

Код $R$ :

options(warn=-1)
library(gmodels)
library(lattice)
library(vegan)
library(grDevices)
library(reshape2)
library(ggplot2)
library(ggthemes)
library(psych)

load("xRegion.RData")

#ФОРМИРОВАНИЕ data.frame ОСНОВНЫЕ СОЦИАЛЬНО-ЭКОНОМИЧЕСКИЕ ПОКАЗАТЕЛИ РЕГИОНОВ РФ 
xreg <-data.frame(P1=x$Pok_001,P2=x$Pok_002,P3=x$Pok_003,P4=x$Pok_004,
P5=x$Pok_005,P6=x$Pok_006,P7=x$Pok_007,P8=x$Pok_008,P9=x$Pok_009,P10=x$Pok_010,
P11=x$Pok_011,P12=x$Pok_012,P13=x$Pok_013,P14=x$Pok_014,P15=x$Pok_015,P16=x$Pok_016,
P17=x$Pok_017,P18=x$Pok_018)
#xreg
#Формирование фактора тип экономики региона
xreg_f <- as.factor(x$Pok_018)

#выбор исследуемого подмножества показателей регионов
w=c(10,11,14,16)
reg.ekon <-xreg[,w]
valreg <-rowSums(reg.ekon)
xreg_nor <-reg.ekon/valreg
Y <- xreg_nor
X=Y
(X <- as.matrix(X))

#Метод le =====================================================

library(Rdimtools)
label=xreg_f

# see the effect bandwidth
out1 = do.llle(X, bandwidth=0.1,K=6,P=3)
out2 = do.llle(X, bandwidth=0.5,K=6,P=3)
out3 = do.llle(X, bandwidth=0.9,K=6,P=3)
# visualize the results

opar <- par(no.readonly=TRUE)
par(mfrow=c(1,3))
plot(out1$Y, col=label, main="bandwidth=0.1")
plot(out2$Y, col=label, main="bandwidth=0.5")
plot(out3$Y, col=label, main="bandwidth=0.9")
par(opar)

#График регионов по классам в двухмерном пространстве LE  

yx.le <-data.frame( le1=out2$Y[,1],le2=out2$Y[,2],Econ=xreg_f)

x11()
g4 <-ggplot(data= yx.le, aes(x=le1,y=le2,colour = Econ),
colors = "Accent")
g4 <-g4 +geom_point(size=6)
g4 <-g4 +labs(title ="Регионы в пространсте двух компонент LE",
x="le1", y="le2",colour = "Классы")
#g4 <-g4 + theme_stata() + scale_colour_stata()
g4 <-g4 + theme_fivethirtyeight() + scale_colour_manual( values = c('darkgoldenrod4', 
'blue','red','green'))
g4

A matrix: 81 × 4 of type dbl

P10	P11	P14	P16
0.1024607357	0.48925564	0.176916377	0.2313672
0.0004702848	0.39228216	0.152835387	0.4544122
0.0070610788	0.63264445	0.041831778	0.3184627
0.0062946660	0.36518432	0.178550653	0.4499704
0.0029042980	0.45653652	0.048061351	0.4924978
0.0036612406	0.77233396	0.040769188	0.1832356
0.0015789528	0.53479642	0.062721044	0.4009036
0.1412626871	0.30145091	0.227111091	0.3301753
0.0061860875	0.66402720	0.104654706	0.2251320
0.0025102570	0.51220351	0.021356685	0.4639295
0.0007240736	0.36139077	0.225483397	0.4124018
0.0029963642	0.54662759	0.102450806	0.3479252
0.0044219285	0.51692122	0.059618262	0.4190386
0.0004073696	0.33198793	0.259307072	0.4082976
0.0016961452	0.53906992	0.068067881	0.3911661
0.0057799303	0.66065459	0.061796314	0.2717692
0.0023675488	0.58979632	0.051511573	0.3563246
0.2517075032	0.36705018	0.013576118	0.3676662
0.5272218355	0.25586073	0.013373691	0.2035437
0.4244979655	0.27268070	0.012134587	0.2906867
0.0008776425	0.76185634	0.032174724	0.2050913
0.0230355529	0.72078980	0.043424792	0.2127499
0.0098563473	0.68108995	0.057009183	0.2520445
0.2124016359	0.40916960	0.004022608	0.3744062
0.0034276470	0.58540020	0.075600275	0.3355719
0.0076922322	0.40696974	0.145005841	0.4403322
0.0196143721	0.30366147	0.126594499	0.5501297
0.0000000000	0.01904875	0.543995445	0.4369558
0.0352299803	0.22928658	0.107467128	0.6280163
0.0185623909	0.36245810	0.135234114	0.4837454
...	...	...	...
0.167042681	0.49841540	0.040367941	0.29417397
0.045323146	0.43000046	0.139309411	0.38536698
0.016284635	0.52403595	0.077257147	0.38242226
0.012267281	0.40167539	0.150432191	0.43562514
0.023469405	0.60376330	0.026051904	0.34671539
0.689281644	0.20286297	0.008220854	0.09963453
0.036633955	0.67476950	0.054015521	0.23458102
0.061176556	0.08104528	0.266281886	0.59149628
0.486146898	0.00878109	0.104101684	0.40097033
0.281106629	0.34792578	0.052480425	0.31848716
0.006186605	0.39759068	0.162362975	0.43385973
0.291003668	0.46942750	0.030307428	0.20926141
0.394878732	0.33598888	0.041555341	0.22757704
0.505473223	0.29954857	0.021596419	0.17338179
0.075183481	0.44251036	0.066501734	0.41580443
0.002749624	0.69173970	0.067695274	0.23781540
0.365852572	0.31607420	0.051640443	0.26643279
0.097501672	0.22863783	0.056512909	0.61734759
0.733874277	0.03391619	0.023521003	0.20868853
0.338709249	0.07432856	0.071566390	0.51539580
0.097758834	0.60492083	0.037476188	0.25984414
0.032009134	0.33473293	0.056708828	0.57654911
0.129754570	0.38570873	0.024285486	0.46025121
0.192388305	0.10900480	0.152317667	0.54628923
0.759687556	0.02743551	0.016557379	0.19631956
0.822158420	0.04875708	0.009013187	0.12007131
0.220680756	0.14168004	0.123307776	0.51433143
0.849457873	0.01327632	0.016787307	0.12047850
0.138838545	0.49228899	0.000560935	0.36831153
0.005931125	0.64552898	0.000000000	0.34853990