Common Components and Specific Weights Analysis

ComDim is the chemometric name for Common Components and Specific Weights Analysis (CCSWA; ACCPS in the French literature). It is an unsupervised method for several quantitative tables measured on the same observations. The method seeks common sample-space dimensions and, for each dimension, a specific weight or salience for every block.

Qannari, Wakeling, Courcoux and MacFie 2000 presented the sequential algorithm later used for this model, in a sensory profiling setting, as an improved estimator for the third model of their 1995 hierarchy of association-matrix models. That 2000 paper does not use the later names CCSWA or ComDim. Qannari, Courcoux and Vigneau 2001 apply the method to preference data under the name Common Components and Specific Weights Analysis. Mazerolles, Hanafi, Dufour, Bertrand and Qannari 2006 present CCSWA as a chemometric method for several tables observed on the same samples. Jouan-Rimbaud Bouveresse and Rutledge 2024 treat CCSWA and ComDim as the same classical unsupervised family.

This card teaches classical unsupervised ComDim / CCSWA. No response $\mathbf{Y}$ enters the extraction. The displayed equations follow Hanafi and Qannari 2008, who give the ACCPS loss, the sequential algorithm attributed to Qannari et al. 2000, and the equivalent maximization of the sum of squared saliences. Tchandao Mangamana, Cariou, Vigneau, Glèlè Kakaï and Qannari 2019 give the inspected unified identity card that distinguishes ComDim from MB-PCA. The 2019 Data Handling in Science and Technology chapter by Cariou, Jouan-Rimbaud Bouveresse, Qannari and Rutledge is the modern ComDim methods chapter. Its PDF body was not opened for this card. No equation is taken from that chapter's metadata.

Let there be $B$ processed blocks ${\mathbf{X}}_b\in {\mathbb{R}}^{N\times P_b}$ for $b=1,\dots ,B$. The blocks share $N$ observations in the same row order. The numbers of variables $P_b$ may differ. Classical ComDim as taught here requires that shared sample mode. It does not cover every multiblock topology.

Code verifies the common row count. It does not silently truncate, reorder, or match samples by position when that count differs.

Hanafi and Qannari work with the matrices of scalar products between observations. Each block may have a different $P_b$, but every ${\mathbf{W}}_b={\mathbf{X}}_b{\mathbf{X}}_b^{\top }$ has size $N\times N$. Common directions are therefore sought in the shared sample space.

After the selected centering, optional within-block scaling, and optional whole-block normalization,

${\mathbf{W}}_b$ is a Gram / cross-product matrix in sample space. It is not called a covariance matrix on this card. A covariance statement would require an explicit divisor and centering convention, which the ACCPS algebra does not insert into ${\mathbf{W}}_b$.

Hanafi and Qannari write the third, ACCPS, model as

$\mathbf{Q}$ has orthonormal columns ${\mathbf{q}}_a$, the common dimensions. ${\mathbf{\Lambda }}_b$ is diagonal, with the specific weights of block $b$ on those dimensions. They note that "common components" is potentially misleading language: a weight may be zero, so a dimension need not underlie every table. A common global direction is therefore not an equal contribution from all blocks.

For a unit-norm common direction $\mathbf{q}$, the specific weight of block $b$ is

Equivalently ${\lambda }_b={\Vert {\mathbf{X}}_b^{\top }\mathbf{q}\Vert }_2^2$. A block does not have one universal ComDim weight. It has a salience ${\lambda }_{ba}$ on each retained dimension. That is the historical meaning of specific weights.

Under the residual-stage convention used here, ${\lambda }_{ba}$ is the block sum of squares along ${\mathbf{q}}_a$. It is not $R^2$ unless that denominator is stated. It is not a variable loading, not a causal weight, and not unique information.

For the first unit-norm common dimension, Hanafi and Qannari minimize

with ${\lambda }_b={\mathbf{q}}^{\top }{\mathbf{W}}_b\mathbf{q}$. They show this is equivalent to

