%% creates toy data for the Variational Garrote
%           http://arxiv.org/abs/1109.0486
%
% p         % training data size
% n         % number of features
% pv        % validation data size
% ns        % non-zero entries (ns<n)
% betah     % output noise (inverse variance)
% correl    % use correlated inputs?
% rho       % correlation level (if correl==true)
% ic        % percent of inversely correlated input
%
% Copyright (C) 2011-, Vicenç Gomez, Bert Kappen
% Last modified: Nov 1, 2013 Vicenç Gomez

p = 100;            % training data size
n = 400;            % number of features
pv = 50;            % validation data size
ns = 10;            % non-zero entries (ns<n)
betah = 2;          % output noise (inverse variance)
correl = 1;         % use correlated inputs?
rho = .8;           % correlation level (if correl==true)
ic = .2;            % percent of inversely correlated input
icn = round(ic*n);

%% ground truth;
w = sparse(1,n);
snonzero = randperm(n);
snonzero(ns+1:end) = [];
w(snonzero) = sign(randn(1,ns));%1:ns;

sigma=sqrt(1/betah);
if correl
    cova = rho.^abs((1:n)'*ones(1,n)-ones(n,1)*(1:n));
else
    cova = eye(n);
end;

%% training set
x = mult_var_gauss(p, cova, correl);

% choose ic% of entries and reverse them
pn=randperm(n);
for pni=1:icn
    pred = pn(pni);
    [xd, idx] = sort(x(pred,:), 'descend');
    x(pred, idx) = sort(x(pred,:), 'ascend');
end
y = w*x+sigma*randn(1,p);
y = y-mean(y);

%% validation set
xv = mult_var_gauss(pv,cova,correl);
yv = w*xv+sigma*randn(1,pv);
yv = yv-mean(yv);
