Introduction to mwTensor

Overview

mwTensor is an R package for Multi-Way Component Analysis (MWCA). It provides solvers for decomposing single or multiple matrices and tensors with flexible choices of matrix factorization algorithms (SVD, NMF, ICA, CX, etc.).

Key functions:

Function	Purpose
`MWCA`	Decompose a single tensor
`CoupledMWCA`	Decompose multiple coupled tensors simultaneously
`checkCoupledMWCA`	Validate parameters before optimization
`initCoupledMWCA`	Generate reproducible initial values
`refineFactor`	One-step decomposition of a factor (experimental)
`MWCAProgram`	Declarative factorization program representation

1. Single Tensor Decomposition with MWCA

MWCA applies a matrix factorization method to each mode of a single tensor.

library(mwTensor)
#> Warning: no DISPLAY variable so Tk is not available

# Create a 3-way tensor (10 x 12 x 8)
set.seed(123)
X <- array(runif(10 * 12 * 8), dim = c(10, 12, 8))

# Default parameters: SVD on each mode, lower dim = 2
params <- defaultMWCAParams(X)
params@algorithms  # decomposition method per mode
#> [1] "mySVD" "mySVD" "mySVD"
params@dims        # target lower dimension per mode
#> [1] 2 2 2

out <- MWCA(params)

# Factor matrices: one per mode
# Each is dims[i] x dim(X)[i]
lapply(out@factors, dim)
#> [[1]]
#> [1]  2 10
#> 
#> [[2]]
#> [1]  2 12
#> 
#> [[3]]
#> [1] 2 8

# Core tensor
dim(out@core)
#> [1] 2 2 2

# Reconstruction error
out@rec_error
#> [1] 8.776407

Choosing different algorithms per mode

params2 <- new("MWCAParams",
    X = X,
    mask = NULL,
    pseudocount = 1e-10,
    algorithms = c("mySVD", "myNMF", "myICA"),
    dims = c(3L, 4L, 2L),
    transpose = FALSE,
    viz = FALSE,
    figdir = NULL)
out2 <- MWCA(params2)
lapply(out2@factors, dim)
#> [[1]]
#> [1]  3 10
#> 
#> [[2]]
#> [1]  4 12
#> 
#> [[3]]
#> [1] 2 8

Available algorithms: mySVD, myALS_SVD, myNMF, myICA, myCX.

2. Coupled Tensor Decomposition with CoupledMWCA

CoupledMWCA decomposes multiple tensors that share dimensions. Shared dimensions get common factor matrices; unshared dimensions get specific factor matrices.

Toy data

Xs <- toyModel("coupled_CP_Easy")
# X1: 20x30 matrix, X2: 30x30x30 tensor, X3: 30x25 matrix
lapply(Xs, dim)
#> $X1
#> [1] 20 30
#> 
#> $X2
#> [1] 30 30 30
#> 
#> $X3
#> [1] 30 25

The coupling structure:

X1 has dimensions I1 (20) and I2 (30)
X2 has dimensions I2 (30), I3 (30), I4 (30)
X3 has dimensions I4 (30) and I5 (25)

I2 is shared between X1 and X2. I4 is shared between X2 and X3.

Setting up the model

# common_model maps: block -> (mode_name -> factor_name)
common_model <- list(
    X1 = list(I1 = "A1", I2 = "A2"),
    X2 = list(I2 = "A2", I3 = "A3", I4 = "A4"),
    X3 = list(I4 = "A4", I5 = "A5"))

# A2 appears in both X1 and X2 -> shared factor for dimension I2
# A4 appears in both X2 and X3 -> shared factor for dimension I4

Running with default parameters

params <- defaultCoupledMWCAParams(Xs, common_model)

# Inspect defaults
params@common_algorithms  # "mySVD" for all
#> $A1
#> [1] "mySVD"
#> 
#> $A2
#> [1] "mySVD"
#> 
#> $A3
#> [1] "mySVD"
#> 
#> $A4
#> [1] "mySVD"
#> 
#> $A5
#> [1] "mySVD"
params@common_dims        # 2 for all
#> $A1
#> [1] 2
#> 
#> $A2
#> [1] 2
#> 
#> $A3
#> [1] 2
#> 
#> $A4
#> [1] 2
#> 
#> $A5
#> [1] 2
params@common_coretype    # "Tucker"
#> [1] "Tucker"

