import os
os.environ["KMP_DUPLICATE_LIB_OK"]="TRUE"
import numpy as np, matplotlib.pyplot as plt, time, warnings
warnings.filterwarnings("ignore")
import torchvision, torchvision.transforms as T
from sklearn.decomposition import PCA
from sklearn.manifold import TSNE
from sklearn.metrics import silhouette_score
import umap
te=torchvision.datasets.FashionMNIST("../Z-ML Unsupervised-Autoencoders/data",train=False,download=True,transform=T.ToTensor())
X=te.data.view(-1,784).float().numpy()/255; y=te.targets.numpy(); classes=te.classes
rng=np.random.default_rng(0); s=rng.choice(len(X),4000,replace=False); Xs,ys=X[s],y[s]
Xp=PCA(50).fit_transform(Xs)                                        # denoise to 50 PCs first (standard)
Zpca=PCA(2).fit_transform(Xs)
t=time.time(); Zt=TSNE(2,perplexity=30,init="pca",random_state=0).fit_transform(Xp); tt=time.time()-t
t=time.time(); Zu=umap.UMAP(n_neighbors=15,min_dist=0.1,random_state=0).fit_transform(Xp); tu=time.time()-t
embs=[("PCA (linear)",Zpca),("t-SNE",Zt),("UMAP",Zu)]
fig,ax=plt.subplots(1,3,figsize=(16,5))
from sklearn.manifold import trustworthiness
from sklearn.neighbors import KNeighborsClassifier, NearestNeighbors
from sklearn.model_selection import cross_val_score
knn=lambda Z: cross_val_score(KNeighborsClassifier(15),Z,ys,cv=5).mean()
hi=NearestNeighbors(n_neighbors=16).fit(Xp).kneighbors(Xp,return_distance=False)[:,1:]
def kept(Z):
    lo=NearestNeighbors(n_neighbors=16).fit(Z).kneighbors(Z,return_distance=False)[:,1:]
    return np.mean([len(set(a)&set(b))/15 for a,b in zip(hi,lo)])
for a,(nm,Z) in zip(ax,embs):
    sc=a.scatter(Z[:,0],Z[:,1],c=ys,cmap="tab10",s=6,alpha=.6); a.set_xticks([]); a.set_yticks([])
    a.set_title(f"{nm}   (15-NN accuracy {knn(Z):.3f})")
cb=fig.colorbar(sc,ax=ax,fraction=0.02,ticks=range(10)); cb.ax.set_yticklabels(classes,fontsize=7)
plt.show()
print(f"{'':16s}{'15-NN accuracy':>16}{'trustworthiness':>17}{'15-NN kept':>13}{'silhouette':>13}")
for nm,Z in embs:
    print(f"  {nm:14s}{knn(Z):>16.3f}{trustworthiness(Xp,Z,n_neighbors=15):>17.3f}{kept(Z):>13.3f}{silhouette_score(Z,ys):>13.3f}")
print("  (10 classes, so 0.100 is chance accuracy. No labels were used to fit any of the three maps.)")
print("\nPCA's linear projection overlaps the classes badly: it keeps only 14% of each point's 50-D neighbourhood, and a")
print("classifier reads the label off it barely better than half the time. t-SNE and UMAP keep three to four times as much")
print("of the local structure and carry the label at ~0.76. For comparison, the autoencoder's 2-D latent code in the")
print("previous notebook reached 0.673 on this same data -- better than PCA, and still short of methods built for the job.")
print("\nNotice that silhouette ranks the three the same way but makes the gap look small (0.15 against -0.04), because it")
print("scores geometry rather than fidelity. Section 4 shows what it does on data with no structure in it at all.")
print(f"\nTiming, since it is usually claimed the other way round: t-SNE {tt:.1f}s, UMAP {tu:.1f}s on these {len(Xs)} points.")

                  15-NN accuracy  trustworthiness   15-NN kept   silhouette

  PCA (linear)             0.520            0.921        0.137       -0.037

  t-SNE                    0.766            0.992        0.519        0.152

  UMAP                     0.750            0.985        0.429        0.137
  (10 classes, so 0.100 is chance accuracy. No labels were used to fit any of the three maps.)

PCA's linear projection overlaps the classes badly: it keeps only 14% of each point's 50-D neighbourhood, and a
classifier reads the label off it barely better than half the time. t-SNE and UMAP keep three to four times as much
of the local structure and carry the label at ~0.76. For comparison, the autoencoder's 2-D latent code in the
previous notebook reached 0.673 on this same data -- better than PCA, and still short of methods built for the job.

