Learning Transformation Groups
and their Invariants

Taco Cohen (10137408)
Supervisor: prof. Max Welling
42 EC

University of Amsterdam

A thesis submitted in conformity with the requirements for the degreeofMSc. in Artificial Intelligence2013

Acknowledgements

I would like to thank Max Welling for his supervision and for the freedom
and encouragement he gave me during the last half year. Despite his busy
schedule, he has always made time to help me when needed.
I would also like to thank Leo Dorst for teaching me many things and forthe many interesting conversations and discussions. Leo has been a great
example for me, and the best teacher I have had.
Thanks to my committee Max, Leo and Joris Mooij for agreeing to readmy thesis on such short notice.
And finally, I wish to thank my parents for their never-ending support.

Abstract

A fundamental problem in vision is that of invariance: how objects areperceived to be essentially the same despite having undergone varioustransformations. When it is known a priori to which transformations a
representation of the image must be invariant, one can try to construct
such an invariant representation by analysis. A potentially more genericsolution is to build systems that learn how to construct an invariant representation by looking at transformed examples. Such a system can belargely agnostic to the representation of the input and the kind of transformation to which invariance is required, so that it can be applied tomany kinds of data.
We approach this problem from a dynamical perspective, by modeling thetransformation group itself. This allows us to not only build an invariant
representation, but also perform pose estimation and understand motion.
We derive an elegant algorithm that learns a so-called maximal toroidal
Lie group from pairs of inputs, where each pair is related by an unkown
transformation from the group. Given one or more pairs of input images,
the model provides a full posterior over elements of the group, describing the inherent uncertainty in the transformation estimate. As it turns
out, the ambiguity in the estimate of the transformation from one image
to itself (the symmetry of the image with respect to the learned group)
provides a maximally informative invariant representation of that image.

Contents
1 Introduction
1.1 The conception of symmetry . . . . . . . . . . . . . . . . . . . . . . .
1.2

1
4

Structure preservation in human knowledge representations . . . . . .

6

1.2.1

Vision . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6

1.2.2
1.2.3

Children’s conception of the world . . . . . . . . . . . . . . . .
Analogy-making . . . . . . . . . . . . . . . . . . . . . . . . . .

7
7

1.2.4

Physics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8

1.3Group structure as an inductive bias . . . . . . . . . . . . . . . . . .

9

1.4
1.5

Why model transformations and invariants? . . . . . . . . . . . . . .
Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

11
13

2 Mathematical Preliminaries

14

2.1

Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

14

2.2

Group theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.1 Elementary definitions . . . . . . . . . . . . . . . . . . . . . .

14
15

2.2.2

Lie groups and Lie algebras . . . . . . . . . . . . . . . . . . .

17

2.2.3

Reduction of a representation . . . . . . . . . . . . . . . . . .

18

2.2.4
2.2.5

Orthogonal & toroidal groups . . . . . . . . . . . . . . . . . .
Relative & absolute invariants, weights and degeneracy . . . .

19
20

3 Representation Learning

21

3.1

Representation learning algorithms . . . . . . . . . . . . . . . . . . .

25

3.2

3.1.1 ICA and ISA . . . . . . . . . . . . . . . . . . . . . . . . . . .
Transformation learning algorithms . . . . . . . . . . . . . . . . . . .

25
26

3.2.1

Anti-symmetric Covariance Analysis . . . . . . . . . . . . . .

26

3.2.2

Gating models . . . . . . . . . . . . . . . . . . . . . . . . . . .

27

3.2.3
3.2.4

Phase-amplitude models . . . . . . . . . . . . . . . . . . . . .
Lie group methods . . . . . . . . . . . . . . . . . . . . . . . .

28
29

i

3.3

Basic principles for representation learning . . . . . . . . . . . . . . .

4 Toroidal Subgroup Analysis
4.1

30
34

The TSA model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.1.1 The likelihood function . . . . . . . . . . . . . . . . . . . . . .

34
34

4.1.2

The prior over phases . . . . . . . . . . . . . . . . . . . . . . .

36

4.1.3

The posterior distribution . . . . . . . . . . . . . . . . . . . .

36

4.1.4 Relation to discrete Fourier transform . . . . . . . . . . . . . .
Learning by expectation-maximization . . . . . . . . . . . . . . . . .

38
39

4.2.1

E-step . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

40

4.2.2

M-step . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

41

5 Low-dimensional subgroups
5.1 The stabilizer subgroup . . . . . . . . . . . . . . . . . . . . . . . . . .

43
43

4.2

5.2

Lie subalgebras and phase wrapping . . . . . . . . . . . . . . . . . . .

45

5.3

Learning one-parameter subgroups of a torus . . . . . . . . . . . . . .

47

5.3.1
5.3.2

47
48

Real-valued weights with explicit phase unwrapping . . . . . .
Integral weights . . . . . . . . . . . . . . . . . . . . . . . . . .

6 Experiments

51

6.1

Random data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

51

6.2
6.3

Learning image filters . . . . . . . . . . . . . . . . . . . . . . . . . . .
Testing rotation invariance . . . . . . . . . . . . . . . . . . . . . . . .

52
53

6.4

Rotation invariant classification . . . . . . . . . . . . . . . . . . . . .

53

6.5

Subgroup learning . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

54

7 Conclusion and Outlook

56

A Derivation of von Mises-Polar Gaussian conjugacy

57

B EM Equations

59

B.1 Expected rotation matrix under von Mises distributed angle . . . . .

59

B.2 Derivation & optimization of EM objective . . . . . . . . . . . . . . .
B.2.1 Optimization with respect to W . . . . . . . . . . . . . . . . .

60
61

Bibliography

63

ii

Chapter 1
Introduction
No man ever steps in the same river twice, for it’s not the same riverand he’s not the same man. – Heraclitus
Not one image captured by our eyes is exactly the same as another one captured
before it, and neither is the visual cortex that interprets it. However, despite thisfact and despite Heraclitus’ philosophical musings, most of us perceive the river andour selves as having some degree of permanence. In other words, we have the ability
to recognize (or confabulate) some invariant essence of a percept, even though thepercept is always changing and never exactly the same.
Which changes are considered essential and which ones are considered irrelevant
depends on context, but in order to make the distinction at all we need to reduce theholistic percept into distinct properties. Consider figure 1.1. We say that the two letters are essentially the same despite a difference in font, and likewise we recognize thesame face with different facial expressions. The properties we use in our description(letter, font, identity, expression) are not explicitly present in the raw representationof the image on our retina, yet they appear completely natural to us. Somehow, wehave decomposed the holistic percept (the entire collection of light intensity values)
into constituent properties so that we may assess the degree and kind of similarity ina meaningful way.
The previous paragraphs could have been the introduction to a thesis in philosophy, psychology, psychophysics or cognitive science, but here we endeavor to understand this phenomenon from a mathematical perspective, and to devise effective
algorithms that learn to see the world in a similar way. That is, we aim to devise learning algorithms that can untangle the different factors of variation, so thatobservations can be represented in terms of distinct properties.

1

Figure 1.1: The same, yet different.
To make the vague notions of “essentially equivalent despite superficial changes”
and “separation into individual properties” precise, we apply two marvelous insights ofmathematician Felix Klein (1849-1925) and mathematician/physicist Hermann Weyl
(1885-1955). Both ideas are deceptively simple, to the point of seeming almost trivial,
but nevertheless have far-reaching consequences and have had a massive impact on
the course of mathematics and natural science.
The principle put forward by Klein in his famous Erlangen Program [Klein, 1893]
was motivated by problems in geometry. In traditional Euclidean geometry, one considers two figures equivalent if they are related by a Euclidean transformation (a rigid
body motion). As we will see in section 1.1, the set of all Euclidean transformations
is one example of a type of mathematical object known as a group. Groups are themathematicians’ way to describe symmetries, and the Euclidean group describes thesymmetries of Euclidean space (this usage of the word “symmetry” will be clarified
later). The key idea is this: instead of taking the Euclidean group as our symmetrygroup, we can take any group and call it a symmetry group. Every choice of symmetry group will induce a different notion of equivalence between figures, and will thus
result in a different kind of geometry. According to Klein, the objects of interest in aparticular geometry are the symmetry transformations and any quantity left invariant
by these transformations. What Klein showed is that most geometries that mathematicians had cared to study (i.e. those that they had found interesting for somereason) could all be obtained by the right choice of symmetry group. Furthermore,
the symmetry groups are related to each other in various ways, and this clarified therelations among the different geometries.
Notice how Klein has given a precise definition of what is to be considered interesting: invariants are worthy of study, other quantities and structures are not. Althoughpaternalistic to the mathematician, such a principle, if representative of our intuitions
about what constitutes interesting structure, is invaluable to those who study naturaland artificial intelligence.

2

Applied to the problem of learning representations of data, Klein’s principle tells
us that the interesting quantities – those that we want to represent explicitly – areinvariant to some group of transformations. We can be a little more specific, too: wewant a maximal invariant. Taken together, the features should identify the input up
to symmetry transformations [Soatto, 2009].
A maximal invariant could be represented in many ways, because any invertible
mapping of a set of features contains the same information and absolute invarianceproperties. As a further constraint then, we will require that each feature has aninvariant meaning that is independent of the other features. But how can we formalize
such a vague notion? To approach this problem we employ a basic principle first
expounded by Hermann Weyl [Weyl, 1939] and later seen to be of use in computer
vision by Kenichi Kanatani [Kanatani, 1990], that states that:Weyl’s principle (intuitive version)
A set of transformed measurement values may rightly be called a single physicalobservable with an invariant meaning if and only if the transformations from a symmetry group do not mix this observable with other observables.
In other words, if two sets of measurements, made at time t and t + 1, are related
by a symmetry transformation, the state of one particular observable (which is afunction of the measurements) at time t + 1 should be completely determined by its
state at time t, and the symmetry transformation that was applied. Information aboutother observables should not be necessary to determine its state, because representing
distinct properties of the observed phenomenon, they transform independently. This
is not to say that different properties cannot interact, but if the transformation is just“a change in perspective” or “a renaming of points in the measurement space” (i.e.
it is a symmetry transformation), then the real phenomenon that is being observedis the same (hence no interactions), and we should be able to tell the next state fromthe current one for each observable independently. The precise definition of Weyl’s
principle stated in terms of the irreducible reduction of a group representation ispresented in section 3.3, and a rudimentary introduction to group theory required tounderstand it is given in chapter 2.
From these basic principles and a small number of additional assumptions, wederive a representation learning algorithm which we call Toroidal Subgroup Analysis.
This algorithm learns from a stream of transforming data to represent that data in
an invariant-equivariant way. That is, the inputs are re-represented in terms of the3

learned features (Weyl’s observables), from which an absolutely invariant representation can easily be computed. At the same time, due to the reduction into independent
observables, transformations in the data can be described explicitly and compactly
in a disjunct pathway. While distinctly geometrical, the algorithm is expressed in thelanguage of probabilistic graphical models, enabling us to make sound probabilistic
inferences about geometrical quantities. When trained on randomly shifting images,
the algorithm learns to perform a discrete Fourier transform of the signal; other typesof transformations result in analogous computations. We then extend the algorithm
to learn subgroups of the learned toroidal group. Using the theory of Lie algebras, wefind that these subgroups correspond to linear subspaces of a space of latent variablesof the model. The resulting representation of transformations in the learned subgroupconsists of a small number of parameters that are easily interpretable.
We will use examples from visual perception as a way to motivate and explain
the TSA algorithm, and the experiments we perform are all done on image patches.
However, we will argue that the principles we adopt are useful for structuring acquired
knowledge in general, because they facilitate generalization and aim at objectivity.
In the following sections of this chapter, we explain what it is exactly that we meanby “symmetry”, “preservation” and “group”, and show how these notions appear toplay a role in structuring human knowledge. Then, we motivate our work on learningtransformation groups and invariant representations from the perspective of machinelearning. In section 1.5, we outline the contributions of this thesis.

1.1

The conception of symmetry

The usage of the word “symmetry” in this thesis may appear somewhat unusual tothose who do not have a background in mathematics or physics. The image thatcomes to mind upon hearing this word could be one of those shown in figure 1.2.
These are indeed examples of symmetry, but the meaning we intend to convey ismore abstract. Before we give the abstract definition of symmetry, let us make someobservations about these concrete geometric symmetries.
In each of the panels in figure 1.2, notice that there are some transformations
that we can apply to the figure so that it remains unchanged. We call these transformations symmetries of the figure. Notice further that there may sometimes be“non-essential” or “superficial” aspects of the figure that are changed by the symmetry transformations. For example, the colored patches on the triangle are permutedby reflections and rotations, even if these transformations leave the triangle shape4

Figure 1.2: Examples of concrete (non-abstract) symmetry.
invariant. For both figures, the symmetry transformations are a subset of the Euclidean transformations of rotation, translation, or reflection (or some combination ofthem).
Obviously, when we have two transformations that both leave the essence of the
figure unchanged, the sequential application of the two will do the same. It is thistrivial observation that forms the basis of the mathematical notion of a group of transformations, from which surprisingly abundant structure emerges. Formally, a group
G is a set of transformations1 together with a product that composes transformations.
To be a group, the set and product must satisfy the following axioms:
1. Closure: if A ∈ G and B ∈ G, then AB ∈ G.
2. Inverse: if A ∈ G then A−1 ∈ G, where AA−1 = I and I is the identity
transformation.
Juxtaposition of two transformations is taken to mean subsequent application, whichcorresponds to matrix multiplication when the transformations are linear and are
represented by matrices in some fixed basis. Normally, the additional axioms of
associativity ((AB)C = A(BC) for all A, B, C ∈ G) and identity (I ∈ G, where
AI = IA = A for all A ∈ G) are included in the definition, but these are superfluous
when it is assumed that the elements of G are transformations.
We can read the definition of a group as a recipe for generating one: given a
number of symmetries of an object, we can generate a whole group of symmetries by
repeatedly multiplying elements and taking inverses, adding them to the group until
we find no more new transformations. So we could say that groups are born fromsymmetries.
1

The elements of a group need not be interpreted as transformations in order to study groups
abstractly, but for now we will not make the distinction between abstract groups and group representations, in order to simplify the exposition.

5

The key idea that leads to the abstract notion of symmetry is the following: insteadof considering the Euclidean group of transformations to be the symmetry group of
space, we can take any group of transformations and call them symmetries. Thismeans that we now have to identify all the objects that can be transformed into eachother by some element of the group, even if this means that objects can be distorted
in non-rigid ways. They are considered “essentially the same” with respect to thechosen group.

1.2

Structure preservation in human knowledge representations

Pure geometry is far from the only application of group theory. In mathematics
more generally, but also in physics and chemistry, symmetries and the groups that
describe them play an important role. While this fact is widely appreciated, it is lessoften recognized that more mundane propositions and non-scientific proto-theories
too, are often of the form “there exists a type of change that leaves some aspect of an
object invariant” or “some structure x is preserved when we do y”. In the next few
sections, we look at examples from vision, developmental psychology, cognitive scienceand physics. We do not claim that all these instances of preservation of structure
must necessarily be modeled using transformation groups, but we do believe thatthe general concept of “preservation of structure” is fundamentally important in our
perception and cognition.

1.2.1

VisionVisual perception is the main area of application in this thesis, and it is there thatthe relevance of the concept of abstract symmetry is most obvious. A lot of researchin computer vision is geared towards constructing or learning image descriptions thatare invariant to some kind of transformation. It has been suggested that the central computational goal of the visual cortex is to learn to discard transformations[Leibo et al., 2011], both generic ones such as translations and rotations, as well asclass-specific or local ones that only apply to a particular region of the input space,such as transformations of facial expressions. We subscribe to this view, but emphasize that besides invariance, equally interesting are the maps to which we build
invariance themselves. In analogy to Klein’s perspective on geometry, we believe thatthe most natural way to define static object classes in vision is to say that they are
the invariants of various transformation groups that act on the visual scene.
61.2.2Children’s conception of the world

