MGM.08 Subspaces, Spans , Dimension, Linear Independence, Basis, Orthonormal Basis
TUTORIALS

Needs["LinearAlgebra`MatrixManipulation`"];
Needs["LinearAlgebra`Orthogonalization`"];
Needs["Graphics`Colors`"];
Off[General::"spell"];
Off[General::"spell1"];
Off[ParametricPlot3D::"ppcom"];
Off[ParametricPlot::"ppcom"];
Off[Plot::"plnr"];
SetOptions[Limit, Analytic -> True];
$Pre = Chop;

CMView = {2.7, 1.6, 1.2};

Vector Drawer

\(\(\( (*\ \(\(:\)\(Date : \ \ Copyright\ 2001\)\), \ 
      Math\ Everywhere, \ \(\(Inc\)\(.\)\)\ 
*) \)\(
\)\( (*
Original\ Conception\ to\ create\ a\ vector\ drawer\ is\ due\ to\ Bill\ Davis . 
Current\ code\ is\ by\ Bruce\ \(\(Carpenter\)\(.\)\)
*) \)\(
\)\( (*\ \ \(\(\(:\)\(Drawbacks : \ \ Since\ the\ vectors\ are\ drawn\ out\ of\)\)\ 

      \ \ \ \ \ context\), \ 
      the\ shapes\ of\ the\ arrowheads\ are\ only\ 
\ \ \ \ \ correct\ with\ AspectRatio\ set\ to\ \(\(Automatic\)\(.\)\)
*) \)\)\)

BeginPackage["Vector3D`","Graphics`Colors`"];

Axes3D::usage = "Axes3D[a, b] creates a Graphics3D object of cartesian axes with x, y, and z running from -a/3 to a, and with axes labels b units beyond the tips of the axes.  Axes3D[a] is Axes3D[a, a/8].";

Perpend::usage = "Perpend[a] returns a unit vector which is perpendicular to the vector a.";

Arrow::usage = "Arrow[a] produces a 2 or 3 dimensional vector from the origin to a. Arrow[a, Tail -> tail] gives a vector from tail to a.";
	
Vector::usage = "Vector[a, b] gives a 2D or 3D vector from a to b.";

ArrowHead::usage = "ArrowHead[a, vec] produces an arrowhead with its tip placed at the point a, pointing in the direction of the vector vec.";

Tail::usage = "Tail -> point puts the tail of the vector at point.";

Norm::usage = "Norm[vec] returns the Euclidean length of vec.";

Normalize::usage = "Normalize[vec] returns a unit vector in the direction of vec, provided vec is nonzero.";

Aperture::usage = "Aperture is the ratio of the radius of the base to the length of the head of a vector.";

TipSize::usage = "TipSize specifies an absolute size for the inner tip length of the head of a vector.";

TipRatio::usage = "TipRatio is the ratio of the inner length of the arrowhead to the length of the shaft of a vector.";

HeadRatio::usage = "HeadRatio is the ratio of the outer length of the arrowhead to the length of the shaft of a vector.";

HeadSize::usage = "HeadSize specifies an absolute size for the head of a vector.";

ScaleFactor::usage = "ScaleFactor specifies the amount to scale a vector in length and may be either a number or a pure function.";

ZeroVectorPointSize::usage = "ZeroVectorPointSize is the size of the point used to represent a zero vector.";

Options[Vector] = Options[Arrow] = Options[ArrowHead] =
{HeadRatio -> 0.18, TipRatio -> 0.14, Aperture -> 0.3, ScaleFactor -> 1, HeadType -> Polygon, ShaftQ -> True, EdgesQ -> True, VectorColor -> RGBColor[0,0,1], ShaftWidth -> 0.005, ZeroVectorPointSize -> 0.01};
		
Begin["`Private`"];

SetOptions[ParametricPlot,AspectRatio -> Automatic];
SetOptions[Plot,AspectRatio -> Automatic];
SetOptions[Graphics,AspectRatio -> Automatic];

Axes3D[u_,v_] := 
      Graphics3D[{
            {Blue,Line[{{-u/3,0,0},{u,0,0}}]},
                    Text["x",{u + v,0,0}],
            {Blue,Line[{{0,-u/3,0},{0,u,0}}]},
                    Text["y",{0, u + v,0}],
            {Blue,Line[{{0,0,-u/3},{0,0,u}}]},
                    Text["z",{0,0,u + v}]}]; 
Axes3D[u_] := Axes3D[u,u/8];
				
Norm[a_?VectorQ] := N[Sqrt[a.a]];

Normalize[a_?VectorQ] := If[a.a == 0, a, a/Norm[a]];

Perpend[a_?VectorQ] := Normalize[{a[[2]],-a[[1]]}]/;Length[a]==2;

Perpend[a_?VectorQ]:=
    Normalize[If[a[[1]]==0,{1,0,0},{a[[2]],-a[[1]],0}]]/;
      Length[a]==3;
      
Base3D = Table[{Cos[2 k Pi/8.], Sin[2 k Pi/8.]},{k,0,8}]; 

Vector[from_?VectorQ, to_?VectorQ, opts___] := Module[
   {sq, sw, ht, hs, hr, tr, ap, sf, vc, zvps, diff = to - from, tip},
   {sq, sw, ht, hs, hr, ts, tr, ap, sf, vc, zvps} = {ShaftQ, ShaftWidth, HeadType, HeadSize, HeadRatio, TipSize, TipRatio, Aperture, ScaleFactor, VectorColor, ZeroVectorPointSize} /. {opts} /. Options[Vector];
   If[from==to, Return[Graphics[{PointSize[zvps],vc,Point[from]}]]];
   diff = If[NumberQ[N[sf]], sf diff, sf[diff]];
   tip = from + diff;
   If[NumberQ[hs], hr = hs/Norm[diff]; tr = 0.75 hr];
   If[NumberQ[ts], tr = ts/Norm[diff]];
   Graphics[{If[sq,{vc,Thickness[sw],Line[{from,tip - tr diff}]},{}],
      With[{trans = ap hr Norm[diff] Perpend[diff]},
      {vc, Thickness[sw], ht[{tip,
                   (tip - hr diff + trans),
                   (tip - tr diff),
                   (tip - hr diff - trans),tip}]}]}]
] /; Length[from]==Length[to]==2

