split 1 column input into 5 column bed file

Question

In f1 below I am trying to split $1 based on each line and create a bed file. If the line is a snp then the : is spilt and the text is $1 the last digit is $2 -1 and $3 and the letter to the left of the > is $4 and the letter to the right is $5. (lines 1 and 2)

If the line has dup then the : is spilt and the text is $1 the last digit is $2 and $3 and the $4 is - and the letter to the right is $5. (line 4)

If the line has del with nothing else then the : is spilt and the text is $1 the last digit is $2 and $3 and the letter to the left of the > is $4 and $5 is -. (line 5)

I am not sure how to format for line 3 (the complex indel)... maybe split $1 on the : and the text is $1 the last digit is $2 and $3 and the letters to the left of the > is $4 and the letters to the right is $5. Thank you :)

f1

chr7:140453145A>T
chr7:140453136A>T 
chr7:140453135_140453136delCAinsTT
chr20:31022287dupA
chr19:13054614delG

desired

chr7    140453144   140453145   A   T 
chr7    140453135   140453136   A   T 
chr7    140453135   140453136   TT  CA
chr20   31022287    31022287    -   A
chr19   13054614    13054614    G   -

awk

awk 'BEGIN{OFS="\t"} {sub(/[^0-9]+$/, "", $1)}' f1

awk 'BEGIN { FS = "[ -]"; OFS="\t" } NF==3 { print $1, $2 - 1, $2, $3 } NF==4 { print $1, $2, $3, $4 }' f1

Answer 1

The following Python script seems to do the job.

#!/usr/bin/env python3

from __future__ import print_function
import re
import sys

for line in sys.stdin:
    if '_' in line:
        match = re.search(r'(\S+):(\d+)_(\d+)del([ACGT]+)ins([ACGT]+)', line)
        print(match.group(1), match.group(2), match.group(3), match.group(5), match.group(4), sep='\t')
    elif 'dup' in line or 'del' in line:
        match = re.search(r'(\S+):(\d+)(del|dup)([ACGT]+)', line)
        if match.group(3) == 'del':
            print(match.group(1), match.group(2), match.group(2), match.group(4), '-', sep='\t')
        else:
            print(match.group(1), match.group(2), match.group(2), '-', match.group(4), sep='\t')
    else:
        match = re.search(r'(\S+):(\d+)([ACGT]+)>([ACGT]+)', line)
        coord = int(match.group(2))
        print(match.group(1), coord - 1, coord, match.group(3), match.group(4), sep='\t')

Invoked like so on the command line.

[standage@lappy ~] $ ./transform < f1
chr7 140453144 140453145 A T
chr7 140453135 140453136 A T
chr7    140453135   140453136   TT  CA
chr20   31022287    31022287    -   A
chr19   13054614    13054614    G   -

Answer 2

The Python script in my first response is heavy on regex matching, which is pretty clunky in Python. I like Python much better than Perl overall, but a throwaway script like this will be clearer and more concise in Perl.

#!/usr/bin/env perl
use strict;

while(<STDIN>) {
    if (m/(\S+):(\d+)_(\d+)del([ACGT]+)ins([ACGT]+)/) {
        print("$1\t$2\t$3\t$5\t$4\n");
    }
    elsif (m/(\S+):(\d+)(del|dup)([ACGT]+)/) {
        if ($3 eq "del") {
            print("$1\t$2\t$2\t$4\t-\n");
        }
        else {
            print("$1\t$2\t$2\t-\t$4\n");
        }
    }
    elsif (m/(\S+):(\d+)([ACGT]+)>([ACGT]+)/) {
        my $coord2 = $2;
        my $coord1 = $2 - 1;
        print("$1\t$coord1\t$coord2\t$3\t$4\n")
    }
}

It's been several years since I wrote Perl on a regular basis, so there are probably ways to make the script even more clear and concise. (Although I have a love/hate relationship with the default/hidden variables, which can make a script more concise but also harder for a newcomer to understand.)

Invoked like so on the command line.

[standage@lappy ~] $ ./transform < f1
chr7 140453144 140453145 A T
chr7 140453135 140453136 A T
chr7    140453135   140453136   TT  CA
chr20   31022287    31022287    -   A
chr19   13054614    13054614    G   -

Answer 3

