WIP: create Universes with Merge in encore

package/CHANGELOG

@@ -20,6 +20,10 @@ mm/dd/yy richardjgowers, kain88-de, lilyminium, p-j-smith, bdice, joaomcteixeira
   * 0.21.0
 
 Fixes
+  * Utility functions in encore now use mda.Merge to create universes
+    instead of mda.Universe (Issue #2576). 
+  * prepare_ensembles_for_convergence_increasing_window() transfers 
+    universes to memory before calling timeseries (Issue #2578)
   * Correct args order of base.AnalysisFromFunction (Issue #2503)  
   * encore.dres() returns dimensionality reduction details instead of a
     reference to itself (Issue #2471)

package/MDAnalysis/analysis/encore/bootstrap.py

@@ -144,15 +144,12 @@ def get_ensemble_bootstrap_samples(ensemble,
     """
 
     ensemble.transfer_to_memory()
+    coords = ensemble.trajectory.timeseries(order='fac')
+    n_frames = coords.shape[0]
 
     ensembles = []
     for i in range(samples):
-        indices = np.random.randint(
-            low=0,
-            high=ensemble.trajectory.timeseries().shape[1],
-            size=ensemble.trajectory.timeseries().shape[1])
-        ensembles.append(
-            mda.Universe(ensemble.filename,
-                        ensemble.trajectory.timeseries(order='fac')[indices,:,:],
-                         format=mda.coordinates.memory.MemoryReader))
+        indices = np.random.randint(low=0, high=n_frames, size=n_frames)
+        ensembles.append(ensemble.copy().load_new(coords[indices]))
+
     return ensembles

package/MDAnalysis/analysis/encore/similarity.py

@@ -661,7 +661,7 @@ def write_output(matrix, base_fname=None, header="", suffix="",
     matrix.square_print(header=header, fname=fname)
 
 
-def prepare_ensembles_for_convergence_increasing_window(ensemble,
+def prepare_ensembles_for_convergence_increasing_window(universe,
                                                         window_size,
                                                         select="name CA"):
     """
@@ -672,7 +672,7 @@ def prepare_ensembles_for_convergence_increasing_window(ensemble,
     Parameters
     ----------
 
-    ensemble : :class:`~MDAnalysis.core.universe.Universe` object
+    universe : :class:`~MDAnalysis.core.universe.Universe` object
         Input ensemble
 
     window_size : int
@@ -684,35 +684,30 @@ def prepare_ensembles_for_convergence_increasing_window(ensemble,
     Returns
     -------
 
-    tmp_ensembles :
+    ensembles :
         The original ensemble is divided into different ensembles, each being
         a window_size-long slice of the original ensemble. The last
         ensemble will be bigger if the length of the input ensemble
         is not exactly divisible by window_size.
 
     """
 
-    ens_size = ensemble.trajectory.timeseries(ensemble.select_atoms(select),
-                                              order='fac').shape[0]
-
-    rest_slices = ens_size // window_size
-    residuals = ens_size % window_size
-    slices_n = [0]
-
-    tmp_ensembles = []
-
-    for rs in range(rest_slices - 1):
-        slices_n.append(slices_n[-1] + window_size)
-    slices_n.append(slices_n[-1] + residuals + window_size)
-
-    for s,sl in enumerate(slices_n[:-1]):
-        tmp_ensembles.append(mda.Universe(
-            ensemble.filename,
-            ensemble.trajectory.timeseries(order='fac')
-            [slices_n[s]:slices_n[s + 1], :, :],
-            format=MemoryReader))
-
-    return tmp_ensembles
+    n_frames = len(universe.trajectory)
+    slice_indices = list(range(0, n_frames, window_size))
+    slice_indices.append(n_frames)
+
+    # gets transferred to memory later anyway
+    universe.transfer_to_memory()
+    ag = universe.select_atoms(select)
+    coordinates = universe.trajectory.timeseries(ag, order='fac')
+    i, j = slice_indices[:2]
+    nu = mda.Merge(ag).load_new(coordinates[i:j])
+    ensembles = [nu]
+    for i, j in zip(slice_indices[1:-1], slice_indices[2:]):
+        nu_ = nu.copy().load_new(coordinates[i:j])
+        ensembles.append(nu_)
+
+    return ensembles
 
 
 def hes(ensembles,
@@ -1102,7 +1097,8 @@ def ces(ensembles,
 
     # Calculate distance matrix if not provided
     if any_method_accept_distance_matrix and not distance_matrix:
-        distance_matrix = get_distance_matrix(merge_universes(ensembles),
+        merged = merge_universes(ensembles, select=select)
+        distance_matrix = get_distance_matrix(merged,
                                               select=select,
                                               ncores=ncores)
     if estimate_error:
@@ -1378,7 +1374,8 @@ def dres(ensembles,
 
     # Calculate distance matrix if not provided
     if any_method_accept_distance_matrix and not distance_matrix:
-        distance_matrix = get_distance_matrix(merge_universes(ensembles),
+        merged = merge_universes(ensembles, select=select)
+        distance_matrix = get_distance_matrix(merged,
                                               select=select,
                                               ncores=ncores)
     if estimate_error:

package/MDAnalysis/analysis/encore/utils.py

@@ -386,25 +386,37 @@ def trm_indices_diag(n):
             yield (i, j)
 
 
-def merge_universes(universes):
+def merge_universes(universes, select='name CA'):
     """
     Merge list of universes into one
 
     Parameters
     ----------
     universes : list of Universe objects
 
+    select: str
+        atoms to select
+
 
     Returns
     ----------
     Universe object
     """
-
-    for universe in universes:
-        universe.transfer_to_memory()
-
-    return mda.Universe(
-        universes[0].filename,
-        np.concatenate(tuple([e.trajectory.timeseries(order='fac') for e in universes]),
-                       axis=0),
-        format=MemoryReader)
+    if not universes:
+        raise ValueError('No universes passed')
+
+    coord_list = []
+    for u in universes:
+        u.transfer_to_memory()
+        ag = u.select_atoms(select)
+        coords_ = u.trajectory.timeseries(asel=ag, order='fac')
+        coord_list.append(coords_)
+
+    try:
+        coords = np.concatenate(coord_list, axis=0)
+    except ValueError:
+        err = 'atom selection "{}" results in mismatching AtomGroups'
+        raise ValueError(err.format(select))
+
+    ag = universes[0].select_atoms(select)
+    return mda.Merge(ag).load_new(coords)

testsuite/MDAnalysisTests/analysis/test_encore.py

@@ -58,17 +58,13 @@ def ens2_template(self):
 
     @pytest.fixture()
     def ens1(self, ens1_template):
-        return mda.Universe(
-            ens1_template.filename,
-            ens1_template.trajectory.timeseries(order='fac'),
-            format=mda.coordinates.memory.MemoryReader)
+        coords = ens1_template.trajectory.timeseries(order='fac')
+        return mda.Merge(ens1_template.atoms).load_new(coords)
 
     @pytest.fixture()
     def ens2(self, ens2_template):
-        return mda.Universe(
-            ens2_template.filename,
-            ens2_template.trajectory.timeseries(order='fac'),
-            format=mda.coordinates.memory.MemoryReader)
+        coords = ens2_template.trajectory.timeseries(order='fac')
+        return mda.Merge(ens2_template.atoms).load_new(coords)
 
     def test_triangular_matrix(self):
         scalar = 2