rTensor2"
In rTensor2: MultiLinear Algebra

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

library(rTensor)
library(rTensor2)

Introduction

rTensor2 is a set of tools to manipulate 3D tensors of data. The package follows the work of Kernfled et al. (2015) https://www.sciencedirect.com/science/article/pii/S0024379515004358. For most operations, the first step is to take a discrete linear transpose along mode 3. The rTensor2 package extends the work of the rTensor package https://www.jstatsoft.org/article/download/v087i10/1270 which also provides 3D tensor tools, but this work was based on the discrete Fourier transform. rTensor2 however, extends this package by allowing tensor manipulation to be done using any one of six discrete transforms.

General principles:

Results are transformation dependent. In other words, multiplying 2 tensors together using a different transformation will give you a different product.
Inverses are not unique. In other words, different transformations yield inverses that are scaled differently. When calculating: $\texttt{tmult}(\texttt{A},\texttt{tInv(A)},\texttt{"tform"})=\texttt{I}$, the identity tensor will have zeros for all lateral slices except the first slice and the first lateral slice will be an identity matrix times some scalar. The scalar will be different depending on which transform was used as well as the number of lateral slices.
If you calculate the inverse tensor and then use the inverse tensor in a product, you must use the same transformation that was used to find the inverse tensor. The same is true for the tensor transpose function.

Tensor Multiplication

The tensor multiplication function ($\texttt{tmult}$) multiplies two tensors. Step one is to perform a discrete transform along mode 3. After transformation, the frontal slices are multiplied to together resulting in a new tensor.

Because the frontal faces are multiplied together, the second mode of the first tensor must match the first mode of the second tensor. As a example, we can multiply a $3 \times 4 \times 2$ tensor by a $4 \times 5 \times 2$ tensor. In addition, mode 3 dimensions must also match. If you are using the discrete wavelet transform, there is an additional requirement that the mode 3 dimension must be a of length $2^n$ where $n>1$ otherwise you must 0 pad the tensor in that dimension.

A <- rand_tensor(modes=c(3,4,2))
B <- rand_tensor(modes=c(4,5,2))
tmult(x=A,y=B,"dct")

Tensor Inverse

The tensor inverse function ($\texttt{tINV}$) performs a discrete transformation down mode 3. Then, each of the frontal slices are inverted. Because the slices must be inverted, the mode 1 and mode 2 dimensions must be equal.

A <-  rand_tensor(modes=c(3,3,2))
Ainv = tINV(A,"dct")
tmult(A,Ainv,"dct") # the result is an identity tensor

Tensor Transpose

The tensor transpose ($\texttt{t_tpose}$) performs a discrete transformation down mode 3. Then, each of the frontal slices are transposed.

A <-  rand_tensor(modes=c(3,5,2))
A
t_tpose(A,"dct")

Eigenvalue Decomposition

The eigenvalue decomposition function ($\texttt{tEIG}$) performs an eigenvalue decomposition on a tensor whose dimensions are $n \times n \times k$ using any discrete transform. The function returns 2 tenors, an $n \times n \times k$ tensor of eigenvectors and a diagonal tensor of eigenvalues which is also $n \times n \times k$. Notice that if we multiply $P D P^{-1}$ we obtain the original tensor.

A = rand_tensor(modes=c(3,3,2))
result = tEIG(A,"dst")
A - tmult(tmult(result$P,result$D,"dst"),tINV(result$P,"dst"),"dst") # zero tensor

LU Decomposition

The LU decomposition factors a 3-mode tensor into a lower triangular tensor and an upper triangular tensor.

Note: In order to perform this operation, the frontal slices need to be square.

A <- rand_tensor(modes=c(3,3,2))
result <- tLU(A,"dht")
A - tmult(result$L,result$U,"dht")

QR Decomposition

The QR decomposition factors a 3-mode tensor into a left singular tensor object tensor and a right singular tensor object.

Note: In order to perform this operation, the frontal slices need to be square.

A <- rand_tensor(modes=c(3,3,2))
result <- tQR(A,"fft")
A - tmult(result$Q,result$R,"fft")

Singular Value Decomposition

As an example, we use singular value decomposition (SVD) for image compression. By factoring a tensor $\cal{A}$ into orthogonal tensors $\cal{U},\cal{S}$ and $\cal{V}$} using the SVD, we can then truncate the tensor to the first $k$ singular values and compare different transforms to the original image.

library(raster)
library(grid)
data(raytrace)

tform = "dst"

A = raytrace$boat
wSVD = tSVD(A,tform)

k = 30 # number of singular values kept
U = wSVD$U
V = wSVD$V
S = wSVD$S
tV = t_tpose(V,tform)

Uk = rand_tensor(modes = c(128, k, 128), drop = FALSE)*0
Sk = rand_tensor(modes = c(k, k, 128), drop = FALSE)*0
Vk = rand_tensor(modes = c(k, 128, 128), drop = FALSE)*0

Uk = U[,1:k,]@data

Sk = S[1:k,1:k,]@data

Vk = tV[1:k,,]@data

Uk = as.tensor(Uk)
Vk = as.tensor(Vk)
Sk = as.tensor(Sk)

X = tmult(Uk, Sk,tform)
X = tmult(X, Vk, tform)

# See how close the compressed image is to the original image
fnorm(X-A)

# feature scale
if (tform=="fft"){
  Xnew=Re(X@data)
} else {
  Xnew = X@data
}
Xnew = X@data
newX = (Xnew-min(Xnew))/(max(Xnew)-min(Xnew))

# View Images

# Compressed image
# grid.raster(newX[,45,])
# title(paste0('Compressed'))
# Original Image
# grid.raster(XT[,45,]@data)