How can I better understand this covariance equation? This is the equation given to me in the lectures:
$$S=\frac1N\sum_{n=1}^N(\mathbf{x}_n-\overline{\mathbf{x}})
(\mathbf{x}_n-\overline{\mathbf{x}})^T,$$
which doesn't make sense to me when I think about it. The $\mathbf{x}_n$ are $D$ dimensional vectors for $D$ features. So subtracting the mean will again result in a vector. Then taking the transpose multiplication will result in a scalar. So as far as I can see, this equation will end up as a scalar instead of a matrix.
How would I be able to correctly interpret this equation? I know there are different equations in matrix form, but I should use this equation.

Let's write this in matrix forms:
$$\overline{\mathbf{x}} =
\left[\matrix{\bar{x}_1\\\bar{x}_2\\\vdots\\\bar{x}_d}\right]
=
\frac{1}{N}\sum_{i=1}^N\mathbf{x}_i
=
\frac{1}{N}\sum_{i=1}^N\left[\matrix{x_{i,1}\\x_{i,2}\\\vdots\\x_{i,d}}\right]
$$
So the (biased) estimate of the covariance matrix is a square matrix, like below:
$$S=\frac1N\sum_{i=1}^N(\mathbf{x}_i-\overline{\mathbf{x}})
(\mathbf{x}_i-\overline{\mathbf{x}})^T
=
\frac1N\sum_{n=1}^N
\left(\left[\matrix{x_{i,1}\\x_{i,2}\\\vdots\\x_{i,d}}\right]-\left[\matrix{\bar{x}_1\\\bar{x}_2\\\vdots\\\bar{x}_d}\right]\right)
\left(\left[\matrix{x_{i,1}\\x_{i,2}\\\vdots\\x_{i,d}}\right]-\left[\matrix{\bar{x}_1\\\bar{x}_2\\\vdots\\\bar{x}_d}\right]\right)^T
=\\
\frac1N\sum_{i=1}^N
\left[\matrix{x_{i,1}-\bar{x}_1\\x_{i,2}-\bar{x}_2\\\vdots\\x_{i,d}-\bar{x}_d}\right]
\left[\matrix{x_{i,1}-\bar{x}_1&x_{i,2}-\bar{x}_2&\cdots&x_{i,d}-\bar{x}_d}\right]=\\
\frac1N\sum_{i=1}^N
\left[\matrix{(x_{i,1}-\bar{x}_1)^2 & (x_{i,1}-\bar{x}_1)(x_{i,2}-\bar{x}_2) & \cdots & (x_{i,1}-\bar{x}_1)(x_{i,d}-\bar{x}_d) \\ (x_{i,1}-\bar{x}_1)(x_{i,2}-\bar{x}_2) & (x_{i,2}-\bar{x}_2)^2 & \cdots & (x_{i,2}-\bar{x}_2)(x_{i,d}-\bar{x}_d) \\ \vdots & \vdots & \ddots & \vdots \\ (x_{i,1}-\bar{x}_1)(x_{i,d}-\bar{x}_d) & (x_{i,2}-\bar{x}_2)(x_{i,d}-\bar{x}_d) & \cdots & (x_{i,d}-\bar{x}_d)^2}\right]
$$
