MATLAB: Do I get different values with ‘fitlmematrix’ with MATLAB R2015a and R2016a

Statistics and Machine Learning Toolbox

If I execute the code from 'testLMM.m' on data loaded from attached MAT-file, the P-values I get with MATLAB R2015a and MATLAB R2016a are very different. What is the reason behind this?

Best Answer

The function 'fitlmematrix' uses Quasi-Newton optimization under the hood. Between MATLAB R2015a and R2016a, we changed the value of the gradient used for performing this optimization for improved performance. This should not affect most computations.
However, in attached data, the input matrices X, Y and Z are ill-conditioned. The order of magnitude of their columns are very different, ranging from 1e01 to 1e05. This results in the divergence in computed results with the new optimization algorithm. 
For all optimization problems, we recommend that you normalize the input matrices in order to get the most accurate results. For example, in this case, the matrices can be normalized as follows:
[rows,~]=size(X);%# X is your matrix
colMax=max(abs(X),[],1);%# take max absolute value to account for negative numbers


normalizedX=X./repmat(colMax,rows,1);
 
[rows,~]=size(Y);%# Y is your matrix
colMax=max(abs(Y),[],1);%# take max absolute value to account for negative numbers
normalizedY=Y./repmat(colMax,rows,1);
 
normalizedZ = cell(1,numel(Z));
for ii = 1:numel(Z)
[rows,~]=size(Z{ii});%# Z is your matrix
colMax=max(abs(Z{ii}),[],1);%# take max absolute value to account for negative numbers
normalizedZ{ii}=Z{ii}./repmat(colMax,rows,1);   
end
 
 
for i=1:13
    lmm = fitlmematrix(normalizedX,normalizedY(:,i),normalizedZ,G,'CovariancePattern',{'Diagonal' 'Diagonal' 'Diagonal' 'Diagonal''Diagonal'},'FitMethod','REML');
    p(i)=coefTest(lmm,H);
end
 
The above code results in the same P-values in both versions.