The knowledge acquired by young children can often be framed in terms of transformations and invariants. These transformations act not on the perceptual arrays
directly (as in vision, for example), but on slightly more abstract mental representations. Consider the following classical experiment by developmental psychologist
Jean Piaget [Piaget, 1952]. The child is shown two identical glasses, filled with thesame amount of liquid. The child is asked whether the two glasses have the sameamount of liquid or if the amount is different. At the ages tested, children are typically able to appreciate that the amount is the same. Then, the experimenter pours
the liquid from one of the containers over to a taller, thinner glass, in plain view of
the child. The child is then asked whether there is still the same amount of liquid in
both glasses, to which a young child will answer “no, there is more water in the tall,thin glass”.
We do not know how this sequence of events is represented in the brain of the child,
but certainly the information about the height of the liquid at each point in time isrepresented somehow. Hence there is a transformation of this mental representation
as the experiment unfolds. As the child develops, it learns to recognize that many
wildly different image sequences all constitute the continuous transformation of thescene called “pouring over” and that this transformation has an invariant: the volumeof the liquid. Piaget invented many more such conservation tasks, showing that manyof the elementary concepts learned by a child (and carried into adulthood) should beviewed in terms of transformations and conserved quantities.

1.2.3

Analogy-making

An analogy is a mapping between two (more or less disparate) cognitive or perceptual
domains, that preserves the structure of the source domain. For example, one canmake the (somewhat naive) analogy that “the atom is like the solar system”: bothhave a central object (the sun ∼ atomic nucleus) and objects orbiting it (planets ∼electrons). The “analogy-function” A maps sun to nucleus and planets to electrons.
The map A radically changes the features of the object: the atomic nucleus is very
different from the sun in terms of size, for instance. But the relational structure ispreserved: orbits(planet, sun) → orbits(A(planet), A(sun))
We can again make an analogy to geometry. In geometry, an object is defined interms of the relations between the points (these relations being distances and angles

7

in Euclidean geometry), and is considered equivalent (or analogous) to a transformed
object if both objects have the same relational structure.
Once an analogical mapping has been established, it can be used for inference
about the target domain. This process is known as analogical inference or relational
generalization [Holyoak and Morrison, 2005]. As an example of analogical inference,
consider the following situation: having already acquired the concept of gravity andfound an analogical mapping between the solar system and the atomic nucleus, one
could conjecture the existence of some force keeping electrons with the nucleus. Although this mode of reasoning can easily fail, it is nevertheless indispensable for generating plausible conjectures. It allows us to generalize in radically non-local ways;
something which machines are currently incapable of.

1.2.4

Physics

In modern physics, practically all laws of nature are derived from symmetry principles.
That is, one starts by assuming that the laws governing the behavior of some physicalsystem are invariant to some group of transformations (symmetries) acting on thestates of the system. From this assumption alone, one can learn a great deal about
the laws themselves, for it is a strong restriction on their form. Furthermore, any
continuous symmetry of a system corresponds to a conserved quantity such as massenergy, momentum, etc – by Noether’s theorem. What matters to us at present isthat physics is yet another example of a knowledge domain where symmetry andinvariance are the main structuring principles.
The grand success of modern physics shows that nature is indeed well describedin terms of symmetries and conserved quantities, but the fact that this body ofknowledge has come about may also inform us about our way of thinking. The
propensity for searching for preserved, unchanging quantities can already be found in
the ideas of Thales – who, according to Aristotle was the first naturalistic philosopher.
He sought to explain phenomena, in particular the nature of matter, without recourse
to gods and deities [Curd, 2012]. In his Metaphysics, Aristotle writes of Thales’ ideas[Aristotle, 1924]:
That of which all things that are consist, the first from which they cometo be, the last into which they are resolved (the substance remaining, butchanging in its modifications), this they say is the element and this theprinciple of things, and therefore they think nothing is either generated or
destroyed, since this sort of entity is always conserved [...] Thales,8the founder of this type of philosophy, says the principle is water(emphasis mine)

The ancient Greeks did not select such a theory from among many other alternatives because it was empirically successful in describing nature (it was not), but
because this is the way we prefer to structure knowledge.

1.3

Group structure as an inductive bias

We have seen some striking examples of symmetry in human knowledge representations. At this point however, it is not quite clear why we see this pattern, and whether
it is something we should strive to replicate in artificial systems or not. The organization of human knowledge could simply be a reflection of the (often symmetric)
structure of the world we perceive, inscribed as it were on the tabula rasa that is thehuman mind. Or it could be that we have prior knowledge of certain symmetries oreven an inductive bias towards group structure in general.
Henri Poincaré wrote on the subject:
The object of geometry is the study of a particular “group”; but the general concept of group pre-exists in our minds, at least potentially. It is
imposed on us not as a form of our sense, but as a form of our understanding.
Only, from among all the possible groups, that must be chosen which willbe, so to speak, the standard to which we shall refer natural phenomena.
In his book Science and Hypothesis [Poincaré, 1904], Poincaré argued that the conceptof a group is innate, and key to reasoning itself. Only, he says, the particular groupthat we typically use to explain the movement of objects around us – the Euclideangroup – is derived from experience. In other words, it is learned ; chosen among manyother possibilities simply because it gives a compact description of some aspects of themovie projected onto our retina, and the impressions arriving at other sensory arrays.
The retina itself provides a very non-homogeneous, high-dimensional representation of
the visual world, and so it is not obvious at all that a representation of the Euclideangroup (of all groups) should act on this space2 . The fact that from this signal weconstruct our conception of 3D space is in itself quite astonishing. According to2

Strictly speaking, the Euclidean group acts on the sensory arrays augmented with latent variables
specifying depth and occlusion, for otherwise transformations are not invertible. We will get backto this topic later.

9

Poincaré, the most plausible explanation is that we have an a priori concept of agroup, the details of which are filled in by experience, using our actions as a sort ofsupervision for the learning process. The mechanism he proposed can be summed up
as follows (quoting again from [Poincaré, 1904]):
1. In the first place, we distinguish two categories of phenomena: The first involuntary, unaccompanied by muscular sensations, and attributed to external objects– they are external changes; the second, of opposite character and attributed
to the movements of our own body, are internal changes.
2. We notice that certain changes of each in these categories may be corrected bya correlative change of the other category.
3. We distinguish among external changes those that have a correlative in the othercategory – which we call displacements; and in the same way we distinguishamong the internal changes those which have a correlative in the first category.
Let us work out an example. We will use xt to indicate the aggregate of perceptual
measurements (pixel values, say) at some time instant t. We could model this (as is
typical in machine learning) as a vector xt ∈ RN . Now we perform an experiment: we
perform a small eye movement denoted g, and observe the resulting percept xt+1 =gxt . We can do this over and over, yielding a data set {(xt , xt+1 ; g t )}t=1...T , fromwhich we can try to learn the functions gxt for all possible g, that predict the changes
in the perceptual space (which we take to be visual space in this example, but we may
substitute auditory space, haptic space, etc. when considering other actions). Noticethat if we have one saccade g and another saccade h, the effect of two successive
saccades hg can also be achieved as one single saccade. Furthermore, for any eye
movement g we can do the reverse movement g −1 . Eye movements are also associative
and we can perform the identity movement, so if we ignore actuator limits3 (we cannot
see into the back of our head), the saccades form a group. As we will see in later
chapters, the group formed by these transformations, together with the projection
function that generates the retinal image from scene parameters, determine a grouprepresentation in the space of pixel values. So we see that in order to predict theeffect of our actions on the world as projected onto our senses – at least in this simple
example – we must learn (or at least somehow develop) a group representation.
3

Such practicalities have not hindered physicists in the successful application of group theory to
the study of crystal structures either, even though crystals are not actually infinite in extent.

10

There is some evidence that this function is indeed learned and not hard-wired.
In a clever experiment, [Cox et al., 2005] modified the identity of visual objects during the period of transient blindness that accompanies eye movements, and showedthat after a training period, this manipulation caused failures in position invarianceof human subjects under normal viewing conditions. It has also been shown that
this manipulation alters the neuronal representation of objects in primate cortex[Li and DiCarlo, 2008], as would be expected.
Both Poincaré and Piaget argued for an inborn preference for group structure
in human cognition. The fact that abundant group structure in our sensorium iscertainly an evolutionary constant supports this idea, and we have presented a very
rough sketch of how such a bias could guide learning. Of course, it is not possible toestablish the claim of innateness of the group concept with any degree of certainty.
One step we can take is to give an existence proof: to show that such a learningalgorithm can indeed work, and do useful things. Our work on Toroidal Subgroup
Analysis can be seen as a small step in this direction.

1.4

Why model transformations and invariants?

The ideas relating to human cognition set forth in the previous sections motivate
our work on learning transformations and invariant representations, but we realizethat this motivation may not be all too compelling to the well-grounded engineer
who is most interested in applications. In this section, we will describe a numberof advantages of the group-theoretical approach as seen from a machine learningperspective, and show how an inductive bias towards group structure is useful.
One of the central problems in machine learning is the classification problem. The
target function f (x) that maps each point x ∈ RN to a label c1 , . . . , cK partitions
the input space RN into non-overlapping sets corresponding to the labels4 . Thispartition defines an equivalence relation ∼, as follows: x ∼ y ⇐⇒ f (x) = f (y);data points x and y are equivalent in some sense (f -equivalent) if they have the samelabel. Every equivalence relation, in turn, corresponds to the action of a group G:x ∼ y ⇐⇒ ∃g ∈ G : y = gx. It is easy to check that this is well defined and satisfies
the axioms reflexivity (because I ∈ G), symmetry (because g ∈ G =⇒ g −1 ∈ G)
and transitivity (because g, h ∈ G =⇒ gh ∈ G) making ∼ an equivalence relation.
So we see that a target function has a symmetry group that consists of all invertible
4
Label noise is easily incorporated in this framework by considering f to be the Bayes-optimal
decision boundary derived from some probability distribution11

maps that permute points within the partitions, but not between them. That is, G
is defined implicitly as the largest group that satisfies ∀g ∈ G, ∀x : f (gx) = f (x).
Is this dual perspective on label functions useful? If we can learn or in some otherway obtain G, we could build an exemplar-based classifier that classifies a new data
point x as class ci iff ∃g ∈ G : gx = x(i) for an exemplar x(i) with f (x(i) ) = ci .
Alternatively, we could construct a representation Φ : RN → RM that is invariant
to the group action: ∀g ∈ G : Φ(gx) = Φ(x). To avoid trivial representations, wewill require that data variation that cannot be explained by the symmetry group isretained: ¬(∃g ∈ G : gx = y) =⇒ Φ(x) 6= Φ(y). In both examples, we can learn ffrom one example per class, once the symmetry group is known; a great reduction insample complexity.
It is quite clear though, that this is just a different view of the same learning
problem, and by itself does not address the fundamental difficulty of learning. The
group G, defined implicitly as the set of all transformations that leave f invariant,can be uncountably infinite and may not admit a convenient parameterization. If
the decision boundary is highly irregular, then its symmetry group will typically becomplicated as well, and so this perspective is no magic bullet. However, as we will
now discuss, learning a symmetry group instead of a label function conveys a numberof advantages, and provides an integrated approach to learning about both static anddynamic aspects of the input.
First, symmetries may be shared between learning tasks. When this is the case,the group structure need only be learned once, and all tasks that share this symmetry
can make use of it. Each of these learning tasks is then simplified, because all thevariability in the data that is caused by the action of the group can be taken out,yielding a problem of reduced sample complexity.
There is increasing interest in learning classes without supervision. This is often
done with clustering algorithms, but it is not quite clear by what principle observations should be grouped together. Due to the closure property inherent in thedefinition, groups form clearly defined “units”, or “wholes”. It is quite clear whena set of transformations form a single unit (a group), but it is not a priori clear at
all when a region of the input space corresponds to a single class. So it may wellbe easier to learn transformation groups without supervision, and have them induce“proto-classes” (the invariants), than it is to use other heuristics to designate certain
regions of the input space as belonging a single class.
A model of the dynamics is not only useful in order to induce static object classes;
it is also of interest in its own right. It is quite clear that (artificial) intelligence12

involves much more than just classification and regression: it involves planning, prediction and action, among many other aspects. All of these involve a dynamic modelof the world. A truly scalable approach to AI must learn such a model from data.
Finally, as we will see in section 3.3, the group theoretical perspective gives usa handle on the problem of “untangling” that we discussed in the beginning of thischapter. Bengio describes its significance, [Bengio, 2013]
Of all the challenges discussed in this paper, this [the disentangling
problem] is probably the most ambitious, and success in solving it will
most likely have far-reaching impact. In addition to the obvious observation that disentangling the underlying factors is almost like pre-solving
any possible task relevant to the observed data, having disentangled representations would also solve other issues, such as the issue of mixing
between modes.

1.5

Contributions

The main conceptual contributions of this work are a group-theoretical analysis oflearning and the statement of a number of basic principles for representation learning(section 3.3). We give a formal definition of an “observable” for application in representation learning. This doctrine is not new in physics [Weyl, 1939], but has neverbeen applied to machine learning. This principle provides a firm theoretical basis forthe untangling problem whose significance we have discussed in section 1.4.
Our main technical contribution, presented in chapter 4, is the Toroidal Subgroup
Analysis algorithm. It applies Weyl’s principle to the toroidal subgroups of SO(n), thelargest compact, connected, Abelian Lie subgroups of the special orthogonal group.
TSA is a probabilistic model that learns to split the input signal in an invariant anda covariant part. As a special case (the 2D translation group), the algorithm learns
to perform a discrete Fourier transform, and hence provides a novel probabilistic
interpretation of the discrete Fourier transform.
TSA relies on a newly derived conjugacy relation between a circularly parameterized Gaussian likelihood and a von Mises prior. Such relations greatly simplify
Bayesian analysis, so this result may be of interest to researchers in circular statistics.
Finally, we use the theory of Lie algebras to derive an extension of the basic
TSA model that aims to learn a low-dimensional subalgebra of a toroidal Lie algebra(chapter 5). Since all compact Abelian Lie groups are the subgroup of some maximal
torus, this work (though incomplete) would complete the compact Abelian case.
13

Chapter 2
Mathematical Preliminaries
Throughout this thesis, we assume that the reader is familiar with linear algebra,
calculus and probability theory. In this chapter, we give a quick introduction togroup theory, which we expect to be less well-known in the field of machine learning.
The notation, definitions and results given in this chapter will form the basis for later
chapters.

2.1

Notation

Scalars x, y, z are written in lower case font, vectors x, y, z in bold lower case font.
Matrices and blades (to be defined shortly) are written as bold upper case romanletters A, B, C. The j-th element of a vector x is denoted xj , the element in row i,column j of matrix A is written Aij . The row i of a matrix is written Wi while thecolumn vector j is written as W·j . We will want to select the elements corresponding
to 2D, basis-aligned subspaces from vectors and matrices, which we denote as follows:
we write xj = (x2j , x2j+1 )T and Wj = (W·2j , W·2j+1 ). For an N × N matrix W,
Wj is an N × 2 matrix. Elements of abstract groups G are denoted by lowercase
letters g, h, . . ., as opposed to concrete matrix or versor representations A, B of theseelements.

2.2

Group theory

In this section we introduce some basic group theoretical concepts that we will needlater on. This introduction is meant to be easy to understand for those not familiar
with advanced mathematics and as such we do not give proofs (for they typicallymake it harder to see the forest for the trees). For the same reasons, we do not try to
make our definitions as general as possible and we do not go far beyond what we will
14

