-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathseq.rb
executable file
·83 lines (64 loc) · 1.47 KB
/
seq.rb
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
78
79
80
81
82
83
#!/usr/bin/env ruby -wKU
# seq.rb tests sequence and quality parser
require "rubygems"
require "positional"
class UnityArray < Array
alias :old_ndx :[]
alias :old_each_with_index :each_with_index
def [](index)
raise StandardError.new("Range error (#{index < 1 ? '<' : '>'})") if index < 1 or index > length
old_ndx(index-1)
end
def each_with_index(&block)
old_each_with_index {|e,i| yield e,i+1}
end
end
class SequenceQuality
class SQ
attr :base, :phred
def initialize(b,p)
@base=b
@phred=p
end
end
attr :bases, :phreds
def initialize(bs, sq)
@bases = bs
@phreds = sq
end
def [](index)
SQ.new @bases[index], @phreds[index]
end
def cutoff(q)
(10**(-q/10)).to_f
end
def mask(len)
str=''
@bases.each_with_index do |b,i|
if @phreds[i] >= len
str << b
else
str << 'X'
end
end
str
end
end
bases='GCTACTGCAAGTTCTAGACT'
qualities='1 8 15 22 60 55 57 56 58 55 60 58 57 59 55 50 44 18 5 6'
bs=UnityArray.new 20
sq=UnityArray.new 20
pb = Positional::MaskedParse.new bs.extend Positional::Decorator::Index
pq = Positional::Parse.new sq.extend Positional::Decorator::Index
bs = pb.parse bases
sq = pq.parse qualities
foo = SequenceQuality.new bs, sq
puts foo[4].base
puts foo[4].phred
puts foo.mask(20)
puts 'XXXACTGCAAGTTCTAGXXX'
puts "Recovery ---"
bf = Positional::MaskedFormat.new foo.bases
qf = Positional::Format.new foo.phreds
puts bf.to_s
puts qf.to_s