Quick start guide

Let's install GCPDecompositions and run our first Generalized CP Tensor Decomposition!

Step 1: Install Julia

Go to https://julialang.org/downloads and install the current stable release.

To check your installation, open up Julia and try a simple calculation like 1+1:

julia> 1 + 12

More info: https://docs.julialang.org/en/v1/manual/getting-started/

Step 2: Install GCPDecompositions

GCPDecompositions can be installed using Julia's excellent builtin package manager.

julia> import Pkg; Pkg.add("GCPDecompositions")

This downloads, installs, and precompiles GCPDecompositions (and all its dependencies). Don't worry if it takes a few minutes to complete.

Tip: Interactive package management with the Pkg REPL mode

Here we used the functional API for the builtin package manager. Pkg also has a very nice interactive interface (called the Pkg REPL) that is built right into the Julia REPL!

Learn more here: https://pkgdocs.julialang.org/v1/getting-started/

Tip: Pkg environments

The package manager has excellent support for creating separate installation environments. We strongly recommend using environments to create isolated and reproducible setups.

Learn more here: https://pkgdocs.julialang.org/v1/environments/

Step 3: Run GCPDecompositions

Let's create a simple three-way (a.k.a. order-three) data tensor that has a rank-one signal plus noise!

julia> dims = (10, 20, 30)(10, 20, 30)
julia> X = ones(dims) + randn(dims); # semicolon suppresses the output

Mathematically, this data tensor can be written as follows:

\[X = \underbrace{ \mathbf{1}_{10} \circ \mathbf{1}_{20} \circ \mathbf{1}_{30} }_{\text{rank-one signal}} + N \in \mathbb{R}^{10 \times 20 \times 30} ,\]

where $\mathbf{1}_{n} \in \mathbb{R}^n$ denotes a vector of $n$ ones, $\circ$ denotes the outer product, and $N_{ijk} \overset{iid}{\sim} \mathcal{N}(0,1)$.

Now, to get a rank $r=1$ GCP decomposition simply load the package and run gcp.

julia> using GCPDecompositions
julia> r = 1 # desired rank1
julia> M = gcp(X, r)10×20×30 CPD{Float64, 3, Vector{Float64}, Matrix{Float64}} with 1 component λ weights: 1-element Vector{Float64}: … U[1] factor matrix: 10×1 Matrix{Float64}: … U[2] factor matrix: 20×1 Matrix{Float64}: … U[3] factor matrix: 30×1 Matrix{Float64}: …

This returns a CPD (short for CP Decomposition) with weights $\lambda$ and factor matrices $U_1$, $U_2$, and $U_3$. Mathematically, this is the following decomposition (read Overview to learn more):

\[M = \sum_{i=1}^{r} \lambda[i] \cdot U_1[:,i] \circ U_2[:,i] \circ U_3[:,i] \in \mathbb{R}^{10 \times 20 \times 30} .\]

We can extract each of these as follows:

julia> M.λ    # to write `λ`, type `\lambda` then hit tab1-element Vector{Float64}:
 77.54481883038551
julia> M.U[1]10×1 Matrix{Float64}: 0.30768881675700266 0.3142778735053032 0.3351005539094988 0.3100851905006007 0.30767494570938175 0.3163050811876798 0.31956523766472333 0.31379215686841844 0.3027509727507048 0.33339678572184
julia> M.U[2]20×1 Matrix{Float64}: 0.23723067173815754 0.22556383549038242 0.23262220194572214 0.21906781260082447 0.19686201437155526 0.2271274564192224 0.20022882185109153 0.23374405055577632 0.22798541795716146 0.1906459345071386 0.2312572692962771 0.21223302185050108 0.21692230384836284 0.25285499559434843 0.22838757431849863 0.24673958948765748 0.22614924145204615 0.219330894677051 0.22895453420439035 0.2076569913991357
julia> M.U[3]30×1 Matrix{Float64}: 0.18525425295352174 0.17448037006039233 0.17684811551802465 0.1798243450274773 0.17527926323098922 0.18446097344814485 0.17829151874293936 0.1897109955348724 0.19278636141530733 0.177977978549376 ⋮ 0.18016173438637337 0.15740854515997166 0.20028917499734672 0.19249450170311042 0.14759277999415193 0.17739504918590063 0.17516983047999943 0.21293013893202548 0.19789399793041665

Let's check how close the factor matrices $U_1$, $U_2$, and $U_3$ (which were estimated from the noisy data) are to the true signal components $\mathbf{1}_{10}$, $\mathbf{1}_{20}$, and $\mathbf{1}_{30}$. We use the angle between the vectors since the scale of each factor matrix isn't meaningful on its own (read Overview to learn more):

julia> using LinearAlgebra: normalize
julia> vecangle(u, v) = acos(normalize(u)'*normalize(v)) # angle between vectors in radiansvecangle (generic function with 1 method)
julia> vecangle(M.U[1][:,1], ones(10))0.03220784964792212
julia> vecangle(M.U[2][:,1], ones(20))0.06877137769482056
julia> vecangle(M.U[3][:,1], ones(30))0.0752818133930448

The decomposition does a pretty good job of extracting the signal from the noise!

The power of Generalized CP Decomposition is that we can fit CP decompositions to data using different losses (i.e., different notions of fit). By default, gcp uses the (conventional) least-squares loss, but we can easily try another! For example, to try non-negative least-squares simply run

julia> M_nonneg = gcp(X, 1; loss = GCPLosses.NonnegativeLeastSquares())10×20×30 CPD{Float64, 3, Vector{Float64}, Matrix{Float64}} with 1 component
λ weights:
1-element Vector{Float64}: …
U[1] factor matrix:
10×1 Matrix{Float64}: …
U[2] factor matrix:
20×1 Matrix{Float64}: …
U[3] factor matrix:
30×1 Matrix{Float64}: …
Congratulations!

Congratulations! You have successfully installed GCPDecompositions and run some Generalized CP decompositions!

Next steps

Ready to learn more?

  • If you are new to tensor decompositions (or to GCP decompositions in particular), check out the Overview page in the manual. Also check out some of the demos!
  • To learn about all the different loss functions you can use or about how to add your own, check out the Loss functions page in the manual.
  • To learn about different constraints you can add, check out the Constraints page in the manual.
  • To learn about different algorithms you can choose, check out the Algorithms page in the manual.

Want to understand the internals and possibly contribute? Check out the developer docs.