#### Note

The ``@`` symbol indicates matrix multiplication on NumPy arrays, and was\n introduced in Python 3.5 / NumPy 1.10. The notation ``plot(*point)`` uses\n Python `argument expansion`_ to \"unpack\" the elements of ``point`` into\n separate positional arguments to the function. In other words,\n ``plot(*point)`` expands to ``plot(3, 2, 5)``.

\n\nNotice that we used matrix multiplication to compute the projection of our\npoint $(3, 2, 5)$ onto the $x, y$ plane:\n\n\\begin{align}\\left[\n \\begin{matrix} 1 & 0 & 0 \\\\ 0 & 1 & 0 \\\\ 0 & 0 & 0 \\end{matrix}\n \\right]\n \\left[ \\begin{matrix} 3 \\\\ 2 \\\\ 5 \\end{matrix} \\right] =\n \\left[ \\begin{matrix} 3 \\\\ 2 \\\\ 0 \\end{matrix} \\right]\\end{align}\n\n...and that applying the projection again to the result just gives back the\nresult again:\n\n\\begin{align}\\left[\n \\begin{matrix} 1 & 0 & 0 \\\\ 0 & 1 & 0 \\\\ 0 & 0 & 0 \\end{matrix}\n \\right]\n \\left[ \\begin{matrix} 3 \\\\ 2 \\\\ 0 \\end{matrix} \\right] =\n \\left[ \\begin{matrix} 3 \\\\ 2 \\\\ 0 \\end{matrix} \\right]\\end{align}\n\nFrom an information perspective, this projection has taken the point\n$x, y, z$ and removed the information about how far in the $z$\ndirection our point was located; all we know now is its position in the\n$x, y$ plane. Moreover, applying our projection matrix to *any point*\nin $x, y, z$ space will reduce it to a corresponding point on the\n$x, y$ plane. The term for this is a *subspace*: the projection matrix\nprojects points in the original space into a *subspace* of lower dimension\nthan the original. The reason our subspace is the $x,y$ plane (instead\nof, say, the $y,z$ plane) is a direct result of the particular values\nin our projection matrix.\n\n\n### Example: projection as noise reduction\n\nAnother way to describe this \"loss of information\" or \"projection into a\nsubspace\" is to say that projection reduces the rank (or \"degrees of\nfreedom\") of the measurement \u2014 here, from 3 dimensions down to 2. On the\nother hand, if you know that measurement component in the $z$ direction\nis just noise due to your measurement method, and all you care about are the\n$x$ and $y$ components, then projecting your 3-dimensional\nmeasurement into the $x, y$ plane could be seen as a form of noise\nreduction.\n\nOf course, it would be very lucky indeed if all the measurement noise were\nconcentrated in the $z$ direction; you could just discard the $z$\ncomponent without bothering to construct a projection matrix or do the matrix\nmultiplication. Suppose instead that in order to take that measurement you\nhad to pull a trigger on a measurement device, and the act of pulling the\ntrigger causes the device to move a little. If you measure how\ntrigger-pulling affects measurement device position, you could then \"correct\"\nyour real measurements to \"project out\" the effect of the trigger pulling.\nHere we'll suppose that the average effect of the trigger is to move the\nmeasurement device by $(3, -1, 1)$:\n\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"trigger_effect = np.array([[3, -1, 1]]).T"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Knowing that, we can compute a plane that is orthogonal to the effect of the\ntrigger (using the fact that a plane through the origin has equation\n$Ax + By + Cz = 0$ given a normal vector $(A, B, C)$), and\nproject our real measurements onto that plane.\n\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# compute the plane orthogonal to trigger_effect\nx, y = np.meshgrid(np.linspace(-1, 5, 61), np.linspace(-1, 5, 61))\nA, B, C = trigger_effect\nz = (-A * x - B * y) / C\n# cut off the plane below z=0 (just to make the plot nicer)\nmask = np.where(z >= 0)\nx = x[mask]\ny = y[mask]\nz = z[mask]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Computing the projection matrix from the ``trigger_effect`` vector is done\nusing [singular value decomposition](svd_) (SVD); interested readers may\nconsult the internet or a linear algebra textbook for details on this method.\nWith the projection matrix in place, we can project our original vector\n$(3, 2, 5)$ to remove the effect of the trigger, and then plot it:\n\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# compute the projection matrix\nU, S, V = svd(trigger_effect, full_matrices=False)\ntrigger_projection_matrix = np.eye(3) - U @ U.T\n\n# project the vector onto the orthogonal plane\nprojected_point = trigger_projection_matrix @ point\nprojected_vector = trigger_projection_matrix @ vector\n\n# plot the trigger effect and its orthogonal plane\nax = setup_3d_axes()\nax.plot_trisurf(x, y, z, color=\"C2\", shade=False, alpha=0.25)\nax.quiver3D(\n *np.concatenate([origin, trigger_effect]).flatten(),\n arrow_length_ratio=0.1,\n color=\"C2\",\n alpha=0.5,\n)\n\n# plot the original vector\nax.plot(*vector, color=\"k\")\nax.plot(*point, color=\"k\", marker=\"o\")\noffset = np.full((3, 1), 0.1)\nax.text(*(point + offset).flat, \"({}, {}, {})\".format(*point.flat), color=\"k\")\n\n# plot the projected vector\nax.plot(*projected_vector, color=\"C0\")\nax.plot(*projected_point, color=\"C0\", marker=\"o\")\noffset = np.full((3, 1), -0.2)\nax.text(\n *(projected_point + offset).flat,\n \"({}, {}, {})\".format(*np.round(projected_point.flat, 2)),\n color=\"C0\",\n horizontalalignment=\"right\",\n)\n\n# add dashed arrow showing projection\narrow_coords = np.concatenate([point, projected_point - point]).flatten()\nax.quiver3D(\n *arrow_coords,\n length=0.96,\n arrow_length_ratio=0.1,\n color=\"C1\",\n linewidth=1,\n linestyle=\"dashed\",\n)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Just as before, the projection matrix will map *any point* in $x, y, z$\nspace onto that plane, and once a point has been projected onto that plane,\napplying the projection again will have no effect. For that reason, it should\nbe clear that although the projected points vary in all three $x$,\n$y$, and $z$ directions, the set of projected points have only\ntwo *effective* dimensions (i.e., they are constrained to a plane).\n\n.. admonition:: Terminology\n :class: sidebar note\n\n In MNE-Python, the matrix used to project a raw signal into a subspace is\n usually called a :term:`projector` or a *projection\n operator* \u2014 these terms are interchangeable with the term *projection\n matrix* used above.\n\nProjections of EEG or MEG signals work in very much the same way: the point\n$x, y, z$ corresponds to the value of each sensor at a single time\npoint, and the projection matrix varies depending on what aspects of the\nsignal (i.e., what kind of noise) you are trying to project out. The only\nreal difference is that instead of a single 3-dimensional point $(x, y,\nz)$ you're dealing with a time series of $N$-dimensional \"points\" (one\nat each sampling time), where $N$ is usually in the tens or hundreds\n(depending on how many sensors your EEG/MEG system has). Fortunately, because\nprojection is a matrix operation, it can be done very quickly even on signals\nwith hundreds of dimensions and tens of thousands of time points.\n\n\n\n## Signal-space projection (SSP)\n\nWe mentioned above that the projection matrix will vary depending on what\nkind of noise you are trying to project away. Signal-space projection (SSP)\n:footcite:`UusitaloIlmoniemi1997` is a way of estimating what that projection\nmatrix should be, by\ncomparing measurements with and without the signal of interest. For example,\nyou can take additional \"empty room\" measurements that record activity at the\nsensors when no subject is present. By looking at the spatial pattern of\nactivity across MEG sensors in an empty room measurement, you can create one\nor more $N$-dimensional vector(s) giving the \"direction(s)\" of\nenvironmental noise in sensor space (analogous to the vector for \"effect of\nthe trigger\" in our example above). SSP is also often used for removing\nheartbeat and eye movement artifacts \u2014 in those cases, instead of empty room\nrecordings the direction of the noise is estimated by detecting the\nartifacts, extracting epochs around them, and averaging. See\n`tut-artifact-ssp` for examples.\n\nOnce you know the noise vectors, you can create a hyperplane that is\northogonal\nto them, and construct a projection matrix to project your experimental\nrecordings onto that hyperplane. In that way, the component of your\nmeasurements associated with environmental noise can be removed. Again, it\nshould be clear that the projection reduces the dimensionality of your data \u2014\nyou'll still have the same number of sensor signals, but they won't all be\n*linearly independent* \u2014 but typically there are tens or hundreds of sensors\nand the noise subspace that you are eliminating has only 3-5 dimensions, so\nthe loss of degrees of freedom is usually not problematic.\n\n\n## Projectors in MNE-Python\n\nIn our example data, `SSP