Quick Perl solution (following on from Daniel's answer saying that this would be "clearer and more concise in Perl"). I was a bit confused with your BED output because you've treated the start/end locations of a single base deletion and insertion different to a SNP. In my example here, the start and end locations for those cases are the same to distinguish from the two-base change for TT->CA [this can be changed]:

cat 10454.txt | \
perl -ne 'if(/^(.*?):([0-9]+)(_([0-9]+))?(del|dup)?([ACGT]+)((ins|>)([ACGT]+))?/){
   print(join("\t",
         $1,
         $2,
         $4?$4:($2),
         ($8 eq "ins")?$9:($5 eq "dup")?"-":$6,
         (($5 eq "dup") || ($8 eq "ins"))?$6:($9)?$9:"-")."\n")
}';
chr7    140453145       140453145       A       T
chr7    140453136       140453136       A       T
chr7    140453135       140453136       TT      CA
chr20   31022287        31022287        -       A
chr19   13054614        13054614        G       -

Regular Expression

This uses enclosing groups to describe the important bits of the sequence. The groups are numbered from 1 upwards, with each number corresponding to a left bracket (in order), capturing the entirety of the sub-expression until the paired right bracket. I'll describe the groups, which are referenced in the subsequent part of the script as $<group_number>:

1. Any number of characters, up to (but not including) the first ':'
2. Any number of the characters '0' to '9' [i.e. numbers]
3. Optionally, a "_", followed by any number of the characters '0' to '9'
4. The numbers from Group [3], i.e. excluding the "_"
5. Optionally, a "dup" or a "del"
6. Any number of characters selected from "A", "C", "G", "T" (i.e. bases)
7. Optionally, an "ins" or ">", followed by any number of bases
8. The "ins" or ">" from Group [7]
9. The bases from Group [7]

Print Code

The print statement prints out a tab-separated list, finishing off with a line break. I'll describe each list element:

1. Group [1]
2. Group [2]
3. If Group [4] exists, then Group [4], otherwise Group [2]
4. If Group [8] is "ins", then Group [9]; [otherwise]
   If Group [5] is "dup", then "-", otherwise Group [6]
5. If Group [5] is "dup", or group [8] is "ins", then Group [6][otherwise]
   If Group [9] exists, then Group [9], otherwise "-"

Answer 4

Here is an R solution, would probably be slower than the suggested Python solution though:

bed <- readLines("bed.txt")

res <- data.frame()

for(i in bed) {

  temp <- strsplit(i, ":") %>% unlist()

  # case N>N

  if(grepl(">", temp[[2]])) {

    begin <- str_extract(temp[[2]], pattern = "^[:digit:]+") %>% as.integer() - 1
    end <- begin + 1
    ref <- str_extract(temp[[2]], pattern = "[A-Z]")
    alt <- str_extract(temp[[2]], pattern = "[A-Z]$")

    temp_res <- data.frame(temp[[1]], begin, end, ref, alt)
    res <- rbind(res, temp_res)
  }

  # case del_ins

  if(grepl("ins", temp[[2]])) {

    begin <- str_extract(temp[[2]], pattern = "^[:digit:]+") %>% as.integer()
    end <- str_extract(temp[[2]], pattern = "_[:digit:]+")
    end <- substring(end, 2, nchar(end)) %>% as.integer()
    ref <- str_extract(temp[[2]], pattern = "[A-Z]+$")
    alt <- str_extract(temp[[2]], pattern = "[A-Z]+")

    temp_res <- data.frame(temp[[1]], begin, end, ref, alt)
    res <- rbind(res, temp_res)
  }


  # case dup

  if(grepl("dup", temp[[2]])) {

    begin <- str_extract(temp[[2]], pattern = "^[:digit:]+") %>% as.integer()
    end <- begin
    ref <- "-"
    alt <- str_extract(temp[[2]], pattern = "[A-Z]$")

    temp_res <- data.frame(temp[[1]], begin, end, ref, alt)
    res <- rbind(res, temp_res)
  }

  # case del only

  if(grepl("del", temp[[2]]) & !grepl("ins", temp[[2]])) {

    begin <- str_extract(temp[[2]], pattern = "^[:digit:]+") %>% as.integer()
    end <- begin
    ref <- str_extract(temp[[2]], pattern = "[A-Z]$")
    alt <- "-"

    temp_res <- data.frame(temp[[1]], begin, end, ref, alt)
    res <- rbind(res, temp_res)
  }

}

print(res)

  temp..1..     begin       end ref alt
1      chr7 140453144 140453145   A   T
2      chr7 140453135 140453136   A   T
3      chr7 140453135 140453136  TT  CA
4     chr20  31022287  31022287   -   A
5     chr19  13054614  13054614   G   -