Vector[from_?VectorQ, to_?VectorQ, opts___] := Module[
   {sq, eq, sw, hs, hr, tr,ap, sf, vc, zvps, diff = to - from, tip},
   {sq, eq, sw, hs, hr, ap, sf, vc, zvps} = {ShaftQ, EdgesQ, ShaftWidth, HeadSize, HeadRatio, Aperture, ScaleFactor, VectorColor, ZeroVectorPointSize} /. {opts} /. Options[Vector];
   If[from==to, Return[Graphics3D[{PointSize[zvps],vc,Point[from]}]]];
   diff = If[NumberQ[N[sf]], sf diff, sf[diff]];
   tip = from + diff;
   If[NumberQ[hs], hr = hs/Norm[diff]];
   tr = hr/2;
   Graphics3D[{If[sq,{vc,Thickness[sw],Line[{from,tip-tr diff}]},{}], If[eq, {}, EdgeForm[]], {SurfaceColor[vc],
      With[{temp = Perpend[diff]},
      Polygon[Append[#,tip]]& /@ Partition[(tip - hr diff + #)&/@((ap hr Norm[diff] Base3D) . {temp,Normalize[Cross[temp,diff]]}),2,1]]}}]
] /; Length[from]==Length[to]==3


Arrow[a_, Tail->b_, opts___] := Vector[b,a+b,opts]

Arrow[a_, opts___] := Vector[Table[0,{Length[a]}],a,opts]

ArrowHead[a_, b_, opts___] := Vector[a-b,a,ShaftQ->False,opts];

End[];

EndPackage[];
Null

Graphics`Colors`GosiaGreen = 
        RGBColor[0,  0.392187,  0];

ThreeAxes[u,v]

ThreeAxes[u_, v_] := 
        Graphics3D[{{Blue, Line[{{-u, 0, 0}, {u, 0, 0}}]}, 
            Text["x", {u + v, 0, 0}], {Blue, 
              Line[{{0, -u, 0}, {0, u, 0}}]}, 
            Text["y", {0, u + v, 0}], {Blue, 
              Line[{{0, 0, 0}, {0, 0, u}}]}, 
            Text["z", {0, 0, u + v}]}];
ThreeAxes[u_] := ThreeAxes[u, (u/8)];
ThreeAxes::"usage" = "ThreeAxes[a,b] makes a standard cartesian \
axis graphics object with x, y, and z running from -a to a, and with axis \
labels b units beyond the tips of the axes.  ThreeAxes[a] is \
ThreeAxes[a,a/8].";

Needs["Statistics`ContinuousDistributions`"];

Clear[matrixplucker];
matrixplucker[expr_, 
        vars_] := Simplify[
          Outer[D, {expr}, 
            vars]][[1]]

T.1) Orthonormal bases are bases consisting of mutually perpendicular unit vectors - they are perpendicular frames.
The advantages and pleasures of having an orthonormal basis.
Getting orthonormal bases via SVD hanger frames for the spanner matrix

Most folks say that a basis for a subspace S of nD is any linearly independent spanning set
           {spanner[1], spanner[2], . . ., 
        spanner[k]}
for S. 
The same folks say that a basis  {spanner[1], spanner[2], . . ., 
        spanner[k]} for a subspace S of nD is an orthonormal basis for S if 
            {spanner[1], spanner[2], . . ., 
        spanner[k]}
is also a perpendicular frame (This is the same as saying.
             spanner[i] . spanner[i] = 1 for i = 1, 2, ... k
and
             spanner[i] . spanner[j] = 0 for i ≠ j.)

T.1.a.i) Life in the orthonormal world is good

List some of the advantages you get when have an orthonormal basis for a subspace S of kD.

Everything gets really easy.
Some examples:

Example 1: Easy calculation of basis coefficients

If  {spanner[1], spanner[2], . . ., 
        spanner[k]} is an orthonormal basis of a subspace S of nD, then every X in S is of the form
         X = (∑_(i = 1)^k) x[i] spanner[i]  
 with  x[i] = X . spanner[i].

Example 2: Easy calculation closest vectors

If  {spanner[1], spanner[2], . . ., 
        spanner[k]} is an orthonormal basis of a subspace S of nD and X is a given vector in nD, then the vector in S closest to X is
         XSclosest = (∑_(i = 1)^k) x[i] spanner[i]  
 with  x[i] = X . spanner[i].

Example 3: Easy calculation Sprojection matrix

If  {spanner[1], spanner[2], . . ., 
        spanner[k]} is an orthonormal basis of a subspace S of nD then you can go directly to the Sprojection matrix
            Sprojection = hanger . aligner
where
        aligner =  {spanner[1], spanner[2], . . ., 
        spanner[k]} 
and
        hanger =  \({spanner[1], spanner[2], \(\(.\)\(\ \)\(.\)\(\ \)\(.\)\), spanner[k]}\^t\)
so that 
          spanner[i] is the ith horizontal row of aligner
and
          spanner[j] is the jth vertical column of hanger.

Life in the orthonormal world is about as good as it gets.

T.1.a.ii) What happens when you don't have an orthonormal basis

Do those nice clean formulas discussed in part i) work for a given basis that is not an orthonormal basis?

T.1.b) Using SVD analysis of the Spanner Matrix to get an orthonormal basis

A subspace S of 6D is defined by the following given spanning set:

Clear[i, spanner];
spanner[1] = {3.43, -0.57, 2.19, -0.94, 1.28};
spanner[2] = {0.97, 1.29, 0.83, 5.23, 1.26};
spanner[3] = {-6.14, 
          3.08, -5.00, -1.97, -0.96};
spanner[4] = {-3.68, 
          1.22, -3.64, -8.14, -0.94};
spanner[5] = {2.71, -2.51, 2.81, 2.91, -0.32};
spanners = Table[spanner[i], {i, 1, 5}]

Use SVD analysis of the corresponding spanner matrix to come up with an orthonormal basis for S.
Use your orthonormal basis to build the Sprojection matrix.

The hangerframe for the Spanner Matrix is a  perpendicular frame (orthonormal basis) that spans S. 
This hangerframe is guaranteed to be an orthonormal basis of S.
Here it is:

SpannerMatrix = Transpose[spanners];
Sdim=Length[SingularValues[SpannerMatrix][[2]]];

Clear[ortho, k];
ortho[k_]:=SingularValues[SpannerMatrix][[1,k]]
Table[ortho[k], {k, 1, Sdim}]

You can use this to make the Sprojection very quickly from rank one matrices

rankOne[k_] := Transpose[{ortho[k]}].{ortho[k]};
sum=(∑_(k = 1)^Sdim)rankOne[k];
MatrixForm[
      sum]

Or you can make the Sprojection matrix directly from SVD:

aligner  =  {ortho[1], ortho[2], 
          ortho[3]};
hanger  =  
        Transpose[{ortho[1], ortho[2], ortho[3]}];
SVDSprojection  =  hanger . aligner;
		MatrixForm[SVDSprojection]

Compare:

Sprojection = 
        SpannerMatrix . PseudoInverse[SpannerMatrix];
MatrixForm[Sprojection]

Nailed it.

T.2)  Vivid Graphics: 
Alien plots coming from projections of highD surfaces onto  three dimensional subspaces

T.2.a.i) Looking at a three dimensional photo of a 4D surface

Here's a function f[s,t]:

Clear[f, s, t];
f[s_, t_]  =  {Cos[s], Sin[s], Cos[t], Sin[t]}

If you could plot this function, the result would be a surface in 4D which some advanced folks call by the name "Clifford Torus." 
Unfortunately, you can't plot this surface, but you can go with any three dimensional subspace S of 4D and plot the projection of this 4D surface onto S.
Just how do you go about doing this?

You do this by going with a spanning set for your three dimensional subspace S or 4D.. 
Here's a sample:

Clear[s, i, spanner];
spanner[1] = {0.327337, -0.121305, 0.39909, 
          0.397503};
spanner[2] = {0.871239, 0.659994, -0.659799, 
          0.497028};
spanner[3] = {0.028996, 0.684016, 
          0.746671, -0.391136};
