Problems with extraregressors

[Bruno L. Giordano]

Hi,

thank you for this nice package! I am running a few analyses with the extraregressors option. I am consistently getting collinearity errors when the final model is computed (see below). No collinearity problem arises during the preceding steps.

Could you please check whether there is an easy fix for this? 

For the moment I am adding motion regressors. However, the goal is to turn off completely the default detrending regressors and substitute them with a DCT basis for compatibility with SPM.

Cheers,

     Bruno

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

PROBLEMS WITH REGRESSORS COLLINEARITY WHEN EXTRAREGRESSORS INCLUDED

HRFMODEL = 'optimize'; WANTPARAMETRIC = 0 or 1;

*** GLMdenoisedata: calculating parametric fits and error estimates. ***

Error using olsmatrix2 (line 56)

Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND =  3.936974e-08.

Error in GLMdenoisedata (line 896)

    mtimescell(olsmatrix2(results.parametric.designmatrix), ...

hrfmodel = 'assume'; wantparametric = 0/1;

*** GLMdenoisedata: calculating denoised data and PC weights. ***

Error using olsmatrix2 (line 56)

Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND =  3.254563e-08.

SAME PROBLEM WITH NO HRF OPTIMIZATION.

HRFMODEL = 'assume'; WANTPARAMETRIC = 0 or 1;

Error in GLMdenoisedata (line 956)

  h = olsmatrix2(cat(2,polymatrix,exmatrix,pcmatrix))*squish(data{p} - modelcomponent,dimdata)';  % parameters x voxels

AFTER COMMENTING LINES 934/986 TO PREVENT COMPUTATION OF FINAL MODEL

HRFMODEL='optimize' or 'assume'; WANTPARAMETRIC=1;

*** GLMdenoisedata: calculating parametric fits and error estimates. ***

895   results.parametric.parameters = ...

Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND =  3.936974e-08.

> In olsmatrix2 (line 55)

  In GLMdenoisedata (line 896)

Error using olsmatrix2 (line 56)

Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND =  3.936974e-08.

Error in GLMdenoisedata (line 896)

    mtimescell(olsmatrix2(results.parametric.designmatrix), ...

[Kendrick Kay]

thank you for this nice package! I am running a few analyses with the extraregressors option. I am consistently getting collinearity errors when the final model is computed (see below). No collinearity problem arises during the preceding steps.

Could you please check whether there is an easy fix for this? 

For the moment I am adding motion regressors. However, the goal is to turn off completely the default detrending regressors and substitute them with a DCT basis for compatibility with SPM.

As I document on the GLMdenoise web site, I generally don't think it's a good idea to always include motion regressors...   In short, doing so goes against the spirit and theory of data-driven noise derivation.

Regarding the detrending regressors, you could replace the polynomials with the DCT.  You can try changing the polydeg input.  Alternatively (depending on what you mean by "compatibility with SPM"), you might have to go in and hack the polynomial creation and use DCT regressors instead.  I have a function constructdctmatrix.m (in github.com/kendrickkay/knkutils/) that might be helpful in this regard.

If I could see your input to GLMdenoisedata, that could help me understand and diagnose your situation

[Bruno L. Giordano]

Hello again,

thanks for the quick reply!

As I document on the GLMdenoise web site, I generally don't think it's a good idea to always include motion regressors...   In short, doing so goes against the spirit and theory of data-driven noise derivation.

Fair point.

I am trying to compare GLMdenoise decoding with what I get from SPM betas. That's why I am trying to add motion and DCT.

The errors I get arise from a correlation between the polynomial regressors and the motion regressors in the design matrix. I tried to orthogonalize the motion regressors (they are correlated with each other) to no effect. The only thing that seems to work is the joint orthogonalization of the polynomial and motion regressors.

If there is nothing easy that can be done about this, I would just drop the motion and take your advice.

Regarding the detrending regressors, you could replace the polynomials with the DCT.  You can try changing the polydeg input.  Alternatively (depending on what you mean by "compatibility with SPM"), you might have to go in and hack the polynomial creation and use DCT regressors instead.  I have a function constructdctmatrix.m (in github.com/kendrickkay/knkutils/) that might be helpful in this regard.

Yes, this is helpful, thanks! 

The easiest approach in terms of code editing is to hack directly constructpolynomialmatrix like I do at the end of the email. One needs to define the number of non-constant DCs ahead (as is, it requires same length for all runs to work properly).

This brings me to another question: is it strictly necessary to normalize the DCT filter to unit length, as you do for the polynomial regressors?

If I could see your input to GLMdenoisedata, that could help me understand and diagnose your situation

I uploaded some data and the calls here: www.psy.gla.ac.uk/~brunog/tmpGLMden.tar.gz

Anyway, I think the most important part of this is the DCT, which I have implemented already.

Cheers,

     Bruno

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

function f = constructpolynomialmatrix(n,degrees)

% do it

if max(degrees)<100%normal behavior

    f = [];

    temp = linspace(-1,1,n)';

    for p=1:length(degrees)

        % construct polynomial

        polyvector = temp .^ degrees(p);

        

        % orthogonalize with respect to earlier polynomials and make unit length

        polyvector = unitlength(projectionmatrix(f)*polyvector);

        

        % record

        f = cat(2,f,polyvector);

    end

else %DCT

    nDCs=max(degrees)-100;

    f=constructdctmatrix(n,0:nDCs); %includes a constant term

end

--

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Bruno L. Giordano, PhD

Institute of Neuroscience and Psychology

58 Hillhead Street, University of Glasgow

Glasgow, G12 8QB, Scotland

T +44 (0) 141 330 5484

Www: http://www.brunolgiordano.net

Email charter: http://www.emailcharter.org/

[Kendrick Kay]

> I am trying to compare GLMdenoise decoding with what I get from SPM betas. That's why I am trying to add motion and DCT.

I see. Empirical tests and checks are always good in my book. It seems you are trying to aim for a fair or controlled comparison. But I guess I would vote for just comparing GLMdenoise vs. the standard SPM approach.

> The errors I get arise from a correlation between the polynomial regressors and the motion regressors in the design matrix. I tried to orthogonalize the motion regressors (they are correlated with each other) to no effect. The only thing that seems to work is the joint orthogonalization of the polynomial and motion regressors.

Orthogonalization is not going to fix the collinearity issue (the regressors, orthogonalized or not, will still span the exact same subspace).

Since it turns out the polynomials and the motion regressors are highly correlated in your dataset, in some sense, that is an argument for just omitting the motion (and the results will empirically be quite similar).

> This brings me to another question: is it strictly necessary to normalize the DCT filter to unit length, as you do for the polynomial regressors?

No, not necessary. In fact, with least-squares estimation, the results will be functionally equivalent. The main advantage is to control your units so you know exactly what you are getting. (There is another, more minor reason to normalize: it makes the numerical stuff more well behaved. But the scaling difference between unnormalized and the normalized version will be very slight, so no big deal)

[Bruno L. Giordano]

Beautiful. Yes, dropping motion makes sense.

Thanks a lot for the prompt help!

    Bruno

[Kendrick Kay]

> The errors I get arise from a correlation between the polynomial regressors and the motion regressors in the design matrix. I tried to orthogonalize the motion regressors (they are correlated with each other) to no effect. The only thing that seems to work is the joint orthogonalization of the polynomial and motion regressors.

Orthogonalization is not going to fix the collinearity issue (the regressors, orthogonalized or not, will still span the exact same subspace).

Just to clarify this:

Suppose you are getting collinearity issues between a set of regressors A = [A1 A2 A3 ...] and B = [B1 B2 B3 ...]

You could orthogonalize the B regressors with respect to themselves e.g. B = [B1 B2-orthogonalize-with-respect-to-B1 B3-orthogongalize-with-respect-to-B1-and-B2 ...], but that's not going to solve the collinearity between the As and the Bs.

However, if you orthogonalize the Bs with respect to both Bs and As, then that will solve the collinearity issue....

But note that if your regressors are super correlated, then be wary that if you orthogonalize, you are very likely amplifying (and including in your regression model) the noise components of your regressors.  And that's probably a bad thing.

Let me know if you want to have a sketch/discussion/Skype about these issues, which might clarify things.

[Bruno L. Giordano]

A follow up that might not surprise you yet please you: 

GLMdenoise is basically showing double the group-level significant encoding voxels I had with a non-denoising SPM-based approach.

Thank you so much: this is a fantastic tool!! I am glad I invested the time necessary to learn how to use it!

Best,

       Bruno