changed defaults for gibbs.Gibbs to g=100, gskip=1

- updated gibbs.Gibbs
- updated tests

Author: orbeckst

CHANGELOG.md

@@ -25,7 +25,13 @@ The rules for this file:
   set the default to 100 (as in the paper) and ensure that the correct value
   is used as Gibbs.gskip (which is relative to the save skip step of Gibbs.g)
   (Issue #48)
-  
+
+### Changed
+* Default kwargs for the skipping in the Gibbs sampler are now
+  gibbs.Gibbs(g=100, gskip=1) (used to be g=50, gskip=2) but for most users
+  gskip for processing data is not important and it makes more sense to focus 
+  on g as the stride at which we sample AND process data (#48, PR #49)
+
 
 ## [1.1.2] - 2025-07-22
 

basicrta/gibbs.py

@@ -116,6 +116,23 @@ class Gibbs(object):
                    directory to load/save results. Allows for multiple cutoffs
                    to be tested in directory containing contacts.
     :type cutoff: float
+    :param g: Gibbs skip parameter for decorrelated samples;
+              only save every `g` samples from full Gibbs sampler chain;
+              default from https://pubs.acs.org/doi/10.1021/acs.jctc.4c01522
+              (NOTE: this value is called *gskip* in cluster.py)
+    :type g: int
+    :param burnin: Burn-in parameter, drop first `burnin` samples as equilibration;
+                   default from https://pubs.acs.org/doi/10.1021/acs.jctc.4c01522
+    :type burnin: int
+    :param gskip: Process data from the subsampled chain (ever `g` samples) at a 
+                  coarser skip interval of `gskip` samples. Thus, in total, samples
+                  are taken at ``g * gskip`` steps from the full chain.
+                  (This is useful for sensitivity analysis where we run the chain with 
+                  a small `g` value and save many samples and then use `gskip` to process
+                  samples at increasingly larger intervals without having to re-run the 
+                  chain.) The default value of 1 means that the samples are processed at
+                  every `g` samples from the full chain.
+    :type gskip: int
 
     EXAMPLE
     -------
@@ -139,7 +156,7 @@ class Gibbs(object):
     """
 
     def __init__(self, times=None, residue=None, loc=0, ncomp=15, niter=110000,
-                 cutoff=None, g=50, burnin=10000, gskip=2):
+                 cutoff=None, g=100, burnin=10000, gskip=1):
         self.times = times
         self.residue = residue
         self.niter = niter

basicrta/tests/test_combine_contacts.py

@@ -265,6 +265,8 @@ def test_gibbs_sampler_integration(self, tmp_path, create_mock_contacts):
             ncomp=2,  # Use 2 components for stability
             niter=1000,  # 1000 steps as requested
             burnin=5,   # 5 burnin steps as requested
+            g=50,
+            gskip=1,
             cutoff=7.0  # Set cutoff for directory creation
         )
 

basicrta/tests/test_gibbs.py

@@ -127,6 +127,16 @@ class TestGibbsSampler:
             'cutoff': 7.0,
             'g': 100
         },
+        {
+            'times': None,  # Will be set from fixture
+            'residue': 'W313',
+            'ncomp': 2,
+            'niter': 1000,
+            'burnin': 5,
+            'cutoff': 7.0,
+            'g': 50,
+            'gskip': 2
+        },        
     ])
     def test_gibbs_run_method(self, tmp_path, synthetic_timeseries, init_kwargs):
         """Test the run() method for Gibbs class with synthetic data."""