STAT 540 Seminar 5

Shaun Jackman — Jan 31, 2013, 9:57 AM

# STAT 540 Seminar 5
# Shaun Jackman
# 2013-01-30

# Brief pause to load the photoRec data and the lattice package
library(lattice)  # if you don't already have this loaded ...
prDat <- read.table("../data/photoRec/GSE4051_data.txt")  # the whole enchilada
str(prDat, max.level = 0)
'data.frame':   29949 obs. of  39 variables:
load("../data/photoRec/GSE4051_design.robj")  # load exp design as 'prDes'
str(prDes)
'data.frame':   39 obs. of  3 variables:
 $ sample  : num  20 21 22 23 16 17 6 24 25 26 ...
 $ devStage: Factor w/ 5 levels "E16","P2","P6",..: 1 1 1 1 1 1 1 2 2 2 ...
 $ gType   : Factor w/ 2 levels "wt","NrlKO": 1 1 1 1 2 2 2 1 1 1 ...

# Write a function to prepare a mini-dataset for a small number of genes
luckyGenes <- c("1419655_at", "1438815_at")

library(reshape2)
prepareData <- function(genes)
  melt(cbind(prDes, t(prDat[genes,])),
    id.vars=colnames(prDes),
    variable.name='gene', value.name='gExp')

jDat <- prepareData(luckyGenes)
str(jDat)
'data.frame':   78 obs. of  5 variables:
 $ sample  : num  20 21 22 23 16 17 6 24 25 26 ...
 $ devStage: Factor w/ 5 levels "E16","P2","P6",..: 1 1 1 1 1 1 1 2 2 2 ...
 $ gType   : Factor w/ 2 levels "wt","NrlKO": 1 1 1 1 2 2 2 1 1 1 ...
 $ gene    : Factor w/ 2 levels "1419655_at","1438815_at": 1 1 1 1 1 1 1 1 1 1 ...
 $ gExp    : num  10.93 10.74 10.67 10.68 9.61 ...
head(jDat)
  sample devStage gType       gene   gExp
1     20      E16    wt 1419655_at 10.930
2     21      E16    wt 1419655_at 10.740
3     22      E16    wt 1419655_at 10.670
4     23      E16    wt 1419655_at 10.680
5     16      E16 NrlKO 1419655_at  9.606
6     17      E16 NrlKO 1419655_at 10.840
tail(jDat)
   sample devStage gType       gene  gExp
73     38  4_weeks    wt 1438815_at 8.211
74     39  4_weeks    wt 1438815_at 8.436
75     11  4_weeks NrlKO 1438815_at 8.465
76     12  4_weeks NrlKO 1438815_at 8.841
77      2  4_weeks NrlKO 1438815_at 8.506
78      9  4_weeks NrlKO 1438815_at 8.952

stripplot(gExp ~ devStage | gene, jDat, group = gType,
  jitter.data = TRUE, auto.key = TRUE,
  type = c("p", "a"), grid = TRUE)

plot of chunk unnamed-chunk-1


# Write a function to stripplot a mini-dataset
makeStripplot <- function(x, ...)
  stripplot(gExp ~ devStage | gene, x, group = gType,
            jitter.data = TRUE, auto.key = TRUE,
            type = c("p", "a"), grid = TRUE,
            ...)
makeStripplot(jDat)

plot of chunk unnamed-chunk-1

makeStripplot(jDat, pch = 17, cex = 3)

plot of chunk unnamed-chunk-1

makeStripplot(newDat <- prepareData('1456341_a_at'))

plot of chunk unnamed-chunk-1

str(newDat)
'data.frame':   39 obs. of  5 variables:
 $ sample  : num  20 21 22 23 16 17 6 24 25 26 ...
 $ devStage: Factor w/ 5 levels "E16","P2","P6",..: 1 1 1 1 1 1 1 2 2 2 ...
 $ gType   : Factor w/ 2 levels "wt","NrlKO": 1 1 1 1 2 2 2 1 1 1 ...
 $ gene    : Factor w/ 1 level "1456341_a_at": 1 1 1 1 1 1 1 1 1 1 ...
 $ gExp    : num  7.04 7.48 7.37 6.94 6.16 ...
