1import os
2import numpy, scipy
3import scipy.io
4import pysces
5from pysces.PyscesModelMap import ModelMap
6from pylab import figure, close, show, ioff, figtext
7from matplotlib.transforms import offset_copy
8import io, string
9from colorsys import hsv_to_rgb as htor
10from matplotlib.figure import Figure as mplfigure
11import warnings
12
13try:
14    input = raw_input  # Py2 compatibility
15except NameError:
16    pass
17
18
19doc = '''
20
21************************************************************************
22*           Welcome to the RateChar add-on module for PySCeS           *
23* Does generic supply-demand parameter scans around internal           *
24* metabolites species of a model and optionally plots the results.     *
25*                                                                      *
26* Copyright(C) J.M Rohwer, C.D. Christensen, J.-H.S. Hofmeyr, 2007-2011*
27* Triple-J Group for Molecular Cell Physiology                         *
28* Stellenbosch University, South Africa                                *
29* RateChar is distributed under the PySCeS (BSD style) licence, see    *
30* LICENCE.txt (supplied with this release) for details                 *
31************************************************************************
32
33'''
34
35
36class RateChar:
37    """Main class for doing generic supply-demand rate characteristics"""
38
39    def __init__(self, mod, min_concrange_factor=100, max_concrange_factor=100):
40        """
41        __init__(mod)
42
43        RateChar Arguments:
44        ==========
45        mod: instantiated PySCeS model
46        """
47
48        print(doc)
49        # default scan parameter values
50        self.min_concrange_factor = min_concrange_factor
51        self.max_concrange_factor = max_concrange_factor
52        self.min_fluxrange_factor = 2
53        self.max_fluxrange_factor = 2
54        self.scan_points = 256
55
56        # initialise model
57        self._out_base_dir = pysces.PyscesModel.OUTPUT_DIR + '/ratechar_out'
58        self._psc_mod = mod
59        self._psc_file = mod.ModelFile
60        self._mod = mod.ModelFile[:-4]
61        self.createModelOutputDir()
62
63    def createModelOutputDir(self):
64        """
65        createModelOutputDir()
66
67        Creates a directory based on the model filename in the output base directory
68        for all model output. For mulitple simulations with different models
69        a directory for each model is created.
70
71        Arguments:
72        =========
73        None
74        """
75        out_base_dir_status = os.path.exists(self._out_base_dir)
76        if not out_base_dir_status:
77            os.mkdir(self._out_base_dir)
78        model_output_dir = self._out_base_dir + '/' + self._mod
79        model_output_dir_status = os.path.exists(model_output_dir)
80        if not model_output_dir_status:
81            os.mkdir(model_output_dir)
82        self._psc_out_dir = model_output_dir
83        # pysces.PyscesModel.output_dir = model_output_dir
84
85    def pscFile2str(self, name):
86        """
87        pscFile2str(name)
88
89        Read psc file and return as string
90
91        Arguments:
92        ==========
93        name - string containing filename
94        """
95        if name[-4:] != '.psc':
96            name += '.psc'
97        F = open(os.path.join(pysces.PyscesModel.MODEL_DIR, name), 'r')
98        fstr = F.read()
99        F.close()
100        return fstr
101
102    def pscMod2str(self, mod):
103        """
104        pscMod2str(name)
105
106        Write PySCeS model out to string
107
108        Arguments:
109        ==========
110        mod - instantiated PySCeS model
111        """
112        F = io.StringIO()
113        mod.showModel(filename=F)
114        fstr = F.getvalue()
115        F.close()
116        return fstr
117
118    def stripFixed(self, fstr):
119        """
120        stripFixed(fstr)
121
122        Take a psc file string and return (Fhead, stripped_fstr)
123        where Fhead is the file header containing the "FIX: " line
124        and stripped_fstr is remainder of file as string
125
126        Arguments:
127        ==========
128        fstr - string representation of psc file
129
130        See also:
131        =========
132        pscFile2str(name)
133        """
134        Fi = io.StringIO()
135        Fi.write(fstr)
136        Fi.seek(0)
137        Fo = io.StringIO()
138        Fhead = None
139        for line in Fi:
140            if line[:4] == "FIX:":
141                Fhead = string.strip(line)
142                Fo.write('\n')
143            else:
144                Fo.write(line)
145        Fo.seek(0)
146        return Fhead, Fo.read()
147
148    def augmentFixString(self, OrigFix, fix):
149        """
150        augmentFixString(OrigFix, fix)
151
152        Add fix to FixString
153
154        Arguments:
155        ==========
156        OrigFix - original FixString
157        fix     - additional species to add to FixString
158        """
159        return OrigFix + ' %s' % fix
160
161    def evalBaseMod(self, show=1, file=0):
162        """
163        evalBaseMod(show=1,file=0)
164
165        Evaluates base model in terms of steady state and MCA.
166        Prints steady state values, elasticities and control coefficients
167        to screen and optionally writes to file.
168
169        Arguments:
170        ==========
171        show (bool) - print steady-state results on screen (default=1)
172        file (bool) - write steady-state results to file (default=0)
173        """
174        basemod = self._psc_mod
175        basemod.mode_elas_deriv = 1
176        basemod.doMca()
177        if show:
178            basemod.showState()
179            basemod.showEvar()
180            basemod.showCC()
181        if file:
182            rFile = open(os.path.join(self._psc_out_dir, self._mod + '_stst.txt'), 'w')
183            basemod.showState(rFile)
184            basemod.showEvar(rFile)
185            basemod.showCC(rFile)
186            rFile.write('\n')
187            rFile.close()
188
189    def setMaxConcRangeFactor(self, max_concrange_factor):
190        """
191        Changes the maximum concentration factor
192        """
193        self.max_concrange_factor = max_concrange_factor
194
195    def setMinConcRangeFactor(self, min_concrange_factor):
196        """
197        Changes the Minimum concentration factor
198        """
199        self.min_concrange_factor = min_concrange_factor
200
201    def setScanPoints(self, scan_points):
202        self.scan_points = scan_points
203
204    def setDefaults(self):
205        self.min_concrange_factor = 100
206        self.scan_points = 256
207        self.max_concrange_factor = 100
208
209    def doRateChar(self, fixed, summary=0, solver=0):
210        """
211        doRateChar(fixed)
212
213        Does rate characteristic analysis around single species (fixed).
214        Stores data in instance of RateCharData class which is attached
215        as attribute 'fixed' to current RateChar instance.
216
217        Arguments:
218        ==========
219        fixed (str)    - species around which to generate rate characteristic
220        summary (bool) - write summary to output file (default 0)
221        solver (int)   - steady-state solver to use: 0 - hybrd, 1 - fintslv, 2 - NLEQ2
222        """
223        assert solver in (0, 1, 2), "\nSolver mode can only be one of (0, 1, 2)!"
224        # if summary:
225        # assert hasattr(self, 'summaryfile'), "\nSummary mode can only be called from doAllRateChar method."
226        # base steady state
227        basemod = self._psc_mod
228        basemod.doMca()
229        assert fixed in basemod.species, "\nInvalid fixed species."
230        basemod_fstr = self.pscMod2str(basemod)
231        basemod.mode_elas_deriv = 1
232        mymap = ModelMap(basemod)
233        OrigFix, Fstr = self.stripFixed(basemod_fstr)
234        print(Fstr)
235        h1 = self.augmentFixString(OrigFix, fixed)
236        model1 = self._mod + '_' + fixed
237        mod = pysces.model(model1, loader="string", fString=h1 + '\n' + Fstr)
238        # mod.doLoad()   # no longer needed as of PySCeS 0.7.0
239        mod.SetQuiet()
240        mod.mode_solver = solver
241        mod.doState()
242        mod.mode_elas_deriv = 1
243        setattr(mod, fixed, getattr(basemod, fixed + '_ss'))
244        mod.doMcaRC()
245        # get fluxes directly linked to fixed species
246        FixSubstrateFor = ['J_' + r for r in getattr(mymap, fixed).isSubstrateOf()]
247        FixProductFor = ['J_' + r for r in getattr(mymap, fixed).isProductOf()]
248        print('Fixed species is Substrate for: ', FixSubstrateFor)
249        print('Fixed species is Product for:   ', FixProductFor)
250        # set scan distances
251        scan_min = getattr(basemod, fixed + '_ss') / self.min_concrange_factor
252        scan_max = getattr(basemod, fixed + '_ss') * self.max_concrange_factor
253        # do scan
254        rng = scipy.logspace(
255            scipy.log10(scan_min), scipy.log10(scan_max), self.scan_points
256        )
257        sc1 = pysces.Scanner(mod)
258        sc1.quietRun = True
259        sc1.printcnt = 50
260        sc1.addScanParameter(fixed, scan_min, scan_max, self.scan_points, log=True)
261        args = [fixed] + FixProductFor + FixSubstrateFor
262        print('Scanner user output: ', args)
263        sc1.addUserOutput(*args)
264        sc1.Run()
265        print(sc1.UserOutputResults)
266
267        if len(sc1.invalid_state_list) != 0:
268            newUserOutputResults = []
269            midpart = 0
270            # the number of results between invalid steady states at the ends
271            resultsperdimension = sc1.UserOutputResults.shape[1]
272            startpoint = 0
273            endpoint = len(sc1.UserOutputResults)
274
275            # start checking for invalid states from the start
276            # move startpoint forward for each NaN
277            # break at first non-NaN result
278            loopdone = False
279            for resultarray in sc1.UserOutputResults:
280                if loopdone:
281                    break
282                for result in resultarray:
283                    if numpy.isnan(result):
284                        startpoint = startpoint + 1
285                        loopdone = False
286
287                        break
288                    else:
289                        loopdone = True
290
291            # print endpoint
292            # checking for invalid states from the end
293            # move endpoint backward for each NaN
294            # break at first non-NaN result
295            loopdone = False
296            for i in range(len(sc1.UserOutputResults)):
297                if loopdone:
298                    break
299                for result in sc1.UserOutputResults[-i - 1]:
300                    if numpy.isnan(result):
301                        endpoint = endpoint - 1
302                        loopdone = False
303
304                        break
305                    else:
306                        loopdone = True
307            # number of results between NaN ends
308            midpart = endpoint - startpoint
309            newUserOutputResults = sc1.UserOutputResults[startpoint:endpoint]
310            newUserOutputResults = scipy.array(newUserOutputResults).reshape(
311                midpart, resultsperdimension
312            )
313            sc1.UserOutputResults = newUserOutputResults
314
315        # capture result by instantiating RateCharData class
316        setattr(
317            self,
318            fixed,
319            RateCharData(
320                self,
321                mod,
322                fixed,
323                args,
324                sc1.UserOutputResults,
325                FixProductFor,
326                FixSubstrateFor,
327            ),
328        )
329        if summary:
330            # test for presence of global summary filename, if present append to it
331            # if absent, write new summary file for single species scan
332            flag = 0
333            if hasattr(self, 'summaryfile'):
334                of = self.summaryfile
335            else:
336                of = open(
337                    os.path.join(
338                        self._psc_out_dir, self._mod + '_' + fixed + '_values.txt'
339                    ),
340                    'w',
341                )
342                flag = 1
343                of.write('# Base model: ' + self._psc_file + '\n')
344                of.write('####################################\n\n')
345            of.write('# Rate Characteristic scan around species: ' + fixed + '\n')
346            of.write('#######################################################\n')
347            of.write('#### Scan Min/Max ####\n')
348            of.write('Scan min: ' + str(scan_min) + '\n')
349            of.write('Scan max: ' + str(scan_max) + '\n')
350            # response coefficients
351            of.write('#### Response coefficients ####\n')
352            for step in (
353                getattr(mymap, fixed).isReagentOf()
354                + getattr(mymap, fixed).isModifierOf()
355            ):
356                of.writelines(
357                    '\\rc{'
358                    + step
359                    + '}{'
360                    + fixed
361                    + '} = '
362                    + str(mod.rc.get(step, fixed))
363                    + '\n'
364                )
365            # elasticities
366            of.write('#### Elasticities ####\n')
367            for step in (
368                getattr(mymap, fixed).isReagentOf()
369                + getattr(mymap, fixed).isModifierOf()
370            ):
371                of.writelines(
372                    '\\ec{'
373                    + step
374                    + '}{'
375                    + fixed
376                    + '} = '
377                    + str(basemod.ec.get(step, fixed))
378                    + '\n'
379                )
380            of.write('#### Partial Response coefficients ####\n')
381            for reaction in getattr(mymap, fixed).isReagentOf():
382                for step in (
383                    getattr(mymap, fixed).isReagentOf()
384                    + getattr(mymap, fixed).isModifierOf()
385                ):
386                    of.writelines(
387                        '\\prc{'
388                        + step
389                        + '}{'
390                        + reaction
391                        + '}{'
392                        + fixed
393                        + '} = '
394                        + str(basemod.ec.get(step, fixed) * mod.cc.get(reaction, step))
395                        + '='
396                        + str(basemod.ec.get(step, fixed))
397                        + '*'
398                        + str(mod.cc.get(reaction, step))
399                        + '\n'
400                    )
401            # steady-state values
402            of.write('#### Steady-state values ####\n')
403            for step in getattr(mymap, fixed).isReagentOf():
404                of.writelines(
405                    '\\J{' + step + '} = ' + str(getattr(basemod, 'J_' + step)) + '\n'
406                )
407            of.writelines(
408                '[' + fixed + '] = ' + str(getattr(basemod, fixed + '_ss')) + '\n'
409            )
410            of.write('\n')
411            if flag:
412                of.close()
413                del of
414                flag = 0
415
416    def doAllRateChar(self, summary=1, plot=0, data=0, solver=0):
417        """
418        doAllRateChar(summary=1, plot=0, data=0, solver=0)
419
420        Does rate characteristic analysis for all variable species of a model
421        by calling doRateChar() for each one.
422
423        Arguments:
424        ==========
425        summary (bool) - write summary to output file (default 1)
426        plot (bool)    - save plots of all rate characteristics to PNG files (default 0)
427        data (bool)    - save data for all graphs to tab-delimited file (default 0)
428        solver (int)   - steady-state solver to use: 0 - hybrd, 1 - fintslv, 2 - NLEQ2
429        """
430        if summary:
431            self.summaryfile = open(
432                os.path.join(self._psc_out_dir, self._mod + '_summary.txt'), 'w'
433            )
434            self.summaryfile.write('# Base model: ' + self._psc_file + '\n')
435            self.summaryfile.write('####################################\n\n')
436        basemod = self._psc_mod
437        for fixed in basemod.species:
438            self.doRateChar(fixed, summary, solver)
439            if plot:
440                fig = eval('self.' + fixed + '.figure()')
441                fig.plotElasRC(file=1)
442                fig.plotPartialRCDemand(file=1)
443                fig.plotPartialRCSupply(file=1)
444
445            if data:
446                eval('self.' + fixed + '.writeDataToFile()')
447        if summary:
448            self.summaryfile.close()
449            del self.summaryfile
450
451    def getRateCharData(self, spec):
452        rcd = getattr(self, spec)
453        return rcd
454
455
456class RateCharData:
457    """Class to store rate characteristic output data """
458
459    def __init__(
460        self, parent, mod, fix, useroutput, scandata, FixProductFor, FixSubstrateFor
461    ):
462        self.parent = parent
463        self.modelfilename = self.parent._psc_file
464        self.mod = mod
465        self.fix = fix
466        self.useroutput = useroutput
467        self.scandata = scandata
468        self.FixProductFor = FixProductFor
469        self.FixSubstrateFor = FixSubstrateFor
470        self.basemod = self.parent._psc_mod
471        self.mymap = ModelMap(self.basemod)
472        self.slope_range_factor = 3.0
473        self.scan_min = (
474            getattr(self.basemod, fix + '_ss') / self.parent.min_concrange_factor
475        )
476        self.scan_max = (
477            getattr(self.basemod, fix + '_ss') * self.parent.max_concrange_factor
478        )
479        self.flux_min = (
480            abs(self.scandata[:, 1:]).min() / self.parent.min_fluxrange_factor
481        )
482        self.flux_max = self.scandata[:, 1:].max() * self.parent.max_fluxrange_factor
483        self.scan_points = len(self.scandata)
484        self.isReagentOf = getattr(self.mymap, self.fix).isReagentOf()
485        self.isModifierOf = getattr(self.mymap, self.fix).isModifierOf()
486        self.isProductOf = getattr(self.mymap, self.fix).isProductOf()
487        self.isSubstrateOf = getattr(self.mymap, self.fix).isSubstrateOf()
488
489    def _wrn(self):
490        warnings.warn(
491            'This is a legacy function, use figure and it\'s  functions instead',
492            stacklevel=2,
493        )
494
495    def plotRateChar(self, file=0):
496        """
497        TAKE NOTE:
498        This function is only included for legacy purposes. All new scripts
499        should use the figure function to create an interactive RCfigure object
500        and use the functions of that object to draw and manipulate plots.
501
502        plotRateChar(file=0)
503
504        Plot the data in the instantiated RateCharData class (fluxes vs.
505        the concentration of the fixed species).
506
507        Arguments:
508        ==========
509        file (bool) - save plot as PNG file (default=0, plots to screen)
510        """
511        self._wrn()
512        fig = self.figure()
513        fig.plotRateChar(file)
514
515    def plotRateCharIndiv(self, file=0):
516        """
517        TAKE NOTE:
518        This function is only included for legacy purposes. All new scripts
519        should use the figure function to create an interactive RCfigure object
520        and use the functions of that object to draw and manipulate plots.
521
522        plotRateCharIndiv(file=0)
523
524        Plot the data in the instantiated RateCharData class (fluxes vs.
525        the concentration of the fixed species for each supply and demand flux
526        in turn).
527
528        Arguments:
529        ==========
530        file (bool) - save plot as PNG file (default=0, plots to screen)
531        """
532        self._wrn()
533        fig = self.figure()
534        fig.plotRateCharIndiv(file)
535
536    def plotElasRC(self, file=0):
537        """
538        TAKE NOTE:
539        This function is only included for legacy purposes. All new scripts
540        should use the figure function to create an interactive RCfigure object
541        and use the functions of that object to draw and manipulate plots.
542
543        plotElasRC(file=0)
544
545        Plot the data in the instantiated RateCharData class (rates,
546        elasticities and response coefficients vs. the concentration of the
547        fixed species).
548
549        Arguments:
550        ==========
551        file (bool) - save plot as PNG file (default=0, plots to screen)
552        """
553        self._wrn()
554        fig = self.figure()
555        fig.plotElasRC(file)
556
557    def plotElasRCIndiv(self, file=0):
558        """
559        TAKE NOTE:
560        This function is only included for legacy purposes. All new scripts
561        should use the figure function to create an interactive RCfigure object
562        and use the functions of that object to draw and manipulate plots.
563
564        plotElasIndivRC(file=0)
565
566        Plot the data in the instantiated RateCharData class (rates,
567        elasticities [sans modifier elasticities] and response coefficients vs.
568        the concentration of the fixed species for each supply and demand
569        reaction in turn).
570
571        Arguments:
572        ==========
573        file (bool) - save plot as PNG file (default=0, plots to screen)
574        """
575        self._wrn()
576        fig = self.figure()
577        fig.plotElasRCIndiv(file)
578
579    def plotPartialRCSupply(self, file=0):
580        """
581        TAKE NOTE:
582        This function is only included for legacy purposes. All new scripts
583        should use the figure function to create an interactive RCfigure object
584        and use the functions of that object to draw and manipulate plots.
585
586        plotPartialRC(file=0)
587
588        Plot the data in the instantiated RateCharData class for supply
589        reactions (fluxes and partial response coefficients vs. the
590        concentration of the fixed species).
591
592        Arguments:
593        ==========
594        file (bool) - save plot as PNG file (default=0, plots to screen)
595        """
596        self._wrn()
597        fig = self.figure()
598        fig.plotPartialRCSupply(file)
599
600    def plotPartialRCDemand(self, file=0):
601        """
602        TAKE NOTE:
603        This function is only included for legacy purposes. All new scripts
604        should use the figure function to create an interactive RCfigure object
605        and use the functions of that object to draw and manipulate plots.
606
607        plotPartialRC(file=0)
608
609        Plot the data in the instantiated RateCharData class for demand
610        reactions (fluxes and partial response coefficients vs. the
611        concentration of the fixed species).
612
613        Arguments:
614        ==========
615        file (bool) - save plot as PNG file (default=0, plots to screen)
616        """
617        self._wrn()
618        fig = self.figure()
619        fig.plotPartialRCDemand(file)
620
621    def plotPartialRCIndiv(self, file=0):
622        """
623        TAKE NOTE:
624        This function is only included for legacy purposes. All new scripts
625        should use the figure function to create an interactive RCfigure object
626        and use the functions of that object to draw and manipulate plots.
627
628        plotPartialRC(file=0)
629
630        Plot the data in the instantiated RateCharData class for each reaction
631        in turn (fluxes and partial response coefficients vs. the
632        concentration of the fixed species).
633
634        Arguments:
635        ==========
636        file (bool) - save plot as PNG file (default=0, plots to screen)
637        """
638        self._wrn()
639        fig = self.figure()
640        fig.plotPartialRCIndiv(file)
641
642    def findRegulatory(self, regTolerance=3):
643
644        """
645        Returns a list of reactions for which the current fixed metabolite is
646        a regulatory metabolite eg. Elasticity = Response Coefficient. This is
647        measured as an angle between the elasticity coefficient and the response
648        coefficient
649
650        Arguments:
651        ==========
652        regTolerance (double) - the tolerance in degrees
653        """
654        regulatoryFor = []
655        for step in getattr(self.mymap, self.fix).isReagentOf():
656            ec_slope = self.basemod.ec.get(step, self.fix)
657            rc_slope = self.mod.rc.get(step, self.fix)
658
659            angleBetween = scipy.Infinity
660            a = scipy.rad2deg(scipy.arctan(ec_slope))
661            b = scipy.rad2deg(scipy.arctan(rc_slope))
662
663            if ec_slope > 0 and rc_slope > 0:
664                angleBetween = abs(a - b)
665            elif ec_slope < 0 and rc_slope < 0:
666                angleBetween = abs(abs(a) - abs(b))
667            elif ec_slope < 0 or rc_slope < 0:
668                angleBetween = abs(a) + abs(b)
669            # print angleBetween
670            if angleBetween < regTolerance:
671                regulatoryFor.append(step)
672        return regulatoryFor
673
674    def writeDataToFile(self):
675
676        """
677        writeDataToFile()
678
679        Write the data in the instantiated RateCharData class (concentration of fixed
680        species and fluxes of all adjoining reactions) to a tab-delimited CSV file.
681
682        Arguments:
683        ==========
684        None
685        """
686        filename = os.path.join(
687            self.parent._psc_out_dir,
688            self.modelfilename[:-4] + '_' + self.fix + '_data.txt',
689        )
690        of = open(filename, 'w')
691        of.write('# Base model: ' + self.modelfilename + '\n')
692        of.write('# Rate Characteristic scan around species: ' + self.fix + '\n')
693        of.write('##########################################\n#')
694        for i in self.useroutput:
695            of.write(i + '\t')
696        of.write('\n')
697        numpy.savetxt(of, self.scandata, fmt='%05.2f', delimiter='\t')
698        # scipy.io.write_array(of, self.scandata, separator='\t', keep_open=0)
699
700    def figure(self, **kwds):
701        """
702        Returns an object of the RcFigure class with the current object as it's
703        parent. This object has a variety of functions to plot the rate
704        characteristic data interactively.
705
706        """
707
708        return RcFigure(self, **kwds)
709
710
711class RcFigure:
712    from pylab import figure, close, show, ioff, figtext
713    from matplotlib.transforms import offset_copy
714    from matplotlib.figure import Figure as mplfigure
715
716    def __init__(
717        self, parent, size=(12.80, 10.24), outformat='png', legendloc=3, lims=0
718    ):
719        ioff()
720        self.legendloc = legendloc
721        self.outformat = outformat
722        self.parent = parent
723        self.figure = figure(num=1, figsize=size)
724        self.figure.clf()
725        self.ax = self.figure.add_subplot(111)
726
727        fix_ss = getattr(self.parent.basemod, self.parent.fix + '_ss')
728        self.ss_line = self.ax.axvline(x=fix_ss, ls=':', color='0.5', lw=1)
729        self.ax.loglog()
730        self.legend = self.ax.legend(loc=self.legendloc)
731        self.title_label = (
732            'Model: ' + self.parent.parent._mod + ' - Rate Characteristic'
733        )
734        self.title = self.ax.set_title(self.title_label)
735        if lims == 0:
736            if numpy.isnan(self.parent.flux_min):
737                self.ax.set_xlim(self.parent.scan_min, self.parent.scan_max)
738            elif numpy.isnan(self.parent.flux_max):
739                self.ax.set_xlim(self.parent.scan_min, self.parent.scan_max)
740            else:
741                self.ax.axis(
742                    [
743                        self.parent.scan_min,
744                        self.parent.scan_max,
745                        self.parent.flux_min,
746                        self.parent.flux_max,
747                    ]
748                )
749        else:
750            self.ax.axis(lims)
751
752        self.lims = self.ax.get_xlim() + self.ax.get_ylim()
753
754        self.huelist = [
755            0.0,  # red
756            0.3333333,  # green
757            0.6666666,  # blue
758            0.0833333,  # orange
759            0.7499999,  # purple
760            0.4166666,  # aqua
761            0.8333333,  # pink
762            0.1666666,  # yellow
763            0.5833333,  # light blue
764            0.2499999,  # lime green
765            0.9166666,
766        ]  # magenta
767
768        self.huedict = self._assign_hue()
769        self.reactions = self._reactions()
770        self.supply_flux, self.demand_flux = self._plot_flux()
771        self.reaction_flux = dict(
772            list(self.supply_flux.items()) + list(self.demand_flux.items())
773        )
774
775        self.supply_total, self.demand_total = self._plot_totals()
776
777        self.supply_elas, self.demand_elas, self.mod_elas = self._plot_elas()
778        self.reaction_elas = dict(
779            list(self.supply_elas.items())
780            + list(self.demand_elas.items())
781            + list(self.mod_elas.items())
782        )
783
784        self.supply_partial_rc, self.demand_partial_rc = self._plot_partial_rc()
785        self.reaction_partial_rc = dict(
786            list(self.supply_partial_rc.items()) + list(self.demand_partial_rc.items())
787        )
788
789        self.supply_rc, self.demand_rc = self._plot_rc()
790        self.reaction_rc = dict(
791            list(self.supply_rc.items()) + list(self.demand_rc.items())
792        )
793
794        self.ax.set_xlabel('$' + self.parent.fix + '$')
795        self.ax.set_ylabel('$Rate$')
796
797        self.legend_status = True
798        self.title_status = True
799        self.clear()
800        self._make_species_dir()
801
802    def _make_species_dir(self):
803        """
804        _make_species_dir()
805
806        Creates a directory for the specific species used in the figure.
807
808        Arguments:
809        =========
810        None
811        """
812
813        out_base_dir_status = os.path.exists(
814            os.path.join(self.parent.parent._psc_out_dir, self.parent.fix)
815        )
816        if not out_base_dir_status:
817            os.mkdir(os.path.join(self.parent.parent._psc_out_dir, self.parent.fix))
818
819    def _reactions(self):
820        reactions = {}
821        for r in self.parent.isSubstrateOf:
822            reactions['J_' + r] = 'Demand Reaction'
823        for r in self.parent.isProductOf:
824            reactions['J_' + r] = 'Supply Reaction'
825        for r in self.parent.isModifierOf:
826            reactions['J_' + r] = 'Modified Reaction'
827        return reactions
828
829    def _assign_hue(self):
830        allreactions = self.parent.isReagentOf + self.parent.isModifierOf
831        h_num = 0
832        while len(allreactions) > len(self.huelist):
833            self.huelist.append(self.huelist[h_num])
834            h_num = h_num + 1
835        huedict = {}
836        for r_num, reaction in enumerate(allreactions):
837            huedict['J_' + reaction] = self.huelist[r_num]
838        return huedict
839
840    def _plot_totals(self):
841
842        supply_total = None
843        demand_total = None
844
845        if len(self.parent.FixProductFor) > 1:
846            res_sup = numpy.zeros(self.parent.scan_points)
847            for flux in range(len(self.parent.FixProductFor)):
848                res_sup += self.parent.scandata[:, flux + 1]
849            supply_total = self.ax.plot(
850                self.parent.scandata[:, 0],
851                res_sup,
852                '--',
853                color='#0000C0',
854                label='$Total$\n$Supply$',
855                linewidth=1,
856            )
857
858        if len(self.parent.FixSubstrateFor) > 1:
859            res_dem = numpy.zeros(self.parent.scan_points)
860            for flux in range(len(self.parent.FixSubstrateFor)):
861                res_dem += self.parent.scandata[
862                    :, flux + len(self.parent.FixProductFor) + 1
863                ]
864            demand_total = self.ax.plot(
865                self.parent.scandata[:, 0],
866                res_dem,
867                '--',
868                color='#00c000',
869                label='$Total$\n$Demand$',
870                linewidth=1,
871            )
872        return supply_total, demand_total
873
874    def _plot_flux(self):
875        supply_flux = {}
876        demand_flux = {}
877
878        # supply plots
879        for column, reaction in enumerate(self.parent.FixProductFor):
880            colr = htor(self.huedict[reaction], 1.0, 1.0)
881            lbl = '$J_{' + reaction[2:] + '}$'
882            supply_flux[reaction] = self.ax.plot(
883                self.parent.scandata[:, 0],
884                self.parent.scandata[:, column + 1],
885                '-',
886                linewidth=1,
887                label=lbl,
888                color=colr,
889            )
890            supply_flux[reaction][0]._permalabel = lbl
891        # demand plots
892        for column, reaction in enumerate(self.parent.FixSubstrateFor):
893            colr = htor(self.huedict[reaction], 1.0, 1.0)
894            lbl = '$J_{' + reaction[2:] + '}$'
895            demand_flux[reaction] = self.ax.plot(
896                self.parent.scandata[:, 0],
897                self.parent.scandata[:, column + len(self.parent.FixProductFor) + 1],
898                '-',
899                linewidth=1,
900                label=lbl,
901                color=colr,
902            )
903            demand_flux[reaction][0]._permalabel = lbl
904
905        return supply_flux, demand_flux
906
907    def _plot_elas(self):
908        supply_elas = {}
909        demand_elas = {}
910        mod_elas = {}
911
912        allreactions = self.parent.isReagentOf + self.parent.isModifierOf
913
914        for step in allreactions:
915            J_ss = getattr(self.parent.basemod, 'J_' + step)
916            fix_ss = getattr(self.parent.basemod, self.parent.fix + '_ss')
917            ec_slope = self.parent.basemod.ec.get(step, self.parent.fix)
918            ec_constant = J_ss / (fix_ss ** ec_slope)
919
920            ydist = scipy.log10(self.lims[3] / self.lims[2])
921            xdist = scipy.log10(self.lims[1] / self.lims[0])
922            xyscale = xdist / ydist
923
924            scale_factor = scipy.cos(scipy.arctan(ec_slope * xyscale))
925            range_min_old = fix_ss / self.parent.slope_range_factor
926            range_max_old = fix_ss * self.parent.slope_range_factor
927            distance = scipy.log10(fix_ss / range_min_old) * scale_factor
928            range_min = fix_ss / (10 ** distance)
929            range_max = fix_ss * (10 ** distance)
930
931            fix_conc = numpy.linspace(range_min, range_max, num=5)
932
933            ec_rate = ec_constant * fix_conc ** (ec_slope)
934
935            colr = htor(self.huedict['J_' + step], 0.4, 1.0)
936            if step in self.parent.isProductOf:
937                the_dict = supply_elas
938            elif step in self.parent.isSubstrateOf:
939                the_dict = demand_elas
940            elif step in self.parent.isModifierOf:
941                the_dict = mod_elas
942
943            lbl = '$\epsilon^{' + step + '}_{' + self.parent.fix + '}$'
944            the_dict['J_' + step] = self.ax.plot(
945                fix_conc, ec_rate, color=colr, linewidth=2, label=lbl
946            )
947            the_dict['J_' + step][0]._permalabel = lbl
948
949        return supply_elas, demand_elas, mod_elas
950
951    def _plot_rc(self):
952        supply_rc = {}
953        demand_rc = {}
954
955        for step in self.parent.isReagentOf:
956            J_ss = getattr(self.parent.basemod, 'J_' + step)
957            fix_ss = getattr(self.parent.basemod, self.parent.fix + '_ss')
958            rc_slope = self.parent.mod.rc.get(step, self.parent.fix)
959            rc_constant = J_ss / (fix_ss ** rc_slope)
960
961            ydist = scipy.log10(self.lims[3] / self.lims[2])
962            xdist = scipy.log10(self.lims[1] / self.lims[0])
963            xyscale = xdist / ydist
964
965            scale_factor = scipy.cos(scipy.arctan(rc_slope * xyscale))
966            range_min_old = fix_ss / self.parent.slope_range_factor
967            range_max_old = fix_ss * self.parent.slope_range_factor
968            distance = scipy.log10(fix_ss / range_min_old) * scale_factor
969            range_min = fix_ss / (10 ** distance)
970            range_max = fix_ss * (10 ** distance)
971
972            fix_conc = numpy.linspace(range_min, range_max, num=5)
973            rc_rate = rc_constant * fix_conc ** (rc_slope)
974
975            colr = htor(self.huedict['J_' + step], 1.0, 0.6)
976
977            if step in self.parent.isProductOf:
978                the_dict = supply_rc
979            elif step in self.parent.isSubstrateOf:
980                the_dict = demand_rc
981            lbl = '$R^{J_{' + step + '}}_{' + self.parent.fix + '}$'
982            the_dict['J_' + step] = self.ax.plot(
983                fix_conc, rc_rate, '--', lw=0.8, color=colr, label=lbl
984            )
985            the_dict['J_' + step][0]._permalabel = lbl
986
987        return supply_rc, demand_rc
988
989    def _plot_partial_rc(self):
990
991        supply_partial_rc = {}
992        demand_partial_rc = {}
993        allreactions = self.parent.isReagentOf + self.parent.isModifierOf
994
995        for reaction in self.parent.isReagentOf:
996            J_ss = getattr(self.parent.basemod, 'J_' + reaction)
997            fix_ss = getattr(self.parent.basemod, self.parent.fix + '_ss')
998            step_dict = {}
999            for step in allreactions:
1000                prc_slope = self.parent.basemod.ec.get(
1001                    step, self.parent.fix
1002                ) * self.parent.mod.cc.get(reaction, step)
1003                prc_constant = J_ss / (fix_ss ** prc_slope)
1004
1005                ydist = scipy.log10(self.lims[3] / self.lims[2])
1006                xdist = scipy.log10(self.lims[1] / self.lims[0])
1007                xyscale = xdist / ydist
1008
1009                scale_factor = scipy.cos(scipy.arctan(prc_slope * xyscale))
1010                range_min_old = fix_ss / self.parent.slope_range_factor
1011                range_max_old = fix_ss * self.parent.slope_range_factor
1012                distance = scipy.log10(fix_ss / range_min_old) * scale_factor
1013                range_min = fix_ss / (10 ** distance)
1014                range_max = fix_ss * (10 ** distance)
1015
1016                fix_conc = numpy.linspace(range_min, range_max, num=5)
1017                prc_rate = prc_constant * fix_conc ** (prc_slope)
1018                # if abs(prc_slope) > 0.01:    # only plot nonzero partial RCs
1019                colr = htor(self.huedict['J_' + step], 0.5, 1.0)
1020                lbl = (
1021                    '$^{'
1022                    + step
1023                    + '}R^{J_{'
1024                    + reaction
1025                    + '}}_{'
1026                    + self.parent.fix
1027                    + '}$'
1028                )
1029                step_dict['J_' + step] = self.ax.plot(
1030                    fix_conc, prc_rate, '-', color=colr, lw=1, label=lbl
1031                )
1032                step_dict['J_' + step][0]._permalabel = lbl
1033                if abs(prc_slope) > 0.01:
1034                    step_dict['J_' + step][0]._nonzero = True
1035                else:
1036                    step_dict['J_' + step][0]._nonzero = False
1037                if reaction in self.parent.isProductOf:
1038                    the_dict = supply_partial_rc
1039                elif reaction in self.parent.isSubstrateOf:
1040                    the_dict = demand_partial_rc
1041                the_dict['J_' + reaction] = step_dict
1042
1043        return supply_partial_rc, demand_partial_rc
1044
1045    def show(self):
1046        """
1047        Shows the plot and updates the legend and title. Calls the methods
1048        update_title and update_legend
1049        """
1050        self.update_title()
1051        self.update_legend()
1052        self.figure.show()
1053
1054    def update_legend(self):
1055        """
1056        Updates the title using the value of legend_status (boolean) for
1057        visibility.
1058
1059        This function is called automatically by save() and show().
1060        """
1061
1062        del self.legend
1063        self.legend = self.ax.legend(loc=self.legendloc)
1064        self.legend.set_visible(self.legend_status)
1065
1066    def update_title(self):
1067        """
1068        Updates the title using the value of title_label for the title label and
1069        title_status (boolean) for visibility.
1070
1071        This function is called automatically by save() and show().
1072        """
1073
1074        del self.title
1075        self.title = self.ax.set_title(self.title_label)
1076        self.title.set_visible(self.title_status)
1077
1078    def save(self, filename):
1079        """
1080        Saves the figure with the given file name (in the current directory when
1081        no location is specified)
1082
1083        Arguments:
1084        ==========
1085        filename    -    The name and location of the file
1086        """
1087
1088        self.update_title()
1089        self.update_legend()
1090        self.figure.savefig(filename + '.' + self.outformat)
1091
1092    def _line_label_vis(self, line, visibility):
1093        if visibility is True:
1094            line._label = line._permalabel
1095        else:
1096            line._label = ''
1097
1098    def set_visible(self, name, affzero=0, **kwds):
1099        """
1100        This function sets the visibility of the lines associated with a certain
1101        reaction. Use the show() method to show the updated figure after this
1102        method.
1103
1104        Arguments:
1105        ==========
1106        name     -    The name of the reaction. This value is in the format
1107                      "J_reactionname", where reactionname is the name of the
1108                      reaction as it appears in the psc file.
1109        affzero  -    If true, partial response coefficient lines with zero
1110                      value slopes will also be affected (default = 0)
1111
1112        Keywords:
1113        ==========
1114        Valid keywords are:
1115        p        -    partial response coefficient
1116        f        -    flux
1117        e        -    elasticity
1118        r        -    response coefficient
1119
1120        For reactions where the fixed metabolite is only a modifier, the only
1121        valid keyword is "e".
1122
1123        Usage example:
1124        fig.set_visible('J_reaction1', l= True, r = False)
1125
1126        This will set the flux line for "reaction1" visible while turning off the
1127        response coefficient line. Elasticity and partial response coefficient
1128        lines are left as they were.
1129
1130        or
1131
1132        fig.set_visible('J_reaction1', f= True, r = False, p = True, showzero = True)
1133
1134        Same as above, but with partial response coefficients also visible, including
1135        those with slopes equal to zero
1136
1137
1138
1139        """
1140        kwrs = {'e': self.reaction_elas, 'r': self.reaction_rc, 'f': self.reaction_flux}
1141        for key in kwrs.keys():
1142            if key in kwds:
1143                line = kwrs[key][name][0]
1144                line.set_visible(kwds[key])
1145                self._line_label_vis(line, kwds[key])
1146        if 'p' in kwds and name in list(self.reaction_partial_rc.keys()):
1147            for key in self.reaction_partial_rc[name].keys():
1148                line = self.reaction_partial_rc[name][key][0]
1149                if affzero:
1150                    line.set_visible(kwds['p'])
1151                    self._line_label_vis(line, kwds['p'])
1152                elif line._nonzero:
1153                    line.set_visible(kwds['p'])
1154                    self._line_label_vis(line, kwds['p'])
1155
1156    def set_type_visible(self, line_set_name, visibility, affzero=0):
1157        """
1158        Sets the visibility of a set of lines. Use the show() method to show the
1159        updated figure after this method.
1160
1161        Arguments:
1162        ==========
1163        line_set_name - The name of the set of lines to enable or disable
1164                        Valid values are:
1165                            supply_rc
1166                            demand_rc
1167                            supply_elas
1168                            demand_elas
1169                            mod_elas
1170                            supply_flux
1171                            demand_flux
1172                            supply_partial_rc
1173                            demand_partial_rc
1174                            supply_total
1175                            demand_total
1176        visibility    - True or False
1177        affzero       - Boolean (default = 0) if set to True, zero slope
1178                        partial response coefficients will also be
1179                        affected
1180        """
1181        values = {
1182            'supply_rc': self.supply_rc,
1183            'demand_rc': self.demand_rc,
1184            'supply_elas': self.supply_elas,
1185            'demand_elas': self.demand_elas,
1186            'mod_elas': self.mod_elas,
1187            'supply_flux': self.supply_flux,
1188            'demand_flux': self.demand_flux,
1189            'supply_partial_rc': self.supply_partial_rc,
1190            'demand_partial_rc': self.demand_partial_rc,
1191            'supply_total': self.supply_total,
1192            'demand_total': self.demand_total,
1193        }
1194        if 'partial' in line_set_name:
1195            for dic in values[line_set_name].keys():
1196                for key in values[line_set_name][dic].keys():
1197                    line = values[line_set_name][dic][key][0]
1198                    if affzero:
1199                        line.set_visible(visibility)
1200                        self._line_label_vis(line, visibility)
1201                    elif line._nonzero:
1202                        line.set_visible(visibility)
1203                        self._line_label_vis(line, visibility)
1204
1205        elif 'total' in line_set_name:
1206            if values[line_set_name]:
1207                values[line_set_name][0].set_visible(visibility)
1208
1209        else:
1210            for key in values[line_set_name].keys():
1211                line = values[line_set_name][key][0]
1212                line.set_visible(visibility)
1213                self._line_label_vis(line, visibility)
1214
1215    def plotRateChar(self, file=0):
1216        """
1217        plotRateChar(file=0)
1218
1219        Plot the data in the instantiated RateCharData class (fluxes vs.
1220        the concentration of the fixed species).
1221
1222        Arguments:
1223        ==========
1224        file (bool) - save plot as PNG file (default=0, plots to screen)
1225        """
1226
1227        filename = os.path.join(
1228            self.parent.parent._psc_out_dir,
1229            self.parent.fix,
1230            self.parent.modelfilename[:-4] + '_' + self.parent.fix + '_flux',
1231        )
1232
1233        self.title_label = (
1234            'Model: ' + self.parent.parent._mod + ' - Rate Characteristic'
1235        )
1236        if self.supply_total:
1237            self.supply_total[0].set_visible(True)
1238        if self.demand_total:
1239            self.demand_total[0].set_visible(True)
1240        self.set_type_visible('supply_flux', True)
1241        self.set_type_visible('demand_flux', True)
1242        self.set_type_visible('supply_elas', False)
1243        self.set_type_visible('demand_elas', False)
1244        self.set_type_visible('mod_elas', False)
1245        self.set_type_visible('supply_rc', False)
1246        self.set_type_visible('demand_rc', False)
1247        self.set_type_visible('supply_partial_rc', False, affzero=True)
1248        self.set_type_visible('demand_partial_rc', False, affzero=True)
1249
1250        if file:
1251            self.save(filename)
1252        else:
1253            self.show()
1254
1255    def clear(self):
1256        if self.supply_total:
1257            self.supply_total[0].set_visible(True)
1258        if self.demand_total:
1259            self.demand_total[0].set_visible(True)
1260
1261        self.set_type_visible('supply_flux', False)
1262        self.set_type_visible('demand_flux', False)
1263        self.set_type_visible('supply_elas', False)
1264        self.set_type_visible('demand_elas', False)
1265        self.set_type_visible('mod_elas', False)
1266        self.set_type_visible('supply_rc', False)
1267        self.set_type_visible('demand_rc', False)
1268        self.set_type_visible('supply_partial_rc', False, affzero=True)
1269        self.set_type_visible('demand_partial_rc', False, affzero=True)
1270
1271    def plotElasRC(self, file=0):
1272        """
1273        plotElasRC(file=0)
1274
1275        Plot the data in the instantiated RateCharData class (rates,
1276        elasticities and response coefficients vs. the concentration of the
1277        fixed species).
1278
1279        Arguments:
1280        ==========
1281        file (bool) - save plot as PNG file (default=0, plots to screen)
1282        """
1283
1284        filename = os.path.join(
1285            self.parent.parent._psc_out_dir,
1286            self.parent.fix,
1287            self.parent.modelfilename[:-4] + '_' + self.parent.fix + '_ElasRC',
1288        )
1289
1290        self.title_label = (
1291            'Model: '
1292            + self.parent.parent._mod
1293            + ' - Response vs. elasticity coefficients'
1294        )
1295        if self.supply_total:
1296            self.supply_total[0].set_visible(True)
1297        if self.demand_total:
1298            self.demand_total[0].set_visible(True)
1299        self.set_type_visible('supply_flux', True)
1300        self.set_type_visible('demand_flux', True)
1301        self.set_type_visible('supply_elas', True)
1302        self.set_type_visible('demand_elas', True)
1303        self.set_type_visible('mod_elas', True)
1304        self.set_type_visible('supply_rc', True)
1305        self.set_type_visible('demand_rc', True)
1306        self.set_type_visible('supply_partial_rc', False, affzero=True)
1307        self.set_type_visible('demand_partial_rc', False, affzero=True)
1308
1309        if file:
1310            self.save(filename)
1311        else:
1312            self.show()
1313
1314    def plotPartialRCSupply(self, file=0):
1315        """
1316        plotPartialRC(file=0)
1317
1318        Plot the data in the instantiated RateCharData class for supply
1319        reactions (fluxes and partial response coefficients vs. the
1320        concentration of the fixed species).
1321
1322        Arguments:
1323        ==========
1324        file (bool) - save plot as PNG file (default=0, plots to screen)
1325        """
1326
1327        filename = os.path.join(
1328            self.parent.parent._psc_out_dir,
1329            self.parent.fix,
1330            self.parent.modelfilename[:-4] + '_' + self.parent.fix + '_partialRCSupply',
1331        )
1332
1333        self.title_label = (
1334            'Model: '
1335            + self.parent.parent._mod
1336            + ' - Partial response coefficients for Supply'
1337        )
1338        if self.supply_total:
1339            self.supply_total[0].set_visible(True)
1340        if self.demand_total:
1341            self.demand_total[0].set_visible(True)
1342        self.set_type_visible('supply_flux', True)
1343        self.set_type_visible('demand_flux', False)
1344        self.set_type_visible('supply_elas', False)
1345        self.set_type_visible('demand_elas', False)
1346        self.set_type_visible('mod_elas', False)
1347        self.set_type_visible('supply_rc', True)
1348        self.set_type_visible('demand_rc', False)
1349        self.set_type_visible('supply_partial_rc', True)
1350        self.set_type_visible('demand_partial_rc', False, affzero=True)
1351
1352        if file:
1353            self.save(filename)
1354        else:
1355            self.show()
1356
1357    def plotPartialRCDemand(self, file=0):
1358        """
1359        plotPartialRC(file=0)
1360
1361        Plot the data in the instantiated RateCharData class for demand
1362        reactions (fluxes and partial response coefficients vs. the
1363        concentration of the fixed species).
1364
1365        Arguments:
1366        ==========
1367        file (bool) - save plot as PNG file (default=0, plots to screen)
1368        """
1369
1370        filename = os.path.join(
1371            self.parent.parent._psc_out_dir,
1372            self.parent.fix,
1373            self.parent.modelfilename[:-4] + '_' + self.parent.fix + '_partialRCDemand',
1374        )
1375
1376        self.title_label = (
1377            'Model: '
1378            + self.parent.parent._mod
1379            + ' - Partial response coefficients for Demand'
1380        )
1381        if self.supply_total:
1382            self.supply_total[0].set_visible(True)
1383        if self.demand_total:
1384            self.demand_total[0].set_visible(True)
1385        self.set_type_visible('supply_flux', False)
1386        self.set_type_visible('demand_flux', True)
1387        self.set_type_visible('supply_elas', False)
1388        self.set_type_visible('demand_elas', False)
1389        self.set_type_visible('mod_elas', False)
1390        self.set_type_visible('supply_rc', False)
1391        self.set_type_visible('demand_rc', True)
1392        self.set_type_visible('supply_partial_rc', False, affzero=True)
1393        self.set_type_visible('demand_partial_rc', True)
1394
1395        if file:
1396            self.save(filename)
1397        else:
1398            self.show()
1399
1400    def plotRateCharIndiv(self, file=0):
1401        """
1402        plotRateCharIndiv(file=0)
1403
1404        Plot the data in the instantiated RateCharData class (fluxes vs.
1405        the concentration of the fixed species for each supply and demand flux
1406        in turn).
1407
1408        Arguments:
1409        ==========
1410        file (bool) - save plot as PNG file (default=0, plots to screen)
1411        """
1412
1413        if self.supply_total:
1414            self.supply_total[0].set_visible(True)
1415        if self.demand_total:
1416            self.demand_total[0].set_visible(True)
1417        self.set_type_visible('supply_flux', False)
1418        self.set_type_visible('demand_flux', False)
1419        self.set_type_visible('supply_elas', False)
1420        self.set_type_visible('demand_elas', False)
1421        self.set_type_visible('mod_elas', False)
1422        self.set_type_visible('supply_rc', False)
1423        self.set_type_visible('demand_rc', False)
1424        self.set_type_visible('supply_partial_rc', False, affzero=True)
1425        self.set_type_visible('demand_partial_rc', False, affzero=True)
1426
1427        for line in self.huedict.keys():
1428            if line in self.mod_elas:
1429                continue
1430            if line not in list(self.reaction_partial_rc.keys()):
1431                continue
1432            if len(self.supply_flux) == 1:
1433                self.set_visible(list(self.supply_flux.keys())[0], f=True)
1434            if len(self.demand_flux) == 1:
1435                self.set_visible(list(self.demand_flux.keys())[0], f=True)
1436
1437            filename = os.path.join(
1438                self.parent.parent._psc_out_dir,
1439                self.parent.fix,
1440                self.parent.modelfilename[:-4]
1441                + '_'
1442                + self.parent.fix
1443                + '_flux_'
1444                + line
1445                + '_',
1446            )
1447
1448            self.title_label = (
1449                'Model: ' + self.parent.parent._mod + ' - Rate Characteristic - ' + line
1450            )
1451            self.set_visible(line, f=True)
1452            if file:
1453                self.save(filename)
1454            else:
1455                self.show()
1456                input("Press Enter to continue...")
1457            self.set_visible(line, f=False)
1458
1459    def plotElasRCIndiv(self, file=0):
1460        """
1461        plotElasIndivRC(file=0)
1462
1463        Plot the data in the instantiated RateCharData class (rates,
1464        elasticities [sans modifier elasticities] and response coefficients vs.
1465        the concentration of the fixed species for each supply and demand
1466        reaction in turn).
1467
1468        Arguments:
1469        ==========
1470        file (bool) - save plot as PNG file (default=0, plots to screen)
1471        """
1472
1473        if self.supply_total:
1474            self.supply_total[0].set_visible(True)
1475        if self.demand_total:
1476            self.demand_total[0].set_visible(True)
1477        self.set_type_visible('supply_flux', False)
1478        self.set_type_visible('demand_flux', False)
1479        self.set_type_visible('supply_elas', False)
1480        self.set_type_visible('demand_elas', False)
1481        self.set_type_visible('mod_elas', False)
1482        self.set_type_visible('supply_rc', False)
1483        self.set_type_visible('demand_rc', False)
1484        self.set_type_visible('supply_partial_rc', False, affzero=True)
1485        self.set_type_visible('demand_partial_rc', False, affzero=True)
1486
1487        for line in self.huedict.keys():
1488            if line in self.mod_elas:
1489                continue
1490            if len(self.supply_flux) == 1:
1491                self.set_visible(list(self.supply_flux.keys())[0], f=True)
1492            if len(self.demand_flux) == 1:
1493                self.set_visible(list(self.demand_flux.keys())[0], f=True)
1494
1495            filename = os.path.join(
1496                self.parent.parent._psc_out_dir,
1497                self.parent.fix,
1498                self.parent.modelfilename[:-4]
1499                + '_'
1500                + self.parent.fix
1501                + '_ElasRC_'
1502                + line
1503                + '_',
1504            )
1505
1506            self.title_label = (
1507                'Model: '
1508                + self.parent.parent._mod
1509                + ' -Response vs. elasticity coefficients - '
1510                + line
1511            )
1512            self.set_visible(line, f=True, e=True, r=True)
1513            if file:
1514                self.save(filename)
1515            else:
1516                self.show()
1517                input("Press Enter to continue...")
1518            self.set_visible(line, f=False, e=False, r=False)
1519
1520    def plotPartialRCIndiv(self, file=0):
1521        """
1522        plotPartialRC(file=0)
1523
1524        Plot the data in the instantiated RateCharData class for each reaction
1525        in turn (fluxes and partial response coefficients vs. the
1526        concentration of the fixed species).
1527
1528        Arguments:
1529        ==========
1530        file (bool) - save plot as PNG file (default=0, plots to screen)
1531        """
1532
1533        if self.supply_total:
1534            self.supply_total[0].set_visible(True)
1535        if self.demand_total:
1536            self.demand_total[0].set_visible(True)
1537        self.set_type_visible('supply_flux', False)
1538        self.set_type_visible('demand_flux', False)
1539        self.set_type_visible('supply_elas', False)
1540        self.set_type_visible('demand_elas', False)
1541        self.set_type_visible('mod_elas', False)
1542        self.set_type_visible('supply_rc', False)
1543        self.set_type_visible('demand_rc', False)
1544        self.set_type_visible('supply_partial_rc', False, affzero=True)
1545        self.set_type_visible('demand_partial_rc', False, affzero=True)
1546
1547        for line in self.huedict.keys():
1548            if line in self.mod_elas:
1549                continue
1550            if line not in list(self.reaction_partial_rc.keys()):
1551                continue
1552            if len(self.supply_flux) == 1:
1553                self.set_visible(list(self.supply_flux.keys())[0], f=True)
1554            if len(self.demand_flux) == 1:
1555                self.set_visible(list(self.demand_flux.keys())[0], f=True)
1556
1557            filename = os.path.join(
1558                self.parent.parent._psc_out_dir,
1559                self.parent.fix,
1560                self.parent.modelfilename[:-4]
1561                + '_'
1562                + self.parent.fix
1563                + '_PartialRC_'
1564                + line
1565                + '_',
1566            )
1567
1568            self.title_label = (
1569                'Model: '
1570                + self.parent.parent._mod
1571                + ' - Partial response coefficients - '
1572                + line
1573            )
1574            self.set_visible(line, f=True, p=True, r=True)
1575            if file:
1576                self.save(filename)
1577            else:
1578                self.show()
1579                input("Press Enter to continue...")
1580            self.set_visible(line, f=False, p=False, r=False)
1581