spm_config_factorial_design

 Configuration file for specification of factorial designs

 This function configures the design matrix (describing the general
 linear model), data specification, and other parameters necessary for
 the statistical analysis. These parameters are saved in a
 configuration file (SPM.mat) in the current directory, and are
 passed on to spm_spm.m (via the Estimate button) which estimates the design. 
 Inference on these estimated parameters is then handled by the SPM 
 results section.

 This function comprises two parts. The first defines 
 the user interface and the second sets up the necessary SPM structures.
 The second part has been largely cannibalised from spm_spm_ui.m which 
 we have retained for developmental continuity. 
 
 It has in common with spm_spm_ui.m its use of the I factor matrix, 
 the H,C,B,G design matrix partitions, the sF,sCFI,CFIforms,sCC,CCforms,
 sGXcalc,sGloNorm,sGMsca option definition variables and use of the 
 functions spm_DesMtx.m, spm_meanby.m and spm_non_sphericity.m.

 It departs from spm_spm_ui.m in that it does not use the design
 definition data structure D. Also, it uses the new SPM.factor field, 
 which for the case of full factorial designs, is used to automatically
 generate contrasts testing for main effects and interactions.

 This function departs from spm_spm_ui.m in that it does not provide 
 the same menu of design options (these were hardcoded in D). Instead it
 provides a number of options for simple designs (1) One-sample t-test,
 (2) Two-sample t-test, (3) Paired t-test and (4) Multiple regression.
 Two facilities are provided for specifying more complicated designs 
 (5) Full-factorial and (6) Flexible-factorial. These should be able to 
 specify all design options (and more) that were available in SPM2.
 For each of these design types one can additionally specify regressors using the
 `covariates' option.

 Options (5) and (6) differ in the 
 efficiency (eg. number of key strokes/button presses) with which a given
 design can be specified. For example, one-way ANOVAs can be specified
 using either option, but (5) is usually more efficient. 

 Full-factorial designs
 ______________________

 This option is best used when you wish to test for all
 main effects and interactions in one-way, two-way or three-way ANOVAs. 

 Design specification proceeds in 2 stages. Firstly, by creating new 
 factors and specifying the 
 number of levels and name for each. Nonsphericity, ANOVA-by-factor (for PET data) 
 and 
 scaling options (for PET data) can also be specified at this stage. Secondly, 
 scans are assigned separately to each cell. This accomodates unbalanced designs.

 For example, if you wish to test for a main effect in the population 
 from which your subjects are drawn
 and have modelled that effect at the first level using K basis functions
 (eg. K=3 informed basis functions) you can use a one-way ANOVA with K-levels. 
 Create a single factor with K levels and then assign the data to each
 cell eg. canonical, temporal derivative and dispersion derivative cells,
 where each cell is assigned scans from multiple subjects.

 SPM will automatically generate the contrasts necessary to test for all 
 main effects and interactions

 Flexible-factorial designs
 __________________________

 In this option the design matrix is created a block at a time. You can 
 decide whether you wish each block to be a main effect or a (two-way) 
 interaction. 
 
 This option is best used for one-way, two-way or 
 three-way ANOVAs but where you do not wish to test for all possible
 main effects and interactions. This is perhaps most useful for PET 
 where there is usually not enough data to test for all possible
 effects. Or for 3-way ANOVAs where you do not wish to test for all 
 of the two-way interactions. A typical example here would be a 
 group-by-drug-by-task analysis where, perhaps, only (i) group-by-drug or 
 (ii) group-by-task interactions are of interest. In this case it is only 
 necessary to have two-blocks in the design matrix - one for each
 interaction. The three-way interaction can then be tested for using a
 contrast that computes the difference between (i) and (ii).
 
 Design specification then proceeds in 3 stages. Firstly, factors
 are created and names specified for each. Nonsphericity, ANOVA-by-factor and 
 scaling options can also be specified at this stage. 

 Secondly, a list of
 scans is produced along with a factor matrix, I. This is an nscan x 4 matrix 
 of factor level indicators (see xX.I below). The first factor must be 
 'replication' but the other factors can be anything. Specification of I and 
 the scan list can be achieved in
 one of two ways (a) the 'Specify All' option allows I
 to be typed in at the user interface or (more likely) loaded in from the matlab
 workspace. All of the scans are then selected in one go. (b) the
 'Subjects' option allows you to enter scans a subject at a time. The
 corresponding experimental conditions (ie. levels of factors) are entered
 at the same time. SPM will then create the factor matrix I. This style of
 interface is similar to that available in SPM2.

 Thirdly, the design matrix is built up a block at a time. Each block 
 can be a main effect or a (two-way) interaction.  
 

 ----------------------------------------------------------------------

 Variables saved in the SPM stucture

 xY.VY         - nScan x 1 struct array of memory mapped images
                 (see spm_vol for definition of the map structure)
 xX            - structure describing design matrix
 xX.I          - nScan x 4 matrix of factor level indicators
                 I(n,i) is the level of factor i corresponding to image n
 xX.sF         - 1x4 cellstr containing the names of the four factors
                 xX.sF{i} is the name of factor i
 xX.X          - design matrix
 xX.xVi        - correlation constraints for non-spericity correction
 xX.iH         - vector of H partition (condition effects) indices,
                 identifying columns of X correspoding to H
 xX.iC         - vector of C partition (covariates of interest) indices
 xX.iB         - vector of B partition (block effects) indices
 xX.iG         - vector of G partition (nuisance variables) indices
 xX.name     - p x 1 cellstr of effect names corresponding to columns
                 of the design matrix
 
 xC            - structure array of covariate details
 xC(i).rc      - raw (as entered) i-th covariate
 xC(i).rcname  - name of this covariate (string)
 xC(i).c       - covariate as appears in design matrix (after any scaling,
                 centering of interactions)
 xC(i).cname   - cellstr containing names for effects corresponding to
                 columns of xC(i).c
 xC(i).iCC     - covariate centering option
 xC(i).iCFI    - covariate by factor interaction option
 xC(i).type    - covariate type: 1=interest, 2=nuisance, 3=global
 xC(i).cols    - columns of design matrix corresponding to xC(i).c
 xC(i).descrip - cellstr containing a description of the covariate
 
 xGX           - structure describing global options and values
 xGX.iGXcalc   - global calculation option used
 xGX.sGXcalc   - string describing global calculation used
 xGX.rg        - raw globals (before scaling and such like)
 xGX.iGMsca    - grand mean scaling option
 xGX.sGMsca    - string describing grand mean scaling
 xGX.GM        - value for grand mean (/proportional) scaling
 xGX.gSF       - global scaling factor (applied to xGX.rg)
 xGX.iGC       - global covariate centering option
 xGX.sGC       - string describing global covariate centering option
 xGX.gc        - center for global covariate
 xGX.iGloNorm  - Global normalisation option
 xGX.sGloNorm  - string describing global normalisation option
 
 xM            - structure describing masking options
 xM.T          - Threshold masking value (-Inf=>None,
                 real=>absolute, complex=>proportional (i.e. times global) )
 xM.TH         - nScan x 1 vector of analysis thresholds, one per image
 xM.I          - Implicit masking (0=>none, 1=>implicit zero/NaN mask)
 xM.VM         - struct array of explicit mask images
                 (empty if no explicit masks)
 xM.xs         - structure describing masking options
                 (format is same as for xsDes described below)
 
 xsDes         - structure of strings describing the design:
                 Fieldnames are essentially topic strings (use "_"'s for
                 spaces), and the field values should be strings or cellstr's
                 of information regarding that topic. spm_DesRep.m
                 uses this structure to produce a printed description
                 of the design, displaying the fieldnames (with "_"'s 
                 converted to spaces) in bold as topics, with
                 the corresponding text to the right% 

_______________________________________________________________________
 Copyright (C) 2005 Wellcome Department of Imaging Neuroscience

0001 function conf = spm_config_factorial_design
0002 % Configuration file for specification of factorial designs
0003 %
0004 % This function configures the design matrix (describing the general
0005 % linear model), data specification, and other parameters necessary for
0006 % the statistical analysis. These parameters are saved in a
0007 % configuration file (SPM.mat) in the current directory, and are
0008 % passed on to spm_spm.m (via the Estimate button) which estimates the design.
0009 % Inference on these estimated parameters is then handled by the SPM
0010 % results section.
0011 %
0012 % This function comprises two parts. The first defines
0013 % the user interface and the second sets up the necessary SPM structures.
0014 % The second part has been largely cannibalised from spm_spm_ui.m which
0015 % we have retained for developmental continuity.
0016 %
0017 % It has in common with spm_spm_ui.m its use of the I factor matrix,
0018 % the H,C,B,G design matrix partitions, the sF,sCFI,CFIforms,sCC,CCforms,
0019 % sGXcalc,sGloNorm,sGMsca option definition variables and use of the
0020 % functions spm_DesMtx.m, spm_meanby.m and spm_non_sphericity.m.
0021 %
0022 % It departs from spm_spm_ui.m in that it does not use the design
0023 % definition data structure D. Also, it uses the new SPM.factor field,
0024 % which for the case of full factorial designs, is used to automatically
0025 % generate contrasts testing for main effects and interactions.
0026 %
0027 % This function departs from spm_spm_ui.m in that it does not provide
0028 % the same menu of design options (these were hardcoded in D). Instead it
0029 % provides a number of options for simple designs (1) One-sample t-test,
0030 % (2) Two-sample t-test, (3) Paired t-test and (4) Multiple regression.
0031 % Two facilities are provided for specifying more complicated designs
0032 % (5) Full-factorial and (6) Flexible-factorial. These should be able to
0033 % specify all design options (and more) that were available in SPM2.
0034 % For each of these design types one can additionally specify regressors using the
0035 % `covariates' option.
0036 %
0037 % Options (5) and (6) differ in the
0038 % efficiency (eg. number of key strokes/button presses) with which a given
0039 % design can be specified. For example, one-way ANOVAs can be specified
0040 % using either option, but (5) is usually more efficient.
0041 %
0042 % Full-factorial designs
0043 % ______________________
0044 %
0045 % This option is best used when you wish to test for all
0046 % main effects and interactions in one-way, two-way or three-way ANOVAs.
0047 %
0048 % Design specification proceeds in 2 stages. Firstly, by creating new
0049 % factors and specifying the
0050 % number of levels and name for each. Nonsphericity, ANOVA-by-factor (for PET data)
0051 % and
0052 % scaling options (for PET data) can also be specified at this stage. Secondly,
0053 % scans are assigned separately to each cell. This accomodates unbalanced designs.
0054 %
0055 % For example, if you wish to test for a main effect in the population
0056 % from which your subjects are drawn
0057 % and have modelled that effect at the first level using K basis functions
0058 % (eg. K=3 informed basis functions) you can use a one-way ANOVA with K-levels.
0059 % Create a single factor with K levels and then assign the data to each
0060 % cell eg. canonical, temporal derivative and dispersion derivative cells,
0061 % where each cell is assigned scans from multiple subjects.
0062 %
0063 % SPM will automatically generate the contrasts necessary to test for all
0064 % main effects and interactions
0065 %
0066 % Flexible-factorial designs
0067 % __________________________
0068 %
0069 % In this option the design matrix is created a block at a time. You can
0070 % decide whether you wish each block to be a main effect or a (two-way)
0071 % interaction.
0072 %
0073 % This option is best used for one-way, two-way or
0074 % three-way ANOVAs but where you do not wish to test for all possible
0075 % main effects and interactions. This is perhaps most useful for PET
0076 % where there is usually not enough data to test for all possible
0077 % effects. Or for 3-way ANOVAs where you do not wish to test for all
0078 % of the two-way interactions. A typical example here would be a
0079 % group-by-drug-by-task analysis where, perhaps, only (i) group-by-drug or
0080 % (ii) group-by-task interactions are of interest. In this case it is only
0081 % necessary to have two-blocks in the design matrix - one for each
0082 % interaction. The three-way interaction can then be tested for using a
0083 % contrast that computes the difference between (i) and (ii).
0084 %
0085 % Design specification then proceeds in 3 stages. Firstly, factors
0086 % are created and names specified for each. Nonsphericity, ANOVA-by-factor and
0087 % scaling options can also be specified at this stage.
0088 %
0089 % Secondly, a list of
0090 % scans is produced along with a factor matrix, I. This is an nscan x 4 matrix
0091 % of factor level indicators (see xX.I below). The first factor must be
0092 % 'replication' but the other factors can be anything. Specification of I and
0093 % the scan list can be achieved in
0094 % one of two ways (a) the 'Specify All' option allows I
0095 % to be typed in at the user interface or (more likely) loaded in from the matlab
0096 % workspace. All of the scans are then selected in one go. (b) the
0097 % 'Subjects' option allows you to enter scans a subject at a time. The
0098 % corresponding experimental conditions (ie. levels of factors) are entered
0099 % at the same time. SPM will then create the factor matrix I. This style of
0100 % interface is similar to that available in SPM2.
0101 %
0102 % Thirdly, the design matrix is built up a block at a time. Each block
0103 % can be a main effect or a (two-way) interaction.
0104 %
0105 %
0106 % ----------------------------------------------------------------------
0107 %
0108 % Variables saved in the SPM stucture
0109 %
0110 % xY.VY         - nScan x 1 struct array of memory mapped images
0111 %                 (see spm_vol for definition of the map structure)
0112 % xX            - structure describing design matrix
0113 % xX.I          - nScan x 4 matrix of factor level indicators
0114 %                 I(n,i) is the level of factor i corresponding to image n
0115 % xX.sF         - 1x4 cellstr containing the names of the four factors
0116 %                 xX.sF{i} is the name of factor i
0117 % xX.X          - design matrix
0118 % xX.xVi        - correlation constraints for non-spericity correction
0119 % xX.iH         - vector of H partition (condition effects) indices,
0120 %                 identifying columns of X correspoding to H
0121 % xX.iC         - vector of C partition (covariates of interest) indices
0122 % xX.iB         - vector of B partition (block effects) indices
0123 % xX.iG         - vector of G partition (nuisance variables) indices
0124 % xX.name     - p x 1 cellstr of effect names corresponding to columns
0125 %                 of the design matrix
0126 %
0127 % xC            - structure array of covariate details
0128 % xC(i).rc      - raw (as entered) i-th covariate
0129 % xC(i).rcname  - name of this covariate (string)
0130 % xC(i).c       - covariate as appears in design matrix (after any scaling,
0131 %                 centering of interactions)
0132 % xC(i).cname   - cellstr containing names for effects corresponding to
0133 %                 columns of xC(i).c
0134 % xC(i).iCC     - covariate centering option
0135 % xC(i).iCFI    - covariate by factor interaction option
0136 % xC(i).type    - covariate type: 1=interest, 2=nuisance, 3=global
0137 % xC(i).cols    - columns of design matrix corresponding to xC(i).c
0138 % xC(i).descrip - cellstr containing a description of the covariate
0139 %
0140 % xGX           - structure describing global options and values
0141 % xGX.iGXcalc   - global calculation option used
0142 % xGX.sGXcalc   - string describing global calculation used
0143 % xGX.rg        - raw globals (before scaling and such like)
0144 % xGX.iGMsca    - grand mean scaling option
0145 % xGX.sGMsca    - string describing grand mean scaling
0146 % xGX.GM        - value for grand mean (/proportional) scaling
0147 % xGX.gSF       - global scaling factor (applied to xGX.rg)
0148 % xGX.iGC       - global covariate centering option
0149 % xGX.sGC       - string describing global covariate centering option
0150 % xGX.gc        - center for global covariate
0151 % xGX.iGloNorm  - Global normalisation option
0152 % xGX.sGloNorm  - string describing global normalisation option
0153 %
0154 % xM            - structure describing masking options
0155 % xM.T          - Threshold masking value (-Inf=>None,
0156 %                 real=>absolute, complex=>proportional (i.e. times global) )
0157 % xM.TH         - nScan x 1 vector of analysis thresholds, one per image
0158 % xM.I          - Implicit masking (0=>none, 1=>implicit zero/NaN mask)
0159 % xM.VM         - struct array of explicit mask images
0160 %                 (empty if no explicit masks)
0161 % xM.xs         - structure describing masking options
0162 %                 (format is same as for xsDes described below)
0163 %
0164 % xsDes         - structure of strings describing the design:
0165 %                 Fieldnames are essentially topic strings (use "_"'s for
0166 %                 spaces), and the field values should be strings or cellstr's
0167 %                 of information regarding that topic. spm_DesRep.m
0168 %                 uses this structure to produce a printed description
0169 %                 of the design, displaying the fieldnames (with "_"'s
0170 %                 converted to spaces) in bold as topics, with
0171 %                 the corresponding text to the right%
0172 %
0173 %_______________________________________________________________________
0174 % Copyright (C) 2005 Wellcome Department of Imaging Neuroscience
0175 
0176 % Will Penny
0177 % $Id: spm_config_factorial_design.m 411 2006-01-20 15:45:10Z john $
0178 
0179 % Define inline types.
0180 %-----------------------------------------------------------------------
0181 
0182 entry = inline(['struct(''type'',''entry'',''name'',name,'...
0183     '''tag'',tag,''strtype'',strtype,''num'',num,''help'',hlp)'],...
0184     'name','tag','strtype','num','hlp');
0185 
0186 files = inline(['struct(''type'',''files'',''name'',name,'...
0187     '''tag'',tag,''filter'',fltr,''num'',num,''help'',hlp)'],...
0188     'name','tag','fltr','num','hlp');
0189 
0190 mnu = inline(['struct(''type'',''menu'',''name'',name,'...
0191     '''tag'',tag,''labels'',{labels},''values'',{values},''help'',hlp)'],...
0192     'name','tag','labels','values','hlp');
0193 
0194 branch = inline(['struct(''type'',''branch'',''name'',name,'...
0195     '''tag'',tag,''val'',{val},''help'',hlp)'],...
0196     'name','tag','val','hlp');
0197 
0198 repeat = inline(['struct(''type'',''repeat'',''name'',name,'...
0199     '''tag'',tag,''values'',{values},''help'',hlp)'],...
0200     'name','tag','values','hlp');
0201 
0202 choice = inline(['struct(''type'',''choice'',''name'',name,'...
0203     '''tag'',tag,''values'',{values},''help'',hlp)'],...
0204     'name','tag','values','hlp');
0205 
0206 %-----------------------------------------------------------------------
0207 
0208 sp_text = ['                                                      ',...
0209       '                                                      '];
0210 
0211 %-------------------------------------------------------------------------
0212 % Covariates for general use
0213 
0214 iCFI    = mnu('Interactions','iCFI',{'None','With Factor 1',...
0215         'With Factor 2','With Factor 3'},{1,2,3,4},'');
0216 iCFI.val={1};
0217 p1 = ['For each covariate you have defined, there is an opportunity to ',...
0218       'create an additional regressor that is the interaction between the ',...
0219       'covariate and a chosen experimental factor. '];
0220 iCFI.help={p1,sp_text};
0221 
0222 iCC    = mnu('Centering','iCC',{'Overall mean','Factor 1 mean',...
0223         'Factor 2 mean','Factor 3 mean','No centering',...
0224         'User specifified value','As implied by ANCOVA',...
0225         'GM'},{1,2,3,4,5,6,7,8},'');
0226 iCC.val={1};
0227 p1 = ['The appropriate centering option is usually the one that ',...
0228       'corresponds to the interaction chosen, and ensures that main ',... 
0229       'effects of the interacting factor aren''t affected by the covariate. ',... 
0230       'You are advised to choose this option, unless you have other ',... 
0231       'modelling considerations. '];
0232 iCC.help={p1,sp_text};
0233 
0234 cname      = entry('Name','cname','s', [1 Inf],'Name of covariate');
0235 
0236 c   = entry('Vector','c','e',[Inf 1],'Vector of covariate values');
0237 
0238 cov  = branch('Covariate','cov',{c,cname,iCFI,iCC},'Covariate');
0239 
0240 cov.help = {'Add a new covariate to your experimental design'};
0241 
0242 covs = repeat('Covariates','covariates',{cov},'');
0243 p1 = ['This option allows for the specification of covariates and ',...
0244       'nuisance variables. Unlike SPM94/5/6, where the design was ',...
0245       'partitioned into effects of interest and nuisance effects ',... 
0246       'for the computation of adjusted data and the F-statistic ',... 
0247       '(which was used to thresh out voxels where there appeared to ',... 
0248       'be no effects of interest), SPM5 does not partition the design ',... 
0249       'in this way. The only remaining distinction between effects of ',... 
0250       'interest (including covariates) and nuisance effects is their ',... 
0251       'location in the design matrix, which we have retained for ',... 
0252       'continuity.  Pre-specified design matrix partitions can be entered. ']; 
0253       
0254 covs.help={p1,sp_text};
0255 
0256 %-------------------------------------------------------------------------
0257 % Covariates for multiple regression
0258 
0259 mcov  = branch('Covariate','mcov',{c,cname},'Covariate');
0260 
0261 mcov.help = {'Add a new covariate to your experimental design'};
0262 
0263 mcovs = repeat('Covariates','covariates',{mcov},'Covariates');
0264 
0265 %-----------------------------------------------------------------------
0266 %  Specify names of factors, numbers of levels and statistical dependencies
0267 
0268 name     = entry('Name','name','s',[1 Inf],'');
0269 name.val = {'Covariate'};
0270 p1=['Enter name of covariate eg. reaction time'];
0271 name.help={p1};
0272 
0273 val      = entry('Value','val','e',[Inf 1],'');
0274 val.help={['Enter the vector of covariate values']};
0275 
0276 variance    = mnu('Variance','variance',{'Equal','Unequal'},{0,1},'');
0277 variance.val={1};
0278 p1=['By default, the measurements in each level are assumed to have unequal variance. '];
0279 p2=['This violates the assumption of ''sphericity'' and is therefore an example of ',...
0280     '''non-sphericity''.'];
0281 p3=['This can occur, for example, in a 2nd-level analysis of variance, one ',...
0282     'contrast may be scaled differently from another.  Another example would ',...
0283     'be the comparison of qualitatively different dependent variables ',...
0284     '(e.g. normals vs. patients).  Different variances ',...
0285     '(heteroscedasticy) induce different error covariance components that ',...
0286     'are estimated using restricted maximum likelihood (see below).'];
0287 p4=['Restricted Maximum Likelihood (REML): The ensuing covariance components ',...
0288      'will be estimated using ReML in spm_spm (assuming the same for all ',...
0289      'responsive voxels) and used to adjust the ',...
0290     'statistics and degrees of freedom during inference. By default spm_spm ',...
0291     'will use weighted least squares to produce Gauss-Markov or Maximum ',...
0292     'likelihood estimators using the non-sphericity structure specified at this ',...
0293     'stage. The components will be found in xX.xVi and enter the estimation ',...
0294     'procedure exactly as the serial correlations in fMRI models.'];
0295 variance.help = {p1,sp_text,p2,sp_text,p3,sp_text,p4,sp_text};
0296 
0297 dept    = mnu('Independence','dept',{'Yes','No'},{0,1},'');
0298 dept.val={0};
0299 p1=['By default, the measurements are assumed to be independent between levels. '];
0300 p2=['If you change this option to allow for dependencies, this will violate ',...
0301     'the assumption of sphericity. It would therefore be an example ',...
0302     'of non-sphericity. One such example would be where you had repeated ',...
0303     'measurements from the same subjects - it may then be the case that, over ',...
0304     'subjects, measure 1 is correlated to measure 2. '];
0305 
0306 dept.help = {p1,sp_text,p2,sp_text,p4,sp_text};
0307 
0308 fname.type    = 'entry';
0309 fname.name    = 'Name';
0310 fname.tag     = 'name';
0311 fname.strtype = 's';
0312 fname.num     = [1 1];
0313 fname.help    = {'Name of factor, eg. ''Repetition'' '};
0314 
0315 levels = entry('Levels','levels','e',[Inf 1],''); 
0316 p1=['Enter number of levels for this factor, eg. 2'];
0317 levels.help ={p1};
0318 
0319 gmsca = mnu('Grand mean scaling','gmsca',...
0320     {'No','Yes'},{0,1},'');
0321 gmsca.val={0};
0322 p0=['This option is only used for PET data.'];
0323 
0324 p1=['Selecting YES will specify ''grand mean scaling by factor'' which could ',...
0325     'be eg. ''grand mean scaling by subject'' if the factor is ''subject''. '];
0326 
0327 p2 =['Since differences between subjects may be due to gain and sensitivity ',...
0328     'effects, AnCova by subject could be combined with "grand mean scaling ',...
0329     'by subject" to obtain a combination of between subject proportional ',...
0330     'scaling and within subject AnCova. '];
0331 gmsca.help={p0,sp_text,p1,sp_text,p2,sp_text};
0332 
0333 ancova = mnu('ANCOVA','ancova',...
0334     {'No','Yes'},{0,1},'');
0335 ancova.val={0};
0336 
0337 p1=['Selecting YES will specify ''ANCOVA-by-factor'' regressors. ',...
0338      'This includes eg. ''Ancova by subject'' or ''Ancova by effect''. ',... 
0339      'These options allow eg. different subjects ',...
0340     'to have different relationships between local and global measurements. '];
0341 ancova.help={p0,sp_text,p1,sp_text};
0342 
0343 factor.type   = 'branch';
0344 factor.name   = 'Factor';
0345 factor.tag    = 'fact';
0346 factor.val    = {fname,levels,dept,variance,gmsca,ancova};
0347 factor.help = {'Add a new factor to your experimental design'};
0348 
0349 factors.type = 'repeat';
0350 factors.name = 'Factors';
0351 factors.tag  = 'factors';
0352 factors.values = {factor};
0353 factors.num  = [1 Inf];
0354 p1 = ['Specify your design a factor at a time. '];
0355 factors.help ={p1,sp_text};
0356 
0357 %-------------------------------------------------------------------------
0358 % Associate each cell in factorial design with a list of images
0359 
0360 scans    = files('Scans','scans','image',[1 Inf],'Select scans');
0361 scans.help = {[...
0362 'Select the images for this cell.  They must all have the same ',...
0363 'image dimensions, orientation, voxel size etc.']};
0364 
0365 levels = entry('Levels','levels','e',[Inf 1],''); 
0366 p1=['Enter a vector or scalar that specifies which cell in the factorial ',...
0367         'design these images belong to. The length of this vector should '...
0368         'correspond to the number of factors in the design'];
0369 p2=['For example, length 2 vectors should be used for two-factor designs ',...
0370         'eg. the vector [2 3] specifies the cell corresponding to the 2nd-level of the first ',...
0371         'factor and the 3rd level of the 2nd factor.'];
0372 levels.help ={p1,sp_text,p2,sp_text};
0373 
0374 icell.type   = 'branch';
0375 icell.name   = 'Cell';
0376 icell.tag    = 'icell';
0377 icell.val    = {levels,scans};
0378 icell.help = {'Enter data for a cell in your design'};
0379 
0380 cells.type = 'repeat';
0381 cells.name = 'Specify cells';
0382 cells.tag  = 'cells';
0383 cells.values = {icell};
0384 cells.num  = [1 Inf];
0385 p1 = ['Enter the scans a cell at a time'];
0386 cells.help={p1,sp_text};
0387 
0388 %-------------------------------------------------------------------------
0389 % Create a design block-by-block
0390 
0391 fac.type   = 'branch';
0392 fac.name   = 'Factor';
0393 fac.tag    = 'fac';
0394 fac.val    = {fname,dept,variance,gmsca,ancova};
0395 p1=['Add a new factor to your design.'];
0396 p2=['If you are using the ''Subjects'' option to specify your scans ',...
0397         'and conditions, you may wish to make use of the following facility. ',...
0398         'There are two reserved words for the names of factors. These are ',...
0399         '''subject'' and ''repl'' (standing for replication). If you use these ',...
0400         'factor names then SPM can automatically create replication and/or ',...
0401         'subject factors without you having to type in an extra entry in the ',...
0402         'condition vector.'];
0403 p3=['For example, if you wish to model Subject and Task effects (two factors), ',...
0404      'under Subjects->Subject->Conditions you can type in simply ',...
0405     '[1 2 1 2] to specify eg. just the ''Task'' factor level. You do not need to ',...
0406     'eg. for the 4th subject enter the matrix [1 4; 2 4; 1 4; 2 4]. '];
0407 
0408 fac.help = {p1,sp_text,p2,sp_text,p3,sp_text};
0409 
0410 facs.type = 'repeat';
0411 facs.name = 'Factors';
0412 facs.tag  = 'facs';
0413 facs.values = {fac};
0414 facs.num  = [1 Inf];
0415 p1=['Specify your design a factor at a time.'];
0416 facs.help={p1,sp_text};
0417 
0418 fnum = entry('Factor number','fnum','e',[Inf 1],''); 
0419 p1=['Enter the number of the factor.'];
0420 fnum.help ={p1};
0421 
0422 fnums = entry('Factor numbers','fnums','e',[2 1],''); 
0423 p1=['Enter the numbers of the factors of this (two-way) interaction.'];
0424 fnums.help ={p1};
0425 
0426 fmain.type   = 'branch';
0427 fmain.name   = 'Main effect';
0428 fmain.tag    = 'fmain';
0429 fmain.val    = {fnum};
0430 fmain.help = {'Add a main effect to your design matrix'};
0431 
0432 inter.type   = 'branch';
0433 inter.name   = 'Interaction';
0434 inter.tag    = 'inter';
0435 inter.val    = {fnums};
0436 inter.help = {'Add an interaction to your design matrix'};
0437 
0438 maininters.type = 'repeat';
0439 maininters.name = 'Main effects & Interactions';
0440 maininters.num  = [1 Inf];
0441 maininters.tag  = 'maininters';
0442 maininters.values = {fmain, inter};
0443 
0444 scans    = files('Scans','scans','image',[1 Inf],'Select scans');
0445 scans.help = {[...
0446 'Select the images to be analysed.  They must all have the same ',...
0447 'image dimensions, orientation, voxel size etc.']};
0448 
0449 imatrix   = entry('Factor matrix','imatrix','e',[Inf Inf],'');
0450 imatrix.val = {0};
0451 
0452 specall.type   = 'branch';
0453 specall.name   = 'Specify all';
0454 specall.tag    = 'specall';
0455 specall.val    = {scans,imatrix};
0456 specall.help = {'Specify (i) all scans in one go and (ii) all conditions using a factor matrix'};
0457 
0458 conds   = entry('Conditions','conds','e',[Inf Inf],'');
0459 
0460 fsubject.type   = 'branch';
0461 fsubject.name   = 'Subject';
0462 fsubject.tag    = 'fsubject';
0463 fsubject.val    = {scans,conds};
0464 fsubject.help = {'Enter data and conditions for a new subject'};
0465 
0466 fsubjects.type = 'repeat';
0467 fsubjects.name = 'Subjects';
0468 fsubjects.tag  = 'fsubjects';
0469 fsubjects.values = {fsubject};
0470 fsubjects.num  = [1 Inf];
0471 
0472 fsuball.type = 'choice';
0473 fsuball.name = 'Specify Subjects or all Scans & Factors';
0474 fsuball.tag  = 'fsuball';
0475 fsuball.values = {fsubjects specall};
0476 
0477 %-------------------------------------------------------------------------
0478 % Two-sample t-test
0479 
0480 scans1    = files('Group 1 scans','scans1','image',[1 Inf],'Select scans');
0481 scans1.help = {[...
0482 'Select the images from sample 1.  They must all have the same ',...
0483 'image dimensions, orientation, voxel size etc.']};
0484 
0485 scans2    = files('Group 2 scans','scans2','image',[1 Inf],'Select scans');
0486 scans2.help = {[...
0487 'Select the images from sample 2.  They must all have the same ',...
0488 'image dimensions, orientation, voxel size etc.']};
0489 
0490 %-------------------------------------------------------------------------
0491 % Paired t-test
0492 
0493 pscans    = files('Scans [1,2]','scans','image',[1 2],'Select scans');
0494 pscans.help = {[...
0495 'Select the pair of images. ']};
0496 
0497 pair.type   = 'branch';
0498 pair.name   = 'Pair';
0499 pair.tag    = 'pair';
0500 pair.val    = {pscans};
0501 pair.help = {'Add a new pair of scans to your experimental design'};
0502 
0503 pairs.type = 'repeat';
0504 pairs.name = 'Pairs';
0505 pairs.tag  = 'pairs';
0506 pairs.num  = [1 Inf];
0507 pairs.values = {pair};
0508 p1 = [' ',...
0509       ' '];
0510 pairs.help ={p1,sp_text};
0511 
0512 %-------------------------------------------------------------------------
0513 % One sample t-test
0514 
0515 t1scans    = files('Scans','scans','image',[1 Inf],'Select scans');
0516 t1scans.help = {[...
0517 'Select the images.  They must all have the same ',...
0518 'image dimensions, orientation, voxel size etc.']};
0519 
0520 %-------------------------------------------------------------------------
0521 % Design menu
0522 
0523 t1 = branch('One-sample t-test','t1',...
0524     {t1scans},'');
0525 
0526 t2 = branch('Two-sample t-test','t2',...
0527     {scans1,scans2,dept,variance,gmsca,ancova},'');
0528 
0529 pt = branch('Paired t-test','pt',...
0530     {pairs,dept,variance,gmsca,ancova},'');
0531 
0532 mreg = branch('Multiple regression','mreg',...
0533     {t1scans,mcovs},'');
0534 
0535 fblock = branch('Flexible factorial','fblock',...
0536     {facs,fsuball,maininters},'');
0537 pb1 = ['Create a design matrix a block at a time by specifying which ',...
0538       'main effects and interactions you wish to be included.'];
0539 
0540 
0541 pb2=['This option is best used for one-way, two-way or ',... 
0542 'three-way ANOVAs but where you do not wish to test for all possible ',...
0543 'main effects and interactions. This is perhaps most useful for PET ',... 
0544 'where there is usually not enough data to test for all possible ',...
0545 'effects. Or for 3-way ANOVAs where you do not wish to test for all ',... 
0546 'of the two-way interactions. A typical example here would be a ',... 
0547 'group-by-drug-by-task analysis where, perhaps, only (i) group-by-drug or ',... 
0548 '(ii) group-by-task interactions are of interest. In this case it is only ',... 
0549 'necessary to have two-blocks in the design matrix - one for each ',...
0550 'interaction. The three-way interaction can then be tested for using a ',...
0551 'contrast that computes the difference between (i) and (ii).'];
0552  
0553 pb3=['Design specification then proceeds in 3 stages. Firstly, factors ',...
0554 'are created and names specified for each. Nonsphericity, ANOVA-by-factor and ',...
0555 'scaling options can also be specified at this stage.']; 
0556 
0557 pb4=['Secondly, a list of ',...
0558 'scans is produced along with a factor matrix, I. This is an nscan x 4 matrix ',... 
0559 'of factor level indicators (see xX.I below). The first factor must be ',... 
0560 '''replication'' but the other factors can be anything. Specification of I and ',...
0561 'the scan list can be achieved in ',...
0562 'one of two ways (a) the ''Specify All'' option allows I ',...
0563 'to be typed in at the user interface or (more likely) loaded in from the matlab ',...
0564 'workspace. All of the scans are then selected in one go. (b) the ',...
0565 '''Subjects'' option allows you to enter scans a subject at a time. The ',...
0566 'corresponding experimental conditions (ie. levels of factors) are entered ',...
0567 'at the same time. SPM will then create the factor matrix I. This style of ',...
0568 'interface is similar to that available in SPM2.'];
0569 
0570 pb5=['Thirdly, the design matrix is built up a block at a time. Each block ',... 
0571 'can be a main effect or a (two-way) interaction. ']; 
0572 fblock.help={pb1,sp_text,pb2,sp_text,pb3,sp_text,pb4,sp_text,pb5,sp_text};
0573 
0574 fd = branch('Full factorial','fd',...
0575     {factors,cells},'');
0576 pfull1=['This option is best used when you wish to test for all ',...
0577 'main effects and interactions in one-way, two-way or three-way ANOVAs. ',...
0578 'Design specification proceeds in 2 stages. Firstly, by creating new ',...
0579 'factors and specifying the ',...
0580 'number of levels and name for each. Nonsphericity, ANOVA-by-factor and ',...
0581 'scaling options can also be specified at this stage. Secondly, scans are ',...
0582 'assigned separately to each cell. This accomodates unbalanced designs.'];
0583 
0584 pfull2=['For example, if you wish to test for a main effect in the population ',... 
0585 'from which your subjects are drawn ',...
0586 'and have modelled that effect at the first level using K basis functions ',...
0587 '(eg. K=3 informed basis functions) you can use a one-way ANOVA with K-levels. ',... 
0588 'Create a single factor with K levels and then assign the data to each ',...
0589 'cell eg. canonical, temporal derivative and dispersion derivative cells, ',...
0590 'where each cell is assigned scans from multiple subjects.'];
0591 
0592 pfull3 = ['SPM will also automatically generate ',...
0593       'the contrasts necessary to test for all main effects and interactions. '];
0594 
0595 fd.help={pfull1,sp_text,pfull2,sp_text,pfull3,sp_text};
0596 
0597 
0598 
0599 des = choice('Design','des',...
0600     {t1,t2,pt,mreg,fd,fblock},'');
0601 
0602 
0603 %-------------------------------------------------------------------------
0604 % Masking options
0605 
0606 im    = mnu('Implicit Mask','im',{'Yes','No'},{1,0},'');
0607 im.val={1};
0608 p1=['An "implicit mask" is a mask implied by a particular voxel ',...
0609     'value. Voxels with this mask value are excluded from the ',...
0610     'analysis. '];
0611 p2=['For image data-types with a representation of NaN ',...
0612     '(see spm_type.m), NaN''s is the implicit mask value, (and ',...
0613     'NaN''s are always masked out). '];
0614 p3=['For image data-types without a representation of NaN, zero is ',...
0615     'the mask value, and the user can choose whether zero voxels ',...
0616     'should be masked out or not.'];
0617 p4=['By default, an implicit mask is used. '];
0618 im.help = {p1,sp_text,p2,sp_text,p3,sp_text,p4,sp_text};
0619 
0620 em = files('Explicit Mask','em','image',[0 1],'');
0621 em.val={''};
0622 em.help = {['Select an explicit mask ']};
0623 p1=['Explicit masks are other images containing (implicit) masks ',...
0624     'that are to be applied to the current analysis.'];
0625 p2=['All voxels with value NaN (for image data-types with a ',...
0626     'representation of NaN), or zero (for other data types) are ',...
0627     'excluded from the analysis. '];
0628 p3=['Explicit mask images can have any orientation and voxel/image ',...
0629     'size. Nearest neighbour interpolation of a mask image is used if ',...
0630     'the voxel centers of the input images do not coincide with that ',...
0631     'of the mask image.'];
0632 em.help = {p1,sp_text,p2,sp_text,p3,sp_text};
0633 
0634 tm_none.type = 'const';
0635 tm_none.name = 'None';
0636 tm_none.tag = 'tm_none';
0637 tm_none.val = {[]};
0638 tm_none.help = {'No threshold masking'};
0639 
0640 athresh = entry('Threshold','athresh','e',[1 1],'');
0641 athresh.val={100};
0642 p1=['Enter the absolute value of the threshold.'];
0643 athresh.help = {p1,sp_text};
0644 
0645 tma = branch('Absolute','tma',...
0646     {athresh},'');
0647 p1=['Images are thresholded at a given value and ',...
0648     'only voxels at which all images exceed the threshold are included. '];
0649 p2=['This option allows you to specify the absolute value of the threshold.'];
0650 tma.help = {p1,sp_text,p2,sp_text};
0651 
0652 p2=['By default, Relative Threshold Masking is turned off. '];
0653 rselect.help = {p1,sp_text,p2,sp_text};
0654 
0655 rthresh = entry('Threshold','rthresh','e',[1 1],'');
0656 rthresh.val={0.8};
0657 p1=['Enter the threshold as a proportion of the global value'];
0658 rthresh.help = {p1,sp_text};
0659 
0660 tmr = branch('Relative','tmr',...
0661     {rthresh},'');
0662 p1=['Images are thresholded at a given value and ',...
0663     'only voxels at which all images exceed the threshold are included. '];
0664 p2=['This option allows you to specify the value of the threshold ',...
0665     'as a proportion of the global value. '];
0666 tmr.help = {p1,sp_text,p2,sp_text};
0667 
0668 tm = choice('Threshold masking','tm',...
0669     {tm_none,tma,tmr},'');
0670 p1=['Images are thresholded at a given value and ',...
0671     'only voxels at which all images exceed the threshold are included. '];
0672 tm.help={p1,sp_text};
0673 
0674 masking = branch('Masking','masking',...
0675     {tm,im,em},'');
0676 p1=['The mask specifies the voxels within the image volume which are to be ',...
0677     'assessed. SPM supports three methods of masking (1) Threshold, ',...
0678     '(2) Implicit and (3) Explicit. The volume analysed ',...
0679     'is the intersection of all masks.'];
0680 masking.help={p1,sp_text};
0681 
0682 %-------------------------------------------------------------------------
0683 % Global calculation
0684 
0685 global_uval = entry('Global values','global_uval','e',[Inf 1],'');
0686 p1=['Enter the vector of global values'];
0687 global_uval.val={0};
0688 global_uval.help={p1,sp_text};
0689 
0690 g_user = branch('User','g_user',{global_uval},'');
0691 g_user.help={'User defined  global effects (enter your own ',...
0692         'vector of global values)'};
0693 
0694 g_mean.type = 'const';
0695 g_mean.name = 'Mean';
0696 g_mean.tag = 'g_mean';
0697 g_mean.val = {[]};
0698 p1=['SPM standard mean voxel value'];
0699 p2=['This defines the global mean via a two-step process. Firstly, the overall ',...
0700         'mean is computed. Voxels with values less than 1/8 of this value are then ',...
0701         'deemed extra-cranial and get masked out. The mean is then recomputed on the ',...
0702         'remaining voxels.'];
0703 g_mean.help = {p1,sp_text,p2,sp_text};
0704 
0705 g_omit.type = 'const';
0706 g_omit.name = 'Omit';
0707 g_omit.tag = 'g_omit';
0708 g_omit.val = {[]};
0709 g_omit.help = {'Omit'};
0710 
0711 globalc = choice('Global calculation','globalc',...
0712     {g_omit,g_user,g_mean},'');
0713 p1=['There are three methods for estimating global effects ',...
0714     '(1) Omit (assumming no other options requiring the global value chosen) ',...
0715     '(2) User defined (enter your own vector of global values) ',...
0716     '(3) Mean: SPM standard mean voxel value (within per image fullmean/8 mask) '];
0717 globalc.help={p0,sp_text,p1,sp_text};
0718 
0719 %-------------------------------------------------------------------------
0720 % Global options
0721 
0722 gmsca_no.type = 'const';
0723 gmsca_no.name = 'No';
0724 gmsca_no.tag  = 'gmsca_no';
0725 gmsca_no.val  = {[]};
0726 gmsca_no.help = {'No overall grand mean scaling'};
0727 
0728 gmscv      = entry('Grand mean scaled value','gmscv','e',[Inf 1],'');
0729 gmscv.val={50};
0730 p1=['The default value of 50, scales the global flow to a physiologically ',...
0731         'realistic value of 50ml/dl/min.'];
0732 gmscv.help={p1,sp_text};
0733 
0734 gmsca_yes=branch('Yes','gmsca_yes',{gmscv},'');
0735 p1 =['Scaling of the overall grand mean simply ',...
0736     'scales all the data by a common factor such that the mean of all the ',...
0737     'global values is the value specified. For qualitative data, this puts ',...
0738     'the data into an intuitively accessible scale without altering the ',...
0739     'statistics. '];
0740 gmsca_yes.help={p1,sp_text};
0741 
0742 gmsca = choice('Overall grand mean scaling','gmsca',...
0743     {gmsca_no,gmsca_yes},'');
0744 
0745 p2=['When proportional scaling global normalisation is used ',...
0746     'each image is separately scaled such that it''s global ',...
0747     'value is that specified (in which case the grand mean is also ',...
0748     'implicitly scaled to that value). When using AnCova or no global ',...
0749     'normalisation, with data from different subjects or sessions, an ',...
0750     'intermediate situation may be appropriate, and you may be given the ',...
0751     'option to scale group, session or subject grand means separately. '];
0752 gmsca.help={p1,sp_text,p2,sp_text};
0753 
0754 
0755 
0756 glonorm = mnu('Normalisation','glonorm',...
0757     {'None','Proportional','ANCOVA'},{1,2,3},'');
0758 glonorm.val={1};
0759 p1 = ['Global nuisance effects are usually ',...
0760     'accounted for either by scaling the images so that they all have the ',...
0761     'same global value (proportional scaling), or by including the global ',...
0762     'covariate as a nuisance effect in the general linear model (AnCova). ',...
0763     'Much has been written on which to use, and when. Basically, since ',...
0764     'proportional scaling also scales the variance term, it is appropriate ',...
0765     'for situations where the global measurement predominantly reflects ',...
0766     'gain or sensitivity. Where variance is constant across the range of ',...
0767     'global values, linear modelling in an AnCova approach has more ',...
0768     'flexibility, since the model is not restricted to a simple ',...
0769     'proportional regression. '];
0770 
0771 p2=['''Ancova by subject'' or ''Ancova by effect'' options are implemented ',...
0772     'using the ANCOVA options provided where each experimental factor ',...
0773     '(eg. subject or effect), is defined. These allow eg. different subjects ',...
0774     'to have different relationships between local and global measurements. '];
0775 
0776 p3 =['Since differences between subjects may be due to gain and sensitivity ',...
0777     'effects, AnCova by subject could be combined with "grand mean scaling ',...
0778     'by subject" (an option also provided where each experimental factor is ',...
0779     'originally defined) to obtain a combination of between subject proportional ',...
0780     'scaling and within subject AnCova. '];
0781 glonorm.help={p1,sp_text,p2,sp_text,p3,sp_text};
0782 
0783 globalm = branch('Global normalisation','globalm',...
0784     {gmsca,glonorm},'');
0785 globalm.help={p0,sp_text,p1,sp_text,p2,sp_text,p3,sp_text};
0786 
0787 %-------------------------------------------------------------------------
0788 % Directory
0789 
0790 cdir = files('Directory','dir','dir',1,'');
0791 cdir.help = {[...
0792 'Select a directory where the SPM.mat file containing the ',...
0793 'specified design matrix will be written.']};
0794 
0795 %-------------------------------------------------------------------------
0796 % Main routine
0797 
0798 conf = branch('Factorial design specification','factorial_design',...
0799     {des,covs,masking,globalc,globalm,cdir},'');
0800 p1=['This interface is used for setting up analyses of PET data. It is also ',...
0801     'used for ''2nd level'' or ''random effects'' analysis which allow ',...
0802     'one to make a population inference. First level models can be used to produce ',...
0803     'appropriate summary data, which can then be used as raw data for a second-level ',...
0804     'analysis. For example, a simple t-test on contrast images from the first-level ',...
0805     'turns out to be a random-effects analysis with random subject effects, inferring ',...
0806     'for the population based on a particular sample of subjects.'];
0807 
0808 p2=['This interface configures the design matrix, describing the general ',...
0809     'linear model, data specification, and other parameters necessary for ',...
0810     'the statistical analysis. These parameters are saved in a ',...
0811     'configuration file (SPM.mat), which can then be ',...
0812     'passed on to spm_spm.m which estimates the design. This is achieved by ',...
0813     'pressing the ''Estimate'' button. Inference on these ',...
0814     'estimated parameters is then handled by the SPM results section. '];
0815 
0816 p3=['A separate interface handles design configuration ',...
0817     'for fMRI time series.'];
0818 
0819 p4=['Various data and parameters need to be supplied to specify the design ',...
0820     '(1) the image files, (2) indicators of the corresponding condition/subject/group ',...
0821     '(2) any covariates, nuisance variables, or design matrix partitions ',...
0822     '(3) the type of global normalisation (if any) ',...
0823     '(4) grand mean scaling options ',...
0824     '(5) thresholds and masks defining the image volume to analyse. ',...
0825     'The interface supports a comprehensive range of options for all these parameters.'];
0826 
0827 conf.help={p1,sp_text,p2,sp_text,p3,sp_text,p4,sp_text};
0828 conf.prog   = @run_stats;
0829 
0830 return;
0831 %-------------------------------------------------------------------------
0832 
0833 function run_stats(job)
0834 
0835 spm_defaults;
0836 
0837 original_dir = pwd;
0838 cd(job.dir{1});
0839     
0840 %-Ask about overwriting files from previous analyses...
0841 %-------------------------------------------------------------------
0842 if exist(fullfile(job.dir{1},'SPM.mat'),'file')
0843     str = {    'Current directory contains existing SPM file:',...
0844             'Continuing will overwrite existing file!'};
0845     if spm_input(str,1,'bd','stop|continue',[1,0],1,mfilename);
0846         fprintf('%-40s: %30s\n\n',...
0847             'Abort...   (existing SPM file)',spm('time'));
0848         return
0849     end
0850 end
0851 
0852 % If we've gotten to this point we're committed to overwriting files.
0853 % Delete them so we don't get stuck in spm_spm
0854 %------------------------------------------------------------------------
0855 files = {'^mask\..{3}$','^ResMS\..{3}$','^RPV\..{3}$',...
0856          '^beta_.{4}\..{3}$','^con_.{4}\..{3}$','^ResI_.{4}\..{3}$',...
0857          '^ess_.{4}\..{3}$', '^spm\w{1}_.{4}\..{3}$'};
0858 
0859 for i=1:length(files)
0860     j = spm_select('List',pwd,files{i});
0861     for k=1:size(j,1)
0862         spm_unlink(deblank(j(k,:)));
0863     end
0864 end
0865 
0866 
0867 %-Option definitions
0868 %-------------------------------------------------------------------
0869 %-Generic factor names
0870 sF = {'sF1','sF2','sF3','sF4'};
0871 
0872 %-Covariate by factor interaction options
0873 sCFI = {'<none>';...                            %-1
0874         'with sF1';'with sF2';'with sF3';'with sF4';...            %-2:5
0875         'with sF2 (within sF4)';'with sF3 (within sF4)'};        %-6,7
0876 
0877 %-DesMtx argument components for covariate by factor interaction options
0878 % (Used for CFI's Covariate Centering (CC), GMscale & Global normalisation)
0879 CFIforms = {    '[]',        'C',    '{}';...            %-1
0880         'I(:,1)',        'FxC',    '{sF{1}}';...            %-2
0881         'I(:,2)',        'FxC',    '{sF{2}}';...            %-3
0882         'I(:,3)',        'FxC',    '{sF{3}}';...            %-4
0883         'I(:,4)',        'FxC',    '{sF{4}}';...            %-5
0884         'I(:,[4,2])',    'FxC',    '{sF{4},sF{2}}';...    %-6
0885         'I(:,[4,3])',    'FxC',    '{sF{4},sF{3}}'    };    %-7
0886 
0887 %-Centre (mean correction) options for covariates & globals            (CC)
0888 % (options 9-12 are for centering of global when using AnCova GloNorm) (GC)
0889 sCC = {        'around overall mean';...                %-1
0890         'around sF1 means';...                    %-2
0891         'around sF2 means';...                    %-3
0892         'around sF3 means';...                    %-4
0893         'around sF4 means';...                    %-5
0894         'around sF2 (within sF4) means';...            %-6
0895         'around sF3 (within sF4) means';...            %-7
0896         '<no centering>';...                    %-8
0897         'around user specified value';...            %-9
0898         '(as implied by AnCova)';...                %-10
0899         'GM';...                        %-11
0900         '(redundant: not doing AnCova)'}';            %-12
0901 %-DesMtx I forms for covariate centering options
0902 CCforms = {'ones(nScan,1)',CFIforms{2:end,1},''}';
0903 
0904 %-Global calculation options                                       (GXcalc)
0905 sGXcalc  = {    'omit';...                        %-1
0906         'user specified';...                    %-2
0907         'mean voxel value (within per image fullmean/8 mask)'};    %-3
0908 
0909 
0910 %-Global normalization options  (GloNorm)
0911 sGloNorm = {    'AnCova';...                        %-1
0912         'AnCova by sF1';...                    %-2
0913         'AnCova by sF2';...                    %-3
0914         'AnCova by sF3';...                    %-4
0915         'AnCova by sF4';...                    %-5
0916         'AnCova by sF2 (within sF4)';...            %-6
0917         'AnCova by sF3 (within sF4)';...            %-7
0918         'proportional scaling';...                %-8
0919         '<no global normalisation>'};                %-9
0920 
0921 
0922 %-Grand mean scaling options                                        (GMsca)
0923 sGMsca = {    'scaling of overall grand mean';...            %-1
0924         'scaling of sF1 grand means';...            %-2
0925         'scaling of sF2 grand means';...            %-3
0926         'scaling of sF3 grand means';...            %-4
0927         'scaling of sF4 grand means';...            %-5
0928         'scaling of sF2 (within sF4) grand means';...        %-6
0929         'scaling of sF3 (within sF4) grand means';...        %-7
0930         '(implicit in PropSca global normalisation)';...    %-8
0931         '<no grand Mean scaling>'    };            %-9
0932 %-NB: Grand mean scaling by subject is redundent for proportional scaling
0933 
0934 % Nonsphericity defaults
0935 xVi.var=zeros(1,4);
0936 xVi.dep=zeros(1,4);
0937 
0938 % Conditions of no interest defaults
0939 B=[];
0940 Bnames={};
0941 
0942 switch strvcat(fieldnames(job.des)),
0943 case 't1',
0944     % One sample t-test
0945     DesName='One sample t-test';
0946     
0947     P=job.des.t1.scans;
0948     n=length(P);
0949     I=[1:n]';
0950     I=[I,ones(n,3)];
0951     
0952     [H,Hnames]=spm_DesMtx(I(:,2),'-','mean');
0953     
0954     SPM.factor(1).name='Group';
0955     SPM.factor(1).levels=1;
0956 
0957 case 't2', 
0958     % Two-sample t-test
0959     DesName='Two-sample t-test';
0960     
0961     P=job.des.t2.scans1;
0962     n1=length(job.des.t2.scans1);
0963     P=[P;job.des.t2.scans2];
0964     n2=length(job.des.t2.scans2);
0965     
0966     I=[];
0967     I=[1:n1]';
0968     I=[I;[1:n2]'];
0969     I=[I,[ones(n1,1);2*ones(n2,1)]];
0970     I=[I,ones(n1+n2,2)];
0971     
0972     [H,Hnames]=spm_DesMtx(I(:,2),'-','Group');
0973     
0974     % Nonsphericity options
0975     xVi.var(2)=job.des.t2.variance;
0976     xVi.dep(2)=job.des.t2.dept;  
0977     
0978     % Names and levels
0979     SPM.factor(1).name='Group';
0980     SPM.factor(1).levels=2;
0981     
0982     % Ancova options
0983     SPM.factor(1).gmsca=job.des.t2.gmsca;
0984     SPM.factor(1).ancova=job.des.t2.ancova;
0985 
0986 case 'pt', 
0987     % Paired t-test
0988     DesName='Paired t-test';
0989     
0990     Npairs=length(job.des.pt.pair);
0991     P=[];
0992     for p=1:Npairs,
0993          P=[P;job.des.pt.pair(p).scans];   
0994     end
0995     
0996     I=ones(Npairs*2,1);
0997     I(:,2)=kron(ones(Npairs,1),[1 2]');
0998     I(:,3)=kron([1:Npairs]',ones(2,1));
0999     I(:,4)=I(:,1);
1000     
1001     [H,Hnames]=spm_DesMtx(I(:,2),'-','Condition');
1002     [B,Bnames]=spm_DesMtx(I(:,3),'-','Subject');
1003     
1004     % Nonsphericity options
1005     xVi.var(2)=job.des.pt.variance;
1006     xVi.dep(2)=job.des.pt.dept;  
1007     
1008     % Names and levels
1009     SPM.factor(1).name='Group';
1010     SPM.factor(1).levels=2;
1011     
1012     % Ancova options
1013     SPM.factor(1).gmsca=job.des.pt.gmsca;
1014     SPM.factor(1).ancova=job.des.pt.ancova;
1015 
1016 case 'mreg', 
1017     % Multiple regression
1018     DesName='Multiple regression';
1019     
1020     P=job.des.mreg.scans;
1021     n=length(P);
1022     I=[1:n]';
1023     I=[I,ones(n,3)];
1024     
1025     % Names and levels
1026     SPM.factor(1).name='';
1027     SPM.factor(1).levels=0;
1028     
1029     H=[];Hnames=[];
1030     [B,Bnames]=spm_DesMtx(I(:,2),'-','mean');
1031     
1032     for i=1:length(job.des.mreg.mcov)
1033         job.cov(i).c=job.des.mreg.mcov(i).c;
1034         job.cov(i).cname=job.des.mreg.mcov(i).cname;
1035         job.cov(i).iCFI=1;
1036         job.cov(i).iCC=1;
1037     end
1038 
1039 case 'fd', 
1040     % Full Factorial Design
1041     DesName='Full factorial';
1042     
1043     [I,P,H,Hnames] = spm_set_factorial_design (job);
1044     
1045     Nfactors=length(job.des.fd.fact);
1046     for i=1:Nfactors,
1047         % Nonsphericity
1048         xVi.var(i+1)=job.des.fd.fact(i).variance;
1049         xVi.dep(i+1)=job.des.fd.fact(i).dept; 
1050         
1051         % Store names and levels
1052         SPM.factor(i).name=job.des.fd.fact(i).name; 
1053         SPM.factor(i).levels=job.des.fd.fact(i).levels; 
1054         
1055         % Ancova options
1056         SPM.factor(i).gmsca=job.des.fd.fact(i).gmsca;
1057         SPM.factor(i).ancova=job.des.fd.fact(i).ancova;
1058     end
1059 
1060 case 'fblock',
1061     % Flexible factorial design
1062     DesName='Flexible factorial';
1063     
1064     if isfield(job.des.fblock.fsuball,'fsubject')
1065         nsub=length(job.des.fblock.fsuball.fsubject);
1066         % Specify design subject-by-subject
1067         P=[];I=[];
1068         subj=[];
1069         for s=1:nsub,
1070             P  = [P; job.des.fblock.fsuball.fsubject(s).scans];
1071             ns = length(job.des.fblock.fsuball.fsubject(s).scans);
1072             cc = job.des.fblock.fsuball.fsubject(s).conds;
1073 
1074             [ccr,ccc] = size(cc);
1075             if ~(ccr==ns) & ~(ccc==ns)
1076                 disp(sprintf('Error for subject %d: conditions not specified for each scan',s));
1077                 return
1078             elseif ccc==ns
1079                 cc=cc';
1080             end
1081             subj=[subj;s*ones(ns,1)];
1082             I = [I; [[1:ns]',cc]];
1083         end
1084         
1085         nf=length(job.des.fblock.fac);
1086         subject_factor=0;
1087         for i=1:nf,
1088             if strcmp(job.des.fblock.fac(i).name,'Repl') | ...
1089                     strcmp(job.des.fblock.fac(i).name,'repl')
1090                     % Copy `replications' column to create explicit `replications' factor
1091                     nI=I(:,1:i);
1092                     nI=[nI,I(:,1)];
1093                     nI=[nI,I(:,i+1:end)];
1094                     I=nI;
1095             end
1096             if strcmp(job.des.fblock.fac(i).name,'Subject') | ...
1097                     strcmp(job.des.fblock.fac(i).name,'subject')
1098                     % Create explicit `subject' factor
1099                     nI=I(:,1:i);
1100                     nI=[nI,subj];
1101                     nI=[nI,I(:,i+1:end)];
1102                     I=nI;
1103                     subject_factor=1;
1104             end
1105             
1106         end
1107         
1108         % Re-order scans and conditions into standard format
1109         % This is to ensure compatibility with how variance components are created
1110         if subject_factor
1111              U=unique(I(:,2:nf+1),'rows');
1112              Un=length(U);
1113              Uc=zeros(Un,1);
1114              r=1;rj=[];
1115              for k=1:Un,
1116                  for j=1:size(I,1),
1117                      match=sum(I(j,2:nf+1)==U(k,:))==nf;
1118                      if match
1119                          Uc(k)=Uc(k)+1;
1120                          Ir(r,:)=[Uc(k),I(j,2:end)];
1121                          r=r+1;
1122                          rj=[rj;j];
1123                      end
1124                  end
1125              end
1126         end
1127         P=P(rj); % -scans
1128         I=Ir;    % -conditions
1129         
1130     else
1131         [ns,nf]=size(job.des.fblock.fsuball.specall.imatrix);
1132         if nf > 4
1133             disp('Error in factorial matrix: number of factors should be less than 5');
1134             return
1135         end
1136         I=job.des.fblock.fsuball.specall.imatrix;
1137         P=job.des.fblock.fsuball.specall.scans;
1138     end
1139     
1140     % Pad out factorial matrix to cover the four canonical factors
1141     [ns,nI]=size(I);
1142     if nI < 4
1143         I = [I, ones(ns,4-nI)];
1144     end
1145 
1146     % Sort main effects and interactions
1147     fmain = struct('fnum',{});
1148     inter = struct('fnums',{});
1149     for k=1:numel(job.des.fblock.maininters)
1150         if isfield(job.des.fblock.maininters{k},'fmain')
1151             fmain(end+1)=job.des.fblock.maininters{k}.fmain;
1152         elseif isfield(job.des.fblock.maininters{k},'inter')
1153             inter(end+1)=job.des.fblock.maininters{k}.inter;
1154         end;
1155     end;
1156 
1157     % Create main effects
1158     H=[];Hnames=[];
1159     nmain=length(fmain);
1160     for f=1:nmain,
1161         fcol=fmain(f).fnum;
1162         fname=job.des.fblock.fac(fcol).name;
1163         
1164         % Augment H partition - explicit factor numbers are 1 lower than in I matrix
1165         [Hf,Hfnames]=spm_DesMtx(I(:,fcol+1),'-',fname);
1166         H=[H,Hf];
1167         Hnames=[Hnames;Hfnames];
1168     end
1169     
1170     % Create interactions
1171     ni=length(inter);
1172     for i=1:ni,
1173         % Get the two factors for this interaction
1174         fnums=inter(i).fnums;
1175         f1=fnums(1);f2=fnums(2);
1176         
1177         % Names
1178         iname{1}=job.des.fblock.fac(f1).name;
1179         iname{2}=job.des.fblock.fac(f2).name;
1180         
1181         % Augment H partition - explicit factor numbers are 1 lower than in I matrix
1182         Isub=[I(:,f1+1),I(:,f2+1)];
1183         [Hf,Hfnames]=spm_DesMtx(Isub,'-',iname);
1184         H=[H,Hf];
1185         Hnames=[Hnames;Hfnames];
1186         
1187     end
1188     
1189     if nmain==0 & ni==0
1190         disp('Error in design specification: You have not specified any main effects or interactions');
1191         return
1192     end
1193     
1194     for i=1:nf,
1195         % Set nonsphericity options
1196         xVi.var(i+1)=job.des.fblock.fac(i).variance;
1197         xVi.dep(i+1)=job.des.fblock.fac(i).dept;
1198         
1199         % Store names and levels
1200         SPM.factor(i).name=job.des.fblock.fac(i).name;
1201         SPM.factor(i).levels=length(unique(I(:,i+1)));
1202         
1203         % Ancova options
1204         SPM.factor(i).gmsca=job.des.fblock.fac(i).gmsca;
1205         SPM.factor(i).ancova=job.des.fblock.fac(i).ancova;
1206     end
1207     
1208     
1209 end
1210 nScan=size(I,1); %-#obs
1211 
1212 % Set up data strucutres for `EEG' non-sphericity routine
1213 nf=length(SPM.factor);nef=nf+1;
1214 If=[I(:,2:nef),I(:,1)];
1215 for i=1:length(SPM.factor), % Other factors
1216     SPM.eeg.Nlevels{i}=SPM.factor(i).levels;    
1217     SPM.eeg.Xind{i}=unique(If(:,1:i),'rows');
1218     SPM.xVi.Qidentical{i}=1-xVi.var(i+1);
1219     SPM.xVi.Qindependent{i}=1-xVi.dep(i+1);
1220 end
1221 SPM.eeg.Nlevels{nef}=length(unique(I(:,1))); % Repetition factor
1222 SPM.eeg.Xind{nef}=unique(If(:,1:nef),'rows');
1223 SPM.xVi.Qidentical{nef}=1;
1224 SPM.xVi.Qindependent{nef}=1;
1225 SPM.eeg.Nfactors=nef;
1226 SPM = spm_eeg_get_vc(SPM);  % Call `EEG' non-sphericity routine
1227 SPM.xVi.I=I;
1228 % Find covariance elements with missing repetitions
1229 ncells=prod([SPM.eeg.Nlevels{1:nf}]);
1230 fullreps=repmat([1:SPM.eeg.Nlevels{nef}]',[ncells,1]);
1231 missing=[];j=1;
1232 for i=1:length(fullreps),
1233     if ~(fullreps(i)==I(j))
1234         missing=[missing i];
1235     else
1236         j=j+1;
1237     end
1238 end
1239 % Remove covariance elements with missing repetitions
1240 for i=1:length(SPM.xVi.Vi)
1241     SPM.xVi.Vi{i}(missing,:)=[];
1242     SPM.xVi.Vi{i}(:,missing)=[];
1243 end
1244 
1245 %-Covariate partition(s): interest (C) & nuisance (G) excluding global
1246 %===================================================================
1247 dstr   = {'covariate','nuisance variable'};
1248 C  = []; Cnames = [];  %-Covariate DesMtx partitions & names
1249 G  = []; Gnames = []; 
1250 
1251 xC = [];                         %-Struct array to hold raw covariates
1252 
1253 % Covariate options:
1254 nc=length(job.cov); % number of covariates
1255 for i=1:nc,
1256     
1257     c      = job.cov(i).c;
1258     cname  = job.cov(i).cname;
1259     rc     = c;                         %-Save covariate value
1260     rcname = cname;                     %-Save covariate name
1261     if job.cov(i).iCFI==1,
1262         iCFI=1;
1263     else
1264         % SPMs internal factor numbers are 1 higher than specified in user
1265         % interface as, internally, the first factor is always `replication'
1266         iCFI=job.cov(i).iCFI+1;
1267     end
1268     switch job.cov(i).iCC,
1269         case 1
1270             iCC=1;
1271         case {2,3,4}
1272             iCC=job.cov(i).iCC+1;
1273         otherwise
1274             iCC=job.cov(i).iCC+3;
1275     end
1276     
1277     %-Centre within factor levels as appropriate
1278     if any(iCC == [1:7]), 
1279         c = c - spm_meanby(c,eval(CCforms{iCC})); 
1280     end
1281     
1282     %-Do any interaction (only for single covariate vectors)
1283     %-----------------------------------------------------------
1284     if iCFI > 1                %-(NB:iCFI=1 if size(c,2)>1)
1285         tI        = [eval(CFIforms{iCFI,1}),c];
1286         tConst    = CFIforms{iCFI,2};
1287         tFnames   = [eval(CFIforms{iCFI,3}),{cname}];
1288         [c,cname] = spm_DesMtx(tI,tConst,tFnames);
1289     elseif size(c,2)>1            %-Design matrix block
1290         [null,cname] = spm_DesMtx(c,'X',cname);
1291     else
1292         cname = {cname};
1293     end
1294     
1295     %-Store raw covariate details in xC struct for reference
1296     %-Pack c into appropriate DesMtx partition
1297     %-----------------------------------------------------------
1298     %-Construct description string for covariate
1299     str = {sprintf('%s',rcname)};
1300     if size(rc,2)>1, str = {sprintf('%s (block of %d covariates)',...
1301                 str{:},size(rc,2))}; end
1302     if iCC < 8, str=[str;{['used centered ',sCC{iCC}]}]; end
1303     if iCFI> 1, str=[str;{['fitted as interaction ',sCFI{iCFI}]}]; end
1304     
1305     tmp       = struct(    'rc',rc,    'rcname',rcname,...
1306         'c',c,        'cname',{cname},...
1307         'iCC',iCC,    'iCFI',iCFI,...
1308         'type',i,...
1309         'cols',[1:size(c,2)] + ...
1310         size([H,C],2) +  ...
1311         size([B,G],2)*(i-1),...
1312         'descrip',{str}                );
1313     if isempty(xC), xC = tmp; else, xC = [xC,tmp]; end
1314     C     = [C,c];
1315     Cnames = [Cnames; cname];
1316     
1317 end    
1318 clear c tI tConst tFnames
1319     
1320 xGX=[];
1321 xM=[];
1322 
1323 
1324 
1325 %===================================================================
1326 % - C O N F I G U R E   D E S I G N -
1327 %===================================================================
1328 
1329 %-Images & image info: Map Y image files and check consistency of
1330 % dimensions and orientation / voxel size
1331 %===================================================================
1332 fprintf('%-40s: ','Mapping files')                               %-#
1333 VY    = spm_vol(char(P));
1334 
1335 %-Check compatability of images (Bombs for single image)
1336 %-------------------------------------------------------------------
1337 spm_check_orientations(VY);
1338 
1339 fprintf('%30s\n','...done')                                      %-#
1340 
1341 %-Global values, scaling and global normalisation
1342 %===================================================================
1343 %-Compute global values
1344 %-------------------------------------------------------------------
1345 switch strvcat(fieldnames(job.globalc))
1346     case 'g_omit',
1347         iGXcalc=1;
1348     case 'g_user',
1349         iGXcalc=2;
1350     case 'g_mean',
1351         iGXcalc=3;
1352 end
1353 
1354 switch job.globalm.glonorm
1355     case 1,
1356         iGloNorm=9;
1357     case 2,
1358         iGloNorm=8;
1359     case 3,
1360         iGloNorm=1;
1361 end
1362 if SPM.factor(1).levels > 1
1363     % Over-ride if factor-specific ANCOVA has been specified
1364     for i=1:length(SPM.factor),
1365         if SPM.factor(i).ancova
1366             iGloNorm=i+2;
1367         end
1368     end
1369 end
1370 
1371 %-Analysis threshold mask
1372 %-------------------------------------------------------------------
1373 %-Work out available options:
1374 % -Inf=>None, real=>absolute, complex=>proportional, (i.e. times global)
1375 M_T = -Inf;
1376 switch strvcat(fieldnames(job.masking.tm)),
1377     case 'tma',
1378         % Absolute
1379         M_T = job.masking.tm.tma.athresh;
1380     case 'tmr',
1381         % Relative
1382         M_T = job.masking.tm.tmr.rthresh*sqrt(-1);
1383         % Need to force calculation of globals
1384         iGXcalc=3;
1385     case 'tm_none'
1386         % None
1387         M_T = -Inf;
1388 end
1389 
1390 if (any(iGloNorm == [1:5]) | iGloNorm==8) & iGXcalc==1
1391     % Over-ride omission of global calculation if we need it
1392     disp(' ');
1393     disp(sprintf('For %s, SPM needs estimates of global activity.',sGloNorm{iGloNorm}));
1394     disp('But you have specified to omit this computation.');
1395     disp('SPM has overridden this omission and will automatically compute ');
1396     disp('globals as the mean value of within brain voxels');
1397     disp(' ');
1398     iGXcalc=3;
1399 end
1400 sGXcalc = sGXcalc{iGXcalc};
1401 
1402 switch iGXcalc, 
1403     case 1
1404         %-Don't compute => no GMsca (iGMsca==9) or GloNorm (iGloNorm==9)
1405         g = [];
1406     case 2
1407         %-User specified globals
1408         g = job.globalc.g_user.global_uval;
1409     case 3 
1410         %-Compute as mean voxel value (within per image fullmean/8 mask)
1411         g = zeros(nScan,1 );
1412         fprintf('%-40s: %30s','Calculating globals',' ')             %-#
1413         for i = 1:nScan
1414             str = sprintf('%3d/%-3d',i,nScan);
1415             fprintf('%s%30s',repmat(sprintf('\b'),1,30),str)%-#
1416             g(i) = spm_global(VY(i));
1417         end
1418         fprintf('%s%30s\n',repmat(sprintf('\b'),1,30),'...done')       %-#
1419     otherwise
1420         error('illegal iGXcalc')
1421 end
1422 rg = g;
1423     
1424 fprintf('%-40s: ','Design configuration')                        %-#
1425 
1426 if iGloNorm==8
1427     iGMsca=8;    %-grand mean scaling implicit in PropSca GloNorm
1428 else
1429     switch strvcat(fieldnames(job.globalm.gmsca))
1430         case 'gmsca_yes',
1431             iGMsca=1;
1432         case 'gmsca_no',
1433             iGMsca=9;
1434     end
1435     if SPM.factor(1).levels > 1
1436         % Over-ride if factor-specific scaling has been specified
1437         for i=1:length(SPM.factor),
1438             if SPM.factor(i).gmsca
1439                 iGMsca=i+2;
1440             end
1441         end
1442     end
1443 end
1444 
1445 %-Value for PropSca / GMsca                                     (GM)
1446 %-------------------------------------------------------------------
1447 switch iGMsca,
1448     case 9                      %-Not scaling (GMsca or PropSca)
1449         GM = 0;                         %-Set GM to zero when not scaling
1450     case 1                                %-Ask user value of GM
1451         GM = job.globalm.gmsca.gmsca_yes.gmscv;
1452         %-If GM is zero then don't GMsca! or PropSca GloNorm
1453         if GM==0, 
1454             iGMsca=9; 
1455             if iGloNorm==8, 
1456                 iGloNorm=9; 
1457             end
1458         end
1459     otherwise
1460         % Grand mean scaling by factor eg. scans are scaled so that the
1461         % mean global value over each level of the factor is set to GM
1462         GM=50;
1463 end
1464 
1465 %-Sort out description strings for GloNorm and GMsca
1466 %-------------------------------------------------------------------
1467 sGloNorm = sGloNorm{iGloNorm};
1468 sGMsca   = sGMsca{iGMsca};
1469 if iGloNorm==8
1470     sGloNorm = sprintf('%s to %-4g',sGloNorm,GM);
1471 elseif iGMsca<8
1472     sGMsca   = sprintf('%s to %-4g',sGMsca,GM);
1473 end
1474 
1475 %-Scaling: compute global scaling factors gSF required to implement
1476 % proportional scaling global normalisation (PropSca) or grand mean
1477 % scaling (GMsca), as specified by iGMsca (& iGloNorm)
1478 %-------------------------------------------------------------------
1479 switch iGMsca, 
1480     case 8
1481         %-Proportional scaling global normalisation
1482         if iGloNorm~=8, error('iGloNorm-iGMsca(8) mismatch for PropSca'), end
1483         gSF    = GM./g;
1484         g      = GM*ones(nScan,1);
1485     case {1,2,3,4,5,6,7}
1486         %-Grand mean scaling according to iGMsca
1487         gSF    = GM./spm_meanby(g,eval(CCforms{iGMsca}));
1488         g      = g.*gSF;
1489     case 9
1490         %-No grand mean scaling
1491         gSF    = ones(nScan,1);
1492     otherwise
1493         error('illegal iGMsca')
1494     end
1495 
1496 
1497 %-Apply gSF to memory-mapped scalefactors to implement scaling
1498 %-------------------------------------------------------------------
1499 for i = 1:nScan
1500     VY(i).pinfo(1:2,:) = VY(i).pinfo(1:2,:)*gSF(i);
1501 end
1502     
1503 %-Global centering (for AnCova GloNorm)                         (GC)
1504 %-If not doing AnCova then GC is irrelevant
1505 if ~any(iGloNorm == [1:7])
1506     iGC = 12;
1507     gc  = [];
1508 else
1509     iGC = 10;
1510     gc = 0;
1511 end
1512         
1513 %-AnCova: Construct global nuisance covariates partition (if AnCova)
1514 %-------------------------------------------------------------------
1515 if any(iGloNorm == [1:7])
1516     
1517     %-Centre global covariate as requested
1518     %---------------------------------------------------------------
1519     switch iGC, case {1,2,3,4,5,6,7}    %-Standard sCC options
1520         gc = spm_meanby(g,eval(CCforms{iGC}));
1521     case 8                    %-No centering
1522         gc = 0;
1523     case 9                    %-User specified centre
1524         %-gc set above
1525     case 10                    %-As implied by AnCova option
1526         gc = spm_meanby(g,eval(CCforms{iGloNorm}));
1527     case 11                    %-Around GM
1528         gc = GM;
1529     otherwise                %-unknown iGC
1530         error('unexpected iGC value')
1531 end
1532 
1533 %-AnCova - add scaled centred global to DesMtx `G' partition
1534 %---------------------------------------------------------------
1535 rcname     = 'global'; 
1536 tI         = [eval(CFIforms{iGloNorm,1}),g - gc];
1537 tConst     = CFIforms{iGloNorm,2};
1538 tFnames    = [eval(CFIforms{iGloNorm,3}),{rcname}];
1539 [f,gnames]  = spm_DesMtx(tI,tConst,tFnames);
1540 clear tI tConst tFnames
1541 
1542 %-Save GX info in xC struct for reference
1543 %---------------------------------------------------------------
1544 str     = {sprintf('%s: %s',dstr{2},rcname)};
1545 if any(iGMsca==[1:7]), str=[str;{['(after ',sGMsca,')']}]; end
1546 if iGC ~= 8, str=[str;{['used centered ',sCC{iGC}]}]; end
1547 if iGloNorm > 1
1548     str=[str;{['fitted as interaction ',sCFI{iGloNorm}]}]; 
1549 end
1550 tmp  = struct(    'rc',rg.*gSF,        'rcname',rcname,...
1551     'c',f,            'cname'    ,{gnames},...
1552     'iCC',iGC,        'iCFI'    ,iGloNorm,...
1553     'type',            3,...
1554     'cols',[1:size(f,2)] + size([H C B G],2),...
1555     'descrip',        {str}        );
1556 
1557 G = [G,f]; Gnames = [Gnames; gnames];
1558 if isempty(xC), xC = tmp; else, xC = [xC,tmp]; end
1559 
1560 
1561 elseif iGloNorm==8 | iGXcalc>1
1562     
1563     %-Globals calculated, but not AnCova: Make a note of globals
1564     %---------------------------------------------------------------
1565     if iGloNorm==8
1566         str = { 'global values: (used for proportional scaling)';...
1567                 '("raw" unscaled globals shown)'};
1568     elseif isfinite(M_T) & ~isreal(M_T)
1569         str = { 'global values: (used to compute analysis threshold)'};
1570     else
1571         str = { 'global values: (computed but not used)'};
1572     end
1573     
1574     rcname ='global';
1575     tmp     = struct(    'rc',rg,    'rcname',rcname,...
1576         'c',{[]},    'cname'    ,{{}},...
1577         'iCC',0,    'iCFI'    ,0,...
1578         'type',        3,...
1579         'cols',        {[]},...
1580         'descrip',    {str}            );
1581     
1582     if isempty(xC), xC = tmp; else, xC = [xC,tmp]; end
1583 end
1584 
1585 
1586 %-Save info on global calculation in xGX structure
1587 %-------------------------------------------------------------------
1588 xGX = struct(...
1589     'iGXcalc',iGXcalc,    'sGXcalc',sGXcalc,    'rg',rg,...
1590     'iGMsca',iGMsca,    'sGMsca',sGMsca,    'GM',GM,'gSF',gSF,...
1591     'iGC',    iGC,        'sGC',    sCC{iGC},    'gc',    gc,...
1592     'iGloNorm',iGloNorm,    'sGloNorm',sGloNorm);
1593     
1594 %-Make a description string
1595 %-------------------------------------------------------------------
1596 if isinf(M_T)
1597     xsM.Analysis_threshold = 'None (-Inf)';
1598 elseif isreal(M_T)
1599     xsM.Analysis_threshold = sprintf('images thresholded at %6g',M_T);
1600 else
1601     xsM.Analysis_threshold = sprintf(['images thresholded at %6g ',...
1602             'times global'],imag(M_T));
1603 end
1604     
1605 %-Construct masking information structure and compute actual analysis
1606 % threshold using scaled globals (rg.*gSF)
1607 %-------------------------------------------------------------------
1608  if isreal(M_T),
1609      M_TH = M_T  * ones(nScan,1);    %-NB: -Inf is real
1610  else,        
1611      M_TH = imag(M_T) * (rg.*gSF); 
1612  end
1613     
1614 %-Implicit masking: Ignore zero voxels in low data-types?
1615 %-------------------------------------------------------------------
1616 % (Implicit mask is NaN in higher data-types.)
1617 type = getfield(spm_vol(P{1,1}),'dt')*[1,0]';
1618 if ~spm_type(type,'nanrep')
1619     M_I = job.masking.im;  % Implicit mask ?
1620     if M_I, 
1621         xsM.Implicit_masking = 'Yes: zero''s treated as missing';
1622     else,   
1623         xsm.Implicit_masking = 'No'; 
1624     end
1625 else
1626     M_I = 1;
1627     xsM.Implicit_masking = 'Yes: NaN''s treated as missing';
1628 end
1629     
1630 %-Explicit masking
1631 %-------------------------------------------------------------------
1632 if isempty(job.masking.em{:})
1633     VM = [];
1634     xsM.Explicit_masking = 'No'; 
1635 else
1636     VM = job.masking.em;
1637     xsM.Explicit_masking = 'Yes';
1638 end
1639 
1640 xM     = struct('T',M_T, 'TH',M_TH, 'I',M_I, 'VM',{VM}, 'xs',xsM);
1641 
1642 %-Construct full design matrix (X), parameter names and structure (xX)
1643 %===================================================================
1644 X      = [H C B G];
1645 tmp    = cumsum([size(H,2), size(C,2), size(B,2), size(G,2)]);
1646 xX     = struct(    'X',        X,...
1647     'iH',        [1:size(H,2)],...
1648     'iC',        [1:size(C,2)] + tmp(1),...
1649     'iB',        [1:size(B,2)] + tmp(2),...
1650     'iG',        [1:size(G,2)] + tmp(3),...
1651     'name',        {[Hnames; Cnames; Bnames; Gnames]},...
1652     'I',        I,...
1653     'sF',        {sF});
1654 
1655 
1656 %-Design description (an nx2 cellstr) - for saving and display
1657 %===================================================================
1658 tmp = {    sprintf('%d condition, +%d covariate, +%d block, +%d nuisance',...
1659         size(H,2),size(C,2),size(B,2),size(G,2));...
1660         sprintf('%d total, having %d degrees of freedom',...
1661         size(X,2),rank(X));...
1662         sprintf('leaving %d degrees of freedom from %d images',...
1663         size(X,1)-rank(X),size(X,1))                };
1664 xsDes = struct(    'Design',            {DesName},...
1665     'Global_calculation',        {sGXcalc},...
1666     'Grand_mean_scaling',        {sGMsca},...
1667     'Global_normalisation',        {sGloNorm},...
1668     'Parameters',            {tmp}            );
1669 
1670 
1671 fprintf('%30s\n','...done')                                      %-#
1672 
1673 %-Assemble SPM structure
1674 %===================================================================
1675 SPM.xY.P    = P;            % filenames
1676 SPM.xY.VY    = VY;            % mapped data
1677 SPM.nscan    = size(xX.X,1); % scan number
1678 SPM.xX        = xX;            % design structure
1679 SPM.xC        = xC;            % covariate structure
1680 SPM.xGX        = xGX;            % global structure
1681 SPM.xM        = xM;            % mask structure
1682 SPM.xsDes    = xsDes;        % description
1683 
1684 % Automatic contrast generation only works for 'Full factorials'
1685 if ~strcmp(DesName,'Full factorial')
1686     % Remove the .factor field to prevent attempted automatic contrast generation
1687     SPM=rmfield(SPM,'factor');
1688 end
1689 
1690 %-Save SPM.mat and set output argument
1691 %-------------------------------------------------------------------
1692 fprintf('%-40s: ','Saving SPM configuration')                    %-#
1693 
1694 if str2num(version('-release'))>=14,
1695     save('SPM', 'SPM', '-V6');
1696 else
1697     save('SPM', 'SPM');
1698 end;
1699 fprintf('%30s\n','...SPM.mat saved')                             %-#
1700 varargout = {SPM};
1701 
1702 %-Display Design report
1703 %===================================================================
1704 fprintf('%-40s: ','Design reporting')                            %-#
1705 fname     = cat(1,{SPM.xY.VY.fname}');
1706 spm_DesRep('DesMtx',SPM.xX,fname,SPM.xsDes)
1707 fprintf('%30s\n','...done')     
1708     
1709 cd(original_dir); % Change back dir
1710 fprintf('Done\n')