docs: def crossValidation(self,seed=0,n_folds=10,fullVector=True,verbose=None,D=None,**keywords):
          """
          Split the dataset in n folds, predict each fold after training the model on all the others
  
          Args:
              seed:        seed
              n_folds:     number of folds to train the model on
              fullVector:  Bolean indicator, if true it stops if no convergence is observed for one of the folds, otherwise goes through and returns a pheno matrix with missing values
              verbose:     if true, prints the fold that is being used for predicitons
              **keywords:  params to pass to the function optimize
          Returns:
              Matrix of phenotype predictions [N,P]
          """
          verbose = dlimix.getVerbose(verbose)
  
          # split samples into training and test
          sp.random.seed(seed)
          r = sp.random.permutation(self.Y.shape[0])
          nfolds = 10
          Icv = sp.floor(((sp.ones((self.Y.shape[0]))*nfolds)*r)/self.Y.shape[0])
  
          RV = {}
          if self.P==1:  RV['var'] = sp.zeros((nfolds,self.n_randEffs))
          else:          RV['var'] = sp.zeros((nfolds,self.P,self.n_randEffs))
  
          Ystar = sp.zeros_like(self.Y)
  
          for fold_j in range(n_folds):
  
              if verbose:
                  print((".. predict fold %d"%fold_j))
  
              Itrain  = Icv!=fold_j
              Itest   = Icv==fold_j
              Ytrain  = self.Y[Itrain,:]
              Ytest   = self.Y[Itest,:]
              vc = VarianceDecomposition(Ytrain)
              vc.setTestSampleSize(Itest.sum())
              for term_i in range(self.n_fixedEffs):
                  F      = self.vd.getFixed(term_i)
                  Ftest  = F[Itest,:]
                  Ftrain = F[Itrain,:]
                  if self.P>1:    A = self.vd.getDesign(term_i)
                  else:           A = None
                  vc.addFixedEffect(F=Ftrain,Ftest=Ftest,A=A)
              for term_i in range(self.n_randEffs):
                  if self.P>1:
                      tct  = self.trait_covar_type[term_i]
                      rank = self.rank[term_i]
                      ftc  = self.fixed_tc[term_i]
                      jitt = self.jitter[term_i]
                      if tct=='lowrank_diag1' or tct=='freeform1':
                          d = D[fold_j,:,term_i]
                      else:
                          d = None
                  else:
                      tct  = None
                      rank = None
                      ftc  = None
                      jitt = None
                      d    = None
                  if term_i==self.noisPos:
                      vc.addRandomEffect(is_noise=True,trait_covar_type=tct,rank=rank,jitter=jitt,fixed_trait_covar=ftc,d=d)
                  else:
                      R = self.vd.getTerm(term_i).getK()
                      Rtrain = R[Itrain,:][:,Itrain]
                      Rcross = R[Itrain,:][:,Itest]
                      vc.addRandomEffect(K=Rtrain,Kcross=Rcross,trait_covar_type=tct,rank=rank,jitter=jitt,fixed_trait_covar=ftc,d=d)
              conv = vc.optimize(verbose=False,**keywords)
              if self.P==1:
                  RV['var'][fold_j,:] = vc.getVarianceComps()[0,:]
              else:
                  RV['var'][fold_j,:,:] = vc.getVarianceComps()
  
              if fullVector:
                  assert conv, 'VarianceDecompositon:: not converged for fold %d. Stopped here' % fold_j
              if conv:
                  Ystar[Itest,:] = vc.predictPhenos()
              else:
                  warnings.warn('not converged for fold %d' % fold_j)
                  Ystar[Itest,:] = sp.nan
  
          return Ystar,RV