Fixation is essential for preserving cellular morphology in biomedical research. However, it may also affect spectra captured in multispectral fluorescence microscopy, impacting molecular interpretations. To investigate fixation effects on tissue, multispectral fluorescence microscopy images of pairs of samples with and without fixation are captured. Each pixel might exhibit overlapping spectra, creating a blind source separation problem approachable with linear unmixing. With multiple excitation wavelengths, unmixing is intuitively extended to tensor factorizations. Yet these approaches are limited by nonlinear effects like attenuation. Further, light exposure during image acquisition introduces subtle Brownian motion between image channels of non-fixed tissue. Finally, hypothesis testing for spectral differences due to fixation is non-trivial as retrieved spectra are paired sequential samples. To these ends, we present three contributions, (1) a novel robust non-negative tensor factorization using the β-divergence and L2,1-norm, which decomposes the data into a low-rank multilinear and group-sparse non-multilinear tensor without making any explicit nonlinear modeling choices or assumptions on noise statistics; (2) a diffeomorphic atlas-based strategy for motion correction; (3) a non-parametric hypothesis testing framework for paired sequential data using functional principal component analysis. PyTorch code for robust non-negative tensor factorization is available at https://github.com/neel-dey/robustNTF.