6

Due to tabix constraints, I need to sort a vcf so that contigs and then positions within contigs occur in numerical order in the vcf. I don't know if the following will sort positions within contigs, but I tried this:

cat in.vcf | vcf-sort -c > out.vcf

but was informed:

Old version of sort command installed, please run without the -c option.

I've read that can cause problems here. Further, my contigs are't numeric, but are something like this: 'Erumb1_s00000030'.

I see two possible solutions. (1) I am on Mac OS Mojave, and if there were a way to install the a new unix-style sort, I could use it. I tried homebrew and didn't find anything. (2) there are likely some bash commands that would do it. If anyone has an idea in either of these directions, I would appreciate it.

Peter Pearman
  • 183
  • 1
  • 6
  • It might be worth posting a question about these tabix constraints. We could be able to help you fix that which would allow you to use bcftools. – terdon Jan 16 '19 at 19:25

2 Answers2

5

You could simply use bcftools sort for it:

$ bcftools sort input.vcf > output.vcf

If you really want to use bash only, you can do this:

$ grep "^#" input.vcf > output.vcf
$ grep -v "^#" input.vcf| sort -k1,1V -k2,2g >> output.vcf

The first command will write the header information to the new vcf file. The second will sort by contig name and position and append the result to the new vcf file.

fin swimmer

finswimmer
  • 1,342
  • 7
  • 15
  • I get this: Writing to /tmp/bcftools-sort.7BfymN [W::vcf_parse] Contig 'Erumb1_s00000030' is not defined in the header. (Quick workaround: index the file with tabix.). So, this is where I started, having to index the .vcf. Looks like the grep approach worked, but I've got to give it a closer look. – Peter Pearman Jan 16 '19 at 17:53
  • 1
    @PeterPearman first check that your sort implementation deals with the -k flag correctly. The command in this answer should work fine for GNU sort (the default on Linux) but I'm not sure it works correctly with BSD sort (default on macOS). The man page I found seems to suggest that the V and g modifiers aren't supported. However, you don't really need the V and you can use n instead of g: sort -k1,1 -k2,2n. – terdon Jan 16 '19 at 18:25
  • 1
    bcftools is very strict about the correct vcf structure. It heavily relies on the header data. This should include contig information and any key used in the INFO and FORMATcolumn. You should always have this information there. bcftools is a very powerful tool. So be strict to your vcf file, that you can use it. – finswimmer Jan 16 '19 at 19:23
  • @terdon is right. V isn't needed. I'm used to it, handling human chromosome names. Doing it that way chr2 is located before chr10. But in vcf file the only conditions are, that the contig names build blocks and inside this blocks the position is ascending. – finswimmer Jan 16 '19 at 19:25
2

Use a subshell. First print the header and then the variant lines, sorting by chromosome and position:

(grep ^"#" input.vcf; grep -v ^"#" input.vcf | sort -k1,1 -k2,2n) > sorted.vcf
terdon
  • 10,071
  • 5
  • 22
  • 48
atongsa
  • 131
  • 4