Added nuisance variables to regression

Author: markallenthornton

demos/fwe_control.m

@@ -85,9 +85,9 @@
 for i = 1:nsim
     imgs = randn(4,4,4,20);
     covariate = [ones(20,1) randn(20,4)];
-    [pcorr_pos,pcorr_neg] = matlab_tfce('regression',2,imgs,[],covariate);
+    [pcorr_pos,pcorr_neg] = matlab_tfce('regression',2,imgs,[],covariate,[],[],[],[],[],[],[1 0 0 1 1]);
     for j = 1:5
         fp(i,j) = sum(pcorr_pos{j}(:)<.05)+sum(pcorr_neg{j}(:)<.05);
     end
 end
-fp/nsim
+sum(fp)/nsim

matlab_tfce.m

@@ -25,10 +25,10 @@
 % corrections, is substantially less conservative due to the fact that
 % it capitalizes on spatial dependency in the data. 
 %
-% [varargout] = matlab_tfce(analysis,tails,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers)
-% [pcorr] = matlab_tfce(analysis,1,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers)
-% [pcorr_pos,pcorr_neg] = matlab_tfce(analysis,2,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers)
-% [pcorr_fac1,pcorr_fac2,pcorr_int] = matlab_tfce(analysis,1,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers)
+% [varargout] = matlab_tfce(analysis,tails,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers,nuisance)
+% [pcorr] = matlab_tfce(analysis,1,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers,nuisance)
+% [pcorr_pos,pcorr_neg] = matlab_tfce(analysis,2,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers,nuisance)
+% [pcorr_fac1,pcorr_fac2,pcorr_int] = matlab_tfce(analysis,1,imgs,imgs2,covariate,nperm,H,E,C,dh,parworkers,nuisance)
 %
 % Arguments:
 %
@@ -89,6 +89,14 @@
 % time consuming so it may be more efficient to do outside the package if
 % one is running several analyses.
 %
+% nuisance -- commandline-only argument for specifying which variables in
+% a multiple regression are nuisance regressors. Accepts a (1 x predictor)
+% binary (0/1 true/false) vector. 1s/trues indicate variables that should
+% be treated as nuisance predictors (vs. variables of interest). These
+% variables will be included in the model, but will not be subjected to
+% TFCE, accelerating calculations. Default is all 0s (i.e. all variables
+% treated as of non-nuisance). Ignorned if analysis not regression.
+%
 % Output: 
 % If tails == 1, a single output image with the same xyz dimensions as imgs
 % consisting of corrected p-values with be returned.
@@ -135,6 +143,7 @@
 if nargin > (fixedargn + 1)
     covariate = varargin{2};
 end
+nuisance = zeros([1 size(covariate,2)]);
 if nargin > (fixedargn + 2)
     if ~isempty(varargin{3})
         nperm = varargin{3};
@@ -165,6 +174,11 @@
         parworkers = varargin{8};
     end
 end
+if nargin > (fixedargn + 8)
+    if ~isempty(varargin{9})
+        nuisance = varargin{9};
+    end
+end
 
 % check that tails are appropriate
 if ~(sum(tails==[1 2]))
@@ -270,6 +284,14 @@
     end
 end
 
+% check nuisance
+if strcmp(analysis,'regression')
+    ns = size(nuisance);
+    if ~isequal(ns,[1 size(covariate,2)])
+        error('nuisance size does not equal predictor number');
+    end
+end
+
 %% analysis calls
 
 % indicate analysis start
@@ -325,9 +347,9 @@
     % multiple regression 
     case 'regression'
         if tails == 1
-            pcorr = matlab_tfce_regression(imgs,covariate,tails,nperm,H,E,C,dh);
+            pcorr = matlab_tfce_regression(imgs,covariate,tails,nperm,H,E,C,dh,nuisance);
         else