head(newDat)
  sample devStage gType         gene  gExp
1     20      E16    wt 1456341_a_at 7.044
2     21      E16    wt 1456341_a_at 7.478
3     22      E16    wt 1456341_a_at 7.374
4     23      E16    wt 1456341_a_at 6.944
5     16      E16 NrlKO 1456341_a_at 6.161
6     17      E16 NrlKO 1456341_a_at 6.931

# Do a two-sample t-test
t.test(gExp ~ devStage,
       subset(newDat, devStage %in% c('P2', '4_weeks')))

    Welch Two Sample t-test

data:  gExp by devStage 
t = -18.84, df = 13.98, p-value = 2.477e-11
alternative hypothesis: true difference in means is not equal to 0 
95 percent confidence interval:
 -4.078 -3.244 
sample estimates:
     mean in group P2 mean in group 4_weeks 
                6.326                 9.987 

# Fit a linear model with a categorical covariate
makeStripplot(mDat <- prepareData('1438786_a_at'))

plot of chunk unnamed-chunk-1

summary(mFit <- lm(gExp ~ devStage, mDat, gType == 'wt'))

Call:
lm(formula = gExp ~ devStage, data = mDat, subset = gType == 
    "wt")

Residuals:
    Min      1Q  Median      3Q     Max 
-1.1565 -0.4400  0.0288  0.4915  1.2065 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)        8.523      0.379   22.50  5.7e-13 ***
devStageP2        -1.450      0.536   -2.71    0.016 *  
devStageP6        -0.107      0.536   -0.20    0.845    
devStageP10       -1.201      0.536   -2.24    0.040 *  
devStage4_weeks    0.081      0.536    0.15    0.882    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.758 on 15 degrees of freedom
Multiple R-squared: 0.497,  Adjusted R-squared: 0.363 
F-statistic: 3.71 on 4 and 15 DF,  p-value: 0.0272 

# Perform inference for a contrast
coef(mFit)
    (Intercept)      devStageP2      devStageP6     devStageP10 
         8.5227         -1.4503         -0.1067         -1.2012 
devStage4_weeks 
         0.0810 
contMat <- rbind(c(0, -1, 0, 1, 0))
obsDiff <- contMat %*% coef(mFit)
sampMeans <- with(subset(mDat, gType == 'wt'),
                  tapply(gExp, devStage, mean))
sampMeans["P2"] - sampMeans["P10"]
    P2 
-0.249 
vcov(mFit)
                (Intercept) devStageP2 devStageP6 devStageP10
(Intercept)          0.1435    -0.1435    -0.1435     -0.1435
devStageP2          -0.1435     0.2870     0.1435      0.1435
devStageP6          -0.1435     0.1435     0.2870      0.1435
devStageP10         -0.1435     0.1435     0.1435      0.2870
devStage4_weeks     -0.1435     0.1435     0.1435      0.1435
                devStage4_weeks
(Intercept)             -0.1435
devStageP2               0.1435
devStageP6               0.1435
devStageP10              0.1435
devStage4_weeks          0.2870
summary(mFit)$coefficients[, "Std. Error"]
    (Intercept)      devStageP2      devStageP6     devStageP10 
         0.3788          0.5357          0.5357          0.5357 
devStage4_weeks 
         0.5357 
sqrt(diag(vcov(mFit)))
    (Intercept)      devStageP2      devStageP6     devStageP10 
         0.3788          0.5357          0.5357          0.5357 
devStage4_weeks 
         0.5357 
estSe <- contMat %*% vcov(mFit) %*% t(contMat)
testStat <- obsDiff/estSe
2 * pt(abs(testStat), df = df.residual(mFit), lower.tail = FALSE)
       [,1]
[1,] 0.3993

# Fit a linear model with two categorical covariates
makeStripplot(oDat <- prepareData("1448690_at"))

