File: ClusterProbability.cpp

package info (click to toggle)
pbseqlib 5.3.1%2Bdfsg-2.1
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 7,136 kB
  • sloc: cpp: 77,246; python: 570; makefile: 312; sh: 111; ansic: 9
file content (44 lines) | stat: -rw-r--r-- 1,178 bytes parent folder | download | duplicates (4)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
#include <alignment/algorithms/anchoring/ClusterProbability.hpp>

#include <cassert>
#include <cmath>

float ComputeAnchorProbability(float pMatch, int minAnchorLength)
{
    assert(minAnchorLength >= 0);
    assert(pMatch < 1.00001 and pMatch > 0);

    int i;
    float totalProbability = 0.0;
    float pMisMatch = 1 - pMatch;
    for (i = 0; i < minAnchorLength; i++) {
        totalProbability += pMatch * pMisMatch;
        pMatch *= pMatch;
    }
    return 1 - totalProbability;
}

float ComputeExpectedAnchorLength(float pMatch, int minAnchorLength, float pAnchor)
{
    int i = 0;
    for (i = 0; i < minAnchorLength; i++) {
        pMatch *= pMatch;
    }
    float pMismatch = 1 - pMatch;
    float pValueEpsilon = 0.000000001;
    float expAnchorLength = 0;

    while (pMatch * pMismatch > pValueEpsilon) {
        expAnchorLength += i * pMatch * pMismatch / pAnchor;
        pMatch *= pMatch;
    }
    return expAnchorLength;
}

float AnchorBasesPerRead(int readLength, float pAnchor) { return pAnchor * readLength; }

float AnchorBasesPerReadSigma(float expAnchorBasesPerRead)
{
    // Assume exponential distribution:
    return sqrt(expAnchorBasesPerRead);
}