Neuroscientific investigation has greatly benefited from the combination of functional Magnetic Resonance Imaging (fMRI) with linearized encoding, which allows to validate and compare computational models of neural activity based on neuroimaging data. In linearized encoding, a multidimensional feature space, usually obtained from a computational model applied to the stimuli, is related to the measured brain activity. This is often done by mapping such space to a dataset (training data, or in-sample), and validating the mapping on a separate dataset (test data, or out-of-sample), to avoid overfitting. When comparing models, the one with the highest explained variance on the test data, as indicated by the coefficient of determination (R2 ), is the one that better reflects the neural computations performed by the brain. An implicit assumption underlying this procedure is that the out-of-sample R2 is an unbiased estimator of the explanatory power of a computational model in the population of stimuli, and can therefore be safely used to compare models. In this work, we show that this is not the case, as the out-of-sample R2 has a negative bias, related to the amount of overfitting in the training data. This phenomenon has dramatic implications for model comparison when models of different dimensionalities are compared. To this aim, we develop an analytical framework that allows us to evaluate and correct biases in both in- and out-of-sample R2, with and without L2 regularization. Our proposed approach yields unbiased estimators of the population R2, thus enabling a valid model comparison. We validate it through illustrative simulations and with an application to a large public fMRI dataset.