need directly in the derivation and motivation of the Toroidal Subgroup Analysis algorithm, so this is far from a complete introduction to group theory. Readers interested
in a more comprehensive treatment are referred to [Weyl, 1939] (a classic but difficult treatment of Lie groups), [J.P. Elliot, 1985] (applications in physics, somewhat
more accessible), [Kanatani, 1990] (applications in computer vision), [Kondor, 2008](applications in machine learning).

2.2.1

Elementary definitions

A group G is a set (finite or infinite), endowed with a product on its elements thatsatisfies the following axioms:
1. a(bc) = (ab)c for all a, b, c ∈ G (Associativity)
2. ∃e ∈ G such that ea = ae = a for all a ∈ G (Identity)
3. If a ∈ G and b ∈ G, then ab ∈ G (Closure)
4. If a ∈ G then a−1 ∈ G, where a−1 a = aa−1 = e (Inverse)
Note that commutativity, ab = ba, is not required. Groups that do satisfy thisadditional constraint are called commutative or Abelian groups. Commutative groups
are much easier to analyze, so they will be the main object of study for us.
As an example of a finite group, consider the set of all permutations on n letters. If
we take multiplication to mean the consecutive application of two permutations, thisset forms a group. The set of all rotations in a fixed rotation plane is an example ofan infinite commutative group. Other (generally non-commutative) examples include
GLn , the group of all invertible linear transformations in n dimensions, the group On
and SOn of orthogonal and special orthogonal transformations (these will be treatedin section 2.2.4).
It is often convenient to specify a finite group in terms of its generators. A
subset S ⊆ G is said to generate G if all elements in G can be obtained by repeatedmultiplication and inversion of the generators. It is often revealing to display afinite group as a graph where each group element is shown as a node and each edge
represents the multiplication by a generator. Such a graph, displayed in figure 2.1, is
known as a Cayley graph.
A subgroup H of G is a subset of G that is also a group. That is, the subset H
must satisfy the group axioms defined above for it to be a subgroup. In figure 2.1,
we can spot a subgroup by considering those sets of nodes that include the identity
15

Figure 2.1: Cayley graphs of the Quaternion group Q4 (left), a group of order 21(middle) and the alternating group A5 describing the symmetries of the icosahedron(right). Nodes represent the elements of the group. Arrows represent the left multiplication by a generator, such as i (blue arrow) and j (red arrow) in Q4 . Figures
courtesy of Nathan Carter, reproduced with permission [Nathan Carter, 2009].
and are connected by only one color of edges. For example, the sets {1, i, −1, −i} and{1, j, −1, −j} are subgroups of the group of quaternions shown in figure 2.1.
A group homomorphism is a mapping φ : G1 → G2 that preserves the groupstructure. More precisely, for any a, b, c ∈ G1 satisfying ab = c it must be the casethat φ(a)φ(b) = φ(c). Notice that the product ab is performed according to the multiplication defined for G1 , while the product φ(a)φ(b) is performed according to themultiplication in G2 . The determinant function is an example of a homomorphismfrom the general linear group (which can be represented by the set of invertible n × nmatrices) to the group of real numbers (excluding zero) under multiplication: fora, b, c ∈ GLn , we have ab = c =⇒ det(a) det(b) = det(c). Geometrically, this saysthat the volume expansions effected by two subsequently applied transformations
multiply. Notice how structure has been preserved in the sense that no equation
involving only products of group elements can be made to fail by applying the homomorphism, while nevertheless we have lost something because φ = det is surjective:
many matrices have the same determinant.
When a homomorphism is invertible (bijective), it is called an ismorphism. In
this case, the groups have exactly the same multiplication structure. This is all
the structure that is intrinsic to the group, so we may consider isomorphic groups
to be equivalent. Still, it may take some work to see that two groups are indeed
isomorphic. For example, the symmetry groups of the cube and the octahedron

16

(obtained by associating each vertex of the cube with a face and vice versa) have thesame symmetry group1 .
A group action is an associative map φ : G × X → X with the property thatφ(e, x) = x for all x ∈ X. Here X is a set, which we will typically take to be
the vector space Rn . A group action is naturally realized in the setting of grouprepresentation theory, where we represent abstract groups with elements a, b, c, . . . bya homomorphic group of matrices A, B, C, . . .. That is, a representation of a group
G is a group homomorphism ρ : G → GL(n) from G to the group of invertible n × nmatrices. In this case, the group action can be defined simply as φ(g, x) = ρ(g)x.
The orbit of an element x ∈ X under the action of G is the set {φ(a, x)|a ∈ G}.
This is an important concept for our applications in representation learning, where
X = Rn is some data space and x ∈ Rn are input vectors. In this setting, theorbits are all the points to which an input x ∈ Rn can be transformed by symmetrytransformations in some group G that we wish to disregard. Our goal then, is tocollapse the orbit into a single point in the representation space (i.e. to construct
an invariant representation) while parameterizing the position on the orbit with acomplementary variable. This is formalized by the notion of an orbit space, denoted
X/G. The points in the orbit space correspond to orbits of G in X.
To make this concrete, consider a 2D space R2 and a representation of the rotationgroup SO(2) that rotates about the origin. A matrix representation of SO(2) is given
by 2 × 2 rotation matrices of the formcos(ϕ) − sin(ϕ)R(ϕ) =sin(ϕ) cos(ϕ)(2.1)

Its action on a vector x ∈ R2 is just R(ϕ)x. The orbits are concentric circles aroundthe origin, and the orbit space R2 /SO(2) is given by the radii of the orbits, i.e. R.

2.2.2

Lie groups and Lie algebras

A Lie group (pronounced Lee-group) is a group that is also a differentiable manifold.
Intuitively, the elements of a Lie-group lie on a smooth surface on which we can takederivatives. To keep the exposition concrete, we will only consider matrix representations of groups from now on. Most Lie groups can be represented by matrices.
The tangent space at the identity of any Lie group gives rise to a structure known
as a Lie algebra. The Lie algebra g of a Lie group G is the tangent space to G at the
1

Interestingly, [Hinton, 1979] has shown that people will confuse these two very different objects
with identical symmetries in a mental imagery task.

17

origin, together with a product on tangent vectors [·, ·] : g × g → g satisfying:[X, aY + bZ] = a[X, Y] + b[X, Z][aX + bY, Z] = a[X, Z] + b[Y, Z][X, Y] = −[Y, X](Bilinearity)(Anti-commutativity)[X, [Y, Z]] + [Z, [X, Y]] + [Y, [Z, X]] = 0(Jacobi Identity)(2.2)
For all scalars a, b and elements X, Y, Z ∈ g. Since we work with matrices, we can be
more concrete: the elements of g are tangent vectors (matrices) to the group manifold
at I, and the Lie bracket can be defined as the commutator: [X, Y] = XY − YX.
The most important fact about Lie groups is that their structure is almost completely determined by the Lie algebra. In other words, the global structure is determined by the infinitesimal structure. Every element of the group can be expressed as
the exponential of an element of the Lie algebra, and the product of two elements ofthe group can be expressed as a sum of elements in the Lie algebra:
exp(X) exp(Y) ≈ exp(X + Y + [X, Y] + . . .)(2.3)

The continuation of this series (which is quite involved) is known as the BakerCampbell-Hausdorff formula. For commutative groups, the commutator (Lie bracket)
is always zero, and hence the matrix exponential behaves exactly like the scalar exponential: exp(X) exp(Y) = exp(X + Y).
A Lie subalgebra h of g is a subspace of g that is closed under the Lie bracket.
Again, things are simpler in the commutative case, for then the commutator is always
zero and hence a subalgebra is just a subspace.

2.2.3

Reduction of a representation

Two representations ρ and π are called equivalent if there exists a single matrix W
such that Wρ(g)W−1 = π(g) for all g ∈ G. In other words, equivalent representations
are related by a change of basis.
Let ρ be a representation of G in Rn . We say that subspace V of Rn is invariantif ρ(g)v ∈ V , for any g ∈ G and v ∈ V . If ρ has a non-trivial invariant subspace, it iscalled reducible (the trivial ones being {0} and Rn ). Otherwise it is called irreducible.
If the orthogonal complement to an invariant subspace V of ρ, V ∗ , is also an invariant
subspace of ρ, then there exists a matrix W that decomposes each element of therepresentation into the same block structure:ρV (g)
0−1Wρ(g)W =0ρV ∗ (g)18(2.4)

where ρV is a representation of G in V and ρV ∗ is a representation of G in V ∗ . This
process is called the reduction of a representation. If it is possible to repeat thisprocess of reduction recursively until we get to the irreducible representations, therepresentation is called fully reducible. All compact groups are fully reducible.

2.2.4

Orthogonal & toroidal groups

The toroidal subgroup analysis algorithm derived in chapter 4 makes use of orthogonal transformations. These are the linear transformations that do not expand or
contract the space on which they act. That is, they leave lengths (and angles) invariant, and this defines them fully. Formulaically, for orthogonal Q, kQxk = kxk
and (Qx)T (Qy) = xT y. A useful computational property of the orthogonal transformations is that they can be inverted for free, because Q−1 = QT . This implies that
QT Q = I, which may also be taken as the defining property.
All orthogonal transformations can be generated by repeated reflection. An evennumber of reflections makes a rotation, an odd number of reflections makes an antirotation or a plain reflection. Rotations have a determinant of +1, while reflections
and anti-rotations have a determinant of −1. The group of all orthogonal transformations in n-dimensions is denoted O(n), while the group of rotations (determinant+1), is denoted SO(n), and called the special orthogonal group.
It can be shown that any orthogonal matrix Q in an even-dimensional space canbe brought into the following canonical form, by an orthogonal change of basis W:


R1


...
WQWT = (2.5)

RM
where each Rj is a 2 × 2 orthogonal matrix. In odd-dimensional spaces, there will bean additional 1 × 1 “block” with value ±1. We can understand this decomposition
easily in three dimensions. For an odd number of reflections, one direction will beflipped, giving a −1 on the singular 1 × 1 block. For an even number of reflections (arotation), the axis is left invariant and so there is a +1 in the 1×1 block. Furthermore,
the orthogonal complement to this direction is a 2D subspace. For a pure rotation(two reflections) or an anti-rotation (three reflections), this 2D space is rotated, while
for a reflection it is rotated by 0 degrees, i.e. R1 = I. Apparently, the concept thatgeneralizes to higher dimensions is not the rotation axis, but the rotation plane. The
axis is just an artifact that appears in odd-dimensional spaces.

19A toroidal subgroup of a compact Lie group G is a compact, connected, commutative subgroup of G. It is easy to see that for the special orthogonal group SO(2M ),toroidal groups take the form of eq. 2.5, where the 2 × 2 blocks are rotation matrices(some of which may be equal to the identity block). The matrix W performs a complete reduction of the toroidal group, by construction. A maximal toroidal subgroupof G is one that is not contained in any other toroidal subgroup of G. We will denotea toroidal subgroup of SO(2M ) in a 2M -dimensional space by TM (W), where W is
the reducing basis as in eq. 2.5. Because a toroidal group is commutative, its Lie
algebra is trivial, i.e. [X, Y] = 0 for all X, Y ∈ t.

2.2.5

Relative & absolute invariants, weights and degeneracy

Let ρ : SO(2) → GL(2M ) be a representation of the rotation group SO(2), parameterized by an angle θ. It can be shown [Kanatani, 1990] that such a representation
can always be reduced into block-diagonal form with 2 × 2 blocks by an orthogonal
matrix W (or equivalently into 1D-complex irreducible representations). Let R(θ)
denote one of these 2D representations. Clearly, we have the following facts:
R(θ)R(θ0 ) = R(θ + θ0 )(2.6)

R(0) = I(2.7)

R(θ + 2π) = R(θ).

(2.8)

R(θ) = exp(nθB),(2.9)

From these, it follows that

where n is an integer and



0 −1
B=
.
1 0

(2.10)The integer n is called the weight of the representation. A vector in the 2D
subspace on which a particular representation R acts, is called an absolute invariant
if n = 0. Otherwise it is called a relative invariant. Notice that while an absoluteinvariant will have the same value after the group acts on it, a relative invariant willchange, but will nevertheless be invariant as a vector (only the basis is changed).
If a representation has multiple irreducible representations of the same weight,
the reduction is not unique. That is, given k representations of the same weight,
any basis for the 2k-dimensional subspace spanned by these representations defines
a valid irreducible reduction. Furthermore, the sum of two relative invariants of thesame weight is another relative invariant of the same weight.
20

Chapter 3
Representation Learning
It is widely recognized that the representation of data is an important factor in
the success of machine learning algorithms [Bengio et al., 2013]. Typically, there issome kind of “nuisance variability” in the raw data whose removal can reduce thecomplexity of the learning task. But even when two representations have the same
information content (i.e. they are in bijection), one can be a more effective substrate
for learning than another (see figure 3.1). Strictly speaking (and as can be seen in thisfigure), the quality of a representation can only be defined with respect to a particular
task and a particular choice of learning algorithm. For example, a representation thatmakes two classes linearly separable is well suited to a linear classifier such as SVM,
while a decision tree classifier may fare poorly on it, because the splits required are notaligned with the coordinate axes (features). Conversely, a representation where theclasses are easily separated by a few splits on each coordinate may be a good choice
for a decision tree classifier but a very bad choice for a linear classifier. Furthermore,
multiple learning tasks can be defined on the same input domain, and different tasks
may suggest different representations of that domain. Nevertheless, it is often possible
and desirable to find one representation that works well for a range of classifiers andlearning tasks in the same input domain.
The observation that representations are important has spurred a great deal ofresearch on various descriptors and feature extraction techniques for visual and audio
data (e.g. [Bruna and Mallat, 2013], [Lowe, 2004], [Davis and Mermelstein, 1980]).
These work well for the tasks and input domains they were designed for, but handcrafting such features is a formidable task, making it infeasible to repeat this for every
input domain. This is especially so if the input domain itself is a more or less abstract(possibly learned) representation of some lower-level perceptual domain, because insuch a situation it may not be easy to understand the input domain analytically. A
black-box that finds useful representations on any input domain would be preferable.
21

Figure 3.1: Two classifiers trained on three different data representations. Left:a somewhat difficult to separate data set consisting of two classes (blue squares,
green circles). Middle: re-representation after performing a Möbius transformation(translation conjugated by circle inversions) to the data. Right: re-representation
after performing an additional rotation. Top: decision boundary (solid line) learned
by a linear SVM. Bottom: decision boundary learned by a decision tree. Clearly, therightmost representation is superior for both classifiers. Other representations (left,
middle) may work for one classifier but not so much for another.
Recently, there have been many attempts at learning representations with desirable qualities, using only unlabeled data. The objectives of this type of work often
go beyond merely learning representations that are useful for classification. According to the website of the International Conference on Learning Representations 2013
[Bengio and Lecun, 2013],
The rapidly developing field of representation learning is concerned withquestions surrounding how we can best learn meaningful and useful representations of data.
Other common objectives of this type of research include: learning representations
that are useful across a range of tasks (transfer learning), learning invariant representations, learning hierarchical representations that are progressively more abstract
and learning representations that explicitly untangle the factors of variation present
in the data. All of these are commendable goals, but at this point seem ill-posed,22

