Helo family!
I am looking for a code to specify arh(1) covariance structure in the PROC IML simulation procedure. I saw the IML code for specifying the R-side ar(1) covariance structure by Gibbs & Kiernan (2020), but I am struggling to figure out how to capture arh(1) instead and generate my simulation samples. Here is the code below by Gibbs & Kiernan (2020), when the R-side covariance structure is ar(1):
proc iml worksize=100; resid=2; rho=0.8; m=5; /* m = number of time points */ cov=j(m,m,1); do i=1 to m; do j=1 to m; cov=[i, j]=resid*rho**abs(i-j); end;
end; mean={0 0 0}; n=200; call randseed(61345); z=RandNormal(n, mean, cov); create AR1Errors from z; append from z; quit;
Kindly assist with the parameters (code) required for arh(1) structure in a similar code as given above. Alternative approaches (codes) are most welcome.
Thank you for your assistance.
The ARH(1) structure is shown in the PROC MIXED documentation.
The structure is the Hadamard product of an outer-product matrix sigma*sigma` and an AR1 matrix (which is a Toeplitz matrix). I think the following IML program will sample normal errors from an ARH1 covariance matrix:
proc iml;
sigma = {8, 16, 4, 8, 1};
rho=0.5;
dim = nrow(sigma);
/* AR1 correlation structure:
https://blogs.sas.com/content/iml/2012/11/05/constructing-common-covariance-structures.html
*/
start AR1Corr(dim, rho);
u = cuprod(j(1,dim-1,rho)); /* cumulative product */
return( toeplitz(1 || u) );
finish;
/* ARH(1) is Hadamard product of outer product sigma*sigma`
and the AR(1) correlation matrix */
OutProd = sigma * sigma`;
AR1 = AR1Corr(dim, rho);
ARH1Cov = OutProd # AR1;
print OutProd, AR1, ARH1Cov;
mean = j(1, dim, 0);
n = 200;
call randseed(61345);
z = RandNormal(n, mean, ARH1Cov);
create ARH1Errors from z[c=('x1':'x5')];
append from z;
quit;
/* visualize */
proc sgscatter data=ARH1Errors;
matrix x1-x5;
run;
The ARH(1) structure is shown in the PROC MIXED documentation.
The structure is the Hadamard product of an outer-product matrix sigma*sigma` and an AR1 matrix (which is a Toeplitz matrix). I think the following IML program will sample normal errors from an ARH1 covariance matrix:
proc iml;
sigma = {8, 16, 4, 8, 1};
rho=0.5;
dim = nrow(sigma);
/* AR1 correlation structure:
https://blogs.sas.com/content/iml/2012/11/05/constructing-common-covariance-structures.html
*/
start AR1Corr(dim, rho);
u = cuprod(j(1,dim-1,rho)); /* cumulative product */
return( toeplitz(1 || u) );
finish;
/* ARH(1) is Hadamard product of outer product sigma*sigma`
and the AR(1) correlation matrix */
OutProd = sigma * sigma`;
AR1 = AR1Corr(dim, rho);
ARH1Cov = OutProd # AR1;
print OutProd, AR1, ARH1Cov;
mean = j(1, dim, 0);
n = 200;
call randseed(61345);
z = RandNormal(n, mean, ARH1Cov);
create ARH1Errors from z[c=('x1':'x5')];
append from z;
quit;
/* visualize */
proc sgscatter data=ARH1Errors;
matrix x1-x5;
run;
Hi Rick,
Thank you so much. You are very right. I got a chance to go through the suggested code and additional reading material. It is quite clear. I thank you.
For a general discussion about how to generate heterogeneous covariance matrices, see "Construct heterogeneous structured covariance matrices in SAS."
SAS Innovate 2025 is scheduled for May 6-9 in Orlando, FL. Sign up to be first to learn about the agenda and registration!
Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.
Find more tutorials on the SAS Users YouTube channel.