1

I don't know what I have in my huge vcf.gz file that looks like this

CHROM       POS         ALT     12345       
1           345632      T       0/1:4,4:8:99:105,0,106
4           032184      C       1/1:46,9:55:99:99,0,1222
6           843290      A       0/1:67,20:87:99:336,0,1641
7           743290      C       0/1:37,20:57:99:336,0,2641
8           329283      T       0/2:99:21:253,0,290:11,10
9           789320      C       2/2:99:21:253,0,290:11,10

And I would like to extract all unique values in the 4th column before the ":". That is in this case:

0/1
1/1
0/2
2/2

Do you have any suggestions?

  • I don't understand the delimiters. Can't be colons otherwise some lines only have 1 value. – user13757 Nov 24 '16 at 13:25
3

With awk, checking that the 4th column is in the right format:

awk 'match($4, /^[0-9]+\/[0-9]+:/) {
       c = substr($4, RSTART, RLENGTH-1)
       if (!seen[c]++) print c
     }'
  • Nice! I'm tempted to counter-suggest $4 ~ /:/ { c = $4 ; sub(/:.*$/, "", c) ; if (!seen[c]) print c }', but I must admit, using using RSTART` and RLENGTH (assigned by match()) is a clever technique. – Nominal Animal Nov 24 '16 at 14:02
  • @NominalAnimal, yes match()+RSTART+RLENGTH is awk's way to get was it being matched by a regexp since it doesn't have capture groups like in perl/sed (except via extensions like gawk's gensub() or RT...). – Stéphane Chazelas Nov 24 '16 at 14:18
0

Let's assume that lines that begin with a letter are comments or ignored.

zcat vcf.gz | awk \
   'BEGIN {
        RS = "[\t\v\f ]*(\r\n|\n\r|\r|\n)" ;
        FS = "[\t\v\f ]+"
    }

    /^[A-Za-z]/ {
        next
    }

    NF >= 4 {
        key = $4 ;
        sub(/:.*$/, "", key) ;
        seen[key]++
    }

    END {
        for (key in seen)
            printf "%s\n", key
    }'

You can write the entire command on one line (remove the \ at end of the first line, though), because I added all the necessary semicolons to do so.

The BEGIN rule sets up universal newline support. All whitespace at end of lines is ignored, and any newline form (CR, LF, CRLF, LFCR) is accepted as a newline. Any number of tabs or spaces are treated as field separators.

The rule /^[A-Za-z]/ applies to all lines that begin with a letter. The next causes them to be ignored.

The next rule applies to all records (lines) with at least four fields. The fourth field is copied to variable key, then everything after the first colon (including the colon) is removed. We use the resulting value as a key to the seen associative array. The value we assign does not really matter, but here, seen[key] will contain the number of times each value was referred to (1 or more).

The END rule is run after all input has been processed. Here, the loop iterates over the keys of seen[] array (in unspecified order), and just prints the keys.

If you wanted to keep the order in the data, or use some specific order for the keys, the above snippet needs slight modifications.

  • Thank you so much. However i did not get the script to work. Im probably missing something. :/ Was I supposed to be abel to copy the whole thing except the \ in the first row after "awk"? Sorry. – Sanna Nov 24 '16 at 15:13
  • @Sanna: You can copy the entire command, as-is, to a Bash command prompt, no changes needed. If you copy the entire command, then remove the newlines (making it all one very long line), then you need to remove the ` also. For me, the command definitely works (when copy-pasted from above to a Bash command prompt); Bash just inserts some >` characters at the beginning of the lines, as usual, to let you know you are inputting a multi-line command. – Nominal Animal Nov 24 '16 at 18:54
0

If perl is okay:

$ perl -lane '($k) = $F[3] =~ m/^([^:]+)/; print $k if !$seen{$k}++ && $. > 1' ip.txt 
0/1
1/1
0/2
2/2
  • ($k) = $F[3] =~ m/^([^:]+)/ get string before : from 4th column
  • Then print if it is unique and not first line (to avoid header)


If the pattern before : has to be necessarily of pattern digits/digits, then change the match to m|^(\d+/\d+):|

0

As a one-liner

gzip -dc input.gz | grep -E "^[0-9]" | cut -f 1 -d : | sort -k 4 -u | awk ' { print $4 } '

Explanation:

* gzip -dc             # decompress and write to standard output   
* grep -E "^[0-9]"     # only look at lines starting with a digit  
* cut -f 1 -d:         # remove the text following the colon  
* sort -k4 -u          # sort on the forth field - uniquely  
* awk ' { print $4 } ' # print the fourth field (that we sorted on)
0

I would use:

grep -o -P '.{1}/.{1}' file

In your case:

zcat vcf.gz | grep -o -P '.{1}/.{1}' 

Edit: for only one ocurrence, add uniq at the end of the line:

zcat vcf.gz | grep -o -P './.' | uniq

or if you want it in order:

zcat vcf.gz | grep -o -P './.' | sort -u

  • .{1} is the same as . – Stéphane Chazelas Nov 24 '16 at 13:45
  • I tried this one. Gave me the answer that I have why many more then 1/1, 0/1 and 1/1 but not in a list only listing unique values ONE time. Never-ending list of all the lines. :) – Sanna Nov 24 '16 at 15:10
  • If you want only one ocurrence for each value, you should put " | uniq ", it would be something like: – Santiago Flores Nov 24 '16 at 18:05
  • Sry, hit enter :| zcat vcf.gz | grep -o -P './.' | uniq (thanks Stephane ;) ) – Santiago Flores Nov 24 '16 at 18:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.