function y = Gsr(X, an, bn, sigmaQ, rlb)
% generate generating the shadow yields as an average of forward shadow rate
% the forward shadow rate is computed using Wu and Xia's method
% the yields are annulized
% INPUT: 
%         X: state variables  Tx3
%         an: Kx1  
%         bn: Kx3  
%         sigmaQ: Kx1 
% OUTPUT:
%         y_hat:  predicted yileds  TxK

T = size(X,1);  
mus = repmat(an', T, 1) + X * bn';  % TxK
z1_temp = (mus-rlb)./ (sigmaQ');    % TxK
z2_temp = rlb + (mus-rlb).*normcdf(z1_temp) + repmat(sigmaQ',T,1).*normpdf(z1_temp);  % TxK
y = z2_temp;





