Latest posts
Going from the binomial to the Poisson distribution
Domagoj Pandza on April 14th, 2018
Mathematics
Show in single article view

In this article, we are going to derive the Poisson distribution from the binomial distribution. If you don't know what the binomial distribution is or how to derive it from "scratch", I wrote an introduction to it in a previous article. As for the Poisson distribution, it is, obviously, also a discrete probability distribution whose shtick is quantifying the likelihood of $k$ events happening in some interval provided you know the average number of such events happening in said interval and that each event is independent of the previous ones. In other words, their timing doesn't affect the timing of the subsequent ones.

One dark example of where the timing of a single event could be considered to affect the timing of future events are passenger aircraft crashes, even though it is often given as an example problem when discussing the Poisson distribution. However, in reality, such events cause global repercussions that affect airlines, plane manufacturers and regulatory bodies in considerable and unpredictable ways. Fleets of aircrafts get grounded, corrections get issued and the issue gets resolved, all of which affects the likelihood of future crashes. If one is pressed for an aviation related problem, a better example that actually obeys the Poisson distribution would be the number of engine shutdowns on an aircraft over a large number of flight hours or cycles.

The boilerplate example is estimating the likelihood of a certain number of floods given statistics accumulated over, say, 10 years. Deviating from the mean on either side will, obviously, adversely affect the likelihood of that number of floods happening during a period of 10 years. All we have to do now is find a way to quantify this concept.

Sending $n$ to $\infty$

It's the right thing to do

As mentioned, we start with the binomial distribution. It tells us the likelihood of having $k$ successes out of $n$ trials. The Poisson distribution tells us the likelihood of $k$ events happening given the average of such events happening during some interval (which is implicit). We know the expected value or the mean of the binomial distribution is equal to $nq$ where n is the number of trials and $q$ is the probability of success of a single trial. The expected value makes sense; the more trials we conduct and the higher the probability of success, the more successes we ought to expect. Thus, another way of expressing the binomial distribution by aliasing $\mu$ as $\lambda$ is:

$$ p_q(k, n) = {n \choose k} q^k(1-q)^{n-k} = {n \choose k} {\left(\frac{\lambda}{n}\right)}^k\left(1-\frac{\lambda}{n}\right)^{n-k} $$

We don't know what the probability $q$ of any one single trial is, but we know an "overall trend" $\lambda$ which we effectively take to be determined over an $\infty$ number of trials, hence to express the Poisson distribution, we allow the number of trials of $n$ go to infinity:

$$ p(k) = \lim_{n\to\infty} p_q(k, n) = \lim_{n\to\infty} {n \choose k} {\left(\frac{\lambda}{n}\right)}^k\left(1-\frac{\lambda}{n}\right)^{n-k} = \lim_{n\to\infty} \frac{n!}{(n-k)!k!} {\left(\frac{\lambda}{n}\right)}^k\left(1-\frac{\lambda}{n}\right)^{n-k} $$

Taking the limit here is not difficult, as we'll show in detail. What creates difficulty for people is the fact that $\lambda = nq$ instigates an indeterminate form if substituted in. However, we make the legitimate assumption that as $n$ goes to infinity, $q$ tends to zero. In other words, we assume $\lambda$ remains constant.

We start off slow, by moving the constant factors outside:

$$ p(k) = \frac{\lambda^k}{k!} \lim_{n\to\infty} \frac{n!}{(n-k)!} {\left(\frac{1}{n}\right)}^k\left(1-\frac{\lambda}{n}\right)^{n-k} $$

Notice that the equation can be rewritten as follows:

$$ p(k) = \frac{\lambda^k}{k!} \lim_{n\to\infty} \left(1-\frac{\lambda}{n}\right)^{n-k} \prod_{i=0}^k { \frac{n - i}{n} } = \frac{\lambda^k}{k!} \lim_{n\to\infty} \left(1-\frac{\lambda}{n}\right)^{n-k} \prod_{i=0}^{k-1} {\left( 1 - \frac{i}{n} \right)}$$

If you're having trouble seeing where that product came from, here it is expanded:

$$ \frac{n!}{(n-k)!} \frac{1}{n^k} = \frac{n}{n} \cdot \frac{n - 1}{n} \cdot \frac{n - 2}{n} \cdots \frac{n - k + 1}{n} = \prod_{i=0}^{k-1} { \frac{n - i}{n} } $$

It is easy to see that the product approaches $1$ when we allow $n$ to tend to infinity, hence:

$$ \lim_{n\to\infty} \prod_{i=0}^{k-1} { \left( 1 - \frac{i}{n} \right) } = 1 \implies p(k) = \frac{\lambda^k}{k!} \lim_{n\to\infty} { \left(1-\frac{\lambda}{n}\right)^{n-k} }$$

We rewrite the expression further chasing simplification:

$$ p(k) = \frac{\lambda^k}{k!} \lim_{n\to\infty} \left(1-\frac{\lambda}{n}\right)^{n} \left(1-\frac{\lambda}{n}\right)^{-k}$$

It's immediately obvious that:

$$ \lim_{n\to\infty} \left(1-\frac{\lambda}{n}\right)^{-k} = 1 \implies p(k) = \frac{\lambda^k}{k!} \lim_{n\to\infty} \left(1-\frac{\lambda}{n}\right)^{n}$$

The remaining factor might look familiar to the way we define $e$ (taking the limit of $(1 + 1/x)^x$ as $x$ goes to infinity). It's not quite there yet, but it stands to reason we can transform it in a way that it can be expressed as a power of $e$:

