-
Notifications
You must be signed in to change notification settings - Fork 0
/
coverage.go
77 lines (65 loc) · 1.78 KB
/
coverage.go
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
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
package main
import (
_ "embed"
"fmt"
"log"
"strconv"
"strings"
)
type Gene struct {
chr string
start int
end int
ensg string
symbol string
}
//go:embed genes.bed
var genesFile string
func getGeneCoverage(targets map[string][]TargetCoverage) map[string]Coverage {
geneLines := strings.Split(string(genesFile), "\n")
geneCoverageMap := make(map[string]Coverage)
for _, line := range geneLines {
if line == "" {
continue
}
fields := strings.Split(line, "\t")
start, err := strconv.Atoi(fields[1])
if err != nil {
log.Fatalf("error parsing start, %v", err)
}
end, err := strconv.Atoi(fields[2])
if err != nil {
log.Fatalf("error parsing end, %v", err)
}
gene := Gene{
chr: fields[0],
start: start,
end: end,
ensg: fields[3],
symbol: fields[4],
}
var overlapCount, totalTargets, totalCovered5x, totalCovered10x, totalCovered20x, totalCovered30x int
for _, target := range targets[gene.chr] {
if !(gene.end < target.start || gene.start > target.end) {
overlapCount++
totalTargets += target.end - target.start
totalCovered5x += target.x5
totalCovered10x += target.x10
totalCovered20x += target.x20
totalCovered30x += target.x30
}
}
geneCoverageMap[gene.ensg] = Coverage{
ENSG: gene.ensg,
Symbol: gene.symbol,
TotalBases: totalTargets,
CoveredBases: totalCovered10x,
BasesCovered5x: totalCovered5x,
BasesCovered10x: totalCovered10x,
BasesCovered20x: totalCovered20x,
BasesCovered30x: totalCovered30x,
}
fmt.Printf("%s\t%s\tCount: %d\tTotal Bases: %d\tCoverage 5x: %d\t10x: %d\t20x: %d\t30x: %d\n", gene.ensg, gene.symbol, overlapCount, totalTargets, totalCovered5x, totalCovered10x, totalCovered20x, totalCovered30x)
}
return geneCoverageMap
}