0
$\begingroup$

In each of the below 2 code snippets, for 20 times, I generate an input matrix $X$ containing 1000 variables and 50 samples, and and an output vector $y$ which is generated from $X$ using a multivariate linear model. And my aim is to perform leave-one-out cross-validation on this data, so predict the response $y$ for each sample using the other samples for training. And my prediction performance metric is the spearman correlation between the true response and predicted response. The 2 code snippets differ in terms of how the prediction is done in the second loop. In the first code, I use correlation coefficients as the variable weights in prediction, while in the second one, I predict using a simple multiple linear regression fit.

    library(foreach)
    library(doParallel)
    registerDoParallel(cores=10)

    randcnt = 20
    varsize = 1000
    samplesize = 100

    cors1 = foreach (r = 1:randcnt, .combine='c') %dopar% {
        set.seed(r*334)
        X = matrix(rnorm(samplesize*varsize), ncol=varsize) # of size #samples x #genes
        set.seed(r*344)
        w = matrix(rnorm(varsize), ncol=1) # of size #genes x 1
        set.seed(r*354)
        eps = matrix(rnorm(samplesize), ncol=1) # of size #samples x 1
        y = X %*% w + eps

        pred = foreach (i = 1:nrow(X), .combine='c') %dopar% {
            print(paste(r, i))
            trainsamp = setdiff(1:nrow(X), i)
            trainX = X[trainsamp,]
            trainy = y[trainsamp]
            testX = X[i,,drop=F]

            w = foreach (i = 1:ncol(trainX), .combine='c') %dopar% {
                cor(trainX[,i], trainy)
            }
            sum(testX*w)
        }
        cor(y, pred, method='spearman')
    }

    library(foreach)
    library(doParallel)
    registerDoParallel(cores=10)

    library(glmnet)

    randcnt = 20
    varsize = 1000
    samplesize = 100

    cors2 = foreach (r = 1:randcnt, .combine='c') %dopar% {
        set.seed(r*334)
        X = matrix(rnorm(samplesize*varsize), ncol=varsize) # of size #samples x #genes
        set.seed(r*344)
        w = matrix(rnorm(varsize), ncol=1) # of size #genes x 1
        set.seed(r*354)
        eps = matrix(rnorm(samplesize), ncol=1) # of size #samples x 1
        y = X %*% w + eps

        pred = foreach (i = 1:nrow(X), .combine='c') %dopar% {
            print(paste(r, i))
            trainsamp = setdiff(1:nrow(X), i)
            trainX = X[trainsamp,]
            trainy = y[trainsamp]
            testX = X[i,,drop=F]

            glmob = glmnet(trainX, trainy, lambda=0)
            predict(glmob, newx=testX)
        }
        cor(y, pred, method='spearman')
    }

Interestingly, I realized that, the correlations from the first method ($cors1$) are oftentimes (18 out of 20 runs) higher than the correlations from the second method ($cors2$). The mean correlation is also higher for the first method. Although the actual predictions from the first method are much bigger than the true correlations in absolute value (as expected), the order of them is better in the correlation-based prediction than in regression-based prediction. This is very confusing to me. What is the reason behind that?

$\endgroup$
  • $\begingroup$ I tried your script, but it says ss not found. Can you fix it? $\endgroup$ – SmallChess Jan 24 '17 at 7:50
  • 2
    $\begingroup$ 1. This isn't a code review site. It would be better to just explain what you're doing and why. 2. Spearman correlation is a terrible way to measure predictive performance. $\endgroup$ – Glen_b Jan 24 '17 at 9:29
  • $\begingroup$ 1. I already explained what I did. I put the code as well so that people can try it to have a better reply. I don't want the code to be reviewed or anything like that. The code is already working. 2. I already know it is a terrible way, but I want to understand why it gives that good ranking of the samples then. $\endgroup$ – user5054 Jan 24 '17 at 18:30

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.