1

I have a tab delimited file with information in it about genetic material. Some information is cut into a smaller tab file with some of the columns extracted, and uniq is used to ensure there are no duplicates. A count is stored which is important later down the pipeline. I want to change the basis of the uniq function from simply being based on the sequence in one field, to a regex from within another field. The field I want to extract from is in this format:

14:50065421-50065521:12397472_t14_w100_x1

but the bit after the 2nd colon changes depending on the file input. I want to use uniq on the basis of the first half, so:

14:50065421-50065521

I've tested the regex '((^[0-9]{0,2}|x|y|MT):[0-9]{0,9}-[0-9]{0,9}:)' and it works on a small sample of the data. I've found some ways using grep and a perl script which both work to remove rows based on the regex but none of them provide the count (which is why uniq is much more ideal). Is there any way to use a regex with uniq? Or is there a better way which will also store the count of those removed?

The current code is this:

cat ${TAB_FILE} | \
    sed -e '1,2d' | \
    cut -f3,4 | \
    sort -k1 -u | \
    sort -k2 | \
    uniq -cf1| \
    sort -rn > t1

cat ${TAB_FILE} | \
    sed -e '1,2d' | \
    awk {'print $3"\t "$6"\t" $7"\t "$4'} | \
    sort -k4 > t2

awk 'FNR==NR {C[$2]=$1;next}FNR==1 \
    {print "Count Chromosome:Positions:QNAME Sequence Exon Transcript_ID"; next}$1 in C \
    {print C[$1], $1, $4, $3, $2}' t1 t2 > t3

cat t3 | awk '{print "//NODECLASS\t\"" $2"_"$1  "\"\t\"Exon " $4 "\"\t\"" $5 "\""}'

From the first step it is column 1 of the cut that I would like to base my regex on, instead of using column 2. Any help would be greatly appreciated, please feel free to ask if you need me to clarify anything.

Example of tab file:

queryHits   subjectHits readname    readSeq geneid  transcriptid    exonnumber  genename    biotype
350851  1   14:50065421-50065521:12397472_t14_w100_x1   CGCTGCCAGCTGCGCGCTCGGGGGAAAAGACGTTGCGCCCCCGCCGACTGCCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCA    ENSG00000165501 ENST00000298288 1   LRR1    protein_coding
350851  5   14:50065421-50065521:12397472_t14_w100_x1   CGCTGCCAGCTGCGCGCTCGGGGGAAAAGACGTTGCGCCCCCGCCGACTGCCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCA    ENSG00000165501 ENST00000318317 1   LRR1    protein_coding
350851  8   14:50065421-50065521:12397472_t14_w100_x1   CGCTGCCAGCTGCGCGCTCGGGGGAAAAGACGTTGCGCCCCCGCCGACTGCCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCA    ENSG00000165501 ENST00000554869 1   LRR1    protein_coding
350852  1   14:50065461-50065561:12655987_t14_w100_x1   CCGCCGACTGCCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCAACCACGGCGAGAACCCGGCCCTGCTAAGGGAGAAGGGAAG    ENSG00000165501 ENST00000298288 1   LRR1    protein_coding
350852  5   14:50065461-50065561:12655987_t14_w100_x1   CCGCCGACTGCCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCAACCACGGCGAGAACCCGGCCCTGCTAAGGGAGAAGGGAAG    ENSG00000165501 ENST00000318317 1   LRR1    protein_coding
350852  8   14:50065461-50065561:12655987_t14_w100_x1   CCGCCGACTGCCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCAACCACGGCGAGAACCCGGCCCTGCTAAGGGAGAAGGGAAG    ENSG00000165501 ENST00000554869 1   LRR1    protein_coding
350853  1   14:50065471-50065571:22679947_t13_w100_x1   CCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCAACCACGGCGAGAACCCGGCCCTGCTAAGGGAGAAGGGAAGCCGTTTCCCG    ENSG00000165501 ENST00000298288 1   LRR1    protein_coding
350853  5   14:50065471-50065571:22679947_t13_w100_x1   CCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCAACCACGGCGAGAACCCGGCCCTGCTAAGGGAGAAGGGAAGCCGTTTCCCG    ENSG00000165501 ENST00000318317 1   LRR1    protein_coding
350853  8   14:50065471-50065571:22679947_t13_w100_x1   CCGGTTTCCCGGGCGCGAGCCCGGATCCAGGTGGTCAGTCCCGGTACGCAACCACGGCGAGAACCCGGCCCTGCTAAGGGAGAAGGGAAGCCGTTTCCCG    ENSG00000165501 ENST00000554869 1   LRR1    protein_coding

