awk working with large files

Question

I have two very large vcf files 2GB and 6GB I want to look for unique combinations of CHROM and POS and output the row that matches. However, because the files a so large my machine always hang and stop processing. Is there are a way to work around the problem of this large file?

I am using this command, taken from an answer to another question:

awk  '{ if(NR==FNR){a[$1$2]=$0}else{if($1$2 in a){print}}}' file1.vcf file2.vcf

file 1

#CHROM  POS ID  REF ALT QUAL
1   10366   rs58108140  G   A   
1   10611   rs189107123 C   G   
1   51954   rs185832753 G   C   
1   13327   rs144762171 G   C

file 2

#CHROM  POS ID  REF ALT QUAL    
1   10366   rs58108140  G   A   
1   51935   rs181754315 C   T   
1   51954   rs185832753 G   C   
1   52058   rs62637813  G   C   
1   52144   rs190291950 T   A

output

1 10366 rs58108140  G   A
1 51954  rs185832753 G   C 

Answer 1

Maybe comm, which is part of GNU coreutils, is more efficient.

comm writes to standard output lines that are common, and lines that are unique, to two input files; a file name of ‘-’ means standard input. Synopsis:

comm [option]… file1 file2

In your case you could sort both files based on the first 3 columns:

file1=file1.vcf
file2=file2.vcf

comm -12 <(awk '{print $1,$2,$3}' $file1 | sort) \
         <(awk '{print $1,$2,$3}' $file2 | sort) > matches


# result:
cat matches
1 10366 rs58108140
1 51954 rs185832753
#CHROM POS ID

This could be of course polished by removing the header. Also if you want the original lines as output, you can grep the matches to the original files:

grep -f matches file1
#CHROM POS ID REF ALT QUAL 
1 10366 rs58108140 G A
1 51954 rs185832753 G C

Edit: In the first version, intermediate temp files were used, but as suggested by terdon in the comments, this can be avoided using Process Substitution allowing the <(..) contstruct.

Answer 2

Sort both files by chromosome and position and read each file one line by line.

Either you iterate to the next line in one file, or go the next in the other, or there's a match that you print out and then iterate both lines.

This is not necessarily an awk approach, but is doable with Perl, Python, or any other scripting language.

In pseudocode:

Sub getCommonLines()
  Sort FileA
  Sort FileB
  Get lineA from FileA
  Get lineB from FileB
  While not EOF in FileA or EOF in FileB:
    If lineA > lineB:
      Get lineB from FileB
    Elseif lineB > lineA:
      Get lineA from FileA
    Else:
      Print lineA or lineB
      Get lineA from FileA
      Get lineB from FileB
End Sub

The use of the > operator is based on lexicographical and numerical ordering provided from sorting files A and B, and comparison should be done only on the fields of interest. If you compare other fields (esp. from VCF files) you may get other fields that give you the wrong or unexpected comparison result.

The term EOF means end-of-file, which means just what it says — you're at the end of the file and there are no more lines to compare between files, so your subroutine is done.

Sorting has a time cost, but memory overhead for modern sorting algorithms is very minimal. You can sort via process substitution, but sorting to an intermediate file can be useful if you have other things you want to do with them, where sorted input is valuable.

This would be a good approach for very large files that you do not edit or update often (or at all), for which a hash table (or "dictionary", or "associative array", all the same terms for the same thing) is impractical because of its large memory overhead.

You can also use a tool like comm with sorted files, as described in another answer/comment. This does basically what my answer describes without you having to understand or implement things. You must remember to provide sorted input, however, as with other Unix tools that do line comparison (uniq, etc.).

Despite this answer being a really long "comment", it is useful to have a basic understanding of how things work under the hood, so that when you use these tools, you know what's going on and why you need to provide inputs in a certain way. It's easy to get a million one-liners off of Stack Exchange sites, but understanding is valuable beyond gold.

Answer 3

Use bedtools intersect. You can select the various combinations of outputs:

https://bedtools.readthedocs.io/en/latest/content/tools/intersect.html