# Diffusion maps

## Algorithm - diffusion maps embeding

#### 2. Normalized affinity matrix $Q$

$P$ 這矩陣可看成是個機率矩陣, 其第 $i$ 個 row 表示從 $x_{i}$ 這個點跳到其他點的機率分佈.

## Implementation in julia

using Plots
using Distances
using LinearAlgebra


# generating the data - a spiral curve
# n: number of sampling
n = 300;

theta = range(0,stop=6*pi,length=n);
r = range(0,stop=1, length=n);
# x and y-coordinates
x = r.*cos.(theta);
y = r.*sin.(theta);
# X is a n-by-2 data matrix
X = [x y];


plot(reshape(X[:,1], 100, 3), reshape(X[:,2], 100, 3), aspect_ratio=:equal, leg=false)


#### 1. Affinity matrix $K$

E = pairwise(Euclidean(), X, dims=1);
sigma = 0.1;
K = exp.(-(E.^2)./(sigma^2));


#### 2. Normalized affinity matrix $Q$

Q = zeros(n,n);
d_sq = zeros(n);
for ii = 1:n
d_sq[ii] = sqrt(sum(K[ii,:]));
end
for ii = 1:n
for jj = 1:n
Q[ii,jj] = K[ii,jj]/(d_sq[ii]*d_sq[jj]);
end
end
U,S,V = svd(Q);


#### 3. 求出 $P$ 矩陣的 eigenvalues 跟 eigenvectors

$P$ 矩陣的 eigenvalues 跟 $Q$ 一樣, eigenvectors 也可以依公式得到:

for ii = 1 : n
V[ii,:] = V[ii,:]/d_sq[ii];
end


p1 = scatter(S[2:11], title="eigenvalues 2:11", leg=false);
p2 = plot(log.(S), title="eigenvalues in log", leg=false);
plot(p1, p2, layout=2, fmt=:png)


#### 4. Define diffusion map $Y$

Y = zeros(n,3);
for ii = 1:3
Y[:,ii] = V[:,ii+1].*S[ii+1];
end


p1 = plot(reshape(Y[:,1], 100, 3), reshape(Y[:,2], 100, 3), title="2D", leg=false)
p2 = plot(reshape(Y[:,1], 100, 3), reshape(Y[:,2], 100, 3), reshape(Y[:,3], 100, 3), title="3D", leg=false)
plot(p1, p2, layout=(1,2), fmt=:png)


## Extension

Diffusion maps 也常被拿來搭配 k-means 做成分群演算法, 算是 spectral clustering 的其中一種.

##### Te-Sheng Lin (林得勝)
###### Associate Professor

The focus of my research is concerned with the development of analytical and computational tools, and further to communicate with scientists from other disciplines to solve engineering problems in practice.