Tchandao Mangamana et al. write the same ComDim objective as maximizing ${\sum }_b{\mathrm{cov}}^2(\mathbf{t},{\mathbf{t}}_b)$ subject to $\Vert \mathbf{t}\Vert =1$, with block components ${\mathbf{t}}_b={\mathbf{W}}_b\mathbf{t}$. With that unit-norm constraint the two statements differ only by the constant convention in their covariance. The educational code maximizes ${\sum }_b{\lambda }_b^2$ directly.

That squared-salience criterion is the algebraic distinction from the MB-PCA / SUM-PCA formulation on the Multiblock PCA card, whose sample-space criterion is a linear sum of captured block inertia. ComDim is not PCA of a concatenated superblock.

Hanafi and Qannari, following Qannari et al. 2000, initialize the weights to 1. The unit vector $\mathbf{q}$ is then the leading eigenvector of ${\sum }_b{\lambda }_b{\mathbf{W}}_b$. The weights are updated by ${\lambda }_b={\mathbf{q}}^{\top }{\mathbf{W}}_b\mathbf{q}$. The steps are repeated. At a stationary point,

The same stationary relation is the weighted summing-up step in Tchandao Mangamana et al.: $\mathbf{t}\propto {\sum }_b{\alpha }_b{\mathbf{t}}_b$ with ${\alpha }_b={\mathbf{t}}^{\top }{\mathbf{t}}_b$. The current MBAnalysis implementation uses this eigen-iteration. The DNRutledge MATLAB repository instead extracts the leading left singular vector of $[\sqrt{{\lambda }_1}{\mathbf{X}}_1\mid \cdots \mid \sqrt{{\lambda }_B}{\mathbf{X}}_B]$, which at a stationary point solves the same symmetric eigenproblem. Those two numerical routes are not different scientific methods. Local-score scaling in that MATLAB code is a different output convention and is not used here.

Hanafi and Qannari prove that the loss decreases monotonically and is bounded below, so the iteration converges. This card does not claim that every start yields a unique global maximizer. The public solution is a stationary iterative solution under the documented start.

Initialization is deterministic: ${\lambda }_b=1$, so the first direction is the leading eigenvector of ${\sum }_b{\mathbf{W}}_b$. That is the Hanafi–Qannari start and the MBAnalysis start. It is not random. It is not the first column of the first table used in comdim_PCA_2020.m. The default numerical tolerance on the relative change of ${\sum }_b{\lambda }_b^2$ is 1e-8, matching MBAnalysis. That default is an implementation parameter, not a scientific law. The DNRutledge code monitors the change in $\mathbf{q}$ instead. This card uses only the relative criterion change.

The extracted $\mathbf{q}$ is the unit-norm global score / common dimension in sample space. Tchandao Mangamana et al. define the associated block component by ${\mathbf{q}}^{(b)}={\mathbf{W}}_b\mathbf{q}$. The natural block loading is ${\mathbf{p}}_b={\mathbf{X}}_b^{\top }\mathbf{q}$, so ${\mathbf{X}}_b{\mathbf{p}}_b={\mathbf{W}}_b\mathbf{q}={\mathbf{q}}^{(b)}$. Loadings and saliences are different quantities.

After a unit global direction is stored, Hanafi and Qannari replace each block by the residual of its orthogonal projection onto that direction:

Cross-product matrices are rebuilt from those residual blocks. The removed variation is the sample-space contribution of the extracted dimension. It is not noise. Under this deflation, successive global directions are orthogonal. Block loadings and block components do not inherit that orthogonality automatically.

If whole-block normalization is applied so that $\Vert {\mathbf{X}}_b{\Vert }_F^2=1$, residual-stage saliences on orthogonal global directions are fractions of that unit block inertia. Without that normalization, a block's initial inertia enters the saliences directly. Variable-wise autoscaling is a different operation and still changes ${\mathbf{W}}_b$.

