\section{Representations for weight states}\label{sec:repweight}
This section discusses border rank and MPS representations for the weight states. For the convenience of the reader, we recall the definition of the (unnormalized) weight states from the main text:
\begin{defi}[Weight states]
Given $a,L,\mathfrak{d}\in\mathds{N}_0$ the \emph{weight state} of weight $a$ on $L$ sites of local dimension $\mathfrak{d}+1$ is:
\begin{align*}
     \ket{\chi_{a,\mathfrak{d},L}}=\sum\limits_{\substack{i_1+i_2+\ldots+i_L=a \\0\leq i_1,\ldots,i_L\leq \mathfrak{d}}}  \ket{i_1,i_2,\ldots,i_L} \in \left( \mathds{C}^{\mathfrak{d}+1}\right)^{\otimes L} .
\end{align*}
Write $\ket{\chi_{a,L}}=\ket{\chi_{a,a,L}}$.
\end{defi}

A state $\ket{\psi}\in\left( \mathds{C}^{d}\right)^{\otimes L}$ admits a representation of tensor rank $r$~\cite{Lan:TensorBook} if there exist product vectors $\ket{u_1},\ldots,\ket{u_r}\in \left( \mathds{C}^{d}\right)^{\otimes L}$ such that
\begin{align}\label{equ:rankdecomposition}
\ket{\psi}=\ket{u_1} + \cdots + \ket{u_r}.
\end{align}

A state $\ket{\psi}\in\left( \mathds{C}^{d}\right)^{\otimes L}$ admits a representation of border rank $r$~\cite{Bini_1980,christandl2018tensor} if there exist product vectors $\ket{u_1(\epsilon)},\ldots,\ket{u_r(\epsilon)}\in \left( \mathds{C}^{d}\right)^{\otimes L}$, depending polynomially on $\epsilon$, and a nonnegative integer $a$ such that
\begin{align}\label{equ:borderrankdecomposition}
\ket{\psi}= \lim _{\epsilon\to0}\frac{1}{\epsilon^a} \left( \ket{u_1(\epsilon)} + \cdots + \ket{u_r(\epsilon)}\right).
\end{align}

The minimal integer $r$ such that $\ket{\psi}$ has a representation of tensor rank (resp. border rank) is called tensor rank (resp. border rank) of $\ket{\psi}$, denoted $\mathrm{R}( \ket{\psi})$ (resp. $\underline{\mathrm{R}}(\ket{\psi})$).

Note that $\ket{\chi_{a,\mathfrak{d},L}}$ can be obtained from $\ket{\chi_{a,L}}$ applying local projection maps $P_{\mathfrak{d}}:\mathds{C}^{a+1}\to \mathds{C}^{\mathfrak{d}+1}$ defined by
\begin{align*}
    P_{\mathfrak{d}}={\textstyle \sum}_{k=0}^{\mathfrak{d}}\ketbra{k}{k}.
\end{align*}
In particular, a border rank $r$ representation of $\ket{\chi_{a,L}}$ immediately provides a border rank $r$ representation of $\ket{\chi_{a,\mathfrak{d},L}}$. Thus, we focus on efficient representations of $\ket{\chi_{a,L}}$.

\begin{prop}[Border rank of weight states]\label{thm:degwweight}
Let $\ket{\chi_{a,L}}\in\left( \mathds{C}^{a+1}\right)^{\otimes L}$ be a weight state. For $\epsilon \in \mathds{C}$, let $\ket{\phi(\epsilon)}\in\mathds{C}^{a+1}$ be the state
\begin{align*}
    \ket{\phi_a(\epsilon)}=\sum\limits_{i=0}^a\epsilon^i\ket{i}.
\end{align*}
Then
\begin{equation}\label{eqn: border rank decomp weight}
    \ket{\chi_{a,L}}= \frac{1}{a!} \lim _{\epsilon \to 0}\frac{1}{\epsilon^a} \left[ \sum_{j=0}^a (-1)^{a-j}{\binom{a}{j}}\ket{\phi(j\epsilon)}^{\otimes L}\right].
\end{equation}
In particular, $\underline{\mathrm{R}}(\ket{\chi_{a,L}}) = a+1$.
\end{prop}
\begin{proof}
The lower bound is immediate via a standard flattening argument, see e.g.~\cite[Section 7.2]{BalBerChrGes:PartiallySymRkW} or~\cite{efremenko94barriers}.
We verify \eqref{eqn: border rank decomp weight} via a direct calculation:
\begin{align*}
\sum_{j=0}^a (-1)^{a-j}{\binom{a}{j}}\ket{\phi(j\epsilon)}^{\otimes L} &= \sum_{j=0}^a (-1)^{a-j} \binom{a}{j} \left[ {\textstyle \sum}_{i=0}^a (j \epsilon)^i \ket{i} \right]^{\otimes L}  \\
&=\sum_{j=0}^a (-1)^{a-j} \binom{a}{j} \left[ {\textstyle \sum}_{i_1 , \dots , i_L } j^{i_1 + \cdots +i_L} \epsilon^{i_1 + \cdots +i_L} \ket{i_1 , \dots , i_L} \right]  \\
&=\sum_{j=0}^a (-1)^{a-j} \binom{a}{j} \sum_{\alpha \geq 0} j^\alpha \epsilon^\alpha \left[ {\textstyle \sum}_{i_1 + \cdots + i_L =\alpha} \ket{i_1 , \dots , i_L} \right]  \\
&= \sum_{\alpha \geq 0} \epsilon^\alpha  \ket{ \chi_{\alpha,a,L}} \left[ {\textstyle \sum}_{j=0}^a  (-1)^{a-j} \binom{a}{j} j^\alpha \right].
\end{align*}
To conclude the proof we show that the coefficients $c_{\alpha,a} = {\textstyle \sum}_{j=0}^a  (-1)^{a-j} \binom{a}{j} j^\alpha$ satisfy $c_{\alpha,a} = 0$ if $\alpha < a$ and $c_{a,a} = a!$. To see this, we use a double induction on $a$ and $\alpha$.

If $(\alpha,a) = (0,0)$ then clearly $c_{\alpha,a}=1$, and if $\alpha = 0$ and $a > 0$, then $c_{\alpha,a} =  {\textstyle \sum}_{j=0}^a  (-1)^{a-j} \binom{a}{j} = (1-1)^a = 0$ by Newton's binomial formula.

Suppose $\alpha >0$; then the $0$-th term of the summation is $0$. We obtain
\begin{align*}
 c_{\alpha,a} & = \sum_{j=1}^a (-1)^{a-j} \binom{a}{j} j^\alpha  \\
 &=\sum_{j=1}^a (-1)^{a-j} \frac{a!}{j!(a-j)!} j^\alpha  \\
 &=a \sum_{j=1}^a (-1)^{(a-1)-(j-1)} \frac{(a-1)!}{(j-1)!(a-j)!} j^{\alpha-1}   \\
 &=a \sum_{j=0}^{a-1} (-1)^{(a-1)-j} \frac{(a-1)!}{j!((a-1)-j)!} (j+1)^{\alpha-1}   \\
 &=a \sum_{j=0}^{a-1} (-1)^{(a-1)-j} \frac{(a-1)!}{j!((a-1)-j)!} \left[ \sum_{\beta=0}^{\alpha-1} \binom{\alpha-1}{\beta} j^\beta  \right] =\\
 &=a \sum_{\beta=0}^{\alpha-1} \binom{\alpha-1}{\beta} \sum_{j=0}^{a-1} (-1)^{(a-1)-j} \binom{a-1}{j} j^\beta = a \sum_{\beta=0}^{\alpha-1} \binom{\alpha-1}{\beta}  c_{\beta, a-1}.
\end{align*}
By the inductive hypothesis, $c_{\beta,a-1} =0$ when $\beta < \alpha-1$ and $c_{a-1,a-1} = (a-1)!$. Since $\beta \leq \alpha-1$, we deduce that the summation above vanishes whenever $\alpha < a$. When $\alpha = a$, the summation reduces to $c_{\alpha, a} = a c_{\alpha-1,a-1} = a(a-1)! = a!$. This concludes the proof.
\end{proof}

A consequence of Proposition \ref{thm:degwweight} is that $\underline{\mathrm{R}}(\ket{\chi_{a,\mathfrak{d},L}}) \leq a+1$, and in particular it is bounded from above independently of $L$. Moreover, all product states in the border rank expression of Proposition \ref{thm:degwweight} are symmetric. Regarding tensor rank, it is a classical fact \cite{Sylv:PrinciplesCalculusForms} that $\mathrm{R}( \ket{\chi_{1,L}}) = L+1$, which implies that $\mathrm{R}(\ket{\chi_{a,L}}) \geq L+1$ grows (at least) linearly in the system size $L$; on the other hand, since the border rank expression for $\ket{\chi_{a,L}}$ has degree $L-a$ in $\epsilon$, a standard interpolation argument (see, e.g., \cite[Thm. 8]{Christandl2017_TR} and \cite[Prop. 6.2]{ChrGesJen:BorderRankNonMult}) shows that $\mathrm{R}(\ket{\chi_{a,L}}) \leq (L-a+1) (a+1)$.

The following result provides a matrix product state representation with open boundary conditions for the weight states. As before, we only deal with the case where $\mathfrak{d}=a$. This representation gives rise to the MPS contraction strategy described in the section \emph{Performing  computations  in  the  ansatz  class} of the main text.

\begin{lemma}\label{lemma: MPS for chi}
Fix $a,L\in\mathds{N}$ and let $P_L$ be the graph consisting of an open chain with $L$ nodes indexed by $\{1 , \dots , L\}$. For $k = 1 , \dots , L$, define a local map $\mathcal{A}^{k}:\mathds{C}^{a+1}\otimes \mathds{C}^{a+1}\to \mathds{C}^{a+1}$ with $\mathcal{A}^{k}=\sum\limits_{j=0}^{a}\ket{j}\otimes A^{(k)}_{j}$ where
 \begin{align*}
 A^{(1)}_j &= \ket{j}\otimes \bra{j}, \\
 A^{(k)}_{j}&=\sum \limits_{i=0}^{a-j}\ket{i}\otimes \bra{i+j} \text{ for }  2\leq k\leq L-1 ,\\
 A^{(L)} _j &=\sum_{i=0}^{a} \ket{a-j}\otimes \bra{i}.
 \end{align*}
Let $\ket{\phi} = ({\textstyle \bigotimes}_{k =1}^L \mathcal{A}_k ) \ket{\Omega^{(P_L)}_{a+1}}$ be the resulting tensor network state. Then $\ket{\phi} = \ket{\chi_{a,L}}$.
\end{lemma}
\begin{proof}
Write $A_j$ for $A_j^{(k)}$ when $k = 2 , \dots , L-1$ and define $A_j = 0$ for $j > a$. It is easy to see that  $A_jA_{j'}=A_{j+j'}$ for all $j,j'$. Therefore, if $i_2+i_3+\ldots+i_{L-1}>a$, we obtain
\begin{align*}
\langle i_1\ldots i_L| \phi\rangle=0.
\end{align*}
Therefore, suppose $i_2+i_3+\ldots+i_{L-1} = b \leq a$, so that $A_{i_2}A_{i_3}\ldots A_{i_{L-1}}=A_b$. Then
\begin{align*}
\langle i_1\ldots i_L| \phi\rangle=1
\end{align*}
if $i_1+i_{L}= a-b$ and $0$ otherwise. This concludes the proof.
\end{proof}

Regarding matrix product states on the ring $C_L$, Lemma \ref{lemma: MPS for chi} implies that $\bTNS^{C_L}_{D,a,\mathfrak{d},d} \subseteq \TNS^{C_L}_{(a+1)D,d}$, relying simply on the fact that $\ket{\Omega^{(C_L)}_{D(a+1)}} =\ket{ \Omega^{(C_L)}_D} \otimes \ket{\Omega^{(P_L)}_{a+1}} \otimes \ket{\Omega^{(e)}_{a+1}}$, where $e$ is the edge which makes the chain $P_L$ into the ring $C_L$.

\section{Separations in efficiency of representations}\label{sec:efficiencyseparations}
As mentioned in the main text, for matrix product states representations on the cycle $C_L$, the maximal separation between the border bond dimension and bond dimension is at most quadratic. For the sake of completeness, we now provide a proof of this folklore result.

\begin{prop}[Separation between border bond and bond dimension]\label{prop: bond quadratic in border bond}
Let $\ket{\psi}\in\left( \mathds{C}^d\right)^{\otimes L}$ be a state with $\bbond^{C_L} (\ket{\psi}) = D$. Then $\bond^{C_L}(\ket{\psi}) \leq D^2$. In particular
\[
 \cTNS^{C_L}_{D} \subseteq \TNS^{C_L}_{D^2}.
\]

\end{prop}
\begin{proof}
Let $P_L$ be the chain obtained from $C_L$ by removing a single edge $e_L$. Then $\TNS^{P_L}_{D^2}$ is a closed set \cite{Landsberg_geometry} and clearly $\TNS^{P_L}_{D^2} \subseteq \TNS_{D^2}^{C_L}$. We show that $\cTNS^{C_L}_D \subseteq \TNS^{P_L}_{D^2}$ and this will conclude the proof.

Note that $\TNS^{P_L}_{D^2}$ can be interpreted as having two copies of maximally entangled stated of bond dimension $D$ on the chain $P_L$. By performing an entanglement swapping protocol, we can then generate a maximally entangled state between sites $1$ and $L$ with one of the maximally entangled states, which shows that we can generate any TNS on the cycle by starting with two copies on the chain.

More formally, let $ T : \mathbb{C}^D\otimes \mathbb{C}^D \to \mathbb{C}$ be the map defined by $T= \sum_1^D \bra{i} \otimes \bra{i}$. Define a family of linear maps $(\mathcal{A}_k : k =1 , \dots , L)$ by $\mathcal{A}^1 = \mathcal{A}^k = \text{id}_{\mathbb{C}^D}$ and $\mathcal{A}^k = T$ for $k =2 , \dots , L-1$. Then $(\bigotimes_k \mathcal{A}^k) \ket{\Omega^{(P_L)} _D }= \ket{\Omega^{(e_L)}_D}$.

