Existence of a unique rotation matrix maximizing the quartimax criterion I am trying to prove that there exists a unique orthonormal matrix $T \in \mathbb{R}^{p \times p}$ (up to permutation and sign changes) maximizing the quartimax criterion as
$$ \underset{T}{\arg\max} \,\, \sum_{i=1}^n \sum_{j=1}^p (AT)_{ij}^4 $$
where $A \in \mathbb{R}^{n \times p}$ has orthonormal columns ($A^TA = I)$ for $n > p$.
I checked that the quartimax criterion is maximized when each row of $AT$ has at most one nonzero element (simple structure) or when $A$ has full rank with $n=p$.
I'm unsure how to proceed for the general case where $A$ does not have that simple structure with $n>p$.
 A: $
\def\bbR#1{{\mathbb R}^{#1}}
\def\o{{\tt1}}\def\p{\partial}
\def\LR#1{\left(#1\right)}
\def\BR#1{\Big(#1\Big)}
\def\bR#1{\big(#1\big)}
\def\trace#1{\operatorname{Tr}\LR{#1}}
\def\skew#1{\operatorname{Skew}\LR{#1}}
\def\cay#1{\operatorname{Cayley}\LR{#1}}
\def\vc#1{\operatorname{vec}\LR{#1}}
\def\qiq{\quad\implies\quad}
\def\grad#1#2{\frac{\p #1}{\p #2}}
\def\hess#1#2#3{\frac{\p^2 #1}{\p #2\,\p #3}}
\def\c#1{\color{red}{#1}}
\def\CLR#1{\c{\LR{#1}}}
\def\fracLR#1#2{\LR{\frac{#1}{#2}}}
\def\gradLR#1#2{\LR{\grad{#1}{#2}}}
$For typing convenience, define the variables
$$\eqalign{
M &= AT \\
M^{\odot 3} &= M\odot M\odot M \\
Q &= 4\,A^TM^{\odot 3} \\
}$$
where the symbol $(\odot)$ denotes Hadamard products/powers.
It will also be convenient to use $J$ to denote the all-ones matrix (which is the same size as $A)\,$ and a colon to denote matrix inner product
$$\eqalign{
A:M &= \sum_{i=1}^n\sum_{j=1}^p A_{ij}M_{ij} \;=\; \trace{A^TM} \\
A:A &= \|A\|^2_F \\
}$$
Introduce an unconstrained matrix variable $U\in\bbR{p\times p}$ and use it to construct $T$
$$\eqalign{
S &= \skew{U} &\doteq \fracLR{U-U^T}2 \\
T &= \cay{S} &\doteq \LR{I+S}^{-1}\LR{I-S} \\
S &= \cay{T} \\
}$$
These definitions have some interesting consequences
$$\eqalign{
\LR{I-S}^{-1}T &= \LR{I+S}^{-1} \\
I+T &= \LR{I+S}^{-1}\BR{(I+\c{S})+(I-\c{S})} \\
    &= 2\;\LR{I+S}^{-1} \\
 dT &= -2\;\LR{I+S}^{-1}\,dS\,\LR{I+S}^{-1} \\
}$$
Now we are ready to calculate the gradient $\bR{{\rm wrt}\:U}$
of the quartimax function
$$\eqalign{
\phi &= J:M^{\odot 4} \\
d\phi &= J:\LR{4M^{\odot 3}\odot dM} \\
 &= {4 M^{\odot 3}}:\LR{A\;dT} \\
 &= \c{Q}:dT \\
 &= -2\,Q:\LR{I+S}^{-1}\,dS\,\LR{I+S}^{-1} \\
 &= -2\,\LR{I-S}^{-1}Q\LR{I-S}^{-1}:\skew{dU} \\
 &= -2\,\skew{\LR{I-S}^{-1}Q\LR{I-S}^{-1}}:{dU} \\
\grad{\phi}{U} 
 &= -2\,\skew{\LR{I-S}^{-1}Q\LR{I-S}^{-1}} \\
 &= {\LR{I+S}^{-1}Q^T\LR{I+S}^{-1}\;-\;\LR{I-S}^{-1}Q\LR{I-S}^{-1}} \\
 &= \LR{I-S}^{-1}\LR{TQ^TT-Q}\LR{I-S}^{-1} \\
}$$
Setting the gradient to zero yields the first order condition (FOC)
$$\eqalign{
Q &= TQ^TT \\
}$$
Since $(Q,T)$ are functions of $U$ this represents a complicated function of $U$ which must be solved for the optimal matrix $U_*$ which will satisfy the FOC.  This will be very difficult.
Conversely, a numerical solution will be very simple. Since you have an exact expression for the gradient $\gradLR{\phi}{U}$ you can employ any gradient ascent algorithm to calculate the optimal $U_*$
Analyzing the numerical solutions for an ensemble of $(n,p,A)$ parameters might reveal something about the structure of $T$.
