我有一个已排序的ID和数字文件(位置)。我需要将第二列中的位置分组为一组中的500。
如果与前一行相比,行的值小于500,则将它们分组为同一组;而如果行的值超过500,则将它们分组为不同的组。
输入文件:
snp00001 200
snp00002 300
snp00003 400
snp00004 500
snp00005 600
snp00006 900
snp00007 1500
snp00008 1800
snp00009 3000
snp00010 3500
snp00011 4000
snp00012 5000
期望输出
snp00001 200 Group1
snp00002 300 Group1
snp00003 400 Group1
snp00004 500 Group1
snp00005 600 Group1
snp00006 900 Group1
snp00007 1500 Group2
snp00008 1800 Group2
snp00009 3000 Group3
snp00010 3500 Group3
snp00011 4000 Group4
snp00012 5000 Group5
额外注意: snp00001 to snp00006将分组为同一组,因为它们之间的范围(snp00002 - snp00001)或(snp00003 - snp00002)或(snp00004 - snp00003) .不到500。
snp00006和snp00007被分组为下一组,因为它们之间的范围(snp00007 - snp00006)超过500。
我试过使用awk,但没有成功。
awk -v step=500 -v OFS='\t' '{if(NR==1 || $2+limit){group++} file="Group"group; print file,$0}' input_file
发布于 2022-04-20 01:35:55
您需要跟踪前一个值,并将当前值与此保存的值进行比较。如果差值超过500,则增加组数。
例如
awk -v group=1 '{if ($2-prev>500) { group++ }} {prev=$2; $3="group" group; print}'
snp00001 200 group1
snp00002 300 group1
snp00003 400 group1
snp00004 500 group1
snp00005 600 group1
snp00006 900 group1
snp00007 1500 group2
snp00008 1800 group2
snp00009 3000 group3
snp00010 3500 group3
snp00011 4000 group3
snp00012 5000 group4
(FWIW,您的9/10/11输出不一致;9->10是500,但不增加组,但10->11也是500,但确实增加组)。
发布于 2022-06-17 12:04:15
使用Raku (以前称为Perl_6)
这是一个有点不同的分组方案,可能被证明是有用的。它使用Raku的~~
smartmatch运算符快速判断某个位置是否在某个范围内(或不位于范围内):
~$ raku -e 'my $i = 1; my $r = 1..500; for lines() {my $a = .words; \
if ($a.[1].Int ~~ $r) {say "$a Group", $i, " ", $r} else { \
repeat { $r+=500 } until ($a.[1].Int ~~ $r); \
say "$a Group", ++$i, " ", $r };}' file
样本输入:
snp00001 200
snp00002 300
snp00003 400
snp00004 500
snp00005 600
snp00006 900
snp00007 1500
snp00008 1800
snp00009 3000
snp00010 3500
snp00011 4000
snp00012 5000
样本输出(从核苷酸1开始每500个核苷酸组SNP):
snp00001 200 Group1 1..500
snp00002 300 Group1 1..500
snp00003 400 Group1 1..500
snp00004 500 Group1 1..500
snp00005 600 Group2 501..1000
snp00006 900 Group2 501..1000
snp00007 1500 Group3 1001..1500
snp00008 1800 Group4 1501..2000
snp00009 3000 Group5 2501..3000
snp00010 3500 Group6 3001..3500
snp00011 4000 Group7 3501..4000
snp00012 5000 Group8 4501..5000
上面的Raku代码声明了一个Group#迭代器$i
和一个初始范围$r
of 1..500
。输入作为lines
,每一行都被分割成(空格分隔的) words
。运行if/else条件:if
第二列~~
智能匹配在$r
范围内,say
行、Group#和range,else
采用500
的$r
范围和repeat
增量,而不是(即until
) ~~
smartmatch成功。然后打印与以前相同的信息,但这次使用的是Group#适当的增量(++$i
)。
上述分组方案的优点是得到的基团都是等间隔长度的,在这种情况下约有500个核苷酸。这种方案可以防止组间隔长度的“膨胀”,就像当多个SNP在一起时可能发生的那样(间隔“扩展”,可能造成“聚类”的错误印象)。
为了使它成为一个更通用的“分组”工具,您可以将范围的右侧抽象为一个变量($m
),以便进行快速分组:
~$ raku -e 'my $i=1; my $m=1000; my $r = 1..$m; for lines() {my $a = .words; if ($a.[1].Int ~~ $r) {say "$a\tGroup$i\t", $r} else { repeat { $r+=$m } until ($a.[1].Int ~~ $r); say "$a\tGroup{++$i}\t", $r };}' file
snp00001 200 Group1 1..1000
snp00002 300 Group1 1..1000
snp00003 400 Group1 1..1000
snp00004 500 Group1 1..1000
snp00005 600 Group1 1..1000
snp00006 900 Group1 1..1000
snp00007 1500 Group2 1001..2000
snp00008 1800 Group2 1001..2000
snp00009 3000 Group3 2001..3000
snp00010 3500 Group4 3001..4000
snp00011 4000 Group4 3001..4000
snp00012 5000 Group5 4001..5000
https://unix.stackexchange.com/questions/699630
复制相似问题