out <- CoupledMWCA(params)

# Common factor matrices
lapply(out@common_factors, dim)
#> $A1
#> [1]  2 20
#> 
#> $A2
#> [1]  2 30
#> 
#> $A3
#> [1]  2 30
#> 
#> $A4
#> [1]  2 30
#> 
#> $A5
#> [1]  2 25

# Convergence
tail(out@rec_error)
#>       95       96       97       98       99      100 
#> 7913.081 7913.081 7913.081 7913.081 7913.081 7913.081

Customizing algorithms and dimensions

params@common_algorithms <- list(
    A1 = "mySVD", A2 = "myNMF", A3 = "myNMF",
    A4 = "mySVD", A5 = "myCX")
params@common_dims <- list(
    A1 = 3, A2 = 3, A3 = 5, A4 = 4, A5 = 4)
params@common_iteration <- list(
    A1 = 5, A2 = 5, A3 = 5, A4 = 5, A5 = 5)

out2 <- CoupledMWCA(params)
lapply(out2@common_factors, dim)
#> $A1
#> [1]  3 20
#> 
#> $A2
#> [1]  3 30
#> 
#> $A3
#> [1]  5 30
#> 
#> $A4
#> [1]  4 30
#> 
#> $A5
#> [1]  4 25

3. Validating Parameters with checkCoupledMWCA

Before running a potentially expensive optimization, validate the parameter object. Unlike CoupledMWCA which stops on the first error, checkCoupledMWCA collects all errors and warnings.

result <- checkCoupledMWCA(params)
result$ok
#> [1] TRUE
result$summary
#> [1] "Validation passed"

# Introduce multiple errors
bad_params <- params
bad_params@common_algorithms$A1 <- "nonExistentAlgo"
bad_params@common_iteration$A2 <- 1.5  # must be integer
bad_params@common_coretype <- "TUCKER"  # must be "Tucker" or "CP"

result <- checkCoupledMWCA(bad_params)
result$ok
#> [1] FALSE
result$errors
#> [1] "nonExistentAlgo is not defined in .GlovalEnv"                  
#> [2] "params@common_iteration must be specified as an integer vector"
#> [3] "params@common_coretype must be 'Tucker' or 'CP'"
result$summary
#> [1] "Validation failed: 3 error(s), 0 warning(s)"

4. Reproducible Initialization with initCoupledMWCA

initCoupledMWCA generates initial factor matrices without running optimization. This is useful for:

Ensuring reproducibility via seeds
Comparing different initialization strategies
Inspecting initial values before committing to a run

params_ok <- defaultCoupledMWCAParams(Xs, common_model)

# Random initialization with seed
init <- initCoupledMWCA(params_ok, seed = 42L)
init@init_policy
#> [1] "random"
lapply(init@common_factors, dim)
#> $A1
#> [1]  2 20
#> 
#> $A2
#> [1]  2 30
#> 
#> $A3
#> [1]  2 30
#> 
#> $A4
#> [1]  2 30
#> 
#> $A5
#> [1]  2 25

Initialization policies

# SVD-based initialization
init_svd <- initCoupledMWCA(params_ok, seed = 42L, init_policy = "svd")

# Non-negative random (safe for NMF)
init_nn <- initCoupledMWCA(params_ok, seed = 42L, init_policy = "nonneg_random")

# Verify non-negativity
all(init_nn@common_factors$A1 >= 0)
#> [1] TRUE

Using initialized values for optimization

A CoupledMWCAInit object can be passed directly to CoupledMWCA. The pre-computed factors are used as the starting point.

# Initialize, then optimize
init <- initCoupledMWCA(params_ok, seed = 42L, init_policy = "svd")
out <- CoupledMWCA(init)

# Same seed -> same result
out_a <- CoupledMWCA(initCoupledMWCA(params_ok, seed = 123L))
out_b <- CoupledMWCA(initCoupledMWCA(params_ok, seed = 123L))
identical(out_a@common_factors$A1, out_b@common_factors$A1)
#> [1] TRUE

5. Factor Refinement with refineFactor (Experimental)

refineFactor applies a one-step matrix factorization to a factor obtained from a previous decomposition. This is a proof-of-concept for “decomposition of a decomposition.”

Mathematical interpretation

Given a factor matrix A (k x n), the refinement decomposes t(A) (n x k):

