Skip to content

tm returns incorrect values on non-pcr preset #25

@phlaster

Description

@phlaster

I stumbled upon this one, when I was trying to find out, how does your pcr kwarg work. Here how it was:

>>> seq = "ACGAGCGAGCGGAGCGGGCG"
>>> import seqfold
>>> seqfold.tm(seq)
72.6
>>> seqfold.tm(seq, pcr=True)
72.6
>>> seqfold.tm(seq, pcr=False)
0.00020166276092434078
>>> seqfold.tm("ACTGTCTATCTACGATCGACGACTAGCTTTCACACACGATCGACGGACTAGCG", pcr=True)
75.6
>>> seqfold.tm("ACTGTCTATCTACGATCGACGACTAGCTTTCACACACGATCGACGGACTAGCG", pcr=False)
0.00020203618743072687

As you can see, with non-pcr preset tm returns a suspiciously not-rounded float. This is caused by:

Early Return in Salt Correction Logic: When R < 0.22, the function returns the salt correction term instead of the Tm.

Another problem here is:

GC% Scaling Error in Salt Correction & misleading docstrings: your docstrings propose gc is a % (so between 0 and 100), then you gc/100 to (presumably) turn it into ratio but it is already in [0,1], as it comes from _gc function! So you underestimate gc by the factor of 2.

Yet another issue is:

Negative concentration may cause logarithm domain error: k = (seq1_conc - seq2_conc / 2.0) * 1e-9 can be negative (e.g., if seq1_conc < seq2_conc/2), causing math.log(k) to fail.

These are is from the top of my head, but I was dealing with some other problems in the other functions as well while translating your code... You can take my Julia implementation to verify yours or ask gpt to translate the code back. I haven't altered function call tree that much, so that shouldn't be hard.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions