Introduzione

Programma: tensori in generale, canonical polyadic decomposition (CP), tucker, tensor train, …

Testi consigliati: Tensor Analysis for Data Science


Lezione 25/02/2026

Cosa sono i tensori e dove trovarli

§ Definizione (Tensore). Un tensore x di ordine d è un array d-dimensionale con entrate in un campo \mathbb{K}.

x \in \mathbb{K}^{n_1 \times \dots \times n_d}

Le sue entrate le indichiamo con x(i_1, \dots, i_d) con i_j \in \{1, \dots, n_j\}. Per d \le 3 usiamo anche la notazione x_{i_1, i_2, i_3} (o simile).

Possiamo osservare che per d=1 e d=2 i tensori coincidono rispettivamente con i vettori e le matrici.

In altri contesti i tensori vengono definiti come oggetti multilineari, ma in questo corso ci concentreremo sulla definizione data sopra ed in generale li tratteremo come oggetti più di “basso livello” rispetto al punto di vista algebrico.

Giusto per completezza, un modo alternativo di definire i tensori è il seguente:

Esempio. Dati V_1, \dots, V_d spazi vettoriali su \mathbb{K}, prendiamo \Phi: V_1 \times \dots \times V_{d-1} \to V_d multi-lineare. Supponiamo di fissare la base \{e_1^{(j)}, \dots, e_{n_j}^{(j)}\} per ogni V_j. Possiamo rappresentare \Phi con il tensore x \in \mathbb{K}^{n_1 \times \dots \times n_d}:

x(i_1, \dots, i_d) = ``\text{componente } i_d \text{ di } \Phi(e_{i_1}^{(1)}, \dots, e_{i_{d-1}}^{(d-1)})"

Si può anche fare usando \Phi: V_1 \times \dots \times V_d \to \mathbb{K} multilineare:

x(i_1, \dots, i_d) = \Phi(e_{i_1}^{(1)}, \dots, e_{i_d}^{(d)})

Nota. Da ora in poi, se non specificato, \mathbb{K} = \mathbb{R}.

Un altro caso in cui in realtà compaiono naturalmente i tensori nello sviluppo di Taylor in più variabili. Ricapitomo brevemente:

Esempio (Sviluppo di Taylor in più variabili). Sia f: \mathbb{R}^d \to \mathbb{R}, x^{(0)} \in \mathbb{R}^d. Sia i = (i_1, \dots, i_d) un multi-indice ed indichiamo:

  • |i| = \sum_{j=1}^d i_j

  • i! = \prod_{j=1}^d i_j!

Lo sviluppo di Taylor p_k(x) di ordine k (attorno a x^{(0)}) è:

p_k(x) = \sum_{|i| \le k} \frac{D^i f(x^{(0)})}{i!} \cdot (x - x^{(0)})^i

ponendo (x - x^{(0)})^i \coloneqq \prod_{j=1}^d (x_j - x_j^{(0)})^{i_j} e le derivate parziali sono espresse come:

D^i f(x^{(0)}) = \frac{\partial^{|i|}}{\partial x_1^{i_1} \dots \partial x_d^{i_d}} f(x^{(0)})

Per k=2, lo sviluppo si riconduce alla forma classica:

p_2(x) = f(x^{(0)}) + \nabla f(x^{(0)})^T (x - x^{(0)}) + \frac{1}{2} (x - x^{(0)})^T H f(x^{(0)}) (x - x^{(0)})

Per generalizzare questo discorso a ordini superiori (k \ge 3, vedremo meglio più avanti quando introdurremo le operazioni sui tensori), è necessario utilizzare i tensori delle derivate parziali:

x(i_1, \dots, i_d) = \frac{\partial^d f}{\partial x_{i_1} \dots \partial x_{i_d}}

Esempio (Equazioni differenziali su griglie). Consideriamo l’equazione di Poisson in d dimensioni:

