code: def _backward_step(self, deltaY, deltaY_sigma2inv_rho_weightT,
                         sigma2_e, weight, mu, mu_Gamma_inv, Gamma_inv,
                         Lambda_0, Lambda_1, H):
          
          n_T = len(Gamma_inv)
          # All the terms with hat before are parameters of posterior
          # distributions of X conditioned on data from all time points,
          # whereas the ones without hat calculated by _forward_step
          # are mean and covariance of posterior of X conditioned on
          # data up to the time point.
          Gamma_inv_hat = [None] * n_T
          mu_Gamma_inv_hat = [None] * n_T
          mu_hat = [None] * n_T
  
          mu_hat[-1] = mu[-1].copy()
          mu_Gamma_inv_hat[-1] = mu_Gamma_inv[-1].copy()
          Gamma_inv_hat[-1] = Gamma_inv[-1].copy()
  
          for t in np.arange(n_T - 2, -1, -1):
              tmp = np.linalg.solve(Gamma_inv_hat[t + 1] - Gamma_inv[t + 1]
                                    + Lambda_1, H)
              Gamma_inv_hat[t] = Gamma_inv[t] + Lambda_0 - np.dot(H.T, tmp)
              mu_Gamma_inv_hat[t] = mu_Gamma_inv[t] \
                  - deltaY_sigma2inv_rho_weightT[t, :] + np.dot(
                      mu_Gamma_inv_hat[t + 1] - mu_Gamma_inv[t + 1]
                      + np.dot(deltaY[t, :] / sigma2_e, weight.T), tmp)
              mu_hat[t] = np.linalg.solve(Gamma_inv_hat[t],
                                          mu_Gamma_inv_hat[t])
          return mu_hat, mu_Gamma_inv_hat, Gamma_inv_hat
docs: backward step for HMM, assuming both the hidden state and noise
              have 1-step dependence on the previous value.