help
help > RE: Clarification on contrasts in CONN 2nd-level multivariate analysis
Oct 12, 2018 07:10 PM | Alfonso Nieto-Castanon - Boston University
RE: Clarification on contrasts in CONN 2nd-level multivariate analysis
Hi Martyn,
Fair enough, you are not convinced until seeing exactly the same results using both approaches, which sounds like healthy (and useful) scepticism so I am going to give it a try. I have taken a closer look at your scripts and these are, I believe, the additional changes there that should make all 7 statistics exactly equal across the three methods (without the need to whiten the data):
1) in the "Main effect of Location" and "Main effect of Texture" sections change the design matrix to include the between-subject factor. Basically, from:
X3 = ones(30,1);
[t3, df3] = spm_ancova(X3,[],Y3,1);
X5 = blkdiag(ones(28,1),ones(32,1));
[t5, df5] = spm_ancova(X5,[],Y5,eye(2));
to:
X3 = blkdiag(ones(14,1),ones(16,1));
[t3, df3] = spm_ancova(X3,[],Y3,[1 1]');
X5 = kron(blkdiag(ones(14,1),ones(16,1)),eye(2));
[t5, df5] = spm_ancova(X5,[]Y5,[1 0 1 0;0 1 0 1]');
2) in the "Main effect of Texture", "Texture x Drink", "Location x Texture", and "Location x Texture x Drink" sections, orthogonalize your contrast matrices (within-subject contrasts) before computing your reduced data Y, basically changing from:
Y5 = kron(eye(30),[1 1 -1 -1 0 0; 0 0 1 1 -1 -1]) * Y;
Y6 = kron(eye(30),[1 1 -1 -1 0 0; 0 0 1 1 -1 -1]) * Y;
Y7 = kron(eye(30),[1 -1 -1 1 0 0; 0 0 1 -1 -1 1]) * Y;
Y8 = kron(eye(30),[1 -1 -1 1 0 0; 0 0 1 -1 -1 1]) * Y;
to:
Y5 = kron(eye(30),orth([1 1 -1 -1 0 0; 0 0 1 1 -1 -1]')') * Y;
Y6 = kron(eye(30),orth([1 1 -1 -1 0 0; 0 0 1 1 -1 -1]')') * Y;
Y7 = kron(eye(30),orth([1 -1 -1 1 0 0; 0 0 1 -1 -1 1]')') * Y;
Y8 = kron(eye(30),orth([1 -1 -1 1 0 0; 0 0 1 -1 -1 1]')') * Y;
That should produce exactly the same F- values (as well as exactly the same degrees of freedom) across the three methods reported in your scripts (partitioned-errrors, two-stage, and SPSS). Please let me know otherwise. The reason why (1) is needed is because your alternative statistics there use a full-factorial model, so even when you are testing the presence of within-subject effects you should still include the between-subject factor if you want to get exactly the same results across the three methods (and then in your contrast specification simply define a contrast that will aggregate across all levels of that factor). The reason why (2) is needed is because none of the statistics here include any form of whitening so they are not invariant to non-rank-reducing linear transformations of your DVs (and orthogonalizing the contrast just makes sure that we are not modifying the scaling of the data that is being entered into spm_ancova, which, again, we need if we want to get exactly the same results across the three methods). To be precise CONN does not use this form of contrast orthogonalization when computing the intermediate volumes as part of CONN's two-stage implementation, which brings the good point that perhaps I should actually have CONN do that (when multiple within-subject contrasts are entered, orthogonalize the within-subject contrast matrix before computing subject-level intermediate volumes). As far as I can tell this should not make a difference in any of the two approaches used by CONN, namely multivariate analyses (since those explicitly compute the data covariance, so they are exactly invariant to non-rank-reducing linear transformations of the data) nor in univariate SPM analyses using ReML covariance modeling (since again the whitening step there will guarantee invariance to these sort of linear transformations), but there might be some subtleties involved, so I will have to think about this a bit more. Any thoughts/comments/suggestions are welcome.
So, at the very least do you think this should be enough to convince you that the two-stage procedure advocated in the Henson&Penny paper is in fact a valid implementation of M-way ANOVAs with partitioned-errors (irrespective of number of levels of a factor)?
Hope this helps, and let me know your thoughts
Alfonso
Originally posted by Martyn McFarquhar:
Fair enough, you are not convinced until seeing exactly the same results using both approaches, which sounds like healthy (and useful) scepticism so I am going to give it a try. I have taken a closer look at your scripts and these are, I believe, the additional changes there that should make all 7 statistics exactly equal across the three methods (without the need to whiten the data):
1) in the "Main effect of Location" and "Main effect of Texture" sections change the design matrix to include the between-subject factor. Basically, from:
X3 = ones(30,1);
[t3, df3] = spm_ancova(X3,[],Y3,1);
X5 = blkdiag(ones(28,1),ones(32,1));
[t5, df5] = spm_ancova(X5,[],Y5,eye(2));
to:
X3 = blkdiag(ones(14,1),ones(16,1));
[t3, df3] = spm_ancova(X3,[],Y3,[1 1]');
X5 = kron(blkdiag(ones(14,1),ones(16,1)),eye(2));
[t5, df5] = spm_ancova(X5,[]Y5,[1 0 1 0;0 1 0 1]');
2) in the "Main effect of Texture", "Texture x Drink", "Location x Texture", and "Location x Texture x Drink" sections, orthogonalize your contrast matrices (within-subject contrasts) before computing your reduced data Y, basically changing from:
Y5 = kron(eye(30),[1 1 -1 -1 0 0; 0 0 1 1 -1 -1]) * Y;
Y6 = kron(eye(30),[1 1 -1 -1 0 0; 0 0 1 1 -1 -1]) * Y;
Y7 = kron(eye(30),[1 -1 -1 1 0 0; 0 0 1 -1 -1 1]) * Y;
Y8 = kron(eye(30),[1 -1 -1 1 0 0; 0 0 1 -1 -1 1]) * Y;
to:
Y5 = kron(eye(30),orth([1 1 -1 -1 0 0; 0 0 1 1 -1 -1]')') * Y;
Y6 = kron(eye(30),orth([1 1 -1 -1 0 0; 0 0 1 1 -1 -1]')') * Y;
Y7 = kron(eye(30),orth([1 -1 -1 1 0 0; 0 0 1 -1 -1 1]')') * Y;
Y8 = kron(eye(30),orth([1 -1 -1 1 0 0; 0 0 1 -1 -1 1]')') * Y;
That should produce exactly the same F- values (as well as exactly the same degrees of freedom) across the three methods reported in your scripts (partitioned-errrors, two-stage, and SPSS). Please let me know otherwise. The reason why (1) is needed is because your alternative statistics there use a full-factorial model, so even when you are testing the presence of within-subject effects you should still include the between-subject factor if you want to get exactly the same results across the three methods (and then in your contrast specification simply define a contrast that will aggregate across all levels of that factor). The reason why (2) is needed is because none of the statistics here include any form of whitening so they are not invariant to non-rank-reducing linear transformations of your DVs (and orthogonalizing the contrast just makes sure that we are not modifying the scaling of the data that is being entered into spm_ancova, which, again, we need if we want to get exactly the same results across the three methods). To be precise CONN does not use this form of contrast orthogonalization when computing the intermediate volumes as part of CONN's two-stage implementation, which brings the good point that perhaps I should actually have CONN do that (when multiple within-subject contrasts are entered, orthogonalize the within-subject contrast matrix before computing subject-level intermediate volumes). As far as I can tell this should not make a difference in any of the two approaches used by CONN, namely multivariate analyses (since those explicitly compute the data covariance, so they are exactly invariant to non-rank-reducing linear transformations of the data) nor in univariate SPM analyses using ReML covariance modeling (since again the whitening step there will guarantee invariance to these sort of linear transformations), but there might be some subtleties involved, so I will have to think about this a bit more. Any thoughts/comments/suggestions are welcome.
So, at the very least do you think this should be enough to convince you that the two-stage procedure advocated in the Henson&Penny paper is in fact a valid implementation of M-way ANOVAs with partitioned-errors (irrespective of number of levels of a factor)?
Hope this helps, and let me know your thoughts
Alfonso
Originally posted by Martyn McFarquhar:
Hi Alfonso,
This is turning into a very interesting conversation!
Apologies for the error in the code, you are indeed correct on all counts. However, I should clarify what I meant as I was simply saying that partitioned-errors are more sensitive than pooled errors and that the two-stage approach, although claiming to implement partitioned errors, does not agree with the partitioned error approaches in other software. My point was just that we need to get the partitioned-errors right.
Even with the changes to the code, the approaches do not agree. I think we should be aiming for better than "close enough" with these methods, and the two-stage approach (based on Will Penny's recommendations) is still not estimating the F-statistic correctly nor the degrees of freedom for the within-subject models with > 2 levels. This can't be due to differences in the treatment of sphericity because the SPSS results I've reported have not adjusted for non-sphericity and nor have the SPM results. As such, the results should be identical.
My thinking is that the non-sphericity correction should not be necessary in order to achieve the "correct" F-statistic when sphericity is assumed. This speaks to a larger issue in terms of the GLM implementation in other packages that do assume sphericity (such as FSL) where the model we use should allow us to achieve the classical results, under assumptions of sphericity. If one wished to assume sphericity then the two-level models (given in the script I sent, not the ones you have indicated) do not actually achieve the classical results under assumptions of sphericity. My aim would be to get the models to agree under assumptions of sphericity and then apply the non-sphericity correction, as I don't think the correction should be seen as an integral part of the model, but rather an adjunct control for liberal/conservative conclusions when the assumptions are not adequately met.
Additionally, although the code you provided does give the same F-statistic it does not appear to give the same degrees of freedom (SPSS-like = 1.9781,55.3860 [which agrees with the Greenhouse-Geisser correction], CONN-like [2-dimensions] = 1.5814, 44.2780, CONN-like [3-dimensions] = 1.9781,55.3860 [again, Greenhouse-Geisser]). And so the 2-dimension model will give a different p-value to the SPSS and 3-dimension model. I suppose the question in then whether this is to be expected or not as to my mind all the models are testing the same hypothesis with the same data and so the conclusions should be invariant to this sort of choice? This is the claim made by the two-stage approach, but as far as I can tell it still does not agree.
Furthermore, even if the degrees of freedom agree with the Greenhouse-Geisser correction, this correction should be applied to the original F-statistic of 1.167, rather than the F-statistic of 1.1006 produced by the code. As such, the SPSS Greenhouse-Geisser corrected results of F(1.98,55.39) = 1.167, p = 0.319 is still not matched by the non-sphericity-corrected results of F(1.98,55.39) = 1.1009, p = 0.340.
I know this seems like splitting hairs, but we need to be convinced that the models we are using are appropriate and I still think the best way to do so is make sure we can achieve equivalency with the classical results in the textbooks (and implemented in other software) before adding on additional corrections (particularly when these corrections are not available in all software). This is easy enough when dealing with 2-levels, but > 2-levels is still a sticking point for me.
Looking forward to hearing your thoughts!
Best wishes,
- Martyn
This is turning into a very interesting conversation!
Apologies for the error in the code, you are indeed correct on all counts. However, I should clarify what I meant as I was simply saying that partitioned-errors are more sensitive than pooled errors and that the two-stage approach, although claiming to implement partitioned errors, does not agree with the partitioned error approaches in other software. My point was just that we need to get the partitioned-errors right.
Even with the changes to the code, the approaches do not agree. I think we should be aiming for better than "close enough" with these methods, and the two-stage approach (based on Will Penny's recommendations) is still not estimating the F-statistic correctly nor the degrees of freedom for the within-subject models with > 2 levels. This can't be due to differences in the treatment of sphericity because the SPSS results I've reported have not adjusted for non-sphericity and nor have the SPM results. As such, the results should be identical.
My thinking is that the non-sphericity correction should not be necessary in order to achieve the "correct" F-statistic when sphericity is assumed. This speaks to a larger issue in terms of the GLM implementation in other packages that do assume sphericity (such as FSL) where the model we use should allow us to achieve the classical results, under assumptions of sphericity. If one wished to assume sphericity then the two-level models (given in the script I sent, not the ones you have indicated) do not actually achieve the classical results under assumptions of sphericity. My aim would be to get the models to agree under assumptions of sphericity and then apply the non-sphericity correction, as I don't think the correction should be seen as an integral part of the model, but rather an adjunct control for liberal/conservative conclusions when the assumptions are not adequately met.
Additionally, although the code you provided does give the same F-statistic it does not appear to give the same degrees of freedom (SPSS-like = 1.9781,55.3860 [which agrees with the Greenhouse-Geisser correction], CONN-like [2-dimensions] = 1.5814, 44.2780, CONN-like [3-dimensions] = 1.9781,55.3860 [again, Greenhouse-Geisser]). And so the 2-dimension model will give a different p-value to the SPSS and 3-dimension model. I suppose the question in then whether this is to be expected or not as to my mind all the models are testing the same hypothesis with the same data and so the conclusions should be invariant to this sort of choice? This is the claim made by the two-stage approach, but as far as I can tell it still does not agree.
Furthermore, even if the degrees of freedom agree with the Greenhouse-Geisser correction, this correction should be applied to the original F-statistic of 1.167, rather than the F-statistic of 1.1006 produced by the code. As such, the SPSS Greenhouse-Geisser corrected results of F(1.98,55.39) = 1.167, p = 0.319 is still not matched by the non-sphericity-corrected results of F(1.98,55.39) = 1.1009, p = 0.340.
I know this seems like splitting hairs, but we need to be convinced that the models we are using are appropriate and I still think the best way to do so is make sure we can achieve equivalency with the classical results in the textbooks (and implemented in other software) before adding on additional corrections (particularly when these corrections are not available in all software). This is easy enough when dealing with 2-levels, but > 2-levels is still a sticking point for me.
Looking forward to hearing your thoughts!
Best wishes,
- Martyn
Threaded View
| Title | Author | Date |
|---|---|---|
| Martyn McFarquhar | Oct 4, 2018 | |
| Alfonso Nieto-Castanon | Oct 5, 2018 | |
| Martyn McFarquhar | Oct 5, 2018 | |
| Alfonso Nieto-Castanon | Oct 5, 2018 | |
| Martyn McFarquhar | Oct 8, 2018 | |
| Alfonso Nieto-Castanon | Oct 8, 2018 | |
| Martyn McFarquhar | Oct 9, 2018 | |
| Alfonso Nieto-Castanon | Oct 9, 2018 | |
| Ali Amad | Oct 11, 2018 | |
| Alfonso Nieto-Castanon | Oct 12, 2018 | |
| Ali Amad | Oct 19, 2018 | |
| Martyn McFarquhar | Oct 11, 2018 | |
| Alfonso Nieto-Castanon | Oct 12, 2018 | |
| Martyn McFarquhar | Oct 12, 2018 | |
| Alfonso Nieto-Castanon | Oct 12, 2018 | |
| Martyn McFarquhar | Oct 15, 2018 | |
| Martyn McFarquhar | Oct 5, 2018 | |