code: def _model(self, beta):
          
  
          # Transform latent variables
          parm = np.array([self.latent_variables.z_list[k].prior.transform(beta[k]) for k in range(beta.shape[0])])
  
          xeps = np.power(self.data-parm[-1],2)
          Y = np.array(self.data[self.max_lag:])
          eps = np.power(Y-parm[-1],2)
          X = np.ones(Y.shape[0])
  
          # ARCH terms
          if self.q != 0:
              for i in range(0,self.q):   
                  X = np.vstack((X,xeps[(self.max_lag-i-1):-i-1]))
              sigma2 = np.matmul(np.transpose(X),parm[0:-self.p-1])
          else:
              sigma2 = np.transpose(X*parm[0])
  
          sigma2 = garch_recursion(parm, sigma2, self.q, self.p, Y.shape[0], self.max_lag)
  
          return np.array(sigma2), Y, eps
docs: Creates the structure of the model
  
          Parameters
          ----------
          beta : np.array
              Contains untransformed starting values for latent variables
  
          Returns
          ----------
          sigma2 : np.array
              Contains the values for the conditional volatility series
  
          Y : np.array
              Contains the length-adjusted time series (accounting for lags)
  
          eps : np.array
              Contains the squared residuals (ARCH terms) for the time series