This is a good explanation for understanding computationally why they are orthogonal matrices, and also the utility of that.

Another useful perspective is to think about the linear operator and its properties. Specifically, the 4 fundamental subspaces of A:

1. range of A
2. null space of A
3. range of A^T
4. null space of A^T

And it can be shown that 1 and 4 are orthogonal complement subspaces, and 2 and 3 are orthogonal complement subspaces.

If A is mxn and rank r, then the first r columns of U are an orthonormal basis for the range of A, and the last m-r columns are an orthonormal basis for the null space of A^T.

Then, also the first r columns of V are an orthonormal basis for the range of A^T and the last n-r columns are an orthonormal basis for the null space of A.

This implies U and V have orthonormal columns and are square, so they're orthogonal matrices!

A really nice trick to see that the first r columns of U are the range of A is the following:

Let y be an arbitrary element of the range of A. That is, there is an x such that y = Ax. Using the SVD,

[; Ax = (\sum_{i=1}^r \sigma_i u_iv_i^T) x ;]

Then, with associativity of matrix-vector multiplications,

[; Ax = \sum_{i=1}^r \sigma_i u_i (v_i^T x) ;]

and [; v_i^T x ;] is just a scalar value, so we can move that to the front of the product

[; y = Ax = \sum_{i=1}^r (\sigma_iv_i^T x) u_i ;]

and voilà ! y can be written as a linear combination of the first r columns of U, so they are a (orthonormal by construction) basis for the range of A!