plot of chunk unnamed-chunk-1

str(oDat)
'data.frame':   39 obs. of  5 variables:
 $ sample  : num  20 21 22 23 16 17 6 24 25 26 ...
 $ devStage: Factor w/ 5 levels "E16","P2","P6",..: 1 1 1 1 1 1 1 2 2 2 ...
 $ gType   : Factor w/ 2 levels "wt","NrlKO": 1 1 1 1 2 2 2 1 1 1 ...
 $ gene    : Factor w/ 1 level "1448690_at": 1 1 1 1 1 1 1 1 1 1 ...
 $ gExp    : num  8.02 9.05 8.71 8.92 6.8 ...
summary(oFitBig <- lm(gExp ~ gType * devStage, oDat))

Call:
lm(formula = gExp ~ gType * devStage, data = oDat)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.1520 -0.5646  0.0278  0.3120  2.2663 

Coefficients:
                           Estimate Std. Error t value Pr(>|t|)    
(Intercept)                  8.6780     0.3987   21.77  < 2e-16 ***
gTypeNrlKO                  -0.8423     0.6090   -1.38  0.17716    
devStageP2                  -1.0290     0.5638   -1.83  0.07830 .  
devStageP6                  -1.9145     0.5638   -3.40  0.00200 ** 
devStageP10                 -2.1932     0.5638   -3.89  0.00054 ***
devStage4_weeks             -2.0820     0.5638   -3.69  0.00091 ***
gTypeNrlKO:devStageP2        0.0698     0.8299    0.08  0.93352    
gTypeNrlKO:devStageP6        0.1653     0.8299    0.20  0.84348    
gTypeNrlKO:devStageP10       0.2258     0.8299    0.27  0.78745    
gTypeNrlKO:devStage4_weeks   0.6461     0.8299    0.78  0.44257    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.797 on 29 degrees of freedom
Multiple R-squared: 0.59,   Adjusted R-squared: 0.463 
F-statistic: 4.64 on 9 and 29 DF,  p-value: 0.000734 
summary(oFitSmall <- lm(gExp ~ gType + devStage, oDat))

Call:
lm(formula = gExp ~ gType + devStage, data = oDat)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.1609 -0.5314  0.0159  0.3075  2.4753 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)        8.580      0.305   28.16  < 2e-16 ***
gTypeNrlKO        -0.614      0.243   -2.53    0.016 *  
devStageP2        -1.010      0.392   -2.57    0.015 *  
devStageP6        -1.848      0.392   -4.71  4.3e-05 ***
devStageP10       -2.097      0.392   -5.34  6.7e-06 ***
devStage4_weeks   -1.775      0.392   -4.52  7.4e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.757 on 33 degrees of freedom
Multiple R-squared: 0.579,  Adjusted R-squared: 0.515 
F-statistic: 9.08 on 5 and 33 DF,  p-value: 1.65e-05 
anova(oFitSmall, oFitBig)
Analysis of Variance Table

Model 1: gExp ~ gType + devStage
Model 2: gExp ~ gType * devStage
  Res.Df  RSS Df Sum of Sq   F Pr(>F)
1     33 18.9                        
2     29 18.4  4     0.497 0.2   0.94

makeStripplot(oDat <- prepareData('1429225_at'))

plot of chunk unnamed-chunk-1

summary(oFitBig <- lm(gExp ~ gType * devStage, oDat))

Call:
lm(formula = gExp ~ gType * devStage, data = oDat)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.5625 -0.1951  0.0268  0.1857  1.1145 

Coefficients:
                           Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   7.313      0.262   27.94   <2e-16 ***
