[R] Add Significance Codes to Data Frame

[Marc Schwartz]

On Jul 14, 2010, at 2:16 PM, darckeen wrote:

> 
> I was hoping that there might be some way to attach significance code like
> the ones from summary.lm to a dataframe.  Anyone know how to do something
> like that.  Here is the function i'd like to add that functionality to:
> 
> 
> add1.coef <- function(model,scope,test="F",p.value=1,order.by.p=FALSE)
> {
> 	num <- length(model$coefficients)
> 	add <- add1(model,scope,test=test)
> 	sub <- subset(add,add$'Pr(F)'<=p.value)
> 	est <- sapply(rownames(sub), function(x) update(model,paste("~ .
> +",x))$coefficients[num+1])
> 	ret <- data.frame(est,sub$'Pr(F)')
> 	
> 	colnames(ret) <- c("Estimate","Pr(F)")
> 	rownames(ret) <- rownames(sub)
> 	ret <- format(ret,digits=1,nsmall=1,scientific=F)
> 	
> 	if (order.by.p) { ret <- ret[order(ret$'Pr(F)'),]}
> 	return(ret)
> }
> 
> fscope <- as.formula(paste("dep.sign.up ~ ", paste(names(lr)[2:10],
> collapse= "+")))
> rslt <- add1.coef(lm(dep.sign.up ~ 1,
> data=lr),fscope,p.value=1,order.by.p=FALSE)
> 
>                   Estimate Pr(F)
> r1.pop.total           0.02  0.09
> r1.pop.household       0.05  0.09
> r1.pop.female          0.03  0.09
> r1.pop.pct.female  14594.39  0.35
> r1.pop.male            0.04  0.08
> r1.pop.pct.male   -13827.51  0.37
> r1.pop.density         0.06  0.09
> r1.age.0.4.pct     14581.65  0.39
> r1.age.5.14.pct     2849.15  0.76



Review the code for print.summary.lm(), in which you will find the use of printCoefmat(), in which you will find the use of symnum(). 

Using lm.D9 from example(lm):

ctl <- c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
trt <- c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
group <- gl(2,10,20, labels=c("Ctl","Trt"))
weight <- c(ctl, trt)
lm.D9 <- lm(weight ~ group)


> printCoefmat(coef(summary(lm.D9)))
            Estimate Std. Error t value  Pr(>|t|)    
(Intercept)  5.03200    0.22022 22.8501 9.547e-15 ***
groupTrt    -0.37100    0.31143 -1.1913     0.249    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 



pv <- coef(summary(lm.D9))[, "Pr(>|t|)"]


> pv
 (Intercept)     groupTrt 
9.547128e-15 2.490232e-01 


Signif <- symnum(pv, corr = FALSE, na = FALSE, 
                 cutpoints = c(0, 0.001, 0.01, 0.05, 0.1, 1), 
                 symbols = c("***", "**", "*", ".", " "))


> Signif
(Intercept)    groupTrt 
        ***             
attr(,"legend")
[1] 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


Note that 'Signif' is:

> str(Signif)
Class 'noquote'  atomic [1:2] ***  
  ..- attr(*, "legend")= chr "0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1"


So you will need to coerce it to a vector before cbind()ing to a data frame:

> as.vector(Signif)
[1] "***" " "  


HTH,

Marc Schwartz