-            [pcorr_pos,pcorr_neg] = matlab_tfce_regression(imgs,covariate,tails,nperm,H,E,C,dh);
+            [pcorr_pos,pcorr_neg] = matlab_tfce_regression(imgs,covariate,tails,nperm,H,E,C,dh,nuisance);
         end
 
     % unrecognized analysis input

matlab_tfce_regression.m

@@ -1,4 +1,4 @@
-function [varargout] = matlab_tfce_regression(imgs,preds,tails,nperm,H,E,C,dh)
+function [varargout] = matlab_tfce_regression(imgs,preds,tails,nperm,H,E,C,dh,nuisance)
 % MATLAB_TFCE_REGRESSION computes TFCE corrected p-values for each
 % parameter in a multiple regression
 %
@@ -14,6 +14,8 @@
 %   -- E extent exponent
 %   -- C connectivity
 %   -- ndh step number for cluster formation
+%   -- nuisance specifies which variables are not of interest
+%   -- omnibusf requests overall model significance test
 %
 %   Output:
 %	If tails == 1:
@@ -112,15 +114,17 @@
     rtstats = bs./serr;
     curexceeds = zeros(nvox,npred);
     for i = 1:npred
-        rbrain=NaN(bsize);
-        rbrain(implicitmask) = rtstats(i,:);
-        rbrain = transform(rbrain,H,E,C,dh);
-        rstats = rbrain(implicitmask);
-        if tails == 2
-            rstats = abs(rstats);
-        end  
-        % compare maxima to t-values and increment as appropriate
-        curexceeds(:,i) = max(rstats) >= tvals(i,:);
+        if ~nuisance(i)
+            rbrain=NaN(bsize);
+            rbrain(implicitmask) = rtstats(i,:);
+            rbrain = transform(rbrain,H,E,C,dh);
+            rstats = rbrain(implicitmask);
+            if tails == 2
+                rstats = abs(rstats);
+            end  
+            % compare maxima to t-values and increment as appropriate
+            curexceeds(:,i) = max(rstats) >= tvals(i,:);
+        end
     end
     exceedances = exceedances + curexceeds;
 
@@ -139,28 +143,38 @@
 end
 
 for i = 1:npred
-    % create corrected p-value image
-    corrected = exceedances(:,i)./nperm;
-    pcorr = ones(bsize);
-    pcorr(implicitmask) = corrected;
+    if ~nuisance(i)
+        % create corrected p-value image
+        corrected = exceedances(:,i)./nperm;
+        pcorr = ones(bsize);
+        pcorr(implicitmask) = corrected;
 
-    % split into positive and negative effects (if needed)
-    if tails == 2
-        btruestat = NaN(bsize);
-        btruestat(implicitmask) = truestat(i,:);
-        pos = btruestat>0;
-        pcorr_pos = pcorr;
-        pcorr_pos(~pos) = 1;
-        pcorr_neg = pcorr;
-        pcorr_neg(pos) = 1;
-    end
-
-    % assign to cormaps
-    if tails == 1
-        cormaps{i} = pcorr;
+        % split into positive and negative effects (if needed)
+        if tails == 2
+            btruestat = NaN(bsize);
+            btruestat(implicitmask) = truestat(i,:);
+            pos = btruestat>0;
+            pcorr_pos = pcorr;
+            pcorr_pos(~pos) = 1;
+            pcorr_neg = pcorr;
+            pcorr_neg(pos) = 1;
+        end
+
+        % assign to cormaps
+        if tails == 1
+            cormaps{i} = pcorr;
+        else
+            cormaps_pos{i} = pcorr_pos;
+            cormaps_neg{i} = pcorr_neg;
+        end
     else
-        cormaps_pos{i} = pcorr_pos;
-        cormaps_neg{i} = pcorr_neg;
+        % assign to cormaps
+        if tails == 1
+            cormaps{i} = NaN;
+        else
+            cormaps_pos{i} = NaN;
+            cormaps_neg{i} = NaN;
+        end
     end
 end