Example of output file where uniq hasn't been used, I want to be able to remove the duplicates such as those on lines 4, 5, and 6, or 8, 9, and 10:

//NODECLASS "Chromosome:Positions:QNAME_Count"  "Exon Exon" "Transcript_ID"

//NODECLASS "14:50067283-50067383:20149917_t14_w100_x1_1"   "Exon 1"    "ENST00000557531"
//NODECLASS "14:50067284-50067366:14257122_t14_w100_x1_2"   "Exon 1"    "ENST00000557531"
//NODECLASS "14:50067285-50067385:2072777_t12_w100_x1_1"    "Exon 1"    "ENST00000557531"
//NODECLASS "14:50074262-50074362:4355312_t12_w100_x1_1"    "Exon 3"    "ENST00000298288"
//NODECLASS "14:50074262-50074362:4355312_t12_w100_x1_1"    "Exon 4"    "ENST00000540712"
//NODECLASS "14:50074262-50074362:4355312_t12_w100_x1_1"    "Exon 4"    "ENST00000554869"
//NODECLASS "14:50067286-50067386:15839225_t12_w100_x1_3"   "Exon 1"    "ENST00000557531"
//NODECLASS "14:50074263-50074363:8914169_t11_w100_x1_1"    "Exon 3"    "ENST00000298288"
//NODECLASS "14:50074263-50074363:8914169_t11_w100_x1_1"    "Exon 4"    "ENST00000540712"
//NODECLASS "14:50074263-50074363:8914169_t11_w100_x1_1"    "Exon 4"    "ENST00000554869"
//NODECLASS "14:50067287-50067387:5439923_t13_w100_x1_1"    "Exon 1"    "ENST00000557531"
//NODECLASS "14:50067287-50067387:14106336_t12_w100_x1_3"   "Exon 1"    "ENST00000557531"
//NODECLASS "14:50074404-50074504:15492363_t11_w100_x1_1"   "Exon 3"    "ENST00000298288"
//NODECLASS "14:50074404-50074504:15492363_t11_w100_x1_1"   "Exon 4"    "ENST00000540712"
//NODECLASS "14:50074404-50074504:15492363_t11_w100_x1_1"   "Exon 4"    "ENST00000554869"
//NODECLASS "14:50074135-50074235:11346262_t11_w100_x1_2"   "Exon 3"    "ENST00000298288"
  • 2
    For people (like me) who aren't in the field, a representative sample of the data would go a long ways. – Jeff Schaller Aug 2 '17 at 16:41
  • 1
    @RachelMiller please share unchanged first input data and expected result, it will help better understanding what you want from which data, thank you – αғsнιη Aug 2 '17 at 16:56
  • Apologies, I've updated with more information now. – Rachel Miller Aug 3 '17 at 10:17
  • @RachelMiller Can you explain how you want the input to be transformed to final output, without these intermediate details? – user218374 Aug 3 '17 at 12:12
  • Basically I want to use uniq to remove any rows which have the same start and end positions (the numbers after the first colon). I found a temporary solution this morning which is to extract these coordinates and then sort | uniq -c on this column instead so only unique rows are printed into the second tab file, and a count of how many were removed is stored as well. This has increased the run time, I think because the writing of the first tab file requires a whole extra column to be extracted/written, so it would be more efficient just to base the uniq on the regex I mentioned. – Rachel Miller Aug 3 '17 at 13:56
0

I am not sure what you want exactly, since you seem overly fixated on your sort/uniq stuff.

However if all you wanted is the pruning of repeating lines, e.g., 4,5,6 and 8,9,10 then you could do this on your last file you are showing which has //NODECLASS lines:

perl -F\" -lane '
   print,next if $. < 3;
   print if ! $h{($F[1] =~ /:(.*?):/)[0]}++;
' NODE_CLASS_file

where the first 2 lines are skipped printed. Whilst for the rest, we look at the 2nd field the number between the : /:(.*?):/ will provide you the number then you need to get into the scalar context by enrobing it in (...)[0] and pass that on as the key to the hash %h. Print the current line only if this key hasn't already been seen.

0

I've managed to find a crude solution whereby I add the start-end positions into the tab file as a separate column. Then when creating the Nodeclass file I extract this column as well, and use sort | uniq -c based on this instead of the sequence column. It seems to be doing what I wanted but it's slightly slower because the script creating the tab file needs to produce a whole extra column! Thanks :)

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.