Notice that silhouette ranks the three the same way but makes the gap look small (0.15 against -0.04), because it
scores geometry rather than fidelity. Section 4 shows what it does on data with no structure in it at all.

Timing, since it is usually claimed the other way round: t-SNE 8.6s, UMAP 13.8s on these 4000 points.

s2=rng.choice(len(X),2000,replace=False); Xp2=PCA(50).fit_transform(X[s2]); y2=y[s2]
fig,ax=plt.subplots(1,3,figsize=(15,4.6))
for a,perp in zip(ax,[5,30,100]):
    Z=TSNE(2,perplexity=perp,init="pca",random_state=0).fit_transform(Xp2)
    a.scatter(Z[:,0],Z[:,1],c=y2,cmap="tab10",s=6,alpha=.6); a.set_xticks([]); a.set_yticks([]); a.set_title(f"perplexity = {perp}")
plt.suptitle("Same data, three perplexities — the embedding is not unique"); plt.tight_layout(); plt.show()
print("Low perplexity shatters clusters into small groups; high perplexity fuses them. Always inspect several settings --")
print("and never read absolute distances or cluster SIZES off a t-SNE/UMAP plot; only local neighbourhood/grouping is trustworthy.")

Low perplexity shatters clusters into small groups; high perplexity fuses them. Always inspect several settings --
and never read absolute distances or cluster SIZES off a t-SNE/UMAP plot; only local neighbourhood/grouping is trustworthy.

import pandas as pd
d=pd.read_csv("stocks_weekly.csv",index_col=0); tickers=list(d.columns)
Rz=(d.values-d.values.mean(0))/d.values.std(0)
sectors={**dict.fromkeys(["AAPL","ADBE","AMD","AVGO","CRM","CSCO","INTC","MSFT","NVDA","ORCL"],"Tech"),
         **dict.fromkeys(["GOOGL","META"],"Comm"),
         **dict.fromkeys(["AXP","BAC","BLK","C","GS","JPM","MS","WFC"],"Financials"),
         **dict.fromkeys(["ABBV","ABT","JNJ","LLY","MRK","PFE","TMO","UNH"],"Healthcare"),
         **dict.fromkeys(["AMZN","HD","LOW","MCD","NKE","TSLA"],"Cons.Disc"),
         **dict.fromkeys(["COST","KO","PEP","PG","WMT"],"Cons.Staples"),
         **dict.fromkeys(["COP","CVX","XOM"],"Energy"),
         **dict.fromkeys(["BA","CAT","GE","HON","LMT","UPS"],"Industrials")}
sec=[sectors.get(t,"Other") for t in tickers]; uniq=sorted(set(sec)); cmap=plt.get_cmap("tab10")
emb=umap.UMAP(n_neighbors=8,min_dist=0.25,random_state=0).fit_transform(Rz.T)
fig,ax=plt.subplots(figsize=(9,6.5))
for k,u in enumerate(uniq):
    m=[i for i in range(len(sec)) if sec[i]==u]; ax.scatter(emb[m,0],emb[m,1],color=cmap(k),s=60,label=u)
for i,t in enumerate(tickers): ax.annotate(t,(emb[i,0],emb[i,1]),fontsize=6.5,alpha=.8)
ax.set_xticks([]); ax.set_yticks([]); ax.set_title("UMAP of 48 stocks by return behaviour — sectors emerge unsupervised"); ax.legend(fontsize=8,loc="best")
plt.tight_layout(); plt.show()
from sklearn.neighbors import NearestNeighbors as _NN
sec_a=np.array(sec)
def sector_nb(Z,k=5):
    nb=_NN(n_neighbors=k+1).fit(Z).kneighbors(Z,return_distance=False)[:,1:]
    return np.mean([np.mean(sec_a[row]==sec_a[i]) for i,row in enumerate(nb)])
