Dear Satra and list,
SPM won't let you estimate a design matrix with empty regressors as all
the other answers have correctly suggested. I will outline a solution
for a conditional design and contrast specification below, but I must
stress that it is complicated and therefore I would therefore strongly
recommend using the batching facilities in SPM2. If you are not
familiar with batching you will find a lot of sample batch files on this
list, e.g. Karl Friston's sample file
http://www.jiscmail.ac.uk/cgi-bin/wa.exe?A2=ind02&L=spm&P=R282882&I=-1
Now, there are 3 constraints to this problem that have to be taken into
account:
(1) There must not be any empty regressors.
(2) For a valid parametric modulation there must be at least 2 trials
(= 2 onset) in the "onset" regressors (i.e. the one where you
specify your onsets)
(3) Contrast specification is tricky (Therefore I recommend a spcific
naming scheme for your regressors (and parametric modulations) which
you can later reference)
Here is some example code for solving problem (1) and (2) ...
In the following I assume that you save your onsets (and parametric
values) in a struct array called "onset". This struct array is setup
like this
onset(1).reg(1).name --- name for regressor 1 in session1
onset(1).reg(1).ons --- onsets for regressor 1 in session1
onset(1).reg(1).pm(1).name --- name for 1st parametric values for regressors 1 in session 1
onset(1).reg(1).pm(1).P --- 1st parametric values for regressors 1 in session 1
onset(1).reg(1).pm(2).name --- name for 2nd parametric values for regressors 1 in session 1
onset(1).reg(1).pm(2).P --- 2nd parametric values for regressors 1 in session 1
onset(1).reg(2).ons --- onsets for regressor 2 in session1
onset(1).reg(2).pm(1).name --- name for 1st parametric values for regressor2 in sesssion1
onset(1).reg(2).pm(1).P --- 1st parametric values for regressor2 in sesssion1
onset(1).reg(2).pm(2).name --- name for 2nd parametric values for regressor2 in sesssion1
onset(1).reg(2).pm(2).P --- 2nd parametric values for regressor2 in sesssion1
...
onset(2).reg(1).name --- name for regressor 1 in session2
onset(2).reg(1).ons --- onsets for regressor 1 in session2
onset(2).reg(1).pm(1).name --- name for 1st parametric values for regressors 1 in session 1
onset(2).reg(1).pm(1).P --- 1st parametric values for regressors 1 in session 1
onset(2).reg(1).pm(2).name --- name for 2nd parametric values for regressors 1 in session 1
onset(2).reg(1).pm(2).P --- 2nd parametric values for regressors 1 in session 1
onset(2).reg(2).ons --- onsets for regressor 2 in session2
onset(2).reg(2).pm(1).name --- name for 1st parametric values for regressor2 in sesssion1
onset(2).reg(2).pm(1).P --- 1st parametric values for regressor2 in sesssion1
onset(2).reg(2).pm(2).name --- name for 2nd parametric values for regressor2 in sesssion1
onset(2).reg(2).pm(2).P --- 2nd parametric values for regressor2 in sesssion1
...
and so on
If you are planning to look for session-specific contrast, you should
specify session-specific names for both regressors and parametric
modulation.
Now, in your batch file when you get to the point of specifying your
regressors (and parametric modulations) then you could something like
the following:
for s = 1:length(SPM.nscan) % loop over sessions
SPM.Sess(s).U = struct;
Rcnt = 1; % counter for the regressors
for r = 1:length(onset(s).reg) % loop over regressors
if ~isempty(onset(s).reg(Rcnt).ons) % test if onset vector is empty
SPM.Sess(s).U(Rcnt).name = {onset(s).reg(Rcnt).name}
SPM.Sess(s).U(Rcnt).ons = onset(s).reg(Ucnt).ons;
SPM.Sess(s).U(Rcnt).dur = 0; % if event-related
Pcnt = 1; % counter for parametric modulations
for p = 1:length(onset(s).reg(Rcnt).pm) % loop over parametric modulations
if length(onset(s).reg(Rcnt).ons > 1 % test if there are more than 1 onsets
SPM.Sess(s).U(Rcnt).P(Pcnt).name = onset(s).reg(Rcnt).pm(Pcnt).name;
SPM.Sess(s).U(Rcnt).P(Pcnt).P = onset(s).reg(Rcnt).pm(Pcnt).P;
SPM.Sess(s).U(Rcnt).P(Pcnt).h = 1; % first order polynomial expansion (linear)
Pcnt = Pcnt + 1;
end
Rcnt = Rcnt + 1;
end
end
% proceed with covariates, assuming there are none
SPM.Sess(s).C.C = [];
SPM.Sess(s).C.name = {};
end
This is it! No further code is necessary for specification of the
regressors. As you can see, it all comes down to saving your onsets and
parametric values in a orderly fashion that can be easily looped over.
If you use a conditional design specification like the one above, then
only those regressors are included that a not empty, and only those
parametric modulations are included for which the regressors has more
than one onset.
The batch file then continues with the specifications of global
normalization etc. and then calls spm_fmri_spm_ui and spm_spm.
Now a solution for problem (2) (contrast specification). Suppose you
want to specify two T-contrast, one "main effect" (where does reg1 show
a positive effect) and one differential contrast (where does reg1 show
a greater effect than reg2). I further assume that you want to average
over sessions and that therefore reg1 and reg2 iare the name for
regressors 1 and 2 IN ALL SESSIONS.
The strategy is to see if the regressors are present and then adjust the
contrast weight accordingly to make then available for a second level
analysis. This adjustment is important because the con-images contain
the sum of the contrast-weighted beta-images (and not the mean). In
matrix notation: con-value = beta*c'. This has important consequences.
Here is an example of the problem: Suppose subject 1 has onset for
regressor 1 in only 2 of 3 sessions, but subject 2 has onset for
regressors 2 in all 3 sessions. If you now put simple [1]-contrasts
weights on the respective columns of the design matrix, then subject 2
will get a much higher con-value because it contains the sum of
3 sessions, as opposed to subject 1 in which the con-value is only the
sum of 2 sessions. Therefore we need to adjust the contrast weights to
deal with this problem.
Here is some code for a batch file that will do this conditional
contrast specification:
% constrast 1: main effect of regressors 1
% define cell array holding contrast names
cTnames = {};
% define a contrast vector of zeros
cT(1,:) = zeros(1,size(SPM.xX.X,2); % every line of cT will contain a T-contrast vector
for c = 1:length(SPM.xX.name) % This cell array contain all the names specified earlier
if ~isempty(strfind(SPM.xX.name{c},onset(1).reg(1).ons))
cT(1,c) = 1;
end
end
% now adjust the contrast weights (if necessary)
if ~isempty(find(cT(1,:)) % regressors 1 is present in at least one session
cTnames{end+1} = 'ME: reg1'; % name of contrast, will show up on the contrast manager
idx = find(cT(1,:); % indices of 1's
l = length(idx); % number of 1's
cT(1,idx) = 1/l;
else
cT = []; % delete cT, if empty
end
This will yield contrast weights of 1/3 for subject 2, and 0.5 for
subject 1 in the example above. Therefore you will obtain properly
scaled con-images that you can test at the second level in a one-sample
t-test.
Now for contrast 2: differential effect of reg1 - reg2
cT(end+1,:) = zeros(1,size(SPM.xX.X,2); % 2nd contrast, but 1st could be empty, therefore cT(end,:)
for c = 1:length(SPM.xX.name)
if ~isempty(strfind(SPM.xX.name{c},onset(1).reg(1).name))
cT(end,c) = 1;
end
if ~isempty(strfind(SPM.xX.name{c},onset(1).reg(2).name))
cT(end,c) = -1;
end
end
% now adjust the contrast weight so that they add up to zero
if ~isempty(find(cT(end,:)==1)) & ~isempty(find(cT(end,:)==-1)) % are there any 1's and -1's in the contrast weights?
cTnames{end+1} = 'reg1-reg2';
pidx = find(cT(end,:)==1); % indices of 1's
pl = length(pidx); % number of 1's
nidx = find(cT(end,:)==-1); % indices of -1's
nl = length(nidx); % number of -1's
cT(end,pidx) = 1/pl; % scale 1's to 1/no. of 1's
cT(end,nidx) = -1/nl; % scale -1's to -1/no. of -1's
end
% finally we put these T-contrast in SPM.xCon and then estimate the
% contrasts
for f = 1:size(cT,1) % loop over all valid T-contrasts
c = cT(f,:); % current (scaled) contrast vector
cname = cTnames{f}; % current contrast name
spm_FcUtil('Set',cname,'T','c',c(:),SPM.xX.xKXs);
end
This applies only to T-contrast, or F-contrasts the case is a bit more
tricky.
I alos want to stress that I copied most of the code from some of my own
batch file and made some appropriate changes, but I have tested the code
above. There might be errors (e.g. typos); if so, please let me know
about it.
I hope this will give you some ideas about implementing a conditional
design specification.
Best,
Jan
On 2005-03-23 (Wed) at 10:01:22 +0100, Neggers, S.F.W. (Bas)
<[log in to unmask]> wrote:
> Hi Satra, List,
>
> I think it would be mathematically impossible to estimate the
> regression coefficients (Betas) for your Data Y (solving B= (X'X)^-1
> x X'Y ) for a design matrix X containing 3 identical regressors, i.e.
> in your case 2 regressors with a constant value for the 'empty'
> condition, plus the intercept. They are perfectly co-linear (=wrong!).
> So you will have to find another way to setup your GLM.
>
> Cheers,
>
> Bas
>
> -----Oorspronkelijk bericht----- Van: SPM (Statistical Parametric
> Mapping) [mailto:[log in to unmask]]Namens Stephen J. Fromm Verzonden:
> woensdag 23 maart 2005 1:40 Aan: [log in to unmask] Onderwerp: Re:
> [SPM] empty conditions in design
>
>
> On Tue, 22 Mar 2005 16:50:12 -0500, Satrajit Ghosh
> <[log in to unmask]> wrote:
>
> >Dear SPMusers,
> >
> >I have a study in which conditions are not balanced across runs.
> >Therefore, run1 may have conditions 1 and 4 and run2 may have
> >conditions 2 and 3. SPM complains if I set the onsets for
> >a condition to be empty. All it should translate to is an empty
> >column in the design matrix, which would be eliminated when the
> >matrix is reduced. The reason for entering empty conditions is that
> >I would like to say each run has all the conditions, but only enter
> >events for the conditions that actually exist in the run.
> >
> >Is this something that's part of the design intentions? Or am I doing
> >something incorrectly?
>
> As far as I can tell, SPM (SPM2 in my case) doesn't let one create
> empty columns.
>
> >
> >Thanks,
> >
> >Satra
>
--
Jan Gläscher Neuroimage Nord
+49-40-42803-7890 (office) Dept. of Neurology, Bldg S10
+49-40-42803-9955 (fax) University Hospital Hamburg-Eppendorf
[log in to unmask] Martinistr. 52
20246 Hamburg
Germany
http://www.uke.uni-hamburg.de/kliniken/neurologie/index_16969.php
---------------------------------------------------------------------
GnuPG/PGP key id: FEC4B55C
fingerprint: 5A36 1EF6 8472 117E 805A F240 3146 A410 FEC4 B55C
---------------------------------------------------------------------
|
