I've placed some R code on Github that is useful for working with higher-order spectral estimators. These estimators are based on the higher-order singular value decomposition of De Lathauwer et. al. [2000] and are useful when your data exhibit tensor-specific structure, such as having approximately low multilinear rank. This code will allow you to: