Feat normtransf optional extrapolation (#528)

* added optional extrapolation to .prep_for_dsi()

* hotfix: remove restrictions on randrealgen

* added ztz approach; tweaked tests

* update tests

* le sigh...

* docstrings

* tweak to nst interpoaltion

* hotfix to update dsi forward run

* correction to quadratic extrapolation

* more updates to pyemu eds

* hack to ensure unique values during nstransform

* add flag to skip nst for cases where all obs values are equal

* cleaner dsi tests

* checkin

* hotfix for handling cases where all obs are equal

* test to catch failing dsi runs

* hotfix to energy threshold

* handle uneven nreal in randrealgen

Author: rhugman

autotest/la_tests.py

@@ -607,7 +607,8 @@ def ends_freyberg_test(tmp_path):
     ends = pyemu.EnDS(pst=pst, sim_ensemble=oe, obscov=os.path.join(test_d, "obs.unc"),predictions=predictions)
 
 
-def ends_freyberg_dsi_test(tmp_path):
+
+def ends_run_freyberg_dsi(tmp_path,nst=False,nst_extrap=None,ztz=False,energy=1.0):
     import pyemu
     import os 
     test_d = "ends_master"
@@ -621,58 +622,44 @@ def ends_freyberg_dsi_test(tmp_path):
     oe_name = pst_name.replace(".pst", ".0.obs.csv")
     oe = pyemu.ObservationEnsemble.from_csv(pst=pst, filename=oe_name).iloc[:100, :]
 
-    
     ends = pyemu.EnDS(pst=pst, sim_ensemble=oe,verbose=True)
     t_d = os.path.join(tmp_path,"dsi_template")
-    ends.prep_for_dsi(t_d=t_d)
+
+    ends.prep_for_dsi(t_d=t_d,
+                      use_ztz=ztz,
+                      apply_normal_score_transform=nst,
+                      nst_extrap=nst_extrap,
+                      energy=energy)
     # copy exe to dsi_template
     #shutil.copy2(os.path.join(test_d,"pestpp-ies.exe"),os.path.join(t_d,"pestpp-ies.exe"))
-
+    filename=os.path.join(t_d,"dsi.0.obs.csv")
+    if os.path.exists(filename):
+        os.remove(filename)
     pst = pyemu.Pst(os.path.join(t_d,"dsi.pst"))
-    pst.control_data.noptmax = 0
+    pst.control_data.noptmax = -1
     pst.write(os.path.join(t_d,"dsi.pst"),version=2)
-    #pyemu.os_utils.run("pestpp-ies dsi.pst",cwd="dsi_template")
-    m_d = os.path.join(tmp_path,"master_dsi")
-    pyemu.os_utils.start_workers(t_d,"pestpp-ies","dsi.pst",num_workers=15,worker_root=tmp_path,
-                                 master_dir=m_d)
+    pyemu.os_utils.run("pestpp-ies dsi.pst",cwd=t_d)
 
-    # run test wtih truncated svd
-    ends.prep_for_dsi(t_d=t_d,truncated_svd=True)
-    # copy exe to dsi_template
-    #shutil.copy2(os.path.join(test_d,"pestpp-ies.exe"),os.path.join(t_d,"pestpp-ies.exe"))
+    #read in the results
+    oe = pyemu.ObservationEnsemble.from_csv(pst=pst, filename=os.path.join(t_d,"dsi.0.obs.csv"))
+    assert oe.shape[0]==50, f"{50-oe.shape} failed runs"
+    phi_vector = oe.phi_vector.sort_values().values
+    assert phi_vector[0] != phi_vector[1],phi_vector
 
-    pst = pyemu.Pst(os.path.join(t_d,"dsi.pst"))
-    pst.control_data.noptmax = 0
-    pst.write(os.path.join(t_d,"dsi.pst"),version=2)
-    pyemu.os_utils.start_workers(t_d,"pestpp-ies","dsi.pst",num_workers=15,worker_root=tmp_path,
-                                 master_dir=m_d)
+def ends_freyberg_dsi_test(tmp_path):
+    ends_run_freyberg_dsi(tmp_path)
 
+def ends_freyberg_dsi_nst_test(tmp_path):
+    ends_run_freyberg_dsi(tmp_path,nst=True,nst_extrap=None)
 
-    # run test wtih normal score transform
-    ends.prep_for_dsi(t_d=t_d,apply_normal_score_transform=True)
-    # copy exe to dsi_template
-    #shutil.copy2(os.path.join(test_d,"pestpp-ies.exe"),os.path.join(t_d,"pestpp-ies.exe"))
-    
-    pst = pyemu.Pst(os.path.join(t_d,"dsi.pst"))
-    pst.control_data.noptmax = 0
-    pst.write(os.path.join(t_d,"dsi.pst"),version=2)
-    pyemu.os_utils.start_workers(t_d,"pestpp-ies","dsi.pst",num_workers=15,worker_root=tmp_path,
-                                 master_dir=m_d)
-    
-    # run test with log-transform
-    pst = pyemu.Pst(pst_name)
-    pst.pestpp_options["predictions"] = predictions
-    pst.observation_data["obstransform"] = "log"
-    ends = pyemu.EnDS(pst=pst, sim_ensemble=oe,verbose=True)
-    ends.prep_for_dsi(t_d=t_d,apply_normal_score_transform=False)
-    # copy exe to dsi_template
-    #shutil.copy2(os.path.join(test_d,"pestpp-ies.exe"),os.path.join(t_d,"pestpp-ies.exe"))
-    
-    pst = pyemu.Pst(os.path.join(t_d,"dsi.pst"))
-    pst.control_data.noptmax = 0
-    pst.write(os.path.join(t_d,"dsi.pst"),version=2)
-    pyemu.os_utils.start_workers(t_d,"pestpp-ies","dsi.pst",num_workers=15,worker_root=tmp_path,
-                                 master_dir=m_d)
+def ends_freyberg_dsi_extrap_test(tmp_path):
+    ends_run_freyberg_dsi(tmp_path,nst=True,nst_extrap='quadratic')
+
+def ends_freyberg_dsi_ztz_test(tmp_path):
+    ends_run_freyberg_dsi(tmp_path,ztz=True)
+
+def ends_freyberg_dsi_svd_test(tmp_path):
+    ends_run_freyberg_dsi(tmp_path,ztz=True,energy=0.999)
 
 
 def plot_freyberg_dsi():
@@ -738,15 +725,28 @@ def dsi_normscoretransform_test():
     oe = pyemu.ObservationEnsemble.from_csv(pst=pst, filename=oe_name).iloc[:100, :]
 
     nstval = randrealgen_optimized(oe.shape[0], 1e-7, 1e4)
+    window_size=3   
+    if oe.shape[0]>40:
+        window_size=5                    
+    if oe.shape[0]>90:
+        window_size=7
+    if oe.shape[0]>200:
+        window_size=9            
     for name in oe.columns:
         print("transforming:",name)
         sorted_values = oe._df.loc[:,name].sort_values().copy()
+        #if all values are the same, skip
+        if sorted_values.iloc[0] == sorted_values.iloc[-1]:
+            print("all values are the same, skipping")
+            continue
+        sorted_values.loc[:] = pyemu.eds.moving_average_with_endpoints(sorted_values.values, window_size)
         transformed_values = np.asarray([normal_score_transform(nstval, sorted_values, value)[0] for value in sorted_values])
         backtransformed_values = np.asarray([inverse_normal_score_transform(nstval, sorted_values, value)[0] for value in transformed_values])
 
         diff = backtransformed_values-sorted_values
         assert max(abs(diff))<1e-7, backtransformed_values
 
+
 if __name__ == "__main__":
     #dsi_normscoretransform_test()
     #ends_freyberg_dev()