/*
******************************************************************
*   (C) Copyright 1999, Peter Lenk. All Rights Reserved.
*  GHBREG1.GSS
*		Generats data for HB Regression Model
*		Y_i 	= X*beta_i + epsilon_i
*		beta_i	= Theta'Z_i + delta_i
*		B		= Z*Theta + Delta
*		epsilon_i is N(0,sigma2*I)
*		delta_i	is N(0,Lambda)
*
*		Common design  matrix X.
*
*****************************************************************
*/
nobs	= 400;		@ Number of subjects 					@
mobs	= 10;		@ Number of observations per subject 	@
sigmat	= 5;		@ True error STD						@

lbd12	= {
			1	.5	-1	.1  0  ,
			0	2	0	3   0   ,
			0	0	3	-1  -.5 ,
			0	0	0	2   0   ,
            0   0   0   0   1  
			};	
lbd12	= lbd12/2;
lambdat	= lbd12'lbd12;
/******************
thetat	= {
			2	-1	-3	4,
			-1	0	-2	3,
			3	2	1	0
			};
***********************/
thetat = { 2 -3 0 4  7 };
rankx	= rows(lambdat);
rankz	= cols(thetat');
xdim	= rankx - 1;
zdim	= rankz - 1;
xdata	= rndn(mobs,xdim);
xnames	= 0 $+ "X" $+ ftocv(seqa(1,1,xdim),2,0);
xmat	= ones(mobs,1)~xdata;
if zdim > 0;
    zdata	= rndn(nobs,zdim);
    znames	= 0 $+ "Z" $+ ftocv(seqa(1,1,zdim),2,0);
    zmat	= ones(nobs,1)~zdata;
else;
    zdata = 0;
    znames = "No Z";
    zmat = ones(nobs,1);
endif;

betat	= zmat*thetat	+ rndn(nobs,rankx)*lbd12;
ydata	= betat*(xmat') + sigmat.*rndn(nobs,mobs);
ynames	= 0 $+ "Y" $+ ftocv(seqa(1,1,mobs),2,0);


create f1 = xdata with ^xnames, 0, 8;
if writer(f1,xdata) /= rows(xdata);
		errorlog "Conversion of XDATA to Gauss File did not work";
endif;
closeall f1;

create f1 = zdata with ^znames, 0, 8;
if writer(f1,zdata) /= rows(zdata);
		errorlog "Conversion of ZDATA to Gauss File did not work";
endif;
closeall f1;

create f1 = ydata with ^ynames, 0, 8;
if writer(f1,ydata) /= rows(ydata);
		errorlog "Conversion of YDATA to Gauss File did not work";
endif;
closeall f1;

save sigmat		= sigmat;
save betat		= betat;
save thetat		= thetat;
save lambdat	= lambdat;

fout = "ydata.xls";
ret = export(ydata,fout,ynames);
fout = "xdata.xls";
ret = export(xdata,fout,xnames);
if zdim > 1;
    fout = "zdata.xls";
    ret = export(zdata,fout,znames);
endif;