EFFICIENT GENERALIZED SPHERICAL CNNS
Year: 2021

Abstract:
Many problems across computer vision and the natural sciences require the analysis of spherical data, for which representations may be learned efficiently by encoding equivariance to rotational symmetries. We present a generalized spherical CNN framework that encompasses various existing approaches and allows them to be leveraged alongside each other. The only existing non-linear spherical CNN layer that is strictly equivariant has complexity OpC2L5q, where C is a measure of representational capacity and L the spherical harmonic bandlimit. Such a high computational cost often prohibits the use of strictly equivariant spherical CNNs. We develop two new strictly equivariant layers with reduced complexity OpCL4q and OpCL3 logLq, making larger, more expressive models computationally feasible. Moreover, we adopt efficient sampling theory to achieve further computational savings. We show that these developments allow the construction of more expressive hybrid models that achieve state-of-the-art accuracy and parameter efficiency on spherical benchmark problems.

text: Many problems across computer vision and the natural sciences require the analysis of spherical data, for which representations may be learned efficiently by encoding equivariance to rotational symmetries. We present a generalized spherical CNN framework that encompasses various existing approaches and allows them to be leveraged alongside each other. The only existing non-linear spherical CNN layer that is strictly equivariant has complexity OpC2L5q, where C is a measure of representational capacity and L the spherical harmonic bandlimit. Such a high computational cost often prohibits the use of strictly equivariant spherical CNNs. We develop two new strictly equivariant layers with reduced complexity OpCL4q and OpCL3 logLq, making larger, more expressive models computationally feasible. Moreover, we adopt efficient sampling theory to achieve further computational savings. We show that these developments allow the construction of more expressive hybrid models that achieve state-of-the-art accuracy and parameter efficiency on spherical benchmark problems.

heading: 1 INTRODUCTION
text: Many fields involve data that live inherently on spherical manifolds, e.g. 360˝ photo and video content in virtual reality and computer vision, the cosmic microwave background radiation from the Big Bang in cosmology, topographic and gravitational maps in planetary sciences, and molecular shape orientations in molecular chemistry, to name just a few. Convolutional neural networks (CNNs) have been tremendously effective for data defined on Euclidean domains, such as the 1D line, 2D plane, or nD volumes, thanks in part to their translation invariance properties. However, these techniques are not effective for data defined on spherical manifolds, which have a very different geometric structure to Euclidean spaces (see Appendix A). To transfer the remarkable success of deep learning to data defined on spherical domains, deep learning techniques defined inherently on the sphere are required. Recently, a number of spherical CNN constructions have been proposed.
Existing CNN constructions on the sphere fall broadly into three categories: fully real (i.e. pixel) space approaches (e.g. Boomsma & Frellsen, 2017; Jiang et al., 2019; Perraudin et al., 2019; Cohen et al., 2019); combined real and harmonic space approaches (Cohen et al., 2018; Esteves et al., 2018; 2020); and fully harmonic space approaches (Kondor et al., 2018). Real space approaches can often be computed efficiently but they necessarily provide an approximate representation of spherical signals and the connection to the underlying continuous symmetries of the sphere is lost. Consequently, such approaches cannot fully capture rotational equivariance. Other constructions take a combined real and harmonic space approach (Cohen et al., 2018; Esteves et al., 2018; 2020), where sampling theorems (Driscoll & Healy, 1994; Kostelec & Rockmore, 2008) are exploited to connect with underlying continuous signal representations to capture the continuous symmetries of the sphere. However, in these approaches non-linear activation functions are computed pointwise in real space, which induces aliasing errors that break strict rotational equivariance. Fully harmonic space spherical CNNs have been constructed by Kondor et al. (2018). A continual connection with underlying continuous signal representations is captured by using harmonic signal representations throughout. Consequently, this is the only approach exhibiting strict rotational equivariance. However, strict equivariance comes at great computational cost, which can often prohibit usage.
˚Corresponding author: jason.mcewen@kagenova.com
In this article we present a generalized framework for CNNs on the sphere (and rotation group), which encompasses and builds on the influential approaches of Cohen et al. (2018), Esteves et al. (2018) and Kondor et al. (2018) and allows them to be leveraged alongside each other. We adopt a harmonic signal representation in order to retain the connection with underlying continuous representations and thus capture all symmetries and geometric properties of the sphere. We construct new fully harmonic (non-linear) spherical layers that are strictly rotationally equivariant, are parameterefficient, and dramatically reduce computational cost compared to similar approaches. This is achieved by a channel-wise structure, constrained generalized convolutions, and an optimized degree mixing set determined by a minimum spanning tree. Furthermore, we adopt efficient sampling theorems on the sphere (McEwen & Wiaux, 2011) and rotation group (McEwen et al., 2015a) to improve efficiency compared to the sampling theorems used in existing approaches (Driscoll & Healy, 1994; Kostelec & Rockmore, 2008). We demonstrate state-of-the-art performance on all spherical benchmark problems considered, both in terms of accuracy and parameter efficiency.

heading: 2 GENERALIZED SPHERICAL CNNS
text: We first overview the theoretical underpinnings of the spherical CNN frameworks introduced by Cohen et al. (2018), Esteves et al. (2018), and Kondor et al. (2018), which make a connection to underlying continuous signals through harmonic representations. For more in-depth treatments of the underlying harmonic analysis we recommend Esteves (2020), Kennedy & Sadeghi (2013) and Gallier & Quaintance (2019). We then present a generalized spherical layer in which these and other existing frameworks are encompassed, allowing existing frameworks to be easily integrated and leveraged alongside each other in hybrid networks.
Throughout the following we consider a network composed of S rotationally equivariant layers Ap1q, ....,ApSq, where the i-th layer Apiq maps an input activation f pi´1q P Hpi´1q onto an output activation f piq PHpiq. We focus on the case where the network input space Hp0q consists of spherical signals (but note that input signals on the rotation group may also be considered).

heading: 2.1 SIGNALS ON THE SPHERE AND ROTATION GROUP
text: Let L2pΩq denote the space of square-integrable functions over domain Ω. A signal f P L2pΩq on the sphere (Ω “ S2) or rotation group (Ω “ SOp3q) can be rotated by ρ P SOp3q by defining the action of rotation on signals by Rρfpωq “ fpρ´1ωq for ω P Ω. An operator A : L2pΩ1q Ñ L2pΩ2q, where Ω1,Ω2 P tS2,SOp3qu, is then equivariant to rotations if RρpApfqq “ ApRρfq for all f P L2pΩ1q and ρ P SOp3q, i.e. rotating the function before application of the operator is equivalent to application of the operator first, followed by a rotation.
A spherical signal f P L2pS2q admits a harmonic representation pf̂0, f̂1, ..., q where f̂ ` P C2``1 are the harmonic coefficients given by the inner product xf, Y `my, where Y `m are the spherical harmonic functions of degree ` and order |m| ď `. Likewise a signal f P L2pSOp3qq on the rotation group admits a harmonic representation pf̂0, f̂1, ...qwhere f̂ ` P Cp2``1qˆp2``1q are the harmonic coefficients with pm,nq-th entry xf,D`mny for integers |m|, |n| ď `, where D` : SOp3q Ñ Cp2``1qˆp2``1q is the unique 2` ` 1 dimensional irreducible group representation of SOp3q on Cp2``1q. The rotation f ÞÑ Rρf of a signal f P L2pΩq can be described in harmonic space by f̂ ` ÞÑ D`pρqf̂ `. A signal on the sphere or rotation group is said to be bandlimited at L if, respectively, xf, Y `my “ 0 or xf,D`mny “ 0 for ` ě L. Furthermore, a signal on the rotation group is said to be azimuthally bandlimited at N if, additionally, xf,D`mny “ 0 for |n| ě N . Bandlimited signals therefore admit finite harmonic representations pf̂0, ..., f̂L´1q. In practice real-world signals can be accurately represented by suitably bandlimited signals; henceforth, we assume signals are bandlimited.