A global dimension ${\mathbf{q}}_a$ is a shared sample-space axis used to represent the blocks. Its salience ${\lambda }_{ba}$ is how much residual-stage inertia of block $b$ lies along that axis under the chosen preprocessing. The block loading and block component show how that block's variables and sample variation relate to the same axis. These three outputs are not interchangeable.

A large salience means substantial block variance along that common direction after preprocessing. It does not mean the best instrument, the most important science, a unique subspace, or a causal driver. A small salience means relatively little of that block's residual inertia is aligned with the axis. It does not mean the block is useless, pure noise, or without other structure.

Both methods are unsupervised multiblock component methods with global and block-level quantities. Tchandao Mangamana et al. give the distinction used here. Both use the link ${\mathbf{t}}_b={\mathbf{W}}_b\mathbf{t}$. MB-PCA, in their Table 1, updates the global component by the unweighted sum of block components and maximizes a linear sum of covariances. ComDim updates by the salience-weighted sum and maximizes the sum of squared covariances. The Open Lab Multiblock PCA card teaches the equivalent CPCA-W / SUM-PCA superblock PCA. Those are different optimization problems. Neither is universally better. ComDim is not another name for MB-PCA.

Tchandao Mangamana et al. state that Hanafi, Kohler and Qannari 2010 showed ComDim to be equivalent to Hierarchical PCA as introduced by Wold, Kettaneh and Tjessem 1996, under the formulation in their Table 5. They also note a different H-PCA version, discussed by Westerhuis, Kourti and MacGregor, that rescales block components to unit length after each update and can have convergence problems. That version is not this card. Equivalence is therefore formulation-specific. This card does not state that every H-PCA implementation is ComDim. GCCA belongs to the other link-function family in their Table 1, based on orthogonal projectors onto block column spaces and correlation criteria. GCCA is not derived here. Correlation and variance criteria are not equivalent.

There is no universal rule such as 95% inertia, a fixed $A$, or one dimension per block. Useful evidence includes salience profiles, residual-stage captured block inertia, stability, and interpretability. If block normalization gives unit total block inertia, the sum of residual-stage saliences over extracted orthogonal dimensions is the captured fraction of that block. Component signs are arbitrary: $\mathbf{q}$ and $-\mathbf{q}$ are the same direction. Saliences are invariant to that flip.

The listing is NumPy. Symmetric eigenproblems use numpy.linalg.eigh. Autoscaling uses the sample standard deviation with denominator $N-1$, matching the Normalization & Scaling card. Whole-block normalization, when selected, divides by the Frobenius norm of the already centered and optionally autoscaled block. sklearn PCA is not ComDim. It may be used only as a single-block check. MBAnalysis ComDim documentation and CRAN 2.2.0 source were inspected. R was not runtime tested. Variable scaling in that package uses $N$ in the standard deviation and is not claimed identical to this listing when autoscale=True.