spanners = Table[spanner[i], {i, 1, 3}]

Next you go from this spanning set for S to an orthonormal basis for S

SpannerMatrix = Transpose[spanners];
Clear[ortho, k];
ortho[
          k_] := SingularValues[SpannerMatrix][[1, k]];
Table[ortho[k], {k, 1, 3}]

To plot the projection of the surface onto S, you plot
	{f[s, t] . ortho[1],  f[s, t] . ortho[2],  f[s, t] . ortho[3]}
and label the corresponding axes accordingly.
Here you go:

Clear[surfaceprojectionplotter, f, s, t];
f[s_, t_]  =  {Cos[s], Sin[s], Cos[t], Sin[t]};

\(\(surfaceprojectionplotter[s_, 
            t_]\  = \ 
	{f[s, t] . ortho[1], f[s, t] . ortho[2], 
            f[s, t] . ortho[3]};\)\(
\)\)
{slow, shigh}  =  {0, 2 π};
{tlow, 
          thigh}  =  {0,  2 π};
\(ParametricPlot3D[
        surfaceprojectionplotter[s, t], {s, slow, shigh}, {t, tlow, thigh}, 
        Axes -> True, 
        AxesLabel -> {"ortho[1]", "ortho[2]", "ortho[3]"}, 

        ViewPoint -> CMView];\)

How does that grab you, sports fans

See what happens when you project onto random three dimensional subspaces of 4D:

Clear[s, i, spanner];
spanner[1] = 
        Table[Random[Real, {-1, 1}], {i, 1, 4}];
spanner[2] = 
        Table[Random[Real, {-1, 1}], {i, 1, 4}];
spanner[3] = 
        Table[Random[Real, {-1, 1}], {i, 1, 4}];
spanners = 
        Table[spanner[i], {i, 1, 3}];
SpannerMatrix = 
        Transpose[spanners];
Clear[ortho, 
        k];
ortho[k_] := SingularValues[SpannerMatrix][[1, k]];

Clear[surfaceprojectionplotter, f, s, t];
f[s_, t_]  =  {Cos[s], Sin[s], Cos[t], Sin[t]};

\(\(surfaceprojectionplotter[s_, 
            t_]\  = \ 
	{f[s, t] . ortho[1], f[s, t] . ortho[2], 
            f[s, t] . ortho[3]};\)\(
\)\)
{slow, shigh}  =  {0, 2 π};
{tlow, 
          thigh}  =  {0,  2 π};
\(ParametricPlot3D[
        surfaceprojectionplotter[s, t], {s, slow, shigh}, {t, tlow, thigh}, 
        Axes -> True, 
        AxesLabel -> {"ortho[1]", "ortho[2]", "ortho[3]"}, 

        ViewPoint -> CMView];\)

Each of these views is a three dimensional slice of a 4D surface.
Almost as good as the Hubble telescope.

T.2.a.ii) Why it's just a little bit interesting

Take the function f[s,t] used to generate the plots in parts i and ii):

f[s, t]

Look at:

f[s, t] . f[s, t]

Apply trig identities:

TrigExpand[f[s, t] . f[s, t]]

This tells you that no matter what s and t are, you are guaranteed that f[s,t] is on the surface of the 4D sphere of radius (Sqrt[2]) centered at {0,0,0,0}.
Why is this fact just a little bit mind-boggling?

This tells you that all the plots produced above are projections of pieces of the 4D sphere of radius (Sqrt[2]) centered at {0,0,0,0} onto  three dimensional subspace of 4D. 
This might make you expect that these projections would look like pieces of spheres.-after all projections of pieces of spheres in 3D onto two dimensional subspaces of 3D look like pieces of circles.
But they don't:

Clear[s, i, spanner];
spanner[1] = 
        Table[Random[Real, {-1, 1}], {i, 1, 4}];
spanner[2] = 
        Table[Random[Real, {-1, 1}], {i, 1, 4}];
spanner[3] = 
        Table[Random[Real, {-1, 1}], {i, 1, 4}];
spanners = 
        Table[spanner[i], {i, 1, 3}];
SpannerMatrix = 
        Transpose[spanners];
Clear[ortho, 
        k];
ortho[k_] := SingularValues[SpannerMatrix][[1, k]];

Clear[surfaceprojectionplotter, f, s, t];
f[s_, t_]  =  {Cos[s], Sin[s], Cos[t], Sin[t]};

\(\(surfaceprojectionplotter[s_, 
            t_]\  = \ 
	{f[s, t] . ortho[1], f[s, t] . ortho[2], 
            f[s, t] . ortho[3]};\)\(
\)\)
{slow, shigh}  =  {0, 2 π};
{tlow, 
          thigh}  =  {0,  2 π};
\(ParametricPlot3D[
        surfaceprojectionplotter[s, t], {s, slow, shigh}, {t, tlow, thigh}, 
        Axes -> True, 
        AxesLabel -> {"ortho[1]", "ortho[2]", "ortho[3]"}, 

        ViewPoint -> CMView];\)

Moral:  4D spheres are alien beasts but they are real.

T.3) The perpendicular complement (S^⊥) of a subspace S of kD is the subspace of kD consisting of all of all Y that are perpendicular to everything in S.
The perpendicular complement (S^⊥) of a subspace S is also the subspace of kD consisting of all hits with   Identity - Sprojection.
How to take a kD vector Y and put it in the form  Y = YS + Y(S^⊥) with YS in S , Y(S^⊥) in (S^⊥) and  
                           ||Y|| = (Sqrt[ ||YS\( || \^2\) +  ||(YS^⊥)\( || \^2\)])
Extending a basis of subspace of kD to a basis of all of kD

T.3.a.i) The perpendicular complement (S^⊥) of a subspace S of kD is the subspace of kD consisting of all hits with 
                  Identity - Sprojection

The perpendicular complement (S^⊥) of a subspace of kD is the subspace of kD consisting of all hits with 
                  Identity - Sprojection. 
A handy spanning set for (S^⊥) is the set of vertical columns of
                  Identity - Sprojection. 
                  
Here's a spanning set for a subspace of 3D:

Clear[s, i, spanner];
spanner[1] = {1.2, 0.9, 0};
spanner[2] = {-1.0, 1.2, -1.5};
spanners = Table[spanner[i], {i, 1, 2}];

And a plot of part of this subspace:

h = 3;
k = 2;
Clear[s, t];
\(\(planeplot = 
        ParametricPlot3D[
          s\ spanner[1] + t\ spanner[2], {s, \(-h\), h}, {t, \(-k\), k}, 
          PlotPoints -> {2, 2}, 
          DisplayFunction -> Identity];\)
\)
\(\(Slabel = 
        Graphics3D[{NavyBlue, 
            Text["S", 
              0.5\ spanner[1] - 
                0.5\ spanner[2]\  + \ {0, 0, 0.3}]}];\)
\)
Splot = 
        Show[planeplot, Slabel, Axes3D[2], ViewPoint -> CMView, 
          PlotRange -> All, Axes -> True, 
          AxesLabel -> {"x", "y", "z"}, Boxed -> False, 
          PlotLabel -> "Piece of S", 
          DisplayFunction -> $DisplayFunction];

