import numpy as np
from sklearn.cluster import KMeans

# Simulate synthetic network traffic data (e.g., [duration, bytes]).
# Three normal clusters and one small attack cluster.
rng = np.random.RandomState(42)
normal1 = rng.normal(loc=[50, 500], scale=[10, 100], size=(500, 2))   # Cluster 1
normal2 = rng.normal(loc=[60, 1500], scale=[8, 200], size=(500, 2))   # Cluster 2
normal3 = rng.normal(loc=[70, 3000], scale=[5, 300], size=(500, 2))   # Cluster 3
attack = rng.normal(loc=[200, 800], scale=[5, 50], size=(50, 2))      # Small attack cluster

X = np.vstack([normal1, normal2, normal3, attack])
# Run K-Means clustering into 4 clusters (we expect it to find the 4 groups)
kmeans = KMeans(n_clusters=4, random_state=0, n_init=10)
labels = kmeans.fit_predict(X)

# Analyze resulting clusters
clusters, counts = np.unique(labels, return_counts=True)
print(f"Cluster labels: {clusters}")
print(f"Cluster sizes: {counts}")
print("Cluster centers (duration, bytes):")
for idx, center in enumerate(kmeans.cluster_centers_):
print(f"  Cluster {idx}: {center}")

In questo esempio, K-Means dovrebbe trovare 4 cluster. Il piccolo cluster di attacco (con una durata insolitamente alta ~200) formerà idealmente il proprio cluster data la sua distanza dai cluster normali. Stampiamo le dimensioni e i centri dei cluster per interpretare i risultati. In uno scenario reale, si potrebbe etichettare il cluster con pochi punti come potenziali anomalie o ispezionare i suoi membri per attività malevole.

Clustering Gerarchico

Il clustering gerarchico costruisce una gerarchia di cluster utilizzando un approccio dal basso verso l'alto (agglomerativo) o un approccio dall'alto verso il basso (divisivo):

