0

I have 2 genetic datasets, both have a column for chromosome position. I am trying to filter file 1 based on each line's chromosome position in comparison to the chromosome positions in file 2. The condition for file 1 line extraction is that only lines that have a chromosome position either more than 5000 larger or less than 5000 smaller than any given chromosome position in all the lines of file 2 is selected. Ultimately I am trying to make sure those genetic variants I extract from file 1 have a greater than 5000 +/- distance than those in file 2.

For example my data looks like:

File 1:

Variant   Chromosome Position
Variant1    14000
Variant2    9000
Variant3    37000
Variant4    21000

File 2:

Variant   Chromosome Position
Variant1    10000
Variant2    20000
Variant3    30000

Expected output (of variants with a greater than 5000 position distance in comparison to each line of file 2):

Variant   Chromosome Position
Variant3    37000

I've found questions extracting data from 2 files, but depending on matching or unqiue. I haven't yet found something applying +/- conditions to solve my question and I'm not sure where else to look/learn. I'm new to linux so I've been trying to construct a command which would have the effect of:

awk 'BEGIN{FS=OFS="\t"} FNR==1{print;next} {if($2>=+5000 | $2<=-5000) print $0}' file1 file2 > newfile

but I'm unsure of the syntax to get this working - any help or guidance would be appreciated.

  • Will there be holes between "ranges" (e.g. Var1 9000 & Var2 21000 -> gap 14001 - 15999)? – RudiC Feb 5 at 20:51
1

Try also (taking "holes" between ranges into consideration)

$ awk '
NR == FNR       {RGMIN[++IX] = $2 - 5000                        # from file2, create min and
                 RGMAX[IX]   = $2 + 5000                        # max range boundaries
                 next
                }
FNR == 1        {print                                          # print file1 header; don't process line further
                 next
                }
                {PR = 1                                         # set PRintout logical variable

                 for (i=2; i<=IX; i++)  PR = PR * ($2 < RGMIN[i] || $2 > RGMAX[i])

                                                                # reset PR to 0 if $2 inside any of the
                                                                # ranges; skip file2's header by starting
                                                                # loop at IX 2

                }
PR                                                              # print if outside ALL ranges
' file2 file1
Variant   Chromosome Position
Variant3    37000
|improve this answer|||||
  • Thank you for this, this seems to have worked. Do you know any resource where I can learn about the functions you've used like RGMIN and RGMAX? I googled for them but can't find any resources – DN1 Feb 6 at 10:18
  • Those two are locally defined array variables: Range min / max boundaries. – RudiC Feb 6 at 21:26

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.