gTypeNrlKO                   -0.260      0.400   -0.65   0.5203    
devStageP2                   -1.158      0.370   -3.13   0.0040 ** 
devStageP6                   -1.250      0.370   -3.38   0.0021 ** 
devStageP10                  -1.072      0.370   -2.90   0.0071 ** 
devStage4_weeks              -0.909      0.370   -2.46   0.0203 *  
gTypeNrlKO:devStageP2         0.280      0.545    0.51   0.6107    
gTypeNrlKO:devStageP6         0.759      0.545    1.39   0.1742    
gTypeNrlKO:devStageP10        1.791      0.545    3.29   0.0026 ** 
gTypeNrlKO:devStage4_weeks    2.239      0.545    4.11   0.0003 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.523 on 29 degrees of freedom
Multiple R-squared: 0.739,  Adjusted R-squared: 0.659 
F-statistic: 9.14 on 9 and 29 DF,  p-value: 2.1e-06 
summary(oFitSmall <- lm(gExp ~ gType + devStage, oDat))

Call:
lm(formula = gExp ~ gType + devStage, data = oDat)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.3578 -0.3652  0.0332  0.4143  1.7121 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)        6.865      0.272   25.22   <2e-16 ***
gTypeNrlKO         0.784      0.217    3.61   0.0010 ** 
devStageP2        -1.093      0.351   -3.12   0.0038 ** 
devStageP6        -0.945      0.351   -2.69   0.0110 *  
devStageP10       -0.251      0.351   -0.71   0.4798    
devStage4_weeks    0.136      0.351    0.39   0.7003    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.677 on 33 degrees of freedom
Multiple R-squared: 0.504,  Adjusted R-squared: 0.429 
F-statistic: 6.71 on 5 and 33 DF,  p-value: 0.000205 
anova(oFitSmall, oFitBig)
Analysis of Variance Table

Model 1: gExp ~ gType + devStage
Model 2: gExp ~ gType * devStage
  Res.Df   RSS Df Sum of Sq    F Pr(>F)    
1     33 15.12                             
2     29  7.95  4      7.17 6.54  7e-04 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

# Ideas for further work

# We wrote functions to prepare and plot data for more than 1 gene. But when we started fitting models and condusting tests, we only worked with 1 gene at a time. Can you use data aggregation strategies from last week to do some of the same work for small sets of genes?
makeStripplot(oDat <- prepareData(c('1448690_at', '1429225_at')))

plot of chunk unnamed-chunk-1

by(oDat, oDat$gene, function(x) {
  print(summary(oFitBig <- lm(gExp ~ gType * devStage, x)))
  print(summary(oFitSmall <- lm(gExp ~ gType + devStage, x)))
  anova(oFitSmall, oFitBig)})

Call:
lm(formula = gExp ~ gType * devStage, data = x)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.1520 -0.5646  0.0278  0.3120  2.2663 

Coefficients:
                           Estimate Std. Error t value Pr(>|t|)    
(Intercept)                  8.6780     0.3987   21.77  < 2e-16 ***
gTypeNrlKO                  -0.8423     0.6090   -1.38  0.17716    
devStageP2                  -1.0290     0.5638   -1.83  0.07830 .  
devStageP6                  -1.9145     0.5638   -3.40  0.00200 ** 
devStageP10                 -2.1932     0.5638   -3.89  0.00054 ***
devStage4_weeks             -2.0820     0.5638   -3.69  0.00091 ***
gTypeNrlKO:devStageP2        0.0698     0.8299    0.08  0.93352    
gTypeNrlKO:devStageP6        0.1653     0.8299    0.20  0.84348    
gTypeNrlKO:devStageP10       0.2258     0.8299    0.27  0.78745    
gTypeNrlKO:devStage4_weeks   0.6461     0.8299    0.78  0.44257    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.797 on 29 degrees of freedom
Multiple R-squared: 0.59,   Adjusted R-squared: 0.463 
F-statistic: 4.64 on 9 and 29 DF,  p-value: 0.000734 