import numpy as np  def _validate_blocks(blocks, label="blocks"):    if not isinstance(blocks, (list, tuple)) or len(blocks) < 1:        raise ValueError("%s must be a non-empty list of 2D arrays." % label)    out = []    n_samples = None    for index, block in enumerate(blocks):        block = np.asarray(block, dtype=float)        if block.ndim != 2:            raise ValueError("Block %d must be a 2D array." % (index + 1))        if min(block.shape) < 1:            raise ValueError("Block %d must have positive dimensions." % (index + 1))        if not np.all(np.isfinite(block)):            raise ValueError("Block %d must contain only finite values." % (index + 1))        if n_samples is None:            n_samples = block.shape[0]        elif block.shape[0] != n_samples:            raise ValueError("All blocks must share the same number of samples.")        out.append(block)    return out  def _leading_unit_eigenvector(matrix):    symmetric = 0.5 * (matrix + matrix.T)    eigenvalues, eigenvectors = np.linalg.eigh(symmetric)    direction = eigenvectors[:, int(np.argmax(eigenvalues))]    norm = np.linalg.norm(direction)    if not np.isfinite(norm) or norm == 0:        raise ValueError("The weighted cross-product matrix has no usable leading direction.")    return direction / norm  def _preprocess_blocks(blocks, center, autoscale, normalize_blocks, tol):    n_samples = blocks[0].shape[0]    n_variables = [block.shape[1] for block in blocks]    if autoscale and n_samples < 2:        raise ValueError("Autoscaling requires at least two observations.")    means = []    variable_scales = []    block_norms = []    processed = []    for block in blocks:        if center:            mean = block.mean(axis=0)        else:            mean = np.zeros(block.shape[1])        centered = block - mean        if autoscale:            scale = centered.std(axis=0, ddof=1)            if np.any(~np.isfinite(scale) | (scale == 0)):                raise ValueError("Scaling cannot be applied to a zero variance variable.")        else:            scale = np.ones(block.shape[1])        scaled = centered / scale        if normalize_blocks:            frobenius = np.linalg.norm(scaled, ord="fro")            if not np.isfinite(frobenius) or frobenius <= tol:                raise ValueError(                    "Block normalization is undefined for a processed block with zero Frobenius norm."                )            scaled = scaled / frobenius        else:            frobenius = 1.0        means.append(mean)        variable_scales.append(scale)        block_norms.append(float(frobenius if normalize_blocks else 1.0))        processed.append(scaled)    return {        "n_samples": n_samples,        "n_variables": n_variables,        "means": means,        "variable_scales": variable_scales,        "block_norms": np.asarray(block_norms, dtype=float),        "processed": processed,    }  def fit_comdim(    blocks,    n_components,    center=True,    autoscale=False,    normalize_blocks=True,    tol=1e-8,    max_iter=200,):    blocks = _validate_blocks(blocks)    n_blocks = len(blocks)    n_components = int(n_components)    if n_components < 1:        raise ValueError("n_components must be a positive integer.")    if max_iter < 1:        raise ValueError("max_iter must be a positive integer.")    prepared = _preprocess_blocks(blocks, center, autoscale, normalize_blocks, 1e-15)    n_samples = prepared["n_samples"]    n_total = sum(prepared["n_variables"])    max_rank = min(n_samples, n_total)    if center:        max_rank = min(n_samples - 1, n_total) if n_samples > 1 else 0    if n_components > max_rank:        raise ValueError("n_components exceeds the rank available after preprocessing.")     residual = [np.array(block, dtype=float, copy=True) for block in prepared["processed"]]    global_scores = np.zeros((n_samples, n_components))    saliences = np.zeros((n_blocks, n_components))    loadings = [np.zeros((n_var, n_components)) for n_var in prepared["n_variables"]]    block_components = [np.zeros((n_samples, n_components)) for _ in range(n_blocks)]    criterion = np.zeros(n_components)    n_iterations = np.zeros(n_components, dtype=int)     for component in range(n_components):        cross_products = [block @ block.T for block in residual]        weights = np.ones(n_blocks)        weighted = sum(cross_products)        direction = _leading_unit_eigenvector(weighted)        previous_fit = 0.0        for block_index, cross_product in enumerate(cross_products):            weights[block_index] = float(direction @ cross_product @ direction)        previous_fit = float(np.sum(weights ** 2))        delta = tol + 1.0        iteration = 0        while delta > tol and iteration < max_iter:            iteration += 1            weighted = sum(                weight * cross_product                for weight, cross_product in zip(weights, cross_products)            )            direction = _leading_unit_eigenvector(weighted)            for block_index, cross_product in enumerate(cross_products):                weights[block_index] = float(direction @ cross_product @ direction)            fit = float(np.sum(weights ** 2))            if previous_fit == 0:                delta = 0.0 if fit == 0 else np.inf            else:                delta = (fit - previous_fit) / previous_fit            previous_fit = fit        if delta > tol:            raise ValueError("ComDim iteration did not meet the selected tolerance.")         global_scores[:, component] = direction        saliences[:, component] = weights        criterion[component] = previous_fit        n_iterations[component] = iteration        projector = np.outer(direction, direction)        for block_index, block in enumerate(residual):            loading = block.T @ direction            local = block @ loading            loadings[block_index][:, component] = loading            block_components[block_index][:, component] = local            residual[block_index] = block - projector @ block     return {        "n_samples": n_samples,        "n_blocks": n_blocks,        "n_variables": prepared["n_variables"],        "n_components": n_components,        "center": center,        "autoscale": autoscale,        "normalize_blocks": normalize_blocks,        "tol": float(tol),        "max_iter": int(max_iter),        "means": prepared["means"],        "variable_scales": prepared["variable_scales"],        "block_norms": prepared["block_norms"],        "processed_blocks": prepared["processed"],        "Q": global_scores,        "saliences": saliences,        "loadings": loadings,        "block_components": block_components,        "criterion": criterion,        "n_iterations": n_iterations,        "block_explained_ss": saliences.copy(),    } 

