Ayman here!!!.  some peices that may be useful:
        self.parcov = parameter covariance matrix
        self.obscov = obseravtion noise covariance matrix
        self.parensmeble = current parameter ensemble
        self.obsensemble = current observation (model output) ensemble

        the ensemble instances have two useful methods:
        Ensemble.covariance_matrix() - get an empirical covariance matrix
        Ensemble.get_deviations() - get an ensemble of deviations around the mean vector

        If you use pyemu.Matrix (and pyemu.Cov) for the linear algebra parts, you don't
        have to worry about the alignment of the matrices (pyemu will dynamically reorder/align
        based in row and/or col names).  The Matrix instance also has a .inv for the inverse
        as well as .s, .u and .v for the SVD components (gets dynamically evaluated if you try to
        access these attributes)



        Following Evensen 2003..

def analysis_evensen(self):
        """Ayman here!!!.  some peices that may be useful:
        self.parcov = parameter covariance matrix
        self.obscov = obseravtion noise covariance matrix
        self.parensmeble = current parameter ensemble
        self.obsensemble = current observation (model output) ensemble

        the ensemble instances have two useful methods:
        Ensemble.covariance_matrix() - get an empirical covariance matrix
        Ensemble.get_deviations() - get an ensemble of deviations around the mean vector

        If you use pyemu.Matrix (and pyemu.Cov) for the linear algebra parts, you don't
        have to worry about the alignment of the matrices (pyemu will dynamically reorder/align
        based in row and/or col names).  The Matrix instance also has a .inv for the inverse
        as well as .s, .u and .v for the SVD components (gets dynamically evaluated if you try to
        access these attributes)



        Following Evensen 2003..
        """

        nz_names = self.pst.nnz_obs_names
        nreals = self.obsensemble.shape[0]
        self.parensemble._transform()
        # nonzero weighted state deviations
        HA_prime = self.obsensemble.get_deviations().loc[:,nz_names].T

        # obs noise pertubations - move to constuctor
        E = (self.obsensemble_0.loc[:,nz_names] - self.pst.observation_data.obsval.loc[nz_names]).T
        #print(E)

        # innovations:  account for any failed runs (use obsensemble index)
        D_prime = (self.obsensemble_0.loc[self.obsensemble.index,nz_names] - self.obsensemble.loc[:,nz_names]).T

        ES = HA_prime.loc[nz_names,E.columns] + E.loc[nz_names,:]
        assert ES.shape == ES.dropna().shape

        ES = pyemu.Matrix.from_dataframe(ES)

        nrmin = min(self.pst.nnz_obs,self.obsensemble.shape[0])
        U,s,v = ES.pseudo_inv_components(maxsing=nrmin,
                                        eigthresh=ES.get_maxsing(self.pst.svd_data.eigthresh),
                                        truncate=True)

        #half assed inverse
        s_inv = s.T
        for i,sval in enumerate(np.diag(s.x)):
           if sval == 0.0:
               break
           s_inv.x[i,i] = 1.0 / (sval * sval)

        X1 = s_inv * U.T
        X1.autoalign = False #since the row/col names don't mean anything for singular components and everything is aligned

        X2 =  X1 * D_prime #these are aligned through the use of nz_names

        X3 = U * X2 #also aligned

        X4 = pyemu.Matrix.from_dataframe(HA_prime.T) * X3

        I = np.identity(X4.shape[0])
        X5 = X4 + I
        #print(X5.x.sum(axis=1))


        # deviations of adj pars
        A_prime = pyemu.Matrix.from_dataframe(self.parensemble.get_deviations().loc[:,self.pst.adj_par_names]).T


        upgrade = (A_prime * X4).to_dataframe().T

        assert upgrade.shape == upgrade.dropna().shape

        upgrade.index = self.parensemble.index
        #print(upgrade)
        parensemble = self.parensemble + upgrade

        parensemble = pyemu.ParameterEnsemble.from_dataframe(df=parensemble,pst=self.pst,istransformed=True)

        assert parensemble.shape == parensemble.dropna().shape
        return parensemble