\left \{ \begin{aligned}
\Delta u &= f && x \in [0, 1]^d \\
u &\equiv 0 && x \in \partial \Omega
\end{aligned} \right.

Discretizzando con le differenze finite su una griglia di passo h = 1/(n+1) data dai punti (i_1 h, \dots, i_d h), l’incognita del problema diventa un vettore di dimensione n^d. Vedremo che in questo contesto è estremamente vantaggioso considerare l’incognita come un tensore:

x \in \mathbb{R}^{n \times \dots \times n}

Esempio (Rappresentazione in basi di funzioni). In generale, se cerchiamo una rappresentazione per una funzione f: \mathbb{R}^d \to \mathbb{R} in una base di funzioni ottenuta come prodotto di basi 1D \{g_j^{(i)}(x_i)\}, possiamo rappresentare f tramite un tensore x \in \mathbb{R}^{n_1 \times \dots \times n_d}:

f(x_1, \dots, x_d) = \sum_{i_1, \dots, i_d} x(i_1, \dots, i_d) \cdot g_{i_1}^{(1)}(x_1) \cdots g_{i_d}^{(d)}(x_d)

Dove g_j^{(i)}(x) è l’i-esimo elemento della base per la variabile x_j.

Osservazione. Una fonte comune di problemi che coinvolgono tensori sono quelli in cui si cerca di approssimare una funzione incognita dove il tensore è definito implicitamente dal problema e (di solito) sarebbe troppo grande per essere rappresentato in maniera esplicita.

Esempio (Tensori multicaratteristica). Tensori che contengono dati multi-caratteristica dove abbiamo accesso esplicito al tensore. Ad esempio, un video di 1000 frame da 640 \times 480 pixel (in scala di grigi) può essere interpretato come un tensore:

x \in \mathbb{R}^{1000 \times 640 \times 480}

Osservazione. Il consumo di risorse cresce estremamente velocemente con l’ordine e le dimensioni del tensore. Ad esempio, una matrice 512 \times 512 è considerata piccola (il calcolo della SVD si fa in meno di un secondo), ma un tensore 512 \times 512 \times 512 occupa già circa 1 GB di memoria.

Sotto-parti di un tensore

§ Definizione (Modi). Un tensore x \in \mathbb{R}^{n_1 \times \dots \times n_d} si dice avere d modi (modes).

§ Definizione (Fibre). I vettori ottenuti da un tensore x fissando tutti gli indici tranne uno sono detti fibre (fibers). In particolare, la fibra del modo k è data da:

x(i_1, \dots, i_{k-1}, :, i_{k+1}, \dots, i_d) \in \mathbb{R}^{n_k}

Nel caso di un tensore di ordine d=3, abbiamo:

  • k=1: column fibers
  • k=2: row fibers
  • k=3: tube fibers

§ Definizione (Fette). Le matrici ottenute fissando tutti gli indici tranne due sono dette fette (slices). Se i due indici liberi di variare sono i primi due, si parla di fette frontali:

x(:, :, i_3, \dots, i_d)

§ Definizione (Hyperslice). Se fissiamo d-k indici, il tensore di dimensione k ottenuto si dice hyperslice.

Ordinamenti per indici multidimensionali

Ci servirà definire un ordinamento su S = \{1, \dots, n_1\} \times \dots \times \{1, \dots, n_d\}, ovvero definire una funzione biunivoca da S all’insieme delle sue linearizzazioni \{1, \dots, \prod_{j=1}^d n_j\}.

§ Definizione (Ordinamento naturale). L’ordinamento naturale è definito come segue:

\mathbb{L}: S \to \left\{1, \dots, \prod_{j=1}^d n_j\right\}
(i_1, \dots, i_d) \mapsto 1 + \sum_{j=1}^d s_j (i_j - 1)

dove s_j = \prod_{h=1}^{j-1} n_h è detto passo dell’ordinamento (con s_1 = 1). La sua inversa \mathbb{L}^{-1} è data da:

\mathbb{L}^{-1}: \left\{1, \dots, \prod_{j=1}^d n_j\right\} \to S
\ell \mapsto \left( 1 + \left\lfloor \frac{(\ell - 1) \pmod{n_1 s_1}}{s_1} \right\rfloor, \dots, 1 + \left\lfloor \frac{(\ell - 1) \pmod{n_d s_d}}{s_d} \right\rfloor \right)

§ Definizione (Reverse ordering). Il reverse ordering è quello definito dalla funzione \mathbb{L}^*:

\mathbb{L}^*: S \to \left\{1, \dots, \prod_{j=1}^d n_j\right\}
(i_1, \dots, i_d) \mapsto 1 + \sum_{j=1}^d s_j^* (i_j - 1)

con s_j^* = \prod_{h=j+1}^{d} n_h (con s_d^* = 1).

Nota. Questi ordinamenti vengono anche detti linearizzazioni e le loro inverse tensorizzazioni.

Quando è necessario scrivere esplicitamente le dimensioni del tensore, scriviamo:

  • \mathbb{L}(i_1, \dots, i_d; n_1, \dots, n_d)
  • \mathbb{L}^{-1}(\ell; n_1, \dots, n_d)

Proprietà. Siano \alpha e \beta indici ottenuti dalla linearizzazione di (i_1, \dots, i_{d'}) e (j_1, \dots, j_d). Definiamo le dimensioni totali M = \prod_{h=1}^{d'} m_h e N = \prod_{h=1}^d n_h. Allora vale:

\mathbb{L}(\alpha, \beta; M, N) = \mathbb{L}((i_1, \dots, i_{d'}, j_1, \dots, j_d); m_1, \dots, m_{d'}, n_1, \dots, n_d)

Dimostrazione. Per induzione. \square

Vettorizzazioni di un tensore

§ Definizione (Vettorizzazione). L’operatore \text{vec} mappa un tensore in un vettore preservandone l’ordine naturale:

\text{vec}: \mathbb{R}^{n_1 \times \dots \times n_d} \to \mathbb{R}^{\prod n_j}
x \mapsto \text{vec}(x) \quad \text{con } (\text{vec}(x))_\ell = x(\mathbb{L}^{-1}(\ell; n_1, \dots, n_d))

Osservazione. I tensori sono memorizzati nell’elaboratore come vettori seguendo l’ordinamento naturale. Di conseguenza, l’operazione \text{vec}(\cdot) non ha alcun costo computazionale in termini di spostamento di dati. Altre forme di vettorizzazione potrebbero risultare più costose.

Matricizzazioni / tensor unwinding

§ Definizione (Modo-k unfolding). Il modo-k unfolding di x \in \mathbb{R}^{n_1 \times \dots \times n_d} è la matrice X_{(k)} \in \mathbb{R}^{n_k \times \prod_{j \neq k} n_j} le cui colonne sono le fibre del modo k, ordinate seguendo l’ordinamento naturale su \{1, \dots, n_1\} \times \dots \times \{1, \dots, n_{k-1}\} \times \{1, \dots, n_{k+1}\} \times \dots \times \{1, \dots, n_d\}.

In componenti:

X_{(k)}(i, j) = x(i_1, \dots, i_{k-1}, i, i_{k+1}, \dots, i_d)

dove j = \mathbb{L}(i_1, \dots, i_{k-1}, i_{k+1}, \dots, i_d; n_1, \dots, n_{k-1}, n_{k+1}, \dots, n_d).

§ Definizione (Reshape). L’operatore reshape lascia immutate le entrate del tensore ma lo ridimensiona cambiando il numero di modi o le loro dimensioni, a patto che il numero totale di elementi rimanga costante.

Dato x \in \mathbb{R}^{n_1 \times \dots \times n_d} e nuovi modi m_1, \dots, m_{d'} tali che \prod_{j=1}^{d'} m_j = \prod_{j=1}^d n_j, definiamo:

\gamma := \text{reshape}(x, m_1 \times \dots \times m_{d'}) \in \mathbb{R}^{m_1 \times \dots \times m_{d'}}

tale che \gamma(i_1, \dots, i_{d'}) = x(j_1, \dots, j_d) con la condizione di uguaglianza degli indici linearizzati:

\mathbb{L}(i_1, \dots, i_{d'}; m_1, \dots, m_{d'}) = \mathbb{L}(j_1, \dots, j_d; n_1, \dots, n_d)

Osservazione. Anche l’operazione di reshape ha costo nullo in termini di movimentazione dati. In particolare, si ha che:

X_{(1)} = \text{reshape}\left(x, n_1 \times \prod_{j=2}^d n_j\right)

da cui segue che \text{vec}(x) = \text{vec}(X_{(1)}), e quindi calcolare X_{(1)} non costa nulla. Strutturalmente X_{(1)} può essere visto come:

X_{(1)} = n_1 \left[ \begin{array}{c|c|c} x(:, \dots, :, 1) & \dots & x(:, \dots, :, n_d) \end{array} \right]

dove le colonne sono le fibre del primo modo.

Osservazione. Formare X_{(n)} con n > 1 ha invece un costo non banale, poiché in generale \text{vec}(X_{(n)}) \neq \text{vec}(x).

Tuttavia, per l’ultimo modo d, si osserva che X_{(d)} ha una struttura tale per cui:

X_{(d)} = n_d \left[ \begin{array}{c} \text{vec}(x(:, \dots, :, 1))^T \\ \hline \vdots \\ \hline \text{vec}(x(:, \dots, :, n_d))^T \end{array} \right]

con n_1 \cdots n_{d-1} colonne. Ne segue che:

X_{(d)}^T = \left[ \text{vec}(x(:, \dots, :, 1)) \mid \dots \mid \text{vec}(x(:, \dots, :, n_d)) \right]

\implies \text{vec}(X_{(d)}^T) = \text{vec}(x). Pertanto, operare su X_{(d)}^T non ha costi di riordinamento in memoria.

§ Definizione (Unwinding generali). Sia S = \{1, \dots, d\} l’insieme dei modi. Si scelgono due sottoinsiemi di indici \mathscr R = \{r_1, \dots, r_\delta\} e \mathscr C = \{c_1, \dots, c_{d-\delta}\} tali che \mathscr R \cup \mathscr C = S e \mathscr R \cap \mathscr C = \emptyset. Si può allora considerare la matricizzazione:

X_{(\mathscr R \times \mathscr C)} \in \mathbb{R}^{M \times N} \quad \text{con } M = \prod_{j \in \mathscr R} n_j, \quad N = \prod_{k \in \mathscr C} n_k

definita da:

X_{(\mathscr R \times \mathscr C)}(\alpha, \beta) = x(i_1, \dots, i_d)

dove \alpha = \mathbb{L}(i_{r_1}, \dots, i_{r_\delta}; n_{r_1}, \dots, n_{r_\delta}) e \beta = \mathbb{L}(i_{c_1}, \dots, i_{c_{d-\delta}}; n_{c_1}, \dots, n_{c_{d-\delta}}).

Osservazione. Il modo-k unfolding è il caso particolare di unwinding generale ottenuto con \mathscr R = \{k\} e \mathscr C = \{1, \dots, k-1, k+1, \dots, d\}. Analogamente, la vettorizzazione \text{vec}(x) corrisponde al caso con \mathscr R = \{1, \dots, d\} e \mathscr C = \emptyset.

Osservazione. Per ogni n \in \{1, \dots, d\}, l’unfolding definito da \mathscr R = \{1, \dots, n\} e \mathscr C = \{n+1, \dots, d\} gode della proprietà:

\text{vec}(X_{(\mathscr R \times \mathscr C)}) = \text{vec}(x)

Tale operazione ha quindi costo nullo in memoria. Questa proprietà sarà fondamentale per la decomposizione tensor train.


Lezione 26/02/2026

Recap della volta scorsa: vec(.) non costa nulla.

Permutazioni di un tensore

Si tratta di una generalizzazione della trasposizione per le matrici.

§ Definizione (Permutazione). Sia x \in \mathbb{R}^{n_1 \times \dots \times n_d} e \pi \in S_d (gruppo simmetrico). Il tensore permutato sarà \mathbb{P}(x, \pi) = \gamma \in \mathbb{R}^{n_{\pi_1} \times \dots \times n_{\pi_d}} tale che:

\gamma(i_{\pi_1}, \dots, i_{\pi_d}) = x(i_1, \dots, i_d)

Osservazione. In generale, applicare una permutazione provoca spostamenti di memoria e quindi ha un costo computazionale non nullo.

Il modo-k unfolding si può scrivere come una permutazione seguita da un reshape:

  1. \pi = (k, 1, \dots, k-1, k+1, \dots, d)
  2. \gamma = \mathbb{P}(x, \pi)
  3. X_{(k)} = \text{reshape}(\gamma, n_k \times \prod_{j \neq k} n_j)

Per un unfolding generale fissiamo \mathscr R = \{r_1, \dots, r_\delta\} e \mathscr C = \{c_1, \dots, c_{d-\delta}\}. Avremo:

  1. \pi = (r_1, \dots, r_\delta, c_1, \dots, c_{d-\delta})
  2. \gamma = \mathbb{P}(x, \pi)
  3. X_{(\mathscr R \times \mathscr C)} = \text{reshape}(\gamma, M \times N) con M = \prod_{i=1}^\delta n_{r_i} e N = \prod_{j=1}^{d-\delta} n_{c_j}

Esempio. In Matlab esiste il comando permute:

T = randn(25, 50, 175, 32);
Y = permute(T, [4 1 2 3]);

Provare a misurare con timeit o tic/toc i costi dei vari unfolding.

§ Definizione (Perfect-shuffle). Si chiama perfect-shuffle k-esimo la matrice di permutazione P_k (di dimensioni opportune) tale che:

P_k \text{vec}(x) = \text{vec}(X_{(k)})

Si osservi che P_1 = I.

Operazioni sui tensori

L’insieme \mathbb{R}^{n_1 \times \dots \times n_d} è uno spazio vettoriale su \mathbb{R} con le operazioni:

  • Somma: (x_1 + x_2)(i_1, \dots, i_d) \coloneqq x_1(i_1, \dots, i_d) + x_2(i_1, \dots, i_d)
  • Prodotto per scalare: (\alpha x)(i_1, \dots, i_d) \coloneqq \alpha \cdot x(i_1, \dots, i_d) per \alpha \in \mathbb{R}
  • Prodotto scalare: \langle x_1, x_2 \rangle \coloneqq \text{vec}(x_1)^T \text{vec}(x_2) = \sum_{i_1, \dots, i_d} x_1(i_1, \dots, i_d) x_2(i_1, \dots, i_d)
  • Norma: \|x\| \coloneqq \|\text{vec}(x)\|_2 = \sqrt{\langle x, x \rangle}

Prodotto esterno tra vettori

Nel caso matriciale, dati u \in \mathbb{R}^{n_1} e v \in \mathbb{R}^{n_2}, il prodotto u v^T è una matrice di rango 1 tale che a_{ij} = u_i v_j. Vogliamo generalizzare questa costruzione ai tensori.

§ Definizione (Prodotto esterno). Dati v^{(1)} \in \mathbb{R}^{n_1}, \dots, v^{(d)} \in \mathbb{R}^{n_d}, definiamo il tensore x = v^{(1)} \circ \dots \circ v^{(d)} \in \mathbb{R}^{n_1 \times \dots \times n_d} tale che:

x(i_1, \dots, i_d) = v_{i_1}^{(1)} \cdots v_{i_d}^{(d)} = \prod_{j=1}^d v_{i_j}^{(j)}

§ Proposizione. Sia x = v^{(1)} \circ \dots \circ v^{(d)}, allora valgono le seguenti proprietà:

  1. \text{vec}(x) = v^{(d)} \otimes \dots \otimes v^{(1)} (Nel caso matriciale A = u v^T \implies \text{vec}(A) = v \otimes u)

  2. X_{(k)} = v^{(k)} \cdot (v^{(d)} \otimes \dots \otimes v^{(k+1)} \otimes v^{(k-1)} \otimes \dots \otimes v^{(1)})^T Ovvero X_{(k)} è una matrice di rango 1 con una struttura di Kronecker nel generatore.

  3. X_{(\mathscr R \times \mathscr C)} = (v^{(r_\delta)} \otimes \dots \otimes v^{(r_1)}) (v^{(c_{d-\delta})} \otimes \dots \otimes v^{(c_1)})^T

Quando abbiamo questa struttura di prodotto esterno, possiamo rappresentare x memorizzando solo i vettori v^{(j)} per j=1, \dots, d. Questo porta a un risparmio enorme in termini di memoria: O(\sum n_j) \ll O(\prod n_j).

Prodotti tensore-matrice

Un tensore può essere visto come l’applicazione di una mappa multilineare \Phi: V_1 \times \dots \times V_d \to \mathbb{R}. Se consideriamo una trasformazione lineare \varphi_k: W_k \to V_k (con \dim W_k = r), possiamo chiederci come cambia la rappresentazione del tensore se applichiamo \varphi_k al modo k-esimo.

§ Definizione (Prodotto nel modo k). Un prodotto nel modo k di un tensore x \in \mathbb{R}^{n_1 \times \dots \times n_d} per una matrice U \in \mathbb{R}^{r \times n_k} è il tensore y = x \times_k U \in \mathbb{R}^{n_1 \times \dots \times r \times \dots \times n_d} definito come:

y(i_1, \dots, i_{k-1}, \alpha, i_{k+1}, \dots, i_d) \coloneqq \sum_{i_k=1}^{n_k} x(i_1, \dots, i_d) U_{\alpha, i_k}

con \alpha \in \{1, \dots, r\}. Equivalentemente, in termini di unfolding:

Y_{(k)} = U X_{(k)}

In altre parole, il prodotto nel modo k corrisponde a moltiplicare a sinistra tutte le fibre del modo k per la matrice U.

Prodotti nei modi come prodotti di Kronecker

Nota (Prodotto di Kronecker). Siano A \in \mathbb{R}^{m \times n} e B \in \mathbb{R}^{p \times q}. Il prodotto di Kronecker è definito a blocchi come:

A \otimes B \coloneqq \begin{bmatrix} a_{11} B & \dots & a_{1n} B \\ \vdots & \ddots & \vdots \\ a_{m1} B & \dots & a_{mn} B \end{bmatrix} \in \mathbb{R}^{mp \times nq}

Se c_{\alpha \beta} è l’entrata di A \otimes B, con \alpha = \mathbb{L}(k, i) e \beta = \mathbb{L}(l, j), allora c_{\alpha \beta} = a_{ij} b_{kl}. Si noti che la linearizzazione qui usata segue la convenzione per cui \mathbb{L}(i, k; m, p) = \mathbb{L}^*(k, i; p, m) = (i-1)p + k.

Vale il seguente risultato classico per le matrici: \text{vec}(AXB) = (B^T \otimes A) \text{vec}(X).

§ Proposizione. Sia x \in \mathbb{R}^{n_1 \times n_2 \times n_3} e siano U \in \mathbb{R}^{m_1 \times n_1}, V \in \mathbb{R}^{m_2 \times n_2}, W \in \mathbb{R}^{m_3 \times n_3}. Allora i seguenti fatti sono equivalenti:

  1. y = x \times_1 U \times_2 V \times_3 W
  2. \text{vec}(y) = (W \otimes V \otimes U) \text{vec}(x)
  3. Y_{(1)} = U X_{(1)} (W \otimes V)^T
  4. Y_{(2)} = V X_{(2)} (W \otimes U)^T
  5. Y_{(3)} = W X_{(3)} (V \otimes U)^T

Dimostrazione.

1 \iff 2: Facciamo vedere che:

\text{vec}(x \times_1 U) = (I_{n_3} \otimes I_{n_2} \otimes U) \text{vec}(x)

Infatti:

\text{vec}(x \times_1 U) = \text{vec}(U X_{(1)}) = \text{vec}(U X_{(1)} I_{n_2 n_3}) = (I_{n_2 n_3} \otimes U) \text{vec}(X_{(1)}) = (I_{n_3} \otimes I_{n_2} \otimes U) \text{vec}(x)

Analogamente per il secondo modo, usando la matrice di permutazione P_2:

\text{vec}(V X_{(2)}) = (I_{n_3} \otimes I_{n_1} \otimes V) \text{vec}(X_{(2)}) = (I_{n_3} \otimes I_{n_1} \otimes V) P_2 \text{vec}(x) = (I_{n_3} \otimes V \otimes I_{n_1}) \text{vec}(x)

Iterando il procedimento si ottiene:

\begin{aligned}
\text{vec}(x \times_1 U \times_2 V \times_3 W) &= (W \otimes I \otimes I) (I \otimes V \otimes I) (I \otimes I \otimes U) \text{vec}(x) \\
&= (W \otimes V \otimes U) \text{vec}(x)
\end{aligned}

2 \iff 4: Sia M = W \otimes U. Consideriamo le entrate di Y_{(2)} per \beta \in \{1, \dots, m_2\} e \hat{\ell} \in \{1, \dots, m_1 m_3\}:

Y_{(2)}(\beta, \hat{\ell}) = \sum_{j=1}^{n_2} \sum_{\ell=1}^{n_1 n_3} V(\beta, j) X_{(2)}(j, \ell) M(\hat{\ell}, \ell)

Siano i, k tali che \mathbb{L}(i, k) = \ell \implies X_{(2)}(j, \ell) = x(i, j, k). Analogamente siano \alpha, \gamma tali che \mathbb{L}(\alpha, \gamma) = \hat{\ell}. Allora M(\hat{\ell}, \ell) = U(\alpha, i) W(\gamma, k). Sostituendo otteniamo:

Y_{(2)}(\beta, \hat{\ell}) = \sum_i \sum_j \sum_k x(i, j, k) U(\alpha, i) V(\beta, j) W(\gamma, k)

che corrisponde alla definizione di x \times_1 U \times_2 V \times_3 W. \square

Nel caso di ordine d > 3 si ha:

  1. y = x \times_1 U_1 \times_2 \dots \times_d U_d
  2. \text{vec}(y) = (U_d \otimes \dots \otimes U_1) \text{vec}(x)
  3. Y_{(n)} = U_n X_{(n)} (U_d \otimes \dots \otimes \widehat{U_n} \otimes \dots \otimes U_1)^T

§ Proposizione. Valgono le seguenti proprietà:

i) x \times_i U \times_j V = x \times_j V \times_i U (per i \neq j) ii) x \times_k U \times_k V = x \times_k (VU)

Costo computazionale dei prodotti nei modi

Calcolare x \times_k U con U \in \mathbb{R}^{r \times n_k} costa come il prodotto tra una matrice r \times n_k e una n_k \times \prod_{j \neq k} n_j, ovvero O(r \prod_{j=1}^d n_j).

Se facciamo i prodotti su due modi diversi, il risultato non cambia ma può cambiare il costo. Quale conviene fare prima?

§ Proposizione. Sia x \in \mathbb{R}^{n_1 \times \dots \times n_d}, U_k \in \mathbb{R}^{m_k \times n_k}, U_\ell \in \mathbb{R}^{m_\ell \times n_\ell}. Consideriamo y = x \times_k U_k \times_\ell U_\ell. Se vale:

\frac{1}{n_k} - \frac{1}{m_k} < \frac{1}{n_\ell} - \frac{1}{m_\ell}

allora (x \times_k U_k) \times_\ell U_\ell è meno costoso dell’altro ordine.

Più in generale, per minimizzare i FLOPS totali, si scelgono i prodotti seguendo l’ordine crescente della quantità \frac{1}{n_j} - \frac{1}{m_j}.

Dimostrazione. Sia N = \prod_{j=1}^d n_j. Il costo di x \times_k U_k \times_\ell U_\ell è dato dal primo prodotto (m_k \times n_k per n_k \times N/n_k) più il secondo (m_\ell \times n_\ell per n_\ell \times (N/n_k) \cdot m_k):

O(m_k N) + O\left(m_\ell \frac{N}{n_k} m_k\right) = O\left(m_k N + m_\ell \frac{N}{n_k} m_k\right)

L’altro ordine avrà costo O(m_\ell N + m_k \frac{N}{n_\ell} m_\ell). Confrontando i due:

\begin{aligned}
m_k N + m_\ell \frac{N}{n_k} m_k &< m_\ell N + m_k \frac{N}{n_\ell} m_\ell \\
\frac{m_k}{m_k m_\ell} + \frac{m_k m_\ell}{n_k m_k m_\ell} &< \frac{m_\ell}{m_k m_\ell} + \frac{m_k m_\ell}{n_\ell m_k m_\ell} \quad (\text{dividendo per } N m_k m_\ell) \\
\frac{1}{m_\ell} + \frac{1}{n_k} &< \frac{1}{m_k} + \frac{1}{n_\ell} \\
\frac{1}{n_k} - \frac{1}{m_k} &< \frac{1}{n_\ell} - \frac{1}{m_\ell}
\end{aligned}

\square


Lezione 04/03/2026

Ottimizzazione computazionale dei prodotti

Abbiamo visto che il prodotto y = x \times_j U_j con x \in \mathbb{R}^{n_1 \times \dots \times n_d} e U_j \in \mathbb{R}^{r \times n_j} può essere calcolato tramite l’unfolding:

(x \times_j U_j)_{(j)} = U_j \cdot X_{(j)}

o tramite vettorizzazione:

\text{vec}(x \times_j U_j) = (I \otimes \dots \otimes I \otimes U_j \otimes I \otimes \dots \otimes I) \text{vec}(x)

Come evitare permutazioni nel calcolo di x \times_j U

Sappiamo che formare X_{(j)} richiede in generale delle permutazioni costose in termini di memoria, tranne in alcuni casi particolari.

Caso j=1: Formare X_{(1)} non comporta permutazioni (è un semplice reshape). Il calcolo può quindi essere ottimizzato come:

  1. X_{(1)} = \text{reshape}(x, n_1 \times \prod_{j=2}^d n_j)
  2. Y = U \cdot X_{(1)}
  3. y = \text{reshape}(Y, r \times n_2 \times \dots \times n_d)

Il costo computazionale è dominato dal prodotto di matrici: O(r \prod_{j=1}^d n_j).

Caso j=d: L’unfolding X_{(d)} richiede permutazioni, ma abbiamo osservato che X_{(d)}^T no (è la vettorizzazione delle fette frontali). Possiamo quindi scrivere:

U \cdot X_{(d)} = (X_{(d)}^T U^T)^T

Quindi il calcolo del prodotto nell’ultimo modo può essere implementato senza permutazioni esplicite del tensore originale:

  1. Si interpreta il tensore come una matrice \prod_{j=1}^{d-1} n_j \times n_d.
  2. Si moltiplica a destra per U^T.
  3. Si effettua il reshape finale.

In sintesi:

(x \times_d U)_{(d)}^T = \text{reshape}(x, \prod_{j=1}^{d-1} n_j \times n_d) \cdot U^T

mantenendo il costo O(r \prod n_j) e minimizzando gli spostamenti di memoria.

Caso j=2, d=3: Consideriamo y = x \times_2 U. In componenti:

y_{i_1, \alpha, i_3} = \sum_{i_2=1}^{n_2} x_{i_1, i_2, i_3} U_{\alpha, i_2}

L’idea è sfruttare il fatto che le fette frontali sono memorizzate in modo contiguo. Possiamo allora operare fetta per fetta:

y(:, :, i_3) = x(:, :, i_3) \cdot U^T

Algoritmicamente:

for j = 1 : n_3
    Y(:, :, j) = X(:, :, j) * U'; % (n1 x n2) * (n2 x r)
end

Caso 1 < j < d, d > 3: Ci riportiamo al caso 3D raggruppando gli indici:

  1. Siano M = \prod_{h=1}^{j-1} n_h e N = \prod_{h=j+1}^d n_h.
  2. Effettuiamo un reshape di x in un tensore di ordine 3: \gamma = \text{reshape}(x, M \times n_j \times N).
  3. Calcoliamo \gamma \times_2 U come nel caso 3D precedente.
  4. Riapplichiamo il reshape finale: y = \text{reshape}(\gamma \times_2 U, n_1 \times \dots \times n_{j-1} \times r \times n_{j+1} \times \dots \times n_d).

Fattorizzazioni: CP (rango basso)

Vogliamo ridurre i costi computazionali sfruttando in generale idee simili all’approssimazione di rango basso nel caso delle matrici.

Nel caso matriciale, X \in \mathbb{R}^{n_1 \times n_2} ha rango r se e solo se X = U V^T con U \in \mathbb{R}^{n_1 \times r} e V \in \mathbb{R}^{n_2 \times r}.

  • Memoria: O((n_1 + n_2) r) < O(n_1 n_2)
  • Accesso: O(r) > O(1)

Se guardiamo alla vettorizzazione: \text{vec}(x) = \text{vec}(U V^T) = \sum_{j=1}^r v_j \otimes u_j.

Domanda. Come possiamo generalizzare ai tensori?

Possiamo chiedere che, per x \in \mathbb{R}^{n_1 \times n_2 \times n_3}:

\text{vec}(x) = \sum_{j=1}^r w_j \otimes v_j \otimes u_j \quad \text{con } u_j \in \mathbb{R}^{n_1}, v_j \in \mathbb{R}^{n_2}, w_j \in \mathbb{R}^{n_3}

ovvero:

x = \sum_{j=1}^r u_j \circ v_j \circ w_j \iff x(i_1, i_2, i_3) = \sum_{j=1}^r (u_j)_{i_1} \cdot (v_j)_{i_2} \cdot (w_j)_{i_3}

§ Definizione (Decomposizione CP). Un tensore x \in \mathbb{R}^{n_1 \times \dots \times n_d} ammette una decomposizione CP (Canonical Polyadic) di rango r se esistono matrici U_1, \dots, U_d con U_j \in \mathbb{R}^{n_j \times r} tali che valgono le seguenti condizioni equivalenti:

  • x = \sum_{h=1}^r U_1^{(h)} \circ \dots \circ U_d^{(h)} \quad con U_j^{(h)} colonna h-esima di U_j

  • \text{vec}(x) = \sum_{h=1}^r U_d^{(h)} \otimes \dots \otimes U_1^{(h)}

  • x(i_1, \dots, i_d) = \sum_{h=1}^r U_1(i_1, h) \cdots U_d(i_d, h)

Per accedere a x(i_1, \dots, i_d) effettuiamo il prodotto entrata per entrata di d vettori (le righe delle U_j) e poi li sommiamo.

Osservazione. I costi di memoria e di accesso se abbiamo una fattorizzazione CP sono:

FormatoMemoriaAccesso
CPO((\sum n_j) r)O(d r)
FULL\prod n_jO(1)

Notazione. Se U_1, \dots, U_d danno una fattorizzazione di x, scriviamo x = \llbracket U_1, \dots, U_d \rrbracket.

§ Definizione (Rango CP). Il rango CP di un tensore x è il minimo r \in \mathbb{N} tale che x ammetta una fattorizzazione CP di rango r.

Esempio (Funzioni separabili). Consideriamo un tensore x le cui entrate sono date da una funzione f(x, y, z) su una griglia. Se f è semiseparabile, ovvero:

f(x, y, z) = \sum_{j=1}^r g_j(x) h_j(y) s_j(z)

allora una fattorizzazione CP di rango r è data da x = \llbracket G, H, S \rrbracket con G_{i \ell} = g_\ell(x_i), H_{j \ell} = h_\ell(y_j), S_{k \ell} = s_\ell(z_k).

Se f non è semiseparabile, si può cercare una \tilde{f} che approssimi f e sia semiseparabile per estrarne una fattorizzazione tale che:

\llbracket G, H, S \rrbracket \approx x

Ad esempio, una serie di Fourier troncata: f(x, y) = \sum_{n=0}^\infty c_n(y) \sin(nx) \approx \sum_{n=0}^r c_n(y) \sin(nx).

(Provare con delle funzioni polinomiali…)

Fattorizzazione CP e unfoldings

Per analizzare gli unfolding di una decomposizione CP, è utile introdurre un nuovo operatore tra matrici.

§ Definizione (Prodotto di Khatri-Rao). Date due matrici A \in \mathbb{R}^{m \times r} e B \in \mathbb{R}^{n \times r}, il prodotto di Khatri-Rao C = A \odot B \in \mathbb{R}^{(mn) \times r} è definito come il prodotto di Kronecker colonna per colonna:

C = [ a_1 \otimes b_1 \mid \dots \mid a_r \otimes b_r ]

In componenti, vale c_{k \ell} = a_{i \ell} b_{j \ell} dove k = \mathbb{L}^*(i, j) = \mathbb{L}(j, i) = (i-1)n + j, con k \in \{1, \dots, mn\}.

§ Proposizione. Se x = \llbracket U_1, \dots, U_d \rrbracket \in \mathbb{R}^{n_1 \times \dots \times n_d}, allora il suo modo-j unfolding è dato da:

X_{(j)} = U_j \cdot (U_d \odot \dots \odot \widehat{U_j} \odot \dots \odot U_1)^T \quad j=1, \dots, d

Dimostrazione. Sia U_j = [ u_1^{(j)} \mid \dots \mid u_r^{(j)} ]. Sappiamo che x = \sum_{h=1}^r u_h^{(1)} \circ \dots \circ u_h^{(d)}. Allora l’unfolding è la somma degli unfolding dei singoli termini di rango 1:

X_{(j)} = \sum_{h=1}^r (u_h^{(1)} \circ \dots \circ u_h^{(d)})_{(j)}

Utilizzando la proprietà dell’unfolding del prodotto esterno:

\begin{aligned}
X_{(j)} &= \sum_{h=1}^r u_h^{(j)} (u_h^{(d)} \otimes \dots \otimes \widehat{u_h^{(j)}} \otimes \dots \otimes u_h^{(1)})^T \\
&= U_j [ u_1^{(d)} \otimes \dots \otimes u_1^{(1)} \mid \dots \mid u_r^{(d)} \otimes \dots \otimes u_r^{(1)} ]^T \\
&= U_j [ U_d \odot \dots \odot \widehat{U_j} \odot \dots \odot U_1 ]^T
\end{aligned}

dove nel prodotto di Khatri-Rao si salta il fattore j-esimo. \square

Rango CP

Trovare il rango CP di un tensore è un problema NP-HARD e ci sono poche proprietà che valgono in generale. Una stima che si può dare è la seguente:

§ Proposizione. Sia x \in \mathbb{R}^{n_1 \times \dots \times n_d} e sia r_j = \text{rank}(X_{(j)}) per j=1, \dots, d. Allora vale:

\max_{j=1, \dots, d} (r_j) \le \text{rank}_{CP}(x) \le \min_{j=1, \dots, d} \left( \prod_{i \neq j} n_i \right)

Dimostrazione. Per quanto riguarda la minorazione, questa è vera perché avere una fattorizzazione di rango CP r implica X_{(j)} = U_j (U_d \odot \dots \odot U_1)^T. Poiché X_{(j)} è il prodotto di due matrici di cui una ha r colonne, il suo rango non può superare r. Ne segue che \text{rank}(X_{(j)}) \le r per ogni j.

Dimostriamo WLOG che \text{rank}_{CP}(X) \le \prod_{i=2}^d r_i: abbiamo che \text{VEC}(X) = \sum_{k=1}^n U_k^{(d)} \otimes \dots \otimes U_k^{(1)} e facciamo vedere che \{ U_k^{(d)} \otimes \dots \otimes U_k^{(2)} \} \subseteq \mathbb{R}^{n_d \cdots n_2} sono indipendenti. Assumiamo per assurdo che non lo siano, in particolare possiamo esprimerne uno come comb. degli altri, senza perdita di generalità supponiamo che

U_1^{(d)} \otimes \dots \otimes U_1^{(2)} = \sum_{k=2}^n c_k (U_k^{(d)} \otimes \dots \otimes U_k^{(2)})

allora

\begin{aligned}
\text{VEC}(X) &= U_1^{(d)} \otimes \dots \otimes U_1^{(2)} \otimes U_1^{(1)} + \sum_{k=2}^n U_k^{(d)} \otimes \dots \otimes U_k^{(1)} \\
&= \sum_{k=2}^n U_k^{(d)} \otimes \dots \otimes U_k^{(2)} \otimes (U_k^{(1)} + c_k U_1^{(1)})
\end{aligned}

\implies X ha rango CP \le n-1 \implies \bot (contraddizione con la minimalità di n).

\implies U_d \otimes \dots \otimes U_2 ha rango n (perché le sue colonne sono \{ U_k^{(d)} \otimes \dots \otimes U_k^{(2)} \}). Dato che U_d \odot \dots \odot U_2 è sottomatrice di U_d \otimes \dots \otimes U_2 si ha n = \text{rank}(U_d \odot \dots \odot U_2) \le \text{rank}(U_d \otimes \dots \otimes U_2) = \prod_{j=2}^d n_j. L’algoritmo si ripete per l’unfolding sul modo j. \square

Osservazione. Segue dal risultato precedente che se \text{rank}_{CP}(X) = 1 \iff \text{rank}(X_{(j)}) = 1 \ \forall j. In generale \text{rank}_{CP}(X) = r \implies \text{rank}(X_{(j)}) \le r \ \forall j.

Osservazione. Se A è di rango 1, a_{ij} \neq 0 allora:

A = \frac{A(:, j) A(i, :)}{a_{ij}}

Nel caso tensoriale (d=3): sia X \in \mathbb{R}^{n_1 \times n_2 \times n_3} di rango 1 allora:

X = \frac{\bar{x} \circ \bar{y} \circ \bar{z}}{X_{i_1, i_2, i_3}^2} \implies X_{i_1, i_2, i_3} \neq 0

con \bar{x} = X(:, i_2, i_3), \bar{y} = X(i_1, :, i_3), \bar{z} = X(i_1, i_2, :).