$$ \lim_{n\to\infty} \left(1-\frac{\lambda}{n}\right)^{n} = \lim_{m\to\infty} \left(\left(1+\frac{1}{m}\right)^{m}\right)^{-\lambda} = e^{-\lambda} \iff m = -\frac{n}{\lambda} $$

Which resolves our distribution to its final form:

$$\Large { p(k) = \frac{\lambda^k}{k!} e^{-\lambda} }$$

While it might not surprise you that the the expected value of this probability distribution is $\lambda$, its variance is also $\lambda$. The moment $e$ crops up in your model, its properties are bound to be interesting.

Deriving the binomial distribution
Domagoj Pandza on April 8th, 2018
Mathematics
Show in single article view

The motivation for writing this article is twofold. First is completeness, as the binomial distribution is an intuitive starting point for elegantly deriving the Poisson distribution which is a topic of a subsequent article, itself a key ingredient of another article. Second is showing how one can intuitively arrive at the binomial distribution as I run into people, from time to time, who have trouble understanding how exactly it works. While its very name might be intimidating to a beginner, it is based on a very simple and intuitive concept, as you'll see in a moment.

At its core, the binomial distribution is a discrete probability distribution which quantifies the likelihood that given $n$ trials, each with a binary outcome, $k$ of them will be a "success". The simplest example would be a series of coin flips, where the sides can be assigned to be either "success" or "failure".

Heads or tails?

Let's throw dice instead

Coins, cards, and dice are the usual suspects when it comes to illustrating concepts in probability theory as they're useful for building intuition. We’ll go with rolling dice just to mix it up a little. Now, here's a question: “If we roll three dice, one by one, what are the odds of getting two sixes?”

I'd like to use coins, but fair coin flips have an associated probability of $0.5$ which makes it equivalent to the probability of failure and thus might make it difficult for people to make some important observations, without introducing magic bias for illustration purposes. So, we stick to fair dice rolls by declaring a roll of six to be "success", "yes", $1$ or whatever you want to call it. Fair dice roll means that each face of the dice has a $1/6$ chance of coming up on top. In our case, the probability of success is thus $1/6$ while the probability of failure is a whopping $1 - 1/6 = 5/6$.

Given the probability of success is $1/6$ for each independent roll, only $1$ out of $6$ times is going to result in rolling a six. On that one out of six times, there's only a $1/6$ chance we roll another six, reducing the odds to $1$ out of $36$.

Why exactly do we multiply probabilities?

In my experience, this is where some folks get lost. They know the rule is to multiply probabilities in the case of independent events, but they don't know why. So, I will endeavor to elaborate further. There are six ways to roll the first dice. And six ways to roll the second one. For each of the six ways to roll the first dice, there are six ways we can roll the second one as they're independent of each other, one roll does not affect the outcome of the other. This brings us to a total of $36$ possibilities, each equally likely, yielding a probability of $1/36$.

On that $1$ out of $36$ times, there is a considerable chance we will not roll a six as the probability of failure is $5/6$. Thus, our chances decrease slightly to a final $5$ times out of $216$, giving us a roughly $2.3\%$ chance to roll two sixes followed by anything else. Two successes, one failure. Effectively, we computed the probability of a single arrangement of two successes and one failure with $p_s = 1/6$ as:

$$ p_sp_s(1 - p_s) = p_s^2(1 - p_s)^{3-2} = \left(\frac{1}{6}\right)^2 \frac{5}{6} = \frac{5}{216}$$

Notice there are three different ways of getting two successes and a failure, in this case, $\{SSF, SFS, FSS\}$. Our computation only deals with the probability of getting any one of those three possible orderings. Thus the probability of getting two successes and one failure is actually three times higher. To generalize this, let $n$ be the number of trials and $k$ be the number of successes, leaving failures implicitly quantified as $n-k$. We need to know how many ways we can order $k$ successes or $n-k$ failures in $n$ trials. Either will do, so we choose what is explicit, the number of successes. If there are $n$ trials, there is $n!$ ways they can be ordered. When you select the "temporal order" of the first trial, there remain only $n -1$ choices for the second placement and so forth. We're only choosing the placement of $k$ trials, the successes. Hence we need to stop the $n!$ factorial after $k$ "iterations", this is done by dividing it by $(n-k)!$ (write it out if you don't see it). Now we have all the ways we can order $k$ successes in the $n$ trial slots. We're only interested in the unique placements of successes (and thus failures which are implicitly in-between) and given that there are only binary outcomes, some of these unique conceptual orderings of trials will end up aliasing one another once their values are "substituted in", we need to remove those.

Let me show you what I mean. In our trivial example with $n=3$ and $k=2$, there is $3! = 6$ ways to order the trials. There is $3! / (3 - 2)! = 6$ ways to order $2$ successes in the $3$ trials. I will use numbers to indicate temporal ordering of the possible configurations: $\{12, 21, 13, 31, 23, 32 \}$. Notice that ordering successes as $12$ and $21$, $13$ and $31$, as well as $23$ and $32$ results in redundancy. Hence, there is only $3$ unique ways we need to account for. There is $k!$ ways to order $k$ elements, hence we finally divide by $k!$. This brings us to the final equation where $q$ is the probability of success, $k$ is the number of successes and $n$ is the total number of trials:

$$ p_q(k, n) = \frac{n!}{(n-k)!k!} q^k(1-q)^{n-k} = {n \choose k} q^k(1-q)^{n-k} $$