Since $\ket{\Omega^{(C_L)}_{D}} = \ket{\Omega^{(P_L)}_{D}} \otimes \ket{\Omega^{(e_L)}_D}$, one immediately has $\ket{\Omega^{(C_L)}_{D} }\in \TNS_{D^2}^{P_L}$ by applying the family $(\mathcal{A}_k : k =1 , \dots , L)$ on one factor of $\ket{\Omega^{(P_L)}_{D^2}} = \ket{{\Omega^{(P_L)}_{D}}}^{\otimes 2}$. Therefore,
\[
 \TNS_D^{C_L} \subseteq \TNS_{D^2}^{P_L}.
\]
Passing to the closures, we obtain $ \cTNS_D^{C_L} \subseteq \cTNS_{D^2}^{P_L}$, and since $\cTNS_{D^2}^{P_L} = \TNS_{D^2}^{P_L}$ is closed we conclude.
\end{proof}

The proof above is purely topological and not constructive: it only ensures the existence of a matrix product states representation of bond dimension $D^2$ under the hypothesis of the existence of a border bond dimension $D$ representation.

The same argument can be generalized to arbitrary graphs as follows. Given a graph $G$, let $T$ be a spanning tree of $G$, i.e., a tree on all vertices of $G$. In particular, if $e$ is an edge of $G$ which is not an edge of $T$, then $T \cup \{e \}$ contains a cycle; the entanglement swapping argument of Proposition \ref{prop: bond quadratic in border bond} shows that $\Omega^{(e)}_D \in \TNS_D^{T}$. In this way, one can obtain the maximally entangled states $\Omega_D^{(e)}$ corresponding to all edges of $G$ which are not edges of $T$.  However, without any assumptions on $G$, one might need a full $\Omega_D^{(T)}$ to generate each edge.

Since $T$ does not contain any cycle, the set $\TNS_D^T$ is a closed set and we conclude
\[
\cTNS^{G}_D \subseteq \TNS^T_{D^\eta} \subseteq \TNS^G_{D^\eta}
\]
where $\eta$ is the number of edges of $G$ which are not edges of $T$. As any tree with $|V|$ vertices has $|E|=|V|-1$ edges, we conclude that  for a graph $G=(V,E)$
\begin{align*}
\cTNS^{G}_D \subseteq \TNS^G_{D^{|E|-|V|+2}}.
\end{align*}

\subsection{Examples of system-size dependent separations in efficiency of representation}\label{sec:analyticalexamples}
In this section, we discuss an example of a state $\ket{\psi_L} \in (\mathbb{C}^3)^{\otimes L}$ for odd $L$ which admits a border bond dimension $2$ translation invariant matrix product state representation with approximation degree $a=1$, yet having translation invariant bond dimension growing with $L$. More precisely,
\begin{align*}
\bbond^{\text{TI-}C_L}(\ket{\psi_L}) &= 2 \text{ (with $a = \mathfrak{d} = 1$)} \\
\bond^{\text{TI-}C_L}(\ket{\psi_L}) & = \Omega \left( \textfrac{L^{1/3}}{\log(L)} \right) .
\end{align*}
Moreover, the state $\ket{\psi_L}$ arises as ground state of a translation invariant local (nearest neighbor) Hamiltonian on the cycle. In particular, we obtain $\ket{\psi_L} \in \bTNS^{\text{TI-}C_L}_{2,1,1,3}$ and $\ket{\psi_L} \notin \TNS^{\text{TI-}C_L}_{D,3}$ if $D$ grows slower than $\Omega( \frac{L^{1/3}}{\log(L)})$.

In the non-translation invariant case, there are no known examples of separations that are larger than a constant factor.

The separation in the translation invariant setting, and in particular the lower bound on $\bond^{\text{TI-}C_L}(\ket{\psi})$ relies on the analogous result for the $W$-state, proved with the same argument as \cite[Corollary 1]{mapsgarcia} via the results of \cite{sanz2010quantum,Rahaman_2020,Michalek_2019_Wielandt}.

\begin{lemma}\label{lem:repwstate}
Let $C_L$ the cycle on $L$ nodes. Then
\begin{align*}
 (a+1)L\geq\bond^{\text{TI-}C_L}(\ket{\chi_{a,L}}) &= \Omega  \left(\frac{L^{1/3}}{\log(L)}\right) ,\\
\bbond^{\text{TI-}C_L} (\ket{\chi_{a,L}}) &\leq a+1.
\end{align*}
\end{lemma}
\begin{proof}
The upper bound $\bbond^{\text{TI-}C_L} (\ket{\chi_{a,L}})\leq a+1$ is obtained directly from the border rank expression of Proposition \ref{thm:degwweight}.  This is because the superposition of $a+1$ symmetric product states can be written as a TI-MPS of bond dimension $a+1$.

The upper bound $\bond^{\text{TI-}C_L}(\ket{\chi_{a,L}}) \leq (a+1)L$ is obtained symmetrizing the non-translation invariant representation obtained in Lemma \ref{lemma: MPS for chi}.

As for the lower bound $\bond^{\text{TI-}C_L}(\ket{\chi_{a,L}}) = \Omega  \left(\frac{L^{1/3}}{\log(L)}\right)$, first note that
the (local and translation invariant) projection map $P^{\otimes L}$ where $P = \ketbra{0}{0} + \ketbra{1}{a}$ sends $\ket{\chi_{a,L}}$ to the $W$-state $\ket{\chi_{1,L}}$. In particular, lower bounds on $\bond^{\text{TI-}C_L}(\ket{\chi_{1,L}})$ are lower bounds for $\ket{\chi_{a,L}}$ as well.

In~\cite[Corollary 1]{mapsgarcia}, the authors related $\bond^{\text{TI-}C_L}(\ket{\chi_{1,L}})$ to quantum Wielandt inequalities~\cite{sanz2010quantum,Rahaman_2020,Michalek_2019_Wielandt}. Although the version of the inequality stated in~\cite[Conjecture 2]{mapsgarcia}, on which \cite[Corollary 1]{mapsgarcia} relies on, is only conjectured, in~\cite{Michalek_2019_Wielandt} the authors show that~\cite[Conjecture 2]{mapsgarcia} holds with $f(D)=\mathcal{O}(D^2\log(D))$.

Reproducing the proof of~\cite[Corollary 1]{mapsgarcia} provides
\begin{align*}
    (D)^3\log(D)^3=\Omega(L),
\end{align*}
whenever $D \geq \bond^{\text{TI-}C_L}(\ket{\chi_{1,L}})$. This implies $ \bond^{\text{TI-}C_L}(\ket{\chi_{1,L}}) =\Omega(e^{\mathcal{W}(L^{1/3})})$, where $\mathcal{W}$ is the Lambert $\mathcal{W}$ (see, e.g., \cite{hoorfar2008inequalities}). Indeed, one can show \cite[(1.1)]{hoorfar2008inequalities} that
\begin{align*}
\mathcal{W}(x)\geq \log(x)-\log \log(x)
\end{align*}
for all $x\geq e$, which implies that $e^{\mathcal{W}(L^{1/3})}=\Omega\left(\frac{L^{1/3}}{\log(L)}\right) $.
\end{proof}

From Lemma \ref{lem:repwstate}, already $\ket{\chi_{a,L}}$ provides an example of a state for which we have a separation.

However, in the main text we presented a local translationally invariant Hamiltonian whose unique translationally invariant ground state presents this behaviour as well. Let us recall its definition and discuss some of its properties. The weight states themselves give separations in the required bond dimension and will be used as stepping stones to prove separations for these examples. But the next examples require a bond dimension $D>1$ to be represented in the ansatz class introduced  class, unlike the weight states.

Let $P :\mathbb{C}^3 \otimes \mathbb{C}^3 \to \mathbb{C}^3 \otimes \mathbb{C}^3$ be the projection defined by
\begin{align}\label{equ:definitonprojection}
P=\ketbra{01}{01}+\ketbra{10}{10}+\ketbra{02}{02}+\ketbra{21}{21};
\end{align}
denote by $P_{i,i+1}$ the linear map $P$ acting on the two factors of $(\mathbb{C}^3)^{\otimes L}$ corresponding to sites $i$ and $i+1$ of the ring $C_L$.
Define $H_{i,i+1} = \text{id} - P_{i,i+1} + \ketbra{2}{2}_i$, where the third summand acts on the $i$-th site. Let $H = \sum_i H_i$: this is a diagonal, translationally invariant Hamiltonian.

We determine its translationally invariant ground state in the case where $L$ is odd.
\begin{prop}
Let $L$ be odd. Then $H$ has a unique translationally invariant ground state given by
\begin{align}\label{equ:TIgroundstate}
\ket{\psi_L}=\frac{1}{\sqrt{L}}\sum\limits_{i=0}^{L-1}S^i\ket{2101010\ldots 10},
\end{align}
where $S$ is the cyclic permutation of the tensor factors, i.e., the shift operator.
\end{prop}
\begin{proof}
Since $H$ is diagonal, every product state in the standard basis is and eigenstate for $H$. Moreover, using the fact that the $P_{i,i+1}$ commute, one obtains that the spectrum of $H$ is contained in the set
\begin{align*}
\left\{\ell+\frac{m}{2L}|\ell,m \in\mathds{N}_{0}\right\}.
\end{align*}
Given a product state in the standard basis, we compute its energy,i.e. the corresponding eigenvalue of $H$.

Note that every state involving two adjacent $\ket{0}$'s or two adjacent $\ket{1}$'s has energy at least $1$, since $1$ is the eigenvalue of the corresponding $H_i$. Moreover, since $L$ is odd, if the state does not involve any $\ket{2}$, there are necessarily adjacent $\ket{0}$'s or $\ket{1}$'s. Every state involving at least two $\ket{2}$'s has energy at least $2 \cdot 1/(2L)$, arising from the $\ketbra{2}{2}$ summands in $H$. There are only two basis state, up to cyclic permutation, involving exactly one $\ket{2}$ and having no adjacent $\ket{0}$'s nor $\ket{1}$'s: they are $\ket{20101 \ldots 01}$ and $\ket {21010 \ldots 10}$. The first has energy at least $1$ arising from $H_{12}$. The state $\ket {21010 \ldots 10}$ has energy $1/(2L)$.

We conclude that $\{ S^k\ket{2101010\ldots 10} : k = 0 , \dots , L-1\}$ span the ground state space of $H$. The only translation invariant state in this span is $\ket{\psi_L}$.
\end{proof}

The following result gives upper and lower bounds for the translation invariant bond and border bond dimension of $\ket{\psi_L}$:
\begin{prop}
Let $\ket{\psi_L} \in (\mathbb{C}^3)^{\otimes L}$ be the state of \eqref{equ:TIgroundstate}. Then
 \begin{align*}
\bbond^{\text{TI-}C_L}(\ket{\psi_L}) &= 2;\\
  \bond^{\text{TI-}C_L}(\ket{\psi_L}) &= \Omega\left( \frac{L^{\frac{1}{3}}}{\log(L)}\right).
 \end{align*}
In fact, $\ket{\psi_L} \in \bTNS^{\text{TI-}C_L}_{2,1,1,3}$ and, for $L$ large enough, $\ket{\psi_L} \notin \TNS^{\text{TI-}C_L}_{D,3}$ if $D$ does not depend on $L$.
\end{prop}
\begin{proof}
Let $\mathcal{Q}: \mathbb{C}^3 \to \mathbb{C}^2$ be the local map defined by
\begin{align*}
\mathcal{Q}=\ketbra{0}{1}+\ketbra{0}{0}+\ketbra{1}{2}.
\end{align*}
Then $\mathcal{Q}^{\otimes L} \ket{\psi_L}=\ket{\chi_{1,L}}$ is the $W$-state. In particular, the lower bound for $\ket{\chi_{1,L}}$ obtained in Lemma~\ref{lem:repwstate} holds for $\bond^{\text{TI-}C_L}(\ket{\psi_L})$ as well.

Now, define $\mathcal{A}(\epsilon): \mathds{C}^{2}\otimes \mathds{C}^{ 2} \to\mathds{C}^3$ by
\begin{align*}
\mathcal{A}(\epsilon)=\ketbra{0}{01}+\ketbra{1}{10}+\epsilon \ketbra{2}{00}.
\end{align*}
One can readily check that
\begin{align*}
\lim_{\epsilon \to 0} \epsilon^{-1} \left( \mathcal{A}(\epsilon)^{\otimes L}\right) \ket{\Omega_2^{(C_L)}}=\ket{\psi_L},
\end{align*}
giving a representation of $\ket{\psi_L}$ of border bond dimension $2$ and approximation degree $1$. In particular, $\ket{\psi_L} \in \bTNS^{C_L}_{2,1,1,3}$ and concludes the proof.
\end{proof}
One can define explicitly a local map $\mathcal{B} : (\mathds{C}^{2}\otimes \mathds{C}^{ 2} )  \otimes \mathbb{C}^2 \to \mathbb{C}^3$ by
\[
(\ketbra{0}{01}+\ketbra{1}{10}) \otimes \bra{0}+\ketbra{2}{00} \otimes \bra{1}
\]
and observe $\ket{\psi_L} = \mathcal{B}^{\otimes L}\ket{ \Omega^{(C_L)}_2} \otimes \ket{\chi_{1,L}}$.

This shows that there are unbounded separations in the bond dimension required to represent states that arise as ground states of local Hamiltonians in the new ansatz class for translationally invariant matrix product state when compared with standard matrix product states.

\subsection{An explicit example of separation between bond and border bond dimension}