chance=np.mean([np.mean(sec_a==u) for u in sec_a])
print(f"Share of a stock's 5 nearest neighbours that share its sector (chance would be {chance:.3f}):")
print(f"   raw 312-week return space (the ceiling) : {sector_nb(Rz.T):.3f}")
print(f"   UMAP, 2-D                               : {sector_nb(emb):.3f}")
print(f"   PCA, 2-D                                : {sector_nb(PCA(2).fit_transform(Rz.T)):.3f}")
_a=[sector_nb(umap.UMAP(n_neighbors=8,min_dist=0.25,random_state=sd).fit_transform(Rz.T)) for sd in range(5)]
print(f"   UMAP across 5 random seeds              : {np.mean(_a):.3f} +/- {np.std(_a):.3f}")
print("\nThe sector structure is in the returns, not in the algorithm: the raw 312-dimensional space already puts a stock")
print("next to its own sector half the time, against a chance rate near 0.15. What UMAP does is carry almost all of that")
print("down to two dimensions -- where PCA's linear projection loses noticeably more of it. That is the honest version of")
print("the claim, and a real one: the value of the map is that it makes visible, at almost no cost in fidelity, structure")
print("that already existed in 312 dimensions where nobody could look at it.")
print(f"\nOne caveat on scale: UMAP is being asked to embed only {len(tickers)} points here, far below the sample sizes it is")
print("designed for. It is stable across seeds on this data, but a 48-point map is a picture to reason from, not evidence.")

Share of a stock's 5 nearest neighbours that share its sector (chance would be 0.147):
   raw 312-week return space (the ceiling) : 0.521
   UMAP, 2-D                               : 0.504
   PCA, 2-D                                : 0.400
   UMAP across 5 random seeds              : 0.507 +/- 0.010

The sector structure is in the returns, not in the algorithm: the raw 312-dimensional space already puts a stock
next to its own sector half the time, against a chance rate near 0.15. What UMAP does is carry almost all of that
down to two dimensions -- where PCA's linear projection loses noticeably more of it. That is the honest version of
the claim, and a real one: the value of the map is that it makes visible, at almost no cost in fidelity, structure
that already existed in 312 dimensions where nobody could look at it.

One caveat on scale: UMAP is being asked to embed only 48 points here, far below the sample sizes it is
designed for. It is stable across seeds on this data, but a 48-point map is a picture to reason from, not evidence.

# The claim to test: do these methods invent clusters? Give them data guaranteed to have none.
noise=rng.normal(size=(2000,50))                                   # one isotropic Gaussian population, no structure
Zt_n=TSNE(2,perplexity=30,init="pca",random_state=0).fit_transform(noise)
Zu_n=umap.UMAP(n_neighbors=15,min_dist=0.1,random_state=0).fit_transform(noise)
fig,ax=plt.subplots(1,2,figsize=(11,4.6))
for a,(nm,Z) in zip(ax,[("t-SNE of pure noise",Zt_n),("UMAP of pure noise",Zu_n)]):
    a.scatter(Z[:,0],Z[:,1],s=5,alpha=.5,color="#2b6cb0"); a.set_xticks([]); a.set_yticks([]); a.set_title(nm)
plt.suptitle("2,000 points drawn from a single Gaussian — every apparent group here is an artefact")
plt.tight_layout(); plt.show()
from sklearn.cluster import KMeans
print("Silhouette of a 5-cluster k-means solution, in each space:")
for nm,Z in [("the original 50-D noise",noise),("its t-SNE embedding",Zt_n),("its UMAP embedding",Zu_n)]:
    print(f"   {nm:26s} {silhouette_score(Z,KMeans(5,n_init=10,random_state=0).fit_predict(Z)):.3f}")
print(f"   {'real Fashion-MNIST t-SNE':26s} {silhouette_score(Zt,ys):.3f}   <- against its TRUE classes, from section 1")
print("\nIn 50 dimensions the noise correctly scores ~0: there is nothing to find. Projected to 2-D by either method it")
print("scores HIGHER than the genuine ten-class structure of Fashion-MNIST does. The clusters in those two panels are")
print("produced entirely by the algorithms, and no amount of staring at the picture would reveal that. This is also why")
print("silhouette cannot be used to judge these embeddings -- it rewards the geometry the method imposes.")

print("\n\nThe two comparative claims usually made for UMAP, measured rather than repeated:")
print("\n1. SPEED.")
for _n in (1000,2000):
    _X=PCA(50).fit_transform(X[rng.choice(len(X),_n,replace=False)])
    _a=time.time(); TSNE(2,perplexity=30,init="pca",random_state=0).fit_transform(_X); _a=time.time()-_a
    _b=time.time(); umap.UMAP(n_neighbors=15,min_dist=0.1,random_state=0).fit_transform(_X); _b=time.time()-_b
    print(f"   n={_n:>5}: t-SNE {_a:>6.1f}s   UMAP {_b:>6.1f}s")
print(f"   n={len(Xs):>5}: t-SNE {tt:>6.1f}s   UMAP {tu:>6.1f}s   (section 1)")
print("   t-SNE is faster at every size tried here. UMAP's advantage is a large-n claim, and these are not large n;")
print("   scikit-learn's Barnes-Hut t-SNE is well optimised in exactly this range.")

