This is a data report on the linear modeling analysis using affymetrix microarrays from RNA isolated from spleen tissue in CC mice at days 2, 4, 7 ,12 post WNV infection. MOre information on linear modeling analysis in R using microarray data can be found here:

http://www.bioconductor.org/packages/devel/bioc/vignettes/limma/inst/doc/usersguide.pdf

# display all load R and bioconductor packages along with their versions
library(limma)

## Warning: package 'limma' was built under R version 3.2.4

sessionInfo()

## R version 3.2.3 (2015-12-10)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 14393)
## 
## locale:
## [1] LC_COLLATE=English_United States.1252 
## [2] LC_CTYPE=English_United States.1252   
## [3] LC_MONETARY=English_United States.1252
## [4] LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.1252    
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] limma_3.26.9
## 
## loaded via a namespace (and not attached):
##  [1] backports_1.0.4 magrittr_1.5    rprojroot_1.1   tools_3.2.3    
##  [5] htmltools_0.3.5 yaml_2.1.14     Rcpp_0.12.8     stringi_1.1.2  
##  [9] rmarkdown_1.3   knitr_1.15.1    stringr_1.1.0   digest_0.6.11  
## [13] evaluate_0.10

library(limma)

# load normalized expression matrix and target file

# https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE91003

norm_matrix <- read.csv(file="C:\\gale_lab\\oas1b_manuscript\\norm_matrix_geo.csv", row.names=1, header=T)

# load target file
targets <- read.csv(file="C:\\gale_lab\\oas1b_manuscript\\target.csv", header=T)

f <- factor(targets$Treatment, levels = unique(targets$Treatment))
design <- model.matrix(~0 + f)
colnames(design) <- levels(f)


cont.matrix <- makeContrasts(x3032x16188_W_Sp_2 - x3032x16188_M_Sp_2, 
x3032x16188_W_Sp_4 - x3032x16188_M_Sp_2, 
x3032x16188_W_Sp_7 - x3032x16188_M_Sp_2, 
x3032x16188_W_Sp_12 - x3032x16188_M_Sp_2, 
x16513x16188_W_Sp_2- x16513x16188_M_Sp_12,
x16513x16188_W_Sp_4- x16513x16188_M_Sp_12,
x16513x16188_W_Sp_7- x16513x16188_M_Sp_12,
x16513x16188_W_Sp_12- x16513x16188_M_Sp_12,
x5489x16557_W_Sp_2 - x5489x16557_M_Sp_2, 
x5489x16557_W_Sp_4 - x5489x16557_M_Sp_2,
x5489x16557_W_Sp_7 - x5489x16557_M_Sp_2, 
x5489x16557_W_Sp_12 - x5489x16557_M_Sp_2, 
x16750x13421_W_Sp_4 - x16750x13421_M_Sp_12,
x16750x13421_W_Sp_7 - x16750x13421_M_Sp_12,
x16750x13421_W_Sp_12 - x16750x13421_M_Sp_12,
x8034x8048_W_Sp_2 - x8034x8048_M_Sp_2, 
x8034x8048_W_Sp_4 - x8034x8048_M_Sp_2, 
x8034x8048_W_Sp_7 - x8034x8048_M_Sp_2, 
x8034x8048_W_Sp_12 - x8034x8048_M_Sp_2, 
x13067x5306_W_Sp_2 - x13067x5306_M_Sp_12,
x13067x5306_W_Sp_4 - x13067x5306_M_Sp_12,
x13067x5306_W_Sp_7 - x13067x5306_M_Sp_12,
x13067x5306_W_Sp_12 - x13067x5306_M_Sp_12,
x16680x8016_W_Sp_2 - x16680x8016_M_Sp_12, 
x16680x8016_W_Sp_4 - x16680x8016_M_Sp_12, 
x16680x8016_W_Sp_7 - x16680x8016_M_Sp_12, 
x16680x8016_W_Sp_12 - x16680x8016_M_Sp_12,levels=design)

#Apply a linear model

fit <- lmFit(norm_matrix, design)

# apply our comparisons via a contrast matrix

fit2 <- contrasts.fit(fit,cont.matrix)

# Apply empirical bayes statistics for differential expression
fit2 <- eBayes(fit2)

# use a threshold of |1.5| fold change and an adjusted p-value of <=.05

results <- decideTests(fit2, lfc=(.58), method="separate", adjust.method="BH", p.value=0.05);

# write DE results to a file 


write.fit(fit2, file="C:\\gale_lab\\oas1b_manuscript\\DE_oas1b.txt", digits=3, method="separate", adjust="BH")

DE_results <- read.table(file="C:\\gale_lab\\oas1b_manuscript\\DE_oas1b.txt", header=T, sep="\t")

#Add gene names to DE list

row.names(DE_results) <- row.names(norm_matrix)

write.csv(DE_results, file="C:\\gale_lab\\oas1b_manuscript\\DE_oas1b.csv")