%EXPORT_TO_R SAS Macro Code

The SAS Analysis blog post 'A macro calls R in SAS for paneled 3d plotting' influenced my macro coding.   The following macro call:

%EXPORT_TO_R(DATA = YOURDATA) 

exports the SAS data set 'YOURDATA' as .csv and produces the R code for setting the R working directory and reading the file.

R-code written by %EXPORT_TO_R :

#  set R working directory
setwd("C:\\Documents and Settings\\wkuuser\\Desktop\\PROJECTS\\Stats Training")
 
#  get data
dat.from.SAS <- read.csv("fromSAS_delete.CSV", header=T)
 
#  check data dimensions
dim(dat.from.SAS)
names(dat.from.SAS)

Created by Pretty R at inside-R.org

Actual SAS code is below.

An application of this can be found here.

One potential drawback is that the file path for where you are exporting the data and the R code is hard coded in the macro. When I'm typically working in SAS and need this capability, I typically have a single 'development' folder for ad hoc support work in R. Ultimately when I've finished the program, I save the final R code in the specific project directory I'm working in.  So it's not a big deal once the macro is set up for the computer you are using.

There are many times where you may be in the SAS environment and require a specific statistical capability (like spatial regression ) that is in R or perhaps you don't have liscence for in SAS (like decision trees  with Enterprise Miner). Being able to quickly move into the R environment is nice. I've worked with submitting R code within the SAS IML environment, & it works great. However, its another system to start and run (outside of base SAS as I prefer the base SAS text editor over the IML environment). Updates to SAS PROC IML are supposed to deliver this capability within the base SAS environment, but I have not been able to get it to work. (I believe there are some timing issues with the SAS version and compatability with recent releases of R).

SAS code for macro:

/***********************************************************

   *  MACRO:      EXPORT_TO_R() - EXPORTS SAS DATA SET TO CSV

   *                            - SETS R WORKING DIRECTORY

   *                            - GENERATES R CODE TO READ DATA

   *  PARAMETERS: DATA   = SAS DATASET FOR EXPORTING

   *             

   *  NOTES: -YOU MUST PHYSICALLY SET THE FILE PATH TO THE DIRECTORY           

   *         -YOU WANT THE DATA TO BE EXPORTED TO

   *         -BECAUSE R IS CASE SENSITIVE, THIS MACRO IS CASE SENSITIVE

   *         

   ***********************************************************/

%MACRO EXPORT_TO_R(DATA =    );

  PROC EXPORT DATA = &DATA OUTFILE = "C:\DOCUMENTS AND SETTINGS\WKUUSER\DESKTOP\PROJECTS\STATS TRAINING\fromSAS_delete.CSV" REPLACE;

  RUN;

 

  PROC SQL;

    CREATE TABLE _TMP0 (STRING CHAR(110));

    INSERT INTO _TMP0

    SET STRING = '#  set R working directory'

       SET STRING = 'setwd("C:\\Documents and Settings\\wkuuser\\Desktop\\PROJECTS\\Stats Training")'

       SET STRING = '                           '

    SET STRING = '#  get data'

    SET STRING = 'dat.from.SAS <- read.csv("fromSAS_delete.CSV", header=T)'

    SET STRING = '           '

       SET STRING = '#  check data dimensions'

       SET STRING = 'dim(dat.from.SAS)'

       SET STRING = 'names(dat.from.SAS)';

    QUIT;

 

  DATA _NULL_;

    SET _TMP0;

    FILE "C:\DOCUMENTS AND SETTINGS\WKUUSER\DESKTOP\PROJECTS\STATS TRAINING\importFromSAS.TXT";

    PUT STRING;

  RUN;

%MEND;

An Intuitive Approach to ROC Curves (with SAS & R)

I developed the following schematic (with annotations) based on supporting documents (link) from the article cited below. The authors used R for their work. The ROC curve in my schematic was output from PROC LOGISTIC in SAS, the scatterplot with marginal histograms was created in R (code below) using the scored data from PROC LOGISTIC exported using my SAS MACRO %EXPORT_TO_R  (link to SAS macro code)

(click to enlarge)

Reference:

Selection of Target Sites for Mobile DNA Integration in the Human Genome
Berry C, Hannenhalli S, Leipzig J, Bushman FD, 2006 Selection of Target Sites for Mobile DNA Integration in the Human Genome. PLoS Comput Biol 2(11): e157. doi:10.1371/journal.pcbi.0020157

quote "The data were analyzed using the R language and environment for statistical computing and graphics "
R code for plot was adapted from code provided via the addicted to R graph gallery : http://addictedtor.free.fr/graphiques/RGraphGallery.php?graph=78

# *------------------------------------------------------------------
# |                
# | import scored logit data from SAS - code generated by SAS MACRO %EXPORT_TO_R
# |  
# |  
# *-----------------------------------------------------------------
 
 
#  set R working directory
setwd("C:\\Documents and Settings\\wkuuser\\Desktop\\PROJECTS\\Stats Training")
 
#  get data
dat.from.SAS <- read.csv("fromSAS_delete.CSV", header=T)
 
#  check data dimensions
dim(dat.from.SAS)
names(dat.from.SAS)
 
 
# *------------------------------------------------------------------
# |                
# |  scatter plot with marginal histograms
# |  
# |  
# *-----------------------------------------------------------------
 
#
# model predicts P(G) so we want these probabilities for each group
#
 
 
# get p(G) data set for the group that is actually green
 
green <- dat.from.SAS[ dat.from.SAS$class=="G",]
dim(green)
 
# get p(G) data set for group that is actually red
 
red <- dat.from.SAS[ dat.from.SAS$class=="R",]
dim(red)
 
# just look at regular histograms for each group
 
 
hist(green$P_G, main = 'histogram for green')
hist(red$P_G, main = 'histogram for red')
 
# in order to do scatter plots n must be the same for each 
# group, randomly sample n = n(green) from red
 
 
# Total number of red observations to match green
N <- 24 
print(N)
 
 
# Randomly arrange the data and select out N size sample for red
# and test set.
 
dat <- red[sample(1:N),]
red.rs <- dat[1:N,]
dim(red.rs)
 
# does the distribution retain original properties? Yes
hist(red.rs$P_G, main = 'histogram for red sample')
 
 
plot(green$P_G, red.rs$P_G) 
 
 
# *------------------------------------------------------------------
# |                
# |  create the marginal plots
# |  
# |  
# *-----------------------------------------------------------------
 
 
 
def.par <- par(no.readonly = TRUE) # save default, for resetting...
 
 
# define histograms
Ghist <- hist(green$P_G,plot=FALSE)
Rhist <- hist(red.rs$P_G, plot=FALSE)
 
top <- max(c(Ghist$counts, Rhist$counts))
Grange <- c(0,1)
Rrange <- c(0,1)
nf <- layout(matrix(c(2,0,1,3),2,2,byrow=TRUE), c(3,1), c(1,3), TRUE)
 
#layout.show(nf)
 
 
par(mar=c(3,3,1,1))
plot(green$P_G, red.rs$P_G, xlim=Grange, ylim=Rrange, xlab="green", ylab="red")
 
par(mar=c(0,3,1,1))
barplot(Ghist$counts, axes=FALSE, ylim=c(0, top), space=0, main = 'green')
 
par(mar=c(3,0,1,1))
barplot(Rhist$counts, axes=FALSE, xlim=c(0, top), space=0, horiz=TRUE, main = 'red')
 
par(def.par)

Created by Pretty R at inside-R.org

SAS Code for Missing Data Analysis

The following code will read your data set and provide a report for missing values for numeric variables. 

/*-------------------------------------------------------*

 |  SIMULATE MISSING DATA

 *-------------------------------------------------------*/

DATA MYDATA;

INPUT Y VAR1 VAR2 VAR3 VAR4 $ VAR5 $;

CARDS;

1      10     33     4  M  A

0      .      21     3  .  B

1      30     .      2  M  C

1      20     76     1  F  A

0      .      24     3  F  .

0      20     22     .  M  A

1      10     .      2  .  A

1      .      49     2  F  C

0      30     .      2  F  B

1      20     59     2  M  B