For small system size and small local dimension, several explicit examples of separation between bond and border bond dimension are known. In \cite{Str:RelativeBilComplMatMult}, Strassen proved that the GHZ state of level three on three parties, $\ket{\mathrm{ghz}_3} =3^{-\frac{1}{2}}\left( \ket{000} + \ket{111} + \ket{222}\right)$, satisfies $\bbond^{C_3} (\ket{\mathrm{ghz}}) = 2$, realized by a degeneration with approximation degree $a=2$; in fact, one can achieve the same result with a degeneration of approximation degree $a=1$ and in addition in \cite{christandl2018tensor} it was proved that $\bond^{C_3} (\ket{\mathrm{ghz}}) = 3$, providing a separation.

In \cite{christandl2018tensor}, a separation was shown for the state $\ket{\lambda} = \sum_{\sigma \in \mathfrak{S}_3} (-1)^\sigma \ket{\sigma(0)\sigma(1)\sigma(2) } + |222\rangle$ as well; here $\mathfrak{S}_3$ denotes the symmetric group on three elements and $(-1)^\sigma$ is the sign of a permutation $\sigma$; in particular $\sum_{\sigma \in \mathfrak{S}_3} (-1)^\sigma \ket{\sigma(0)\sigma(1)\sigma(2) }$ is the unique (up to scaling) alternating tensor in $\mathbb{C}^3 \otimes \mathbb{C}^3 \otimes \mathbb{C}^3 $. It was proved that $\bbond^{C_3} (\ket{\lambda}) =2$ and $\bond^{C_3} = 3$.

For these examples, comparing the numerical results in the standard tensor network ansatz class with the ones in the new ansatz class, we observe improvements both in the distance and in the convergence. However, in this case, the class $\bTNS^{C_3} _{2,1,1,3}$ is dense in $\mathbb{C}^3 \otimes \mathbb{C}^3 \otimes \mathbb{C}^3$, so the improvements occur for a trivial reason, as the optimization procedure can just ``walk straight'' to the global minimum.

In order to see an example where the optimization procedure is not trivial, we consider the state
\begin{equation}\label{eqn: tensor degen 3x3 mamu}
\begin{aligned}
 \ket{T} &=\sqrt{17}^{-1}( \ket{005}+\ket{016}+\ket{040}+\ket{126}+\ket{160}+\ket{227}+\ket{251}+\ket{262}+\ket{338}\\
 &+\ket{373}+\ket{384}+\ket{430}+\ket{501}+\ket{632}+\ket{703}+\ket{714}+\ket{824}) \in \mathbb{C}^9 \otimes \mathbb{C}^9 \otimes \mathbb{C}^9,
 \end{aligned}
\end{equation}
introduced in the \emph{Numerical examples} section of the main text. Lemma \ref{lemma: bond and border bond of T} below shows $\ket{T} \notin \TNS^{C_3}_{3,9}$ and $\ket{T} \in \bTNS^{C_3}_{3,1,9}$. Moreover, the set $\bar{\bTNS}^{C_3}_{3,1,9}$ is a proper subvariety of $\mathbb{C}^9 \otimes \mathbb{C}^9 \otimes \mathbb{C}^9$, therefore the new ansatz class is non-trivial.

\begin{lemma}\label{lemma: bond and border bond of T}
Let $\ket{T}$ be the state defined in \eqref{eqn: tensor degen 3x3 mamu}. Then
\begin{align*}
\bond^{C_3}{(\ket{T})} &\geq 4 , \\
\bbond^{C_3}{(\ket{T})}& = 3 \text{ and the degeneration has approximation degree $a = 1$}.
\end{align*}
In particular, $\ket{T} \notin \TNS^{C_3}_{3,9}$ and $\ket{T} \in \bTNS^{C_3}_{3,1,9}$.
\end{lemma}
\begin{proof}
Assume by contradiction $\bond^{C_3}(\ket{T}) \leq 3$ and observe that $\ket{\Omega^{(C_3)}_3} \in \mathbb{C}^9 \otimes \mathbb{C}^9 \otimes \mathbb{C}^9$, where each $\mathbb{C}^9$ has a local structure $\mathbb{C}^3 \otimes \mathbb{C}^3$. For $i=1,2,3$, let $\mathcal{A}_i : (\mathbb{C}^3 \otimes \mathbb{C}^3) \to \mathbb{C}^9$ be the linear maps such that $(\bigotimes \mathcal{A}^i) (\ket{\Omega^{(C_3)}_3}) = \ket{T}$. A direct calculation shows that there is no proper subspace $U \subsetneq \mathbb{C}^9$ such that $\ket{T} \in U \otimes \mathbb{C}^9 \otimes \mathbb{C}^9$ or $\ket{T} \in \mathbb{C}^9 \otimes U \otimes \mathbb{C}^9$ or $\ket{T} \in \mathbb{C}^9 \otimes \mathbb{C}^9 \otimes U$. Thus, the maps $\mathcal{A}^i$ are invertible, and $\ket{T}$ is a point of the orbit of $\ket{\Omega^{(C_3)}_3}$ under the action of $GL_9 \times GL_9 \times GL_9$. As a consequence, the stabilizer of $\ket{T}$ and the one of $\ket{\Omega^{(C_3)}_3}$ are conjugate and in particular they have the same dimension. But this is not the case: a direct calculation shows the stabilizer of $\ket{\Omega^{(C_3)}_3}$ is $26$-dimensional whether the stabilizer of $\ket{T}$ is $28$-dimensional. This shows the lower bound $\bond^{C_3} (\ket{T}) \geq 4$.

We provide a degeneration of $\ket{\Omega^{(C_3)}_3}$ with approximation degree $a=1$ which realizes $\ket{T}$. In other words, we determine three linear maps $\mathcal{A}^{i}(\epsilon) : \mathbb{C}^3 \otimes \mathbb{C}^3 \to \mathbb{C}^9$ such that $\lim _{\epsilon \to 0} ( \frac{1}{\epsilon} \bigotimes _i \mathcal{A}^i(\epsilon)) \ket{\Omega^{(C_3)}_3} = \ket{T}$. Write $\mathcal{A}^i(\epsilon) = \mathcal{A}_0^{i} +  \epsilon \mathcal{A}_1^i$ and express $\mathcal{A}_0^i$ and $\mathcal{A}_1^i$ as matrices with entries in $\mathbb{C}^9$:
\begin{align*}
 \mathcal{A}_0^1 = \mathcal{A}_0^2 = \left(\begin{array}{ccc}
                          \ket{0} & \ket{1} & \\
                          & \ket{2} & \\
                          & & \ket{3}
                         \end{array} \right)
                         & \quad
 \mathcal{A}_0^3 = \left(\begin{array}{ccc}
                         & \ket{1}   & \ket{3} \\
                          & & \ket{4} \\
                          \ket{0} & \ket{2} &
                         \end{array} \right) \\
 \mathcal{A}_1^1 = \mathcal{A}_1^2 = \left(\begin{array}{ccc}
                         & & \ket{4} \\
                         \ket{5} & & \ket{6} \\
                          \ket{7} & \ket{8} &
                         \end{array} \right)
&          \quad
\mathcal{A}_1^3 = \left(\begin{array}{ccc}
                          \ket{5} &  & \\
                          \ket{6} & \ket{7} & \\
                          & & \ket{8}
                         \end{array} \right)
                         & \quad
 \end{align*}
A direct calculation shows that $\lim _{\epsilon \to 0} ( \frac{1}{\epsilon} \bigotimes _i \mathcal{A}^i(\epsilon)) \ket{\Omega^{(C_3)}_3} = \ket{T}$.
\end{proof}

\section{Approximate representations}\label{sec:approximaterep}
As we already showed, there exist states that arise as ground states of local Hamiltonians and admit a more efficient \emph{exact} representation in $\textrm{bTNS}$. However, by definition, if a state lies on the boundary of the TNS variety, it can be approximated arbitrarily well by ``usual" TNS, albeit with local maps whose entries diverge as we increase the precision.

Most rigorous results available in the literature only ascertain that tensor network states of a polynomial bond dimension \emph{approximate} ground states well~\cite{Verstraete_2006}. Thus, it is of high importance to investigate the interplay between the speed at which the local maps diverge and how well we can approximate the target state by the restrictions arising from degenerations. This issue that was also raised in~\cite{Hackbusch_2019} recently. In other words, it is important to understand how small one should pick $\epsilon$ so that $\ket{\phi(\epsilon)}\simeq \ket{\phi(0)}$ for the task at hand. First, recall that to estimate most quantities of physical interest up to an additive error $\delta>0$,  such as the energy or two-point correlation functions of the ground state, it usually suffices to estimate states to a precision of $\mathcal{O}(\delta L^{-\alpha})$ in the Euclidean norm, where again $L$ is the system size and $\alpha\in\mathds{N}_0$. This is because these quantities are expectation values of observables with operator norm  $\mathcal{O}(L^{\alpha})$.

With that in mind, we are ready to discuss how to pick $\epsilon$ for a degeneration to ensure the desired precision. To gain some intuition about how the local maps behave in the regime of interest, we consider the case of the $W$-state, before considering the general setting.

Consider the family $\ket{\psi(\epsilon)}$ of TI-MPS on a cycle of size $L$ with bond dimension $2$ and matrices given by:
\begin{align}\label{equ:MPSdefining}
    A_0(\epsilon)= \left(\begin{matrix}1&0\\0& e^{\frac{\pi i}{L}}\end{matrix}\right), \quad A_1(\epsilon)=\left(\begin{matrix}\epsilon&0\\0&0\end{matrix}\right).
\end{align}
One can readily check that
\begin{align*}
    \lim\limits_{\epsilon\to 0}\epsilon^{-1}\ket{\psi(\epsilon)}=\ket{W},
\end{align*}
where $\ket{W}$ is the unnormalized $W$-state:
\begin{align*}
\ket{W}=\ket{100\ldots 0}+\ket{010\ldots 0}+\ldots+\ket{000\ldots 1}.
\end{align*}
Moreover, the transfer matrix $\mathcal{E}(\epsilon)$ of this family of MPS is given by:
\begin{align*}
   \mathcal{E}(\epsilon)= \left(\begin{matrix}1+\epsilon^2&0&0&0\\0& e^{\frac{\pi i}{L}}&0&0\\0&0& e^{\frac{\pi i}{L}}&0\\0&0&0&1\end{matrix}\right).
\end{align*}
Thus, we have that:
\begin{align*}
    \|\psi(\epsilon)\|^2=\text{tr}{\mathcal{E}(\epsilon)^L}=(1+\epsilon^2)^L-1.
\end{align*}
It is also straightforward to compute the overlap of $\ket{\psi(\epsilon)}$ with the $W$-state, as:
\begin{align*}
    \scalar{W}{\psi(\epsilon)}=L\text{tr}{A_0^{L-1}A_1}=L\epsilon.
\end{align*}
Combining these expressions and normalizing all involved states we have:
\begin{align}\label{equ:overlapWdeg}
    \left\|\frac{\ket{\psi(\epsilon)}}{\|\ket{\psi(\epsilon)}\|}-\frac{\ket{W}}{\sqrt{L}}\right\|^2=2\lb1-\frac{L\epsilon}{\left( L\left( \left( 1+\epsilon^2\right)^L-1\right)\rb^{\frac{1}{2}}}\right)=2\frac{\left( \left( 1+\epsilon^2\right)^L-1\right)^{\frac{1}{2}}-\sqrt{L}\epsilon}{\left( \left( 1+\epsilon^2\right)^L-1\right)^{\frac{1}{2}}}
\end{align}
Performing a Taylor expansion, we see that
\begin{align*}
2\frac{\left( \left( 1+\epsilon^2\right)^L-1\right)^{\frac{1}{2}}-\sqrt{L}\epsilon}{\left( \left( 1+\epsilon^2\right)^L-1\right)^{\frac{1}{2}}}=\mathcal{O}(L\epsilon^2).
\end{align*}
Thus, in order to ensure that
\begin{align*}
 \left\|\frac{\ket{\psi(\epsilon)}}{\|\ket{\psi(\epsilon)}\|}-\frac{\ket{W}}{\sqrt{L}}\right\|^2\leq \frac{\delta^2}{L^{2\alpha}}
\end{align*}
for some $\alpha\in \mathds{N}_0$, we need to pick $\epsilon=\mathcal{O}(\delta L^{-\frac{2\alpha+1}{2}})$.

Inserting this bound in the definitions of the matrices in Eq. \eqref{equ:MPSdefining} and normalizing the underlying MPS, we see that in order to obtain an approximation up to an error $\delta L^{-\alpha}$ to the W-state in Euclidean norm we need to pick the matrices with entries scaling like:
\begin{align*}
        A_0= \left(\frac{e^{\log(L)\frac{2\alpha+1}{2}}}{\delta}\right)^{\frac{1}{L}}\left(\begin{matrix}1&0\\0& e^{\frac{\pi i}{L}}\end{matrix}\right), \quad A_1=\left(\frac{e^{-\log(L)\frac{2\alpha+1}{2}}}{\delta}\right)\left(\begin{matrix}1&0\\0&0\end{matrix}\right).
\end{align*}
Note that $A_0$ has entries of constant order, while $A_1$ has entries of inverse polynomial order for $\alpha$ fixed. Moreover, the norm of the underlying MPS will be of order $\delta L^{-\alpha}$, which implies that trying to approximate the $W$-state with this family of MPS leads to ill-conditioned MPS. Setting the vanishing entries to $0$ leads to a state that is \emph{orthogonal} to the state of interest.
We expect that this behavior is universal, i.e., that \emph{any} approximation of states on the boundary of the set of TNS of a given bond dimension up to inverse polynomial in system size  precision leads to tensors with entries that are polynomially small and ill-conditioned.
Unfortunately, we still lack the tools to prove such statements for any approximation of a degeneration and can only prove it for some given degeneration, as we did for the $W$ state:
\begin{lemma}[Approximate representation from degeneration]\label{lem:approximate_rep}
Let $\mathcal{A}^k(\epsilon):\left( \mathds{C}^{D}\right)^{\otimes k_v}\to \mathds{C}^d$ be a family of local maps whose entries depend polynomially on $\epsilon$ defining a degeneration to $\ket{\phi_0}$ with $\|\phi_0\|=1$. Denote the expansion of
\begin{align*}
\ket{\psi(\epsilon)}=\left(\bigotimes_{v\in V}\mathcal{A}^v(\epsilon)\right)\left(\bigotimes\limits_{y\in E} \ket{\Omega^y}   \right)
\end{align*}
as a polynomial in terms of $\epsilon$ by:
\begin{align}\label{equ:expansionhigherorder}
\ket{\psi(\epsilon)}=\epsilon^{a}\ket{\phi_0}+\sum\limits_{l=1}^{e}\epsilon^{a+l}\ket{\phi_l}
\end{align}
for some $a\in\mathds{N}_0$ and set
\begin{align*}
\tau(\epsilon)=\sum\limits_{l=1}^{e}\epsilon^{l}\|\ket{\phi_l}\|,
\end{align*}
Then, for $\tau(\epsilon)\leq\frac{1}{2}$ we have:
\begin{align}\label{equ:firstapprox}
\left\|\frac{\ket{\psi(\epsilon)}}{\|\ket{\psi(\epsilon)}\|}-\ket{\phi_0}\right\|\leq 4\tau(\epsilon)
\end{align}
\end{lemma}
\begin{proof}

