sas - 使用 PROC IML 处理矩阵和样本

Question

我正在尝试从某些分布中抽取随机样本，如下所示：

我的代码运行但数字看起来很奇怪。所以我不确定出了什么问题，也许是一些运营商。元素非常大。我的尝试：

C_hat=(((x`)*x)**(-1))*((x`)*z);
S=((z-x*c_hat)`)*((z-x*c_hat));

*draw sigma;
sigma = shape(RandWishart(1, 513 - 3 - 2,s**(-1)),4,4);
*draw vec(c);
vec_c_hat= colvec(c_hat`); *vectorization of c_hat;
call randseed(4321);
vec_c = RandNormal(1,vec_c_hat,(sigma`)@(((x`)*x)**(-1)));
c = shape(vec_c,4,4);
print c;

Answer 1

由于您没有提供数据或参考，因此很难猜测您的“奇怪”和“极大”数字是否正确。但是，该程序看起来大部分是正确的，因此请检查您的数据。

您的程序的一个小问题是您正在使用 SHAPE 函数将 vec_c 向量重塑为矩阵。您应该使用 SHAPECOL 函数（或转置结果）。

以下程序使用随 SAS 分发的 Sashelp.Cars 数据来初始化 X 和 Z 矩阵。该程序计算一个随机矩阵 C，它接近于数据的逆叉积矩阵。我还添加了一些中间计算和注释。此版本在 Sashelp.Cars 数据上按预期工作：

proc iml;
use sashelp.cars;
read all var {weight wheelbase enginesize horsepower} into X;
read all var {mpg_city mpg_highway} into Z;
close;

*normal equations and covariance;
xpx = x`*x;
invxpx = inv(xpx);
C_hat = invxpx*(x`*z);
r = z-x*c_hat;
S = r`*r;

*draw sigma;
call randseed(4321);
DF = nrow(X)-ncol(X)-2;
W = RandWishart(1, DF, inv(S));  /* 1 x (p*p) row vector */
sigma = shape(W, sqrt(ncol(W))); /* reshape to p x p matrix */
*draw vec(c);
vec_c_hat = colvec(c_hat`);      /* stack columns of c_hat */
vec_c = RandNormal(1, vec_c_hat, sigma@invxpx);
c = shapecol(vec_c, nrow(C_hat), ncol(C_hat));   /* reshape w/ SHAPECOL */
print C_hat, c;