diff --git a/mv_classify.m b/mv_classify.m
index b320048..1552020 100644
--- a/mv_classify.m
+++ b/mv_classify.m
@@ -382,6 +382,11 @@
             sz_Xtrain = size(Xtrain);
             sz_Xtest = size(Xtest);
             
+            % JM ADDITION: NFEAT MIGHT HAVE CHANGED DURING PREPROCESSING
+            nfeat = [sz_Xtrain ones(1, numel(cfg.dimension_names) - ndims(sz_Xtrain))];
+            nfeat = nfeat(feature_dim);
+            if isempty(nfeat), nfeat = 1; end
+
             for ix = dim_loop                       % ---- search dimensions ----
                                 
                 % Training data for current search position
diff --git a/optimization/LedoitWolfEstimate.m b/optimization/LedoitWolfEstimate.m
index b90f9ae..91b48b9 100644
--- a/optimization/LedoitWolfEstimate.m
+++ b/optimization/LedoitWolfEstimate.m
@@ -23,7 +23,7 @@
 
 % remove mean
 if n > 1
-    X = X - repmat(mean(X,1), n, 1);
+    X = X - mean(X,1);
 else
     X = X - mean(X);
 end
diff --git a/preprocess/mv_preprocess.m b/preprocess/mv_preprocess.m
index f27b940..f4dc4b9 100644
--- a/preprocess/mv_preprocess.m
+++ b/preprocess/mv_preprocess.m
@@ -26,6 +26,22 @@
     % call preprocessing function
     [cfg.preprocess_param{pp}, X, clabel] = cfg.preprocess_fun{pp}(cfg.preprocess_param{pp}, X, clabel);
     
+    % the current preprocessing step may affect downstream preprocessing, e.g.
+    % when undersampling and subsequent confounding class specific mean subtraction
+    if cfg.preprocess_param{pp}.is_train_set && isfield(cfg.preprocess_param{pp}, 'keep_idx_train')
+      for p = (pp+1):numel(cfg.preprocess_fun)
+        for fn = cfg.preprocess_param{p}.select_data
+          cfg.preprocess_param{p}.([fn{1} '_train']) = cfg.preprocess_param{p}.([fn{1} '_train'])(cfg.preprocess_param{pp}.keep_idx_train,:,:,:);
+        end
+      end
+    elseif ~cfg.preprocess_param{pp}.is_train_set && isfield(cfg.preprocess_param{pp}, 'keep_idx_test')
+      for p = (pp+1):numel(cfg.preprocess_fun)
+        for fn = cfg.preprocess_param{p}.select_data
+          cfg.preprocess_param{p}.([fn{1} '_test'])  = cfg.preprocess_param{p}.([fn{1} '_test'])(cfg.preprocess_param{pp}.keep_idx_test,:,:,:);
+        end
+      end
+    end
+    
     % swap between train/test set
     cfg.preprocess_param{pp}.is_train_set = 1 - cfg.preprocess_param{pp}.is_train_set;
 
diff --git a/preprocess/mv_preprocess_pca.m b/preprocess/mv_preprocess_pca.m
index a4a8960..c3c7d74 100644
--- a/preprocess/mv_preprocess_pca.m
+++ b/preprocess/mv_preprocess_pca.m
@@ -16,7 +16,9 @@
 % .feature_dimension - which dimension codes the features (eg the channels)
 %                     (default 2)
 % .normalize        - if 1, all PCs are scaled to have variance 1 (default
-%                     1). This only works if X has a rank of at least n. 
+%                     1). This only works if X has a rank of at least n.
+% .omitrows         - if 1, rows with nans are omitted from the covariance
+%                     computation (default 0)
 %
 % Note: features x features covariance matrices are calculated across the
 % target dimension. A covariance matrix is calculated for every element of
@@ -70,11 +72,19 @@
     % to flip the matrix for the covariance calculation
     if f < t
         for ix = dim_loop