or even “not-posed”, because we have no formal definition for such intuitive ideas as“meaningfulness”, “abstractness” and “untangling”.
Despite the paucity of formal concepts, the field has been very successful, achieving
the stated goal of “usefulness” in a spectacular way. Many benchmarks for perceptualclassification tasks are now dominated by approaches that learn several layers ofinternal representations of the data. On what is arguably the most challenging image
recognition data set, ImageNet, [Krizhevsky et al., 2012] trained a deep, 60-million
parameter convolutional neural network, winning the competition by a large margin(however, this approach did not involve unsupervised learning). Using a multi-layer
network trained by Independent Subspace Analysis (section 3.1.1), [Le et al., 2011]
made significant improvements on the KTH, Hollywood2, UCF and YouTube action
recognition data sets. In speech recognition, large gains were reported on the TIMIT
data set using a deep belief network [Mohamed et al., 2012]. See [Hinton et al., 2012]
for an overview of deep neural networks in speech recognition. The architectureof [Le et al., 2012], based on sparse autoencoders was able to learn fairly high-levelobject classes such as faces and cats using only unsupervised image data.
Let us look at the goals of representation learning in a little more detail. We(and [Bengio and Lecun, 2013]) want a meaningful representation, but what makes
a representation meaningful? As a first approximation, we can consider the degree
to which the representation makes explicit such concepts as feel natural to a humanobserver. If the model produces units that are uniquely sensitive to a region of pixelspace that a human would label cat, or face (e.g. [Le et al., 2012]), then it haslearned a meaningful representations because these appear to be completely naturalobject classes. We can try out various models and algorithms, and keep the ones thatlearn representations that “make sense” to us. To actually guide learning towards
such representations however, we need a formal definition of what it is we seek. The
regions of pixel-space that correspond to natural classes are not arbitrary, so we shouldseek to characterize them in a way that does not refer to human judgment or intuition
about what are the natural classes. Clearly, human judgment itself requires such an
objective principle, although supervision (i.e. deferring to other humans’ judgment)
can play some role.
The related notion of “untangling” is used in two related but distinct meanings.
The first, as used in [DiCarlo and Cox, 2007] and [DiCarlo et al., 2012], involves flattening of object manifolds. The pixel space manifold generated by continuously changing scene parameters such as object position, lighting position, lighting intensity, etc.,
is highly convoluted, and appears tangled when projected onto three dimensions (see23figure 3 in [DiCarlo and Cox, 2007]). The authors of this paper suggest that thecentral problem in object recognition is the flattening of these manifolds.
Bengio uses the word in a slightly different way, to mean the isolation of distinct factors of variation at each point on the object manifold [Bengio et al., 2013].
For example, some images can be factored (disentangled) into subject identity, pose,
position, camera position, etc. It was noted that current representation learning algorithms seem to do some amount of disentangling, but (quoting from [Bengio, 2013]):
Although these observations [of successful disentangling] are encouraging,
we do not yet have a clear understanding as to why some representation
algorithms tend to move towards more disentangled representations, andthere are other experimental observations suggesting that this is far from
sufficient.
These two uses can be clarified by considering the notion of a manifold in a little
more detail. A (presentation of a) differentiable manifold is a set of differentiable,
invertible maps called charts φα that map some open subset Uα of a flat parameter
space Rn to a curved subset of Rm . In the machine learning setting, Rm is where we
make observations, Rn is a space of latent or hidden variables. Flattening a manifold,in the sense of [DiCarlo and Cox, 2007], means to learn the map(s) φ−1
α , so that we
can move around the data manifold by controlling latent variables in a linear space
Rn . Untangling in the sense of [Bengio et al., 2013] is more specific, in that it should
also prefer certain bases for this linear space to other bases. That is, we prefer arepresentation where we can assign some familiar names to the coordinate axes; weprefer a basis (width, height, x-pos, y-pos)T over a linear transformation (“linear
tangling”) of these, even though such a transformed space is equally flat and equally
informative about the percept. This is essentially the same issue as that of learning
meaningful representations that we discussed before, only applied to a factorial representation (object properties) instead of an object class based representation (cat,face, ...). We can restate the fundamental question in the language of manifolds:
what principle can we use to express a preference for one basis over another?
In his review, [Bengio, 2013] lists several priors (or rather, inductive biases) that
may help solve the problem of disentangling. These include: smoothness, multiple
explanatory factors, hierarchical organization, natural clustering, sparsity and several others (see section 6.2 of this paper for details). In our view, only the priorof sparsity and the related notion of statistical independence gets at the core of the24

problem as described above – the others seem to be generally useful ideas in statistical modeling, but not aimed specifically at disentangling. Two variables x, y arestatistically independent when p(x, y) = p(x)p(y) so that p(x|y) = p(x); variable yprovides no information about x and vice versa. Thus these variables may be said toconstitute different factors of variation (literally, in the sense that the joint density
factorizes), and this is the basis for the Independent Component and Independent
Subspace Analysis algorithms, which we will discuss in more detail in section 3.1.1.
However, we do not believe that independence is ultimately the right way toformalize disentangling. The reason is that many distinct factors of variation (thenatural properties we use to describe objects), are in fact correlated. We want toseparate factors that can in principle be varied independently, even if in a particular
data set (the totality of experience of an intelligent agent, for instance) these factors
are correlated and hence not independent. We think statistical independence is areasonably good proxy for a better notion of separation, which we will present insection 3.3.
In the rest of this chapter, we review a number of representation learning algorithms, highlighting similarities and differences with the Toroidal Subgroup Analysis
algorithm that will be presented in chapter 4. We distinguish three representation
learning paradigms: supervised, unsupervised and transformation learning. Classicalneural networks can be viewed as learning an internal representation in a supervisedmanner, and will not be reviewed here. In section 3.1 we will review more recent
work on unsupervised representation learning, and in section 3.2, we will look attransformation learning algorithms. Finally, in section 3.3 we present a number ofbasic principles for representation learning that are informed by group theory.

3.1

Representation learning algorithms

3.1.1

ICA and ISA

Independent Component Analysis (ICA) refers to a class of techniques for separatinga signal into independent components. The main equation of the ICA model is
x = As,(3.1)

where x is the input signal, A is the mixing matrix that is to be learned, and s isa vector of source coefficients. The goal of learning is to obtain coefficients sj thatare as independent as possible. Typical objectives include minimization of mutual
information and maximization of kurtosis [Hyvärinen and Oja, 2000].
25

If the signal conforms reasonably well with the model assumptions, such as when
the signal consists of overlapping speech signals, ICA is very effective in separating
the sources. However, it is not always possible to find independent components by
taking a linear combination of the input variables, so that the found sources show
some residual dependencies. This has led to the development of Independent Subspace Analysis (ISA) algorithms [Cardoso, 1998], where dependencies within groupsof source variables (subspaces) are allowed but dependencies between groups are discouraged. When trained on natural images using 2D subspaces, this algorithm finds
subspaces that are locally invariant under image shifts [Hyvarinen and Hoyer, 2000].

3.2

Transformation learning algorithms

3.2.1

Anti-symmetric Covariance Analysis

In a very interesting but largely overlooked paper, [Bethge et al., 2007] provide ananalysis of unsupervised learning of invariant representations, and derive an antisymmetric version of canonical correlation analysis for this purpose. Their algorithm
is based on the observation that for wide-sense stationary processes, the optimal linear
predictor is given by the time-lagged covariance matrix. Furthermore, they show thatfor whitened data, this matrix is orthogonal.
The analysis is simplest in the noiseless case. Consider the time-series xt , t =1 . . . T generated by xt = Rxt−1 . In their experiments, the data x are vectorized
image patches and the matrix R is a representation of an image rotation operator in
the space of pixels. Let X and Y be matrices with column t equal to xt−1 and xt ,respectively. Assuming whiteness of the instantaneous covariance, the time-lagged
covariance is equal to the optimal linear predictor;
1
YXT
T
1
= Y(RT Y)T
T


1
T
=
YY R
T
=R

Ct,t−1 =(3.2)

This matrix, estimated either as the time-lagged covariance or simply by leastsquares, is then decomposed to yield invariant filter pairs. As discussed in section
2.2.4, an orthogonal matrix can always be reduced into block-diagonal form with 2×2
rotation blocks by an orthogonal change of basis. A pair of basis vectors corresponding26

to one of these 2 × 2 blocks spans an invariant subspace of the transformation R, sothat the norm of the projection of a signal onto this subspace provides an invariant
representation of the signal. [Bethge et al., 2007] show that this basis matrix W canbe as the eigenvectors of the squared anti-symmetric part of R.
The ACA algorithm assumes that there is a single transformation R relating
subsequent observations, possibly with added noise. In the case of image rotations,this means that R represents a rotation by only one particular angle. Even thoughimage rotations by any angle can be represented by a single basis matrix W withvarying block-diagonal matrices, the ACA algorithm cannot learn W in this scenario.
This limitation is a significant barrier to application of this learning algorithm in realworld settings, where the stimulus is not tightly controlled. In the TSA algorithm,
we weaken this assumption to the assumption that the subsequent observations are
related by an unknown and variable element of a compact commutative group.

3.2.2

Gating models

A more flexible class of models are the so-called gating models, such as the gatedBoltzmann machine [Memisevic and Hinton, 2010], gated autoencoder [Memisevic, 2011],
synchrony-kmeans [Konda et al., 2013] and the closely related spatiotemporal energy
models [Adelson and Bergen, 1985]. For concreteness and brevity, we focus on thegated autoencoder. This model encodes the transformation between two input vectors by computing products of filter responses. One such factor signals the presenceof one feature (filter) in the first image and another in the second image. If the filter
pairs learn to span the shared invariant 2D subspaces of a group of transformations,then using appropriate (possibly learned) pooling matrices, the factor responses canbe turned into transformation-invariant representations of image content and contentinvariant representations of the transformation.
While effective and easy to train, gating models are parameterized inefficiently.
A single filter pair can detect one particular angle in one particular 2D subspace.
Therefore, to effectively model transformations from a commutative group, the modelhas to tile each invariant subspace with many filter pairs, effectively discretizing thetransformations.
Gating models were initially motivated by the wish to incorporate multiplicative
or “gating” interactions in neural networks; the relation with the decomposition oforthogonal transformations into rotations in invariant 2D subspaces was found later[Memisevic, 2012]. The gating perspective greatly increases the expressive power ofneural network models using a neurally realistic operation (taking the product of two
27

scalars), but we feel that this is too low a level of analysis. Useful mathematicalabstractions such as Lie groups exist, and we believe using them can guide the development of new transformation learning models as well as help us interpret existingmodels and the representations they produce.

3.2.3

Phase-amplitude models

Motivated by findings in natural image statistics, [Cadieu and Olshausen, 2012] proposed a model for learning representations of form and motion. The first finding isthat responses of pairs of filters learned by sparse coding of natural images tend toshow circularly symmetric kurtotic distributions. The fact that these dependencies
occur even though the prior in such models is factorial, suggests that a better model
should incorporate dependencies among pairs of filter responses. In particular, thecircular symmetry of these distributions suggests that a polar-coordinate representation is the right way to describe the signal in each 2D subspace spanned by filter
pairs. The second finding is that the log-amplitudes tend to show linear dependencies.
That is, if a signal is projected onto two learned 2D subspaces, the logarithms of the
projection amplitudes show linear dependencies.
The model projects the image patch onto a set of complex filters, each of which isdecomposed into an amplitude and phase. This is similar to the TSA model presented
in chapter 4, although the basis of [Cadieu and Olshausen, 2012] is not orthogonal
as in TSA, and the phase variables represent the absolute phase, i.e. phase relative
to an arbitrary real axis, instead of using phase offsets between subsequent patches.
This first layer is trained by optimizing an objective function that encourages sparseand slowly changing amplitudes, and phases that allow the signal to be reconstructed
properly.
Once trained, the log amplitudes and phase differences of the data are computed,
and used as input data for a second layer. The log-amplitudes are modeled using asparse and slow objective, similar to the first layer. The phase differences are modeledby von Mises distributions whose means are linear functions of a number of sparseand slow components.
While the model is specified in probabilistic language, only point estimates for thelatent variables are obtained, in contrast to full posteriors we obtain in our model.
Furthermore, the values for the latent variables are obtained by a gradient descentprocedure as opposed to a one-shot computation, precluding real-time applications.

28

3.2.4

Lie group methods

Several researchers have looked at the problem of learning Lie groups, the earliestof which is [Rao and Ruderman, 1999], later improved in [Miao and Rao, 2007]. The
basic idea in each case is to learn a set of generator matrices Gj , that are linearly
combined using a set of latent variables sj and exponentiated to obtain the transformations relating observed vectors x and y:!
y = exp

X

sj Gj x.

(3.3)

j

The model of [Miao and Rao, 2007] is a probabilistic model trained using EM. In
the E-step, the expected complete-data log likelihood with respect to the posterior
over sj is computed. In the M-step, this function is optimized with respect to the generators Gj (among other parameters), where the matrix exponential is approximatedlinearly. Once trained, the model can generate arbitrarily large transformations bycomputing the matrix exponential of an arbitrary linear combination of the generators, but training itself can only be done using infinitesimally small transformations
due to the linear approximation.
[Sohl-Dickstein et al., 2010] avoid having to compute a full matrix exponential, byusing an eigendecomposition:

exp As = exp WΛW−1 s = W exp (Λs)W−1 .

(3.4)

where A = WΛW−1 denotes the eigen decomposition of A into a complex matrix
W and a complex diagonal matrix Λ. The matrix exponential of a diagonal matrix
Λs is simply the elementwise exponential of the diagonal entries. This is much fasterthan computing a full matrix exponential as in the model of [Miao and Rao, 2007],
but precludes a simple linear combination of generators in the Lie algebra, without
recombining A = WΛW−1 . Instead, [Sohl-Dickstein et al., 2010] use a multiplicative
combination of one-parameter groups in a fixed order:!
y=Y

Wj exp (Λj sj )Wj−1

x.

(3.5)

j

The parameters and latent variables are obtained by gradient descent on an objective
based on 3.5, with additional smoothing and regularization terms.

29

3.3

Basic principles for representation learning

In this chapter’s introduction, we looked at the goals of representation learning, andsaw that there are currently no solid foundations on which to build a theory for thisfield. With a clear picture of the goals of- and current approaches to representation
learning, and a rudimentary understanding of group theory we are now ready to lay
down some basic principles which will form the basis of our work in later chapters,
and which we hope others will find useful.
The object classes we recognize do not partition the perceptual space arbitrarily.
If we are to build representation learning algorithms that can discover these classes
without much supervision, we should try to understand the structure of this partitioning. As we have seen in the introduction, all equivalence relations (partitionings
of a set, classes) correspond to a group action. If we can learn the effect of our own
movements (or more generally “actions”) on the input space, we have defined theobject classes implicitly:
The most natural proto-classes are defined by a group G and its action on a data
space X, as the points in the orbit space X/G.
Typically, we take take X = RN to be some feature space consisting of pixels, optical flow vectors, the output of a representation learning algorithm, etc. –see [Kanatani, 1990] for many examples of group representations in such feature
spaces. As discussed in section 1.4, when the orbits are aligned with the actualclass-boundaries, a representation in terms of orbits can drastically reduce the sample complexity, because we can learn the label of an entire orbit from one or a few
examples.
Besides volitional actions, we can consider actions that are unobserved or externally generated. In the next chapter, we show that assuming the changes are caused
by a toroidal group acting linearly on the input space, we can learn the group representation without being told which transformation was responsible for each observedchange. In general we will want to learn a group action that can explain as manyobserved changes as possible using as few generators as possible. That is, we want thegroup to act transitively on high-density regions of input space (typically manifolds)
so that we will often have some element of the group that can explain an observedchange. At the same time, the group should not contain elements that can transform a high-density observation into a very low-density observation (at least for smalltransformations).
30