Call:
lm(formula = gExp ~ gType + devStage, data = x)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.1609 -0.5314  0.0159  0.3075  2.4753 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)        8.580      0.305   28.16  < 2e-16 ***
gTypeNrlKO        -0.614      0.243   -2.53    0.016 *  
devStageP2        -1.010      0.392   -2.57    0.015 *  
devStageP6        -1.848      0.392   -4.71  4.3e-05 ***
devStageP10       -2.097      0.392   -5.34  6.7e-06 ***
devStage4_weeks   -1.775      0.392   -4.52  7.4e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.757 on 33 degrees of freedom
Multiple R-squared: 0.579,  Adjusted R-squared: 0.515 
F-statistic: 9.08 on 5 and 33 DF,  p-value: 1.65e-05 


Call:
lm(formula = gExp ~ gType * devStage, data = x)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.5625 -0.1951  0.0268  0.1857  1.1145 

Coefficients:
                           Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   7.313      0.262   27.94   <2e-16 ***
gTypeNrlKO                   -0.260      0.400   -0.65   0.5203    
devStageP2                   -1.158      0.370   -3.13   0.0040 ** 
devStageP6                   -1.250      0.370   -3.38   0.0021 ** 
devStageP10                  -1.072      0.370   -2.90   0.0071 ** 
devStage4_weeks              -0.909      0.370   -2.46   0.0203 *  
gTypeNrlKO:devStageP2         0.280      0.545    0.51   0.6107    
gTypeNrlKO:devStageP6         0.759      0.545    1.39   0.1742    
gTypeNrlKO:devStageP10        1.791      0.545    3.29   0.0026 ** 
gTypeNrlKO:devStage4_weeks    2.239      0.545    4.11   0.0003 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.523 on 29 degrees of freedom
Multiple R-squared: 0.739,  Adjusted R-squared: 0.659 
F-statistic: 9.14 on 9 and 29 DF,  p-value: 2.1e-06 


Call:
lm(formula = gExp ~ gType + devStage, data = x)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.3578 -0.3652  0.0332  0.4143  1.7121 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)        6.865      0.272   25.22   <2e-16 ***
gTypeNrlKO         0.784      0.217    3.61   0.0010 ** 
devStageP2        -1.093      0.351   -3.12   0.0038 ** 
devStageP6        -0.945      0.351   -2.69   0.0110 *  
devStageP10       -0.251      0.351   -0.71   0.4798    
devStage4_weeks    0.136      0.351    0.39   0.7003    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.677 on 33 degrees of freedom
Multiple R-squared: 0.504,  Adjusted R-squared: 0.429 
F-statistic: 6.71 on 5 and 33 DF,  p-value: 0.000205 
oDat$gene: 1448690_at
Analysis of Variance Table

Model 1: gExp ~ gType + devStage
Model 2: gExp ~ gType * devStage
  Res.Df  RSS Df Sum of Sq   F Pr(>F)
1     33 18.9                        
2     29 18.4  4     0.497 0.2   0.94
-------------------------------------------------------- 
oDat$gene: 1429225_at
Analysis of Variance Table

Model 1: gExp ~ gType + devStage
Model 2: gExp ~ gType * devStage
  Res.Df   RSS Df Sum of Sq    F Pr(>F)    
1     33 15.12                             
2     29  7.95  4      7.17 6.54  7e-04 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

# In lecture we also experimented with a quantitative version of devStage, which we called age. This opens the door to modelling with a quantitative covariate. Can you fit linear and quadratic models to the expression data for one or several genes?
makeStripplot(oDat <- prepareData(c('1448690_at', '1429225_at')))

plot of chunk unnamed-chunk-1

devStageToDays <- c('E2'=-2, 'P2'=2, 'P6'=6, 'P10'=10, '4_weeks'=28)
oDat <- prepareData(c('1448690_at', '1429225_at'))
oDat$age <- devStageToDays[oDat$devStage]
str(oDat)
'data.frame':   78 obs. of  6 variables:
 $ sample  : num  20 21 22 23 16 17 6 24 25 26 ...
 $ devStage: Factor w/ 5 levels "E16","P2","P6",..: 1 1 1 1 1 1 1 2 2 2 ...
 $ gType   : Factor w/ 2 levels "wt","NrlKO": 1 1 1 1 2 2 2 1 1 1 ...
 $ gene    : Factor w/ 2 levels "1448690_at","1429225_at": 1 1 1 1 1 1 1 1 1 1 ...
 $ gExp    : num  8.02 9.05 8.71 8.92 6.8 ...
 $ age     : num  -2 -2 -2 -2 -2 -2 -2 2 2 2 ...