The MATLAB listing implements the same Hanafi–Qannari / MBAnalysis eigen-iteration with base eig on the symmetrized weighted sum of cross-products. It does not call comdim_PCA_2020. That author repository was inspected: default whole-block normalization after column centering, iteration until the change in $\mathbf{q}$ is small, and rescaled local scores ${\mathbf{X}}_b{\mathbf{p}}_b/({\mathbf{p}}_b^{\top }{\mathbf{p}}_b)$. Those local scores are not compared as raw numbers with ${\mathbf{W}}_b\mathbf{q}$. Signs of eigenvectors may differ across languages. Align signs before comparing global scores. Saliences are sign-invariant.

function model = fit_comdim(blocks, nComponents, center, autoscale, normalizeBlocks, tol, maxIter)    if nargin < 3, center = true; end    if nargin < 4, autoscale = false; end    if nargin < 5, normalizeBlocks = true; end    if nargin < 6, tol = 1e-8; end    if nargin < 7, maxIter = 200; end    blocks = validateBlocks_comdim(blocks, "blocks");    nBlocks = numel(blocks);    nSamples = size(blocks{1}, 1);    nVariables = zeros(1, nBlocks);    for b = 1:nBlocks        nVariables(b) = size(blocks{b}, 2);    end    nComponents = double(nComponents);    if nComponents < 1 || nComponents ~= floor(nComponents)        error('nComponents must be a positive integer.');    end    if autoscale && nSamples < 2        error('Autoscaling requires at least two observations.');    end    nTotal = sum(nVariables);    maxRank = min(nSamples, nTotal);    if center        if nSamples > 1            maxRank = min(nSamples - 1, nTotal);        else            maxRank = 0;        end    end    if nComponents > maxRank        error('nComponents exceeds the rank available after preprocessing.');    end     means = cell(1, nBlocks);    variableScales = cell(1, nBlocks);    blockNorms = ones(1, nBlocks);    processed = cell(1, nBlocks);    residual = cell(1, nBlocks);    for b = 1:nBlocks        if center            means{b} = mean(blocks{b}, 1);        else            means{b} = zeros(1, nVariables(b));        end        centered = blocks{b} - means{b};        if autoscale            scale = std(centered, 0, 1);            if any(~isfinite(scale) | scale == 0)                error('Scaling cannot be applied to a zero variance variable.');            end        else            scale = ones(1, nVariables(b));        end        variableScales{b} = scale;        scaled = centered ./ scale;        if normalizeBlocks            frobenius = norm(scaled, "fro");            if ~isfinite(frobenius) || frobenius <= 1e-15                error('Block normalization is undefined for a processed block with zero Frobenius norm.');            end            blockNorms(b) = frobenius;            scaled = scaled / frobenius;        end        processed{b} = scaled;        residual{b} = scaled;    end     Q = zeros(nSamples, nComponents);    saliences = zeros(nBlocks, nComponents);    loadings = cell(1, nBlocks);    blockComponents = cell(1, nBlocks);    for b = 1:nBlocks        loadings{b} = zeros(nVariables(b), nComponents);        blockComponents{b} = zeros(nSamples, nComponents);    end    criterion = zeros(1, nComponents);    nIterations = zeros(1, nComponents);     for a = 1:nComponents        W = cell(1, nBlocks);        for b = 1:nBlocks            W{b} = residual{b} * residual{b}';        end        weights = ones(nBlocks, 1);        weighted = zeros(nSamples, nSamples);        for b = 1:nBlocks            weighted = weighted + W{b};        end        q = leadingUnitEigenvector_comdim(weighted);        for b = 1:nBlocks            weights(b) = q' * W{b} * q;        end        previousFit = sum(weights.^2);        delta = tol + 1;        iteration = 0;        while delta > tol && iteration < maxIter            iteration = iteration + 1;            weighted = zeros(nSamples, nSamples);            for b = 1:nBlocks                weighted = weighted + weights(b) * W{b};            end            q = leadingUnitEigenvector_comdim(weighted);            for b = 1:nBlocks                weights(b) = q' * W{b} * q;            end            fit = sum(weights.^2);            if previousFit == 0                if fit == 0                    delta = 0;                else                    delta = Inf;                end            else                delta = (fit - previousFit) / previousFit;            end            previousFit = fit;        end        if delta > tol            error('ComDim iteration did not meet the selected tolerance.');        end        Q(:, a) = q;        saliences(:, a) = weights;        criterion(a) = previousFit;        nIterations(a) = iteration;        for b = 1:nBlocks            loading = residual{b}' * q;            local = residual{b} * loading;            loadings{b}(:, a) = loading;            blockComponents{b}(:, a) = local;            residual{b} = residual{b} - q * (q' * residual{b});        end    end     model = struct();    model.nSamples = nSamples;    model.nBlocks = nBlocks;    model.nVariables = nVariables;    model.nComponents = nComponents;    model.center = logical(center);    model.autoscale = logical(autoscale);    model.normalizeBlocks = logical(normalizeBlocks);    model.tol = tol;    model.maxIter = maxIter;    model.means = means;    model.variableScales = variableScales;    model.blockNorms = blockNorms;    model.processedBlocks = processed;    model.Q = Q;    model.saliences = saliences;    model.loadings = loadings;    model.blockComponents = blockComponents;    model.criterion = criterion;    model.nIterations = nIterations;    model.blockExplainedSs = saliences;end function blocks = validateBlocks_comdim(blocks, label)    if nargin < 2, label = "blocks"; end    if ~iscell(blocks) || isempty(blocks)        error('%s must be a non-empty cell array of 2D arrays.', label);    end    nSamples = [];    for b = 1:numel(blocks)        block = double(blocks{b});        if ~ismatrix(block) || any(size(block) < 1)            error('Block %d must be a 2D array with positive dimensions.', b);        end        if ~all(isfinite(block), "all")            error('Block %d must contain only finite values.', b);        end        if isempty(nSamples)            nSamples = size(block, 1);        elseif size(block, 1) ~= nSamples            error('All blocks must share the same number of samples.');        end        blocks{b} = block;    endend function q = leadingUnitEigenvector_comdim(matrix)    symmetric = 0.5 * (matrix + matrix');    [vectors, values] = eig(symmetric);    [~, index] = max(real(diag(values)));    q = real(vectors(:, index));    nrm = norm(q);    if ~isfinite(nrm) || nrm == 0        error('The weighted cross-product matrix has no usable leading direction.');    end    q = q / nrm;end 

