| | 114 | * **Remove batch effects** |
| | 115 | * With [[http://www.bioconductor.org/packages/release/bioc/html/edgeR.html|edgeR]], the treatments can be adjusted for differences between the batch with addictive model formula by design = model.matrix( ~Batch+Treatment). For example, with target.txt like this: |
| | 116 | |
| | 117 | |
| | 118 | ||= Sample =||= genotype =||= date =|| |
| | 119 | || KO.1 || KO || old || |
| | 120 | || KO.2 || KO || old || |
| | 121 | || WT.1 || WT || old || |
| | 122 | || KO.3 || KO || new || |
| | 123 | || KO.4 || KO || new || |
| | 124 | || WT.2 || WT || new || |
| | 125 | |
| | 126 | |
| | 127 | Sample R codes for reducing the date effect (old/new) are |
| | 128 | |
| | 129 | {{{ |
| | 130 | target <- read.delim("target.txt", header=T) |
| | 131 | genotype <- factor(target$genotype, levels=c("WT", "KO")) |
| | 132 | mydate <- factor(target$date, levels=c("old", "new")) |
| | 133 | Design <-model.matrix(~mydate+genotype) |
| | 134 | colnames(Design) |
| | 135 | # [1] "(Intercept)" "mydatenew" "genotypeKO" |
| | 136 | |
| | 137 | Y2 <- calcNormFactors(Y) # Y is the DGEList |
| | 138 | Y2 <- estimateGLMCommonDisp(Y2, Design, verbose=TRUE) |
| | 139 | Y2 <- estimateGLMTrendedDisp(Y2, Design) |
| | 140 | Y2 <- estimateGLMTagwiseDisp(Y2, Design) |
| | 141 | Fit <- glmFit(Y2, Design) |
| | 142 | Lrt <- glmLRT(Fit, coef="genotypeKO") # coef is pointed to the column name in Design |
| | 143 | }}} |
| | 144 | |
| | 145 | Lrt$table has the log2FC(Fold change) and P-value. |
| | 146 | |
| | 147 | Detail information can be found in [[http://www.bioconductor.org/packages/release/bioc/html/edgeR.html|edgeR]] |
| | 148 | |
| | 149 | |