makeXYplot <- function(x, ...)
  xyplot(gExp ~ age | gene, x, group = gType,
            jitter.data = TRUE, auto.key = TRUE,
            type = c("p", "a"), grid = TRUE,
            ...)
makeXYplot(oDat)

plot of chunk unnamed-chunk-1


by(oDat, oDat$gene, function(x) {
  print(summary(oFitBig <- lm(gExp ~ gType + age, x)))
  print(summary(oFitSmall <- lm(gExp ~ gType + age + I(age^2), x)))
  anova(oFitSmall, oFitBig)})

Call:
lm(formula = gExp ~ gType + age, data = x)

Residuals:
   Min     1Q Median     3Q    Max 
-1.287 -0.610 -0.279  0.449  2.920 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   7.6135     0.2490   30.58   <2e-16 ***
gTypeNrlKO   -0.6608     0.3050   -2.17   0.0370 *  
age          -0.0431     0.0147   -2.94   0.0058 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.952 on 36 degrees of freedom
Multiple R-squared: 0.275,  Adjusted R-squared: 0.235 
F-statistic: 6.83 on 2 and 36 DF,  p-value: 0.00305 


Call:
lm(formula = gExp ~ gType + age + I(age^2), data = x)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.1092 -0.5313  0.0055  0.3882  2.4675 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  8.03332    0.21162   37.96  < 2e-16 ***
gTypeNrlKO  -0.61699    0.23762   -2.60    0.014 *  
age         -0.23522    0.04053   -5.80  1.4e-06 ***
I(age^2)     0.00685    0.00139    4.94  1.9e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.741 on 35 degrees of freedom
Multiple R-squared: 0.573,  Adjusted R-squared: 0.536 
F-statistic: 15.7 on 3 and 35 DF,  p-value: 1.27e-06 


Call:
lm(formula = gExp ~ gType + age, data = x)

Residuals:
   Min     1Q Median     3Q    Max 
-1.088 -0.613 -0.116  0.454  1.889 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   6.2171     0.2058   30.21   <2e-16 ***
gTypeNrlKO    0.7469     0.2521    2.96   0.0054 ** 
age           0.0247     0.0121    2.04   0.0489 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.787 on 36 degrees of freedom
Multiple R-squared: 0.269,  Adjusted R-squared: 0.229 
F-statistic: 6.64 on 2 and 36 DF,  p-value: 0.00352 


Call:
lm(formula = gExp ~ gType + age + I(age^2), data = x)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.2352 -0.5075 -0.0413  0.4185  1.6772 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  6.36851    0.21880   29.11   <2e-16 ***
gTypeNrlKO   0.76265    0.24567    3.10   0.0038 ** 
age         -0.04453    0.04191   -1.06   0.2953    
I(age^2)     0.00247    0.00143    1.72   0.0938 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 0.766 on 35 degrees of freedom
Multiple R-squared: 0.326,  Adjusted R-squared: 0.269 
F-statistic: 5.66 on 3 and 35 DF,  p-value: 0.00287 
oDat$gene: 1448690_at
Analysis of Variance Table

Model 1: gExp ~ gType + age + I(age^2)
Model 2: gExp ~ gType + age
  Res.Df  RSS Df Sum of Sq    F  Pr(>F)    
1     35 19.2                              
2     36 32.6 -1     -13.4 24.4 1.9e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 
-------------------------------------------------------- 
oDat$gene: 1429225_at
Analysis of Variance Table

Model 1: gExp ~ gType + age + I(age^2)
Model 2: gExp ~ gType + age
  Res.Df  RSS Df Sum of Sq    F Pr(>F)  
1     35 20.5                           
2     36 22.3 -1     -1.74 2.97  0.094 .
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1