print("\n2. GLOBAL STRUCTURE. Rank correlation between pairwise distances in 50-D and in each 2-D map:")
from scipy.spatial.distance import pdist
from scipy.stats import spearmanr
_s=rng.choice(len(Xp),1500,replace=False); _dh=pdist(Xp[_s])
for nm,Z in [("PCA-2",Zpca),("t-SNE (PCA init)",Zt),("UMAP",Zu)]:
    print(f"   {nm:18s} {spearmanr(_dh,pdist(Z[_s])).statistic:+.3f}")
print("   UMAP preserves global geometry LESS well than t-SNE here, not more. The usual comparison is against t-SNE")
print("   started from a random layout; initialised from PCA, as it is throughout this notebook, t-SNE keeps the coarse")
print("   arrangement it was handed. And PCA beats both, which is the point -- it is the method that optimises for it.")

Silhouette of a 5-cluster k-means solution, in each space:
   the original 50-D noise    0.014

   its t-SNE embedding        0.306
   its UMAP embedding         0.334

   real Fashion-MNIST t-SNE   0.152   <- against its TRUE classes, from section 1

In 50 dimensions the noise correctly scores ~0: there is nothing to find. Projected to 2-D by either method it
scores HIGHER than the genuine ten-class structure of Fashion-MNIST does. The clusters in those two panels are
produced entirely by the algorithms, and no amount of staring at the picture would reveal that. This is also why
silhouette cannot be used to judge these embeddings -- it rewards the geometry the method imposes.


The two comparative claims usually made for UMAP, measured rather than repeated:

1. SPEED.

   n= 1000: t-SNE    0.9s   UMAP    1.1s

   n= 2000: t-SNE    2.4s   UMAP    3.2s
   n= 4000: t-SNE    8.6s   UMAP   13.8s   (section 1)
   t-SNE is faster at every size tried here. UMAP's advantage is a large-n claim, and these are not large n;
   scikit-learn's Barnes-Hut t-SNE is well optimised in exactly this range.

2. GLOBAL STRUCTURE. Rank correlation between pairwise distances in 50-D and in each 2-D map:

   PCA-2              +0.891

   t-SNE (PCA init)   +0.695

   UMAP               +0.598
   UMAP preserves global geometry LESS well than t-SNE here, not more. The usual comparison is against t-SNE
   started from a random layout; initialised from PCA, as it is throughout this notebook, t-SNE keeps the coarse
   arrangement it was handed. And PCA beats both, which is the point -- it is the method that optimises for it.

print("QUICK REFERENCE -- which dimensionality-reduction tool when:")
print("  PCA / factor analysis : compression + interpretable factors; distances & variance meaningful; linear")
print("  Autoencoder           : nonlinear compression; reusable codes; scales to new data")
print("  t-SNE                 : 2-D visualization of clusters; local structure only; no new-point transform")
print("  UMAP                  : 2-D visualization; local structure only; CAN transform new points")
print("\nRule of thumb: PCA/AE to COMPRESS and feed a model; t-SNE/UMAP to SEE the data -- never read distances off the latter.")
print("On the two claims usually used to choose UMAP over t-SNE, the measurements above say: not at this scale. The one")
print("that does separate them is that UMAP can embed new points and vanilla t-SNE cannot -- which matters whenever the")
print("map has to be applied to data arriving later, and is a difference in kind rather than in degree.")

QUICK REFERENCE -- which dimensionality-reduction tool when:
  PCA / factor analysis : compression + interpretable factors; distances & variance meaningful; linear
  Autoencoder           : nonlinear compression; reusable codes; scales to new data
  t-SNE                 : 2-D visualization of clusters; local structure only; no new-point transform
  UMAP                  : 2-D visualization; local structure only; CAN transform new points

Rule of thumb: PCA/AE to COMPRESS and feed a model; t-SNE/UMAP to SEE the data -- never read distances off the latter.
On the two claims usually used to choose UMAP over t-SNE, the measurements above say: not at this scale. The one
that does separate them is that UMAP can embed new points and vanilla t-SNE cannot -- which matters whenever the
map has to be applied to data arriving later, and is a difference in kind rather than in degree.

Unsupervised Learning IV — t-SNE & UMAP¶

Manifold learning for visualization, vs PCA and autoencoders¶

1. The headline — four embeddings of Fashion-MNIST¶

2. How t-SNE behaves — the perplexity knob¶

3. A finance application — mapping the market by sector¶

4. Caveats — how to read (and not read) these plots¶

5. Summary — and the end of the Unsupervised subsection¶