Throw some of the members of (S^⊥) into this plot and describe what you see.

Set up the Sprojection matrix:

SpannerMatrix  =  Transpose[spanners];
Sprojection  =  
        SpannerMatrix . 
          PseudoInverse[SpannerMatrix];
MatrixForm[Sprojection]

(S^⊥) consists of all hits with Identity - Sprojection. 
See some::

a =  
        3;
point[k_]  :=  {Random[Real, {-a, a}], 
          Random[Real, {-a, a}], Random[Real, {-a, a}]};
randomhits =Table[Graphics3D[{PointSize[0.02],CadmiumOrange,Point[(IdentityMatrix[3]  -  Sprojection).point[k]]}],{k, 1, 300}];
label = Graphics3D[{Black,Text["\!\(S\^\*
StyleBox[\"⊥\",\nFontSize->16]\)",
(IdentityMatrix[3]  -  
                        Sprojection) . (1.1 {a, a, a})]}];
both =Show[Splot,randomhits,label,PlotLabel->"S and \!\(S\^\*
StyleBox[\"⊥\",\nFontSize->16]\)"];

Look at everything from a view point parallel to S:

Show[both, ViewPoint -> 12 spanner[1], Axes -> False];

The vectors in (S^⊥) are perpendicular to everything in S.
If you want to see why this is automatic, then go on.

T.3.a.ii)  If Y is in (S^⊥),  then Y is perpendicular to everything in S

Explain this:
Given a subspace S of kD and a vector X in (S^⊥), then it's automatic that X  is perpendicular to everything in S.

Go with an X in (S^⊥). This means X is a hit with Identity - Sprojection.
In other words, there is guaranteed to be a Y in kD with
             X = (Identity - Sprojection).Y = Y - Sprojection.Y
Remember this: Sprojection.Y is the vector in S that is closest to Y. 
And because perpendicular distance is the closest distance,
             X = Y - Sprojection.Y 
is perpendicular to everything in S.

T.3.a.iii)  You can decompose Y = YS + Y(S^⊥) with YS in S , Y(S^⊥) in (S^⊥) and  
                           ||Y|| = (Sqrt[ ||YS\( || \^2\) +  ||(YS^⊥)\( || \^2\)])

Explain this:
Given a subspace S of kD and given any kD vector Y, you can decompose 
              Y = YS + Y(S^⊥) 
with YS in S and with Y(S^⊥) in (S^⊥) and with

||Y|| = (Sqrt[ ||YS\( || \^2\) +  ||(YS^⊥)\( || \^2\)]).

Go with your Y and set
           YS = Sprojection.Y 
and set
           Y(S^⊥) = (Identity - Sprojection).Y  = Y - Sprojection.Y.
It' s now automatic that 
           Y = YS +  Y(S^⊥) 
and that YS is in S and Y(S^⊥) is in (S^⊥).
And:

\(\(\(||\)\(Y . Y\)\( || \^2\)\)\) = Y.Y = (YS  + (YS^ ⊥) ). (YS  +(YS^ ⊥) )

= (YS.YS) + (YS.(YS^ ⊥) ) + ((YS^ ⊥) .YS) + ((YS^ ⊥) .(YS^ ⊥) )
             = (YS.YS)  + 0 + 0 + ((YS^ ⊥) .(YS^ ⊥) )

= \(\(\(||\)\(YS\)\( || \^2\)\(\ \)\(+\ \ \(\(||\)\(YS\^ ⊥ \)\( || \^2\)\)\)\)\).
So

||Y|| = \(\@\(\(\ \)\(\(||\)\(YS\)\( || \^2\)\(\ \)\(+\ \ \(\(||\)\(YS\^ ⊥ \)\( || \^2\)\)\)\)\)\)

T.3.a.iv)  If Y is perpendicular to everything in S, then Y is in  (S^⊥)

Explain this:
Go with a subspace S of kD. 
If a kD vector Y is perpendicular to everything in S, then it's automatic that  Y is in (S^⊥).

Go with a Y that is perpendicular to everything in S.
Decompose
                Y = YS + Y(S^⊥) 
with YS in S and Y(S^⊥) in (S^⊥).

The goal is to explain why YS = {0,0,. . .,0}, forcing Y = Y(S^⊥)  and putting Y into (S^⊥).

To this end, take the equation
               Y = YS + Y(S^⊥).
Dot YS through this equation to get
          YS.Y = YS.YS + YS.(YS^ ⊥) .
This is the same as
             0 =  YS.YS  + 0

This results in 
            0 =  YS.YS   = \(\(\(||\)\(YS\)\( || \^2\)\)\)
and leads to the inescapable conclusion that
            YS = {0,0,. . . ,0}.

T.3.a.v)    (S^⊥) consists exactly of all Y that are perpendicular to everything in S

Explain this:
Go with a subspace S of kD. You are guaranteed that (S^⊥) consists exactly of all Y that are perpendicular to everything in S.

Put PS equal to the set of vectors in kD that are perpendicular to everything in S. 

Part ii) above says that  if Y is in (S^⊥),  then Y is perpendicular to everything in S.
So    (S^⊥)⊂ PS.

Part iv) above says that  if Y is perpendicular to everything in S, then Y is in  (S^⊥).
So    PS ⊂ (S^⊥).
The conclusion:
         PS = (S^⊥).

T.3.a.vi) The projection matrix for Perp[S] is (Identity - Sprojection)

Given a subspace S of kD,you can express the projection matrix for Perp[S] as 
           Perp[S] projection  =  IdentityMatrix - Sprojection.
Explain why this works.

Take any vector X in kD
Put 
           Y = (IdentityMatrix - Sprojection) . X
To explain why
           Perp[S] projection = (IdentityMatrix - Sprojection),
 you've got to explain why
->Y is in Perp[S] 
and why
-> Y is closer to X than any other member of Perp[S].

Explanation of why Y is in Perp[S]

Because Sprojection.X is the closest member of S to X and because perpendicular distance is the shortest distance.
           Y = (IdentityMatrix - Sprojection) . X   =    X - 
            Sprojection . X
is perpendicular to everything in S. This tells you that Y is Perp[S].

Explanation of why Y is closer to X than any other member of Perp[S].

Because perpendicular distance is the shortest distance, this is the same as explaining why X - Y is perpendicular to everything in Perp[S]. 
To do this, take any Z in Perp[S]  and look at this little calculation:
            Z.(X - Y) = Z.( X - (X - Sprojection.X)) = Z.(Sprojection.X) = 0
because Z is In Perp[S] and Sprojection.X is in S and everything in Perp[S] is perpendicular to everything in S.

This tells you that if you take any member of Perp[S], then (X - Y) is perpendicular to it. In other words,
 X - Y is perpendicular to everything in Perp[S].  
 Done.

T.3.b.i)  Using the (S^⊥) idea to extend a given basis for a subspace  S of  kD to a basis for all of kD

Here's a basis (linearly independent spanning set) for a subspace S of 5D

Clear[s, i, spanner];
spanner[1]=Table[Random[Real, {-2, 2}], {j, 1, 5}];
spanner[2]=Table[Random[Real, {-2, 2}], {j, 1, 5}];
spanner[3]=Table[Random[Real, {-2, 2}], {j, 1, 5}];
spanners=Table[spanner[i],{i,1,3}]