1      20     76     1  .  A

0      .      24     3  F  C

0      20     22     .  F  C

1      10     .      2  M  .

1      .      49     2  M  C

0      30     .      2  F  A

;

RUN;

/*-------------------------------------------------------*

 | ANALYZE AND REPORT MISSING CHARACTER VARIABLES

 *-------------------------------------------------------*/

PROC FREQ DATA = mydata;

TABLES _CHARACTER_ / MISSING;

RUN;

/*-------------------------------------------------------*

 | ANALYZE AND REPORT MISSING NUMERIC VARIABLES

 *-------------------------------------------------------*/

PROC MEANS DATA = MYDATA NMISS;

VAR  VAR1 VAR2 VAR3 ; /* ENTER NUMERIC VARIABLES OF INTEREST*/

OUTPUT OUT=T (DROP=_TYPE_ _FREQ_) NMISS=/AUTONAME;

RUN;

PROC TRANSPOSE DATA = T PREFIX=NMISS OUT=S1;

VAR   _NUMERIC_;

RUN;

DATA S2;

SET S1;

PMISS = NMISS1/16*100; /*DENOMINATOR = TOTAL N IN DATA SET*/

RUN;

PROC PRINT DATA = S2;

RUN;

(output below)

*--CHARACTER VARIABLES--*

VAR4	Frequency	Percent	Cumulative Frequency	Cumulative Percent
	3	18.75	3	18.75
F	7	43.75	10	62.50
M	6	37.50	16	100.00

VAR5	Frequency	Percent	Cumulative Frequency	Cumulative Percent
	2	12.50	2	12.50
A	6	37.50	8	50.00
B	3	18.75	11	68.75
C	5	31.25	16	100.00

 *----NUMERIC VARIABLES ---*

Obs	_NAME_	NMISS1	PMISS
1	VAR1_NMiss	5	31.25
2	VAR2_NMiss	5	31.25
3	VAR3_NMiss	2	12.50

Decision Tree Mechanics with R and SAS

(click to enlarge)

(click to enlarge)

(click to enlarge)

(click to enlarge)

(click to enlarge)

(click to enlarge)
(click to enlarge)

# *------------------------------------------------------------------
# | PROGRAM NAME: R_tree_basic 
# | DATE:4/26/11   
# | CREATED BY: Matt Bogard 
# | PROJECT FILE:P:\R  Code References\Data Mining_R              
# *----------------------------------------------------------------
# | PURPOSE: demo of basic decision tree mechanics               
# |
# *------------------------------------------------------------------
 
rm(list=ls()) # get rid of any existing data 
ls() # view open data sets
 
setwd('/Users/wkuuser/Desktop/R Data Sets') # mac 
setwd("P:\\R  Code References\\R Data") # windows
 
library(rpart) # install rpart decision tree library
 
# *------------------------------------------------------------------
# | get data            
# *-----------------------------------------------------------------
 
dat1 <-  read.csv("basicTree.csv", na.strings=c(".", "NA", "", "?"), encoding="UTF-8")
plot( dat1$x2, dat1$x1, col = dat1$class) # plot data space
 
# fit decision tree
 
(r <- rpart(class ~ x1 + x2, data = dat1)) 
 
plot(r)
text(r)
 
library(rattle) # data mining package
drawTreeNodes(r) # for more detailed tree plot supported by rattle
 
# *------------------------------------------------------------------
# |  
# | 
# | chi square test - 1st split 
# |   
# |         
# *-----------------------------------------------------------------
 
# create a categorical for the first cutoff for x2 > 6.5
 
dat1$cutoff <- (ifelse (dat1$x2 >= 6.5, "x2 >=6.5 ", "x2 < 6.5"))
 
# library(MASS) # required for cross tabulation 
 
tab1 <- table(dat1$cutoff,dat1$class) # cross tabulation
print(tab1) 
 
Xsq <-chisq.test(tab1, correct = FALSE)# chi-squared test for independence
print(Xsq)
 
print(Xsq$exp) # print expected values
 
# *------------------------------------------------------------------
# |  chi square test - 1st split - choose an arbitrarily higher and 
# |  lower split value and compare to optimal split chosen by tree            
# *-----------------------------------------------------------------
 
