code: def _msge_with_gradient_underdetermined(data, delta, xvschema, skipstep, p):
      
      t, m, l = data.shape
      d = None
      j, k = 0, 0
      nt = np.ceil(t / skipstep)
      for trainset, testset in xvschema(t, skipstep):
  
          a, b = _construct_var_eqns(atleast_3d(data[trainset, :, :]), p)
          c, d = _construct_var_eqns(atleast_3d(data[testset, :, :]), p)
  
          e = sp.linalg.inv(np.eye(a.shape[0]) * delta ** 2 + a.dot(a.T))
  
          cc = c.transpose().dot(c)
  
          be = b.transpose().dot(e)
          bee = be.dot(e)
          bea = be.dot(a)
          beea = bee.dot(a)
          beacc = bea.dot(cc)
          dc = d.transpose().dot(c)
  
          j += np.sum(beacc * bea - 2 * bea * dc) + np.sum(d ** 2)
          k += np.sum(beea * dc - beacc * beea) * 4 * delta
  
      return j / (nt * d.size), k / (nt * d.size)
docs: Calculate mean squared generalization error and its gradient for
      underdetermined equation system.