By the triangle inequality we have:
\begin{align}\label{equ:unnormalized2}
\left\|\frac{\ket{\psi(\epsilon)}}{\|\psi(\epsilon)\|}-\ket{\phi_0}\right\|\leq\left\|\frac{\ket{\psi(\epsilon)}}{\epsilon^a}-\ket{\phi_0}\right\|+\left\|\frac{\ket{\psi(\epsilon)}}{\|\ket{\psi(\epsilon)}\|}-\frac{\ket{\psi(\epsilon)}}{\epsilon^a}\right\|
\end{align}
Let us analyze the first term in the expression above. Expanding the polynomials we see that:
\begin{align}\label{equ:firsttriangle}
\left\|\frac{\ket{\psi(\epsilon)}}{\epsilon^a}-\ket{\phi_0}\right\|=\left\|\sum\limits_{l=1}^{e}\epsilon^{l}\ket{\phi_l}\right\|\leq \tau(\epsilon)
\end{align}
by a triangle inequality.
Let us now analyze the second term. We have:
\begin{align}\label{equ:normalized}
&\left\|\frac{\ket{\psi(\epsilon)}}{\|\ket{\psi(\epsilon)}\|}-\frac{\ket{\psi(\epsilon)}}{\epsilon^a}\right\|=\left|\frac{1}{\|\ket{\psi(\epsilon)}\|}-\frac{1}{\epsilon^a}\right|\|\ket{\psi(\epsilon)}\|=\left|\frac{1}{\epsilon^{-a}\|\ket{\psi(\epsilon)}\|}-1\right|\epsilon^{-a}\|\ket{\psi(\epsilon)}\|
\end{align}
Again by a triangle inequality and the assumption that $\|\phi_0\|=1$, it follows from eq.~\eqref{equ:expansionhigherorder} that:
\begin{align*}
\|\epsilon^{-a}\psi(\epsilon)\|\leq 1+\tau(\epsilon).
\end{align*}

Moreover, by a reverse triangle inequality
\begin{align*}
\left|\epsilon^{-a}\|\ket{\psi(\epsilon)}\|-\|\ket{\phi_0}\|\right|=\left|\epsilon^{-a}\|\ket{\psi(\epsilon)}\|-1\right|\leq\|\epsilon^{-a}\ket{\psi(\epsilon)}-\ket{\phi_0}\|\leq \tau(\epsilon)
\end{align*}
and we assumed that $\|\phi_0\|=1$.
This yields
\begin{align*}
\epsilon^{-a}\|\psi(\epsilon)\|\geq 1-\tau(\epsilon)
\end{align*}
and we obtain
\begin{align*}
\left|\frac{1}{\epsilon^{-a}\|\ket{\psi(\epsilon)}\|}-1\right|\leq\frac{\tau(\epsilon)}{1-\tau(\epsilon)}.
\end{align*}
Together with Eq.~\eqref{equ:normalized} we conclude that
\begin{align}\label{equ:boundonnormalized}
\left\|\frac{\ket{\psi(\epsilon)}}{\|\ket{\psi(\epsilon)}\|}-\frac{\ket{\psi(\epsilon)}}{\epsilon^a}\right\|\leq \frac{\tau(\epsilon)(1+\tau(\epsilon))}{1-\tau(\epsilon)}\leq 3\tau(\epsilon),
\end{align}
as we assumed that $\tau(\epsilon)\leq 1/2$.
Inserting the bounds of Equation~\eqref{equ:boundonnormalized} and~\eqref{equ:firsttriangle} into~\eqref{equ:unnormalized2} we obtain the claim in~\eqref{equ:firstapprox}.

\end{proof}
Thus, we see that in order to obtain an approximation up to an error $\mathcal{O}(\delta L^{-\alpha})$, it suffices to pick $\epsilon_0=\frac{\delta}{L^{\alpha}M}$, where $M=\max \|\phi_l\|$. Indeed, for this choice of parameters, we have $\tau(\epsilon_0)=\mathcal{O}(\delta L^{-\alpha})$. Moreover, we have:
\begin{align*}
\|\ket{\psi(\epsilon_0)}\|=\mathcal{O}\left(\frac{\delta^a}{L^{a\alpha}M^a}\right).
\end{align*}
and we expect that the condition number of the underlying maps scales in a similar way. This implies, for instance, that in order to obtain an estimate of the energy of the normalized state $\ket{\psi(\epsilon_0)}/\|\psi(\epsilon_0)\|$, it will be necessary to perform the multiplication $\scalar{\psi(\epsilon_0)}{H\psi(\epsilon_0)}\times\|\psi(\epsilon_0)\|^{-2}$, which entails multiplying two numbers of order $L^{-2a}\delta^{2a}$ and $L^{2a}\delta^{-2a}$, respectively.
For $a$ constant this requires a polynomial precision in the computation of expectation values and in case $a$ is at least linear in system size, this requires exponential precision.

We conclude that approximating states arising from degenerations by restrictions requires a great amount of precision when computing expectation values and will lead to numerical instabilities as the system size or degree of the degenerations increases. On the other hand, our methods remain stable.

This is even more pronounced in the case of tensor networks that are hard to contract exactly. In those cases, it is not feasible to contract the restrictions to the precision we just established are necessary. However, as highlighted in Section~\ref{sec:stableinterpolation}, approximate contractions yield good approximations in bTNS.

However, note that in principle some of the low order terms of $\tau(\epsilon)$ may vanish as well, which would lead to better convergence estimates than the ones discussed above. Let us illustrate this phenomenon in the case of the weight states:

\subsection{Approximate representations for the weight states}

Let us analyze more closely the trade-off between the rank of border rank representations of weight states and the speed at which they converge to the state. As noted in~\cite{Hackbusch_2019}, the representation given in~\eqref{equ:MPSdefining} for the W-state is suboptimal w.r.t. to the convergence speed. Indeed, one can check that:
\begin{align}\label{equ:betterrepW}
\frac{1}{2\epsilon}\left(\lb \ket{0}+\epsilon \ket{1}\right)^{\otimes L}-\left( \ket{0}-\epsilon \ket{1}\right)^{\otimes L}\right)=\ket{W}+\mathcal{O}(\epsilon^2),
\end{align}
while the degeneration given in~\eqref{equ:MPSdefining} only converges up to $\mathcal{O}(\epsilon)$.
This improved representation is closely connected with the fact that we can also obtain the weight states as:
\begin{align}\label{equ:representationderivative}
\ket{\chi_{a,L}}=a!\frac{d^a}{d \epsilon^a}\ket{\Gamma_a(\epsilon)}\big |_{\epsilon=0},
\end{align}
where
\begin{align*}
\ket{\Gamma_a(\epsilon)}=(\ket{0}+\epsilon \ket{1}+\ldots+\epsilon^{a}\ket{a})^{\otimes L}.
\end{align*}
For instance, the degeneration given in Eq.~\eqref{equ:betterrepW} arises from the symmetric difference formula for the derivative of a function:
\begin{align*}
\frac{f(x+\epsilon)-f(x-\epsilon)}{2\epsilon}=f'(x)+\mathcal{O}(\epsilon^2).
\end{align*}
Evaluating this formula for $\ket{\Gamma_a(\epsilon)}$ we obtain this improved degeneration. More generally, so-called central finite difference formulas~\cite[Section 10.10]{Quarteroni_2006} allow us to estimate the $a$-th derivative of a function up to an error of order $\mathcal{O}(\epsilon^{2(k+1)})$ by evaluating the function at $2\lfloor \frac{a+1}{2}\rfloor+2k$ points. This immediately yields:
\begin{lemma}\label{lem:betterdeg}
For $a,L\in \mathds{N}$ and $k\in\mathds{N}_0$ the weight states $\ket{\chi_{a,L}}$ admit a representation of symmetric border rank  $2\lfloor \frac{a+1}{2}\rfloor+2k$, $\ket{\chi_{a,L,k}(\epsilon)} $,  s.t.
\begin{align*}
\|\epsilon^{-a}\ket{\chi_{a,L,k}(\epsilon)} -\ket{\chi_{a,L}}\|=\mathcal{O}(\epsilon^{2(1+k)}).
\end{align*}
\end{lemma}
\begin{proof}
This immediately follows from the central finite difference formulas together with the representation of the weight states given in~\eqref{equ:representationderivative}.
\end{proof}
Moreover, it is possible to explicitly compute central difference formulas for a desired order of the derivative and accuracy by solving a system of linear equations~\cite{Quarteroni_2006}.
One corollary of Lemma~\ref{lem:betterdeg} is that, even though the weight states do not admit exact system-size independent TI representations as matrix product states, the accuracy with which we can approximate them increases exponentially with the bond dimension.

In principle these representations can be used to obtain more efficient approximate contractions schemes for states in bTNS. Indeed, by tensoring with a state with $\epsilon^{-a}\ket{\chi_{a,L,k}(\epsilon)}$ for $k$ large enough instead of $\ket{\chi_{a,L}}$ for small enough $\epsilon$ will lead to approximate contraction schemes that do not depend on $L$ and do not blow up the bond dimension of the underlying states. But it remains unclear at this stage how numerically stable these methods are.

