☰

두 글자로만 구성된 줄과 한 글자로 시작하고 특정 길이를 충족하는 줄의 경우 grep/awk/sed

tag-icon • tag-icon linux text-processing awk sed grep
두 글자로만 구성된 줄과 한 글자로 시작하고 특정 길이를 충족하는 줄의 경우 grep/awk/sed

DNA 염기서열 분석 데이터에 대해 잘 아는 사람이 얼마나 되는지 모르겠지만 이것이 내 파일의 일부인 경우(">"로 시작하는 줄은 ID이고 문자로 시작하는 줄은 DNA 염기서열입니다):

>NB501013:9:HJJ75BGXX:4:13609:24076:18015/2
GGGGGGGAAAAAAA
>NB501013:9:HJJ75BGXX:4:21602:19346:16945/2
CTCGTCGCATCACAAAGGGAT
>NB501013:9:HJJ75BGXX:3:11407:17650:13229/2
CCGCGGGCCGGTGCGGGGGTTTTTTTGTTTTTTTGGTTACAACGGGTGGG
>NB501013:9:HJJ75BGXX:3:13509:1817:13239/2
CAGCCC
>NB501013:9:HJJ75BGXX:4:22611:20567:13384/2
GAATA

이 줄을 제거하고 싶습니다: GGGGGGGAAAAAAA

시퀀싱 ID와 함께(이 작업을 수행하는 데 사용할 수 있다는 것을 알고 있습니다 grep -B1). 그런데 두 글자로만 구성된 줄을 삭제하는 방법을 아는 사람이 있나요?

또한 5글자보다 짧은 시퀀스의 경우 해당 ID와 해당 ID를 제거하고 싶습니다. 모든 ID가 길기 때문에 특정 길이보다 긴 줄을 단순히 grep할 수 없으므로 어떻게든 다음 grep -v으로 시작하는 코드를 사용해야 합니다. 문자(따라서 ">"로 시작하지 않는 줄)이며 특정 길이보다 깁니다.

따라서 내 예제 출력은 다음과 같습니다.

>NB501013:9:HJJ75BGXX:4:21602:19346:16945/2
CTCGTCGCATCACAAAGGGAT
>NB501013:9:HJJ75BGXX:3:11407:17650:13229/2
CCGCGGGCCGGTGCGGGGGTTTTTTTGTTTTTTTGGTTACAACGGGTGGG
>NB501013:9:HJJ75BGXX:3:13509:1817:13239/2
CAGCCC

답변1

erl 호환 gexp 모듈을 grep사용해 보십시오 :PCRE

  • 두 글자 조합 삭제:

    pcregrep -Mv '>.*\n([ACGT])\1*([ACGT])\2*(\1|\2)*$' file

    산출:

    >NB501013:9:HJJ75BGXX:4:21602:19346:16945/2
CTCGTCGCATCACAAAGGGAT
>NB501013:9:HJJ75BGXX:3:11407:17650:13229/2
CCGCGGGCCGGTGCGGGGGTTTTTTTGTTTTTTTGGTTACAACGGGTGGG
>NB501013:9:HJJ75BGXX:3:13509:1817:13239/2
CAGCCC
>NB501013:9:HJJ75BGXX:4:22611:20567:13384/2
GAATA

  • 5자 이하의 조합을 삭제하세요.

     pcregrep -Mv '>.*\n[ACGT]{1,5}$' file

    산출:

    >NB501013:9:HJJ75BGXX:4:13609:24076:18015/2
GGGGGGGAAAAAAA
>NB501013:9:HJJ75BGXX:4:21602:19346:16945/2
CTCGTCGCATCACAAAGGGAT
>NB501013:9:HJJ75BGXX:3:11407:17650:13229/2
CCGCGGGCCGGTGCGGGGGTTTTTTTGTTTTTTTGGTTACAACGGGTGGG
>NB501013:9:HJJ75BGXX:3:13509:1817:13239/2
CAGCCC

답변2

#!/usr/bin/env perl
#
# Usage: thisscriptname < someinputfile

use strict;
use warnings;

while (1) {
  exit if eof;
  # rash assumption there are always pairs of ID and sequence lines
  # NOTE these contain a newline, so many need chomp() depending
  # on what you do with them...
  my $id = readline;
  my $seq = readline;

  # calculate unique sequence letters via hash (is there also a U
  # or something? been a few decades since AP bio...)
  my %chars;
  $chars{$_}++ for $seq =~ m/[ATGC]/g;

  # business logic time!
  if (keys %chars > 2 and length $seq >= 5) {
    print $id;
    print $seq;
  }
}

답변3

파일을 뒤집고 DNA 서열을 테스트하는 것을 고려할 수 있습니다. 테스트에 통과하면 이 줄을 무시하고다음철사:

tac file |
  awk '!/^>/ && (length($1) < 5 || $1 == "GGGGGGGAAAAAAA") {getline; next} 1' |
  tac

관련 정보