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.

Sign up for GitHub

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

Report position in read for reads in a pileup #260

Closed
wants to merge 2 commits into from
Closed

Report position in read for reads in a pileup #260

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
1a04afc
Report alignment position for pileup
Jun 22, 2022
995cc1d
Merge remote-tracking branch 'upstream/master'
Jun 22, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
17 changes: 10 additions & 7 deletions src/cljam/algo/pileup.clj
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -48,17 +48,18 @@
93)))
short-array))))

(deftype PosBase [^char base indel])
(deftype PosBase [^char base indel ^int pos])

(defn- seqs-at-ref
[idx ^String s]
(->> idx
(map (fn [[op x xs]]
(let [c (if (number? x) (.charAt s x) x)]
(let [c (if (number? x) (.charAt s x) x)
rel-pos (if (number? x) (inc x) -1)]
(case op
:m (->PosBase c nil)
:d (->PosBase c xs)
:i (->PosBase c (subs s (first xs) (last xs)))))))
:m (->PosBase c nil rel-pos)
:d (->PosBase c xs rel-pos)
:i (->PosBase c (subs s (first xs) (last xs)) rel-pos)))))
object-array))

(defn index-cigar
Expand Down Expand Up @@ -92,7 +93,8 @@
^PosBase pb (aget ^objects (:seqs-at-ref aln) relative-pos)
base (.base pb)
indel (.indel pb)
deletion? (number? indel)]
qpos (.pos pb)
deletion? (= \* base)]
(plpio/->PileupBase
(zero? relative-pos)
(.mapq aln)
Expand All @@ -101,7 +103,8 @@
(flag/reversed? (.flag aln))
(= ref-pos (.end aln))
(when-not deletion? indel)
(when deletion? indel)
(when deletion? indel) ;; TODO this is incorrect IMO - D in the cigar refers to the base _after_ this in the ref, so this field is set for the wrong column!
qpos
(.qname aln)
(dissoc aln :seqs-at-ref :quals-at-ref))))

Expand Down
5 changes: 3 additions & 2 deletions src/cljam/io/pileup.clj
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@
^boolean end?
insertion ;; String?
deletion ;; int?
^int pos
qname ;; String?
alignment]) ;; SAMAlignment?

Expand Down Expand Up @@ -68,7 +69,7 @@
(when (zero? upper-base-int)
(throw (ex-info (format "Invalid character %s in %s" base column) {:column column :base base})))
(if (= len (+ i (if mapq 3 1)))
(persistent! (conj! results (PileupBase. start? mapq upper-base -1 reverse? false nil nil nil nil)))
(persistent! (conj! results (PileupBase. start? mapq upper-base -1 reverse? false nil nil -1 nil nil)))
(let [x (.charAt column (+ i (if mapq 3 1)))
ins? (= x \+)
del? (= x \-)
Expand All @@ -91,7 +92,7 @@
(if indel-num (+ 1 indel-num-chars indel-num) 0)
(if end? 1 0)))]
(recur next-pos
(conj! results (PileupBase. start? mapq upper-base -1 reverse? end? (when ins? indel-seq) (when del? indel-num) nil nil))))))
(conj! results (PileupBase. start? mapq upper-base -1 reverse? end? (when ins? indel-seq) (when del? indel-num) -1 nil nil))))))
(persistent! results)))))

(defn- parse-pileup-line
Expand Down
25 changes: 13 additions & 12 deletions test/cljam/io/pileup_test.clj
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,7 @@
:end? true,
:insertion nil,
:deletion nil,
:pos -1,
:qname nil,
:alignment nil}]})

Expand All @@ -50,7 +51,7 @@
(deftest parse-bases-col
(testing "without-ref"
(are [?in ?out]
(= (mapv #(into {:qname nil :alignment nil} %) ?out)
(= (mapv #(into {:pos -1, :qname nil :alignment nil} %) ?out)
(mapv rec->map (#'plpio/parse-bases-col nil ?in)))
"" []
"A" [{:start? false, :mapq nil, :base \A, :qual -1, :reverse? false, :end? false, :insertion nil, :deletion nil}]
Expand Down Expand Up @@ -78,7 +79,7 @@
{:start? true, :mapq 34, :base \A, :qual -1, :reverse? true, :end? false, :insertion nil, :deletion nil}]))
(testing "with-ref"
(are [?in ?out]
(= (mapv #(into {:qname nil :alignment nil} %) ?out)
(= (mapv #(into {:pos -1, :qname nil :alignment nil} %) ?out)
(mapv rec->map (#'plpio/parse-bases-col \T ?in)))
"" []
"A" [{:start? false, :mapq nil, :base \A, :qual -1, :reverse? false, :end? false, :insertion nil, :deletion nil}]
Expand Down Expand Up @@ -117,26 +118,26 @@
"chr1\t10\tN\t0\t\t" {:rname "chr1", :pos 10, :ref \N, :count 0,
:pile []}
"chr1\t10\tN\t1\tA\tI" {:rname "chr1", :pos 10, :ref \N, :count 1,
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? false, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}]}
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? false, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}]}
"chr1\t10\tN\t2\taA\tIB" {:rname "chr1", :pos 10, :ref \N, :count 2,
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}]}
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}]}
"chr1\t10\tN\t2\taA-2NN$\tIB" {:rname "chr1", :pos 10, :ref \N, :count 2,
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? true, :insertion nil, :deletion 2, :qname nil, :alignment nil}]}))
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? true, :insertion nil, :deletion 2, :pos -1, :qname nil, :alignment nil}]}))
(testing "with-ref"
(are [?in ?out]
(= ?out (rec->map (#'plpio/parse-pileup-line ?in)))
"chr1\t10\tA\t0\t\t" {:rname "chr1", :pos 10, :ref \A, :count 0,
:pile []}
"chr1\t10\ta\t1\t.\tI" {:rname "chr1", :pos 10, :ref \a, :count 1,
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? false, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}]}
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? false, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}]}
"chr1\t10\tA\t2\t,.\tIB" {:rname "chr1", :pos 10, :ref \A, :count 2,
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}]}
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}]}
"chr1\t10\tA\t2\t,.-2CA$\tIB" {:rname "chr1", :pos 10, :ref \A, :count 2,
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? true, :insertion nil, :deletion 2, :qname nil, :alignment nil}]})))
:pile [{:start? false, :mapq nil, :base \A, :qual 40, :reverse? true, :end? false, :insertion nil, :deletion nil, :pos -1, :qname nil, :alignment nil}
{:start? false, :mapq nil, :base \A, :qual 33, :reverse? false, :end? true, :insertion nil, :deletion 2, :pos -1, :qname nil, :alignment nil}]})))

(deftest reader
(with-open [r (plpio/reader test-pileup-file)]
Expand Down