Check to see that these spanners are a basis of S:

SpannerMatrix  =  Transpose[spanners];
rank  =  
      Length[SingularValues[SpannerMatrix][[2]]]

Come up with two 5D vectors X[1] and X[2] so that 
when you put 
            spanner[4] = X[1] and spanner[5] = X[2], 
then
          {spanner[1], spanner[2], spanner[3], spanner[4], spanner[5]} 
is a basis (linearly independent spanning set) for all of 5D.

Make the Sprojection matrix:

SpannerMatrix  =  Transpose[spanners];
Sprojection  =  
        SpannerMatrix . 
          PseudoInverse[SpannerMatrix];
MatrixForm[Sprojection]

You can extend the original spanning set of S by throwing in a basis of (S^ ⊥) .

Because  (S^ ⊥)  is all hits with Identity - Sprojection, you can get a basis for (S^ ⊥)  by looking at stretch factors and the hanger frame for  Identity - Sprojection:

A = IdentityMatrix[5]  -  Sprojection;
stretches = SingularValues[A][[2]]

Clear[hangerframe];
{hangerframe[1], hangerframe[2]} =SingularValues[A][[1]]

Put  spanner[4] = hangerframe[1] and spanner[5] = hangerframe[2] to get your basis of all of 5D:

spanner[4] = hangerframe[1];
spanner[5] = hangerframe[2];
extendedspanners =Table[spanner[k], {k, 1, 5}]

Check to see that {spanner[1],spanner[2],spanner[3],spanner[4],spanner[5]} spans all five dimensions of 5D:

extendedSpannerMatrix  =  
        Transpose[extendedspanners];
rank = Length[SingularValues[extendedSpannerMatrix][[2]]]

Yep, {spanner[1],spanner[2],spanner[3],spanner[4],spanner[5]} spans all five dimensions of 5D.
And you're out of here.

T.4)  The null space of a matrix A is the subspace of hitdimD of vectors hit onto {0,0,,...,0} by A.
The SVD aligner frame for A spans a subspace S  of hitdimD. 
The null space of A is the perpendicular complement (S^ ⊥) of this subspace. 
This fact tells you how to come up with an orthonormal basis for the null space of A

T.4.a.i)  The null space of A is the subspace of hitdimD of vectors hit onto {0,0,,...,0} by A

Here's a matrix:

Clear[a];
a[i_, j_]  :=  Random[Real, {-5, 5}]
A  = Table[a[i, j], {i, 1, 3}, {j, 1, 7}];
MatrixForm[A]

When you ask Mathematica for the null space of A, you get:

NullSpace[A]

What are these vectors?

Look again:

Clear[spanner];
{spanner[1], spanner[2], spanner[3], spanner[4]}  = 
      NullSpace[A]

Mathematica delivers a basis for the subspace of 7D consisting of the Y's in 7D that are squashed to  {0,0,. . .,0} by a hit with A.  
Try this out for cleared coefficients x[j]:

Expand[A . ((∑_(j = 1)^4) x[j] spanner[j])]

This subspace of 7D is called the Null Space of A.

T.4.a.ii) The SVD aligner frame for A spans a subspace S  of hitdimD. The null space of A is the perpendicular complement  (S^ ⊥) of this subspace. This fact tells you how to come up with an orthonormal basis for the null space of A

Here's a new matrix A:

hangdim  =  4;
hitdim  =  7;
Clear[a];
a[i_, j_] := Random[Real, {-5, 5}]
A = 
        Table[a[i, j], {i, 1, hangdim}, {j, 1, 
            hitdim}];
MatrixForm[A]

Use SVD and perpendicular complements to come up with an orthonormal basis for the null space of A.

Calculate the SVD alignerframe consisting of vectors with corresponding positive stretch factors:

rankA =  
        Length[SingularValues[A][[2]]];
Clear[alignerframe];
alignerframe[
          k_]  :=  SingularValues[A][[3, 
            k]];
Sspanners  = Table[alignerframe[k], {k, 1, rankA}]

The set {alignerframe[1],alignerframe[2],. . ., alignerframe[rankA]} spans a subspace S of hitdimD.
The truth is that the null space of A is the perpendicular complement (S^ ⊥)  of S.

On way of explaining why the null space of A is the perpendicular complement (S^ ⊥)  of S is to explain  the following two facts:

Fact 1  Any X in (S^ ⊥)  is also in the null space of A.

Go with an X in (S^ ⊥)  and look at

A.X = (∑_(k = 1)^rankA)stretch[k] (X . alignerframe[k]) hangerframe[k]       ( with stretch[k] > 0)

= (∑_(k = 1)^rankA)stretch[k] ( 0 )  hangerframe[
            k]

= (∑_(k = 1)^rankA){0, 0, . . ., 0}

=  {0 ,0, . . .,0}.

This explains why any X in (S^ ⊥)  is also in the null space of A.

Fact 2  Any X in the null space of A is also in (S^ ⊥) .

Go with an X in the null space of A. 
This gives:

{0,0,. . .,0} = A.X = (∑_(k = 1)^rankA)stretch[k] (X . alignerframe[k]) hangerframe[k].
So

(∑_(k = 1)^rankA)stretch[k] (X . alignerframe[k]) hangerframe[k] = {0,0,. . .,0}
Take a p with 1≤  p ≤ rankA and dot both sides with hangerframe[p] to get

(∑_(k = 1)^rankA)stretch[k] (X . alignerframe[k]) hangerframe[p].hangerframe[k] = hangerframe[p].{0,0,. . .,0} = 0.

This gives:
                stretch[p] (X.alignerframe[p])  = 0.

The result:
                (X.alignerframe[p])  = 0.

And because p is any integer with 1≤ p ≤ rankA, this tells you that X is perpendicular to each of the spanning vectors 
           {alignerframe[1],alignerframe[2],. . .,alignerframe[rankA]} of S. 
Consequently X is perpendicular to everything in S. 
This signals that X is in (S^ ⊥) .
This explains why any X in the null space of A is also in (S^ ⊥)

To go after an orthonormal basis for S (perpendicular frame spanning S) remember that (S^ ⊥)  consists of all hits with Identity - Sprojection:

SpannerMatrix  =  
        Transpose[Sspanners];
SProjection  =  
        SpannerMatrix . 
          PseudoInverse[SpannerMatrix];
B  =  IdentityMatrix[hitdim]  - SProjection; 
MatrixForm[
      B ]

So you can get an orthonormal basis for the null space of A  by going with the SVD hangerframe of 
                   B = Identity - Sprojection.
Here you go:

rankB =  Length[SingularValues[B][[2]]];
Clear[hangerframe];
hangerframe[k_]  := SingularValues[B][[1, k]];
nullspacespanners  = Table[hangerframe[k], {k, 1, rankB}]

There you go.
Try it out:

Clear[x, j];
Expand[
      A . (∑_(k = 1)^rankB) x[k] hangerframe[k]]

