[R] sas to r

[Frank E Harrell Jr]

Here is an improvement to what I sent earlier today.

g <- function(y) {
   s <- apply(y, 2,
              function(z) {
                z <- z[!is.na(z)]
                n <- length(z)
                if(n==0) c(NA,NA,NA,0) else
                if(n==1) c(z, NA,NA,1) else {
                  m <- mean(z)
                  s <- sd(z)
                  c(Mean=m, CV=s/m, SE=s/sqrt(n), N=n)
                }
              })
   w <- as.vector(s)
   names(w) <-  as.vector(outer(rownames(s), colnames(s), paste, sep=''))
   w
}
library(Hmisc)
dat2 <-  with(dat1,
               summarize(cbind(BdA, sla),
                         llist(wshed, site, species),
                         g,
                         subset=species %in% c('b','c','p'),
                         stat.name=NULL)
               )

options(digits=3)
dat2  # is a data frame

    wshed site species MeanBdA CVBdA SEBdA NBdA Meansla  CVsla SEsla Nsla
1      1    A       b   100.5 0.133  4.23   10     195 0.1813 11.20   10
2      1    A       c    99.7 0.101  3.17   10     206 0.1024  6.68   10
3      1    A       p   101.4 0.239  7.65   10     188 0.1580  9.39   10
4      1    B       b   109.9 0.118  4.09   10     203 0.1433  9.21   10
5      1    B       c    98.4 0.193  6.01   10     221 0.1250  8.72   10
6      1    B       p   102.9 0.216  7.03   10     203 0.1446  9.29   10
7      2    A       b    95.8 0.241  7.31   10     195 0.2011 12.40   10
8      2    A       c    98.7 0.194  6.04   10     207 0.1274  8.33   10
9      2    A       p   102.2 0.217  7.01   10     191 0.1709 10.31   10
10     2    B       b    97.8 0.235  7.27   10     191 0.2079 12.58   10
11     2    B       c   100.9 0.164  5.24   10     194 0.0987  6.07   10
12     2    B       p   103.0 0.144  4.69   10     209 0.0769  5.35    9

# Another approach, but does casewise deletion of NAs and computes all
# possible marginal summaries
dat2 <- summary(cbind(BdA, sla) ~ wshed + site + species,
              data=dat1, fun=g, method='cross')
print.data.frame(dat2)  # after Sept 2004 just say dat2
# dat2 is similar to but not really a data frame, with a matrix of 
statistics S
# Remove method='cross' to just get all one-way summaries

-- 
Frank E Harrell Jr   Professor and Chair           School of Medicine
                      Department of Biostatistics   Vanderbilt University