The ultimate completion of this program (learning to explain every observedchange, and considering two inputs equivalent when they are in one orbit) leads
to a system that will see every input as one and the same thing. However, nothing
prevents us from performing analysis (“to break up into parts”), i.e. to consider only
subgroups at a time, and to consider the classes induced by these subgroups.
When we observe a change x → y, and seek to explain it using some elementof a group g ∈ G, there can be many candidates that are consistent with this observations. Typically, the set of candidates forms a group. A special case of this iswhen we ask which transformations g ∈ G can transform x into itself. These are thesymmetries of x with respect to G. It is clearly desirable to know the ambiguities in
our transformation estimate, and to describe this set of possibilities (using a posterior
distribution over the group elements, for example). So we have,
A representation of data x ∈ X should include a characterization of the symmetries of x with respect to the learned group G.
A concrete example that shows why this is useful can be found in the widely
studied problem of estimating camera motion and 3D scene structure from point
correspondences. As shown in [Ma et al., 1999], typical camera motions allow us toreconstruct the scene geometry and camera motion only up to subgroups of the Euclidean group. The utility of characterizing the ambiguity is obvious in this example,
but it is no less useful in the case of a learned group representation acting on highdimensional observations instead of 3D points.
Up until now, we have ignored the issue of invertibility. Many changes of practicalrelevance are not strictly invertible when one considers only the observation at one
instant in time. For example, by translating an image, pixels fall of on one side andthese cannot be recovered through an inverse translation. While this can be seen asa limitation of the proposed approach, we believe it is actually a strength. We canalways assume that non-invertibility implies incomplete information; it signals theneed to introduce latent variables or invoke some memory mechanism to make thetransformation on the union of observed and latent variables invertible. Our notionof Euclidean space is an example of this. We never directly observe Euclidean space,in the sense that sensory measurement values do not come in the form of points in a
3D space. However, if we construct a more abstract representation that does involvesome kind of 3D-description and remember past observations, any Euclidean motion
in a scene becomes invertible on this augmented space (and thus an understanding31

of object permanence is achieved). This idea is taken to the extreme in the field ofphysics, where all quantities of interest are inferred, and all transitions are invertible
at the elementary level.
We can take this idea of molding the representation Φ(x) to conform to priorassumptions of group structure one step further. As we have seen in chapter 2, anabstract (Lie) group has infinitely many representations, so we can pick one with
useful computational properties – as long as we can learn a way to embed x into a
representation space Φ(x) that transforms according to the group representation of
choice. One important computational requirement is that we can compute inverses
cheaply, which suggests that we seek a data representation that makes the grouprepresentation orthogonal. All compact groups (and more generally the amenablegroups) are unitarizable, and hence orthogonalizable (by a procedure described in
section 11.4.1 of [Doran and Lasenby, 2003]). Orthogonal transformations leave invariant a (not necessarily Euclidean) metric in the representation space, and so they
are sometimes called isometries.
In the representation space, symmetries should become isometries.
We have taken this maxim from [Dorst et al., 2007], where it is applied to model
Euclidean and conformal geometry. Besides free invertibility, this approach has theimportant advantage that all metric and non-metric geometric relations between
the represented objects in this geometry are covariantly preserved by the symmetry group (as long as these relations can be expressed in geometric algebra; see[Dorst et al., 2007] chapter 17 for details).
The next problem is untangling. The group theoretic point of view provides aclear principle for disentangling, which we have called Weyl’s principle [Weyl, 1939]
after [Kanatani, 1990]:
A k-tuple zj of latent variables can rightly be called a single observable if andonly if it defines an inequivalent irreducible representation of a symmetry group of
the phenomenon being observed.
By “symmetry group of the phenomenon” we are referring to the previous principles: it can be anything we wish, but for practical purposes we use a group that canoften be invoked to explain observed transformation sequences, and that appears toa symmetry of one or more functions we are interested in. By “defining an irrep” we32mean that the vectors zj (for j = 1 . . . J) are obtained by a linear change of variablesfrom the measurements x ∈ RN , such that all elements of the group only transformthe elements of each zj among themselves. See section 2.2.3 for details on irreducible
representations.
More generally, the goal can be to separate the variables into sets that transform
independently (as an ordered set), using some class of non-linear functions. This canin principle be done, as was shown by [Kanatani, 1990] for camera rotation reconstruction problems. Kanatani achieved this non-linear untangling by analysis, andnot by a learning algorithm. To our knowledge non-linear untangling by something
like Weyl’s principle has never been tried in machine learning.
Weyl’s principle also provides a new way to understand and construct convolutional neural networks [LeCun and Bottou, 1998], [Serre et al., 2007]. This classical architecture is characterized by two features: convolution and pooling. Several
learned filters are convolved with the image, which is equivalent to saying the weights
of filters at each location in the image are shared. After a layer of filtering, thenetwork locally pools several filter response values at each location, by squaring andsumming, or by taking maxima. This can be done within one feature map or between
feature maps.
As was discussed in section 2.2.5, the sum of two relative invariants is a relative
invariant if and only if they define equivalent representations, i.e. if they have thesame weight. So if zj = Wj T x (for some matrix of filters W) defines an irreducible
representation, we may pool the values in zj at nearby locations. We should neversum-pool the values of inequivalent representations, because the result will not be an
invariant. However, we can perform max-pooling among absolute invariants computedfrom inequivalent representations, because the result is again an invariant.

33

Chapter 4
Toroidal Subgroup Analysis
In this chapter we describe a model and an algorithm for learning toroidal groups
from examples. These groups are a natural place to start developing a catalog oflearnable group structures, because their defining properties (compact, connected,
commutative) make them easy to deal with, both analytically and computationally.

4.1

The TSA model

4.1.1

The likelihood function

The TSA model assumes that observed pairs of vectors {x(i) , y(i) }i=1...N are each related by some unknown transformation from a toroidal subgroup TM (W) of SO(2M ).
Recall (section 2.2.4) that every transformation Q in a toroidal group TM (W) can be
written as Q = WR(ϕ)WT , where Q and W are 2M × 2M orthogonal matrices, andR(ϕ) is block-diagonal with M rotation blocks of size 2 × 2. The matrix W performs
the reduction of TM (W) into irreducible representations, which are the 2 × 2 blocksof R(ϕ). We will parameterize the group in reduced form, by considering W as theparameters to be learned. The toroidal group is determined by W, while a particulartransformation from the group is identified with the phase-vector ϕ(i) . This gives thefollowing data model (suppressing indices i from now on):y = WR(ϕ)WT x + ,

(4.1)

where ϕ = (ϕ1 , . . . , ϕM )T is a latent vector of phase-variables (one vector per datapair) and  is isotropic Gaussian noise with mean zero and variance σn2 . It followsthaty ∼ N (WR(ϕ)WT x, σn2 ).

34

(4.2)

v[j]
W2j

ϕj

u[j]

W2j−1

Figure 4.1: Likelihood function p(vj |uj ).
By design, the transformation Q is more easily described in the W-basis. We
define u(i) = WT x(i) and v(i) = WT y(i) , so that equation 4.2 can be written asv ∼ N (R(ϕ)u, σn2 ).

(4.3)

In general, obtaining the pdf of a transformed random vector involves scaling by a
Jacobian determinant factor (see [DeGroot and Schervish, 2002], section 3.9), but inour case this reduces to | det(W)| = 1.
The distribution over v has diagonal covariance (for it is again isotropic), so it
factors into M independent distributions over each invariant subspace:p(v|u) =MY
j=1

p(vj |uj ) =MYN (vj |R(ϕj )uj , σn2 ),

(4.4)

j=1

where we use a bold index j to denote the part of u or v in the j-th invariant 2D
subspace (as opposed to simply taking the j-th coordinate vj ) and R(ϕj ) is a 2 × 2
rotation matrix that rotates by angle ϕj . The distribution p(vj |uj ) is shown in figure
4.1.

35

0.8

0
0.5

0.7

1.0
2.0

0.6

4.0
0.5

0.4

0.3

0.2

0.1

−π

0.0

−π/2

0

π/2

π

Figure 4.2: The von Mises distribution for µ = π/3 and various values of κ (shownin the legend).

4.1.2

The prior over phases

For each data pair, we have M latent phase variables ϕj , j = 1 . . . M that describe therotation in each invariant subspace. Since these are periodic variables, we use the vonMises distribution ([Mardia and Jupp, 1999, Mardia, 1975]) which we denote by VM.
The von Mises distribution is the maximum entropy distribution for a periodic variable with given circular mean and circular variance [Jammalamadaka and Sengupta, 2001].
It is given byp(ϕj |µ, κ) = VM(ϕj |µ, κ) =eκ cos(ϕj −µ),2πI0 (κ)(4.5)

where I0 is the order-0 modified Bessel function of the first kind [Abramowitz and Stegun, 1965],µ is the preferred angle and κ is a precision parameter. We have plotted the von Mises
distribution for µ = −π/3 and various values of κ in figure 4.2.
Q
We assume (for now) that ϕj are independent, i.e. p(ϕ) = j p(ϕj ). The conditional independence structure is shown in figure 4.3.

4.1.3

The posterior distribution

For EM learning (section 4.2), we need to evaluate the posterior distribution overthe latent variables ϕj given the data x, y. The posterior distribution in invariant
36

µ

κ

σn

ϕ

y

W

M
x
N

Figure 4.3: TSA as a probabilistic graph. Shaded nodes indicate observed variables,
unshaded ones denote latent variables. The plates denote repetition of the variable
for each subspace j = 1, . . . , M and data pair i = 1, . . . , D. It is also possible to haveseparate variables κj and µj for each subspace j.
subspace j isp(ϕj |uj , vj ) = R 2π

p(vj |uj , ϕj )p(ϕj )p(vj |uj , ϕj )p(ϕj )dϕj
1
2
∝ exp κ cos(ϕj − µ) − 2 kvj − R(ϕj )uj k .
2σn
0

(4.6)

While it is not obvious from the form of eq. 4.6, the posterior is again a von Mises

distribution. That is, p(ϕj |uj , vj ) ∝ exp κ0j cos(ϕj − µ0j ) . As far as we know, thisconjugacy relation has not been published before. Conjugacy relations greatly simplify sequential Bayesian inference, because the form of the distribution is invariant
under Bayesian updating and therefore remains tractable indefinitely.
We list the update equations for µ and κ below. Because the computation issomewhat tedious, the details are left for Appendix A. Leaving out index j forreadability, writing v for vj ,s2κσn2 (vT R(µ)u) + kvk2 kuk2
κ0 =+ κ2 ,
σn4
(4.7)


2
(v
u
−
v
u
)
+
κ
σ
sin(µ)
1221nµ0 = tan−1.
(v1 u1 + v2 u2 ) + κ σn2 cos(µ)
We can gain some insight into eq. 4.7 by using the outer product v ∧ u fromgeometric algebra [Dorst et al., 2007]. When applied to vectors, this product constructs an element called a 2-blade that represents the weighted subspace spanned by
its arguments. The weight of a blade is analogous to the magnitude of a vector, andcan be retrieved by geometric multiplication by a unit-weight 2-blade of the sameorientation, here given by I = W2j ∧ W2j+1 , the unit 2-blade representing span of theunit vectors W2j , W2j+1 .
37

The coordinates of a 2-blade v ∧ u in an orthonormal frame can be computed as
the upper triangle of vuT − uvT [Bethge et al., 2007]. In the 2D case, this is a single
scalar. The numerator and denominator in the equation for µ0 then become(v1 u2 − v2 u1 ) = I(v ∧ u)

= kvkkuk sin(θ)(v1 u1 + v2 u2 ) = (v · u)

= kvkkuk cos(θ).

(4.8)

Equation 4.8 shows the symmetry in eq. 4.7, and makes it clear that to computethe posterior mean, we should convert the prior mean µ and observed angle θ tounit vectors (cos(µ), sin(µ))T and (cos(θ), sin(θ))T , compute a weighted average of
those vectors and convert the result back to an angle using the arctangent function.1
The averaging weights are the strength of our prior belief in µ, the precision κ, andthe strength of belief we should have in θ, which is kukkvkσn−2 . This geometricallyintuitive result follows simply from consistent application of Bayes rule, as shown inAppendix A.
The update equation for the precision parameter κ captures the propagation of
measurement noise  through the angle computation process. When the vectors u
and v in eq 4.7 are short relative to the noise variance σn2 , or when the prior µ does
not fit well to the observed data, the posterior precision will be small. Another way
to read the κ-update in equation 4.7 is that we must divide the data vectors u andv by the standard deviation of the noise: p = u/σn and q = v/σn . The update thenbecomes
κ0 =p
2κ qT R(µ)p + kqk2 kpk2 + κ2 ,(4.9)

This obviates the need for excluding vectors uj , vj whose norm is below somearbitrary threshold, as is necessary in the model of Cadieu & Olshausen to avoid
instabilities [Cadieu and Olshausen, 2012]. TSA doe not exclude any µj , but instead
provides each one with a confidence κj .

4.1.4

Relation to discrete Fourier transformIf we make no prior assumptions, i.e. if we set κ = 0, we can simplify even further:κ0j = kqj kkpj k.

(4.10)Let us assume the filters W are sinusoids. As we will demonstrate in chapter
6, these emerge when TSA is trained on image shifts. Taking the input vectors toconstitute discrete samplings of a 1D signal, we can describe the discretely sampled
1

In practice, the atan2 function is preferred.

38

sinusoid filters analytically as follows. Let Wj have columns equal to the real andimaginary parts of ω −jk (for k = 1 . . . N ), where ω = e2πi/N is the N -th primitive rootof unity. The DFT of x is given by:
Xj =

N
−1
X

xn ω −jn .

(4.11)

n=0

So the real and imaginary part of Xj are given by the two components of Wj T x forthe first N/2 frequencies. The other half are complex conjugates of the first when x
is real, and so they carry no information.
Now suppose we are interested in the transformation taking some arbitrary fixed
vector e = ( √12 , . . . , √12 )T to p = WT x. The vector e is to be interpreted in the
W-basis, given by sinusoids. By eq. 4.10, we find κ0j = kpj k = kWj T xk = |Xj | (themodulus of the Fourier transform). Summarizing: using a flat prior and a Fourier
basis for W, the TSA representation of the transformation from e to x in terms of µj
and κj is equal to the DFT representation of x in terms of complex phase arg(Xj ) andamplitude |Xj |. So the model we have specified provides a probabilistic interpretationof the DFT, and generalizes it.

4.2

Learning by expectation-maximization

The model defined in the previous section contains a latent vector of angles ϕ that
relates each input pair x, y. We saw that given some estimate of W, we can computethe posterior over ϕ. As we will see in this section, given the posterior over ϕ, we can
optimize the expected likelihood of the data with respect to W. Situations like thisare a common occurrence, and are typically solved by the Expectation-Maximization(EM) algorithm [Dempster et al., 1977].
The EM algorithm finds a maximum-likelihood estimate of the parameters in thepresence of latent variables, by iterating through an E-step and an M-step. In the
E-step, the posterior over latent variables p(ϕ|x, y) given the current estimate of theparameters W is evaluated. In the M-step, the expected complete-data log-likelihood
under the posterior over latent variables is optimized with respect to the parameters.
The details for the TSA models are given below.
We work with a data set {(x(i) , y(i) )}i=1...D of pairs of vectors. We will sometimes
store the observed and latent variables column-wise in matrices X, Y, Φ. We assume
the pairs i, i0 are independent from one another.

39

4.2.1

E-stepIn the E-step, we compute the following expectation, which is to serve as an objectivein the M-step:
Q(W, W0 ) = E [ln p(Y, Φ|X, W)]Φ|X,Y,W0X=
E[ln p(y(i) |ϕ(i) , x(i) , W)p(ϕ(i) |µ, κ)]ϕ(i) |x(i) ,y(i) ,W0 .