T.5) The sum  S1 + S2 of two subspaces of kD is a subspace of kD.
The intersection  S1⋂ S2 of two subspaces of kD is a subspace of kD.
How to find orthonormal bases for them.
No matter what Y in kD you go with, you are guaranteed that ||Sprojection.Y|| ≤ ||Y||.
Saying ||Sprojection.Y|| = ||Y|| for a given Y in kD  is the same as saying Y is in S.
Saying ||Sprojection.Y|| < ||Y||  for a given Y in kD  is the same as saying Y is  not in S.
Saying that a given Y in kD is in S1⋂ S2 is the same as saying
                            ||S2Projection.S1Projection.Y|| = ||Y||

T.5.a.i)  Adding two subspaces of kD

Given two subspaces S1 and S2 of kD, most folks calculate
                  S1 + S2 
by taking everything in S1 and adding this to everything in S2.
Confirm that  S1 + S2 is a subspace of kD by explaining how to get a spanning set for 
                 S1 + S2 
from a given spanning sets S1 and a given spanning set for S2.

To do this you've got to start with a spanning sets for S1 and S2 and exhibit a spanning set for S1 + S2.

 The subspace S1 is specified by a spanning set

spanners1 ={spanner1[1],spanner1[2],. . .,spanner1[p]}.

The members of S1 look like this:

(∑_(i = 1)^p)x[i] spanner1[i]

where the x[i]'s are real numbers.

The subspace S1 is specified by a spanning set

spanners2 ={spanner2[1],spanner2[2],. . .,spanner2[q]}.

The members of S2 look like this:

(∑_(j = 1)^q)y[j] spanner2[j]

where the y[j]'s are real numbers.

You get S1 + S2 by taking everything in S1 and adding this to everything in S2.

In other words S1 + S2 is everything of this form:
         (∑_(i = 1)^p)x[i] spanner1[i]+ (∑_(j = 1)^q)y[j] spanner2[j].
This tells you that a spanning set for  S1 + S2  is the spanning set you get by throwing 
spanners1 and spanners2 into one big spanning set. 
In fancier terms, a spanning set for S1 + S2 is:
            spanners1 ⋃ spanners2

T.5.a.ii)  Calculating the dimension of S1 + S2 and getting an orthonormal basis for S1 + S2

A subspace S1 of 8D is specified by the following spanning set:

Clear[i, 
        spanner1];
spanner1[1] = {0.0, -0.3, -0.2, 0.2, 
          0.4, 0.4, 0.5, 0.9};
spanner1[2] = {0.9, 0.8, 0.7, 0.4, 
          0.0, -0.1, 0.0, -0.7};
spanner1[
          3] = {-0.5, -0.3, 0.9, 0.4, -0.3, 0.2, -0.5, 
          0.4};
spanner1[4] = {1.4, 0.8, -0.4, 0.2, 0.7, 0.1, 
          1. , -0.2};
spanners1=Table[spanner1[i], {i, 1, 4}]

Another subspace S2 of 8D  is specified by the following spanning set:

Clear[j, 
        spanner2];
spanner2[1] = {0.2, -0.5, -0.3, 0.4, 
          0.4, -0.8, 0, 0};
spanner2[2] = {0.4, 0.2, 1.4, 1. , 
          0.1, 0.5, 0, 0.6};
spanners2=Table[spanner2[j], {j, 1, 2}]

Calculate:
▪ the dimension of S1,
▪ the dimension of S2 and 
▪ the dimension of S1 + S2.
Exhibit an orthonormal basis of S1 + S2

The dimension of S1 is:

SpannerMatrix1=Transpose[spanners1];
rank=Length[SingularValues[SpannerMatrix1][[2]]];
S1dim = rank

The dimension of S2 is:

SpannerMatrix2 = Transpose[spanners2];
rank=Length[SingularValues[SpannerMatrix2][[2]]];
S2dim = rank

To calculate the dimension of S1 + S2, you have to come up with a spanning set. According to part i) immediately above, you just take spanners1 and throw in spanners2 to get 
        spanners1 ⋃ spanners2 
as a spanning set for S1 + S2:

S1plusS2spanners =spanners1⋃spanners2

Now you can make the spanner matrix for S1 + S2 and use it to calculate the dimension of S1 + S2:

SpannerMatrix12 = Transpose[S1plusS2spanners ];
rank = Length[SingularValues[SpannerMatrix12][[2]]];
S1plusS2dim  =  rank

The dimension of S1 + S2  is 4.
An orthonormal basis for S1 + S2 comes from the hanger frame that SpannerMatrix12 uses to hang its hits:

hangerframe = SingularValues[SpannerMatrix12][[1]]

There you go.

T.5.b.i)  Using the S1 Projection and the S2 Projection matrices to come up with an orthonormal basis of S1⋂ S2  and to get the  S1⋂ S2 Projection matrix

A subspace S1 of 7D is specified by the following spanning set:

Clear[i, 
        spanner1];
spanner1[1] = {0.0, -0.3, -0.2, 0.2, 
          0.4, 0.4, 0.5};
spanner1[2] = {0.9, 0.8, 0.7, 0.4, 
          0.0, -0.1, 0.0};
spanner1[3] = {-0.5, -0.3, 
          0.9, 0.4, -0.3, 0.2, -0.5};
spanner1[4] = {1.4, 
          0.8, -0.4, 0.2, 0.7, 0.1, 1. };
spanners1=Table[spanner1[i], {i, 1, 4}]

Another subspace S2 of 7D  is specified by the following spanning set:

Clear[j, 
        spanner2];
spanner2[1] = {0.2, -0.5, -0.3, 0.4, 
          0.4, -0.8, 0};
spanner2[2] = {0.4, 0.2, 1.4, 1. , 
          0.1, 0.5, 0};
spanner2[
          3]  =  {-1.4, -1.1, 0.2, 0, -0.3, 
          0.3, -0.5};
spanners2=Table[spanner2[j], {j, 1, 3}]

Come up with an orthonormal basis for S1 ⋂ S2 and use the result to calculate the dimension of S1 ⋂ S2.
Come up with the  S1⋂ S2Projection matrix.

Calculate the S1projection matrix:

SpannerMatrix1=Transpose[spanners1];
S1projection = SpannerMatrix1.PseudoInverse[SpannerMatrix1];
MatrixForm[S1projection]

Calculate the S2projection matrix:

SpannerMatrix2=Transpose[spanners2];
S2projection = SpannerMatrix2.PseudoInverse[SpannerMatrix2];
MatrixForm[S2projection]

Multiply out S1projection.S2projection

product =S2projection.S1projection;
MatrixForm[product]

This is not necessarily the projection matrix for S1 ⋂ S2, but you can extract the  projection matrix for S1 ⋂ S2 from S2projection.S1projection.
Here's how:

Look at the stretch factors of S2projection.S1projection:

stretches = SingularValues[product][[2]]

Two of them equal 1. This is your signal that the dimension of S1 ⋂ S2 is 2. 
To get an orthonormal basis for  S1 ⋂ S2, you go with the first two alignerframe vectors for
               product  = S2projection.S1projection:

Clear[alignerframe];
alignerframe[k_] := SingularValues[S2projection.S1projection][[3,k]] ;
orthonormal = Table[alignerframe[k],{k,1,2}]

