I have a gene expression data-set with log2-transformed expression values (no NAs) for 495 genes for 59 samples for which values of a continuous response variable (r) are also known (no NAs). I want to use leave-one-out cross validation to test if r of the test sample can be predicted from the sample's gene expression.
For this, I intend to use the samr R package for Significance Analysis of Microarrays to identify significant genes associated with r in the training set of samples. Then, I want to generate a linear model using the significant genes as variables, which will then be used to predict r of the test sample. I have tried the following code to begin with, but when I generate the model and examine it, I see many NAs in the model summary, which makes me suspect that I am doing something wrong.
Can someone tell me what I might be doing wrong?
Secondly, I will appreciate any comment on the use of nperms (in SAM) with a value of 100. Is it too low for an expression data-set for 495 genes.
# rVals with the r values is read as a vector from a row of a table for phenotypic data read from a tab-delimited file with sample-names as column names and phenotype features as row-names # geneVals is the log2-transformed gene expression data-set read as a matrix from a tab-delimited file with sample-names as column names and gene-names as row-names # Perform SAM with FDR of 5% and obtain list of significant genes sam <- SAM(x=geneVals, y=rVals, resp.type=c("Quantitative"), testStatistic=c("standard"), regression.method=c("standard"), logged2=TRUE, fdr.output=0.05, eigengene.number=1, knn.neighbors=10, nperms=100, genenames=as.vector(rownames(geneVals))) sigGenes <- rbind(sam$siggenes.table$genes.up, sam$siggenes.table$genes.lo) # Generate linear model toModel <- data.frame(t(rbind(rVals, geneVals)), check.names=FALSE) myModel <- lm(toModel[c('r', sigGenes[,c("Gene ID")])]) # Examine model summary(myModel) ...output... Call: lm(formula = toModel[c("rVals", sigGenes[, c("Gene ID")])]) Residuals: ALL 59 residuals are 0: no residual degrees of freedom! Coefficients: (58 not defined because of singularities) Estimate Std. Error t value Pr(>|t|) (Intercept) -18.29363 NA NA NA `let-7e` -1.70545 NA NA NA `miR-125a-5p` 2.43177 NA NA NA `miR-151-5p` 2.67439 NA NA NA ...