Innovation representation using gram-schmidt construction MATLAB

Question

A zero mean random vector y consisting of three random variables has a 3x3 covariance matrix given by R. I need to determine the innovation representation of y i.e y= B*E using gram schmidt construction where E is the vector of uncorrelated components and B is a lower triangular matrix.

I have visited numerous pages giving tutorials on how to do this in MATLAB but only when my input is a matrix with independent column vectors. Here, I have been given covariance matrix and I can't seem to make a connection on how to utilize this matrix and implement it on MATLAB.

Following is the code I got from MathWorks site created by Reza Ahmadzadeh:

function v = GramSchmidt(v)


k = size(v,2);


for ii = 1:1:k
    v(:,ii) = v(:,ii) / norm(v(:,ii));
    for jj = ii+1:1:k
        v(:,jj) = v(:,jj) - proj(v(:,ii),v(:,jj));
    end
end

    function w = proj(u,v)
        % This function projects vector v on vector u
        w = (dot(v,u) / dot(u,u)) * u;
    end

end

I don't have such an input matrix. Or is it just that I am not able to understand this code?

Any help would be deeply appreciated. This is the first time I am working on such a project on MATLAB and a little help would really make me understand this concept much better.

Answer 1

I downloaded the function you refer to here and I started messing around with it. Actually, the v argument it requires is nothing but the random vector to which the Gram-Schmidt algorithm is applied, your y. Proof:

% Create sample data and check it's correlation...
y = randi(10,10,3);
y = y - repmat(mean(y),10,1);
corr(y)

% Compute the Gram-Schmidt using a well known formulation...
[Q,R] = GramSchmidt_Standard(y);

% Compute the Gram-Schmidt using the Reza formulation...
UNK = GramSchmidt_Reza(y);

% Q and UNK are identical...
Q
UNK

function [Q,R] = GramSchmidt_Standard(y)
    [m,n] = size(y);
    Q = zeros(m,n);
    R = zeros(n,n);

    for j = 1:n
        v = y(:,j);

        for i = 1:j-1
            R(i,j) = Q(:,i).' * y(:,j);
            v = v - R(i,j) * Q(:,i);
        end

        R(j,j) = norm(v);
        Q(:,j) = v / R(j,j);
    end
end

function v = GramSchmidt_Reza(v)
    function w = proj(u,v)
        w = (dot(v,u) / dot(u,u)) * u;
    end

    k = size(v,2);

    for ii = 1:1:k
        v(:,ii) = v(:,ii) / norm(v(:,ii));

        for jj = ii+1:1:k
            v(:,jj) = v(:,jj) - proj(v(:,ii),v(:,jj));
        end
    end
end

So you can pick the function you prefer and proceed with your computations.