The projection matrix  for S1⋂ S2 is:

aligner  =  {alignerframe[1], alignerframe[2]};
hanger  =  
        Transpose[ aligner];
S1intersectS2Projection = hanger . aligner;
MatrixForm[S1intersectS2Projection]

T.5.b.ii) Why did that work?

In the last problem,you calculated an  SVD aligner frame for the product  S2Projection.S1Projection.
And then you got your orthonormal basis for S1⋂ S2 by taking all the aligner frame vectors with corresponding SVD stretch factors equal to 1.
Why did that do the job?
Why will the same thing work in other situations?

Go with subspaces S1 and S2 of kD and a given kD vector Y.
It so happens that saying that Y is in S1⋂ S2 is the same as saying that
          ||S2Projection.S1Projection.Y|| = ||Y||.

That's why you go with the aligner frame vectors for the product  S2Projection.S1Projection that have corresponding stretch factors equal to 1.

T.5.b.iii) Explanation of the logical sequence that says for subspaces S. S1 and S2 of kD: 
            1 No matter what Y in kD you go with, you are guaranteed that ||Sprojection.Y|| ≤ ||Y||.
             2  Saying ||Sprojection.Y|| = ||Y|| for a given Y in kD  is the same as saying Y is in S.
             3  Saying ||Sprojection.Y|| < ||Y||  for a given Y in kD  is the same as saying Y is  not in S.
      ↦ 4  Saying that a given Y in kD is in S1⋂ S2 is the same as saying
                            ||S2Projection.S1Projection.Y|| = ||Y||

Here is the logical sequence for an explanation of why when you are given two subspaces S1 and S2 of kD, then saying that a given Y in kD is in S1⋂ S2 is the same as saying
                            ||S2Projection.S1Projection.Y|| = ||Y||.

For subspaces S, S1 and S2 of kD:

1 No matter what Y in kD you go with, you are guaranteed that ||Sprojection.Y|| ≤ ||Y||.

Explanation:
Take Y break it into Y = Y1 + Y2  with Y1 in S and Y2 perpendicular to S. (So that SProjectionY = Y1).
Now note that
                 ||SProjectionY|| = ||Y1|| ≤ \(\@\(\(||\)\(Y1\)\( || \^2\)\(+\(\(||\)\(Y2\)\( || \^2\)\)\)\)\) = ||Y||.

2  Saying ||Sprojection.Y|| = ||Y|| for a given Y in kD  is the same as saying Y is in S.

Explanation:
Take Y break it into Y = Y1 + Y2  with Y1 in S and Y2 perpendicular to S. (So that SProjectionY = Y1).
It's automatic that ||SProjectionY|| = ||Y1||. 

||SProjectionY|| = ||Y1|| = ||Y|| = \(\@\(\(||\)\(Y1\)\( || \^2\)\(+\(\(||\)\(Y2\)\( || \^2\)\)\)\)\)  ⟺ Y2 = {0,0,...,0} ⟺ Y is in S.

3 Saying ||Sprojection.Y|| < ||Y|| for a given Y in kD is the same as saying Y is not in S.

According to 1 , you know that ||Sprojection.Y|| ≤ ||Y||. 
2 tells you that saying  ||Sprojection.Y|| = ||Y|| for a given Y in kD  is the same as saying Y is in S.
So saying ||Sprojection.Y|| < ||Y|| is the same as saying Y is not in S.

↦ 4  Given two subspaces S1 and S2 of kD, then saying that a given Y in kD is in S1⋂ S2 is the same as saying
                            ||S2Projection.S1Projection.Y|| = ||Y||.

First half of the explanation:
 Assume Y is in S1 ⋂ S2 and explain why ||S2Projection.S1Projection.Y||  =  ||Y||.
The assumption that Y is in S1⋂ S2,  guarantees that S1Projection.Y = Y and  S2Projection.Y = Y.   
So if Y is in both S1⋂ S2, then
              S2Projection.S1Projection.Y  = S2Projection.Y = Y. 
And so
               ||S2Projection.S1Projection.Y||  =  ||Y||. 
              
Second half of the explanation::
Assume ||S2Projection.S1Projection.Y||  =  ||Y|| and explain why Y is in S1⋂ S2.

According to  3, if Y is not in S1, then ||S1projection.Y|| < ||Y||.
But according to 1  
      ||Y|| = ||S2Projection.S1Projection.Y||  ≤ ||S1Projection.Y||  < || Y ||.  
The upshot:
The assumption that Y is not in S1 leads to the contradiction  ||Y|| <  ||Y||.
This explains why  Y is in S1. 

The fact that Y is in S1 guarantees that S1Projection.Y = Y.
So according to  1 ,
                ||Y||  = ||S2Projection.S1Projection.Y||  = ||S2Projection.Y||.
This gives 
                ||S2Projection.Y|| =  ||Y||.
Now 2  steps in to tell you that Y is also in S2. 


The conclusion:
Saying ||S2Projection.S1Projection.Y||  =  ||Y|| is the same as saying  Y is in S1⋂ S2.
There you go.

T.6) Getting an orthonormal basis from a given basis via the Gram-Schmidt process

T.6.a.i) Using the GramSchmidt instruction

A subspace S of 6D is defined by the following given spanning set:

Clear[i, spanner];
spanner[1] = {0.43, -0.57, 2.19, -0.94, 1.28};
spanner[2] = {0.97, -3.29, -0.83, 5.23, 1.26};
spanner[3] = {-6.14, 3.08, 
          0.00, -1.97, -0.96};
spanner[4] = {-3.68, 
          1.22, -3.64, -2.14, -0.94};
spanners = Table[spanner[i], {i, 1, 4}]

Check the dimension of S:

SpannerMatrix = Transpose[spanners];

Sdim=Length[SingularValues[SpannerMatrix][[2]]]

This tells you that the dimension of S is the same as the number of vectors in the spanning set. 
And this signals that the given spanning set for S is also a basis (linearly independent spanning set) for S.
Old matrix hands know that you can get an orthonormal basis (perpendicular frame spanning S) for S by applying what lots of folks like to call the "Gram-Schmidt Orthogonalization Process:"

Needs["LinearAlgebra`Orthogonalization`"];
\(orthospanners = \ GramSchmidt[\ spanners]
\)

What are these vectors?

Take another look:

orthospanners =  GramSchmidt[ spanners]

These four vectors give you an orthonormal basis for the subspace S spanned by the original basis. They were produced through a process widely known as the Gram-Schmidt process. 

If you want to know what the Gram-Schmidt process is, go on.

T.6.a.ii) How the Gram-Schmidt process works

How does the Gram-Schmidt process take in a basis (linearly independent spanning set)
        {spanner[1], spanner[2], spanner[3], . . ., 
        spanner[k]}
for a subspace S of kD and spit out an orthonormal basis (perpendicular frame) for the same subspace S?

Here's how it goes:
Go with the given basis   
      {spanner[1], spanner[2], spanner[3], . . ., 
        spanner[k]}
for  a subspace S of nD.

Gram-Schmidt  Step1:

Put Y[1] = spanner[1]. 
and put
     X[1] = Y[1]/(Sqrt[Y[1] . Y[1]]).

