Seminar 6: limma versus lm

shannonerdelyi — Feb 12, 2014, 2:18 PM

#####################################################################
# STAT 540: Seminar 6
# High Volume Linear Models
# Author: Shannon Erdelyi
# Date: 12-02-2014
#####################################################################

# source code with
#   read 'dat' and 'des'
#   function 'miniDF(genes)' to subset the data
dir <- "/Users/shannonerdelyi/Dropbox/UBC/W2014/STAT 540"
source(paste0(dir, "/Seminars/photoRec_Functions.R"))

#####################################################################
# Compare results from limma to lm on probe 1440645_at 
#####################################################################

# get wild type data only for probe 1440645_at
wt <- subset(miniDF(rownames(dat)), gType=="wt" & gene=="1440645_at")
head(wt)
         sidChar sidNum devStage gType  gExp       gene
681799 Sample_20     20      E16    wt 5.449 1440645_at
681800 Sample_21     21      E16    wt 5.417 1440645_at
681801 Sample_22     22      E16    wt 5.514 1440645_at
681802 Sample_23     23      E16    wt 5.353 1440645_at
681806 Sample_24     24       P2    wt 5.740 1440645_at
681807 Sample_25     25       P2    wt 5.797 1440645_at

# linear model
mod <- lm(gExp ~ devStage, wt)
summary(mod)

Call:
lm(formula = gExp ~ devStage, data = wt)

Residuals:
    Min      1Q  Median      3Q     Max 
-0.1425 -0.0833 -0.0258  0.0768  0.3237 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)       5.4333     0.0675   80.48  < 2e-16 ***
devStageP2        0.3990     0.0955    4.18  0.00081 ***
devStageP6        0.1952     0.0955    2.05  0.05881 .  
devStageP10       0.9200     0.0955    9.64  8.1e-08 ***
devStage4_weeks   3.9613     0.0955   41.49  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.135 on 15 degrees of freedom
Multiple R-squared: 0.994,  Adjusted R-squared: 0.992 
F-statistic:  589 on 4 and 15 DF,  p-value: 2.73e-16 
# the parameter estimates are the same 
# the F statistic is larger (589 > 425.4)