Wickramasinghe et al. BMC Genomics 2012, 13:45 http://www.biomedcentral.com/1471-2164/13/45 280 ratio, gel electrophoresis and capillary electrophor- esis with an Experion bio-analyzer (Bio-Rad, Hercules, CA). RNA sequencing and data analysis Gene expression analysis (three samples per lactation stage) was conducted by Illumina RNA sequencing (RNA-Seq) technology. Messenger RNA was isolated and purified using an RNA-Seq sample preparation kit (Illumina, San Diego, CA). Then mRNA was fragmented to approximately 200 bp fragments and first and second strand cDNA were synthesized, followed by end repair and adapter ligation. The fragments were purified and sequenced at the UC Davis Genome Center DNA Tech- nologies Core Facility using the Illumina Genome Ana- lyzer (GAII). Short sequence reads of 36-40 bp were assembled and analyzed in RNA-Seq and expression analysis application of CLC Genomics Workbench 3.7 (CLC Bio, Aarhus, Denmark). The bovine genome Btau 4.0 http://www.ncbi.nlm.nih.gov/genome/guide/cow/ index.html was utilized as the reference genome for the assembly. The following criteria were used to filter the unique sequence reads: minimum length fraction of 0.9; minimum similarity fraction of 0.8; maximum number of two mismatches. Data were normalized by calculating the reads per kilobase per million mapped reads (RPKM = total exon reads/mapped reads in millions × exon length in kb) [18] for each gene and annotated with NCBI bovine genome assembly (27,368 unique genes). The experiment was carried out in two steps with day 15 and day 250 samples collected from the same cows and day 90 samples collected from different cows. Initial model checking was conducted to test for independence of samples collected at two time points from the same cow. Using R, a model was fit with animal and stage of lactation as two (non-interacting) factors on the expres- sion of 27,368 unique genes in day 15 and day 250 sam- ples. P values were obtained for the animal effect and the distribution of the p values was plotted. The major- ity of the genes had reasonably uniform p values between 0.4-0.5, and there was no significant animal effect on the analysis. Because of this uniform distribu- tion of p values and the 235-day interval between sam- ples collected from the same cow, samples were assumed to be independent of each other. T-tests and ANOVA were performed on log2-transformed data (0.5 was added to each number before log transformation to deal with zero counts) to identify the genes with signifi- cant changes in expression (p ≤ 0.05, FDR q ≤ 0.3) between the stages of lactations. GO annotation Coding sequences of genes with high expression in each stage of lactation and genes that showed statistically Page 12 of 14 significant changes between the lactation stages were obtained from the ENSEMBL biomart martview applica- tion http://www.ensembl.org/biomart/martview. These sequences were imported to the Blast2GO program to perform the blastx, mapping and GO annotation [28]. Statistical assessment of annotation differences between lactation stages were performed using all expressed genes as background and Fisher’s Exact Test for multi- ple test correction in Blast2GO. Pathway analysis MetaCore pathway analysis by GeneGo, a Thomson Reuters business, was used to identify the significant Gene GO pathways and Gene GO metabolic networks in genes with high expression and statistically significant changes in expression in each stage of lactation. This is calculated using a built in function of MetaCore soft- ware that uses a variation of the Fisher’s exact test adjusted for multiple sample testing using the Benja- mini-Hochberg FDR analysis. Human orthologs were used for the analysis. Additional material Additional file 1: Figure S1: Biological process Blast2GO annotation results of highly expressed genes at three stages of lactation. Biological process GO terms are shown in the x-axis. Percentages of genes belonging to each category are shown in y-axis. For all the three stages, the highest numbers of biological process GO-terms were found for cellular and metabolic process. Additional file 2: Table S1. Top GeneGo pathways identified in the genes with ubiquitous expression in MSC. Additional file 3: Table S2. Top GeneGo pathways identified in the genes with statistically significant changes in expression between transition lactation MSC and peak lactation MSC. Additional file 4: Table S3. Top GeneGo pathways identified in the genes with statistically significant changes in expression between transition lactation MSC and late lactation MSC. Additional file 5: Figure S2: Comparison of biological process GO terms for the genes with statistically significant changes in expression between day 90 and day 250. Biological process GO terms are shown in the x-axis. Percentages of genes belonging to each category are shown in y-axis. In both stages there is enrichment of cellular process and metabolic process GO terms. Day 250 milk had comparatively higher number of GO terms for cellular process, biological regulation, localization, signaling, cellular component organization, cell proliferation, death, immune system process, growth and multi-organism process. Additional file 6: Table S4. Top GeneGo pathways identified in the genes with statistically significant changes in expression between peak lactation MSC and late lactation MSC. Acknowledgements The authors thank Doug Gisi and Sharlie Cunningham at the UC Davis dairy facility for their valuable help in sample selection and managing the dairy herd, and Charlie Nicolet at UC Davis genome center for his excellent technical expertise to perform Illumina GAII sequencing. We thank Tom Famula for his help in writing an R script used in the statistical analysis. S. Wickramasinghe was supported by the Austin Eugene Lyons Graduate