% InfoRetrievalQR.m
% Example from M. Berry, Z. Drmac and E. Jessup
% "Matrices, Vector Spaces, and Information Retrieval"
% SIAM Review, (41)#2, 335-362 (1999).

A0 = [1 0 0 1 0;
      1 0 1 1 1;
      1 0 0 1 0;
      0 0 0 1 0;
      0 1 0 1 1;
      0 0 0 1 0];
%    at=A0(:,4); A0(:,4)=A0(:,5); A0(:,5)=at;
[t,d] = size(A0);
n = sqrt(sum(A0.*A0));
% Column normalized A
A = A0.*repmat(1./n,t,1);

% 1st query vector
q1 = [1; 0; 1; 0; 0; 0];
q1n = q1/norm(q1);  % normalize q1
% cos theta as a row vector
c1 = q1n.'*A,

% 2nd query vector
q2 = [1; 0; 0; 0; 0; 0];
q2n = q2/norm(q2);  % normalize q2
% cos theta as a row vector
c2 = q2n.'*A,

% QR factorization with column pivoting
[Q,R,P] = qr(A),
%[Q,R,P] = qr(A,0),
%[Q,R] = qr(A),
f=find(abs(diag(R))>eps); rk = f(end),
QA = Q(:,1:rk),
RA = R(1:rk,:),
% the reduced cosines should be the same as before
c1r =  q1n.'*QA*RA,
c1r =  q2n.'*QA*RA,

% Use a lower ranking (rank = 3) approximation for A
rkNew = rk - 1;
QAnew = Q(:,1:rkNew),
RAnew = R(1:rkNew,:),
Anew = QAnew*RAnew;
nNew = sqrt(sum(Anew.*Anew));
% Column normalized A
Anew = Anew.*repmat(1./nNew,t,1);

c1New =  q1n.'*Anew,
c1New =  q2n.'*Anew,

% Try using an even lower ranking (rank = 2) approximation for A
rkNew2 = rk - 2;
QAnew2 = Q(:,1:rkNew2),
RAnew2 = R(1:rkNew2,:),
Anew2 = QAnew2*RAnew2;
nNew2 = sqrt(sum(Anew2.*Anew2));
% Column normalized A
Anew2 = Anew2.*repmat(1./nNew2,t,1);

c1New2 =  q1n.'*Anew2,
c1New2 =  q2n.'*Anew2,