P-ComDim is a supervised extension built for a response block. It is not derived here and is not classical ComDim. Path-ComDim (Cariou, Qannari, Rutledge and Vigneau 2018) extends ComDim to blocks with an a priori pattern of directed relations. That paper's path structure is a model assumption. It is not treated here as a proof of causality. AComDim applies ComDim to ANOVA-style factor and residual tables. The 2024 review describes ComDim-ICA, ComDim-PLS, and ComDim-OPLS as variants that replace or adapt the internal decomposition. None of those algorithms is implemented on this card.

1. 1.

   Qannari, E. M., Wakeling, I., Courcoux, P., & MacFie, H. J. H. (2000). Defining the underlying sensory dimensions. Food Quality and Preference, 11(1-2), 151-154.

   doi:10.1016/S0950-3293(99)00069-5
2. 2.

   Qannari, E. M., Courcoux, P., & Vigneau, E. (2001). Common components and specific weights analysis performed on preference data. Food Quality and Preference, 12(5-7), 365-368.

   doi:10.1016/S0950-3293(01)00026-X
3. 3.

   Mazerolles, G., Hanafi, M., Dufour, E., Bertrand, D., & Qannari, E. M. (2006). Common components and specific weights analysis: A chemometric method for dealing with complexity of food products. Chemometrics and Intelligent Laboratory Systems, 81(1), 41-49.

   doi:10.1016/j.chemolab.2005.09.004