heading: 2.2 CONVOLUTION ON THE SPHERE AND ROTATION GROUP
text: A standard definition of convolution between two signals f, ψ P L2pΩq on either the sphere (Ω “ S2) or rotation group (Ω “ SOp3q) is given by
pf ‹ ψqpρq “ xf,Rρψy “ ż
Ω
dµpωqfpωqψ˚pρ´1ωq, (1)
where dµpωq denotes the Haar measure on Ω and ¨˚ complex conjugation (e.g. Wandelt & Górski, 2001; McEwen et al., 2007; 2013; 2015b; 2018; Cohen et al., 2018; Esteves et al., 2018). In particular, the convolution satisfies
ppRρfq ‹ ψqpρ1q “ xRρf,Rρ1ψy “ xf,Rρ´1ρ1ψy “ pRρpf ‹ ψqqpρ1q (2) and is therefore a rotationally equivariant linear operation, which we shall denote by Lpψq. The convolution of bandlimited signals can be computed exactly and efficiently in harmonic space as {pf ‹ ψq` “ f̂ `ψ̂`˚, ` “ 0, ..., L´ 1, (3) which for each degree ` is a vector outer product for signals on the sphere and a matrix product for signals on the rotation group (see Appendix B for further details). Convolving in this manner results in signals on the rotation group (for inputs on both the sphere and rotation group). However, in the spherical case, if the filter is invariant to azimuthal rotations the resultant convolved signal may be interpreted as a signal on the sphere (see Appendix B).

heading: 2.3 GENERALIZED SIGNAL REPRESENTATIONS
text: The harmonic representations and convolutions described above have proven useful for describing rotationally equivariant linear operators Lpψq. Cohen et al. (2018) and Esteves et al. (2018) define spherical CNNs that sequentially apply this operator, with intermediary representations taking the form of signals on SOp3q and S2 respectively. Alternatively, for intermediary representations we now consider the more general space of signals introduced by Kondor et al. (2018), to which the aforementioned notions of rotation and convolution naturally extend.
In describing the generalization we first note from Section 2.1 that all bandlimited signals on the sphere and rotation group can be represented as a set of variable length vectors of the form
f “ tf̂ `t P C2``1 : ` “ 0, .., L´ 1; t “ 1, ..., τ `fu, (4) where τ `f “ 1 for signals on the sphere and τ `f “ minp2`` 1, 2N ´ 1q for signals on the rotation group. The generalization is to let FL be the space of all such sets of variable length vectors, with τf unrestricted. This more general space contains the spaces of bandlimited signals on the sphere and rotation group as strict subspaces. For a generalized signal f P FL we adopt the terminology of Kondor et al. (2018) by referring to f̂ `t as the t-th fragment of degree ` and to τf “ pτ0f , ..., τ L´1 f q, specifying the number of fragments for each degree, as the type of f . The action of rotations upon FL can be naturally extended from their action upon L2pS2q and L2pSOp3qq. For f P FL we define the rotation operator f ÞÑ Rρf by f̂ `t ÞÑ D`pρqf̂ `t , allowing us to extend the usual notion of equivariance to operators A : FL Ñ FL.

heading: 2.4 GENERALIZED CONVOLUTIONS
text: The convolution described by Equation 1 provides a learnable linear operator Lpψq that satisfies the desired property of equivariance. Nevertheless, given the generalized interpretation of signals on S2 and SOp3q as signals in FL, the notion of convolution can also be generalized (Kondor et al., 2018). In order to linearly and equivariantly transform a signal f P FL of type τf into a new signal f ˚ψ P FL of any desired type τpf˚ψq, we may specify a filter ψ “ tψ̂` P Cτ ` fˆτ ` pf˚ψq : ` “ 0, ..., L ´ 1u, which in general is not an element of FL, and define a transformation f ÞÑ f ˚ ψ by
{pf ˚ ψq ` t “ τ`f ÿ
t1“1 f̂ `t1 ψ̂ ` t,t1 , ` “ 0, ..., L´ 1; t “ 1, ..., τ `pf˚ψq. (5)
The degree-` fragments of the transformed signal pf ˚ ψq are simply linear combinations of the degree-` fragments of f , with no mixing between degrees. Equation 3 shows that this is precisely the form taken by convolution on the sphere and rotation group. In fact Kondor & Trivedi (2018) show that all equivariant linear operations take this general form; the standard notion of convolution is just a special case. One benefit to the generalized notion is that the filter ψ is not forced to occupy the same domain as the signal f , thus allowing control over the type τpf˚ψq of the transformed signal. We use LpψqG to denote this generalized convolutional operator.

heading: 2.5 NON-LINEAR ACTIVATION OPERATORS
text: For FL to be a useful representational space, it must be possible to not only linearly but also nonlinearly transform its elements in an equivariant manner. However, equivariance and non-linearity is not enough. Equivariant linear operators cannot mix information corresponding to different degrees. Therefore it is of crucial importance that degree mixing is achieved by the non-linear operator.

heading: 2.5.1 POINTWISE ACTIVATIONS
text: When the type τf of f P FL permits an interpretation as a signal on S2 or SOp3q we may perform an inverse harmonic transform to map the function onto a sample-based representation (e.g. Driscoll & Healy, 1994; McEwen & Wiaux, 2011; Kostelec & Rockmore, 2008; McEwen et al., 2015a). A non-linear function σ : C Ñ C may then be applied pointwise, i.e. separately to each sample, before performing a harmonic transform to return to a representation in FL. We denote the corresponding non-linear operator as Nσpfq “ FpσpF´1pfqqq, where F represents the harmonic (i.e. Fourier) transform on S2 or SOp3q. The computational cost of the non-linear operator is dominated by the harmonic transforms. While costly, fast algorithms can be leveraged (see Appendix A). While inverse and forward harmonic transforms on S2 or SOp3q that are based on a sampling theory maintain perfect equivariance for bandlimited signals, the pointwise application of σ (most commonly ReLU) is only equivariant in the continuous limit L Ñ 8. For any finite bandlimit L, aliasing effects are introduced such that equivariance becomes approximate only, as shown by the following experiments.
We consider 100 random rotations ρ P SOp3q, for each of 100 random signal-filter pairs pf, ψq, and compute the mean equivariance error dpApRρfq,RρpAfqq for operator A, where dpf, gq “ }f ´ g}{}f} is the relative distance between signals. For convolutions the equivariance error is 4.4ˆ 10´7 for signals on S2 and 5.3 ˆ 10´7 for signals on SOp3q (achieving floating point precision). By comparison the equivariance error for a pointwise ReLU is 0.34 for signals on S2 and 0.37 for signals on SOp3q. Only approximate equivariance is achieved for the ReLU since the nonlinear operation spreads information to higher degrees that are not captured at the original bandlimit, resulting in aliasing. To demonstrate this point we reduce aliasing error by oversampling the realspace signal. When oversampling by 2ˆ or 8ˆ for signals on SOp3q the equivariance error of the ReLU is reduced to 0.10 and 0.01, respectively. See Appendix D for further experimental details.
Despite the high cost of repeated harmonic transforms and imperfect equivariance, this is nevertheless the approach adopted by Cohen et al. (2018), Esteves et al. (2018) and others, who find empirically that such models maintain a reasonable degree of equivariance.

heading: 2.5.2 TENSOR-PRODUCT ACTIVATIONS
text: In order to define a strictly equivariant non-linear operation that can be applied to a signal f P FL of any type τf the decomposability of tensor products between group representations may be leveraged, as first considered by Thomas et al. (2018) in the context of neural networks.
Given two group representations D`1 and D`2 of SOp3q on C2`1`1 and C2`2`1 respectively, the tensor-product group representation D`1 bD`2 of SOp3q on C2`1`1 b C2`2`1 is defined such that pD`1 bD`2qpρq “ D`1pρq bD`2pρq for all ρ P SOp3q. Decomposing D`1 bD`2 into a direct sum of irreducible group representations then constitutes finding a change of basis for C2`1`1 b C2`2`1 such that pD`1 bD`2qpρq is block diagonal, where for each ` there is a block equal to D`pρq. The necessary change of basis for û`1 b v̂`2 P C2`1`1 b C2`2`1 is given by
pû`1 b v̂`2q`m “ `1 ÿ
m1“´`1
`2 ÿ
m2“´`2
C`1,`2,`m1,m2,mû `1 m1 v̂ `2 m2 , (6)
where C`1,`2,`m1,m2,m P C denote Clebsch-Gordan coefficients whose symmetry properties are such that pû`1 b v̂`2q`m is non-zero only for |`1´ `2| ď ` ď `1` `2. The use of Equation 6 arises naturally in quantum mechanics when coupling angular momenta.
This property is useful since if f̂ `1 P C2`1`1 and f̂ `2 P C2`2`1 are two fragments that are equivariant with respect to rotations of the network input, then a rotation of ρ applied to the network input results
in f̂ `1 b f̂ `2 transforming as
rf̂ `1 b f̂ `2s` ÞÑ rpD`1pρqf̂ `1q b pD`2pρqf̂ `2qs` (7) “ rpD`1 bD`2qpρqpf̂ `1 b f̂ `2qs` (8) “ D`pρqrf̂ `1 b f̂ `2s`, (9)
where the final equality follows by block diagonality with respect to the chosen basis. Therefore, if fragments f̂ `1 and f̂ `2 are equivariant with repsect to rotations of the network input, then so is the fragment pC`1,`2,`qJpf̂ `1 b f̂ `2q P C2``1, where we have written Equation 6 more compactly. We now describe how Kondor et al. (2018) use this fact to define equivariant non-linear transformations of elements in FL.
A signal f “ tf̂ `t P C2``1 : ` “ 0, .., L ´ 1; t “ 1, ..., τ `fu P FL may be equivariantly and non-linearly transformed by an operator Nb : FL Ñ FL defined as Nbpfq “ tpC`1,`2,`qJpf̂ `1t1 b f̂ `2 t2 q : ` “ 0, ..., L´1; p`1, `2q P P ` L; t1 “ 0, ..., τ `1 f ; t2 “ 0, ..., τ
`2 f u, (10)
where for each degree ` P t0, ..., L´ 1u the set
P`L “ tp`1, `2q P t0, ..., L´ 1u2 : |`1 ´ `2| ď ` ď `1 ` `2u (11) is defined in order to avoid the computation of trivially equivariant all-zero fragments. We make the dependence on P`L explicit since we redefine it in Section 3.
Unlike the pointwise activations discussed in the previous section this operator is strictly equivariant, with a mean relative equivariance error of 5.0 ˆ 10´7 (see Appendix D). Note that applying this operator to signals on the sphere or rotation group results in generalized signals that are no longer on the sphere or rotation group. This is the rationale for the generalization to FL: to unlock the ability to introduce non-linearity in a strictly equivariant manner. Note, however, that g “ Nbpfq has type τg “ pτ0g , ..., τL´1g q where τ `g “ ř
p`1,`2qPP`L τ `1f τ `2 f and therefore application of this non-
linear operator results in a drastic expansion in representation size, which is problematic.

heading: 2.6 GENERALIZED SPHERICAL CNNS
text: Equipped with operators to both linearly and non-linearly transform elements of FL, with the latter also performing degree mixing, we may consider a network with representation spaces Hp0q “ ... “ HpSq “ FL. We consider the s-th layer of the network to take the form of a triple Apsq “ pL1,N ,L2q such that Apsqpf ps´1qq “ L2pN pL1pf ps´1qqqq, where L1,L2 : FL Ñ FL are linear operators and N : FL Ñ FL is a non-linear activation operator. The approaches of Cohen et al. (2018) and Esteves et al. (2018) are encompassed in this framework as Apsq “ pLpψq,Nσ, Iq, where I denotes the identity operator and filters ψ may be defined to encode real-space properties such as localization (see Appendix C). The framework of Kondor et al. (2018) is also encompassed as Apsq “ pI,Nb,LpψqG q. Here the generalized convolution L pψq G comes last to counteract the representation-expanding effect of the tensor-product activation and prevent it from compounding as signals pass through the network. Appendix E lends intuition regarding relationships that may be captured by tensor-product activations followed by generalized convolutions.
For any intermediary representation f piq P FL we may transition from equivariance with respect to the network input to invariance by discarding all but the scalar-valued fragments corresponding to ` “ 0 (equivalent to average pooling for signals on the sphere and rotation group). Finally, note that within this general framework we are free to consider hybrid approaches whereby layers proposed by Cohen et al. (2018); Esteves et al. (2018); Kondor et al. (2018) and others, and those presented subsequently, can be leveraged alongside each other within a single model.

heading: 3 EFFICIENT GENERALIZED SPHERICAL CNNS
text: Existing approaches to spherical convolutional layers that are encompassed within the above framework are computationally demanding. They require the evaluation of costly harmonic transforms
on the sphere and rotation group. Furthermore, the only strictly rotationally equivariant non-linear layer is that of Kondor et al. (2018), which has an even greater computational cost, scaling with the fifth power of bandlimit — thereby limiting spatial resolution — and quadratically with the number of fragments per degree — thereby limiting representational capacity. This often prohibits the use of strictly equivariant spherical networks.
In this section we introduce a channel-wise structure, constrained generalized convolutions, and an optimized degree mixing set in order to construct new strictly equivariant layers that exhibit much improved scaling properties and parameter efficiency. Furthermore, we adopt efficient sampling theory on the sphere and rotation group to achieve additional computational savings.

heading: 3.1 EFFICIENT GENERALIZED SPHERICAL LAYERS
text: For an activation f P FL the value τ̄f “ 1L řL´1 `“1 τ ` f represents a resolution-independent proxy for its representational capacity. Kondor et al. (2018) consider the separate fragments contained within f to subsume the traditional role of separate channels and therefore control the capacity of intermediary network representations through specification of τf . This is problematic because, whereas activation functions usually act on each channel separately and therefore have a cost that scales linearly with representational capacity (usually controlled by the number of channels), for the activation function Nb not only does the cost scale quadratically with representational capacity τ̄f , but so too does the size of Nbpfq. This feeds forward the quadratic dependence to the cost of, and number of parameters required by, the proceeding generalized convolution.
More specifically, note that computation of g “ Nbpfq requires the computation of řL´1 `“0 τ ` g fragments, where τ `g “ ř
p`1,`2qPP`L τ `1f τ `2 f . The size of P`L is OpL`q for each ` and therefore the ex-
panded representation has size řL´1 `“0 τ ` g , of order Opτ̄2fL3q. By exploiting the sparsity of ClebschGordan coefficients (C`1,`2,`m1,m2,m “ 0 if m1 `m2 ‰ m) each fragment pC `1,`2,`qJpf̂ `1 b f̂ `2q can be computed in Op`minp`1, `2qq. Hence, the total cost of computing all necessary fragments has complexity OpC2L5q, where C “ τ̄f captures representational capacity.

heading: 3.1.1 CHANNEL-WISE TENSOR-PRODUCT ACTIVATIONS
text: As is more standard for CNNs we maintain a separate channels axis, with network activations taking the form pf1, ..., fKq P FKL where fi P FL all share the same type τf . The non-linearity Nb may then be applied to each channel separately at a cost that is reduced byK-times relative to its application on a single channel with the same total number of fragments. This saving arises since for each ` we need only compute K ř
p`1,`2qPP`L τ `1f τ `2 f fragments rather than ř p`1,`2qPP`L pKτ `1f qpKτ `2 f q.
Figure 1 visualizes this reduction for the case K “ 3. Note, however, that for practical applications K „ 100 is more typical. The K-times reduction in cost is therefore substantial and allows for intermediary activations with orders of magnitude more representational capacity.
By introducing this multi-channel approach and using C “ K rather than C “ τ̄f to control representational capacity, we reduce the complexity of Nb with respect to representational capacity from OpC2q to OpCq.

heading: 3.1.2 CONSTRAINED GENERALIZED CONVOLUTION
text: Although much reduced, for a signal f P FKinL the channel-wise application of Nb still results in a drastically expanded representation g “ Nbpfq, to which a representation-contracting generalized convolution must be applied in order to project onto a new activation g1 “ LpψqG pgq P FKoutL of the desired type τg1 and number of channels Kout. However, under our multi-channel structure computational and parameter efficiency can be improved significantly by decomposing LpψqG into three separate linear operators, Lpψ1qG1 , L pψ2q G2 and L pψ3q G3 .
The first, Lpψ1qG1 , acts uniformly across channels, performing a linear projection down onto the desired type, and should be interpreted as a learned extension of Nb which undoes the drastic expansion. The second, Lpψ2qG2 , then acts channel-wise, taking linear combinations of the (contracted number of) fragments within each channel. The third, Lpψ3qG3 , acts across channels, taking linear combinations to learn new features. More concretely, the three filters are of the form ψ1 “ tψ̂`1 P C τ`gˆτ ` g1 : ` “ 0, ..., L´ 1u, ψ2 “ tψ̂`,k2 P C τ` g1ˆτ ` g1 : ` “ 0, ..., L ´ 1; k “ 1, ...,Kinu and ψ3 “ tψ̂`3 P CKinˆKout : ` “ 0, ..., L ´ 1u, rather than a single filter of the form ψ “ tψ̂` P CKinˆτ ` gˆKoutˆτ ` g1 : ` “ 0, ..., L´ 1u, leading to a large reduction in the number of parameters as τ `g is invariably very large.
By applying the first step uniformly across channels we minimize the parametric dependence on the expanded representation and allow new features to be subsequently learned much more efficiently. Together the second and third steps can be seen as analogous to depthwise separable convolutions often used in planar convolutional networks.

heading: 3.1.3 OPTIMIZED DEGREE MIXING SETS
text: We now consider approaches to reduce the OpL5q complexity with respect to spatial resolution L. In the definition of Nb each element of P`L independently defines an equivariant fragment. Therefore a restricted Nb in which only a subset of P`L is used for each degree ` still defines a strictly equivariant operator, while reducing computational complexity. In order to make savings whilst remaining at resolution L it is necessary to consider subsets of P`L that scale better than OpL2q. The challenge is to find such subsets that do not hamper the ability of the resulting operator to inject non-linearity and mix information corresponding to different degrees `.
Whilst various subsetting approaches are possible, the following argument motivates an approach that we have found to be particularly effective. If p`1, `3q P P`L, then representational space is designated to capture the relationship between `1 and `3-degree information. However, if resources have been designated already to capture the relationship between `1 and `2-degrees, as well as between `2 and `3-degrees, then some notion of the relationship between `1 and `3-degrees has been captured already. Consequently, it is unnecessary to designate further resources for this purpose.
More generally, consider the graphG`L “ pNL,P`Lqwith nodes NL “ t0, ..., L´1u and edges P`L. A restricted tensor-product activation can be constructed by using a subset of P`L that corresponds to a subgraph ofG`L. The subgraph ofG ` L captures some notion of the relationship between incoming `1 and `2-degree information if it contains a path between nodes `1 and `2. Therefore we are interested in subgraphs for which there exists a path between any two nodes if there exists such a path in the original graph, guaranteeing that any degree-mixing relationship captured by the original graph is also captured by the subgraph.
The smallest subgraph satisfying this property is a minimum spanning tree (MST) of G`L. The set of edges corresponding to any MST has at most L elements and we choose to consider its union with the set of loop-edges in G`L (of the form p`1, `1q), which proved particularly important for injecting non-linearity. We denote the resulting set as P̄`L and note that it satisfies |P̄`L| ď 2L. Therefore the tensor-product activation N̄b corresponding to Equation 11 with P`L replaced by P̄`L has reduced spatial complexity OpL4q. Given that many minimal spanning trees of the unweighted graphG`L exist for each `, we select the ones that minimize the cost of the resulting activation N̄b by assigning to each edge p`1, `2q inG`L a weight equal to the cost of computing pC`1,`2,`qJpf̂ `1b f̂ `2q and selecting the MST of the weighted graph.
An example of P`L and a MST subset P̄`L is shown in Figure 2, where the dashed line in Figure 2c shows the general form of the MST. Using this as a principled starting point we consider the further reduced MST (RMST) subset P̃`L corresponding to centering the MST at the edge p`, `q and retaining only the edges that fall a distance of 2i away on the dotted line for some i P N. We use Ñb to denote the corresponding operator and note that it has further reduced spatial complexity of OpL3 logLq. We demonstrate in Section 4 that networks that make use of the MST tensor-product activation achieve state-of-the-art performance. Replacing the MST with RMST activation results in a small but insignificant degradation in performance, which is offset by the reduced computational cost.

heading: 3.1.4 REDUCTION IN COMPUTATIONAL AND MEMORY FOOTPRINTS
text: The three modifications proposed in Sections 3.1.1 to 3.1.3 are motivated by improved scaling properties. Importantly, they also translate to large reductions in the computational and memory cost of strictly equivariant layers in practice, as detailed in Appendix F. Even at a modest bandlimit of L “ 64 and relatively small number of channels K “ 4, for example, the modifications together lead to a 51-times reduction in the number of flops required for computations and 16-times reduction in the amount of memory required to store representations, weights and gradients for training.

heading: 3.2 EFFICIENT SAMPLING THEORY
text: By adopting sampling theorems on the sphere we provide access to underlying continuous signal representations that fully capture the symmetries and geometric properties of the sphere, and allow standard convolutions to be computed exactly and efficiently through their harmonic representations, as discussed in greater detail in Appendices A and B. We adopt the efficient sampling theorems on sphere and rotation group of McEwen & Wiaux (2011) and McEwen et al. (2015a), respectively, which reduce the Nyquist rate by a factor of two compared to those of Driscoll & Healy (1994) and Kostelec & Rockmore (2008), which have been adopted in other spherical CNN constructions (e.g. Cohen et al., 2018; Kondor et al., 2018; Esteves et al., 2018; 2020). The sampling theorems adopted are equipped with fast algorithms to compute harmonic transforms, with complexity OpL3q for the sphere and OpL4q for the rotation group. When imposing an azimuthal bandlimitN ! L, the complexity of transforms on the rotation group can be reduced to OpNL3q, which we often exploit in our standard (non-generalized) convolutional layers.

heading: 4 EXPERIMENTS
text: Using our efficient generalized spherical CNN framework (implemented in the fourpiAI 1 code) we construct networks that we apply to numerous spherical benchmark problems. We achieve stateof-the-art performance, demonstrating enhanced equivariance without compromising representational capacity or parameter efficiency. In all experiments we use a similar architecture, consisting of 2–3 standard convolutional layers (e.g. S2 or SOp3q convolutions proceeded by ReLUs), followed by 2–3 of our efficient generalized layers. We adopt the efficient sampling theory described in Section 3.2 and encode localization of spatial filters as discussed in Appendix C. Full experimental details may be found in Appendix G.
1https://www.kagenova.com/products/fourpiAI/
Table 2: Test root mean squared (RMS) error for QM7 regression problem
RMS Params Montavon et al. (2012) 5.96 - Cohen et al. (2018) 8.47 1.4M Kondor et al. (2018) 7.97 ą1.1M Ours (MST) 3.16 337k Ours (RMST) 3.46 335k

heading: 4.1 ROTATED MNIST ON THE SPHERE
text: We consider the now standard benchmark problem of classifying MNIST digits projected onto the sphere. Three experimental modes NR/NR, R/R and NR/R are considered, indicating whether the training/test sets have been randomly rotated (R) or not (NR). Results are presented in Table 1, which shows that we closely match the prior state-of-the-art performance obtained by Esteves et al. (2020) on the NR/NR and R/R modes, whilst outperforming all previous spherical CNNs on the NR/R mode, demonstrating the increased degree of equivariance achieved by our model.
Results are shown for models using both the MST-based and RMST-based mixing sets within the tensor-product activation. The results obtained when using the full sets P`L are very similar to those obtained when using the MST-based sets (e.g. full sets achieved an accuracy of 99.39 for R/R).

heading: 4.2 ATOMIZATION ENERGY PREDICTION
text: We consider the problem of regressing the atomization energy of molecules given the molecule’s Coulomb matrix and the positions of the atoms in space, using the QM7 dataset (Blum & Reymond, 2009; Rupp et al., 2012). Results are presented in Table 2, which shows that we dramatically outperform other approaches, whilst using significantly fewer parameters.

heading: 4.3 3D SHAPE RETRIEVAL
text: We consider the 3D shape retrieval problem on the SHREC’17 (Savva et al., 2017) competition dataset, containing 51k 3D object meshes. We follow the pre-processing step of Cohen et al. (2018), where several spherical projections of each mesh are computed, and use the official SHREC’17 data splits. Results are presented in Table 3 for the standard SHREC precision and recall metrics, which shows that we achieve state-of-the-art performance compared to other spherical CNN approaches, achieving the highest three of five performance metrics, whilst using significantly fewer parameters.

heading: 5 CONCLUSIONS
text: We have presented a generalized framework for CNNs on the sphere that encompasses various existing approaches. We developed new efficient layers to be used as primary building blocks in this framework by introducing a channel-wise structure, constrained generalized convolutions, and optimized degree mixing sets determined by minimum spanning trees. These new efficient layers exhibit strict rotational equivariance, without compromising on representational capacity or parameter efficiency. When combined with the flexibility of the generalized framework to leverage the strengths of alternative layers, powerful hybrid models can be constructed. On all spherical benchmark problems considered we achieve state-of-the-art performance, both in terms of accuracy and parameter efficiency. In future work we intend to improve the scalability of our generalized framework further still. In particular, we plan to introduce additional highly scalable layers, for example by extending scattering transforms (Mallat, 2012) to the sphere, to further realize the potential of deep learning on a host of new applications where spherical data are prevalent.

heading: A REPRESENTATIONS OF SIGNALS ON THE SPHERE AND ROTATION GROUP
text: To provide further context for the discussion presented in the introduction and to elucidate the properties of different sampling theory on the sphere and rotation group, we concisely review representations of signals on the sphere and rotation group.
A.1 DISCRETIZATION
It is well-known that a completely regular point distribution on the sphere does in general not exist (e.g. Tegmark, 1996). Consequently, while a variety of spherical discretization schemes exists (e.g. icosahedron, HEALPix, graph, and other representations), it is not possible to discretize (i.e. to sample or pixelize) the sphere in a manner that is invariant to rotations, i.e. a discrete sampling of rotations of the samples on the sphere will in general not map onto the same set of sample positions. This differs to the Euclidean setting and has important implications when constructing convolution operators on the sphere, which clearly are a critical component of CNNs.
Since convolution operators are in general built using a translation operator – equivalently a rotation operator when on the sphere – it is thus not possible to construct a convolution operator directly on a discretized representation of the sphere that captures all of the symmetries of the underlying spherical manifold. While approximate discrete representations can be considered, and are nevertheless useful, such representations cannot capture all underlying spherical symmetries.
A.2 SAMPLING THEORY
Alternative representations, however, can capture all underlying spherical symmetries. Sampling theories on the sphere (e.g. Driscoll & Healy, 1994; McEwen & Wiaux, 2011) provide a mechanism to capture all information content of an underlying continuous function on the sphere from a finite set of samples (and similarly on the rotation group; Kostelec & Rockmore 2008; McEwen et al. 2015a). A sampling theory on the sphere is equivalent to a cubature (i.e. quadrature) rule for the exact integration of a bandlimited functions on the sphere. While optimal cubature on the sphere remains an open problem, the most efficient sampling theory on the sphere and rotation group is that developed by McEwen & Wiaux (2011) and McEwen et al. (2015a), respectively.
On a compact manifold like the sphere (and rotation group), harmonic (i.e. Fourier) space is discrete. Hence, a finite set of harmonic coefficients captures all information content of an underlying continuous bandlimited signal. Since such a representation provides access to the underlying continuous signal, all symmetries and geometric properties of the sphere are captured perfectly. Such representations have been employed extensively in the construction of wavelet transforms on the sphere, where the use of sampling theorems on the sphere and rotation group yield wavelet transforms of discretized continuous signals that are theoretically exact (e.g. McEwen et al., 2013; 2015b; 2018). Harmonic signal representations have also been exploited in spherical CNNs to access all underlying spherical symmetries and develop equivariance network layers (Cohen et al., 2018; Kondor et al., 2018; Esteves et al., 2018; 2020).
A.3 EXACT AND EFFICIENT COMPUTATION
Signals on the sphere f P L2pS2q may be decomposed into their harmonic representations as
fpωq “ 8 ÿ
`“0
ÿ̀
m“´` f `mY ` mpωq, (12)
where their spherical harmonic coefficients are given by
f `m “ xf, Y `my “ ż
S2 dµpωqfpωqY `˚m pωq, (13)
for ω P S2. Similarly, signals on the rotation group g P L2pSOp3qq may be decomposed into their harmonic representations as
gpρq “ 8 ÿ
`“0
2`` 1 8π2 ÿ̀
m“´`
ÿ̀
n“´` g`mnD `˚ mnpρq (14)
where their harmonic (Wigner) coefficients are given by
g`mn “ xg,D`˚mny “ ż
SOp3q dµpρqgpρqD`mnpρq, (15)
for ρ P SOp3q. Note that we adopt the convention where the conjugate of the Wigner D-function is used in Equation 14 since this leads to a convenient harmonic representation when considering convolutions (cf. McEwen et al., 2015a; 2018).
As mentioned above, sampling theory pertains to strategies to capture all of the information content of band limited signals from a finite set of samples. Since the harmonic space of the sphere and rotation group is discrete, this is equivalent to an exact quadrature rule for the computation of harmonic coefficients by Equation 13 and Equation 15 from sampled signals.
The canonical equiangular sampling theory on the sphere was that developed by Driscoll & Healy (1994), and subsequently extended to the rotation group by Kostelec & Rockmore (2008). More recently, novel sampling theorems on the sphere and rotation group were developed by McEwen & Wiaux (2011) and McEwen et al. (2015a), respectively, that reduce the Nyquist rate by a factor of two. Previous CNN constructions on the sphere (e.g. Cohen et al., 2018; Kondor et al., 2018; Esteves et al., 2018; 2020) have adopted the more well-known sampling theories of Driscoll & Healy (1994) and Kostelec & Rockmore (2008). In contrast, we adopt the more efficient sampling theories of McEwen & Wiaux (2011) and McEwen et al. (2015a) to provide additional efficiency savings, implemented in the open source ssht 2 and so3 3 software packages (we also make use of a TensorFlow implementation of these algorithms in our private tensossht 4 code). Note also that the sampling schemes associated with the theory of McEwen & Wiaux (2011) (and other minor variants implemented in ssht ) align more closely with the one-to-two aspect ratio of common spherical data, such as 360˝ photos and videos.
All of the sampling theories discussed are equipped with fast algorithms to compute harmonic transforms, with complexity OpL3q for transforms on the sphere (Driscoll & Healy, 1994; McEwen & Wiaux, 2011) and complexity OpL4q for transforms on the rotation group (Kostelec & Rockmore, 2008; McEwen et al., 2015a). Note that algorithms that achieve slightly lower complexity have been developed (Driscoll & Healy, 1994; Healy et al., 2003; Kostelec & Rockmore, 2008) but these are known to suffer stability issues (Healy et al., 2003; Kostelec & Rockmore, 2008). By imposing an azimuthally bandlimit N , where typically N ! L, the complexity of transforms on the rotation group can be reduced to OpNL3q (McEwen et al., 2015a), which we exploit in our networks. These fast algorithms to compute harmonic transforms on the sphere and rotation group can be leveraged to yield the exact and efficient computation of convolutions through their harmonic representations (see Appendix B). By computing convolutions in harmonic space, pixelization and quadrature errors are avoided and computational complexity is reduced to the cost of the respective harmonic transforms.

heading: B CONVOLUTION ON THE SPHERE AND ROTATION GROUP
text: For completeness we make explicit the standard (non-generalized) convolution operations on the sphere and rotation group that we adopt. The general form of convolution for signals f P L2pΩq either on the sphere (Ω “ S2) or rotation group (Ω “ SOp3q) is specified by Equation 1, with harmonic representation given by Equation 3. Here we provide specific expressions for the convolution for a variety of cases, describe the normalization constants that arise and may be absorbed into learnable filters, and derive the corresponding harmonic forms. In practice all convolutions are computed in harmonic space since the computation is then exact, avoiding pixelisation or quadrature errors, and efficient when fast algorithms to compute harmonic transforms are exploited (see Appendix A).
2http://www.spinsht.org/ 3http://www.sothree.org/ 4Available on request from https://www.kagenova.com/.
B.1 CONVOLUTION ON THE SPHERE
Given two spherical signals f, ψ P L2pS2q their convolution, which in general is a signal on the rotation group, may be decomposed as
pf ‹ ψqpρq “ xf,Rρψy (16)
“ ż
S2 dΩpωqfpωqψ˚pρ´1ωq (17)
“ ÿ
`m
ÿ
`1m1
ÿ
n
f `mD `1˚ m1npρqψ` 1˚ n
ż
S2 dΩpωqY `mpωqY `
1˚ m1 pωq (18)
“ ÿ
`m
ÿ
`1m1
ÿ
n
f `mD `1˚ m1npρqψ` 1˚ n δ``1δmm1 (19)
“ ÿ
`mn
` f `mψ `˚ n ˘ D`˚mnpρq, (20)
yielding harmonic coefficients
pf ˚ ψq`mn “ 8π2
2`` 1f ` mψ `˚ n . (21)
The constants 8π2{p2`` 1q may be absorbed into learnable parameters.
B.2 CONVOLUTION ON THE SPHERE WITH AXISYMMETRIC FILTERS
When convolving a spherical signal f P L2pS2q with an axisymmetric spherical filter ψ P L2pS2q that is invariant to azimuthal rotations, the resultant pf ‹ ψq may be interpreted as a signal on the sphere. To see this note that an axisymmetric filter ψ has harmonic coefficients ψ`n “ ψ`0δn0 that are non-zero only for m “ 0. Denoting rotations by their zyz-Euler angles ρ “ pα, β, γq and substituting into Equation 20 we see that the convolution may be decomposed as
pf ‹ ψqpα, β, γq “ ÿ
`mn
` f `mψ `˚ 0 δn0 ˘ D`˚mnpα, β, γq (22)
“ ÿ
`m
f `mψ `˚ 0 D `˚ m0pα, β, 0q (23)
“ ÿ
`m
f `mψ `˚ 0
c
4π
2`` 1Y ` mpβ, αq. (24)
We may therefore interpret pf ‹ ψq as a signal on the sphere with spherical harmonic coefficients
pf ‹ ψq`m “ c 4π
2`` 1f ` mψ `˚ 0 . (25)
The constants a 4π{p2`` 1q may be absorbed into learnable parameters.
B.3 CONVOLUTION ON THE ROTATION GROUP
Given two signals f, ψ P L2pSOp3qq on the rotation group their convolution may then be decomposed as
pf ‹ ψqpρq “ xf,Rρψy (26)
“ ż
SOp3q dµpρ1qfpρ1qψ˚pρ´1ρ1q (27)
“ ż
SOp3q dµpρ1q
„
ÿ
`
2`` 1 8π2 ÿ
mn
f `mnD `˚ mnpρ1q
„
ÿ
`1
2`1 ` 1 8π2 ÿ
m1n1
ψ` 1˚ m1n1D `1 m1n1pρ´1ρq 
(28)
“ ÿ
`
2`` 1 8π2 ÿ
mn
f `mn ÿ
`1
2`1 ` 1 8π2 ÿ
m1n1
ψ` 1˚ m1n1
ż
SOp3q dµpρ1qD`˚mnpρ1qD`
1
m1n1pρ´1ρ1q
(29)
“ ÿ
`
2`` 1 8π2 ÿ
mn
f `mn ÿ
`1
2`1 ` 1 8π2 ÿ
m1n1
ψ` 1˚ m1n1
ż
SOp3q dµpρ1qD`˚mnpρ1q
ÿ
k
D` 1˚ km1pρqD` 1 kn1pρ1q
(30)
“ ÿ
`
2`` 1 8π2 ÿ
mn
f `mn ÿ
`1
2`1 ` 1 8π2 ÿ
m1n1
ψ` 1˚ m1n1 ÿ
k
D` 1˚ km1pρq
8π2
2`` 1δ`` 1δmkδnn1 (31)
“ ÿ
`mm1
2`` 1 8π2 D`˚mm1pρq ˆ ÿ
n
f `mnψ `˚ m1n
˙
, (32)
where for Equation 30 we make use of the relation (e.g. Marinucci & Peccati, 2011; McEwen et al., 2018)
D`mnpρ´1ρ1q “ ÿ
k
D`˚kmpρqD`knpρ1q. (33)
This decomposition yields harmonic coefficients
pf ˚ ψq`mn “ ÿ
m1
f `mm1ψ `˚ nm1 . (34)

heading: C FILTERS ON THE SPHERE AND ROTATION GROUP
text: When defining filters we look to encode desirable real-space properties, such as locality and regularity. However, in practice considerable computation may be saved by defining the filters in harmonic space and saving the cost of harmonic transforming ahead of harmonic space convolutions. We describe here how filters motivated by their real space properties may be defined directly in harmonic space.
C.1 DIRAC DELTA FILTERS ON THE SPHERE
Spherical filters may be constructed as a weighted sum of Dirac delta functions on the sphere. This construction is useful as the harmonic representation has an analytic form that may be computed efficiently. Furthermore, various real space properties can be encoded through sensible placement of the Dirac delta functions.
The spherical Dirac delta function δω1 centered at ω1 “ pθ1, φ1q P S2 is defined as
δω1pωq “ 1
sin θ δRpcos θ ´ cos θ1qδRpφ´ φ1q, (35)
where δR is the familiar Dirac delta function on the reals centered at 0. The Dirac delta on the sphere may be represented in harmonic space by
pδω1q`m “ Y `˚m pω1q “ N `mP `mpcos θ1qe´imφ 1 , (36)
which follows form the sifting property of the Dirac delta, and where Y `m denote the spherical harmonic functions, P `mpxq are associated Legendre functions and
N `m “
d
2`` 1 4π pl ´mq! pl `mq! (37)
is a normalizing constant.
This representation may then be used to define a filter ψ P L2pS2q as a weighted sum of spherical Dirac delta functions, with weights wij assigned to Dirac delta functions centered at points tpθi, φjq : i “ 1, ..., Nθ; j “ 1, ..., Nφu. The associated harmonic space representation is given by
ψ`m “ ÿ
i,j
wijN ` mP ` mpcos θiqe´imφj (38)
“ ÿ
i
N `mP ` mpcos θiq
ÿ
j
wije ´imφj , (39)
where fast Fourier transforms may be leveraged to compute the inner sum if the Dirac deltas are spaced evenly azimuthally (e.g. if φj “ 2πj{Nφ). Alternative arbitrary samplings can of course be considered if useful for a problem at hand.
When defining filters in this manner one should be careful not to over-parametrize by assigning more weights than needed to define a filter at the harmonic bandlimit of the signal with which we wish to convolve. For example, if the filter is to be convolved with a signal bandlimited at L then a maximum of 2L ´ 1 Dirac deltas should be placed along each ring of constant θ. One may also choose to interpolate the weights from a smaller number of learnable parameters acting as anchor points, allowing higher resolution filters to be defined with fewer learnable parameters.
C.2 DIRAC DELTA FILTERS ON THE ROTATION GROUP
Similarly a Dirac delta function δρ1 on the rotation group SOp3q centered at position ρ1 “ pα1, β1, γ1q P SOp3q is defined as
δρ1pρq “ 1
sinβ δRpα´ α1qδRpcosβ ´ cosβ1qδRpγ ´ γ1q, (40)
with harmonic form pδρ1q`mn “ D`mnpρ1q “ e´imα 1 d`mnpβ1qe´inγ 1 , (41)
where d`mn are Wigner (small) d-matrices.
The filter ψ P L2pSOp3qq corresponding to a weighted sum of Dirac deltas with weights wijk assigned to Dirac delta functions centered at points tpαi, βj , γkq : i “ 1, ..., Nα; j “ 1, ..., Nβ ; k “ 1, ..., Nγu has harmonic form
ψ`mn “ ÿ
i,j,k
wijk e ´imαjd`mnpβiqe´inγk (42)
“ ÿ
j
d`mnpβiq ÿ
i
e´imαj ÿ
k
wijke ´inγk , (43)
where again fast Fourier transforms may be leveraged to compute the inner two sums assuming the Dirac deltas are spaced evenly in α and γ. The outer sums of Equation 39 and Equation 43 can also be computed by fast Fourier transforms by decomposing the Wigner d-matrices into their Fourier representation (cf. Trapani & Navaza, 2006; McEwen & Wiaux, 2011). One should again be careful not to over-parametrize.

heading: D EQUIVARIANCE TESTS
text: To test rotational equivariance of operators we consider Nf “ 100 random signals tfiu Nf i“1 in L2pΩ1q with harmonic coefficients sampled from the standard normal distribution and Nρ “ 100
random rotations tρju Nρ j“1 sampled uniformly on SOp3q. In order to measure the extent to which an operator A : L2pΩ1q Ñ L2pΩ2q is equivariant we evaluate the mean relative error
dpApRρjfiq,Rρj pAfiqq “ 1
Nf
1
Nρ
Nf ÿ
i“1
Nρ ÿ
j“1
}ApRρjfiq ´Rρj pAfiqq} }ApRρjfiq}
(44)
resulting from pre-rotation of the signal, followed by application of A, as opposed to post-rotation after application of A, where the operator norm } ¨ } is defined using the inner product x¨, ¨yL2pΩ2q. Table 4 presents the mean relative equivariance errors computed. We consider the three standard convolutions described in Appendix B (with a random filter ψi for each signal fi, generated in the same manner as fi), the pointwise ReLU activation described in Section 2.5.1 for signals on the sphere (Ω1 “ S2) and rotation group (Ω1 “ SOp3q), and the composition of tensor-product activation with a generalized convolution, described in Sections 2.5.2 and 2.4, respectively. We follow the tensor-product activation with a generalized convolution in order to project down onto the sphere to allow the same notion of error to be adopted as for the other operators. For consistency with the context in which we leverage these operators, all experiments are performed using singleprecision arithmetic.
We see that all three standard notions of convolution and the composition of the tensor-product activation and generalized convolution are all strictly equivariant to floating point machine precision, with errors on the order of 10´7. The pointwise ReLU operator is not strictly equivariant, with a mean relative error of 0.37 for signals on the rotation group and 0.34 for signals on the sphere. These errors reduce when the signals are oversampled before application of the ReLU, indicating that the error is due to aliasing induced by the spreading of information to higher degrees not captured at the original bandlimit. For example, for the pointwise ReLU operator on the rotation group oversampling by factors of 2ˆ, 4ˆ and 8ˆ results in a reduction in the mean relative equivariance error from 0.37 to 0.098, 0.032 and 0.0096, respectively.

heading: E CONNECTION BETWEEN THE TENSOR PRODUCT ACTIVATION AND POINTWISE SQUARING
text: To provide some intuition on the manner in which the tensor-product based activation introduces non-linearity into representations we describe its relationship to pointwise squaring for signals on the sphere. Here we consider the operator N : L2pS2q Ñ L2pS2q satisfying pN fqpxq “ f2pxq for all x P S2, which differs subtly to Nσ with σpxq “ x2 (using notation from Section 2.5.1), which corresponds to obtaining a sample-based representation at a finite bandlimit L and applying the squaring at the sample positions. For the special case σpxq “ x2 we can directly compute the harmonic representation corresponding to the equivariance-preserving continuous limit LÑ8.
Given a spherical signal f P L2pS2q with generalized representation f “ tf̂ `0 P C2``1 : ` “ 0, .., L´ 1u, the generalized representation of the signal f2 P L2pS2q is given as
f2 “ t ÿ
p`1,`2qPP`,L
pG`1,`2,`qJpf̂ `10 b f̂ `2 0 q : ` “ 0, 1, ..., L´ 1u, (45)
where G`1,`2,` P Cp2`1`1qˆp2`2`1qˆp2``1q are Gaunt coefficients defined as
Gj1j2j3m1m2m3 “ ż
S2 dµpωqY j1m1pωqY j2 m2pωqY j3˚ m3 pωq. (46)
Gaunt coefficients are related to the Clebsch-Gordan coefficients by
Gj1j2j3m1m2m3 “ w j1j2j3Cj1j2j3m1m2m3 , (47)
where wj1j2j3 “ p´1qm3 b p2j1`1qp2j2`1q 4πp2j3`1q C j1j2j3 0 0 0 . Therefore, the continuous squaring operation corresponds to passing f through a tensor-product activation Nb followed by a generalized convolution back down onto the sphere (single fragment per degree) with weight assigned to the p`1, `2q P P`L fragment in degree-` given by wl1l2l3 . This demonstrates that activations that are learnable within our framework can have very simple realspace interpretations. Even when confining outputs to the sphere we found it to be beneficial to allow the down-projection to be learnable rather than enforcing the weights given above for pointwise squaring. Learned activations will remain quadratic, however, given that output fragments are linear combinations of products between input fragments.

heading: F COMPARISON OF COMPUTATIONAL AND MEMORY FOOTPRINTS
text: We perform a quantitative analysis of the computational cost and memory requirements of our proposed layers, and comparisons to prior approaches, to demonstrate how the complexity savings made through our proposals in Section 3 translate to tangible efficiency improvements.
We consider the simplest comparison, between a multi-channel generalized signal f “ pf1, ..., fKq P FKL where each channel fi has type τfi “ p1, ..., 1q (and therefore corresponds to a signal on the sphere) and a uni-channel signal g P FL of type τg “ pK, ...,Kq. The setting corresponding to signal f captures the efficiency improvements of our proposed layers, while the setting corresponding to g represents the case without these improvements. Notice that the total number of fragments is the same in both f and g. We compare the number of floating point operations (FLOPs) and amount of memory required to perform a tensor-product based activation followed by a generalized convolution projecting back down onto a signal of the same type as the input. When applied to f , MST mixing sets P̄`L (Section 3.1.3) and constrained generalized convolutions (Section 3.1.2) are
used. When applied to g full mixing sets P`Land unconstrained generalized convolutions are used, as in Kondor et al. (2018). Considered are the costs for a single training instance (batch size of 1).
Figure 3 shows the computational costs and memory requirements, in terms of floating point operations and megabytes respectively, for K “ 4 and various spatial bandlimits L. We adopt the convention whereby complex addition and multiplication require 2 and 6 floating point operations respectively. At low bandlimits we see the saving arising from the channel-wise structure. Note that the saving illustrated here is relatively small since K “ 4 for these experiments (so that the L scaling is apparent), whereas in practice typically K „ 100. The saving then increases linearly in response to increases in the bandlimit of the input, as expected given the OpL5q and OpL4q spatial complexities. We see that even in this simple case, with a relatively small number of channels (K “ 4), both the computational and memory footprints are reduced by orders of magnitude. At a bandlimit of L “ 128 the computational cost is 101-times reduced and the memory requirement is 29-times reduced.

heading: G ADDITIONAL INFORMATION ON EXPERIMENTS
text: G.1 ROTATED MNIST ON THE SPHERE
For our MNIST experiments we used a hybrid model with the architecture shown in Figure 4. The first block includes a directional convolution on the sphere that lifts the spherical input (τ `
fp0q “ 1)
onto the rotation group (τ ` fp1q “ minp2` ` 1, 2N1 ´ 1q). The second block includes a convolution on the rotation group, hence its input and output both live on the rotation group. We then apply a restricted generalized convolution to map to type τ `
fp3q “ rτmax{
? 2`` 1s, where τmax“5. The same
type is used for the following three channel-wise tensor-product activations and two restricted generalized convolutions until the final restricted generalized convolution maps down to a rotationally invariant representation (τ `
fp5q “ δ`0). As is traditional in convolution networks we gradually de-
crease the resolution, with pL0, L1, L2, L3, L4, L5q “ p20, 10, 10, 6, 3, 1q, and increase the number of channels, with pK0,K1,K2,K3,K4,K5q “ p1, 20, 22, 24, 26, 28q. We proceed these convolutional layers with a single dense layer of size 30, sandwiched between two dropout layers (keep probability 0.5), and then fully connect to the output of size 10.
We train the network for 50 epochs on batches of size 32, using the Adam optimizer (Kingma & Ba, 2015) with a decaying learning rate starting at 0.001. For the restricted generalized convolutions we follow the approach of Kondor et al. (2018) by using L1 regularization (regularization strength 10´5) and applying a restricted batch normalization across fragments, where the fragments are only scaled by their average and not translated (to preserve equivariance).
G.2 ATOMIZATION ENERGY PREDICTION
When regressing the atomization energy of molecules there are two inputs to the model: the number of atoms of each element contained in the molecule; and spherical cross-sections of the potential energy around each atom. We adopt the high-level QM7-specific architecture of Cohen et al. (2018) which contains a spherical CNN as a sub-model, for which we substitute our own. This results in an overall model that is invariant to both rotations of the molecule around each constituent atom and to permutations of the ordering of the atoms.
The first (non-spherical) input is mapped onto a scalar output using a multi-layer perceptron (MLP) with three hidden layers of sizes 100, 100 and 10 (and ReLU activations). The second input, multiple spherical cross sections for each atom, are separately projected using a shared spherical CNN (of architecture described below) onto lower dimensional vectors of size 64. The mean vector is then taken across atoms (ensuring invariance w.r.t. permutations of the atoms) and mapped onto a scalar output using an MLP with a single hidden layer of size 512 (with a ReLU activation). The predicted energy is then taken to be the sum of the two scalar outputs.
As a starting point we train the first MLP to regress the atomization energies alone (achieving RMS „20), before pairing it with the spherical model (and its connected MLP). We then train the joint model for 60 epochs, again with the Adam optimizer, a decaying learning rate (starting at 2.5 ˆ 10´4), regularizing the efficient generalized layers with L2 regularization (strength 2.5ˆ 10´6) and batch sizes of 32.
For the spherical component we again adopt the convolutional architecture shown in Figure 4 except with one fewer efficient generalized layer. We use bandlimits of pL0, L1, L2, L3, L4q“p10, 6, 6, 3, 1q, channels of pK0,K1,K2,K3,K4q“p5, 16, 24, 32, 40q and τmax “ 6. One minor difference is that this time we include a skip connection between the ` “ 0 components of the fourth and fifth layer. We proceed the convolutional layers with two dense layers of size p256, 64q and use batch normalization between each layer.
G.3 3D SHAPE RETRIEVAL
To project the 3D meshes of the SHREC’17 data onto spherical representations (bandlimited at L “ 128) we adopt the preprocessing approach of Cohen et al. (2018) and augment the data with random rotations and translations.
We construct a model with an architecture that is again similar to that described in Appendix G.1 but with an additional axisymmetric convolutional layer prepended to the start of the network and one fewer efficient generalized layers. We use bandlimits pL0, L1, L2, L3, L4, L5q “ p128, 32, 16, 16, 6, 1q, channels pK0,K1,K2,K3,K4,K5q “ p6, 20, 30, 40, 60, 70q and τmax “ 6 for the efficient generalized layers. The convolutional layers are followed by a dense layer of size 128 which is fully connected to the output (of size 55).
We again train with the Adam optimizer, a decaying learning rate (starting at 5 ˆ 10´4) and batch sizes of 8, this time until performance on the validation set showed no improvement for at least 4 epochs (36 epochs in total). We perform batch normalization between convolutional layers and dropout preceding the dense layer. We regularize the efficient generalized layers with L2 regularization (strength 10´5). When testing our model we average the output probabilities over 15 augmentations of the data.

Review summaries:
The paper introduces a framework for computationally efficient and exactly rotation-equivariant spherical CNNs. The work most closely resembles the Fourier space method of Kondor et al., but improves on it in a number of ways: firstly, a channel-wise structure is introduced for the tensor product nonlinearities, which avoids the degree blowup of this operation while still allowing mixing between different harmonic degrees. Secondly, computational complexity of linear layers is reduced by factorizing it into three operators, two of which operate similar to depthwise-separable convolutions and one of which acts uniformly across channels. Thirdly, an optimized sparse degree mixing set is proposed, based on a minimum spanning tree. Finally, a more efficient sampling theorem is used that reduces the Nyquist rate by a factor of two compared to the ones used in previous works on spherical CNNs.

This paper introduces a generalized spherical convolution operation that is strictly equivariant to rotation. The authors show that the spherical convolution operations introduced in prior works can be encompassed by the proposed approach. Because spherical convolutions introduce significant computational overhead, the authors also introduce an array of methods that reduce the computational cost while maintaining the model accuracy. Experiment results on multiple benchmark datasets show that the proposed approach outperforms the alternative approaches while having less number of parameters.