1. Agglomerativo (Dal Basso Verso l'Alto): Iniziare con ogni punto dati come un cluster separato e unire iterativamente i cluster più vicini fino a quando rimane un singolo cluster o viene soddisfatta una condizione di arresto.
2. Divisivo (Dall'Alto Verso il Basso): Iniziare con tutti i punti dati in un singolo cluster e dividere iterativamente i cluster fino a quando ogni punto dati è il proprio cluster o viene soddisfatta una condizione di arresto.

Il clustering agglomerativo richiede una definizione della distanza inter-cluster e un criterio di collegamento per decidere quali cluster unire. I metodi di collegamento comuni includono il collegamento singolo (distanza dei punti più vicini tra due cluster), il collegamento completo (distanza dei punti più lontani), il collegamento medio, ecc., e la metrica di distanza è spesso euclidea. La scelta del collegamento influisce sulla forma dei cluster prodotti. Non è necessario specificare in anticipo il numero di cluster K; è possibile "tagliare" il dendrogramma a un livello scelto per ottenere il numero desiderato di cluster.

Il clustering gerarchico produce un dendrogramma, una struttura ad albero che mostra le relazioni tra i cluster a diversi livelli di granularità. Il dendrogramma può essere tagliato a un livello desiderato per ottenere un numero specifico di cluster.

Assunzioni e Limitazioni

Il clustering gerarchico non assume una forma particolare del cluster e può catturare cluster annidati. È utile per scoprire tassonomie o relazioni tra gruppi (ad esempio, raggruppare il malware per sottogruppi familiari). È deterministico (nessun problema di inizializzazione casuale). Un vantaggio chiave è il dendrogramma, che fornisce informazioni sulla struttura di clustering dei dati a tutte le scale – gli analisti della sicurezza possono decidere un taglio appropriato per identificare cluster significativi. Tuttavia, è computazionalmente costoso (tipicamente $O(n^2)$ tempo o peggio per implementazioni naive) e non fattibile per dataset molto grandi. È anche una procedura avido – una volta che una fusione o una divisione è stata effettuata, non può essere annullata, il che può portare a cluster subottimali se si verifica un errore all'inizio. Gli outlier possono anche influenzare alcune strategie di collegamento (il collegamento singolo può causare l'effetto "chaining" in cui i cluster si collegano tramite outlier).

from sklearn.cluster import AgglomerativeClustering
from scipy.cluster.hierarchy import linkage, dendrogram

# Perform agglomerative clustering (bottom-up) on the data
agg = AgglomerativeClustering(n_clusters=None, distance_threshold=0, linkage='ward')
# distance_threshold=0 gives the full tree without cutting (we can cut manually)
agg.fit(X)

print(f"Number of merge steps: {agg.n_clusters_ - 1}")  # should equal number of points - 1
# Create a dendrogram using SciPy for visualization (optional)
Z = linkage(X, method='ward')
# Normally, you would plot the dendrogram. Here we'll just compute cluster labels for a chosen cut:
clusters_3 = AgglomerativeClustering(n_clusters=3, linkage='ward').fit_predict(X)
print(f"Labels with 3 clusters: {np.unique(clusters_3)}")
print(f"Cluster sizes for 3 clusters: {np.bincount(clusters_3)}")

DBSCAN (Clustering Spaziale Basato sulla Densità delle Applicazioni con Rumore)

DBSCAN è un algoritmo di clustering basato sulla densità che raggruppa insieme punti che sono strettamente vicini, mentre contrassegna i punti nelle regioni a bassa densità come outlier. È particolarmente utile per dataset con densità variabili e forme non sferiche.

DBSCAN funziona definendo due parametri:

• Epsilon (ε): La distanza massima tra due punti per essere considerati parte dello stesso cluster.
• MinPts: Il numero minimo di punti richiesti per formare una regione densa (punto centrale).

DBSCAN identifica punti centrali, punti di confine e punti di rumore:

• Punto Centrale: Un punto con almeno MinPts vicini entro una distanza ε.
• Punto di Confine: Un punto che si trova entro una distanza ε da un punto centrale ma ha meno di MinPts vicini.
• Punto di Rumore: Un punto che non è né un punto centrale né un punto di confine.

Il clustering procede scegliendo un punto centrale non visitato, contrassegnandolo come un nuovo cluster, quindi aggiungendo ricorsivamente tutti i punti raggiungibili per densità da esso (punti centrali e i loro vicini, ecc.). I punti di confine vengono aggiunti al cluster di un punto centrale vicino. Dopo aver espanso tutti i punti raggiungibili, DBSCAN passa a un altro punto centrale non visitato per avviare un nuovo cluster. I punti non raggiunti da alcun punto centrale rimangono etichettati come rumore.

Assunzioni e Limitazioni

Assunzioni e Punti di Forza:: DBSCAN non assume cluster sferici – può trovare cluster di forma arbitraria (anche cluster a catena o adiacenti). Determina automaticamente il numero di cluster in base alla densità dei dati e può identificare efficacemente gli outlier come rumore. Questo lo rende potente per dati reali con forme irregolari e rumore. È robusto agli outlier (a differenza di K-Means, che li costringe nei cluster). Funziona bene quando i cluster hanno densità approssimativamente uniforme.

Limitazioni: Le prestazioni di DBSCAN dipendono dalla scelta di valori appropriati per ε e MinPts. Potrebbe avere difficoltà con dati che presentano densità variabili – un singolo ε non può adattarsi sia a cluster densi che sparsi. Se ε è troppo piccolo, etichetta la maggior parte dei punti come rumore; se è troppo grande, i cluster potrebbero fondersi in modo errato. Inoltre, DBSCAN può essere inefficiente su dataset molto grandi (naivamente $O(n^2)$, anche se l'indicizzazione spaziale può aiutare). Negli spazi delle caratteristiche ad alta dimensione, il concetto di "distanza entro ε" potrebbe diventare meno significativo (la maledizione della dimensionalità), e DBSCAN potrebbe necessitare di una sintonizzazione attenta dei parametri o potrebbe non riuscire a trovare cluster intuitivi. Nonostante ciò, estensioni come HDBSCAN affrontano alcuni problemi (come la densità variabile).

from sklearn.cluster import DBSCAN

# Generate synthetic data: 2 normal clusters and 5 outlier points
cluster1 = rng.normal(loc=[100, 1000], scale=[5, 100], size=(100, 2))
cluster2 = rng.normal(loc=[120, 2000], scale=[5, 100], size=(100, 2))
outliers = rng.uniform(low=[50, 50], high=[180, 3000], size=(5, 2))  # scattered anomalies
data = np.vstack([cluster1, cluster2, outliers])

# Run DBSCAN with chosen eps and MinPts
eps = 15.0   # radius for neighborhood
min_pts = 5  # minimum neighbors to form a dense region
db = DBSCAN(eps=eps, min_samples=min_pts).fit(data)
labels = db.labels_  # cluster labels (-1 for noise)

# Analyze clusters and noise
num_clusters = len(set(labels) - {-1})
num_noise = np.sum(labels == -1)
print(f"DBSCAN found {num_clusters} clusters and {num_noise} noise points")
print("Cluster labels for first 10 points:", labels[:10])

Analisi delle Componenti Principali (PCA)

La PCA è una tecnica per la riduzione della dimensionalità che trova un nuovo insieme di assi ortogonali (componenti principali) che catturano la massima varianza nei dati. In termini semplici, la PCA ruota e proietta i dati su un nuovo sistema di coordinate in modo che la prima componente principale (PC1) spieghi la massima varianza possibile, la seconda PC (PC2) spieghi la massima varianza ortogonale a PC1, e così via. Matematicamente, la PCA calcola gli autovettori della matrice di covarianza dei dati – questi autovettori sono le direzioni delle componenti principali, e i corrispondenti autovalori indicano la quantità di varianza spiegata da ciascuno. È spesso utilizzata per l'estrazione delle caratteristiche, la visualizzazione e la riduzione del rumore.

Nota che questo è utile se le dimensioni del dataset contengono dipendenze o correlazioni lineari significative.

La PCA funziona identificando le componenti principali dei dati, che sono le direzioni di massima varianza. I passaggi coinvolti nella PCA sono:

1. Standardizzazione: Centrare i dati sottraendo la media e scalando a varianza unitaria.
2. Matrice di Covarianza: Calcolare la matrice di covarianza dei dati standardizzati per comprendere le relazioni tra le caratteristiche.
3. Decomposizione degli Autovalori: Eseguire la decomposizione degli autovalori sulla matrice di covarianza per ottenere gli autovalori e gli autovettori.
4. Selezionare le Componenti Principali: Ordinare gli autovalori in ordine decrescente e selezionare i primi K autovettori corrispondenti ai più grandi autovalori. Questi autovettori formano il nuovo spazio delle caratteristiche.
5. Trasformare i Dati: Proiettare i dati originali sul nuovo spazio delle caratteristiche utilizzando le componenti principali selezionate. La PCA è ampiamente utilizzata per la visualizzazione dei dati, la riduzione del rumore e come passo di preprocessing per altri algoritmi di machine learning. Aiuta a ridurre la dimensionalità dei dati mantenendo la sua struttura essenziale.

Autovalori e Autovettori

Un autovalore è uno scalare che indica la quantità di varianza catturata dal suo corrispondente autovettore. Un autovettore rappresenta una direzione nello spazio delle caratteristiche lungo la quale i dati variano di più.

Immagina che A sia una matrice quadrata e v sia un vettore non nullo tale che: A * v = λ * v dove:

• A è una matrice quadrata come [ [1, 2], [2, 1]] (ad esempio, matrice di covarianza)
• v è un autovettore (ad esempio, [1, 1])

Allora, A * v = [ [1, 2], [2, 1]] * [1, 1] = [3, 3] che sarà l'autovalore λ moltiplicato per l'autovettore v, rendendo l'autovalore λ = 3.

Autovalori e Autovettori nella PCA

Spieghiamo questo con un esempio. Immagina di avere un dataset con molte immagini in scala di grigi di volti di 100x100 pixel. Ogni pixel può essere considerato una caratteristica, quindi hai 10.000 caratteristiche per immagine (o un vettore di 10000 componenti per immagine). Se vuoi ridurre la dimensionalità di questo dataset utilizzando la PCA, seguiresti questi passaggi:

1. Standardizzazione: Centrare i dati sottraendo la media di ciascuna caratteristica (pixel) dal dataset.
2. Matrice di Covarianza: Calcolare la matrice di covarianza dei dati standardizzati, che cattura come le caratteristiche (pixel) variano insieme.

• Nota che la covarianza tra due variabili (pixel in questo caso) indica quanto cambiano insieme, quindi l'idea qui è scoprire quali pixel tendono ad aumentare o diminuire insieme con una relazione lineare.
• Ad esempio, se il pixel 1 e il pixel 2 tendono ad aumentare insieme, la covarianza tra di loro sarà positiva.
• La matrice di covarianza sarà una matrice 10.000x10.000 in cui ogni voce rappresenta la covarianza tra due pixel.

3. Risolvi l'equazione degli autovalori: L'equazione degli autovalori da risolvere è C * v = λ * v dove C è la matrice di covarianza, v è l'autovettore e λ è l'autovalore. Può essere risolta utilizzando metodi come:

• Decomposizione degli Autovalori: Eseguire la decomposizione degli autovalori sulla matrice di covarianza per ottenere gli autovalori e gli autovettori.
• Decomposizione ai Valori Singolari (SVD): In alternativa, puoi utilizzare SVD per decomporre la matrice dei dati in valori e vettori singolari, che possono anche fornire le componenti principali.

4. Selezionare le Componenti Principali: Ordinare gli autovalori in ordine decrescente e selezionare i primi K autovettori corrispondenti ai più grandi autovalori. Questi autovettori rappresentano le direzioni di massima varianza nei dati.

Assunzioni e Limitazioni

La PCA assume che gli assi principali di varianza siano significativi – è un metodo lineare, quindi cattura correlazioni lineari nei dati. È non supervisionato poiché utilizza solo la covarianza delle caratteristiche. I vantaggi della PCA includono la riduzione del rumore (componenti a bassa varianza spesso corrispondono a rumore) e la decorrelazione delle caratteristiche. È computazionalmente efficiente per dimensioni moderatamente elevate ed è spesso un utile passo di preprocessing per altri algoritmi (per mitigare la maledizione della dimensionalità). Una limitazione è che la PCA è limitata a relazioni lineari – non catturerà strutture complesse non lineari (mentre autoencoder o t-SNE potrebbero). Inoltre, le componenti PCA possono essere difficili da interpretare in termini di caratteristiche originali (sono combinazioni di caratteristiche originali). Nella cybersecurity, bisogna essere cauti: un attacco che causa solo un cambiamento sottile in una caratteristica a bassa varianza potrebbe non apparire nelle prime PC (poiché la PCA dà priorità alla varianza, non necessariamente all'"interessantezza").

from sklearn.decomposition import PCA

# Create synthetic 4D data (3 clusters similar to before, but add correlated features)
# Base features: duration, bytes (as before)
base_data = np.vstack([normal1, normal2, normal3])  # 1500 points from earlier normal clusters
# Add two more features correlated with existing ones, e.g. packets = bytes/50 + noise, errors = duration/10 + noise
packets = base_data[:, 1] / 50 + rng.normal(scale=0.5, size=len(base_data))
errors = base_data[:, 0] / 10 + rng.normal(scale=0.5, size=len(base_data))
data_4d = np.column_stack([base_data[:, 0], base_data[:, 1], packets, errors])

# Apply PCA to reduce 4D data to 2D
pca = PCA(n_components=2)
data_2d = pca.fit_transform(data_4d)
print("Explained variance ratio of 2 components:", pca.explained_variance_ratio_)
print("Original shape:", data_4d.shape, "Reduced shape:", data_2d.shape)
# We can examine a few transformed points
print("First 5 data points in PCA space:\n", data_2d[:5])

Modelli di Miscele Gaussiane (GMM)

Un Modello di Miscele Gaussiane assume che i dati siano generati da una miscela di diverse distribuzioni Gaussiane (normali) con parametri sconosciuti. In sostanza, è un modello di clustering probabilistico: cerca di assegnare dolcemente ogni punto a uno dei K componenti Gaussiani. Ogni componente gaussiano k ha un vettore medio (μ_k), una matrice di covarianza (Σ_k) e un peso di miscelazione (π_k) che rappresenta quanto è prevalente quel cluster. A differenza di K-Means che fa assegnazioni "dure", GMM dà a ogni punto una probabilità di appartenere a ciascun cluster.

L'adattamento di GMM viene tipicamente eseguito tramite l'algoritmo di Massimizzazione delle Aspettative (EM):

• Inizializzazione: Iniziare con stime iniziali per le medie, le covarianze e i coefficienti di miscelazione (o utilizzare i risultati di K-Means come punto di partenza).

• E-step (Aspettativa): Date le attuali parametri, calcolare la responsabilità di ciascun cluster per ogni punto: essenzialmente r_nk = P(z_k | x_n) dove z_k è la variabile latente che indica l'appartenenza al cluster per il punto x_n. Questo viene fatto usando il teorema di Bayes, dove calcoliamo la probabilità posteriore di ciascun punto appartenente a ciascun cluster in base ai parametri attuali. Le responsabilità vengono calcolate come:

r_{nk} = \frac{\pi_k \mathcal{N}(x_n | \mu_k, \Sigma_k)}{\sum_{j=1}^{K} \pi_j \mathcal{N}(x_n | \mu_j, \Sigma_j)}

dove:

• ( \pi_k ) è il coefficiente di miscelazione per il cluster k (probabilità a priori del cluster k),

• ( \mathcal{N}(x_n | \mu_k, \Sigma_k) ) è la funzione di densità di probabilità gaussiana per il punto ( x_n ) dato la media ( \mu_k ) e la covarianza ( \Sigma_k ).

• M-step (Massimizzazione): Aggiornare i parametri utilizzando le responsabilità calcolate nell'E-step:

• Aggiornare ogni media μ_k come la media ponderata dei punti, dove i pesi sono le responsabilità.

• Aggiornare ogni covarianza Σ_k come la covarianza ponderata dei punti assegnati al cluster k.

• Aggiornare i coefficienti di miscelazione π_k come la responsabilità media per il cluster k.

• Iterare i passi E e M fino alla convergenza (i parametri si stabilizzano o il miglioramento della verosimiglianza è al di sotto di una soglia).

Il risultato è un insieme di distribuzioni gaussiane che modellano collettivamente la distribuzione complessiva dei dati. Possiamo utilizzare il GMM adattato per raggruppare assegnando a ciascun punto il gaussiano con la massima probabilità, o mantenere le probabilità per l'incertezza. Si può anche valutare la verosimiglianza di nuovi punti per vedere se si adattano al modello (utile per la rilevazione di anomalie).

Assunzioni e Limitazioni

GMM è una generalizzazione di K-Means che incorpora la covarianza, quindi i cluster possono essere ellissoidali (non solo sferici). Gestisce cluster di dimensioni e forme diverse se la covarianza è completa. Il clustering morbido è un vantaggio quando i confini dei cluster sono sfocati – ad esempio, nella cybersecurity, un evento potrebbe avere tratti di più tipi di attacco; GMM può riflettere quell'incertezza con probabilità. GMM fornisce anche una stima di densità probabilistica dei dati, utile per rilevare outlier (punti con bassa probabilità sotto tutti i componenti della miscela).

D'altra parte, GMM richiede di specificare il numero di componenti K (anche se si possono utilizzare criteri come BIC/AIC per selezionarlo). EM può a volte convergere lentamente o a un ottimo locale, quindi l'inizializzazione è importante (spesso si esegue EM più volte). Se i dati non seguono effettivamente una miscela di gaussiane, il modello potrebbe non adattarsi bene. C'è anche il rischio che una gaussiana si restringa per coprire solo un outlier (anche se la regolarizzazione o i limiti minimi di covarianza possono mitigare ciò).

from sklearn.mixture import GaussianMixture

# Fit a GMM with 3 components to the normal traffic data
gmm = GaussianMixture(n_components=3, covariance_type='full', random_state=0)
gmm.fit(base_data)  # using the 1500 normal data points from PCA example

# Print the learned Gaussian parameters
print("GMM means:\n", gmm.means_)
print("GMM covariance matrices:\n", gmm.covariances_)

# Take a sample attack-like point and evaluate it
sample_attack = np.array([[200, 800]])  # an outlier similar to earlier attack cluster
probs = gmm.predict_proba(sample_attack)
log_likelihood = gmm.score_samples(sample_attack)
print("Cluster membership probabilities for sample attack:", probs)
print("Log-likelihood of sample attack under GMM:", log_likelihood)

from sklearn.ensemble import IsolationForest

# Combine normal and attack test data from autoencoder example
X_test_if = test_data  # (120 x 2 array with 100 normal and 20 attack points)
# Train Isolation Forest (unsupervised) on the test set itself for demo (in practice train on known normal)
iso_forest = IsolationForest(n_estimators=100, contamination=0.15, random_state=0)
iso_forest.fit(X_test_if)
# Predict anomalies (-1 for anomaly, 1 for normal)
preds = iso_forest.predict(X_test_if)
anomaly_scores = iso_forest.decision_function(X_test_if)  # the higher, the more normal
print("Isolation Forest predicted labels (first 20):", preds[:20])
print("Number of anomalies detected:", np.sum(preds == -1))
print("Example anomaly scores (lower means more anomalous):", anomaly_scores[:5])

# 1 ─────────────────────────────────────────────────────────────────────
#    Create synthetic 4-D dataset
#      • Three clusters of “normal” traffic (duration, bytes)
#      • Two correlated features: packets & errors
#      • Five outlier points to simulate suspicious traffic
# ──────────────────────────────────────────────────────────────────────
import numpy as np
import matplotlib.pyplot as plt
from sklearn.manifold import TSNE
from sklearn.preprocessing import StandardScaler

rng = np.random.RandomState(42)

# Base (duration, bytes) clusters
normal1 = rng.normal(loc=[50, 500],  scale=[10, 100], size=(500, 2))
normal2 = rng.normal(loc=[60, 1500], scale=[8,  200], size=(500, 2))
normal3 = rng.normal(loc=[70, 3000], scale=[5,  300], size=(500, 2))

base_data = np.vstack([normal1, normal2, normal3])       # (1500, 2)

# Correlated features
packets = base_data[:, 1] / 50 + rng.normal(scale=0.5, size=len(base_data))
errors  = base_data[:, 0] / 10 + rng.normal(scale=0.5, size=len(base_data))

data_4d = np.column_stack([base_data, packets, errors])  # (1500, 4)

# Outlier / attack points
outliers_4d = np.column_stack([
rng.normal(250, 1, size=5),     # extreme duration
rng.normal(1000, 1, size=5),    # moderate bytes
rng.normal(5, 1, size=5),       # very low packets
rng.normal(25, 1, size=5)       # high errors
])

data_viz = np.vstack([data_4d, outliers_4d])             # (1505, 4)

# 2 ─────────────────────────────────────────────────────────────────────
#    Standardize features (recommended for t-SNE)
# ──────────────────────────────────────────────────────────────────────
scaler = StandardScaler()
data_scaled = scaler.fit_transform(data_viz)

# 3 ─────────────────────────────────────────────────────────────────────
#    Run t-SNE to project 4-D → 2-D
# ──────────────────────────────────────────────────────────────────────
tsne = TSNE(
n_components=2,
perplexity=30,
learning_rate='auto',
init='pca',
random_state=0
)
data_2d = tsne.fit_transform(data_scaled)
print("t-SNE output shape:", data_2d.shape)  # (1505, 2)

# 4 ─────────────────────────────────────────────────────────────────────
#    Visualize: normal traffic vs. outliers
# ──────────────────────────────────────────────────────────────────────
plt.figure(figsize=(8, 6))
plt.scatter(
data_2d[:-5, 0], data_2d[:-5, 1],
label="Normal traffic",
alpha=0.6,
s=10
)
plt.scatter(
data_2d[-5:, 0], data_2d[-5:, 1],
label="Outliers / attacks",
alpha=0.9,
s=40,
marker="X",
edgecolor='k'
)

plt.title("t-SNE Projection of Synthetic Network Traffic")
plt.xlabel("t-SNE component 1")
plt.ylabel("t-SNE component 2")
plt.legend()
plt.tight_layout()
plt.show()