Repository: chrovis/cljam
Pull request #166 (2019-05-23 04:46:48)

Normalize variants in a VCF file
~~**This PR is depending on #164**~~

#### Summary
Add `cljam.io.vcf.util.normalize/normalize` which converts variants into canonical forms.
cf. `bcftools norm -m-`

See the test cases for typical conversion examples.

```
ref0: TATTTTTTTTTTTT TC  TTTTT
alt1:                TTCC     <- not right-trimmed
alt2:                TTC      <- not right-trimmed

ref0: TATTTTTTTTTTTT T  CTTTTT
alt1:                TTC
alt2:                TT       <- not left-aligned

ref0: T A TTTTTTTTTTTT T  CTTTTT
alt1:                  TTC
alt2:   AT
```

**input** multiallelic, not left-aligned, not right-trimmed
https://github.com/chrovis/cljam/blob/41881832c20f96c25d91027a73cab09a34378950/test/cljam/io/vcf/util/normalize_test.clj#L328-L330

**output** biallelic, left-aligned, right-trimmed
https://github.com/chrovis/cljam/blob/41881832c20f96c25d91027a73cab09a34378950/test/cljam/io/vcf/util/normalize_test.clj#L313-L318

#### Tests

- `lein check` 🆗
- `lein test :all` 🆗
- `lein cloverage` 🆗
- `lein eastwood` 🆗


codecov[bot]: # [Codecov](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=h1) Report
> Merging [#166](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=desc) into [master](https://codecov.io/gh/chrovis/cljam/commit/7b5d627de556181632afb1051c6c4a70555431d4?src=pr&el=desc) will **increase** coverage by `0.38%`.
> The diff coverage is `97.38%`.

[![Impacted file tree graph](https://codecov.io/gh/chrovis/cljam/pull/166/graphs/tree.svg?width=650&token=hrCuiz6NWa&height=150&src=pr)](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=tree)

```diff
@@            Coverage Diff             @@
##           master     #166      +/-   ##
==========================================
+ Coverage    86.1%   86.48%   +0.38%     
==========================================
  Files          71       72       +1     
  Lines        5419     5610     +191     
  Branches      489      493       +4     
==========================================
+ Hits         4666     4852     +186     
- Misses        264      265       +1     
- Partials      489      493       +4
```


| [Impacted Files](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=tree) | Coverage Δ | |
|---|---|---|
| [src/cljam/io/vcf/util/normalize.clj](https://codecov.io/gh/chrovis/cljam/pull/166/diff?src=pr&el=tree#diff-c3JjL2NsamFtL2lvL3ZjZi91dGlsL25vcm1hbGl6ZS5jbGo=) | `96.92% <96.92%> (ø)` | |
| [src/cljam/io/vcf/util.clj](https://codecov.io/gh/chrovis/cljam/pull/166/diff?src=pr&el=tree#diff-c3JjL2NsamFtL2lvL3ZjZi91dGlsLmNsag==) | `92.43% <98.36%> (+1.9%)` | :arrow_up: |

------

[Continue to review full report at Codecov](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=continue).
> **Legend** - [Click here to learn more](https://docs.codecov.io/docs/codecov-delta)
> `Δ = absolute <relative> (impact)`, `ø = not affected`, `? = missing data`
> Powered by [Codecov](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=footer). Last update [7b5d627...4188183](https://codecov.io/gh/chrovis/cljam/pull/166?src=pr&el=lastupdated). Read the [comment docs](https://docs.codecov.io/docs/pull-request-comments).


Patch:
diff --git a/src/cljam/io/vcf/util/check.clj b/src/cljam/io/vcf/util/check.clj
new file mode 100644
index 0000000..c51c2eb
--- /dev/null
+++ b/src/cljam/io/vcf/util/check.clj
@@ -0,0 +1,10 @@
+(ns cljam.io.vcf.util.check
+  (:require [cljam.io.sequence :as io-seq]))
+
+(defn same-ref?
+  "Checks if a REF allele matches the reference sequence."
+  [seq-reader {:keys [chr ^long pos ^String ref]}]
+  (let [ref-len (.length ref)
+        ref-region {:chr chr :start pos :end (dec (+ pos ref-len))}
+        ref-seq (io-seq/read-sequence seq-reader ref-region {:mask? false})]
+    (and ref-seq (.equalsIgnoreCase ^String ref-seq ref))))
diff --git a/src/cljam/io/vcf/util/normalize.clj b/src/cljam/io/vcf/util/normalize.clj
new file mode 100644
index 0000000..ec32bf3
--- /dev/null
+++ b/src/cljam/io/vcf/util/normalize.clj
@@ -0,0 +1,155 @@
+(ns cljam.io.vcf.util.normalize
+  (:require [clojure.string :as cstr]
+            [cljam.io.vcf.util :as vcf-util]
+            [cljam.io.sequence :as io-seq]))
+
+(defn- biallelic-map
+  [m id->number ^long i ^long n-allele]
+  (let [ploidy (or (some-> m :GT vcf-util/parse-genotype count) 2)]
+    (into
+     {}
+     (map
+      (fn [[k v]]
+        (let [n (id->number k)]
+          [k (cond
+               (= k :GT) (vcf-util/biallelic-genotype v (inc i))
+               (= n "A") [(nth v i)]
+               (= n "R") [(first v) (nth v (inc i))]
+               (= n "G") (vcf-util/biallelic-coll ploidy n-allele (inc i) v)
+               :else v)])))
+     m)))
+
+(defn- get-id->number
+  [coll]
+  (into {} (map (juxt (comp keyword :id) :number)) coll))
+
+(defn make-splitter
+  "Returns a function that splits a multiallelic variant into a sequence of
+  biallelic variants."
+  [{:keys [info format]} header]
+  (let [info-id->number (get-id->number info)
+        format-id->number (get-id->number format)
+        samples (map keyword (drop 9 header))]
+    (fn split [{:keys [alt] :as v}]
+      (let [n (count alt)]
+        (if (<= n 1)
+          [v]
+          (->> alt
+               (map-indexed
+                (fn sample-ith-alt [i allele]
+                  (-> v
+                      (assoc :alt [allele])
+                      (update :info biallelic-map info-id->number i n)
+                      (into
+                       (map
+                        (fn [s]
+                          [s (biallelic-map (v s) format-id->number i n)]))
+                       samples))))))))))
+
+(defn- regions-before [chr ^long pos ^long window]
+  (->> (range pos 1 (- window))
+       (map
+        (fn [^long e]
+          {:chr chr,
+           :start (Math/max 1 (- e window)),
+           :end (dec e)}))))
+
+(defn- char-equals-ignore-case?
+  ([_]
+   true)
+  ([x y]
+   (= (Character/toUpperCase ^char x)
+      (Character/toUpperCase ^char y)))
+  ([x y & z]
+   (and (char-equals-ignore-case? x y)
+        (apply char-equals-ignore-case? y z))))
+
+(defn- trim-right
+  "Trims all duplicated allele sequences from right. If reached the start,
+  reference sequences of length `window` will be read from `seq-reader` and
+  be prepended to the alleles."
+  [seq-reader {:keys [chr ^long pos ref alt]
+               {:keys [END]} :info :as v} ^long window]
+  (let [alleles (cons ref alt)
+        min-end (long (dec (+ pos (apply min (map count alleles)))))
+        regs (regions-before chr pos window)
+        ref-seqs (sequence ;; unchunk
+                  (mapcat
+                   (fn [r]
+                     (-> seq-reader
+                         (io-seq/read-sequence r {:mask? false})
+                         cstr/reverse)))
+                  regs)
+        allele-seqs (map #(concat (cstr/reverse %) ref-seqs) alleles)
+        matched-length (->> allele-seqs
+                            (apply map char-equals-ignore-case?)
+                            (take-while true?)
+                            count
+                            long
+                            (Math/min (dec min-end)))]
+    (if (pos? matched-length)
+      (let [pos' (or (some (fn [{:keys [start end]}]
+                             (when (<= start (- min-end matched-length) end)
+                               start)) regs)
+                     pos)
+            [ref' & alt'] (map
+                           #(->> %2
+                                 (take (inc (- (dec (+ pos (count %1))) pos')))
+                                 (drop matched-length)
+                                 reverse
+                                 cstr/join)
+                           alleles
+                           allele-seqs)]
+        (cond-> (assoc v :pos pos' :ref ref' :alt alt')
+          (integer? END) (assoc-in [:info :END] (dec (+ pos' (count ref'))))))
+      v)))
+
+(defn- trim-left
+  "Trims all duplicated allele sequences from left."
+  [{:keys [ref alt] :as v}]
+  (let [alleles (cons ref alt)
+        n-trim-left (->> alleles
+                         (apply map char-equals-ignore-case?)
+                         (take-while true?)
+                         count)
+        min-len (long (apply min (map count alleles)))
+        n-trim-left (min (dec min-len) n-trim-left)]
+    (if (pos? n-trim-left)
+      (-> v
+          (update :pos + n-trim-left)
+          (update :ref subs n-trim-left)
+          (assoc :alt (map #(subs % n-trim-left) alt)))
+      v)))
+
+(defn realign
+  "Left-align and trim REF/ALT alleles."
+  ([seq-reader variant]
+   (realign seq-reader variant {}))
+  ([seq-reader
+    {:keys [ref alt] :as variant}
+    {:keys [window] :or {window 100}}]
+   (cond
+     (some->> (seq alt)
+              (every? (comp #{:snv :mnv :insertion :deletion :complex}
+                            :type
+                            (partial vcf-util/inspect-allele ref))))
+     (trim-left (trim-right seq-reader variant window))
+
+     (and ref (nil? (seq alt)))
+     (cond-> (update variant :ref subs 0 1)
+       (integer? (:END (:info variant)))
+       (assoc-in [:info :END] (:pos variant)))
+
+     :else variant)))
+
+(defn normalize
+  "Splits multiallelic variants and realigns all alleles. Returns a lazy
+  sequence of normalized variants. Note that the result may be unsorted due to
+  the realignment. Returns a transducer when no variant is provided. The
+  transducer does not guarantee thread safety of `ref-reader`."
+  ([ref-reader meta-info header]
+   (comp
+    (mapcat (make-splitter meta-info header))
+    (map (partial realign ref-reader))))
+  ([ref-reader meta-info header variants]
+   (sequence (normalize ref-reader meta-info header) variants)))


Test patch:
diff --git a/test/cljam/io/vcf/util/check_test.clj b/test/cljam/io/vcf/util/check_test.clj
new file mode 100644
index 0000000..b6f911e
--- /dev/null
+++ b/test/cljam/io/vcf/util/check_test.clj
@@ -0,0 +1,17 @@
+(ns cljam.io.vcf.util.check-test
+  (:require [clojure.test :refer :all]
+            [cljam.io.protocols :as protocols]
+            [cljam.io.vcf.util.check :as check]))
+
+(defn- seq-reader [s]
+  (reify protocols/ISequenceReader
+    (read-sequence [_ {:keys [^long start ^long end]} _]
+      (subs s (dec start) end))))
+
+(deftest same-ref?
+  (are [?seq ?variant ?expected]
+       (= ?expected (check/same-ref? (seq-reader ?seq) ?variant))
+    "A" {:pos 1, :ref "A"} true
+    "ATT" {:pos 2, :ref "tt"} true
+    "TTT" {:pos 3, :ref "A"} false
+    "ATGC" {:pos 4, :ref "N"} false))
diff --git a/test/cljam/io/vcf/util/normalize_test.clj b/test/cljam/io/vcf/util/normalize_test.clj
new file mode 100644
index 0000000..a418682
--- /dev/null
+++ b/test/cljam/io/vcf/util/normalize_test.clj
@@ -0,0 +1,330 @@
+(ns cljam.io.vcf.util.normalize-test
+  (:require
+   [clojure.test :refer :all]
+   [cljam.test-common :refer :all]
+   [cljam.io.vcf.util.normalize :as norm]
+   [cljam.io.protocols :as protocols]
+   [cljam.io.sequence :as io-seq]
+   [cljam.util.chromosome :as chr]))
+
+(deftest make-splitter
+  (let [meta-info {:info [{:id "XF" :number 0}
+                          {:id "DP" :number 1}
+                          {:id "XT" :number 2}
+                          {:id "AC" :number "A"}
+                          {:id "XR" :number "R"}
+                          {:id "XG" :number "G"}]
+                   :format [{:id "XF" :number 0}
+                            {:id "DP" :number 1}
+                            {:id "XT" :number 2}
+                            {:id "GT" :number 1}
+                            {:id "AC" :number "A"}
+                            {:id "XR" :number "R"}
+                            {:id "XG" :number "G"}]}
+        header ["CHROM" "POS" "ID" "REF" "ALT" "QUAL" "FILTER" "INFO"
+                "FORMAT" "SAMPLE01" "SAMPLE02"]]
+    (are [?variant ?expected]
+         (= ?expected ((norm/make-splitter meta-info header) ?variant))
+
+      {:chr "1", :pos 1, :ref "T", :alt ["A"],
+       :id nil, :qual 100.0, :filter [:PASS],
+       :info {:XF :exists, :DP 10, :XT [1 2],
+              :AC [2], :XR [100.0 200.0], :XG [0 1 2]}
+       :FORMAT [:XF :DP :XT :GT :AC :XR :XG],
+       :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20],
+                  :GT "0/1", :AC [40], :XR [2 3], :XG [10 11 12]},
+       :SAMPLE02 {:XF :exists, :DP 300, :XT [30 40],
+                  :GT "0/0/1", :AC [80], :XR [1 7], :XG [0 1 2]}}
+
+      [{:chr "1", :pos 1, :ref "T", :alt ["A"],
+        :id nil, :qual 100.0, :filter [:PASS],
+        :info {:XF :exists, :DP 10, :XT [1 2],
+               :AC [2], :XR [100.0 200.0], :XG [0 1 2]}
+        :FORMAT [:XF :DP :XT :GT :AC :XR :XG],
+        :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20],
+                   :GT "0/1", :AC [40], :XR [2 3], :XG [10 11 12]},
+        :SAMPLE02 {:XF :exists, :DP 300, :XT [30 40],
+                   :GT "0/0/1", :AC [80], :XR [1 7], :XG [0 1 2]}}]
+
+      ;; -----
+
+      {:chr "2", :pos 20, :ref "T", :alt ["A", "C"],
+       :id nil, :qual 100.0, :filter [:PASS],
+       :info {:XF :exists, :DP 10, :XT [1 2],
+              :AC [2 4], :XR [100.0 200.0 300.0], :XG [0 1 2 3 4 5]}
+       :FORMAT [:XF :DP :XT :GT :AC :XR :XG],
+       :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20], :GT "0/2",
+                  :AC [40 50], :XR [2 3 4], :XG [10 11 12 13 14 15]},
+       :SAMPLE02 {:XF :exists, :DP 300, :XT [30 40], :GT "0/1/2",
+                  :AC [80 90], :XR [1 7 9], :XG [0 1 2 3 4 5 6 7 8 9]}}
+
+      [{:chr "2", :pos 20, :ref "T", :alt ["A"],
+        :id nil, :qual 100.0, :filter [:PASS],
+        :info {:XF :exists, :DP 10, :XT [1 2],
+               :AC [2], :XR [100.0 200.0], :XG [0 1 2]}
+        :FORMAT [:XF :DP :XT :GT :AC :XR :XG],
+        :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20], :GT "0/0",
+                   :AC [40], :XR [2 3], :XG [10 11 12]},
+        :SAMPLE02 {:XF :exists, :DP 300, :XT [30 40], :GT "0/1/0",
+                   :AC [80], :XR [1 7], :XG [0 1 2 3]}}
+       {:chr "2", :pos 20, :ref "T", :alt ["C"],
+        :id nil, :qual 100.0, :filter [:PASS],
+        :info {:XF :exists, :DP 10, :XT [1 2],
+               :AC [4], :XR [100.0 300.0], :XG [0 3 5]}
+        :FORMAT [:XF :DP :XT :GT :AC :XR :XG],
+        :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20], :GT "0/1",
+                   :AC [50], :XR [2 4], :XG [10 13 15]},
+        :SAMPLE02 {:XF :exists, :DP 300, :XT [30 40], :GT "0/0/1",
+                   :AC [90], :XR [1 9], :XG [0 4 7 9]}}]
+
+      ;; -----
+
+      {:chr "3", :pos 300, :ref "T", :alt ["A", "C"],
+       :id nil, :qual 100.0, :filter [:PASS],
+       :info {}
+       :FORMAT [:XF :DP :XT:AC :XR :XG],
+       :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20],
+                  :AC [40 50], :XR [2 3 4], :XG [10 11 12 13 14 15]},
+       :SAMPLE02 {:XF :exists, :DP 20, :XT [10 20],
+                  :AC [40 50], :XR [2 3 4], :XG [10 11 12 13 14 15]}}
+      [{:chr "3", :pos 300, :ref "T", :alt ["A"],
+        :id nil, :qual 100.0, :filter [:PASS],
+        :info {}
+        :FORMAT [:XF :DP :XT:AC :XR :XG],
+        :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20],
+                   :AC [40], :XR [2 3], :XG [10 11 12]},
+        :SAMPLE02 {:XF :exists, :DP 20, :XT [10 20],
+                   :AC [40], :XR [2 3], :XG [10 11 12]}}
+       {:chr "3", :pos 300, :ref "T", :alt ["C"],
+        :id nil, :qual 100.0, :filter [:PASS],
+        :info {}
+        :FORMAT [:XF :DP :XT:AC :XR :XG],
+        :SAMPLE01 {:XF :exists, :DP 20, :XT [10 20],
+                   :AC [50], :XR [2 4], :XG [10 13 15]},
+        :SAMPLE02 {:XF :exists, :DP 20, :XT [10 20],
+                   :AC [50], :XR [2 4], :XG [10 13 15]}}])))
+
+(defn- seq-reader [s]
+  (reify protocols/ISequenceReader
+    (read-sequence [_ {:keys [^long start ^long end]} _]
+      (subs s (dec start) end))))
+
+(deftest trim-right
+  (are [?seq ?pos ?ref ?alts ?expected]
+       (= ?expected
+          (@#'norm/trim-right
+           (seq-reader ?seq)
+           {:pos ?pos, :ref ?ref :alt ?alts} 3))
+    "AATGACCGACCGACCTTGA"
+    11 "CGA" ["C"]
+    {:pos 11, :ref "CGA", :alt ["C"]}
+
+    "AATGACCGACCGACCTTGA"
+    11 "CGA" ["CA"]
+    {:pos 11, :ref "CG", :alt ["C"]}
+
+    "AATGACCGACCGACCTTGA"
+    11 "CGA" ["CA" "CGAA"]
+    {:pos 11, :ref "CG", :alt ["C" "CGA"]}
+
+    "AATGACCGACCGACCTTGA"
+    11 "CGAC" ["TCGAC" "ACGAC" "GCGAC"]
+    {:pos 8, :ref "GAC", :alt ["GACT" "GACA" "GACG"]}
+
+    "ATGCACTCCGTTGCATCCCCCTG"
+    17 "C" ["CC"]
+    {:pos 14, :ref "CAT", :alt ["CATC"]}
+
+    "ATGCACTCCGTTGCATCCCCCTG"
+    18 "C" ["CC"]
+    {:pos 15, :ref "AT", :alt ["ATC"]}
+
+    "ATGCACTCCGTTGCATCCCCCTG"
+    19 "C" ["CC"]
+    {:pos 16, :ref "T", :alt ["TC"]}
+
+    "ATGCACTCCGTTGCATCCCCCTG"
+    20 "C" ["CC"]
+    {:pos 14, :ref "CAT", :alt ["CATC"]}
+
+    "ATGCACTCCGTTGCATCCCCCTG"
+    21 "C" ["CC"]
+    {:pos 15, :ref "AT", :alt ["ATC"]}
+
+    "ATGATGATGATG"
+    11 "TG" ["TATG" "TATGATG"]
+    {:pos 8, :ref "TG", :alt ["TGAT" "TGATATG"]}
+
+    "ACCCCCCCCC"
+    9 "CC" ["C"]
+    {:pos 1, :ref "AC", :alt ["A"]}
+
+    "ACCCCCCCCC"
+    8 "CCC" ["C" "CC" "CCCC"]
+    {:pos 1, :ref "ACC", :alt ["A" "AC" "ACCC"]}
+
+    "TACCCCCCCCC"
+    9 "CCC" ["C" "CC" "CCCC"]
+    {:pos 1, :ref "TACC", :alt ["TA" "TAC" "TACCC"]}
+
+    "GTACCCCCCCCC"
+    10 "CCC" ["C" "CC" "CCCC"]
+    {:pos 1, :ref "GTACC", :alt ["GTA" "GTAC" "GTACCC"]}
+
+    "CGTACCCCCCCCC"
+    11 "CCC" ["C" "CC" "CCCC"]
+    {:pos 2, :ref "GTACC", :alt ["GTA" "GTAC" "GTACCC"]}
+
+    "TCGTACCCCCCCCC"
+    12 "CCC" ["C" "CC" "CCCC"]
+    {:pos 3, :ref "GTACC", :alt ["GTA" "GTAC" "GTACCC"]}
+
+    "CCCCC"
+    5 "C" ["CC"]
+    {:pos 1, :ref "C", :alt ["CC"]}))
+
+(deftest realign
+  (testing "not affected"
+    (are [?seq ?variant ?expected]
+         (= ?expected (norm/realign (seq-reader ?seq) ?variant))
+      "G"
+      {:pos 1, :ref "G", :alt ["T"]}
+      {:pos 1, :ref "G", :alt ["T"]}
+
+      "TGGG"
+      {:pos 1, :ref "TGGG", :alt ["TAC" "TG" "TGGGG", "AC"], :info {:END 4}}
+      {:pos 1, :ref "TGGG", :alt ["TAC" "TG" "TGGGG", "AC"], :info {:END 4}}
+
+      ""
+      {:pos 10, :ref "A", :alt nil}
+      {:pos 10, :ref "A", :alt nil}
+
+      "NNNNNNAC"
+      {:pos 7, :ref "A", :alt ["a", "AT"], :info {:END 7}}
+      {:pos 7, :ref "A", :alt ["a", "AT"], :info {:END 7}}
+
+      "CAT"
+      {:pos 1, :ref "C", :alt ["CC"]}
+      {:pos 1, :ref "C", :alt ["CC"]}))
+
+  (testing "realigned"
+    (are [?seq ?variant ?expected]
+         (= ?expected (norm/realign (seq-reader ?seq) ?variant))
+      ""
+      {:pos 10, :ref "AA", :alt nil, :info {:END 11}}
+      {:pos 10, :ref "A", :alt nil, :info {:END 10}}
+
+      "AT"
+      {:pos 1, :ref "AT", :alt nil, :info {:END 2}}
+      {:pos 1, :ref "A", :alt nil, :info {:END 1}}
+
+      "ATGC"
+      {:pos 2, :ref "TGC", :alt nil, :info {:END 4}}
+      {:pos 2, :ref "T", :alt nil, :info {:END 2}}
+
+      "AT"
+      {:pos 1, :ref "AT", :alt ["TT"], :info {:END 2}}
+      {:pos 1, :ref "A", :alt ["T"], :info {:END 1}}
+
+      "AT"
+      {:pos 1, :ref "AT", :alt ["AA"], :info {:END 2}}
+      {:pos 2, :ref "T", :alt ["A"], :info {:END 2}}
+
+      "ATTTTTTTTTTTTT"
+      {:pos 1, :ref "ATTTTTTTTTTTTT", :alt ["ATTTTTTTTTTTTTTT"], :info {:END 15}}
+      {:pos 1, :ref "A", :alt ["ATT"], :info {:END 1}}
+
+      "GATG"
+      {:pos 1, :ref "GATG", :alt ["GACT"]}
+      {:pos 3, :ref "TG", :alt ["CT"]}
+
+      "TAAACCCTAAA"
+      {:pos 1, :ref "TAAACCCTAAA", :alt ["TAA" "TAACCCTAAA"]}
+      {:pos 1, :ref "TAAACCCTA", :alt ["T" "TAACCCTA"]}
+
+      "CACA"
+      {:pos 1, :ref "CACA", :alt ["CAC"]}
+      {:pos 3, :ref "CA", :alt ["C"]}
+
+      "CAT"
+      {:pos 1, :ref "CAT", :alt ["AT"], :info {:END 3}}
+      {:pos 1, :ref "CA", :alt ["A"], :info {:END 2}}
+
+      "CAC"
+      {:pos 1, :ref "CAC", :alt ["CC", "CTC"]}
+      {:pos 1, :ref "CA", :alt ["C", "CT"]}
+
+      "GTTTTTTTTTTTTT"
+      {:pos 1, :ref "GTTTTTTTTTTTTT",
+       :alt ["GTTTTTTTTTTTTTTC", "GTTTTTTTTTTTTTT"]}
+      {:pos 14, :ref "T", :alt ["TTC", "TT"]}
+
+      "ACTCTGGGGGGGGGGGGGGGGGGGGAGGGGA"
+      {:pos 22, :ref "G", :alt ["GAG"]}
+      {:pos 21, :ref "G", :alt ["GGA"]}
+
+      "ACTCTGGGGGGGGGGGGGGGGGGGGAGGGGA"
+      {:pos 22, :ref "G", :alt ["GAG", "GCG"]}
+      {:pos 21, :ref "G", :alt ["GGA", "GGC"]}))
+
+  (testing "skipped"
+    (are [?seq ?variant ?expected]
+         (= ?expected (norm/realign (seq-reader ?seq) ?variant))
+      ""
+      {}
+      {}
+
+      "G" ;; ref
+      {:pos 1, :ref "G", :alt ["G"], :info {:END 1}}
+      {:pos 1, :ref "G", :alt ["G"], :info {:END 1}}
+
+      "ATAAAAAAAAA" ;; ref
+      {:pos 10, :ref "AA", :alt ["AA"], :info {:END 11}}
+      {:pos 10, :ref "AA", :alt ["AA"], :info {:END 11}}
+
+      "G" ;; spanning deletion
+      {:pos 1, :ref "G", :alt ["G" "*"], :info {:END 2}}
+      {:pos 1, :ref "G", :alt ["G" "*"], :info {:END 2}}
+
+      "C" ;; INS
+      {:pos 1, :ref "C", :alt ["C<ctg1>"]}
+      {:pos 1, :ref "C", :alt ["C<ctg1>"]}
+
+      "NC" ;; BND
+      {:pos 2, :ref "C", :alt ["CTT[1:123["]}
+      {:pos 2, :ref "C", :alt ["CTT[1:123["]}
+
+      "NNC" ;; Symbolic
+      {:pos 3, :ref "C", :alt ["<DUP>"]}
+      {:pos 3, :ref "C", :alt ["<DUP>"]})))
+
+(deftest realign-real-seq
+  (are [?variant ?realigned]
+       (= ?realigned
+          (with-open [r (io-seq/reader medium-twobit-file)]
+            (norm/realign r ?variant)))
+    {:chr "chr1", :pos 30658, :ref "T", :alt ["TT"]}
+    {:chr "chr1", :pos 30639, :ref "C", :alt ["CT"]}
+
+    {:chr "chr1", :pos 30640, :ref "T", :alt ["TCCT"]}
+    {:chr "chr1", :pos 30636, :ref "T", :alt ["TTCC"]}))
+
+(deftest normalize
+  (is (= [{:chr "1", :pos 2, :id nil, :ref "A", :alt ["AT"],
+           :qual nil, :filter "PASS", :info {:AC [1]},
+           :format [:GT], :SAMPLE01 {:GT "0/1"} :SAMPLE02 {:GT "0/0"}}
+          {:chr "1", :pos 15, :id nil, :ref "T", :alt ["TTC"],
+           :qual nil, :filter "PASS", :info {:AC [2]},
+           :format [:GT], :SAMPLE01 {:GT "1/0"} :SAMPLE02 {:GT "0/1"}}]
+         (sort-by
+          (juxt (comp chr/chromosome-order-key :chr) :pos)
+          (norm/normalize
+           (seq-reader "TATTTTTTTTTTTTTCTTTTT")
+           {:filter [{:id "PASS"}],
+            :info [{:id "AC", :type "Integer", :number "A"}],
+            :format [{:id "GT", :type "String", :number 1}]}
+           ["CHROM" "POS" "ID" "REF" "ALT" "QUAL" "FILTER" "INFO"
+            "FORMAT" "SAMPLE01" "SAMPLE02"]
+           [{:chr "1", :pos 15, :id nil, :ref "TC", :alt ["TTCC" "TTC"],
+             :qual nil, :filter "PASS", :info {:AC [2 1]},
+             :format [:GT], :SAMPLE01 {:GT "1/2"} :SAMPLE02 {:GT "0/1"}}])))))