The program in this section consists of four parts:
Read the data into SAS/IML vectors.
Transfer the data to R.
Call R functions to analyze the data.
Transfer the results of the analysis into SAS/IML vectors.
Read the data. The following statements read the Weight
and Height
variables from the Sashelp.Class
data set into SAS/IML vectors with the same names:
proc iml; use Sashelp.Class; read all var {Weight Height}; close Sashelp.Class;
Transfer the data to R. The following statements run the ExportMatrixToR subroutine in order to transfer data from a SAS/IML matrix into an R matrix. The names of the corresponding R vectors that
contain the data are w
and h
.
/* send matrices to R */ call ExportMatrixToR(Weight, "w"); call ExportMatrixToR(Height, "h");
Call R functions to perform some analysis. The SUBMIT statement with the R option is used to send statements to R. Comments in R begin with a hash mark (#, also called a number sign or a pound sign).
submit / R; Model <- lm(w ~ h, na.action="na.exclude") # a ParamEst <- coef(Model) # b Pred <- fitted(Model) Resid <- residuals(Model) endsubmit;
The R program consists of the following steps:
The lm
function computes a linear model of w
as a function of h
. The na.action=
option specifies how the model handles missing values (which in R are represented by NA). In particular, the na.exclude
option specifies that the lm
function should not omit observations with missing values from residual and predicted values. This option makes it easier
to merge the R results with the original data when the data contain missing values.
Various information is retrieved from the linear model and placed into R vectors named ParamEst
, Pred
, and Resid
.
Transfer the data from R. The ImportMatrixFromR subroutine transfers the ParamEst
vector from R into a SAS/IML vector named pe
. This vector is printed by the SAS/IML PRINT statement. The predicted values (Pred
) and residual values (Resid
) can be transferred similarly. The parameter estimates are used to compute the predicted values for a series of hypothetical
heights, as shown in Figure 11.3.
call ImportMatrixFromR(pe, "ParamEst"); print pe[r={"Intercept" "Height"}]; ht = T( do(55, 70, 5) ); A = j(nrow(ht),1,1) || ht; pred_wt = A * pe; print ht pred_wt;
Figure 11.3: Results from an R Analysis
pe | |
---|---|
Intercept | -143.0269 |
Height | 3.8990303 |
ht | pred_wt |
---|---|
55 | 71.419746 |
60 | 90.914898 |
65 | 110.41005 |
70 | 129.9052 |
You cannot directly transfer the contents of the Model
object. Instead, various R functions are used to extract portions of the Model
object, and those simpler pieces are transferred.