% November 16, 2020
%
% This function illustrates a simple randomized algorithm for estimating
% the norm of a given vector x. It relies on the fact that if g is a
% Gaussian vector of the same dimensions as x, then
%
%    Expectation of norm(x)^2 = Expectation of ip(g,x)^2
%
% where ip(g,x) is the inner product of g and x.

function main_gaussian_embedding

single_estimate
average_estimates
use_a_matrix

return

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function single_estimate

n    = 5;

%%% Set an arbitrary starting vector.
x    = randn(n,1);
x(1) = 2;
x(3) = -3;
fprintf(1,'||x||^2    = %10.5f\n',norm(x)^2)


%%% Compute randomized estimates:
while 1
  fprintf(1,'Estimate   = %10.5f\n',(randn(1,n)*x)^2)
  keyboard
end

return

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function average_estimates

n    = 5;
nexp = 10000;

%%% Set an arbitrary starting vector.
x    = randn(n,1);
x(1) = 2;
x(3) = -3;

%%% Estimate the norm of the vector.
normest = zeros(1,nexp);
for i = 1:nexp
  G          = randn(1,n);
  normest(i) = (G*x)^2;
end

fprintf(1,'||x||^2               = %10.5f\n',norm(x)^2)
fprintf(1,'Mean of all estimates = %10.5f\n',mean(normest))

figure(1)
hold off
histogram(normest,linspace(0,3*norm(x)^2,min([round(nexp/20),30])))
hold on
plot(norm(x)^2,0,'r.','MarkerSize',40)
hold off

keyboard
return

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function use_a_matrix

n    = 5;
d    = 100;
nexp = 10000;

%%% Set an arbitrary starting vector.
x    = randn(n,1);
x(1) = 2;
x(3) = -3;

%%% Estimate its length
normest = zeros(1,nexp);
for i = 1:nexp
  G          = (1/sqrt(d))*randn(d,n);
  normest(i) = norm(G*x)^2;
end

fprintf(1,'||x||^2               = %10.5f\n',norm(x)^2)
fprintf(1,'Mean of all estimates = %10.5f\n',mean(normest))

figure(1)
hold off
histogram(normest,linspace(0,3*norm(x)^2,min([round(nexp/20),30])))
hold on
plot(norm(x)^2,0,'r.','MarkerSize',40)
hold off

keyboard

return