(4.12)iwhere W is the optimization variable and W0 is our current parameter estimate,which is fixed.
Only the likelihood term p(y(i) |ϕ(i) , x(i) , W) in eq. 4.12 depends on W, so we canP(i) (i)(i)instead maximize2
i E[ln p(y |ϕ , x , W)]ϕ(i) |x(i) ,y(i) ,W0 . Taking logarithms andwriting out the expectation, we find this to be equal to
#
"
Z
Y
1 X
(i)
(i)
− 2
VM(ϕj |µj , κj ) ky(i) − WR(ϕ)WT x(i) k2 dϕ,(4.13)
2σn i ϕ∈TM j
where we integrate over the whole torus TM = {(ϕ1 , . . . , ϕM ) | ϕj ∈ [0, 2π]}. NoQ(i)(i)
tice that the product of von Mises distributions j VM(ϕj |µj , κj ) is equal top(ϕ|x(i) , y(i) , W0 ), although the dependencies on x, y and W0 are implicit in the
(i)

(i)

use of posterior parameters κj and µj .
As shown in Appendix B.2, we can remove terms that do not depend on W andsimplify, yielding the following maximization objective:#Z 2π "
(i)
(i)
X
exp
(κ
cos(ϕ
−
µ
))
j
j
j
y(i)T Wj
R(ϕj ) dϕj Wj T x(i)(i)
2πI0 (κj )0i,j=X I1 (κ(i)
j )(i)i,j I0 (κj )(4.14)(i)y(i)T Wj R(µj )Wj T x(i)where Wj denotes columns 2j − 1 and 2j of W spanning subspace j, and R(ϕj ) is a
2 × 2 matrix that rotates by an angle of ϕj .
The result is rather intuitive. For each data pair i and invariant subspace j, we(i)rotate the part of x(i) that lies in subspace j using the posterior mean angle µj andthen compute a dot product with the part of y(i) that lies in subspace j. This will
be at a maximum when the rotated x is fully aligned with y, which is what we want(i)to achieve. The ratio of Bessel functions

I1 (κj )(i)

I0 (κj )2, plotted in figure 4.4, acts to weigh

The term “E-step” is somewhat of a misnomer for us, because we do not evaluate the true
expectation but instead evaluate an objective that is proportional to it.

40

1.0

0.8

I1 (κ)/I0 (κ)

0.6

I0 (κ) ∗ e−|κ|I1 (κ) ∗ e−|κ|

0.4

0.2

0.0
0

4

2

κ

6

8

10

Figure 4.4: Ratio of Bessel functions I1 (κ)/I0 (κ), computed from the exponentiallyscaled Bessel functions.
the terms so that those terms with low uncertainty (high κ) are weighted higher
than those with high uncertainty (low κ). In practice, it is important for numericalstability to use the exponentially-scaled Bessel functions (available in most scientificcomputing libraries) to compute this ratio:
I1 (κ) exp (−|κ|)I1 (κ)=.
I0 (κ)
I0 (κ) exp (−|κ|)

4.2.2(4.15)

M-stepTo optimize the objective in eq. 4.14, we employ a simple gradient ascent scheme.
The gradient is(i)i
X I1 (κj ) h
d
(i)
(i) T
0
(i) (i)T(i) (i)TQ(Wj , W ) =x y Wj R(µj ) + y x Wj R(µj ) .
(i)
2
dWj
i I0 (κj )σn(4.16)
Clearly, Wj = 0 will make the gradient zero, but this is typically not a maximum
and does not take into account the orthogonality constraint on W. So we optimize
Q by gradient ascent, and use a simple SVD procedure to enforce the orthogonality
constraint. After each gradient update, we compute USVT = W and then set
W := UVT . More sophisticated methods for optimization over the orthogonal groupexist [Edelman et al., 1998, Plumbley, 2004], but our approach is simple and works
well in practice. Alternatively, we can normalize the length of each column of W andorthogonalize column 2j against column 2j + 1 by a single step of Gram-Schmidt.
41

This approximate procedure appears to work reasonably well; the subspaces learn tobecome near-orthogonal although this constraint is not strictly enforced.
In practice, we do not completely optimize the objective in each M-step, butinstead perform a single gradient-based update per M-step. Such a procedure has thesame convergence guarantees as regular EM, and works well in practice.

42

Chapter 5
Low-dimensional subgroups
It is often the case that the transformations that act on the data can be described with
far fewer parameters than the M angles used in the basic TSA model of the previous
chapter. For example, when rotations act on a data set consisting of image patches,
the transformations can be described with a single parameter: the rotation angle.
When (cyclic) translations act on the data, we should need only two parameters:
the amount of shift in each direction. If the data set consists of images of humanfaces with changing facial expressions, the transformations could be described byparameters like changes in smiling, opening-closing of the eyes, raisedness of theeyebrows, etc [Susskind et al., 2011].
These examples show that we are typically interested in groups that are of lower
dimensionality than the M -dimensional Toroidal groups. Since every compact, connected commutative Lie-group is a subgroup of a maximal toroidal group, we canfirst learn the toroidal group and then search for the appropriate low-dimensional
subgroup. The subspace angles ϕj are coordinates in the Lie algebra of TM (W), so
instead of looking for subgroups, we can look for subalgebras. The advantage of thisapproach is that we can seek linear subspaces in the Lie algebra instead of non-linear
manifolds in the space of matrices.

5.1

The stabilizer subgroup

For a data point x, we can ask the question: what is the set of transformations fromthe learned toroidal group TM (W) that will leave x unchanged? This set is a groupcalled the stabilizer subgroup of x in TM (W),stab(x) = {Q ∈ TM (W) | Qx = x}.

43

(5.1)

For a given x, which transformations are in stab(x)? Since we are working witha toroidal group, the answer is pretty simple. When the component of x in subspace
j is zero, W·jT x = uj = 0, then we can rotate it by any angle and the result willbe the same (still 0). Conversely, when x has a non-zero projection uj , rotation by
any non-zero (mod 2π) angle will change the vector. So we see that a transformation
Q = WR(ϕ)WT is in stab(x) iff ϕj = 0 whenever uj 6= 0. In practice, we will
be mostly concerned with transformations Q for which the constraint Qx = x is
approximately satisfied, i.e. those that rotate in subspace j only when kuj k and theposterior precision κ0j are small.
The relevance of the stabilizer subgroup to our discussion is that we can only
determine the transformation in TM (W) that takes x to y modulo the stabilizer
subgroup. That is, if Qx = y for Q ∈ TM (W) then also QSx = y for S ∈ stab(x).
The stabilizer subgroup expresses the ambiguity in the transformation estimate. One
nice property of the TSA model is that it can naturally incorporate new evidence todiminish the (simultaneous) stabilizer subgroup down to {e}. That is, given a setof correspondence pairs {x(i) , y(i) }i=1...D , we can use Bayesian updating to get the
optimal toroidal rotation relating the x’s to the y’s (assuming they are all relatedby the same transformation, and the torus determined by W is fixed). In a more
traditional estimation setting, one would always need at least M correspondence
pairs in order to estimate the M angles, but by virtue of the conjugacy relation we
can easily incorporate new information sequentially and give a best guess at any point
in time, while being explicit about the remaining uncertainty.
The phenomenon of low-norm subspace projections appears in bilinear gating
models too, where it has been called the subspace aperture problem by Memisevic [Memisevic, 2012, Memisevic and Exarchakis, 2013]. In this perspective, the lownorm subspaces pose a problem because they make it hard to estimate the “true”
transformation from only one correspondence pair. The solution offered is to create
a new layer in the network that pools over multiple factors, thus sharing informationand yielding a representation of the transformation that is supposed to be independent of image content. But can this goal be achieved at all? What is the “true”
transformation that is to be represented, and can it always be retrieved from one
correspondence pair?
The group theoretical perspective provides clarity: if we are estimating an elementof an M -parameter group such as TM (W), it cannot generally be done given only one
correspondence pair, no matter what clever trick is used. We need M pairs. In orderto deal properly with the case where we have insufficient information to determine
44

the transformation uniquely with high precision, we should describe the uncertainty(both due to noise and the symmetry of the data) by a posterior distribution. If
the transformation group under consideration contains fewer parameters, we need
fewer data, so we should seek to diminish the number of parameters. Apparently thetoroidal group is not the real symmetry group of interest; it is merely the centralizer
in SO(N ) of the group we really are interested in. If the number of parameters to be
estimated equals one (as is the case for image rotations), this single parameter can be
estimated from one correspondence pair. Whether it is actually represented as a single
number or as a distributed representation as in the gating models is not relevant.
Already for two-parameter groups such as x, y-translation of an image, it becomes
impossible to estimate the transformation parameters from a single correspondence
pair. This leads to the well-known aperture problem, which says that for images
with translational symmetry, the x, y-translation parameters cannot be determineduniquely.
In the section 5.3, we discuss how subgroups can be learned in order to reduce theambiguity in the transformation estimates as obtained from a single correspondence
pair. However, a non-trivial stabilizer is not only a problem: it also provides valuable
information about the signal. We can ask, “what are the transformations in TM (W)that take x to x?” The answer is stab(x), characterized by κ0 inferred from the pair(x, x). So even if we are really interested in a subgroup G/TM (W), the symmetries ofx with respect to the centralizer of G in SO(N ) (a maximal torus, TM (W)), provide
a maximal invariant to transformations from G. For example, the equilateral triangle
in figure 1.2, call it x, has three rotational symmetries. If we rotate the image, ouruncertainty about the transformation that took x to x remains the same: it is aninvariant representation of the triangle.

5.2

Lie subalgebras and phase wrapping

When the subspace angles ϕj are dependent, as they must be when the correspondence
pairs are related by transformations from a non-trivial subgroup of the toroidal group,
we can try to predict the angles in underdetermined subspaces from those where kuj kis large. What model should we use for these dependencies? A cursory look at thedata would suggest that we need a fairly complex model, and indeed all related workwhere similar problems are encountered use high-dimensional filter banks to model
the relations between subspace angles [Cadieu and Olshausen, 2012, Memisevic, 2013](see section 3.2.2 and 3.2.3 for details).
45

6

5

5

4

4
z

y

6

3

6
5
4
3 z
2
1

3

2

2

1

1

1
2
x

3
4
5

0
0

1

2

3
x

4

5

6

0
0

1

2

3
x

4

5

6

6

1

2

3

4

5

6

y

Figure 5.1: Example of phase-wrapping. Points on the dark green line correspond tovectors ϕ in three dimensions. The x, y, z axes correspond to ϕ1 , ϕ2 , ϕ3 . The slopeof the line is (1, 2, 3), corresponding to the values of the weights a1 , a2 , a3 (see text).
Left, middle: projection onto the x-y and x-z planes, right: 3D x-y-z space.
Can we do better? Remember that while any group can serve as a group of symmetries, it must indeed be a group; requesting invariance to some arbitrary set oftransformations does not make sense. So while we are free to model high density
regions of ϕ-space, the set of transformations corresponding to this region must forma group. As we have seen in section 2.2.2, subgroups H of a toroidal group TM (W)correspond to linear subspaces h of the Lie algebra t of TM (W). Since the vector ϕ corresponds to an element on the Lie algebra t, linear subspaces of ϕ-space correspond toP
subgroups of TM (W) through the exponential map Q = exp j ϕj Bj = WR(ϕ)WT ,T
T
where Bj is the generator of rotations in subspace j: Bj = W·2j W·2j+1−W·2j+1 W·2j.
There is however one major complication. For a fixed h, any ϕ ∈ h corresponds toan element in H through the exponential map, but doing bottom up inference of ϕfrom correspondence pairs y = Qx will not necessarily yield ϕ ∈ h. The reason is that
an element of h may correspond to an angle ϕj that is outside [0, 2π], but bottom-up
inference will always yield an angle in that range. Consider a 1-D subalgebra of t,P
generated by an element A ∈ t: A = j aj Bj . The elements of the subgroup are P
given by H = exp sA = exp s j aj Bj . Whenever saj > 2π, the inferred phase ϕjwraps, and so the observed vectors ϕ will not lie on a line. Figure 5.1 shows a visual
example in 3 dimensions.
Roughly speaking, current algorithms such as [Memisevic, 2012] deal with thepiecewise linear space of angles using a large filter bank. A more parameter-efficientapproach that will also yield a more comprehensible representation of the transformation, is to learn the weights aj and introduce a new random variable corresponding to
the group parameter s. Recall from section 2.2.5, that the weights must be integers
46

in the case of a representation of SO(2), where s ∈ [0, 2π]. If we let the weights bereal numbers, the learned subgroup can become non-compact and non-periodic bychoosing weights that are irrational. (Of course, truly irrational numbers do not exist
in a finite machine, but nearly-irrational weights can still yield impractically large
periods.)
We have experimented with various approaches for learning the weights, bothin the real-valued case (which may be easier to deal with even if we know we are
only interested in integral weights) and integral case. The first approach we tried(section 5.3.1) is to explicitly represent the unwrapped phase variables, who in turnare connected to a group parameter. Another approach is to connect the means of thewrapped phases directly to a group parameter. We present this approach in section
5.3.2.
So far, these attempts have been moderately but not completely successful. Learning the weights reliably and correctly turns out to be surprisingly difficult. The maindifficulty is in the exponentially large search space. When each weight is allowed totake values 1, . . . , K, there are K M possible joint assignments to the weights a1 , . . . aM
for a single 1-parameter subgroup. For a 16 × 16 image patch and K = 10, this gives
12810 = 270 possibilities. For the continuous case, we can use gradient descent in the
128 dimensional space, but this space appears to have many local optima, due to theperiodicities of the objective function.
So far, we have only worked with 1-parameter subgroups, but it is conceptuallystraightforward to generalize to higher dimensional subgroups.

5.3

Learning one-parameter subgroups of a torus

5.3.1

Real-valued weights with explicit phase unwrapping

We introduce a new vector of real-valued parameters a = (a1 , . . . , aM )T that will
represent the weights of the irreducible representations. Next, we introduce a setof unwrapped phases ϕ̄j that determine the mean of the von Mises distribution overwrapped phases ϕj : p(ϕj ) = VM(ϕj |ϕ̄j , κ). The unwrapped phases are assumed to benormally distributed with mean determined by a and a latent variable s corresponding
to the group parameter: p(ϕ̄|s, aj , σ̄ 2 ) = N (ϕ̄|aj s, σ̄ 2 ). Since the group parameter
need not be periodic, we use a Gaussian prior on it: p(s) = N (s|µs , σs2 ).
To perform EM learning, we need to optimize an expectation with respect to theposterior p(ϕ, ϕ̄, s|x, y, a, W). This is probably not possible analytically, so we use a

47

3.0
p̃(ϕ̄)

envelope

2.5
2.0
1.5
1.0
0.5
0.0
−10

−5

5

0

ϕ̄

10

Figure 5.2: The unnormalized conditional distribution of ϕ̄ given its neighbors,
wrapped by a scaled Gaussian.
sampling approximationL1X
ln p(ϕ(l) , ϕ̄(l) , s(l) , y|x, Θ),E[ln p(ϕ, ϕ̄, s, y|x, Θ)]ϕ,ϕ̄,s|x,y,Θ0 ≈
L l=1

(5.2)

where Θ represents all the parameters in the model and the samples l are drawn
from the posterior over latent variables given data and old parameter estimates. Thisalgorithm is known as Monte Carlo EM (MC-EM), [Levine and Casella, 2001].
An easy way to obtain samples from the posterior is to use Gibbs sampling. Thatis, we iteratively sample from the conditionals of each variable given its neighbors intheir current state. The conditional p(ϕ|ϕ̄, x, y) is a von Mises distribution derived
using the update equations in eq. 4.7, using µ = ϕ̄j . The conditional p(s|ϕ̄) is a
Gaussian and,p(ϕ̄j |ϕj , s) ∝ exp (κ cos(ϕj − ϕ̄j ) −1(ϕ̄j − aj s)2 ).
2σ̄ 2

(5.3)

Because the periodic part can be bounded, this distribution can be enveloped by ascaled Gaussian, as shown in figure 5.2. This allows us to perform rejection samplingof ϕ̄.

