88from tqdm import tqdm
99from MDAnalysis .analysis .base import Results
1010from basicrta .util import confidence_interval
11+ import multiprocessing
1112from multiprocessing import Pool , Lock
1213import MDAnalysis as mda
1314from basicrta import istarmap
@@ -47,8 +48,6 @@ def run(self, run_resids=None, g=100):
4748 :param run_resids: Resid(s) for which to run a Gibbs sampler.
4849 :type run_resids: int or list, optional
4950 """
50- from basicrta .util import run_residue
51-
5251 with open (self .contacts , 'r+b' ) as f :
5352 contacts = pickle .load (f )
5453
@@ -96,6 +95,38 @@ def run(self, run_resids=None, g=100):
9695 pass
9796
9897
98+ def run_residue (residue , time , proc , ncomp , niter , cutoff , g , from_combined = False ):
99+ """Run Gibbs sampler for a single residue.
100+
101+ :param residue: Residue name
102+ :type residue: str
103+ :param time: Residence times data
104+ :type time: array-like
105+ :param proc: Process number for progress bar positioning
106+ :type proc: int
107+ :param ncomp: Number of mixture components
108+ :type ncomp: int
109+ :param niter: Number of iterations
110+ :type niter: int
111+ :param cutoff: Cutoff value used in contact analysis
112+ :type cutoff: float
113+ :param g: Gibbs skip parameter
114+ :type g: int
115+ :param from_combined: Whether data comes from combined contacts
116+ :type from_combined: bool
117+ """
118+ x = np .array (time )
119+ if len (x ) != 0 :
120+ try :
121+ proc = int (multiprocessing .current_process ().name .split ('-' )[- 1 ])
122+ except ValueError :
123+ proc = 1
124+
125+ gib = Gibbs (times = x , residue = residue , loc = proc , ncomp = ncomp , niter = niter , cutoff = cutoff , g = g )
126+ gib ._from_combined_contacts = from_combined
127+ gib .run ()
128+
129+
99130class Gibbs (object ):
100131 r"""Gibbs sampler to estimate parameters of an exponential mixture for a set
101132 of data. Results are stored in :class:`gibbs.results`, which uses
@@ -301,8 +332,8 @@ def cluster(self, method="GaussianMixture", **kwargs):
301332 pindicator [tmpind , mapinds [i ]] += 1
302333
303334 pindicator = (pindicator .T / pindicator .sum (axis = 1 )).T
304- setattr ( self .processed_results , ' indicator' , pindicator )
305- setattr ( self .processed_results , ' labels' , all_labels )
335+ self .processed_results . indicator = pindicator
336+ self .processed_results . labels = all_labels
306337
307338 def process_gibbs (self , show = True ):
308339 r"""
@@ -329,9 +360,8 @@ def process_gibbs(self, show=True):
329360
330361 self .cluster (n_init = 117 , n_components = lmode )
331362 labels , presorts = mixture_and_plot (self , show = show )
332- setattr (self .processed_results , 'labels' , labels )
333- setattr (self .processed_results , 'indicator' ,
334- self .processed_results .indicator [:, presorts ])
363+ self .processed_results .labels = labels
364+ self .processed_results .indicator = self .processed_results .indicator [:, presorts ]
335365
336366 attrs = ["weights" , "rates" , "ncomp" , "residue" , "iteration" , "niter" ]
337367 values = [fweights , frates , lmode , self .residue , indices , self .niter ]
@@ -364,7 +394,7 @@ def _sample_indicator(self):
364394 # sample indicator
365395 s = np .argmax (rng .multinomial (1 , z ), axis = 1 )
366396 indicator [i ] = s
367- setattr ( self , ' indicator' , indicator )
397+ self . indicator = indicator
368398 return indicator [burnin_ind ::self .gskip ]
369399
370400 def save (self ):
@@ -719,8 +749,8 @@ def _estimate_params(self):
719749 params = np .array ([[wh [1 ][np .argmax (wh [0 ])], rh [1 ][np .argmax (rh [0 ])]]
720750 for wh , rh in zip (whists , rhists )])
721751
722- setattr ( rp , ' parameters' , params )
723- setattr ( rp , ' intervals' , np .array ([wbounds , rbounds ]) )
752+ rp . parameters = params
753+ rp . intervals = np .array ([wbounds , rbounds ])
724754
725755 def estimate_tau (self ):
726756 r"""
0 commit comments