You might have recognized "n-choose-k" from combinatorics, ways of choosing $k$ successes out of $n$ trials which, while it matches the value, doesn't quite relay the semantics. And if you're wondering about the name, originally, it was Bernoulli who proved that $p_q(k,n)$ we derived in this article is equal to the value of the $k$th term of the polynomial expansion of the binomial power $(q + (1-q))^n$. Careful, don't collapse it. Obviously, it needs to resolve to $1$, we're talking about probabilities.

The distribution itself is inherently discrete which makes defining the cumulative distribution function rather trivial. As the cumulative distribution function tells us the probability that a random variable will end up being less than or equal to $k$, it merely amounts to adding up discrete probabilities up to $k$:

$$ P_q(k, n) = \sum_{i=0}^{k} { {n \choose i} q^i(1-q)^{n-i} } $$

What value to expect?

A quick derivation

Like all distributions, we need to discuss the expected value of a random variable drawn from it. In this case, it simply comes down to the average of all the possible values weighted by their probability. You already know how the average is defined from childhood, but if you're wondering why it is defined the way it is, I wrote an article on the subject which gives it an intuitive treatment. As expected, we start with:

$$ \mu = \frac{ \sum_{i=0}^{n} { x_i w_i } } {\sum_{i=0}^{n} { w_i }} = \sum_{i=0}^{n} { x_i p(x_i) } $$

Substituting the binomial distribution function for $p(x)$ and $k$ for $x_i$ which can neatly alias the iterator $i$:

$$ \mu = \sum_{k=0}^{n} { k \space {n \choose k} q^k(1-q)^{n-k} } = \sum_{k=0}^{n} { k \space \frac{n!}{(n-k)!k!} q^k(1-q)^{n-k} } $$

From here, we need to simplify. If we could get rid of the leftmost $k$ in a manner which would make the binomial theorem applicable, we'd be able to completely obliterate the sum. To this end, I'd like to point out the following holds:

$$ {n \choose k} = \frac{n}{k} {n - 1 \choose k - 1} = \frac{n}{k} \frac{(n-1)!}{(k-1)!((n-1) - (k-1))!} \space \iff \space k {n \choose k} = n {n - 1 \choose k - 1} $$

Given this relationship, we can transform the $k$ out of the sum. Furthermore, if it isn't immediately obvious, $(n-1) - (k-1) = n - k$. Here's the result after our first round of transformations:

$$ \mu = \sum_{k=0}^{n} { k \space {n \choose k} q^k(1-q)^{n-k} } = n\sum_{k=1}^{n} { {n - 1 \choose k - 1} q^k(1-q)^{(n-1) - (k-1)} } $$

As $k - 1$ needs to be larger than or equal to $0$, we set the lower limit of the sum to $1$. To make the binomial theorem applicable, we finalize by extracting a single factor of $q$. Then, we just transform the range to make it more readable by using $j = k - 1$ and $m = n - 1$:

$$\begin{split} \mu &= n\sum_{k=1}^{n} { {n - 1 \choose k - 1} q^k(1-q)^{(n-1) - (k-1)} } \\ & = nq\sum_{k=1}^{n} { {n - 1 \choose k - 1} q^{k-1}(1-q)^{(n-1) - (k-1)} } \\ & = nq\sum_{j=0}^{m} { {m \choose j} q^{j}(1-q)^{m - j} } \\ \end{split}$$

And now we just apply it:

$$ \mu = nq\sum_{j=0}^{m} { {m \choose j} q^{j}(1-q)^{m - j} } = nq(q + (1-q))^m = nq $$

This is my preferred way of arriving at the mean, other approaches feel either circular (depending on the knowledge of the mean in the first place) or are just convoluted. The same approach works for the variance, so I will just tell you it's equal to $nq(1-q)$ should you decide to work it out yourself. Remember, variance is defined as the expected value of the squared deviation of a random variable drawn from the distribution to the mean of the distribution. Simplified, it comes down to:

$$ E[(X - E[X])^2] = E[X^2 - 2XE[X] + E[X]^2] = E[X^2] - 2E[X]^2 + E[X]^2 = E[X^2] - E[X]^2$$

You already know the squared expected value as it is the mean squared, $(nq)^2$, all you have to do is re-do the dance we did before to calculate the expected value of the squared random variable, meaning you calculate the sum of squares of the possible values weighted by the associated probabilities. Enjoy.

Arriving at quaternions from "scratch"
Domagoj Pandza on March 21st, 2018
Mathematics
Show in single article view

A lot of people use quaternions to represent object orientations and to perform rotations in a dependable way, without actually understanding how they work. This is peculiar as the discovery of quaternions is often cited as perhaps the most documented discovery in the history of mathematics. Unfortunately, Hamilton gave the world a poor first impression because he had trouble understanding quaternions. A decade-long campaign of misinterpretation and confusion later, quaternions $\mathbb{H}$ were whisked away to the plane of oblivion for more than a century, at least as far the mainstream public was concerned. But not before its carcass was picked clean by the likes of Gibbs and Heaviside who independently brought about vector algebra and calculus in the modern sense, among other things.