# higher x2 split value
 
dat1$h1 <- (ifelse (dat1$x2  >= 7.5, "x2 >=7.5 ", "x2 < 7.5"))
tab2 <- table(dat1$h1,dat1$class) # cross tabulation
print(tab2) 
 
Xsq <-chisq.test(tab2, correct = FALSE)# chi-squared test for independence
print(Xsq) # chi square value is lower
 
# lower x2 split value
 
dat1$l1 <- (ifelse (dat1$x2  >= 5.5, "x2 >=5.5 ", "x2 < 5.5"))
tab3 <- table(dat1$l1,dat1$class) # cross tabulation
print(tab3) 
 
Xsq <-chisq.test(tab3, correct = FALSE)# chi-squared test for independence
print(Xsq) # chi square value is lower
 
# look at current dat1 data set summary
dim(dat1)
names(dat1)
 
# *------------------------------------------------------------------
# |  
# |
# | chi square test - 2nd split   
# |
# |          
# *-----------------------------------------------------------------
 
# *------------------------------------------------------------------
# | get data            
# *-----------------------------------------------------------------
 
# to get the data in the 2nd split we have to first subset or partition the 
# data where x2 >= 6.5, hence the partition in the 'recursive partitioning' 
# algorithm used by decision trees) 
 
dat2 <- dat1[dat1$x2 >= 6.5,]
dim(dat2) # N = 44 which is correct, recall print(tab1)
 
# *------------------------------------------------------------------
# | create a categorical for the second cutoff for x1 >= 4.5           
# *-----------------------------------------------------------------
 
dat2$cutoff <- (ifelse (dat2$x1 >= 4.5, "x1 >=4.5 ", "x1 < 4.5"))
 
# *------------------------------------------------------------------
# | cross tab & chi square test       
# *-----------------------------------------------------------------
 
tab4 <- table(dat2$cutoff,dat2$class) # cross tabulation
print(tab4) 
 
Xsq <-chisq.test(tab4, correct = FALSE)# chi-squared test for independence
print(Xsq) # X-squared = 44, df = 1, p-value = 3.284e-11
 
# *------------------------------------------------------------------
# |  chi square test - 2nd split - choose an arbitrarily higher and 
# |  lower split value and compare to optimal split chosen by tree            
# *-----------------------------------------------------------------
 
# higher x1 split value
 
dat2$h2 <- (ifelse (dat2$x1 >= 5.5, "x1 >=5.5 ", "x1 < 5.5"))
tab5 <- table(dat2$h2,dat2$class) # cross tabulation
print(tab5) 
 
Xsq <-chisq.test(tab5, correct = FALSE)# chi-squared test for independence
print(Xsq) # X-squared  is lower
 
# lower x1 split value
 
dat2$l2 <- (ifelse (dat2$x1 >= 3.5, "x1 >=3.5 ", "x1 < 3.5"))
tab6 <- table(dat2$l2,dat2$class) # cross tabulation
print(tab6) 
 
Xsq <-chisq.test(tab6, correct = FALSE)# chi-squared test for independence
print(Xsq) # X-squared  is lower
 
# look at data in dat2
dim(dat2)
names(dat2)
 
# note dat2 inherits the variables h1 and l1 (cutoffs for x2) from dat1, but they are not 
# relevant to our analysis in the 2nd split (which focuses on cutoffs for x1)
 
 
# *------------------------------------------------------------------
# |  
# | export data sets to SAS to repeat analysis using a SAS data set,
# | base SAS, and Enterprise Miner
# |   
# |          
# *-----------------------------------------------------------------
 
 
# export data sets to SAS to repeat analysis using a SAS data set,
# base SAS, and Enterprise Miner
 
library(foreign)
# export data and SAS code for analyzing dat1
write.foreign(dat1, "dat1.txt", "dat1_read.sas",   package="SAS") 
# export data and SAS code for analyzing dat2
write.foreign(dat2, "dat2.txt", "basicTreedat2_read.sas",   package="SAS")  
 
 
 
 
 
 

Created by Pretty R at inside-R.org