t(A) ~ U * V

where U is (n x dim) and V is (dim x k). This gives:

sub_factors = t(U): new lower-rank basis vectors (dim x n)
coef = t(V): coefficients (k x dim)
A ~ coef %*% sub_factors

set.seed(42)
X <- matrix(runif(20 * 30), nrow = 20, ncol = 30)
params_mat <- defaultMWCAParams(X)
params_mat@dims <- c(5L, 5L)
fit <- MWCA(params_mat)

# Factor 1 is 5 x 20. Refine it down to dim=2.
ref <- refineFactor(fit, 1L, algorithm = "mySVD", dim = 2L)
dim(ref@sub_factors)  # 2 x 20
#> [1]  2 20
dim(ref@coef)         # 5 x 2
#> [1] 5 2

# Approximation quality
max(abs(ref@source_factor - ref@coef %*% ref@sub_factors))
#> [1] 0.5069393

Refinement on CoupledMWCA results

params_c <- defaultCoupledMWCAParams(Xs, common_model)
params_c@common_dims <- list(A1=3, A2=3, A3=3, A4=3, A5=3)
fit_c <- CoupledMWCA(params_c)

# Refine common factor A2 (3 x 30) with NMF
ref_A2 <- refineFactor(fit_c, "A2", algorithm = "myNMF", dim = 2L)
dim(ref_A2@sub_factors)  # 2 x 30
#> [1]  2 30
dim(ref_A2@coef)         # 3 x 2
#> [1] 3 2

6. Factorization Programs with MWCAProgram (Experimental)

MWCAProgram provides a declarative way to describe factorization problems, including one-step factor refinements. This is an internal representation layer designed for future recursive decomposition support.

Defining a program

prog <- MWCAProgram(
    blocks = list(
        X1 = MWCAProgramBlock(
            modes = c("I1", "I2"),
            factor_map = c(I1 = "A1", I2 = "A2")),
        X2 = MWCAProgramBlock(
            modes = c("I2", "I3", "I4"),
            factor_map = c(I2 = "A2", I3 = "A3", I4 = "A4")),
        X3 = MWCAProgramBlock(
            modes = c("I4", "I5"),
            factor_map = c(I4 = "A4", I5 = "A5"))),
    factors = list(
        A1 = MWCAProgramFactor(mode = "I1", dim = 3),
        A2 = MWCAProgramFactor(mode = "I2", dim = 3),
        A3 = MWCAProgramFactor(mode = "I3", dim = 5, algorithm = "myNMF"),
        A4 = MWCAProgramFactor(mode = "I4", dim = 4),
        A5 = MWCAProgramFactor(mode = "I5", dim = 4)))
print(prog)
#> MWCAProgram
#>   Blocks:      3
#>   Factors:     5
#>   Refinements: 0
#>   Block details:
#>     X1 [common]: modes=I1,I2 -> factors={A1, A2}
#>     X2 [common]: modes=I2,I3,I4 -> factors={A2, A3, A4}
#>     X3 [common]: modes=I4,I5 -> factors={A4, A5}
#>   Factor details:
#>     A1 [common]: mode=I1, dim=3, status=decomposed, algo=mySVD
#>     A2 [common]: mode=I2, dim=3, status=decomposed, algo=mySVD
#>     A3 [common]: mode=I3, dim=5, status=decomposed, algo=myNMF
#>     A4 [common]: mode=I4, dim=4, status=decomposed, algo=mySVD
#>     A5 [common]: mode=I5, dim=4, status=decomposed, algo=mySVD

Validation

v <- validateMWCAProgram(prog)
v$ok
#> [1] TRUE
v$summary
#> [1] "Program is valid"

Compiling to solver parameters

compiled <- compileMWCAProgram(prog, Xs)
is(compiled, "CoupledMWCAParams")
#> [1] TRUE
compiled@common_dims
#> $A1
#> [1] 3
#> 
#> $A2
#> [1] 3
#> 
#> $A3
#> [1] 5
#> 
#> $A4
#> [1] 4
#> 
#> $A5
#> [1] 4
compiled@common_algorithms$A3  # "myNMF" as specified
#> [1] "myNMF"

Programs with refinement

prog_ref <- MWCAProgram(
    blocks = prog$blocks,
    factors = prog$factors,
    refinements = list(
        R1 = MWCAProgramRefinement(
            source_factor = "A2",
            algorithm = "mySVD",
            dim = 2)))