Gram-Schmidt   Step2:

Put 
      Y[2] = 
        spanner[2] - (component of spanner[2] in the direction of X[
                1]).
 Put 
      X[2] = Y[2]/(Sqrt[Y[2] . Y[2]]) .
Note:
       X[1] and X[2] are unit vectors,
      X[1] . X[2] = 0 ,
 and 
    {X[1], X[2]} spans the subspace spanned by {spanner[1], spanner[2]}

Gram-Schmidt  Step3:

Put

\(Y[3] = \(\(spanner[3]\)\(
\)
        \(		\)\(-\ \((component\ of\ spanner[3]\ in\ the\ direction\ of\ X[
                1])\)\)\(\ \)\)\).

- (component of spanner[3] in the direction of X[
              2]).

and put 
      X[3] = Y[3]/(Sqrt[Y[3] . Y[3]]).
Note:Note:
       X[1], X[2] and X[3] are all unit vectors,
     
       X[1] . X[2] = 0 
       X[1] . X[3] = 0 
       X[2] . X[3] = 0 
 and 
{X[1], X[2], X[3]} spans the subspace spanned by {spanner[1], spanner[2], spanner[3]}.

Gram-Schmidt  Step j:

Put 
      Y[j] = 
        spanner[j] - (∑_(i = 1)^(j - 1))component of spanner[
                j] in the direction of X[i].

      X[j] = Y[j]/(Sqrt[Y[j] . Y[j]]) .
 Note
      X[1], X[2],. . . , X[j] are all unit vectors,
      
      X[p] . X[i] = 0 for i,p = 1,2,...,j-1  with i ≠ p.
and
     {X[1], X[2], X[3], . . ., X[j]} 
 spans the subspace spanned by 
     {spanner[1], spanner[2], spanner[3],  ... , spanner[j]}

Starting with
 given spanning set   
      {spanner[1], spanner[2], spanner[3], . . ., 
        spanner[k]},
the process stops when  j = k and produces vectors 
         {X[1], X[2],  ... , X[k]} 
with the guarantees that:
1) Each X[i] is a unit vector .
2) And for each stage j, 
       {X[1], X[2],  ... , X[j]} 
spans the same subspace as
       {spanner[1], spanner[2], . . ., 
          spanner[j]}.
 
 This tells you that
         {X[1], X[2],  ... , X[k]}  
  spans the same subspace as
         {spanner[1], spanner[2], spanner[3], . . ., 
        spanner[k]} 
So {X[1], X[2], . . ., X[k]}  is guaranteed to be an
 orthonormal basis of S (a perpendicular frame that frames up S).

T.6.a.iii) Trying it out

Here is the linearly independent spanning set from part i)

Clear[i, spanner];
spanner[1] = {0.43, -0.57, 2.19, -0.94, 1.28};
spanner[2] = {0.97, -3.29, -0.83, 5.23, 1.26};
spanner[3] = {-6.14, 3.08, 
          0.00, -1.97, -0.96};
spanner[4] = {-3.68, 
          1.22, -3.64, -2.14, -0.94};
spanners = Table[spanner[i], {i, 1, 4}]

An orthonormal basis for the same subspace is:

orthospanners =  GramSchmidt[ spanners]

Illustrate how the  Gram-Schmidt process produces this orthonormal basis for S:

Copy, paste and edit:
Go with the given basis   
      {spanner[1], spanner[2], spanner[3], sapnner[4]}
for  a subspace S of 6D.

Gram-Schmidt Step1:

Put Y[1] = spanner[1]. 
and put
     X[1] = Y[1]/(Sqrt[Y[1] . Y[1]]):

Clear[X, Y];
Y[1] = spanner[1];
X[1] = Y[1]/(Sqrt[Y[1] . Y[1]])

Gram-Schmidt  Step2:

Put 
      Y[2] = 
        spanner[2] - (component of spanner[2] in the direction of X[
                1]).
 Put 
      X[2] = Y[2]/(Sqrt[Y[2] . Y[2]])

Y[2]=spanner[2]-(spanner[2] . X[1]) X[1];
X[2] = Y[2]/(Sqrt[Y[2] . Y[2]])

Gram-Schmidt  Step3:

Put

\(Y[3] = \(\(spanner[3]\)\(
\)
        \(		\)\(-\ \((component\ of\ spanner[3]\ in\ the\ direction\ of\ X[
                1])\)\)\(\ \)\)\).

- (component of spanner[3] in the direction of X[
              2]).

and put 
      X[3] = Y[3]/(Sqrt[Y[3] . Y[3]]).

Y[3]=spanner[3]-(spanner[3] . X[1] )X[1]- (spanner[3] . X[2]) X[2];
X[3] = Y[3]/(Sqrt[Y[3] . Y[3]])

Gram-Schmidt Step 4:

Put

\(Y[4] = \(\(spanner[4]\)\(
\)
        \(		\)\(-\ \((component\ of\ spanner[4]\ in\ the\ direction\ of\ X[
                1])\)\)\(\ \)\)\).

- (component of spanner[4] in the direction of X[
              2]).

- (component of spanner[4] in the direction of X[
              3]).

and put 
      X[4] = Y[4]/(Sqrt[Y[4] . Y[4]]).

Clear[i];
Y[4]=spanner[4]-(∑_(i = 1)^3)(spanner[4] . X[i])  X[i];
X[4] = Y[4]/(Sqrt[Y[4] . Y[4]])

Compare:

{X[1], X[2], X[3], X[4]}

orthospanners =  GramSchmidt[ spanners]

It will work this way every time.

T.6.b.i) SVD versus Gram-Schmidt: They both give you an orthonormal basis. 
Usually they don't give the same basis

Given a linearly independent spanning set for a subspace S of nD,
you have your choice of an orthonormal basis S via 
1) By going with the GramSchmidt process (as above) 
or 
2) By going with the SVD hanger frame for the associated SpannerMatrix (as in the Basics)
Does it matter which you use?.

Try out the two on this spanning set for a subspace S of 6D:

Clear[i, spanner];
spanner[1] = {0.43, -0.57, 2.19, -0.94, 1.28, 
          1.00};
spanner[2] = {0.97, -3.29, -0.83, 5.23, 1.26, 
          0.98};
spanner[3] = {-6.14, 3.08, 0.00, -1.97, -0.96, 
          0.75};
spanner[4] = {-3.68, 
          1.22, -3.64, -2.14, -0.94, -1.26};
spanner[5] = {1.23, 1.78, 6.90, -3.77, -0.86, 
          2.53};
spanners = Table[spanner[i], {i, 1, 5}]

SpannerMatrix  =  
        Transpose[spanners];
SVDorthonormal =SingularValues[SpannerMatrix ][[1]]

GramSchmidtorthonormal= GramSchmidt[ spanners]

The results are different but each gives a fully legitimate orthonormal basis for S.

If you are away from the computer, you might prefer Gram-Schmidt because, in principle, it can be done by hand.

T.6.b.ii) Down side of Gram Schmidt

What's the down side of Gram-Schmidt?

The down side of Gram-Schmidt is that before you use it, you have to know that the given spanning set is linearly independent (is a basis). There is no such problem when you go with SVD.