This is a frequently asked question (see SO, Quora, Reddit). I
must admit that I haven't fully understood the generalized Brackprop of convolution
on 3D volumes, especially when using the chain of calculus without using
computational graphs. But I will give it a try.
I will use mostly the notation of the DL book (scalar $a$, matrix
$\boldsymbol{M}$, $\boldsymbol{\mathsf{M}}$ etc.)
and for some quantities the notation of Nielsen's book, e.g., $\delta^l_{...}$
for the neuron error at layer $l$ and I will consider the convolutional layer
to be a linear layer, i.e., the output is just the result of the convolution
without the non-linearity $z^l_{...}$. This is not relevant since it does not
influence the dimensionality of the operations. The convolution operations are
one-based index.
Furthermore, we will consider that the error at the output of the convolutional
layer is known, so
there is no need to consider any other layer. Further, to unclutter notation,
we will consider one instance batches and we will omit the batch dimension. Finally, we will assume padding zero and stride one (sorry, but otherwise, it becomes way too complicated).
Now, we have:
$l$ is the $l$th layer where $l=1$ is the first layer and $l=\ell$ is the last layer.
The Input tensor $\boldsymbol{\mathsf{X}} \in \mathbb{R}^{n \times
n \times
n}$ with dimensions width $w$, height $h$ and three channels $c$.
The Kernel tensor $\boldsymbol{\mathsf{K}} \in \mathbb{R}^{n \times
n \times
n \times n}$ with dimensions width $m$, height $n$, three channels $c$,
and $k$
kernels that will produce $k$ feature maps.
Every kernel $k$ in $\boldsymbol{\mathsf{K}}_{:,:,:,k}$ has weights
$w_{...}$.
The derivative of the loss function $J$ w.r.t. an input activation
$z^l_{x,y,z}$, say, is defined as
\begin{equation}
\delta^l_{x,y,z} = \frac{\partial J}{\partial \; z^l_{x,y,z}}
\end{equation}
First, we need to restate the equations in the corresponding form:
The 3D convolution acts such as one element of the output
$\mathsf{M}{\,}^{l}_{i,j,k}$ is the scalar resulting from multiplying and
adding the corresponding elements of the input with one kernel which weights
are
$w^l_{m,n,c,k}$. I.e.,
\begin{align}
{\big( (\boldsymbol{\mathsf{I}} * \boldsymbol{\mathsf{K}})(i,j,k)
\big)}^l
&= \sum_{c}
\sum_{m} \sum_{n} x^{l - 1}_{i + m - 1,j + n - 1,c} \cdot w^l_{m,n,c,k}
\label{eq:cross_correlation_3d_1}\\
z^l_{i,j,k} &= {\big( (\boldsymbol{\mathsf{I}} *
\boldsymbol{\mathsf{K}})(i,j,k) \big)}^l + b^l\\
\mathsf{M}{\,}^{l}_{i,j,k} &=
z^l_{i,j,k}.\label{eq:cross_correlation_3d_3}
\end{align}
Here $i,j,k$ are the output iterators.
The strategy we will follow is the one explained by Nielsen with some priority
changes, i.e., first we will derive the backpropagation rule
($\delta^l_{\cdot,\cdot,\cdot}$ in term of $\delta^{l +
1}_{\cdot,\cdot,\cdot}$).
Second, derive the error of the weights (derivative of the loss w.r.t. the
kernel weights) in terms of $\delta^{l + 1}_{\cdot,\cdot,\cdot}$.
The 2D case derivation suggests that, in the end, we will have some kind of
convolution of the layer error with the kernels.
Let us start with the error of one element at layer $\delta^l_{m,n,u}$ w.r.t.
to the error at layer $\delta^{l + 1}_{i,j,k}$. In order to apply the chain
rule, we
have to identify which output elements $z^{l + 1}_{i,j,k}$ are affected by the
input element
$z^l_{m,n,u}$, so that we can add the partial derivatives. One important
observation is that one input element affects all feature maps (albeit not all
elements), so we must sum over all output channels $k$. Another (trivial?)
observation is that the multi-channel $\boldsymbol{\mathsf{K}}$ does not move
across the depth dimension,
therefore the index does not contain offsets. After having considered the
channels, we must further consider the effect of the input element
$z^l_{m,n,u}$ in one
specific output channel $z^l_{:,:,k}$. Using the 2D analogy, we can say that
given a
specific
channels $k$, the input element affects some elements in that channel, namely
the
elements fixed at the position $m,n$, and that are the results of the kernel
having been slit by by $-a, -b$ (input perspective, we introduce auxiliary
indices $a$ and $b$). Finally, we know that all
kernels have the same dimensionality (they are packed in one 4D-tensor), thus,
we can use the same iterators $a,b$ for all kernels. These observations lead to
\begin{align}
\delta^l_{m,n,u} &= \sum_{k} \sum_a \sum_b \frac{\partial
J}{\partial \; z^{l +
1}_{m - a + 1,n - b + 1,k}}
\cdot \frac{\partial z^{l + 1}_{m - a + 1,n - b + 1,k}}{\partial
z^l_{m,n,u}}\label{eq:delta1}\tag{1}\\
&= \sum_{k} \sum_a \sum_b \delta^{l +
1}_{m - a + 1,n - b + 1,k}
\cdot \frac{\partial z^{l + 1}_{m - a + 1,n - b + 1,k}}{\partial
z^l_{m,n,u}}.\label{eq:delta2}\tag{2}
\end{align}
That was one of the most difficult steps, where I could have committed errors.
Feel free to correct me.
Here, we substituted, as usual, the error $\frac{\partial
J}{\partial \; z^{l + 1}_{m - a + 1,n - b + 1,k}}$ by its definition
$\delta^{l +
1}_{m - a + 1,n - b + 1,k}$.
You know where this is heading.
We now focus on the second term on the right and expand it. Here, again we
introduce auxiliary indices $p, q, r$, and use $x^l_{...} = z^l_{...}$ (as
mentioned, we deal with linear layers).
\begin{align}
\frac{\partial z^{l + 1}_{m - a + 1,n - b + 1,k}}{\partial z^l_{m,n,u}} &=
\frac{\partial}{\partial z^l_{m,n,u}} \Big(
\sum_{p} \sum_{q} \sum_{r} z^{l}_{(m - a + 1) + p - 1,(n - b + 1) + q - 1,r}
\cdot w^{l +
1}_{p,q,r,k}\Big)\label{eq1}\tag{3}\\
&=
\frac{\partial}{\partial z^l_{m,n,u}} \Big(
\sum_{p} \sum_{q} \sum_{r} z^{l}_{m - a + p,n - b + q,r}
\cdot w^{l + 1}_{p,q,r,k}\Big)\label{eq2}\tag{4}\\
&=
\frac{\partial}{\partial z^l_{m,n,u}} \Big( z^{l}_{m,n,u}
\cdot w^{l + 1}_{a,b,u,k}\Big)\label{eq3}\tag{5}\\
&= w^{l + 1}_{a,b,u,k}.\label{eq4}\tag{6}
\end{align}
In Eq. \ref{eq1} we expanded as mentioned above, then in Eq. \ref{eq2} we
reduced the indices as $(m - a + 1) + p - 1 = m - a + p, (n - b + 1) + q - 1 =
n - b + q$. In Eq. \ref{eq3} we observe that all partial derivatives w.r.t.
$z^l_{m,n,u}$, except when $a = p, b = q, r = u$. Accordingly, the weights must
have indices $a,b,u,k$. Then, plugging Eq. \ref{eq4} in Eq.\ref{eq:delta2}, we
have:
\begin{align}
\delta^l_{m,n,u} &= \sum_{k} \sum_a \sum_b \delta^{l +
1}_{m - a + 1,n - b + 1,k}
\cdot w^{l + 1}_{a,b,u,k}.\label{eq:final1}\tag{7}
\end{align}
The second most difficult part is how to interpret Eq.\ref{eq:final1}. As
mentioned above, it is clear that this is some kind of convolution between
the error at the end of the convolutional layer (the error it receives from the
upper layer) with the kernel tensor. One thing we can do is transform
(rotate?) the kernel to rearrange the indices, i.e.,
\begin{align}
\delta^l_{m,n,u} &= \sum_{k} \sum_a \sum_b \delta^{l +
1}_{m - a + 1,n - b + 1,k}
\cdot w^{l + 1}_{a,b,k,u}.\label{eq:final2}\tag{8}
\end{align}
In Eq.\ref{eq:final2}, we can see that we have to collect all partial
derivatives corresponding to elements in all channels $k$ having been generated
by channel $u$. This is somehow not intuitive and difficult to picture (at
least for me) because
all input channels are combined/merged to create a feature map. So, how can we
untangle a feature map to collect specific the input elements of one channel
that were melted to form the feature map? It turns out that one can if you
think in algebraic terms. The output is a tensor, and therefore, provided one
has the right indices, the operation can be conducted. This standard 3D
convolution may surely be expressed as a 3D cross-correlation. Ideas are
welcome.
An interesting concrete example is given by Aggarwal 2018 on page 355:
In order to understand the transposition above, consider a situation
in which we use 20 filters on the 3-channel RGB volume in order to
create an output volume of depth 20. While backpropagating, we will
need to take a gradient volume of depth 20 and transform to a
gradient volume of depth 3. Therefore, we need to create 3 filters for
backpropagation, each of which is for the red, green, and blue colors.
We pull out the 20 spatial slices from the 20 filters that are
applied to the red color invert them using the the approach of Figure
8.7, and then create a single 20-depth filter for backpropagating gradients with respect to the red slice. Similar approaches are used
for the green and blue slices. The transposition and inversion in
Equation 8.3 corresponds to these operations.
The derivation of the derivatives of the cost function w.r.t. the kernel
weights and bias are left as exercise ;)
