michellab / Sire

Sire Molecular Simulations Framework
http://siremol.org
GNU General Public License v3.0
95 stars 26 forks source link

FreeEnergyAnalysis.py Error #309

Closed lefnelso closed 4 years ago

lefnelso commented 4 years ago

After an FEP run I am trying to obtain the MBAR and TI values but the following error keeps occurring (I am using python 3.7 and openmm 7.3.1):

Traceback (most recent call last): File "/home/laurennelson/sire-combRules.app/share/Sire/scripts/analyse_freenrg.py", line 621, in free_energy_obj.run_mbar(test_overlap) File "/home/laurennelson/sire-combRules.app/lib/python3.7/site-packages/Sire/Tools/FreeEnergyAnalysis.py", line 85, in run_mbar (deltaF_ij, dDeltaF_ij, theta_ij) = MBAR_obj.getFreeEnergyDifferences() ValueError: not enough values to unpack (expected 3, got 2)

Any thoughts?

lefnelso commented 4 years ago

The analysis can be run on other people's computers too ( @SofiaBariami ) so the output files are fine.

jmichel80 commented 4 years ago

Can you check what version of pymbar you are using ? Maybe @ppxasjsm can comment, but my immediate thought is that this could be due to an API change. Share more details on how you compiled /home/laurennelson/sire-combRules.app/

ppxasjsm commented 4 years ago

There has been some changes to MBAR recently. Version information would be very helpful. I will likely have to update the way we wrap MBAR once they've finished with the most recent version of the revamp. What is the current version you get?

lohedges commented 4 years ago

This is from the latest mbar.py file:

  def getFreeEnergyDifferences(self, compute_uncertainty=True,
    uncertainty_method=None, warning_cutoff=1.0e-10, return_theta=False, return_dict=False):
        """Get the dimensionless free energy differences and uncertainties
           among all thermodynamic states.

        Parameters
        ----------
        compute_uncertainty : bool, optional
            If False, the uncertainties will not be computed (default: True)
        uncertainty_method : string, optional
            Choice of method used to compute asymptotic covariance method,
            or None to use default.  See help for computeAsymptoticCovarianceMatrix()
            for more information on various methods. (default: svd)
        warning_cutoff : float, optional
            Warn if squared-uncertainty is negative and larger in magnitude
            than this number (default: 1.0e-10)
        return_theta : bool, optional
            Whether or not to return the theta matrix.  Can be useful for complicated differences.
        return_dict: bool, default False
            If true, returns are in a dictionary otherwise a tuple is returned
        Returns
        -------
        'Delta_f' : np.ndarray, float, shape=(K, K)
            Deltaf_ij[i,j] is the estimated free energy difference
            If return_dict, key is 'Delta_f'
        'dDelta_f' : np.ndarray, float, shape=(K, K)
            If compute_uncertainty==True,
            dDeltaf_ij[i,j] is the estimated statistical uncertainty
            (one standard deviation) in Deltaf_ij[i,j].  Otherwise not included.
            If return_dict, key is 'dDelta_f'
        'Theta' : np.ndarray, float, shape=(K, K)
            The theta_matrix if return_theta==True, otherwise not included.
            If return_dict, key is 'Theta'

It looks like the last two return arguments are now optional. To get dDelta_f you need to pass compute_uncertainty=True and for Theta you need return_theta=True.

It should be easy to wrap the function call in a try-except block, or check the version of pymbar and use an if-else block to use the correct function and return signature.

lohedges commented 4 years ago

Sorry, looking again at the docstring the default is compute_uncertainty=True so you only need to add return_theta=True. This is why you are seeing two return arguments rather than three.

Cheers.

lohedges commented 4 years ago

Bizarrely it looks like that return argument has been optional for some time...

Could you try replacing line 85 of /home/laurennelson/sire-combRules.app/lib/python3.7/site-packages/Sire/Tools/FreeEnergyAnalysis.py with the following:

        try:
            (deltaF_ij, dDeltaF_ij, theta_ij) = MBAR_obj.getFreeEnergyDifferences()
        except:
            (deltaF_ij, dDeltaF_ij, theta_ij) = MBAR_obj.getFreeEnergyDifferences(return_theta=True)

Let me know if it fixes your problem and I'll push the changes.

ppxasjsm commented 4 years ago

It might be that the conda channel was only updated recently?

lohedges commented 4 years ago

Possibly. I did download a two year old archive from their conda-forge channel and the return argument was optional even then.

lefnelso commented 4 years ago

Hi all,

Thanks for your help this fix sorted the issue: try: (deltaF_ij, dDeltaF_ij, theta_ij) = MBAR_obj.getFreeEnergyDifferences() except: (deltaF_ij, dDeltaF_ij, theta_ij) = MBAR_obj.getFreeEnergyDifferences(return_theta=True)