5.3.2

Integral weights

Since we are interested in compact subgroups of the torus, it is arguably better toconstrain the aj to be integer valued. We will also connect the group parameter s
directly to the ϕ, making explicit phase-unwrapping unnecessary.

48

The prior on ϕj in the basic TSA model is replaced by a von Mises distribution
whose mean is determined by aj and s:p(ϕj |s, aj , κj ) = VM(ϕj |aj s, κj ).

(5.4)

The prior on s must now also be a periodic distribution, so we use a von Mises:p(s|ν, γ) = VM(s|ν, γ).
For EM learning, we again require samples from the posterior p(ϕ, s|x, y, Θ). It
is possible to integrate out ϕ, but the resulting density involves a Bessel function of
a complicated function of s, making rejection sampling difficult. We have not tried it
yet, but importance sampling may work well.
We now describe a Gibbs sampling scheme for this model. Suppressing parameters,
we have p(ϕ|x, y, s) ∝ p(ϕ|s)p(y|x, ϕ), which is a von Mises times a polar-Gaussian,and hence is again a von Mises with parameters µ0 , κ0 (as derived in appendix A).
The other conditional required for Gibbs sampling is!
p(s|ϕ) ∝ p(ϕ|s)p(s) ∝ exp γ cos(s − ν) +Xκj cos(ϕj − aj s) .

(5.5)

j

In general, there can be multiple aj with the same integer value, and we can groupthese together using the formula for a linear combinationA).