This article aims to show how quaternions can be derived and why they manage rotations in such an efficient, reliable and elegant manner. It also sets the stage for follow up articles on geometric algebra which will more clearly illustrate the useful properties quaternions have (as they're a subalgebra of it). Not to mention open your eyes to see just how badly one man's confusion in 1843. gave mathematical physics a limp we're yet to get rid of.

A quick rundown of \(\mathbb{C}\)

Complex numbers are the algebraic basis for arriving at quaternions, as their properties motivate the development of quaternions. As you should already know, complex numbers arise when one desires to solve the canonical \(x^2 + 1 = 0\) by defining an indeterminate \(i\) as the solution. Complex numbers have a real and imaginary part, both of which are real numbers. The imaginary part enumerates the indeterminate and when combined with the real part, it forms an algebraic extension of real numbers and enables solutions of polynomials such as the one that motivated the extension. We can write a complex number as:

$$ z = x1 + yi $$

They may look a bit peculiar to people, but they're as "legit" as real numbers are. Here's a simple example, we will square the expression above and procure yet another complex number:

$$ z^2 = (x + yi)^2 = x^2 + 2xyi + y^2i^2 = x^2 - y^2 + 2xyi$$

Much like the set of real numbers \(\mathbb{R}\), \(\mathbb{C}\) is a field which, among other terms, predicates being closed under multiplication. Just like one can plot real numbers on the number line, complex numbers may be plotted in the complex plane which can be visualized as two real number lines setup perpendicular to one another. One plots the real part on the horizontal "axis" and the imaginary part on the vertical one.

The complex plane is not a complex vector space.

Some people erroneously conflate complex vector spaces and the complex plane. \(1\) and \(i\) are not to be thought of as linearly independent bases. As mentioned, the set of complex numbers is a field and can merrily replace real numbers as scalars of a vector space, in which case it becomes a complex vector space. Complex numbers themselves are not vectors. Complex numbers can be constructed as ordered pairs with special operations defined on them. But an ordered pair is not a vector. Historically, vector spaces came after quaternions which came after complex numbers.

Like real numbers, complex numbers define absolute value to coincide with the distance from the origin. Unlike real numbers, though, we have to avail ourselves of some geometric insights to compute the absolute value or modulus. Thus, the absolute value of \(z\) is equal to \(\sqrt{x^2 + y^2}\). While we're here, we may as well introduce the complex conjugate by showing how one can compute the absolute value with its help:

$$ \sqrt{x^2 + y^2} = \sqrt{x^2 - y^2i^2 + xyi - xyi} = \sqrt{(x + yi)(x - yi)} = \sqrt{zz^*} $$

In the complex plane, all complex numbers can be thought of as directed segments from the origin thanks to the associated absolute value / modulus. The behavior or basic addition and subtraction may be reminiscent of what you know of vector algebra, which should give you a clear historical perspective where that came from. However, complex numbers also have interesting geometric interpretations when it comes to multiplying and dividing them. To best illustrate this without actually needing to whip out graphs, we will transform to polar form and further compact it by using Euler's relation:

$$ z = r(\cos\theta + i\sin\theta) = re^{i\theta} $$

To transform to polar form, \(r\) is just the absolute value of \(z\), \(|z|\), the distance to the origin. You recover \(\theta\) with \(\tan^{-1}(y/(r+x))\), which is the angle to the horizontal "axis". Now, I want you to observe what happens when you multiply two complex numbers:

$$ z_1z_2 = r_1e^{i\alpha} r_2e^{i\beta} = r_1r_2e^{i(\alpha + \beta)} $$

Interestingly, the angles from the horizontal axis are being added together while the radii are multiplied together. As we already mentioned, the radius is the absolute value or modulus, the distance to the origin. The fact that the modulus of the product is equal to the product of the two moduli was called the law of the moduli by Hamilton and was one of the requirements when he was searching for quaternions. This rather trivial-looking law is of importance, as we can now talk of a complex number $R$ with $r=1$ which has the ability to rotate $z$ by an angle $\theta$:

$$ Rz = zR = e^{i\theta} re^{i\alpha} = re^{i(\alpha + \theta)} = z' $$

This was a profound insight, one sorely known to Hamilton. Keep in mind the commutativity of the complex product, as it is a notable feature courtesy of the commutativity of the addition of real numbers. This result may also encourage one to think about the division of two directed segments:

$$ \frac{z_1}{z_2} = \frac{r_1}{r_2}e^{i(\alpha - \beta)} $$

When two directed segments are parallel (the angle is $0$), the result is the positive ratio of their lengths (and thereby norms). When they're anti-parallel (angle of $\pi$), the exponential degenerates to $e^{i\pi} = \cos\pi + i\sin\pi = -1$ which intuitively signs the ratio of the lengths to be negative. The last option is that they have some general angle of separation, at which point division encodes both angular and length information. In case of unit directed segments, their division can be used to recover the angle between them. Hamilton desperately wanted to generalize this to the third dimension, but given that you have never heard of vector division in this sense in classic vector algebra, you may guess on his success. To be honest, the last two concepts discussed, rotating by a single product and the idea of the quotient of two directed segments, brought the entire ship down through a case of chronic underestimation and shortsightedness. The dogmatic perserverence that the same must hold in three dimensions instead of assuming that they're just special cases of a much more general construct. But more on that later.

Going from \(\mathbb{C}\) to \(\mathbb{H}\)

While the geometric interpretation was never far from Hamilton's mind, the way he arrived at quaternions was purely algebraic. Which incidentally may be one of the problems, because Hamilton and his followers tried to shoehorn the geometric interpretation of complex numbers after the fact without giving proper weight to the whole "adding a new dimension" thing. So, we retread the steps, as promised, but for brevity, we don't make all the wrong turns as part of our derivation.

Hamilton did the obvious expansion of complex theory, he added yet another imaginary part and indeterminate $j$ which squared to $-1$ and was presumed perpendicular to both the real and imaginary axis (when interpreting geometrically). Trouble catches up to us once we try to multiply two generic triplets in this new algebra we're crafting and attempt to see whether it is closed under multiplication like \(\mathbb{C}\) is, but you might guess where this is going:

$$\begin{split} q_1q_2 & = (x_1 + y_1i + z_1j)(x_2 + y_2i + z_2j) \\ & = x_1x_2 - y_1y_2 - z_1z_2 \\ & + (x_1y_2 + x_2y_1)i \\ & + (x_1z_2 + x_2z_1)j \\ & + (y_1z_2 - y_2z_1)ij \\ \end{split}$$

As expected, a higher grade term emerges and is not going away. To close the algebra, we make the safe assumption $ij = k$ and add $k$ into the mix, which Hamilton himself was dumbfounded by, describing it as a paradox:

$$\begin{split} q_1q_2 & = (x_1 + y_1i + z_1j + w_1k)(x_2 + y_2i + z_2j + w_2k) \\ & = x_1x_2 - y_1y_2 - z_1z_2 - w_1w_2 \\ & + (x_1y_2 + x_2y_1)i + (x_1z_2 + x_2z_1)j + (x_1w_2 + x_2w_1)k \\ & + (y_1z_2 - y_2z_1)ij + (y_1w_2 - y_2w_1)ik + (z_1w_2 - z_2w_1)jk \\ \end{split}$$

Given that we established $ij = -ji$, key being the anticommutative property, we assert $ik = -ki$ and $jk = -kj$ as well. We also hold that $i^2 = j^2 = k^2 = -1$ as we're trying to generalize complex numbers in some way. $ij = k$ imposes $ik = iij = -j$, $jk = -ikk = i$ and $ijk = kk = -1$:

$$\begin{split} q_1q_2 & = (x_1 + y_1i + z_1j + w_1k)(x_2 + y_2i + z_2j + w_2k) \\ & = x_1x_2 - y_1y_2 - z_1z_2 - w_1w_2 \\ & + (x_1y_2 + x_2y_1 + z_1w_2 - z_2w_1)i \\ & + (x_1z_2 + x_2z_1 + y_2w_1 - y_1w_2)j \\ & + (x_1w_2 + x_2w_1 + y_1z_2 - y_2z_1)k \\ \end{split}$$

All of this is a pretty big deal, as in the case of triplets there is nothing you can do to close the algebra or, in other words, make $ij$ go away. Once we uncover the exponential form in the next few sections, we will show that the law of moduli holds, along with some information on why it matters. Though, nothing is preventing you from bruteforcing it right now. In case you forgot what the law of moduli was, it's $|q_1||q_2| = |q_1q_2|$. Simple, but extremely potent.

Keeping up with the Hamiltons

Mankind divided over vectors divided

Hamilton knew what he wrapped his hands around the day it dawned on him to shift to the fourth dimension. Historical documents suggest he immediately started searching for the geometric interpretation of his discovery, getting basic rotations going by the end of the first day. There were three main notions from complex numbers that he wanted to generalize. The ability of complex numbers to represent directed segments, the ability of complex numbers to rotate those segments via a single complex product and the notion that one could recover useful information akin to the dot product through the concept of vector division. In fact, Hamilton was so obsessed with vector division that he actually defined the quaternion in terms of it: A quaternion is the quotient of two vectors, in the general case.

Hamilton introduced the terms vector and scalar

As you might know, in the first half of the 19th century, there was no formal notion of vectors. There were many informal approaches which were used to obtain useful results, concepts of directed segments and arrows can be observed in works from the 17th century, such as in the works of Sir Isaac Newton. But these constructs were not called vectors.

Hamilton based the term vector on the Latin verb vehere meaning "to carry". Similarly, he was the one who came up with the term scalar. The name vector was introduced to augment Hamilton's developing quaternion theory, where he considered the imaginary part to form a directed segment (now referred to as a vector), whereas the real part was considered the scalar. Furthermore, in Hamilton's quaternion algebra, a real number could be expressed as a real quaternion with a zero imaginary vector $\mathbf{0}$, $(a, \mathbf{0})$. Contrasting this is a pure quaternion which consists of a non-existent scalar part and the imaginary vector $(0, \mathbf{v})$.

That obviously isn't the case today, so what happened? A geometric interpretation. In common vernacular, everything is subject to interpretation and there are many views that humans accept as valid even though they're most certainly expressing a subset of the whole truth. Due to Hamilton's heavy bias towards the geometric interpretation of complex numbers and making strong, unsubstantiated assumptions about the behavior of quaternions, he was wrong on many levels. This is a bit disappointing as a physicist should know that many things which are true in two dimensions are outright false or false in the general case in three dimensions or more. One example is the relationship between vector quantities of angular momentum and angular velocity in modern physics. In two dimensions, they share direction. In three dimensions, they can be offset by a considerable angle giving rise to what is sometimes referred to as gyro motion.

Geometric interperation of \(\mathbb{H}\)

Back on track!

The crux of the problem here is that Hamilton failed to recognize the two-dimensional interpretation of complex numbers as the special case, not the norm. Initial surveying of the geometric interpretation of quaternions was done manually, by trying to map out where the result gets placed in three dimensional space (by treating the imaginary part as a directed segment or vector, as positions/directions versus rotations). We will do something a bit different. We recognize how obvious the geometric interpretation of complex numbers becomes in exponential form, so that might be a good start for us. There are multiple ways of arriving at an exponential form, by writing out the series or alternatively going back to Euler's relation with a few tricks up our sleeve. First, verify for yourself that the following statements are all true:

$$|q|^2 = |s|^2 + |\mathbf{v}|^2$$

$$q = s + \mathbf{v} = |q|\Bigg( \frac{|s|}{|q|} + \frac{|\mathbf{v}|}{|q|} \frac{\mathbf{v}}{|\mathbf{v}|} \Bigg)$$

$$ \begin{split} \mathbf{\hat{v}}^2 &= (v_xi + v_yj + v_zk)(v_xi + v_yj + v_zk) \\ &= \mathbf{\hat{v}} \times \mathbf{\hat{v}} - \mathbf{\hat{v}} \cdot \mathbf{\hat{v}} = \mathbf{0} - 1 = -1 \end{split} $$

Transform the modulus expression by dividing both sides by $|q|^2$:

$$1 = \Bigg( \frac{|s|}{|q|} \Bigg) ^2 + \Bigg( \frac{|\mathbf{v}|}{|q|} \Bigg) ^2 $$

Apply the trigonometric identity $\cos^2\theta + \sin^2\theta = 1$ by assuming there exists $\theta$ which satisfies:

$$ \frac{|s|}{|q|} = \cos\theta \qquad \frac{|\mathbf{v}|}{|q|} = \sin\theta $$

As \(\mathbf{\hat{v}}^2 = -1\), the Euler relation holds here as well:

$$q = s + \mathbf{v} = |q|( \cos\theta + \mathbf{\hat{v}} \sin\theta) = |q|e^{\mathbf{\hat{v}}\theta} $$

Satisfying the law of moduli

With the exponential form in hand, it becomes easy to show $|q_1||q_2| = |q_1q_2|$:

$$|q_1||q_2| = \sqrt{q_1q_1^*q_2q_2^*} = \sqrt{|q_1|^2|q_2|^2e^{\mathbf{\hat{v}}_1(\alpha - \alpha) + \mathbf{\hat{v}}_2(\beta - \beta)}} = \sqrt{q_1q_2(q_1q_2)^*} = |q_1q_2|$$

Only 4 normed division algebras satisfy this rule (and thus exist), $\mathbb{R}$, $\mathbb{C}$, $\mathbb{H}$ and an algebra that follows naturally from quaternions, octonions $\mathbb{O}$.

It is immediately visible that a compact exponential form will not trivially reveal how exactly we can achieve rotations through multiplication due to the added dimensionality. Thus we will not make use of Euler's relation here and just proceed in the expanded form. We will attempt, much like Hamilton, to multiply a quaternion representing the rotation \(r\) with a quaternion \(p\) representing the position we're interested in rotating. Hopefully, by expanding, we will determine which part of the rotation quaternion defines the plane of rotation and whether or not an actual proper rotation takes place.

$$rp = |p|( \cos\theta + \mathbf{\hat{r}}_v \sin\theta)( \cos\alpha + \mathbf{\hat{p}}_v \sin\alpha) $$

The position quaternion has no scalar component at all, which in this form needs to be eliminated by setting \(\alpha = \pi/2\). And then we just:

$$ \begin{split} rp &= |p|( \cos\theta + \mathbf{\hat{r}}_v \sin\theta)\mathbf{\hat{p}}_v \\ &= |p|( \mathbf{\hat{p}}_v\cos\theta + \mathbf{\hat{r}}_v\mathbf{\hat{p}}_v \sin\theta) \\ &= |p|( \mathbf{\hat{p}}_v\cos\theta + (\mathbf{\hat{r}}_v \times \mathbf{\hat{p}}_v - \mathbf{\hat{r}}_v \cdot \mathbf{\hat{p}}_v) \sin\theta) \end{split} $$

It follows that the only way to have a rotated pure quaternion result is iff \( \mathbf{\hat{r}}_v \cdot \mathbf{\hat{p}}_v = 0 \), which is to say that the point needs to be in the rotation plane. This is not very useful. Hamilton was very excited about this, though. He thought he had finally generalized the geometric interpretation of complex numbers to three dimensions with quaternions and some help from the "paradoxical" fourth dimension. This little escapade of ours hasn't been a total bust, though. We now know for sure that \(\mathbf{\hat{r}}_v = \mathbf{\hat{n}}\), it is the rotation plane's normal. Our next step is to explore this further. It is reasonable to expect arbitrary rotations are possible and that the result will be a pure quaternion with its norm preserved. To this end, we must find a way to extrapolate the scalar part from a quaternion. Since we don't have much to go on, removing the scalar part is a good start. This is fairly easy if you remember your conjugation:

$$q + q^* = 2w = 2S(q) \iff S(q) = \frac{q + q^*}{2}$$

Seems rather trivial, but it's useful nonetheless, as we can posit requirements such as that whatever the missing piece of our product, we desire $S(\hat{r}p\hat{q}) = 0$. We also know that proper rotations preserve norm, so $|\hat{r}p\hat{q}| = |p|$. We can constrain $\hat{q}$ by:

$$ S(\hat{r}p) = S(p\hat{q}) \iff S(\hat{r}p\hat{q}) = S(\hat{r}p) - S(p\hat{q}) = 0 \iff \hat{q} = \hat{r}^*$$

Conjugates negate the imaginary vector part leaving the scalar part untouched. That's how we formulated the scalar extraction function in the first place. Finally, let's check the norm requirement knowing that we're working with unit quaternions and that $|q| = |q^*|$:

$$|\hat{r}p\hat{r}^*| = |\hat{r}||p||\hat{r}^*| = |p|$$

This allows us to propose the following transform:

$$ \Large { p' = \hat{r}p\hat{r}^* = \hat{r}p\hat{r}^{-1} }$$

Right off the bat, one can guess there's a double transform by \(\theta\) going on here. However, to get the details, we have to work it out:

$$ \begin{split} \mathbf{p'} = \hat{r}p\hat{r}^* &= |p|[\cos\theta + \mathbf{\hat{n}}\sin\theta]\mathbf{\hat{p}_v}[\cos\theta - \mathbf{\hat{n}}\sin\theta] \\ &= |p|[ \mathbf{\hat{p}_v}\cos\theta + \mathbf{\hat{n}}\mathbf{\hat{p}_v} \sin\theta] [\cos\theta - \mathbf{\hat{n}}\sin\theta] \\ &= |p|[ \mathbf{\hat{p}_v}\cos^2\theta + (\mathbf{\hat{n}}\mathbf{\hat{p}_v} - \mathbf{\hat{p}_v}\mathbf{\hat{n}})\sin\theta\cos\theta - \mathbf{\hat{n}}\mathbf{\hat{p}_v}\mathbf{\hat{n}} \sin^2\theta] \\ &= |p|[ \mathbf{\hat{p}_v}\cos^2\theta + (\mathbf{\hat{n}} \times \mathbf{\hat{p}_v} - \mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v} + \mathbf{\hat{n}} \times \mathbf{\hat{p}_v} + \mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\sin\theta\cos\theta - (\mathbf{\hat{n}} \times \mathbf{\hat{p}_v} - \mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}} \sin^2\theta] \\ &= |p|[ \mathbf{\hat{p}_v}\cos^2\theta + 2\mathbf{\hat{n}} \times \mathbf{\hat{p}_v}\sin\theta\cos\theta - (\mathbf{\hat{n}} \times \mathbf{\hat{p}_v} \times \mathbf{\hat{n}} - (\mathbf{\hat{n}} \times \mathbf{\hat{p}_v}) \cdot \mathbf{\hat{n}} - (\mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}})\sin^2\theta] \\ &= |p|[ \mathbf{\hat{p}_v}\cos^2\theta + 2\mathbf{\hat{n}} \times \mathbf{\hat{p}_v}\sin\theta\cos\theta - (\mathbf{\hat{p}_v} (\mathbf{\hat{n}} \cdot \mathbf{\hat{n}}) - 2(\mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}})\sin^2\theta] \\ &= |p|[ \mathbf{\hat{p}_v}(\cos^2\theta - \sin^2\theta) + 2\mathbf{\hat{n}} \times \mathbf{\hat{p}_v}\sin\theta\cos\theta + 2(\mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}}\sin^2\theta] \\ \end{split}$$

At this stage, we can confirm that the scalar part has been completely eliminated, but we still need to simplify the expression to understand the semantics of quaternion rotation. The following trigonometric identities should help us do so, which you may already know:

$$ \cos(2\theta) = \cos^2\theta - \sin^2\theta \qquad \sin(2\theta) = 2\sin\theta\cos\theta \qquad \sin^2(\theta) = \frac{1 - \cos(2\theta)}{2} $$

And we can finally unravel the product:

$$ \begin{split} \mathbf{p'} = \hat{r}p\hat{r}^* &= |p|[\cos\theta + \mathbf{\hat{n}}\sin\theta]\mathbf{\hat{p}_v}[\cos\theta - \mathbf{\hat{n}}\sin\theta] \\ &= |p|[ \mathbf{\hat{p}_v}(\cos^2\theta - \sin^2\theta) + 2\mathbf{\hat{n}} \times \mathbf{\hat{p}_v}\sin\theta\cos\theta + 2(\mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}}\sin^2\theta] \\ &= |p|[ (\mathbf{\hat{p}_v} - (\mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}})\cos(2\theta) + \mathbf{\hat{n}} \times \mathbf{\hat{p}_v}\sin(2\theta) + (\mathbf{\hat{n}} \cdot \mathbf{\hat{p}_v})\mathbf{\hat{n}}] \\ &= (\mathbf{p_v} - (\mathbf{\hat{n}} \cdot \mathbf{p_v})\mathbf{\hat{n}})\cos(2\theta) + \mathbf{\hat{n}} \times \mathbf{p_v}\sin(2\theta) + (\mathbf{\hat{n}} \cdot \mathbf{p_v})\mathbf{\hat{n}} \end{split}$$

If you're having trouble seeing what's going, maybe the following substitutions can help:

$$ \mathbf{p'} = \hat{r}p\hat{r}^* = (\mathbf{p_v} - \mathbf{p_{\parallel}})\cos(2\theta) + \mathbf{p_{\perp}}\sin(2\theta) + \mathbf{p_{\parallel}} $$

Disregarding the factor doubling the angle, this is the exact same conical transformation you arrive at geometrically in linear algebra. In fact, this was its first appearance. Double angle is easily mitigated by supplying the rotation quaternion with half the intended angle. The beauty of quaternions is that they're extremely compact and implicit with respect to other representations, they're tolerant of precision loss and you can trivially recover the plane of rotation and the angle involved. But most of all, quaternion rotation doesn't suffer from what people annoyingly keep referring as the gimbal lock. We're talking about potential singularities caused by mapping deficiencies of the chosen representation.

Relationship with rotation group $SO(3)$

I will attempt to keep this as simple as possible, which will most definitely come at the expense of accuracy. I recognize that a lot of people don't have a background in abstract algebra and topology, and quite frankly, they don't even want one. In the future, I might make a series taking a deeper dive into this, but for now, let's just show what makes quaternions so well suited for representing rotations.

The three-dimensional rotation group $SO(3)$ is the group of all rotations around the origin of \(\mathbb{R}^3\) under composition. In other words, composing rotations will result in yet another rotation. The order in which you compose rotations matters, hence the group is non-abelian. The group of all rotations under composition sounds very abstract, because it is. It doesn't really describe how you can make a rotation, that's the job of a representation which conceptually maps to the group of all rotations with varying degrees of success. Representations always come first in terms of "discovery". The most common representation this rotation group is identified with are \(3\times3\) orthonormal rotation matrices which form a group under matrix multiplication and encode rotations with a discernible rotation axis and angle. But the one that helped formalize the various properties of the group was actually the discovery of the quaternion group.

To illustrate this in a lower dimension, we discussed the geometric interpretation of the multiplication of complex numbers at the start of this article and how complex numbers of absolute value \(1\) can represent rotations, encoding the angle as shown in the compact polar form \(e^{i\theta} = \cos\theta + i\sin\theta\). Strictly speaking, complex numbers of absolute value $1$ form a group under multiplication called \(\mathbf{U}(1)\). Clearly the set of all such complex numbers forms a unit circle in the complex plane, which we tag as \(S^1\). Hence, it stands to reason it can have an isomorphism (an invertible mapping) to the group of all \(2\times2\) orthonormal rotation matrices which can be identified with \(SO(2)\):

$$ \Large { e^{i\theta} \Rightarrow \begin{pmatrix} \cos\theta & -\sin\theta \\ \sin\theta & \cos\theta \end{pmatrix} } $$

Thus, we prove the isomorphism. For three dimensions, things are not so simple. As you saw when we derived the geometric interpretation of the group of unit quaternions under multiplication. The rotation axis is no longer implicit. Quaternions, and the time you surely spent constructing \(3\times3\) orthonormal rotation matrices, showed us that we need both a rotation axis and an angle. And speaking of unit quaternions, what kind of topology are we dealing with? Well, if you take a look at how their norm is computed, all the unit quaternions define a sphere in four dimensional space, a 3-sphere $S^3$ which happens to be three dimensional. This is to facilitate both encoding the rotation axis and angle as we derived in previous sections.

Following are some observations about rotations by an angle around an axis. A rotation around axis \(\mathbf{n}\) with angle \(\theta\) is the same rotation as the one around \(-\mathbf{n}\) with angle \(-\theta\). The two rotation axes lie on the same line. Additionally, when the rotation angle is \(\pi\), the sign of the axis doesn't change the result, it's the same rotation. Thus, since the set of all unit quaternions is encoding all possible axes and angles, this also means it has redundant elements which map to the same rotation. We'll get back to this in a jiffy.

For the following remark, imagine the \(2\)-sphere to more easily understand what's going on. If you have a line going through a sphere, it will intersect it in two places. This much is clear. When this line goes through the center of the sphere, the two intersections are said to be antipodal to each other. This is useful information because of what we've already discussed, namely \(R(\mathbf{n}, \theta) = R(\mathbf{-n}, -\theta)\). Given the 3-sphere formed by unit quaternions, these mutually redundant points are on opposite ends, they're antipodal points and need to be taken care of if we're going to show that the group of quaternions can be mapped to the group of all rotations. This is done with a simple mapping of one point to the other, which in topological freestyle lingo is called "gluing" them together. This leaves us in the topological space of the group of all rotations, which happens to be the real projective space, \(\mathbb{RP}^3\).

In case you're not familiar, the real projective space is, really, the mathematical take on something you might be very familiar with — perspective. Chronologically, the study of real projective spaces was inspired by painters, after they figured out how to paint their surroundings in a way that doesn't look like a doodle of a Bohemian peasant from 1403. In other words, when they figured out perspective. It comes down to having the set of all lines that pass through the origin in \(\mathbb{R}^n\) be points in \(\mathbb{RP}^{n-1}\). For example, when you have a 3-sphere in \(\mathbb{R}^4\) centered on the origin (much like the one we've been describing), in \(\mathbb{RP}^3\) antipodal points are the same one, as they sit on the line passing through the origin. In fact, this is one way to define a real projective space. By identifying antipodal points on a 3-sphere. And in turn, we can surmise that the group of all rotations \(SO(3)\) is diffeomorphic to \(\mathbb{RP}^3\), a fact not obvious without the quaternion group. Diffeomorphism is an isomorphism between smooth manifolds. In simple terms, it's a well-behaved reversible mapping between the topological space of the group of all 3D rotations to the real projective space.

We can now conclude there is a 2:1 mapping between the quaternion and rotation group because \( \hat{q}p\hat{q}^* = (-\hat{q})p(-\hat{q}^*) \) (remember the antipodal points), making the quaternion group a double cover of \(SO(3)\). More specifically, a surjective homomorphism. This means any change in the target space can be replicated in the source space. Which for you means no more annoyed users who had their spaceship camera controls freak out because they disagreed with your definition of which way is up.