function theta_draws=sample_post_commcoeff_varspec_reg(Z,F,Omegas,theta_bar,V_theta,ndraws)


%  Sample from the posterior of a unrestricted multivariate regression,
%  with hetero or homoskedastic innovations, common coefficients and
%  variable-specific regressors.
%  The draw is extracted from the posterior implemented following the strategy
%  described in Gelman, Carlin, Stern, Dunson, Vehtari, and Rubin (2013)
%
%  Specification is:
%
%  Z_t   = F_t * theta     + u_t, u_t multnormal(0,Omega_t)
%  prior for theta is multinormal ( theta_bar, V_theta)
%
%  Z_t is nx1, F_t  is nxm , theta is mx1
%  theta_bar is mx1, V_theta is mxm
%
%  sample line: theta_draws=sample_post_commcoeff_varspec_reg(Z,F,Omegas,theta_bar,V_theta)
%
%  INPUTS:   Z is nxT ,  F is nxmxT  
%  Omegas is nxnxT (heteroskedastic) or nxn  (homoskedastic)
%            theta_bar is mx1, V_theta is mxm   (prior moments)
%  OUTPUTS:  theta_draws is mx1
%            or m x ndraws if there is a sixth argument for the number of draws

% Copyright (c) 2013-2018, Francesco Corsello, All Rights Reserved.

b_plurimae=not(nargin<6);
[n,m,T]=size(F);


% evaluating if heteroskedastic or homoskedastic

ndim_om=numel(size(Omegas));

if ndim_om>2 && ndim_om<4
    Omega_diam = spblkdiag(Omegas);
elseif ndim_om>=2 && ndim_om<3
    Omega_diam =spblkdiag(repmat(Omegas,[1,1,T]) );
else
    error('Matrix Omegas is not 2-d or 3-d matrix, please specify it correctly')
end

%stack in rows and prepare covariance matrix
Z_vstack=reshape(Z,n*T,1);
F_vstack=reshape( permute( F, [ 1 3 2 ] ), [ n*T m ] );

if b_plurimae
    theta_draws=draw_stacked_regression_posterior(Z_vstack,F_vstack,Omega_diam,theta_bar, V_theta,ndraws);
else
    theta_draws=draw_stacked_regression_posterior(Z_vstack,F_vstack,Omega_diam,theta_bar, V_theta);
end

