Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

bin value are not unique #477

Open
shawnhsueh opened this issue Oct 30, 2022 · 3 comments
Open

bin value are not unique #477

shawnhsueh opened this issue Oct 30, 2022 · 3 comments

Comments

@shawnhsueh
Copy link

shawnhsueh commented Oct 30, 2022

I found that in line 544-548 of pymbar/fes.py you wrote

                # how do we label the bins? if N-dimensional:
                # bins[0] + bins[1]*bin_length[1]**1 + bins[2]*bin_length[2]**2
                sample_label[n] = int(
                    np.sum([bin_n[n][d] * bin_length[d] ** d for d in range(dims)])
                )

However, this setting can not give unique bin values. For example, given bin_length=[100,2], bins=(30,1) and bins=(32,0) would have the same bin value of 32 (i.e. 30*100**0+1*2**1=32*100**0+0*2**1=32).

Non-unique bin values would cause under-counting of histogram_data["bin_order"] in line 566-576 of pymbar/fes.py

        if b == 0:
            bin_order = {}
            i = 0
            for bv in bin_label.values():
                if bv not in bin_order:
                    bin_order[bv] = i
                    i += 1
            histogram_data["bin_order"] = bin_order
            histogram_data["bin_label"] = bin_label
        else:
            bin_order = self.histogram_data["bin_order"]

I would suggest to correct line 546-548 of pymbar/fes.py to be

                sample_label[n] = int(
                    np.sum([bin_n[n][d] * np.sum(bin_length[:d+1]) ** d for d in range(dims)])
                )
@mikemhenry
Copy link
Contributor

@shawnhsueh would you be willing to make this change and submit a PR?

Thanks!

@mikemhenry
Copy link
Contributor

We can then discuss the proposed change there. Normally it is good to talk about the plan for the change before writing the code, but this seems like a small enough fix to jump straight to a PR (pull request).

@shawnhsueh
Copy link
Author

OK. I've created a pull request #478.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants