Skip to content

Add support for generating CSI file from bgzip compressed VCF. - #190

Merged
athos merged 1 commit into
masterfrom
feature/csi-writer-for-vcf
Apr 17, 2020
Merged

athos merged 1 commit into
masterfrom
feature/csi-writer-for-vcf

Conversation

@niyarin

@niyarin niyarin commented Feb 26, 2020

Copy link
Copy Markdown
Contributor

Add support for generating CSI file from VCF.

@niyarin
niyarin requested a review from alumi as a code owner February 26, 2020 11:38
@niyarin
niyarin requested a review from a team February 26, 2020 11:38
@alumi
alumi requested review from a team and athos and removed request for a team February 27, 2020 09:20
Comment thread src/cljam/algo/vcf_indexer.clj Outdated
Comment thread src/cljam/io/vcf/reader.clj Outdated
Comment thread test/cljam/io/util/bin_test.clj Outdated
Comment thread src/cljam/io/util/bin.clj Outdated
Comment thread src/cljam/io/vcf.clj Outdated
Comment thread src/cljam/io/vcf.clj Outdated
Comment thread src/cljam/io/vcf/reader.clj Outdated
Comment thread src/cljam/io/vcf/reader.clj Outdated
Comment thread src/cljam/io/vcf.clj Outdated
Comment thread src/cljam/io/util/bin.clj Outdated

@athos athos left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for working on another feature! I just added some comments on relatively superficial matters.

Comment thread src/cljam/io/util/bin.clj Outdated
Comment thread src/cljam/io/util/bin.clj Outdated
Comment thread src/cljam/io/vcf.clj Outdated
Comment thread src/cljam/algo/vcf_indexer.clj Outdated
@codecov

codecov Bot commented Mar 2, 2020 •

Copy link
Copy Markdown

Codecov Report

Merging #190 into master will increase coverage by 0.14%.
The diff coverage is 92.22%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master     #190      +/-   ##
==========================================
+ Coverage   86.67%   86.81%   +0.14%     
==========================================
  Files          76       77       +1     
  Lines        6057     6205     +148     
  Branches      502      516      +14     
==========================================
+ Hits         5250     5387     +137     
+ Misses        305      302       -3     
- Partials      502      516      +14     
Impacted Files Coverage Δ
src/cljam/io/protocols.clj 100.00% <ø> (ø)
src/cljam/io/csi.clj 91.46% <89.06%> (-8.54%) ⬇️
src/cljam/io/util/bin.clj 94.33% <93.33%> (-3.17%) ⬇️
src/cljam/algo/vcf_indexer.clj 100.00% <100.00%> (ø)
src/cljam/io/bam_index/writer.clj 84.37% <100.00%> (+0.08%) ⬆️
src/cljam/io/sam/util.clj 92.30% <100.00%> (+11.22%) ⬆️
src/cljam/io/vcf.clj 96.00% <100.00%> (+0.16%) ⬆️
src/cljam/io/vcf/reader.clj 88.75% <100.00%> (+2.08%) ⬆️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 6d961d7...2444372. Read the comment docs.

@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch from 9ae9ce7 to 803b2d5 Compare March 2, 2020 06:59
Comment thread src/cljam/io/vcf.clj Outdated
Comment thread src/cljam/io/util/bin.clj Outdated
Comment thread test/cljam/io/util/bin_test.clj Outdated
Comment thread test/cljam/algo/vcf_indexer_test.clj Outdated
@niyarin

niyarin commented Mar 2, 2020

Copy link
Copy Markdown
Contributor Author

Thank you for reviewing.
I fixed bad expressions and added some test cases about vcf not ascending order chroms.

@alumi alumi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the updates! I added a few more comments.

Comment thread src/cljam/io/csi.clj Outdated
Comment thread src/cljam/algo/vcf_indexer.clj Outdated
Comment thread src/cljam/io/csi.clj Outdated
Comment thread src/cljam/io/vcf/reader.clj Outdated
Comment thread src/cljam/io/csi.clj Outdated
Comment thread test/cljam/algo/vcf_indexer_test.clj
@niyarin

niyarin commented Mar 16, 2020

Copy link
Copy Markdown
Contributor Author

Sorry for late, I fixed.

@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch from cf1ce51 to cddf8c2 Compare March 16, 2020 01:26

@alumi alumi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks like lein test :all is failing 🔥

@niyarin

niyarin commented Mar 17, 2020

Copy link
Copy Markdown
Contributor Author

I'm so sorry,I will fix test cases same as others.

Comment thread src/cljam/io/csi.clj Outdated
Comment thread test/cljam/algo/vcf_indexer_test.clj Outdated
Comment thread test/cljam/io/csi_test.clj Outdated

@athos athos left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Added some more trivial comments.

Comment thread src/cljam/io/csi.clj Outdated
[bin (->> offsets
(map #(chunk/->Chunk (:file-beg %) (:file-end %)))
reverse)]))
sort

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
sort
Comment thread src/cljam/io/csi.clj Outdated
[offsets shift depth]
(let [chr-offsets (->> (partition-by :chr offsets)
(map #(vector (:chr (first %)) %))
sort

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
sort
Comment thread src/cljam/io/csi.clj Outdated
from offsets {:file-beg :file-end :beg :end :chr }"
[offsets shift depth]
(let [chr-offsets (->> (partition-by :chr offsets)
(map #(vector (:chr (first %)) %))

@athos athos Mar 18, 2020 •

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These two lines seem to be equivalent to (group-by :chr offsets) in this case.

Comment thread src/cljam/algo/vcf_indexer.clj Outdated

(defn create-index
"Creates a CSI index file from the VCF file."
[in-bgziped-vcf out-csi {:keys [shift depth] :or {shift 14 depth 5}}]

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

s/bgziped/bgzipped/

@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch 3 times, most recently from d55f8f2 to 79451d1 Compare March 23, 2020 07:18
@niyarin

niyarin commented Mar 23, 2020

Copy link
Copy Markdown
Contributor Author

Thank you for pointing.

I made the following changes.

  • Loffset in CSI decided not to match the implementation of hts_lib.

    • Specifically, the loffset of the key that does not exist during the beg-end of vcf came closer to the target data.
  • Added an operation to move smaller chunks to higher binning in bidx calculation,same as hts_lib.

@alumi alumi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Run formatter and linters every time before you push your commits
  • It seems like generated CSI files are not compatible with bcftools
    • Please make sure that both cljam.io.vcf.read-variants-randomly and bcftools view works as well
  • Add a test case for VCF files with skipped contigs as in #190 (comment)
@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch from 79451d1 to 3bf6899 Compare March 24, 2020 09:47
@niyarin

niyarin commented Mar 24, 2020

Copy link
Copy Markdown
Contributor Author

Sorry, I applied linter.

@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch from e655e42 to 4fd65bd Compare March 25, 2020 02:29
@niyarin

niyarin commented Mar 25, 2020

Copy link
Copy Markdown
Contributor Author

Bcftools (hts_lib) will generate an error if some tabix data (such as format,col_seq , names and so on) is not inserted into CSI aux.
So I fixed CSI for portability.

@alumi alumi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks! The indexing looks good. Added some more trivial comments.

Comment thread src/cljam/io/csi.clj Outdated
Comment thread src/cljam/io/csi.clj Outdated
Comment thread src/cljam/io/vcf/reader.clj Outdated
Comment thread test/cljam/algo/vcf_indexer_test.clj Outdated
Comment thread src/cljam/io/csi.clj Outdated
Comment thread src/cljam/io/csi.clj
Comment thread src/cljam/io/csi.clj Outdated
Comment thread src/cljam/io/util/bin.clj
@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch 2 times, most recently from 5d3a3aa to 985fed3 Compare April 6, 2020 02:39
@niyarin

niyarin commented Apr 6, 2020

Copy link
Copy Markdown
Contributor Author

Thank you for pointing.I fixed them.
I think that arguments of reg->bin is [beg,end) (one-based, half-close-half-open).

@alumi alumi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for the changes. I've added some suggestions around reg->bin. Would you check them, please?

Comment thread src/cljam/io/sam/util.clj Outdated
@@ -68,4 +48,4 @@
[aln]
(let [beg (dec (:pos aln))

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
(let [beg (dec (:pos aln))
(let [beg (:pos aln)
Comment thread src/cljam/io/bam_index/writer.clj Outdated
Comment on lines +76 to +77
(let [bin (util-bin/reg->bin (.pos aln) (inc (.end aln))
linear-index-shift linear-index-depth)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
(let [bin (util-bin/reg->bin (.pos aln) (inc (.end aln))
linear-index-shift linear-index-depth)
(let [bin (util-bin/reg->bin
(.pos aln) (.end aln) linear-index-shift linear-index-depth)
Comment thread src/cljam/io/util/bin.clj Outdated
Comment on lines +74 to +75
beg (Math/min max-pos (Math/max 0 (dec beg)))
end (Math/min max-pos (Math/max 0 (dec end)))]

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
beg (Math/min max-pos (Math/max 0 (dec beg)))
end (Math/min max-pos (Math/max 0 (dec end)))]
beg (dec (Math/max 1 (Math/min max-pos beg)))
end (dec (Math/max 1 (Math/min max-pos end)))]
Comment thread test/cljam/io/util/bin_test.clj Outdated
Comment on lines +61 to +66
(deftest reg->bin-test
(is (= (util-bin/reg->bin 1 1 14 6) 37449))
(is (= (util-bin/reg->bin 1 1 14 7) 299593))
(is (= (util-bin/reg->bin 1 32769 15 6) 4681))
(is (= (util-bin/reg->bin 1 2 0 2) 1))
(is (= (util-bin/reg->bin 240877561 240877568 14 6) 52150)))

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
(deftest reg->bin-test
(is (= (util-bin/reg->bin 1 1 14 6) 37449))
(is (= (util-bin/reg->bin 1 1 14 7) 299593))
(is (= (util-bin/reg->bin 1 32769 15 6) 4681))
(is (= (util-bin/reg->bin 1 2 0 2) 1))
(is (= (util-bin/reg->bin 240877561 240877568 14 6) 52150)))
(deftest reg->bin-test
(testing "BAI-compatible"
(are [?start ?end ?bin]
(= ?bin (util-bin/reg->bin ?start ?end 14 5))
0 0 4681
0 1 4681
1 1 4681
1 16384 4681
16384 16384 4681
16384 16385 585
16385 16385 4682
131072 131072 4688
131072 131073 73
8388608 8388609 1
67108864 67108865 0
536870912 536870912 37448
536870912 536870913 37448))
(testing "1-by-1"
(are [?start ?end ?bin]
(= ?bin (util-bin/reg->bin ?start ?end 0 2))
0 1 9
1 1 9
2 2 10
1 2 1
1 8 1
1 9 0
8 9 0
9 9 17
63 64 8
64 64 72
64 64 72))
(testing "various min-shfits and depths"
(is (= (util-bin/reg->bin 1 1 14 6) 37449))
(is (= (util-bin/reg->bin 1 1 14 7) 299593))
(is (= (util-bin/reg->bin 1 32769 15 6) 4681))
(is (= (util-bin/reg->bin 240877561 240877568 14 6) 52150))))
@niyarin

niyarin commented Apr 10, 2020

Copy link
Copy Markdown
Contributor Author

Thank you,I fixed.

@alumi alumi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM 👍

I found a bug in the random reading of VCF files. Could you fix it in another PR?

;; with test-resources/vcf/test-chr-skipped.vcf.gz.csi
(with-open [r (vcf/reader "test-resources/vcf/test-chr-skipped.vcf.gz")]
        (doall (vcf/read-variants-randomly r {:chr "chr3"} {})))
Execution error (NullPointerException) at cljam.io.csi.CSI/get_chunks (csi.clj:19).
@alumi

alumi commented Apr 14, 2020

Copy link
Copy Markdown
Member

@niyarin Please don't forget to squash the commits before merging this branch 😀

@athos athos left a comment •

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@niyarin Thanks for your persistent effort to improve the code for this feature! 🙏Now the implementation works fine 👍

Just one more thing, can you add a mention to CSI, as well as Tabix, in the docstring of cljam.io.vcf/read-variants-randomly as another supported index format? The current description looks to me like it only supports the Tabix format. Also, you may want to avoid using the name tabix-data in the body of that function.

Once you make a commit for the above point and squash all the commits into one (or a few), I will merge the patch soon.

@niyarin
niyarin force-pushed the feature/csi-writer-for-vcf branch from 2c71be5 to 2444372 Compare April 16, 2020 10:32
@niyarin

niyarin commented Apr 16, 2020 •

Copy link
Copy Markdown
Contributor Author

Thank you,I added a comment that VCF can also use CSI and renamed 'tabix-data' to 'index-data'.

@athos
athos merged commit fa2fa54 into master Apr 17, 2020
@athos
athos deleted the feature/csi-writer-for-vcf branch April 17, 2020 02:01
@alumi

alumi commented Apr 17, 2020

Copy link
Copy Markdown
Member

🎉 🎉 🎉

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

3 participants