4. 4.

   Hanafi, M., Mazerolles, G., Dufour, E., & Qannari, E. M. (2006). Common components and specific weight analysis and multiple co-inertia analysis applied to the coupling of several measurement techniques. Journal of Chemometrics, 20(5), 172-183.

   doi:10.1002/cem.988
5. 5.

   Hanafi, M., & Qannari, E. M. (2008). Nouvelles propriétés de l'analyse en composantes communes et poids spécifiques. Journal de la Société Française de Statistique, 149(2), 75-96.

6. 6.

   Hanafi, M., Kohler, A., & Qannari, E. M. (2010). Shedding new light on hierarchical principal component analysis. Journal of Chemometrics, 24(11-12), 703-709.

   doi:10.1002/cem.1334
7. 7.

   Cariou, V., Qannari, E. M., Rutledge, D. N., & Vigneau, E. (2018). ComDim: From multiblock data analysis to path modeling. Food Quality and Preference, 67, 27-34.

   doi:10.1016/j.foodqual.2017.02.012
8. 8.

   Cariou, V., Jouan-Rimbaud Bouveresse, D., Qannari, E. M., & Rutledge, D. N. (2019). ComDim methods for the analysis of multiblock data in a data fusion perspective. In Cocchi, M. (Ed.), Data Fusion Methodology and Applications. Data Handling in Science and Technology, 31, 179-204.

   doi:10.1016/B978-0-444-63984-4.00007-7
9. 9.

   Tchandao Mangamana, E., Cariou, V., Vigneau, E., Glèlè Kakaï, R. L., & Qannari, E. M. (2019). Unsupervised multiblock data analysis: A unified approach and extensions. Chemometrics and Intelligent Laboratory Systems, 194, 103856.

   doi:10.1016/j.chemolab.2019.103856
10. 10.

    Jouan-Rimbaud Bouveresse, D., & Rutledge, D. N. (2024). A synthetic review of some recent extensions of ComDim. Journal of Chemometrics, 38(5), e3454.

    doi:10.1002/cem.3454
11. 11.

    Wold, S., Kettaneh-Wold, N., & Tjessem, K. (1996). Hierarchical multiblock PLS and PC models for easier model interpretation and as an alternative to variable selection. Journal of Chemometrics, 10(5-6), 463-482.

    doi:10.1002/(SICI)1099-128X(199609)10:5/6<463::AID-CEM445>3.0.CO;2-L
12. 12.

    Mahieu, B., Tchandao Mangamana, E., Vigneau, E., & Cariou, V. (2026). MBAnalysis: Multiblock exploratory and predictive data analysis. R package version 2.2.0, function ComDim. CRAN.