result <- executeMWCAProgram(prog_ref, Xs)
is(result$fit, "CoupledMWCAResult")
#> [1] TRUE
length(result$refinements)  # 1
#> [1] 1
dim(result$refinements$R1@sub_factors)
#> [1]  2 30

Factor status: decomposed, fixed, frozen

prog_status <- MWCAProgram(
    blocks = list(
        X1 = MWCAProgramBlock(
            modes = c("I1", "I2"),
            factor_map = c(I1 = "A1", I2 = "A2")),
        X2 = MWCAProgramBlock(
            modes = c("I2", "I3"),
            factor_map = c(I2 = "A2", I3 = "A3"))),
    factors = list(
        A1 = MWCAProgramFactor(mode = "I1", dim = 3, status = "decomposed"),
        A2 = MWCAProgramFactor(mode = "I2", dim = 3, status = "fixed"),
        A3 = MWCAProgramFactor(mode = "I3", dim = 3, status = "frozen")))
# decomposed: updated by solver
# fixed:      initialized but not updated (fix=TRUE)
# frozen:     identity-like, not decomposed (decomp=FALSE)
validateMWCAProgram(prog_status)$ok
#> [1] TRUE

Session Info

sessionInfo()
#> R version 4.4.3 (2025-02-28)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Rocky Linux 9.5 (Blue Onyx)
#> 
#> Matrix products: default
#> BLAS:   /opt/R/4.4.3/lib64/R/lib/libRblas.so 
#> LAPACK: /opt/R/4.4.3/lib64/R/lib/libRlapack.so;  LAPACK version 3.12.0
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=C              
#>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
#>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
#>  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
#> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
#> 
#> time zone: Asia/Tokyo
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] mwTensor_1.2.2
#> 
#> loaded via a namespace (and not attached):
#>  [1] dotCall64_1.2       gtable_0.3.6        ellipse_0.5.0      
#>  [4] spam_2.11-3         xfun_0.57           bslib_0.10.0       
#>  [7] ggplot2_4.0.3       tagcloud_0.7.0      ggrepel_0.9.8      
#> [10] lattice_0.22-6      vctrs_0.7.3         tools_4.4.3        
#> [13] generics_0.1.4      parallel_4.4.3      tibble_3.3.1       
#> [16] rARPACK_0.11-0      pkgconfig_2.0.3     Matrix_1.7-5       
#> [19] RColorBrewer_1.1-3  S7_0.2.2            mixOmics_6.30.0    
#> [22] groupICA_0.1.1      lifecycle_1.0.5     compiler_4.4.3     
#> [25] farver_2.1.2        stringr_1.6.0       fields_17.3        
#> [28] iTensor_1.0.2       codetools_0.2-20    misc3d_0.9-2       
#> [31] einsum_0.1.2        htmltools_0.5.9     maps_3.4.3         
#> [34] sass_0.4.10         yaml_2.3.12         pillar_1.11.1      
#> [37] jquerylib_0.1.4     tidyr_1.3.2         MASS_7.3-65        
#> [40] BiocParallel_1.40.2 cachem_1.1.0        nlme_3.1-167       
#> [43] RSpectra_0.16-2     nnTensor_1.3.0      tidyselect_1.2.1   
#> [46] digest_0.6.39       stringi_1.8.7       dplyr_1.2.1        
#> [49] reshape2_1.4.4      purrr_1.2.2         splines_4.4.3      
#> [52] fastmap_1.2.0       grid_4.4.3          cli_3.6.6          
#> [55] rTensor_1.4.9       ccTensor_1.0.3      magrittr_2.0.5     
#> [58] corpcor_1.6.10      scales_1.4.0        plot3D_1.4.2       
#> [61] rmarkdown_2.31      matrixStats_1.4.1   igraph_2.1.4       
#> [64] otel_0.2.0          gridExtra_2.3       evaluate_1.0.5     
#> [67] knitr_1.51          tcltk_4.4.3         viridisLite_0.4.3  
#> [70] mgcv_1.9-1          jointDiag_0.4       rlang_1.2.0        
#> [73] Rcpp_1.1.1-1.1      glue_1.8.1          geigen_2.3         
#> [76] jsonlite_2.0.0      R6_2.6.1            plyr_1.8.9