+          if pparam.omitrows
+            C = C + cov(squeeze(X(ix{:},:))', 'omitrows');
+          else
             C = C + cov(squeeze(X(ix{:},:))');
+          end
         end
     else
         for ix = dim_loop
+          if pparam.omitrows
+            C = C + cov(squeeze(X(ix{:},:)), 'omitrows');
+          else
             C = C + cov(squeeze(X(ix{:},:)) );
+          end
         end
     end
     
diff --git a/preprocess/mv_preprocess_undersample.m b/preprocess/mv_preprocess_undersample.m
index 3c178d3..339988a 100644
--- a/preprocess/mv_preprocess_undersample.m
+++ b/preprocess/mv_preprocess_undersample.m
@@ -20,7 +20,9 @@
 % since it does not introduce any dependencies between samples.
 
 if pparam.is_train_set || pparam.undersample_test_set
-    
+    % ensure column
+    clabel = clabel(:);
+
     sd = sort(pparam.sample_dimension(:))';
     uc = unique(clabel);
     nclasses = numel(uc);
@@ -35,19 +37,26 @@
     
     % undersample the majority class(es)
     rm_samples = abs(N - min(N));
+    ix_rm = zeros(0,1);
     for cc=1:nclasses
         if rm_samples(cc)>0
             ix_this_class = find(clabel == uc(cc));
-            ix_rm = randperm( numel(ix_this_class), rm_samples(cc));
-
-            % Remove samples from all sample dimensions
-            for rm_dim=sd
-                s_dim = s;
-                s_dim(rm_dim) = {ix_this_class(ix_rm)};
-                X(s_dim{:})= [];
-            end
-            clabel(ix_this_class(ix_rm))= [];
+            ix_rm = [ix_rm; ix_this_class(randperm( numel(ix_this_class), rm_samples(cc)))];
         end
     end
+
+    % Remove samples from all sample dimensions
+    for rm_dim=sd
+        s_dim = s;
+        s_dim(rm_dim) = {ix_rm};
+        X(s_dim{:})= [];
+    end
+    clabel(ix_rm)= [];
+ 
+    if pparam.is_train_set
+      pparam.keep_idx_train = setdiff(1:sum(N), ix_rm);
+    elseif pparam.undersample_test_set
+      pparam.keep_idx_test = setdiff(1:sum(N), ix_rm);
+    end
 end
 
diff --git a/preprocess/mv_preprocess_zscore.m b/preprocess/mv_preprocess_zscore.m
index 4f073e5..6f7b88d 100644
--- a/preprocess/mv_preprocess_zscore.m
+++ b/preprocess/mv_preprocess_zscore.m
@@ -15,6 +15,8 @@
 % .dimension        - which dimension(s) of the data matrix will be used
 %                     for z-scoring. Typically the dimension representing
 %                     the samples (default 1)
+% .omitnan          - boolean flag to omit nans or not (default 0). If the
+%                     data contains nans, it makes sense to use 1
 %
 % Nested preprocessing: For train data, preprocess_param.mean and 
 % preprocess_param.standard_deviation are calculated
@@ -22,8 +24,13 @@
 % obtained from the train data are used to scale the test data.
 
 if pparam.is_train_set
+  if pparam.omitnan
+    pparam.mean = mean(X, pparam.dimension, 'omitnan');
+    pparam.standard_deviation = std(X, [], pparam.dimension, 'omitnan');    
+  else
     pparam.mean = mean(X, pparam.dimension);
     pparam.standard_deviation = std(X, [], pparam.dimension);    
+  end
 end
 
 % Remove mean from data and divide by standard deviation
diff --git a/utils/mv_calculate_performance.m b/utils/mv_calculate_performance.m
index b73df78..d2467a3 100644
--- a/utils/mv_calculate_performance.m
+++ b/utils/mv_calculate_performance.m
@@ -279,7 +279,7 @@
         % https://en.wikipedia.org/wiki/Student%27s_t-test#Equal_or_unequal_sample_sizes.2C_equal_variance
         perf = cell(sz_output);
         
-         % Aggregate across samples, for each class separately
+        % Aggregate across samples, for each class separately
         if nextra == 1
             % Get means
             M1 = cellfun( @(cfo,lab) nanmean(cfo(lab==1,:,:,:,:,:,:,:,:,:),1), model_output, y, 'Un',0);
diff --git a/utils/mv_check_inputs.m b/utils/mv_check_inputs.m
index 6ed20dc..59d376e 100644
--- a/utils/mv_check_inputs.m
+++ b/utils/mv_check_inputs.m
@@ -22,13 +22,15 @@
 end
 
 if iscell(clabel) || ~all(ismember(clabel,1:n_classes))
-    warning('clabel should be a vector consisting of integers 1 (class 1), 2 (class 2), 3 (class 3) and so on. Relabelling them accordingly.');
+    clabelorig = clabel;
+    warning('clabel should almost always be a vector consisting of integers 1 (class 1), 2 (class 2), 3 (class 3) and so on. Relabelling them accordingly.');
     newlabel = nan(numel(clabel), 1);
     for i = 1:n_classes
         newlabel(ismember(clabel, u(i))) = i; % set to 1:nth classes
     end
     clabel = newlabel;
     if has_second_dataset
+        clabel2orig = clabel2;
         newlabel = nan(numel(clabel2), 1);
         for i = 1:n_classes
             newlabel(ismember(clabel2, u(i))) = i;
@@ -270,3 +272,10 @@
 
 %% cfg: set defaults for classifier hyperparameter
 cfg.hyperparameter = mv_get_hyperparameter(cfg.classifier, cfg.hyperparameter);
+if isfield(cfg.hyperparameter, 'relabel_design') && ~cfg.hyperparameter.relabel_design
+    warning('overruling the relabelling of the design, switching back to the original');
+    clabel = clabelorig;
+    if has_second_dataset
+      clabel2 = clabel2orig;
+    end
+end
diff --git a/utils/mv_get_hyperparameter.m b/utils/mv_get_hyperparameter.m
index 5ae78b2..55a93be 100644
--- a/utils/mv_get_hyperparameter.m
+++ b/utils/mv_get_hyperparameter.m
@@ -167,6 +167,23 @@
         mv_set_default(param,'lambda_y','auto');
         mv_set_default(param,'form', 'auto');
 
-
-    otherwise, error('Unknown model ''%s''',model)
+ 
+  otherwise
+        % check for the existence of a triplet of functions:
+        %   train_<model>
+        %   test_<model>
+        %   params_<model>
+        % if these conditions are met then it may be something that has
+        % been developed at home, let's give it a shot
+        if exist(sprintf('train_%s', model), 'file') == 2 && ...
+           exist(sprintf('test_%s', model), 'file') ==2 && ...
+           exist(sprintf('params_%s', model), 'file') == 2
+          try
+            param = feval(sprintf('params_%s', model), param);
+          catch
+            error('Failed to get default parameters for model %s', model);
+          end
+        else
+           error('Unknown model ''%s''',model)
+        end
 end
diff --git a/utils/mv_get_preprocess_param.m b/utils/mv_get_preprocess_param.m
index d568ddd..f7135ff 100644
--- a/utils/mv_get_preprocess_param.m
+++ b/utils/mv_get_preprocess_param.m
@@ -87,6 +87,7 @@
         mv_set_default(preprocess_param,'target_dimension',3);
         mv_set_default(preprocess_param,'normalize',1);
         mv_set_default(preprocess_param,'select_data',[]);
+        mv_set_default(preprocess_param,'omitrows',0);
     
     case 'impute_nan'
         mv_set_default(preprocess_param,'is_train_set',1);
@@ -119,4 +120,5 @@
         mv_set_default(preprocess_param,'is_train_set',1);
         mv_set_default(preprocess_param,'dimension',1);     
         mv_set_default(preprocess_param,'select_data',[]);
+        mv_set_default(preprocess_param,'omitnan',0);
 end
\ No newline at end of file