\section{Stable Interpolation}\label{sec:stableinterpolation}
In this section, we discuss the stability result related to the border rank contraction strategy presented in the main text to compute the expectation values of an observable for a state in the new ansatz class. Given a state $\ket{\phi} \in \bTNS^{G}_{D,a,\mathfrak{d},d}$, the border rank contraction strategy considers the state $\ket{\phi(\epsilon)} = \ket{\Omega^{(G)}_D} \otimes \ket{\chi_{a,\mathfrak{d},L}(\epsilon)}$ where $\ket{\chi_{a,\mathfrak{d},L}(\epsilon)}$ is a border rank expression for $\ket{\chi_{a,\mathfrak{d},L}}$ (e.g., the one arising from Proposition \ref{thm:degwweight}. Given an observable $O : (\mathbb{C}^d)^{\otimes L} \to (\mathbb{C}^d)^{\otimes L}$, define
\begin{align}\label{equ:expectationvalueeps}
p(\epsilon) = \epsilon^{-2a}\scalar{\phi(\bar{\epsilon})}{O\phi(\epsilon)};
\end{align}
then $p(0) = \scalar{\phi}{O\phi}$ is the desired expectation value.

Rather than evaluating $p(0)$, in the border rank strategy we recover the value of $p(0)$ via Lagrange interpolation, by evaluating $p(\epsilon)$ at $2(L\mathfrak{d}-a) +1$ points. Expanding the polynomial, we observed that the evaluation of $p(\epsilon)$ can be performed in terms of $(a+1)^2$ overlaps of states in $\TNS^{G}_{D,d}$.

In the main text, we only considered exact contraction schemes to compute such overlaps. However, beyond the case of 1D systems, exact contraction of tensor network states is known to be hard~\cite{PhysRevLett.98.140506,Haferkamp_2020} and, in practice, approximate contraction schemes are used~\cite{Orus_practical}. Thus, it is of paramount importance for the applicability of the border rank contraction scheme that the interpolation step is robust: in other words, even if the values of $p(\epsilon)$ at the $2(L\mathfrak{d} -a)$ interpolation points are only exact up to an error $\delta$, one should still be able to recover the value at $0$ up to an error $\mathcal{O}(\delta)$.

Choosing the interpolation points to be roots of unity achieves this result.

\begin{lemma}[Stable polynomial interpolation]
Let $p:\mathds{C}\to\mathds{C}$ be a polynomial of degree $k$ and $z_j=e^{\frac{2\pi ij}{k}}$ for $0\leq j \leq k-1$ be $k$-th roots of unity. Let $c_0 , \dots , c_{k-1} \in\mathds{C}$ be $k$ points such that there exists $\delta > 0$ with
\begin{align}\label{equ:localapprox}
    \left|p(z_j)-c_j\right|\leq \delta.
\end{align}
Then:
\begin{align*}
    \left|p(0) -\frac{1}{k}{\textstyle \sum}_0^{k-1} c_j\right|\leq\delta.
\end{align*}
\end{lemma}
\begin{proof}
Via Lagrange interpolation
\begin{align*}
    p(0)=\frac{1}{k}{\textstyle \sum}_0^{k-1} p(z_j).
\end{align*}
Therefore
\begin{align*}
 \left|p(0) - \frac{1}{k}{\textstyle \sum}_0^{k-1} c_j\right| = \frac{1}{k}\left|{\textstyle \sum}_0^{k-1} ( p(z_j) - c_j ) \right| \leq \frac{1}{k}{\textstyle \sum}_0^{k-1}\bigl|p(z_j) - c_j\bigr| \leq \delta
\end{align*}
by a triangle inequality.
\end{proof}

This result shows that, choosing the interpolation points to be the $2aL$-th roots of unity, if the approximate contraction scheme only provides the values of $z_j^{-2a}p(z_j)$ up to an error $\delta$.
Note that $z_j^{-2a}$ is a complex number of modulus $1$, so  $c_j$ satisfying $|p(z_j)-c_j|\leq\delta$ also satisfies $|z_j^{-2a}p(z_j)-c_jz_j^{-2a}|$.
The resulting approximate value of $p(0)$ is up to an error $\delta$ as well.

\bibliographystyle{plain}


\section{Introduction}

Convolutional Neural Networks (ConvNets) have achieved phenomenal success in computer vision~\cite{krizhevsky2012imagenet,simonyan2014very,szegedy2015going,he2015delving,ioffe2015batch,he2016deep,xie2017aggregated,huang2017densely}.
While shallow ConvNets with a few layers have existed for decades~\cite{lecun1989backpropagation, denker1989neural, lecun1990handwritten, lecun1998gradient}, it is only until recently that networks with hundreds or even thousands of layers can be effectively trained.
Such \emph{deep} ConvNets are able to learn sophisticated decision rules for complex practical data, therefore are usually indispensable for obtaining state-of-the-art performance.

\begin{figure}[t]
\centering
\includegraphics[width=\linewidth]{figs/teaser.pdf}
\label{fig:net}
\vspace{-1.2em}
\caption{A basic block of an Isometric Network (ISONet). Each block contains only convolution and nonlinear activation layers, with an identity skip connections (for R-ISONet only). The convolution is initialized as the (Kronecker) delta kernel and regulated to be (near) orthogonal during training. The activation is Shifted ReLU (SReLU) $\phi(\cdot)$ with a learnable parameter $b$ for obtaining a balance between nonlinearity and isometry. The figure shows three examples of the SReLU with $b = 1, -1$ and $-3$.
}
\end{figure}

Training \emph{deep} ConvNets is inherently difficult~\cite{pascanu2013difficulty,glorot2010understanding}.
Despite widespread interests in this problem, it has not been made possible until the development of proper weight initialization~\cite{glorot2010understanding,he2015delving}, feature map normalization~\cite{ioffe2015batch},
and residual learning techniques~\cite{he2016deep,srivastava2015highway}.
Following these initial works, a wide variety of network architectural components including novel nonlinear activation, weight initialization, weight regularization, and residual learning techniques have been purposed in the literature (see Section~\ref{sec:previouswork} for an overview).
Such techniques are motivated and justified from diverse perspectives, such as prevention of dead neurons \cite{maas2013rectifier}, promotion of self-normalization \cite{klambauer2017self}, reduction of filter redundancy \cite{wang2019orthogonal}, improvement of generalization~\cite{jia2019orthogonal}, to name a few.

The abundance of existing architectural components and the diversity of their design principles posit ConvNets architectural design a difficult task.
After all, which combination of components should one use for their networks?
Such a challenge motivates us to pose the following question:
\emph{Is there a central guiding principle for training very deep ConvNets?}

\textbf{Isometric learning.}
We show that the \emph{isometric} property, where each layer of the network preserves the inner product for both forward and backward propagation, plays a central role in training deep ConvNets.
To illustrate this, we take a minimalist approach and show that a vanilla deep network (i.e., composed of interleaved convolution and nonlinear layers only) can be trained as long as both the convolution and nonlinear activation are close to an isometry.

Specifically, we design the Isometric Networks (ISONets) where the convolution layer is simply initialized as {\em the identity} and is regularized to be {\em near orthogonal} during training, and the nonlinear layer is the {\em Shifted ReLU} (SReLU) obtained by shifting regular ReLU towards being an identity around the origin with a learnable parameter (see Figure~\ref{fig:net}).
We perform extensive experiments with so designed simple networks on image classification using the ImageNet dataset~\cite{deng2009imagenet,russakovsky2015imagenet}.
Our results show that an ISONet with more than $100$ layers is trainable and can achieve surprisingly competitive performance.
To the best of our knowledge, this is the best performing \emph{vanilla} ConvNet for such tasks.

From a practical perspective, isometric learning provides an important guiding principle for addressing challenges in network architectural design.
One such challenge is that the commonly used normalization layers in modern networks
\cite{ioffe2015batch,ba2016layer,wu2018group}
require certain statistical independence assumptions to hold and large enough batch size or channel number for precise estimation of such statistics.
This drawback significantly limits their applications to robust learning \cite{sun2019test}, contrastive learning \cite{chen2020simple}, implicit models \cite{bai2020multiscale}, object detection and so on.

To address such a challenge, we develop a Residual ISONet (R-ISONet) that does not contain any normalization layer.
In particular,
we introduce skip connection to ISONet as it helps to better promote isometry of the network~\cite{tarnowski2019dynamical}.
We evaluate the performance of R-ISONet for object detection and instance segmentation on the COCO dataset~\cite{lin2014microsoft}.
In these applications, the batch size is very small due to the high resolution of the images, so that batch normalization becomes ineffective.
Our experiment shows that R-ISONet obtains better performance than standard ResNet.
Comparing with existing techniques for addressing the small-batch challenge~\cite{luo2018differentiable,wu2018group,singh2019filter}, R-ISONet has the benefit that it does not suffer from a slowdown during inference.

Isometry and many closely related notions have been implicitly and explicitly studied and utilized in many previous works for improving deep networks. In particular, techniques to promote isometric property have been widely practiced in the literature. However, they are often used together with many other design techniques and none of the previous work has clearly demonstrated that the isometry property {\em alone} ensures surprisingly strong performance for deep networks. At the end (Section \ref{sec:previouswork}), we discuss how isometric learning offers a simple but unified framework for understanding numerous seemingly different and diverse ideas in the literature, supporting the hypothesis that isometric learning offers a central guiding principle for designing and learning deep networks.

\section{Isometric Networks}
\label{sec:isonet}

In this section, we introduce the isometric principle and present the isometric network (ISONet).
In Section~\ref{sec:prelim}, we formally introduce the isometric principle in a vanilla net as enforcing both the convolution and nonlinear activation layers to be close to an \emph{isometry}.
Subsequently, Section~\ref{sec:convolution} derives the notion of isometry for convolution, and explains how it can be enforced at initialization and in the training process.
In Section~\ref{sec:activation}, we discuss principles for designing nonlinear activation functions in isometric learning and argue that one needs to strike a balance between isometry and nonlinearity, which cannot be achieved at the same time.
Finally, in Section~\ref{sec:residual}, skip connections  (or residual structure) can be naturally introduced to ISONet to further improve isometry, which leads to the R-ISONet.

\subsection{Isometric Learning}
\label{sec:prelim}

We first develop a vanilla network that is composed of interleaved convolution\footnote{We consider convolution for 2D signals (e.g., images), though our analysis trivially generalizes to arbitrary dimensional signals.} and nonlinear activation layers,
i.e.,
\begin{equation}
    \boldsymbol{x}^\ell = \phi(\boldsymbol{y}^\ell), ~~\boldsymbol{y}^\ell = \mathcal{A}^\ell \boldsymbol{x}^{\ell-1}, ~~~\ell=1, \ldots, L,\label{eq:forward}
\end{equation}
where $\mathcal{A}^\ell: \Re^{C^{\ell-1} \times H\times W} \to \Re^{C^{\ell} \times H\times W}$ denotes a convolution operator, $\phi(\cdot)$ denotes a point-wise nonlinear activation function, and $\boldsymbol{x}^0 \in \Re^{C^0 \times H \times W}$ is the input signal.
Assuming a squared loss given by Loss$=\frac{1}{2}\|\boldsymbol{z} - \boldsymbol{x}^L\|_2^2$, the backward propagation dynamic at $\boldsymbol{x}^0$ is given by
\begin{equation}\label{eq:backward}
    \frac{\partial \text{Loss}}{\partial \boldsymbol{x}^0} = (\mathcal{A}^1)^* \mathcal{D}^1 \cdots (\mathcal{A}^L)^* \mathcal{D}^L (\boldsymbol{z} - \boldsymbol{x}^L),
\end{equation}
where $(\mathcal{A}^\ell)^*: \Re^{ C^{\ell} \times H\times W} \to \Re^{C^{\ell-1} \times H \times W}$ is the adjoint operator of $\mathcal{A}^\ell$, and $\mathcal{D}^\ell: \Re^{C^{\ell} \times H \times W} \to \Re^{C^{\ell} \times H \times W}$ operates point-wise by multiplying the $(c, h, w)$-th entry of the input by $\phi'(\boldsymbol{y}_{c,h,w}^\ell)$.
For the operator $\mathcal{A}^\ell$ and $\phi(\cdot)$ in the forward dynamic and $(\mathcal{A}^\ell)^*$ and $\mathcal{D}^\ell$ in the backward dynamic, we may define the notion of \emph{isometry} as follows:
\begin{definition}[Isometry]
\label{def:isometry}
A map $\mathcal{A}: \Re^C \to \Re^M$ is called an isometry if
\begin{equation}
    \langle \mathcal{A} \boldsymbol{x}, \mathcal{A} \boldsymbol{x}' \rangle = \langle \boldsymbol{x}, \boldsymbol{x}' \rangle, ~~\forall \{\boldsymbol{x}, \boldsymbol{x}'\} \subseteq \Re^C.
\end{equation}
\end{definition}

Our ISONet is designed to maintain that all transformations ($\mathcal{A}, \mathcal{A}^*, \phi$, and $\mathcal{D}$) in the forward and the backward dynamics are close to an \emph{isometry}.
In the following, we develop techniques for enforcing isometry in convolution and nonlinear activation layers.
As we show in Section~\ref{sec:experiments}, combining these components enables effective training of ISONet with more than $100$ layers on image classification tasks.

\subsection{Isometry in Convolution}
\label{sec:convolution}

In this section, we show how to impose the property of isometry for convolutional layers in ConvNets, and explain how it can be achieved at initialization and through network training.

\textbf{Notation.}
For computation concerning convolution, we treat a 2D signal $\boldsymbol{\xi} \in \Re^{H \times W}$ as a function defined on the discrete domain $[1, \ldots, H] \times [1, \ldots, W]$ and further extended to the domain $\bb Z \times \bb Z$ by padding zeros.
Similarly, we regard any 2D kernel $\boldsymbol{\alpha} \in \Re^{k\times k}$ with $k = 2k_0 + 1$ for some integer $k_0$ as a function defined on the discrete domain $[k_0-1, \ldots, 0, \ldots, k_0+1] \times [k_0-1, \ldots, 0, \ldots, k_0+1]$ and extended to $\bb Z \times \bb Z$ by padding zeros.
$\boldsymbol{\xi}[i, j]$ represents the value of the function $\boldsymbol{\xi}$ at the coordinate $(i, j)$.

Given any $\boldsymbol{\xi} \in \Re^{H \times W}$ and $\boldsymbol{\alpha} \in \Re^{k \times k}$, the convolution of $\boldsymbol{\xi}$ and $\boldsymbol{\alpha}$ is defined as
\begin{equation}\label{eq:def-conv}
    (\boldsymbol{\alpha} * \boldsymbol{\xi})[i, j] = \sum_{p=-k_0}^{k_0} \sum_{q=-k_0}^{k_0} \boldsymbol{\xi}[i-p, j-q] \cdot \boldsymbol{\alpha}[p, q],
\end{equation}
and the correlation of $\boldsymbol{\xi}$ and $\boldsymbol{\alpha}$ is defined as
\begin{equation}\label{eq:def-corr}
    (\boldsymbol{\alpha} \star \boldsymbol{\xi})[i, j] = \sum_{p=-k_0}^{k_0} \sum_{q=-k_0}^{k_0} \boldsymbol{\xi}[i+p, j+q] \cdot\boldsymbol{\alpha}[p, q].
\end{equation}

In contrast to the conventions in signal processing, the convolution layers (i.e. $\{\mathcal{A}^\ell\}_{\ell=1}^L$ in \eqref{eq:forward}) in modern deep learning frameworks actually perform multi-channel \textit{correlation} operations that map a $C$-channel 2D signal to an $M$-channel 2D signal. Let $\boldsymbol{x} = (\boldsymbol{\xi}_1, \ldots, \boldsymbol{\xi}_C)\in \Re^{C\times H\times W}$ be the input signal where $\boldsymbol{\xi}_c \in \Re^{H \times W}$ for each $c \in \{1, \ldots, C\}$, and let
\begin{equation}\label{eq:kernel}
\resizebox{.9\hsize}{!}{$
    \boldsymbol{A} =
    \begin{pmatrix}
    \boldsymbol{\alpha}_{11} & \boldsymbol{\alpha}_{12} & \boldsymbol{\alpha}_{13} & \dots & \boldsymbol{\alpha}_{1C} \\
    \boldsymbol{\alpha}_{21} & \boldsymbol{\alpha}_{22} & \boldsymbol{\alpha}_{23} & \dots & \boldsymbol{\alpha}_{2C} \\
    \vdots & \vdots & \vdots & \ddots & \vdots \\
    \boldsymbol{\alpha}_{M1} & \boldsymbol{\alpha}_{M2} & \boldsymbol{\alpha}_{M3} & \dots & \boldsymbol{\alpha}_{MC}
\end{pmatrix}
\in \Re^{M\times C \times k \times k},$
}
\end{equation}
be the convolution kernel where each $\boldsymbol{\alpha}_{m,c}$ is a kernel of size $k \times k$ for $c\in\{1, \ldots, C\}$ and $m \in \{1, \ldots, M\}$.
The convolution operator $\mathcal{A}$ associated with $\boldsymbol{A}$ is given by
\begin{equation}\label{eq:conv_A}
    \mathcal{A} \boldsymbol{x} := \sum_{c=1}^C \big(\boldsymbol{\alpha}_{1c} \star \boldsymbol{\xi}_c, \ldots, \boldsymbol{\alpha}_{Mc} \star\boldsymbol{\xi}_c\big) \in \Re^{M \times H \times W}.
\end{equation}

The adjoint of $\mathcal{A}$ that appears in the backward dynamics \eqref{eq:backward}, denoted as $\mathcal{A}^*$, is a mapping from $\Re^{M \times H\times W}$ to $\Re^{C\times H \times W}$:
\begin{equation}\label{eq:conv_At}
    \mathcal{A}^* \boldsymbol{y} := \sum_{m=1}^M \big(\boldsymbol{\alpha}_{m1} * \boldsymbol{\eta}_m, \ldots, \boldsymbol{\alpha}_{mC} * \boldsymbol{\eta}_m\big) \in \Re^{C \times H \times W},
\end{equation}
where
$\boldsymbol{y} = (\boldsymbol{\eta}_1, \ldots, \boldsymbol{\eta}_M)\in \Re^{M \times H \times W}$ and $\boldsymbol{\eta}_m \in \Re^{H \times W}$.

We first study under what conditions of kernel $\boldsymbol{A}$ such that $\mathcal{A}$ and $\mathcal{A}^*$ are isometric. Starting from Definition~\ref{def:isometry} we arrive at the following theorem (see appendix for the proof).

\begin{theorem}\label{thm:convolution_isometry}
Given a convolution kernel $\boldsymbol{A}\in \Re^{M\times C \times k \times k}$ in \eqref{eq:kernel}, the operator $\mathcal{A}$ is an isometry if and only if
\begin{equation}\label{eq:isometry-conv}
    \sum_{m=1}^M \boldsymbol{\alpha}_{mc} \star \boldsymbol{\alpha}_{mc'} =
    \begin{cases}
        \mathbf{\delta} & ~\text{if}~c = c',\\
        \textbf{0} & ~\text{otherwise},\\
    \end{cases}
\end{equation}
and the operator $\mathcal{A}^*$ is an isometry if and only if
\begin{equation}\label{eq:coisometry-conv}
    \sum_{c=1}^C \boldsymbol{\alpha}_{mc} \star \boldsymbol{\alpha}_{m'c} =
    \begin{cases}
    \mathbf{\delta} & ~\text{if}~m = m',\\
    \textbf{0} & ~\text{otherwise}.\\
    \end{cases}
\end{equation}
In above, $\mathbf{\delta}$ is the Kronecker delta function defined on $\bb Z \times \bb Z$ that takes value $1$ at coordinate $(0, 0)$ and $0$ otherwise.
\end{theorem}

For the case of $M = C$, the isometry of $\mathcal{A}$ is equivalent to the isometry of $\mathcal{A}^*$\footnote{When $M=C$, the isometry of $\mathcal{A}$ implies that it is surjective, therefore unitary, hence $\mathcal{A}^*$ is also an isometry.}.
We refer to $\boldsymbol{A}$ as an \emph{orthogonal convolution} kernel if $M=C$ and either $\mathcal{A}$ or $\mathcal{A}^*$ is isometric.

On the other hand, if $M \ne C$ then $\mathcal{A}$ and $\mathcal{A}^*$ cannot be both isometric.
Nonetheless, we can still enforce isometry of $\mathcal{A}$ (in the case $M\ge C$) or the isometry of $\mathcal{A}^*$ (in the case $M\le C$) by the conditions \eqref{eq:isometry-conv} and \eqref{eq:coisometry-conv}, respectively.

\vspace{1em}
\begin{example}[Delta kernel]
\label{eg:delta}
    We refer to a convolution kernel $\boldsymbol{A} \in \Re^{ C \times M \times k \times k}$ as the (Kronecker) delta kernel, denoted as $\mathbf{\delta}^{C\times M \times k \times k}$, if its entries indexed by $(i, i, k_0, k_0)$ (assuming $k = 2k_0+1$) for all $i \in \{1, \ldots, \min(C, M)\}$ are set to $1$ and all other entries are set to zero.
    If $M \ge C$ (resp., $M \le C$), then the operator $\mathcal{A}$ (resp., $\mathcal{A}^*$) associated with a delta kernel is an isometry.
    Moreover, if $M = C$ then a delta kernel is orthogonal.
\end{example}

\textbf{Isometry at initialization.} As shown in Example~\ref{eg:delta}, we may initialize all convolution kernels to be isometric by setting them to be the delta kernel, which we refer to as the Delta initialization.
Other orthogonal initialization, such as the Delta Orthogonal initialization \cite{xiao2018dynamical}, are also plausible choices. We empirically find that Delta initialization works well and often outperforms Delta Orthogonal initialization. Such initialization is also commonly used in initializing Recurrent Neural Networks~\cite{le2015simple}.

\textbf{Isometry during training.} Isometry at initialization does not guarantee that isometry will be preserved through the training process.
In addition to Delta initialization, we enforce isometry by penalizing the difference between the left and right hand sides of \eqref{eq:isometry-conv} (or \eqref{eq:coisometry-conv}).
As we shown below, this can be easily implemented via modern deep learning packages.

In particular, given an input that contains a batch of $N$ multi-channel signals $\boldsymbol{X} := (\boldsymbol{x}_1, \ldots, \boldsymbol{x}_N) \in \Re^{N \times C \times H \times W}$ where each $\boldsymbol{x}_n \in \Re^{C \times H \times W}$ for each $n \in \{1, \ldots, N\}$, and a convolution kernel $\boldsymbol{A} \in \Re^{M\times C \times k \times k}$ defined as in \eqref{eq:kernel}, a convolution function in a typical deep learning framework, denoted as $\text{Conv}(\boldsymbol{A}, \boldsymbol{X})$, performs correlation of $\boldsymbol{A}$ on each of the $N$ signals and stack the resulting $N$ multi-channel signals to a tensor:
\begin{equation}\label{eq:conv-operator}
    \text{Conv}(\boldsymbol{A}, \boldsymbol{X}) = (\mathcal{A} \boldsymbol{x}_1, \ldots, \mathcal{A} \boldsymbol{x}_N).
\end{equation}
By replacing $\boldsymbol{X}$ in \eqref{eq:conv-operator} with the kernel $\boldsymbol{A}$ itself and using \eqref{eq:conv_A}, we obtain
\begin{equation}
\begin{split}\label{eq:convAA}
    & \text{Conv}(\boldsymbol{A}, \boldsymbol{A}) \\
    =& \; \big( \mathcal{A}(\boldsymbol{\alpha}_{11}, \ldots, \boldsymbol{\alpha}_{1C}), \ldots, \mathcal{A}(\boldsymbol{\alpha}_{M1}, \ldots, \boldsymbol{\alpha}_{MC}) \big)\\
    =&
    \begin{small}
    \begin{pmatrix}
        \sum_{c=1}^C \boldsymbol{\alpha}_{1c}\star \boldsymbol{\alpha}_{1c} & \cdots & \sum_{c=1}^C \boldsymbol{\alpha}_{Mc}\star \boldsymbol{\alpha}_{1c}\\
        \vdots & \ddots & \vdots \\
        \sum_{c=1}^C \boldsymbol{\alpha}_{1c}\star \boldsymbol{\alpha}_{Mc} & \cdots & \sum_{c=1}^C \boldsymbol{\alpha}_{Mc}\star \boldsymbol{\alpha}_{Mc}
    \end{pmatrix}.
    \end{small}
\end{split}
\end{equation}
Comparing \eqref{eq:convAA} with the condition in \eqref{eq:coisometry-conv}, we see that the operator $\mathcal{A}^*$ is an isometry if and only if Conv$(\boldsymbol{A}, \boldsymbol{A})$ is a delta kernel $\mathbf{\delta}^{M \times M \times k \times k}$.
Similarly, $\mathcal{A}$ is an isometry if and only if Conv$(\boldsymbol{A}^\top, \boldsymbol{A}^\top) = \mathbf{\delta}^{C \times C \times k \times k}$, where $\boldsymbol{A}^\top$ denotes a kernel with the first two dimension of $\boldsymbol{A}$ transposed.
From these results, we can enforce isometry by adding one of the following regularization terms to the objective function:
\begin{align}\label{eq:orth-reg}
    L(\boldsymbol{A}) &= \frac{\gamma}{2} \|\text{Conv}(\boldsymbol{A}, \boldsymbol{A}) - \mathbf{\delta}^{M \times M \times k \times k}\|_F^2, ~~\text{or}\\
    L(\boldsymbol{A}^\top) &= \frac{\gamma}{2} \|\text{Conv}(\boldsymbol{A}^\top, \boldsymbol{A}^\top) - \mathbf{\delta}^{C \times C \times k \times k}\|_F^2,
\end{align}
where $\gamma$ is a regularization coefficient. We use $L(\boldsymbol{A})$ when $C > M$ and $L(\boldsymbol{A}^\top)$ otherwise.

\subsection{Isometry in Nonlinear Activation}
\label{sec:activation}

The rectified linear unit (ReLU) is one of the most popular activation functions for deep learning applications \cite{nair2010rectified}.
Defined entry-wise on the input as $\phi(x) = \max(0, x)$, ReLU is an identity map for nonnegative input, but completely removes all negative component.
For an input signal that is normalized to zero mean and unit variance, ReLU is far from being an isometry.
How can we develop an isometric nonlinear activation layer?

\begin{figure}[t]
\centering
    \begin{minipage}[t]{0.3\linewidth}
    \centering
    \includegraphics[width=\linewidth]{figs/srelu/linear.pdf}
    \scriptsize{Purely isometric \\ Linear \\ \vspace{1em} a) Linear}
    \end{minipage}
    ~
    \begin{minipage}[t]{0.3\linewidth}
    \centering
    \includegraphics[width=\linewidth]{figs/srelu/srelu.pdf}
    \scriptsize{Partially isometric \\ Weakly nonlinear \\ \vspace{1em} b) SReLU}
    \end{minipage}
    ~
    \begin{minipage}[t]{0.3\linewidth}
    \centering
    \includegraphics[width=\linewidth]{figs/srelu/relu.pdf}
    \scriptsize{Not isometric \\ Strongly nonlinear \\ \vspace{1em} c) ReLU}
    \end{minipage}
    \caption{Illustration of isometry in nonlinear activation}
    \label{fig:srelu}
\end{figure}

Unfortunately, isometry is innately at odds with nonlinearity.
By Mazur-Ulam theorem, any surjective isometry between two normed spaces over $\Re$ is linear up to a translation \cite{nica2013mazur}.
Therefore, isometry can never be achieved if any nonlinearity is present. On the other hand, a purely isometric network without nonlinearity is not interesting either, since the entire network becomes a linear transformation which has very limited learning ability.
The design of the nonlinear activation function requires striking a good {\em trade-off between isometry and nonlinearity}.

We advocate using a Shifted ReLU (SReLU) for achieving such a trade-off. SReLU is obtained by interpolating regular ReLU (for obtaining nonlinearity) and the identity (for obtaining isometry).
Given a signal $\boldsymbol{y} \in \Re^{C\times H \times W}$ where $C$ is the number of channels, SReLU is defined point-wise as
\begin{equation}\label{eq:SReLU}
    \phi_{b}(y) = \max(y, b),
\end{equation}
where $b$ is a parameter that is shared for each channel of $\boldsymbol{y}$.
SReLU becomes the regular ReLU if $b = 0$ and the identity if $b \to -\infty$, therefore can be interpreted as a trade-off between nonlinearity and isometry  (see Figure~\ref{fig:srelu}).
In terms of backward dynamic (see \eqref{eq:backward}), the operator $\mathcal{D}$ associated with $\phi_{b}(\cdot)$ is given as $\mathcal{D}(x) = 1_{y \ge b} \cdot x$, which is an identity map if $y$ is larger than $b$.

In our experiments, we find it beneficial to initialize $b$ to be negative values, which indicates the importance of such a trade-off between ReLU and identity.
For simplicity, in all our experiments we initialize $b$ to be $-1$ and optimize it in the training process.

\subsection{Isometry in Residual Structure}
\label{sec:residual}

We may further improve the isometric property of each layer of ISONet with a residual component \cite{he2016deep}, which we call a residual ISONet (R-ISONet).
The core idea is to introduce a skip connection to ISONet, so that the network learns a \emph{residual} component that is added onto the path of signal propagation (see Figure~\ref{fig:net}).
Such a network architecture is automatically near isometric if the residual component is small enough relative to the identity shortcut~\cite{tarnowski2019dynamical}.
Motivated by such an observation, we add a scalar multiplier $s$ at the end of each residual branch and initialize $s$ to be zero ($s$ is shared within each channel).
Similar design ideas have previously been explored to improve ResNet~\cite{szegedy2017inception,goyal2017accurate, zhang2019fixup}.

Since skip connection helps to enforce isometry, R-ISONet exhibits better performance than ISONet.
For image classification on ImageNet, R-ISONet obtains almost as good performance as regular ResNet, even though R-ISONet has no normalization layers.
Since no normalization layer is required, R-ISONet is particularly suited for applications where normalization layers are ineffective.
As we show in the next section, R-ISONet obtains better performance than regular ResNet for object detection tasks on COCO dataset.

\section{Experiments}
\label{sec:experiments}

\subsection{Experimental Setup}

\textbf{Setup.}
We test the performance of ISONets on the ILSVRC-2012 image classification dataset~\cite{deng2009imagenet, russakovsky2015imagenet}. The training set contains $\app1.28$ million images from $1,\!000$ categories. Evaluation is performed on the validation set which contains $50,\!000$ images. 1-Crop, Top-1 accuracy is reported.

\textbf{Network architectures.}
The design of network architecture follows from the design of ResNet \cite{he2016deep, goyal2017accurate}, with the only difference being that we change the residual blocks to the block shown in Figure~\ref{fig:net}. In particular, the number of channels and layers in each stage remain the same as ResNet. We find that (R-)ISONet tends to overfit since BatchNorm is removed. Therefore we add one dropout layer right before the final classifier. The dropout probability is $0.4$ for R-ISONet and $0.1$ for ISONet.

\textbf{Implementation details.}
The hyperparameter settings follow from prior works~\cite{he2016deep, goyal2017accurate}. All models are trained using SGD with weight decay $0.0001$, momentum $0.9$ and mini-batch size $256$. The initial learning rate is 0.1 for
R-ISONet and 0.02 for ISONet. The models are trained for 100 epochs with the learning rate subsequently divided by a factor of $10$ at the $30$th, $60$th and $90$th epochs. To stabilize training in the early stage, we perform linear scheduled warmup~\cite{goyal2017accurate} for the first $5$ epochs for both our methods and baseline methods. The orthogonal regularization coefficient $\gamma$ in \eqref{eq:orth-reg} (when used) is 0.0001 except for ISONet-101 where a stronger regularization (0.0003) is needed. For other baselines without normalization or skip connections, we choose the largest learning rate that makes them converge.

\subsection{ISONet: Isometry Enables Training Vanilla Nets}

In this section, we verify through extensive experiments that isometry is a central principle for training deep networks.
In particular, we show that more than $100$-layer ISONet can be effectively trained with neither normalization layers nor skip connections.
In addition, we demonstrate through ablation study that all of the isometric components in ISONet are necessary for obtaining good performance.

\begin{table}[t]
    \centering
    \small
    \setlength{\tabcolsep}{4pt}
    \renewcommand{\arraystretch}{1.5}
    \begin{tabular}{c|c|ccc|c}
    & Method & SReLU & {\begin{mytabular}[1.2]{c} Delta \\ Init. \end{mytabular}} & \begin{mytabular}[1.2]{c} Ortho. \\ Reg. \end{mytabular} & \begin{mytabular}[1.2]{c} Top-1 \\ Acc. (\%) \end{mytabular}  \\
    \hline
    \hline
    (a) & \demph{ResNet} & & & & \demph{73.29}\\
    \hline
    (b) & Vanilla & & & & 63.09 \\
    \hline
    (c) & & & \checkmark & \checkmark & 46.83 \\
    (d) & & \checkmark & & & 67.35 \\
    (e) & & \checkmark & & \checkmark & 68.50 \\
    (f) & & \checkmark & \checkmark & & 68.55 \\
    \hline
    (g) & ISONet & \checkmark & \checkmark & \checkmark & \textbf{70.45} \\
    \end{tabular}
    \caption{\textbf{Isometric learning (with SReLU, Delta initialization and orthogonal regularization) enables training ISONets on ImageNet without BatchNorm and skip connection. }
    (a) Regular $34$-layer ResNet with ReLU activation, Kaiming initialization and no orthogonal regularization. (b) Same as ResNet but without BatchNorm and skip connection. (g) Our ISONet with the same backbone as Vanilla. (c-f) Ablation for (g) that shows that the combination of all three isometric learning components is necessary for effectively training ISONet.
    }
    \label{tab:abl_iso}
\end{table}

\begin{figure*}[t]
\centering
\subfigure{\includegraphics[width=0.24\linewidth]{figs/bias/bias1.pdf}}
\subfigure{\includegraphics[width=0.24\linewidth]{figs/bias/bias2.pdf}}
\subfigure{\includegraphics[width=0.24\linewidth]{figs/bias/bias3.pdf}}
\subfigure{\includegraphics[width=0.24\linewidth]{figs/bias/bias4.pdf}}
\vspace{-8pt}
\caption{Histogram of the parameter $b$ in SReLU computed for each of the four stages of a 34-layer ISONet.
}\label{fig:bias}
\end{figure*}

\begin{table}[t]
    \centering
    \small
    \setlength{\tabcolsep}{3.5pt}
    \renewcommand{\arraystretch}{1.6}
    \begin{tabular}{c|ccccc}
    Orthogonal coefficient $\gamma$ & 0 & $1e^{-5}$ & $3e^{-5}$ & $1e^{-4}$ & $3e^{-4}$ \\
    \hline
    \hline
    Top-1 Accuracy (\%) & 68.55 & 70.08 & 70.44 & \textbf{70.45} & 69.32 \\
    \end{tabular}
    \caption{\textbf{Effect of orthogonal coefficients for 34-layer ISONets}. The performance is not sensitive to the coefficient in a wide range.}
    \label{tab:abl_coeff}
    \vspace{-0.5em}
\end{table}

\begin{table}[t]
    \centering
    \small
    \setlength{\tabcolsep}{4pt}
    \renewcommand{\arraystretch}{1.4}
    \begin{tabular}{c|c|c|c|c}
    & Activation & Adaptive $b$? & ~~Initial $b$~~ & Top-1 Acc. (\%) \\
    \hline
    \hline
    (a) & ReLU & & - & 46.83 \\
    \hline
    (b) & \multirow{3}{*}{SReLU} & & -1.0 & 70.27 \\
    (c) & & \checkmark & -1.0 & 70.45 \\
    (d) & & \checkmark & -2.0 & 70.04 \\
    \end{tabular}
    \caption{\textbf{Effect of parameter $b$ in SReLU for 34-layer ISONets. }}
    \label{tab:abl_act}
    \vspace{-0.5em}
\end{table}

\textbf{Ablation study of ISONet components.}
To demonstrate how each isometric component affects the trainability of deep neural networks, we perform experiments with $34$-layer vanilla neural networks and report results in Table~\ref{tab:abl_iso}.
In Table~\ref{tab:abl_iso} (a), we show the result of regular ResNet.
If the normalization layers and skip connections in ResNet are removed, the network performance drops by more than $10$ points as shown in Table~\ref{tab:abl_iso} (b).
This shows a plain vanilla network cannot be easily optimized.
In fact, $34$-layer is the deepest vanilla network we found trainable in our experiments.
In contrast, ISONet with neither normalization layers nor skip connections can be effectively trained and obtain competitive $70.45\%$ Top-1 accuracy, as shown in Table~\ref{tab:abl_iso} (g).
To our best knowledge, this is the only vanilla network that achieves $>70\%$ accuracy on ImageNet.

We further show that all the three isometric components (i.e., Delta initialization, orthogonal regularization and SReLU) are indispensable for training deep ISONet.
Firstly, enforcing isometry in convolution alone produces \raise.17ex\hbox{$\scriptstyle\sim$}$46\%$ Top-1 Accuracy, which is far lower than ISONet (see Table~\ref{tab:abl_iso} (c)).
The reason is that orthogonal convolution is not suitable for vanilla networks with ReLU activation functions, since the energy may decrease exponentially through multiple layers of the network.
Secondly, enforcing isometry in activation functions alone cannot achieve good performance, as shown in Table~\ref{tab:abl_iso} (d).
Thirdly, for convolution layers, isometry must be enforced both at initialization as well as during training.
This can be seen in Table~\ref{tab:abl_iso} (e) and (f), where the performance of ISONet drops by $2$ points if either Delta initialization or orthogonal regularization is removed.
This demonstrates the critical role of isometry during the entire training process.

\textbf{Analysis for orthogonal regularization.}
To understand the effect of the orthogonal regularization term in \eqref{eq:orth-reg}, we conduct experiments with different values of the regularization coefficient $\gamma$ and report Top-1 accuracy in Table~\ref{tab:abl_coeff}.
It is shown that our method is not sensitive to the choice of $\gamma$ in a wide range (from $1e^{-5}$ to $1e^{-4}$), whereas too large ($3e^{-4}$) or no regularization will hurt the performance.

\textbf{Analysis for SReLU.}
To verify the design idea of SReLU, we perform an ablation study on the effect of parameter $b$ in \eqref{eq:SReLU} and report the results in Table~\ref{tab:abl_act}.
Comparing (b) and (c), a trainable threshold in SReLU is better than a fixed one.
Comparing (c) and (d), changing the initial value of $b$ from $-1$ to $-2$ reduces the accuracy by less than $0.5$ percent.
We also plot the histogram of the learned $b$ in SReLU computed within each of the four stages (corresponding to spatial resolutions $\{56, 28, 14, 7\}$) of ISONet and report the results in Figure~\ref{fig:bias}.
In all cases, the histogram is bell-shaped and concentrates near $-1$. This explains why a fixed $b= -1$ produces a good performance as shown in Table~\ref{tab:abl_act} (b).

\textbf{Training very deep ISONets. }
To demonstrate the effectiveness of isometric principle for training deep networks, we report performance of ISONet with increasing number of layers and report the results in
Table~\ref{tab:depth}.
Firstly, we observe that a standard vanilla network, either without (Table~\ref{tab:depth} (c)) or with (Table~\ref{tab:depth} (d)) BatchNorm, cannot converge when networks are extremely deep.
In contrast, ISONet maintains a competitive accuracy when the depth is beyond $100$ layers.
Adding dropout to ISONet further improves the performance by \raise.17ex\hbox{$\scriptstyle\sim$} 0.5 points. On the other hand, the performance ISONet is not as good as ResNet.
In the following, we show that when combined with skip connection, the performance of our R-ISONet is on par with ResNet.

\begin{table}[t]
    \centering
    \small
    \setlength{\tabcolsep}{4pt}
    \renewcommand{\arraystretch}{1.5}
    \begin{tabular}{c|c|c|cccc}
    & \multirow{2}{*}{Method} & \multirow{2}{*}{\footnotesize{Dropout}} & \multicolumn{4}{|c}{Top-1 Accuracy (\%)} \\
    & & & d=18 & d=34 & d=50 & d=101 \\
    \hline
    \hline
    (a) & \multirow{2}{*}{\demph{ResNet}} &  & \demph{69.67} & \demph{73.29} & \demph{76.60} & \demph{77.37} \\
    (b) & & \demph{\checkmark} & \demph{68.91} & \demph{73.35} & \demph{76.40} & \demph{77.99} \\
    \hline
    (c) & Vanilla & & 65.67 & 63.09 & \footnotesize{N/A} & \footnotesize{N/A} \\
    (d) & Vanilla+BN & & 68.98 & 69.43 & 70.00 & \footnotesize{N/A} \\
    \hline
    (e) & \multirow{2}{*}{ISONet} & & 67.94 & 70.45 & 70.73 & 70.38 \\
    (f) & & \checkmark & 68.10 & 70.90 & 71.20 & 71.01
    \end{tabular}
    \caption{\centering\textbf{ISONet with varying depth $d \in \{18, 34, 50, 101\}$ can be effectively trained on ImageNet.} (a, b) Regular ResNet with or without dropout.
    (c) Same as ResNet but without BatchNorm and skip connection. (d) Same as ResNet but without skip connection.
    (e, f) Our ISONet with the same backbone as Vanilla.
    N/A stands for not converging or obtaining an accuracy lower than $40$.
    }
    \vspace{-0.5em}
    \label{tab:depth}
\end{table}

\begin{table}[h!]
    \centering
    \small
    \setlength{\tabcolsep}{4pt}
    \renewcommand{\arraystretch}{1.5}
    \begin{tabular}{c|c|c|cccc}
    & \multirow{2}{*}{Method} & \multirow{2}{*}{\footnotesize{Dropout}} & \multicolumn{4}{|c}{Top-1 Accuracy (\%)} \\
    & & & d=18 & d=34 & d=50 & d=101 \\
    \hline
    \hline
    (a) & \multirow{2}{*}{ResNet} & & 69.67 & 73.29 & 76.60 & 77.37 \\
    (b) & & \checkmark & 68.91 & 73.35 & 76.40 & 77.99 \\
    \hline
    (c) & R-Vanilla & & 65.66 & \footnotesize{N/A} & \footnotesize{N/A} & \footnotesize{N/A} \\
    \hline
    (d) & Fixup & & 68.63 & 71.28 & 72.40 & 73.18 \\
    (e) & Fixup++ & & 67.37 & 72.56 & 76.00 & 76.17 \\
    \hline
    (f) & \multirow{2}{*}{R-ISONet} & & 69.06 & 72.17 & 74.20 & 75.44 \\
    (g) & & \checkmark & 69.17 & 73.43 & 76.18 & 77.08
    \end{tabular}
    \caption{\textbf{Performance of R-ISONet with varying depth $d \in \{18, 34, 50, 101\}$ on ImageNet. }
    All networks do not have BatchNorm except (a, b).
    (c) ResNet with all BatchNorm removed. (d, e) Same as (c), but with Fixup initialization \cite{zhang2019fixup}. (f, g) Our R-ISONet performs comparable to ResNet in (a, b) and better than all networks in (c, d, e).
    N/A stands for not converging or obtaining an accuracy lower than $40$. Fixup results are reproduced from running the code in \cite{zhang2019fixup} with $100$ training epochs. Fixup++ denotes Fixup with Mixup data augmentation.
    }
    \label{tab:depth2}
\end{table}

\subsection{R-ISONet: a Practical Network without Normalization Layers}
\label{sec:exp-R-ISONet}

In this section, we evaluate the performance of R-ISONet for ImageNet classification and show that it is comparable to regular ResNet, despite the fact that R-ISONet does not contain BatchNorm. The results are reported in Table~\ref{tab:depth2}.

From Table~\ref{tab:depth2} (c), ResNet without BatchNorm cannot converge for $d\ge34$. In contrast, our R-ISONet (Table~\ref{tab:depth2} (f)) can be effectively trained for varying network depth.
In addition, R-ISONet is better than previous methods without BatchNorm such as Fixup~\cite{zhang2019fixup} (Table~\ref{tab:depth2} (d)).

Since R-ISONet is prone to overfitting due to the lack of BatchNorm, we add dropout layer right before the final classifier and report the results in Table~\ref{tab:depth2} (g).
The results show that R-ISONet is comparable to ResNet with dropout (Table~\ref{tab:depth2} (b)) and is better than Fixup with Mixup regularization~\cite{zhang2017mixup}.

The superior performance that R-ISONet obtains even without any normalization layers makes it a favorable choice for tasks where statistics in BatchNorm cannot be precisely estimated.
In Section~\ref{sec:exp-object}, we demonstrate such advantage of R-ISONet for object detection and instance segmentation tasks.

\subsection{Transfer Learning on Object Detection}
\label{sec:exp-object}

We further evaluate our method for object detection and instance segmentation tasks on COCO dataset~\cite{lin2014microsoft}.
The dataset contains $115$k training images and $5$k validation images. We use the mAP (over different IoU threshold) metric to evaluate the performance of the model, and the ``dilated conv5'' variant of Faster RCNN~\cite{ren2015faster} as our detector.
The backbone of the network includes all the convolution layers in R-ISONet.
The RoIAlign operator is added on top of the last convolution with the last stage having convolution with stride $1$ and dilation $2$.
Two fully-connected layers are then used to predict the bounding box scores and regression output using this feature.
This design follows prior works~\cite{wu2019detectron2, dai2017deformable, he2019momentum}.
The hyper-parameter settings of both ResNet and R-ISONet follows detectron2~\cite{wu2019detectron2}. We train our model for 90k iterations. The learning rate is initialized to be $0.01$ and divide by $10$ at 60k and 80k iterations.
The training is performed on 8 GPUs, each of which holds 2 images. To keep a fair comparison with standard protocols in object detection, dropout is not added.

The results are reported in Table~\ref{tab:abl_det}.
Although the classification accuracy of our R-ISONet is lower than that of ResNet with the same depth, the detection and instance segmentation performance of R-ISONet is better. This demonstrates that our model has better feature transfer abilities and can mitigate the disadvantages introduced by BatchNorm.

\begin{table}[t]
    \centering
    \small
    \setlength{\tabcolsep}{5pt}
    \renewcommand{\arraystretch}{1.5}
    \begin{tabular}{c|c|cc}
    & Methods & mAP$^{\text{bbox}}$ & mAP$^{\text{mask}}$ \\
    \hline
    \hline
    \multirow{2}{*}{34 layer} & ResNet & 35.0 & 32.2 \\
    & R-ISONet & \textbf{36.2} & \textbf{33.0} \\
    \hline
    \multirow{2}{*}{50 layer} & ResNet  & 37.0 & 33.9 \\
    & R-ISONet & \textbf{37.3} & \textbf{34.4} \\
    \end{tabular}
    \caption{\textbf{Performance of R-ISONet with varying depth for object detection on COCO dataset. } R-ISONet \textit{outperforms} standard ResNet, indicating that our model has better transfer ability.
    }
    \label{tab:abl_det}
    \vspace{-0.5em}
\end{table}

\section{Related Work and Discussion}
\label{sec:previouswork}

In this paper, we contend that isometry is a central design principle that enables effective training and learning of deep neural networks. In the literature, numerous instantiations of this principle have been explicitly or implicitly suggested, studied, and exploited, often as an additional heuristic or regularization, for  improving existing network training or performance. Now supported by the strong empirical evidence presented earlier about isometric learning, we reexamine these ideas in the literature from the unified perspective of isometry.

Arguably, the notion of isometry was first explored in the context of weight initialization (Section \ref{sec:c_init}) and then in diverse contexts such as weight regularization (Section~\ref{sec:c_regu}), design of nonlinear activation (Section~\ref{sec:c_activation}), and training techniques for residual networks (Section~\ref{sec:c_residual}).

\subsection{In Context of Weight Initialization}
\label{sec:c_init}

Early works on weight initialization are based on the principle that the variance of the signal maintains a constant level as it propagates forward or backward through multiple layers \cite{lecun2012efficient,glorot2010understanding}.
A popular method that provides such a guarantee is the Kaiming Initialization \cite{he2015delving} which derives a proper scaled Gaussian initialization for weights in vanilla convolutional networks with ReLU activation functions.
Derivations for general activation functions beyond ReLU is more difficult, but nonetheless can be achieved by working in the regime where the network is infinitely wide using mean-field theory \cite{poole2016exponential,schoenholz2016deep}.
From the perspective of isometry, the aforementioned works guarantee that the average of the squared singular values of the input-output Jacobian matrix is close to $1$.
However, this \emph{does not mean that the Jacobian is an isometry}, which requires that \emph{all} the singular values concentrate at $1$.
In fact, with the common practice of Gaussian weight initialization, isometry can never be achieved regardless of the choice of activation functions \cite{pennington2018emergence}.

To address such an issue, orthogonal weight initialization has been extensively studied in the past few years.
For linear networks with arbitrary depth, orthogonality of each composing layers trivially leads to an isometry of the input-output Jacobian, and the benefit over Gaussian initialization in terms of training efficiency has been empirically observed \cite{saxe2013exact} and theoretically justified \cite{hu2020provable}.
For deep nonlinear networks, isometry may also be achieved if the network works in a local regime that the nonlinearity becomes approximately linear \cite{pennington2018emergence}, and empirical good performance of orthogonal initialization is also observed \cite{mishkin2015all} particularly when combined with proper scaling.
This line of work culminates at a recent work \cite{xiao2018dynamical} which shows that orthogonal initialization enables effective training of ConvNets with $10,\!000$ layers.
Nonetheless, the performance of the network in \citet{xiao2018dynamical} is far below the state-of-the-art networks, perhaps due to the fact that isometry beyond the initialization point is not guaranteed.

The Delta initialization adopted in our method is a particular case of orthogonal initialization for convolution kernels. Despite the existence of other orthogonal initialization \cite{xiao2018dynamical}, we advocate the Delta initialization due to its simplicity and good empirical performance in our evaluation for visual recognition tasks.

\subsection{In Context of Weight Regularization}
\label{sec:c_regu}

A plethora of works has explored the idea of regulating the convolution operators in a ConvNet to be orthogonal in the training process \cite{harandi2016generalized,jia2017improving,cisse2017parseval,bansal2018can,zhang2019approximated,li2019efficient,huang2020controllable}.
It is also found that orthogonal regularization helps with training GANs \cite{brock2018large,liu2019oogan} and RNNs \cite{arjovsky2016unitary,lezcano2019cheap}.
Despite widely observed performance improvement, the explanation for the effectiveness of orthogonal regularization is rather diverse: it has been justified from the perspective of
alleviating gradient vanishing or exploding \cite{xie2017all},
stabilizing distribution of activation over layers~\cite{huang2018orthogonal},
improving generalization~\cite{jia2019orthogonal},
and so on.

In our framework, the benefit of orthogonal regularization lies in that it enforces the network to be close to an isometry.
To make our argument precise, we first distinguish between the two related concepts of orthogonal weights and orthogonal convolution.
When enforcing orthogonality of convolution kernels in convolutional neural networks, all of the works mentioned above are based on flattening a 4D kernel into a 2D matrix and imposing orthogonality on the matrix.
We refer to such a method as imposing \emph{weight orthogonality}, which is not the same as the orthogonality of convolution as discussed in Section~\ref{sec:convolution}.
For example, in \citet{bansal2018can} a kernel $\boldsymbol{A} \in \Re^{M \times C \times k \times k}$ is reshaped into a matrix of shape $M \times (C \times k \times k)$ and is enforced to be row-orthogonal. That is,
\begin{equation}\label{eq:coisometry-conv-flatten}
    \sum_{c=1}^C \langle \text{vec}(\boldsymbol{\alpha}_{mc}), \text{vec}(\boldsymbol{\alpha}_{m'c}) \rangle = \begin{cases}
    1 & ~\text{if}~m = m',\\
    0 & ~\text{otherwise},\\
    \end{cases}
\end{equation}
where vec$(\boldsymbol{\alpha}_{mc})$ denotes a vector of length $k \times k$ obtained from flattening $\boldsymbol{\alpha}_{mc}$.
Comparing \eqref{eq:coisometry-conv-flatten} with \eqref{eq:coisometry-conv}, it is clear that such a regularization is necessary, but not sufficient, for the operator $\mathcal{A}^*$ to be isometric.
Therefore imposing weight orthogonality as in previous works can be interpreted as partially enforcing network isometry as well.

To the best of our knowledge, the only works that have derived orthogonality for convolution kernels are \citet{xiao2018dynamical} and \citet{wang2019orthogonal}.
\citet{xiao2018dynamical} provides a means of constructing orthogonal convolution kernels, but it can only represent a subset of all orthogonal convolutions \cite{li2019preventing}.
In addition, while the method in \citet{xiao2018dynamical} can be used to generate a random orthogonal convolution, it does not provide a means of enforcing orthogonality in the training process.
We become aware of a very recent work \cite{wang2019orthogonal} during the preparation of this paper, which derives the notion of orthogonal convolution that is equivalent to ours.
However, the objective of our work is different from theirs. Instead of improving an existing network by adding an additional regularization, our work aims to show that isometry is the most central principle to design effective deep vanilla networks. And as shown in Table~\ref{tab:abl_iso}, using orthogonal regularization alone is not enough for training deep vanilla networks.
Meanwhile, the derivation in \cite{wang2019orthogonal} is based on expressing convolution as matrix-vector multiplication using a doubly block-Toeplitz matrix and using the notion of orthogonality for matrices.
Such a derivation is limited to 2D convolution while generalization to higher dimensional convolution may become very cumbersome.
Moreover, the definition is restricted to discrete time signals and does not adapt to continuous time signals.
In contrast, our definition of orthogonal convolution can be extended for higher-dimensional convolution and for continuous times signals by properly redefining the operator ``*'' in \eqref{eq:def-conv}, therefore may bear broader interests.

\subsection{In Context of Nonlinear Activation}
\label{sec:c_activation}

Many of the important variants of nonlinear activation functions developed over the past few years can be interpreted as obtaining closer proximity to isometry.
Early works on improving ReLU such as Leaky ReLU \cite{maas2013rectifier}, Parametric ReLU \cite{he2015delving} and Randomized ReLU \cite{xu2015empirical}, which are generically defined as
\begin{equation}\label{eq:LReLU}
    f(x) =
    \begin{cases}
      x & \text{if}~~ x \ge 0, \\
      \alpha x  & \text{otherwise},
    \end{cases}
\end{equation}
with the parameter $\alpha >0$ being fixed, learned and randomly chosen, respectively.
Motivation for this family of activation functions comes from the dead neuron issue of ReLU, which states that negative input values give rise to zero gradients that prevents effective training of the network.

The Exponential Linear Unit (ELU) \cite{clevert2015fast} and its variants \cite{trottier2017parametric,klambauer2017self}, defined as
\begin{equation}\label{eq:ELU}
    f(x) =
    \begin{cases}
      \gamma x & \text{if}~~ x \ge 0, \\
      \alpha (\exp(x/\beta)-1)  & \text{otherwise,}
    \end{cases}
\end{equation}
with $\alpha, \beta$ and $\gamma$ being trainable or fixed parameters, is another family of activation functions that may solve the dead neuron issue.
Another important motivation for ELU is that it pushes the mean of the activation closer to zero, therefore alleviates the issue with regular ReLU that the bias of the input signal is always shifted towards positivity.
Such benefit is further confirmed by a theoretical study \cite{klambauer2017self} which shows that a zero mean and unit variance signal converges towards being zero mean and unit variance after propagating over multiple layers.
Empirically, a very recent work \cite{huang2019sndcnn} demonstrates that a 50-layer neural network without batch normalization and skip connection may be effectively trained for speech recognition tasks by using ELU instead of ReLU.

The SReLU that we adopt in the isometric network is closely related to the ELU family.
To the best of our knowledge, a variant of SReLU (in which the parameter is fixed and not optimizable) first appears in \citet{clevert2015fast}, the paper that proposes ELU, where it is shown to outperform ReLU and Leaky ReLU by a large margin and has similar performance as ELU.
Besides, SReLU is also adopted in \cite{qiu2018frelu} as an important baseline.
Recently, SReLU is rediscovered in \cite{xiang2017effects,singh2019filter} for alleviating the bias shifting issue.

We argue that aside from the benefits claimed by the original papers, the activation function families \eqref{eq:LReLU}, \eqref{eq:ELU} as well as the SReLU are advantageous over the regular ReLU as they bring the network closer to an isometry.
At an intuitive level, both \eqref{eq:LReLU} (with $0 < \alpha < 1$), \eqref{eq:ELU} (with $\gamma=1$ and $0 < \alpha < 1$ and $\beta=1$) and SReLU lie between ReLU and identity, therefore may be considered as a tradeoff between obtaining nonlinearity and isometry.
From a theoretical perspective, the analysis from \cite{pennington2018emergence} reveals that SReLU combined with random orthogonal initialization can achieve dynamic isometry, while regular ReLU cannot.

\subsection{In Context of Residual Learning}
\label{sec:c_residual}

Ever since the inception of residual learning \cite{he2016deep}, there has been no lack of effort in further enhancing its performance. In particular, there is a line of work that demonstrates the benefit of reducing the energy on the residual branch of a ResNet. This is studied from different perspectives such as dynamic isometry~\cite{tarnowski2019dynamical,qiu2018frelu}, training failure modes~\cite{taki2017deep,hanin2018start,balduzzi2017shattered}. Such design is also empirically applied in different ways such as BatchNorm initialization~\cite{goyal2017accurate}, adding a small scalar in the residual branch~\cite{zagoruyko2017diracnets,szegedy2017inception}, and 0-initialized convolution~\cite{zhang2019fixup}, all of which brings the network closer to an isometry. With a careful design along this line of study \cite{zhang2019fixup}, deep ResNet can be effectively trained to obtain competitive performance on visual recognition tasks, even without the help of BatchNorm layers. This line of research as well as our simple R-ISONet verifies the effectiveness of isometry in the design the network structures of residual learning.

\section{Conclusion}
In this paper, we have demonstrated through a principle-guided design and strong empirical evidence why isometry is likely to be the main key property that enables effective learning of deep networks and ensures high performance on real visual recognition tasks. With this design principle, one may achieve competitive performance with much simplified networks and eased training. We also argue that isometric learning provides a unified principle that helps explain numerous  ideas, heuristics and regularizations scattered in the literature that exploit this property and are found effective. We believe that the isometry principle may help people design or discover new simple network operators and architectures with much-improved performance in the future.

\section*{Acknowledgements}
The authors acknowledge support from Tsinghua-Berkeley Shenzhen Institute Research Fund.
Haozhi is supported in part by DARPA Machine Common Sense. Xiaolong is supported in part by DARPA Learning with Less Labels.
We thank Yaodong Yu and Yichao Zhou for insightful discussions on orthogonal convolution. We would also like to thank the members of BAIR for fruitful discussions and comments.

\clearpage