docs: def _perturbation(self):
          """
          Internal function for parameter initialization
          Returns Gaussian perturbation
          """
          if self.P>1:
              scales = []
              for term_i in range(self.n_randEffs):
                  _scales = sp.randn(self.diag[term_i].shape[0])
                  if self.jitter[term_i]>0:
                      _scales  = sp.concatenate((_scales,sp.zeros(1)))
                  scales.append(_scales)
              scales = sp.concatenate(scales)
          else:
              scales = sp.randn(self.vd.getNumberScales())
          return scales