P of cosines (see appendix
F
This yields a distribution of the form p(s|ϕ) ∝ exp
j=1 ρj cos(j(s − δj )) , which
is known as a generalized von Mises distribution [Gatto and Jammalamadaka, 2007,
Gatto, 2008].
We can again use rejection sampling to sample from this unnormalized
distribu
PFtion, because it can be bounded by a constant function exp
j=1 ρj . If we wantto sample from the same distribution
 specified in the form of eq. 5.5, the boundingPM
constant becomes exp γ + j=1 κj . Rejection sampling then proceeds as follows:
1. Sample a proposal sample s0 ∼ U(0, 2π).
PM
2. Sample u0 ∼ U(0, exp γ + j=1 κj ). The pair (s0 , u0 ) is distributed uniformly
under the uniform proposal distribution.
3. Let p̃(s) denote the unnormalized distribution from which we wish to sample.
If u0 > p̃(s0 ), reject the sample and repeat from step 1. Otherwise accept s0 as
a sample.

49

(l)

Having obtained a set of samples {ϕ(l) , sl=1...L }, we optimize the expected complete
data log-likelihood,
L
L
1X
1X
1
ln p(ϕ(l) , s(l) , y|x, Θ) =− 2 ky − WR(ϕ(l) )WT xk2 − 2M ln σn (5.6)
L l=1
L l=1 2σn
M
X

(l)

(5.7)

+ γ cos(s(l) − ν) − ln 2πI0 (γ).

(5.8)

+

κj cos(ϕj − aj s(l) ) − ln 2πI0 (κj )j=1This objective must be summed over the whole data set i = 1 . . . D. The objective
is easily optimized by gradient descent with respect to σn , κj , ν, γ. The gradientwith respect to W is the same as derived in the previous chapter, except that it is
now summed over l. Optimization with respect to aj can be done by exhaustively
checking every value in some range aj = 0 . . . F . However, in practice, when κj islarge the samples of ϕj tend to be closely aligned with aj s so that the optimal valueof aj is not changed, and when κj is small the values of aj fluctuate randomly. We
have various alternatives, such as making aj a random variable and sample it (either
by itself or jointly with s), but have so far not found an algorithm that consistently
yields good results. We hope to improve this in the future, but for now we will only
perform experiments with the continuous-weight version of the algorithm presented
in section 5.3.1.

50

Chapter 6
Experiments
6.1

Random data

To check the correctness of the basic TSA algorithm described in chapter 4, we
generated random 18-dimensional vectors x(i) ∼ N (0, I), i = 1 . . . 10000 and rotatedeach one by a random matrix of the form(i)R(ϕ1 )
..
R(i) = 
.
(i)



,
(i)
R(ϕ9 )(6.1)(i)

where ϕ(i) ∼ U(0, 2π) and R(ϕj ) is a 2 × 2 rotation matrix that rotates by angle ϕj .
The matrices R(i) are elements of the torus T9 (I) since W = I. We initialized W ata random orthogonal matrix and trained the TSA model until convergence. Figure
6.1 shows the learned matrix W with high values corresponding to light pixels andlow values corresponding to dark pixels. As expected, there are two types of nonidentifiability in the matrix W. First, the blocks are permuted in a random orderand second, because there is no inherent ordering of 2D subspaces. Secondly, each2 × 2 block is orthogonal but not necessarily equal to the identity I2×2 , because only
the subspace but not the basis can be identified. Finally, notice that in each block,either the diagonal or anti-diagonal values are equal, while the other (anti-diagonal,
resp. diagonal) values are opposite. This corresponds to a pure rotation and an antirotation (also known as a roto-reflection), which is a reflection (which is orthogonal)
followed by a rotation.

51

0

5

10

15
0

5

10

15

Figure 6.1: Retrieved block-matrix W. Lighter pixels correspond to higher values.

6.2

Learning image filters

Next, we trained the model on randomly generated image patches, that were rotatedor cyclicly shifted by random amounts. To be precise, a random patch was sampledfrom an isotropic Gaussian, x(i) ∼ N (0, I) and then the corresponding patch y(i) was
generated from x(i) by rotation by a uniformly random angle θ ∼ U(0, 2π) or shiftedby sx , sy ∼ U(0, w) where w = 16 is the width of the patch. The result is shown infigure 6.2.

Figure 6.2: Filters learned from random rotations and shifts.

52pixel-space distance

−π

−21 π

0

1
2

π

π

Figure 6.3: Blue, dashed line: distance of a reference patch (center) to rotated copies
of it. Green, solid line: distance between a reference patch and the back-rotation by
WR(µ0 )WT of a rotated copy (where µ0 is the vector of inferred subspace angles).

6.3

Testing rotation invariance

To test the theoretically derived invariance properties of the representation, we performed the following experiment. We trained the basic TSA model on randomly
rotated image patches, as in section 6.2. Then, we selected 10.000 reference patchesfrom the MNIST data set of handwritten digits [LeCun and Bottou, 1998], scaled
them to 16 × 16 pixels, and rotated each of them by an angle that increased from −π
to π in 100 steps. For each one of the 100 angles θ, we computed the average distance
between the reference patches x(i) and their rotated copies y(i,θ) (i = 1 . . . 10000).
Next, we computed the posterior mean angle µ(i,θ) for each pair x(i) , y(i,θ) , and computed z(i,θ) = WR(µ(i,θ) )WT x(i) . Figure 6.3 shows the average distance kx(i) − y(i,θ) k
and kx(i) − z(i,θ) k for each angle θ.

6.4

Rotation invariant classification

We performed an experiment to see if the learned features are useful for classification. We classified rotated handwritten digits from the MNIST-rot data set[Larochelle et al., 2007], which were rescaled from 28 × 28 to 16 × 16, because we
did not have time to train a full basis W of size 784 × 784. This data set consistsof 10 classes of handwritten digits, which were rotated by random angles θ ∈ [0, pi].
The data set has 12000 training examples and 50000 test examples.
53

Figure 6.4: Classification results obtained on rescaled MNIST-rot, using the raw pixel
representation and the invariant κ representation.
Using a basis W trained as in section 6.2, we obtained an invariant representationfrom each digit x(i) by computing the posterior precision vector κ(i) of the transformation from x(i) to itself. The prior means and precisions were set to 0. As explained
in section 5.1, the vector κ(i) obtained in this way describes the symmetries of theimage with respect to the toroidal group. Since a transformation from the learned
group (which includes all image rotations) does not change the symmetries of the
figure, this representation is invariant.
Figure 6.4 shows classification results for various training set sizes and two data
representations: the invariant representation κ and the raw pixel values. These wereobtained using a random forest classifier with 100 trees. The invariant representation
performs better for every training set size, but the benefit is most pronounced forsmall training sets. It should be noted that our results are lower than those reported
in [Larochelle et al., 2007], both for the invariant and raw pixel representations, which
is most likely due to having rescaled the images.

6.5

Subgroup learning

We tested the subgroup-learning algorithm presented in section 5.3.1 on image rotations. The 2D rotation group SO(2), being compact, connected and commutative, is54Group parameter s
−3

−2

−1

1
0
angle (rad)2

3

Figure 6.5: Inferred s-value for various rotation angles. Each line corresponds to adigit image.
a subgroup of some maximal torus (which is in SO(N ), where N is the dimensionalityof the pixel space).
We trained the algorithm on random rotations of random 16 × 16 image patches,
reduced by PCA to 64 dimensions1 . Then, we chose one MNIST digit from eachclass 09, and rotated it by 100 angles in [−π, π]e. Figure 6.5 shows one line for each
digit, with the angle between the reference digit and the rotated copy shown on thex-axis, and the sampled value of s after 10 steps of Gibbs sampling on the y-axis.
Notice that the lines are largely overlapping, showing that the same representation(s-value) is used to represent the same rotation on different images. Furthermore,
notice that there is still a phase-wrap at π2 and − π2 , which shows that the model has
failed to make use of the lowest-frequency filters. This means that it will have thesame internal representation of two different rotations, and therefore will not be able
to apply this transformation properly to a new image. Other than that, the model
has linearized the 1-parameter group; it has discovered the natural concept “rotation
angle”.

1

The training data will indeed lie in a low-dimensional subspace, but this will not be so for new
randomly drawn images. However, as the results indicate it does work well enough on digit images.

55

Chapter 7
Conclusion and Outlook
We have presented an analysis of group theoretical methods in representation learning,
and a model and algorithm for learning representations toroidal groups. Our workon learning subgroups of a maximal toroidal group has so far been only moderately
successful, and we hope to improve this in the future. Doing so would make itpossible to create convolutional networks whose architecture is informed by grouprepresentation theory, thus opening the door for application of our work to more
serious problems such as object and action recognition in full-size images and videos.
A sample of other interesting topics for future work:
1. Modeling non-commutative groups, possibly using toroidal groups as building
blocks.
2. In our work so far we have assumed that the pairs (x(i) , y(i) ) are all independent,
but it is also possible to consider sequence models with the same conditional
independencies as HMM/LDS models but using emission probabilities definedby TSA.
3. Build sparse versions of TSA using a kurtotic prior on κ and an overcomplete
dictionary W.
4. A full Bayesian analysis, i.e. using priors over the parameters and samplingthem from their posterior (taking into account orthogonality constraints).

56

Appendix A
Derivation of von Mises-Polar
Gaussian conjugacy
Assume a von Mises prior on ϕ:p(ϕ|µ, κ) = VM(ϕ|µ, κ) =eκ cos(ϕ−µ)
2πI0 (κ)

(A.1)
(A.2)

where I0 is the order-0 modified Bessel function of the first kind. Furthermore, let u
and v be 2D vectors and assume a circularly parameterized Gaussian distribution onv:
p(v|u, ϕ, σ) = N (v|R(ϕ)u, σ 2 )
where

(A.3)




cos(ϕ) − sin(ϕ)R(ϕ) =sin(ϕ) cos(ϕ).

(A.4)

The posterior distribution in is given by Bayes rulep(v|u, ϕ)p(ϕ)
p(ϕ|u, v) = R 2πp(v|u, ϕ)p(ϕ)dϕ
0


1
2
∝ exp κ cos(ϕ − µ) − 2 kv − R(ϕ)uk
2σ

(A.5)

kv − R(ϕ)uk2 = vT v + uT u − 2vT R(ϕ)u.

(A.6)

Notice that

The terms in eq. A.6 that do not depend on ϕ can be factored out of the normalizing
integral of eq. A.5 as well as the numerator, and so they cancel. This gives us theposterior:p(ϕ|u, v) ∝ eκ cos(ϕ−µ)+σ57−2 vT R(ϕ)u(A.7)

It turns out this is again a von Mises distribution, although it is not obvious in thecurrent form. To get the update equations for µ and κ, we write out the simplified
exponent of in terms of sines, as follows:κ cos(ϕ − µ) + σ −2 vT R(ϕ)u= κ cos(ϕ − µ) + σ −2 (v1 cos(ϕ)u1 + v2 sin(ϕ)u1 − v1 sin(ϕ)u2 + v2 cos(ϕ)u2 )


1
1
1
−2
v1 u1 sin(ϕ + π) + v2 u1 sin(ϕ) − v1 u2 sin(ϕ) + v2 u2 sin(ϕ + π) .
= κ sin(ϕ − µ + π) + σ
2
2
2
(A.8)
Now we apply the following identity:Xai sin(ϕ + δi ) = a sin(ϕ + δ),

(A.9)

i

where a2 =P

ij ai aj cos(δi − δj ) andPai sin(δi )tan δ = P i,
i ai cos(δi )(A.10)

making the substitution and simplifying, we find the parameters of the von Mises
posterior:r2κσ 2 (vT R(µ)u) + kuk2 kvk2
+ κ2 ,4σ
(v1 u2 − v2 u1 )σ −2 + κ sin(µ)
0−1µ = tan.
(v1 u1 + v2 u2 )σ −2 + κ cos(µ)
κ0 =

58

(A.11)

Appendix B
EM Equations
B.1

Expected rotation matrix under von Mises distributed angle

Here we derive the formula for the expected 2 × 2 rotation matrix under a von Mises
distributed rotation angle. This result will be useful later.
Let ϕ be a scalar angle, andcos(ϕ) − sin(ϕ)R(ϕ) =sin(ϕ) cos(ϕ)
We show that

Z 2π
0I1 (κ)
exp (κ cos(ϕ − µ))R(ϕ) dϕ =R(µ).
2πI0 (κ)
I0 (κ)

59

(B.1)

(B.2)

First, substitute u = ϕ − µ:Z 2πexp (κ cos(u))
R(u + µ) du
2πI0 (κ)
0

Z 2π
1
exp (κ cos(u))R(u) du R(µ)=
2πI0 (κ) 0
Z 2π Z 2π 
1=exp (κ cos(u)) cos(u)du −exp (κ cos(u)) sin(u)du B R(µ)
2πI0 (κ)
00 Z 2π Z 2πd1exp (κ cos(u))d cos(u) B R(µ)
exp (κ cos(u))du +=2πI0 (κ)
dκ00 Z 2π1d=exp (κ cos(u))du R(µ)
2πI0 (κ) dκ 02π1exp (κ cos(u))BR(µ)+
2πI0 (κ)κ
01 d=
I0 (κ)R(µ) + 0
I0 (κ) dκ
I1 (κ)=R(µ),
I0 (κ)(B.3)
Where we defined


0 1
B=
(B.4)
−1 0
We also used the standard results:
Z 2πexp(κ cos(u))du = 2πI0 (κ)
0d
I0 (κ) = I1 (κ).
dκ

B.2

(B.5)

Derivation & optimization of EM objective

We now derive the equations required for EM learning in the basic TSA model.
Let X, Y, Φ, M , K be matrices containing x(i) , y(i) , ϕ(i) , µ(i) , κ(i) in the i-th
column. Let Wj denote the two columns in W spanning subspace j.
At step t of the algorithm, do:
(i)
(i)
(i)
E step: compute the parameters µj and κi of the posterior p(ϕj |x(i) , y(i) , Wj (t−1) , µ, κ),using the update equations derived in appendix A.
M step: we first define Q as the expected complete data log-likelihood under theposterior over ϕ, and then optimize it with respect to W.

60

Q(W, W(t−1) ) = E[ln p(Y|Φ, X, W)p(Φ|µ, κ)]X=E[ln p(y(i) |ϕ(i) , x(i) , W) + ln p(ϕ(i) |µ, κ)]

(B.6)

i

where the expectation is taken with respect to the posterior of ϕ, which depends on
W(t−1) .

B.2.1

Optimization with respect to W

Only the likelihood term has a dependence on W, so we leave out the prior term ofequation B.6:
Z(i) (i)(i)E[ln p(y |ϕ , x , W)] =p(ϕ(i) |µ(i) κ(i) ) ln p(y(i) |ϕ(i) , x(i) , W)dϕ(i)ϕ(i) ∈TM (W)(B.7)
Writing out the defining conditionals we find B.7 to equal
#
"
Z
(i)
(i)
(i)
Y exp (κ(i)− WR(ϕ(i) )WT x(i) k2 (i)j cos(ϕj − µj )) ky
dϕ .
−
(i)
2σn2
2πI0 (κj )ϕ(i) ∈TM (W)j(B.8)
We expand the squared norm to find the ϕ-dependent terms (suppressing iindices):ky − WR(ϕ)WT xk2 = yT y − 2yT WR(ϕ)WT x+ xT WR(ϕ)T WT WR(ϕ)WT x(B.9)= kyk2 + kxk2 − 2yT WR(ϕ)WT x
The last step assumes that W is orthogonal, and so we must enforce this constraint
later on. Since the first two terms kyk2 and kxk2 do not depend on ϕ, they can betaken out of the expectation. Since they don’t depend on W either, they do not affectthe position of the optimum and so we leave them out.
Since R(ϕ(i) ) is block diagonal, we can split the bilinear form into a sum of M
smaller bilinear forms:
yT WR(ϕ)WT x =XyT Wj R(ϕj )WjT x,(B.10)jwhere each term in the sum only depends on one of the M angles ϕj for data point i.
Leaving out the constant terms kxk2 and kyk2 and the constant factor 2σ12 out ofneq. B.8, we are left to maximize

61

"

XZ
i

=

ϕ(i) ∈TM (W)

ϕ(i) ∈TM (W)

#"
X

(i)
2πI0 (κj )

j

"

XZ
i,j

(i)
(i)
Y exp (κ(i)
j cos(ϕj − µj ))#(i)y(i)T Wj R(ϕj )Wj T x(i) dϕ(i)jY exp (κ(i) cos(ϕ(i) − µ(i) ))
l
l
l
(i)
2πI0 (κl )
l

#



y

(i)T

(i)
Wj R(ϕj )Wj T x(i)



dϕ(i)

#
"
(i)
(i)


X Z 2π exp (κ(i)
cos(ϕ
−
µ
))
j
j
j
(i)
(i)
T (i)
(i)T
dϕj
=
y
W
R(ϕ
)W
x
j
j
j
(i)
2πI0 (κj )
0
i,j
"
#
Z 2π
(i)
(i)
(i)
X
exp (κj cos(ϕj − µj ))
(i)
(i)
T
=
yi W jR(ϕj ) dϕj Wj T x(i)(i)
2πI0 (κj )0i,j=X I1 (κ(i)
j )(i)i,j I0 (κj )(i)y(i)T Wj R(µj )Wj T x(i) .
(B.11)The last step uses the result derived in section B.1.
Finally we evaluate the gradient:i
X I1 (κj(i) ) h
d
(i)T
(i)
(t−1)
(i) (i)T(i) (i)T
.
+
y
x
W
R(µ
x
y
W
R(µ
Q(Wj , Wj
)=
j
j
j
j
(i)
2
dWj
i 2I0 (κj )σn(B.12)(we reintroduced the factor 2σn2 that we left out before) Keep in mind that feasible
solutions must satisfy orthogonality, so we optimize this by projected gradient descent.

62

Bibliography[Abramowitz and Stegun, 1965] Abramowitz, M. and Stegun, I. A. (1965). Handbookof Mathematical Functions. Dover Pub.
[Adelson and Bergen, 1985] Adelson, E. H. and Bergen, J. R. (1985). Spatiotemporalenergy models for the perception of motion. J Opt Soc Am A, 2(2):284–299.
[Aristotle, 1924] Aristotle (1924). Metaphysics. Clarendon Press.
[Bengio, 2013] Bengio, Y. (2013). Deep learning of representations: Looking forward.
arXiv preprint arXiv:1305.0445.
[Bengio et al., 2013] Bengio, Y., Courville, A., and Vincent, P. (2013). Representation Learning: A Review and New Perspectives. IEEE transactions on patternanalysis and machine intelligence, pages 1–30.
[Bengio and Lecun, 2013] Bengio, Y. and Lecun, Y. (2013). International Conference
on Learning Representations.
[Bethge et al., 2007] Bethge, M., Gerwinn, S., and Macke, J. H. (2007). Unsupervised
learning of a steerable basis for invariant image representations. Proceedings of
SPIE Human Vision and Electronic Imaging XII (EI105), pages 64920C–64920C–
12.
[Bruna and Mallat, 2013] Bruna, J. and Mallat, S. (2013). Invariant scattering convolution networks. IEEE transactions on pattern analysis and machine intelligence,35(8):1872–86.
[Cadieu and Olshausen, 2012] Cadieu, C. F. and Olshausen, B. a. (2012). Learningintermediate-level representations of form and motion from natural movies. Neuralcomputation, 24(4):827–66.

63

[Cardoso, 1998] Cardoso, J. (1998). Multidimensional independent component analysis. Proceedings of the 1998 IEEE International Conference on Acoustics, Speechand Signal Processing, 4.
[Cox et al., 2005] Cox, D. D., Meier, P., Oertelt, N., and DiCarlo, J. J. (2005). ’Breaking’ position-invariant object recognition. Nature Neuroscience, 8(9):1145–1147.
[Curd, 2012] Curd, P. (2012). Presocratic Philosophy. The Stanford Encyclopedia of
Philosophy (Winter 2012 Edition), Edward N. Zalta (ed.),.
[Davis and Mermelstein, 1980] Davis, S. and Mermelstein, P. (1980). Comparisonof parametric representations for monosyllabic word recognition in continuously
spoken sentences. Acoustics, Speech and Signal Processing, IEEE Transactions on,
28(4).
[DeGroot and Schervish, 2002] DeGroot, M. H. and Schervish, M. J. (2002). Probability and Statistics. Addison Wesley, 3rd edition.
[Dempster et al., 1977] Dempster, A., Laird, N., and Rubin, D. (1977). Maximum
likelihood from incomplete data via the EM algorithm. Journal of the Royal Statistical Society. Series B (Methodological), 39(1):1–38.
[DiCarlo and Cox, 2007] DiCarlo, J. J. and Cox, D. D. (2007). Untangling invariant
object recognition. Trends in cognitive sciences, 11(8):333–41.
[DiCarlo et al., 2012] DiCarlo, J. J., Zoccolan, D., and Rust, N. C. (2012). How does
the brain solve visual object recognition? Neuron, 73(3):415–34.
[Doran and Lasenby, 2003] Doran, C. and Lasenby, A. (2003). Geometric algebra forphysicists. Cambridge University Press.
[Dorst et al., 2007] Dorst, L., Fontijne, D., and Mann, S. (2007). Geometric Algebra for Computer Science: An Object-Oriented Approach to Geometry. Morgan
Kaufmann Publishers Inc., San Francisco, CA, USA.
[Edelman et al., 1998] Edelman, A., Arias, T. a., and Smith, S. T. (1998). The
Geometry of Algorithms with Orthogonality Constraints. SIAM Journal on Matrix
Analysis and Applications, 20(2):303–353.
[Gatto, 2008] Gatto, R. (2008). Some computational aspects of the generalized vonMises distribution. Statistics and Computing, 18(3):321–331.
64

[Gatto and Jammalamadaka, 2007] Gatto, R. and Jammalamadaka, S. (2007). The
generalized von Mises distribution. Statistical Methodology, 4(3):341–353.
[Hinton, 1979] Hinton, G. (1979). Some Demonstrations of the Effects of Structural
Descriptions in Mental Imagery*. Cognitive Science, 3(3):231–250.
[Hinton et al., 2012] Hinton, G., Deng, L., Yu, D., Dahl, G. E., Mohamed, A.-r.,
Jaitly, N., Senior, A., Vanhoucke, V., Nguyen, P., Sainath, T. N., and Kingsbury,
B. (2012). Deep Neural Networks for Acoustic Modeling in Speech Recognition.
IEEE Signal Processing Magazine, 29(6):82–97.
[Holyoak and Morrison, 2005] Holyoak, K. J. and Morrison, R., editors (2005). The
Cambridge Handbook of Thinking and Reasoning. Cambridge University Press,
Cambridge, UK.
[Hyvarinen and Hoyer, 2000] Hyvarinen, A. and Hoyer, P. (2000). Emergence ofPhase- and Shift-Invariant Features by Decomposition of Natural Images into Independent Feature Subspaces. Neural Computation, 1720:1705–1720.
[Hyvärinen and Oja, 2000] Hyvärinen, A. and Oja, E. (2000). Independent component analysis: algorithms and applications. Neural networks : the official journalof the International Neural Network Society, 13(4-5):411–30.
[Jammalamadaka and Sengupta, 2001] Jammalamadaka, S. R. and Sengupta, A.
(2001). Topics in Circular Statistics. World Scientific Pub Co Inc, har/dskt edition.
[J.P. Elliot, 1985] J.P. Elliot, P. D. (1985). Symmetry in Physics: Principles andSimple Applications, Volume 1. Oxford University Press, USA.
[Kanatani, 1990] Kanatani, K. (1990). Group Theoretical Methods in Image Understanding. Springer-Verlag New York, Inc., Secaucus, NJ, USA.
[Klein, 1893] Klein, F. (1893). A comparative review of recent researches in geometry.
Bulletin of the American Mathematical Society, 1.
[Konda et al., 2013] Konda, K., Memisevic, R., and Michalski, V. (2013).
role of spatio-temporal synchrony in the encoding of motion.

The

arXiv preprint

arXiv:1306.3162, pages 1–9.
[Kondor, 2008] Kondor, R. (2008). Group theoretical methods in machine learning
Risi Kondor. PhD thesis, Columbia University.
65

[Krizhevsky et al., 2012] Krizhevsky, A., Sutskever, I., and Hinton, G. (2012). ImageNet classification with deep convolutional neural networks. Advances in Neural
Information Processing Systems, pages 1–9.
[Larochelle et al., 2007] Larochelle, H., Erhan, D., Courville, A., Bergstra, J., andBengio, Y. (2007). An empirical evaluation of deep architectures on problems withmany factors of variation. Proceedings of the 24th International Conference on
Machine Learning (ICML’07).
[Le et al., 2012] Le, Q., Ranzato, M., Monga, R., Devin, M., Kai, C., Corrado, G. S.,
Dean, J., and Ng, A. Y. (2012). Building high-level features using large scale unsupervised learning. Proceedings of the 29th International Conference on Machine
Learning.
[Le et al., 2011] Le, Q. V., Zou, W. Y., Yeung, S. Y., and Ng, A. Y. (2011). Learninghierarchical invariant spatio-temporal features for action recognition with independent subspace analysis. In Computer Vision and Pattern Recognition (CVPR),
2011 IEEE Conference on, pages 3361–3368. Comput. Sci. Dept., Stanford Univ.,
Stanford, CA, USA, IEEE.
[LeCun and Bottou, 1998] LeCun, Y. and Bottou, L. (1998). Gradient-based learning
applied to document recognition. Proceedings of the IEEE.
[Leibo et al., 2011] Leibo, J., Mutch, J., and Poggio, T. (2011). Learning to discounttransformations as the computational goal of visual cortex. Nature Precedings,
pages 1–3.
[Levine and Casella, 2001] Levine, R. and Casella, G. (2001). Monte Carlo EM Algorithm. Journal of Computational & Graphical Statistics, 10(3):422–439.
[Li and DiCarlo, 2008] Li, N. and DiCarlo, J. J. (2008). Unsupervised Natural Experience Rapidly Alters Invariant Object Representation in Visual Cortex. Science,
321(5895):1502–1507.
[Lowe, 2004] Lowe, D. G. (2004). Distinctive Image Features from Scale-Invariant
Keypoints. International Journal of Computer Vision, 60(2):91–110.
[Ma et al., 1999] Ma, Y., Soatto, S., Kosecka, J., and Sastry, S. (1999). Euclidean
reconstruction and reprojection up to subgroups. In Computer Vision, 1999. The
Proceedings of the Seventh IEEE International Conference on, volume 2, pages
773—-780 vol.2. California Univ., Berkeley, CA, IEEE.
66

[Mardia, 1975] Mardia, K. (1975). Statistics of directional data. Journal of the Royal
Statistical Society. Series B (Methodological), 37(3):349–393.
[Mardia and Jupp, 1999] Mardia, K. V. and Jupp, P. E. (1999). Directional Statistics.
John Wiley & Sons, 1 edition.
[Memisevic, 2011] Memisevic, R. (2011). Gradient-based learning of higher-order
image features. 2011 International Conference on Computer Vision, pages 1591–
1598.
[Memisevic, 2012] Memisevic, R. (2012). On multi-view feature learning. International Conference on Machine Learning.
[Memisevic, 2013] Memisevic, R. (2013). Learning to Relate Images. IEEE transactions on pattern analysis and machine intelligence, pages 1–19.
[Memisevic and Exarchakis, 2013] Memisevic, R. and Exarchakis, G. (2013). Learning invariant features by harnessing the aperture problem. International Conference
on Machine Learning, 28.
[Memisevic and Hinton, 2010] Memisevic, R. and Hinton, G. E. (2010). Learning
to Represent Spatial Transformations with Factored Higher-Order Boltzmann Machines. Neural Computation, 22(6):1473–1492.
[Miao and Rao, 2007] Miao, X. and Rao, R. P. N. (2007). Learning the Lie groupsof visual invariance. Neural computation, 19(10):2665–93.
[Mohamed et al., 2012] Mohamed, A.-r., Dahl, G., and Hinton, G. (2012). Acoustic
modeling using deep belief networks. IEEE transactions on Audio, Speech, andLanguage Processing, 20(1):14–22.
[Nathan Carter, 2009] Nathan Carter (2009). Visual Group Theory website.
[Piaget, 1952] Piaget, J. (1952). The Origins of Intelligence in Children. International
University Press, New York.
[Plumbley, 2004] Plumbley, M. (2004). Lie group methods for optimization withorthogonality constraints. Proceedings of the International Conference on Independent Component Analysis and Blind Signal Separation.
[Poincaré, 1904] Poincaré, H. (1904). Science and hypothesis. The Walter Scott
Publishing Co.
67

[Rao and Ruderman, 1999] Rao, R. and Ruderman, D. (1999). Learning Lie groups
for invariant visual perception. Advances in neural information processing systems,
816:810–816.
[Serre et al., 2007] Serre, T., Wolf, L., Bileschi, S., Riesenhuber, M., and Poggio, T.
(2007). Robust object recognition with cortex-like mechanisms. IEEE transactions
on pattern analysis and machine intelligence, 29(3):411–26.
[Soatto, 2009] Soatto, S. (2009). Actionable information in vision. Machine learningfor computer vision, pages 2138–2145.
[Sohl-Dickstein et al., 2010] Sohl-Dickstein, J., Wang, J., and Olshausen, B.
(2010).

An unsupervised algorithm for learning lie group transformations.

http://arxiv.org/abs/1001.1027.
[Susskind et al., 2011] Susskind, J., Hinton, G., Memisevic, R., and Pollefeys, M.
(2011). Modeling the joint density of two images under a variety of transformations.
Cvpr 2011, pages 2793–2800.
[Weyl, 1939] Weyl, H. (1939). The classical groups: their invariants and representations. Princeton University Press.

68