Extrair dados citados e rotulados de uma determinada coluna

1

Eu tenho um grande arquivo GTF , como abaixo:

 # ./stringtie -p 4 -G /home/humangenome_hg19/homo_gtf_file.gtf -o strAD1_as/transcripts.gtf -l strAD1 /home/software/star-2.5.2b/bin/Linux_x86_64/mapA1Aligned.sortedByCoord.out.bam                               
# StringTie version 1.3.2d                              
1   StringTie   transcript  30267   31109   1000    +   .   gene_id "strAD1.1"; transcript_id "strAD1.1.1"; reference_id "ENST00000469289"; ref_gene_id "ENSG00000243485"; ref_gene_name "MIR1302-10"; cov "0.028725"; FPKM "0.053510"; TPM "0.109957";
1   StringTie   exon    30267   30667   1000    +   .   gene_id "strAD1.1"; transcript_id "strAD1.1.1"; exon_number "1"; reference_id "ENST00000469289"; ref_gene_id "ENSG00000243485"; ref_gene_name "MIR1302-10"; cov "0.014218";
1   StringTie   exon    30976   31109   1000    +   .   gene_id "strAD1.1"; transcript_id "strAD1.1.1"; exon_number "2"; reference_id "ENST00000469289"; ref_gene_id "ENSG00000243485"; ref_gene_name "MIR1302-10"; cov "0.072139";

Eu quero ter a 9ª coluna com apenas gene_id , transcript_id , reference_id e ref_gene_id . Eles estão na 9ª coluna e separados por espaço (as próprias colunas são separadas por TAB). Você poderia por favor me ajudar como eu posso tal coluna com um comando simples no Linux? Eu não quero usar o Excel para isso.

    
por Mary 13.05.2017 / 09:37

2 respostas

3

Idealmente, como os dados estão no formato GTF, deve-se usar um analisador GTF para analisá-los. Atualmente, não tenho esse analisador ou biblioteca de análise instalada, portanto, minha solução é baseada apenas nos dados que você forneceu na pergunta.

Para extrair a 9ª coluna:

$ cut -f 9 data.gtf
gene_id "strAD1.1"; transcript_id "strAD1.1.1"; reference_id "ENST00000469289"; ref_gene_id "ENSG00000243485"; ref_gene_name "MIR1302-10"; cov "0.028725"; FPKM "0.053510"; TPM "0.109957";
gene_id "strAD1.1"; transcript_id "strAD1.1.1"; exon_number "1"; reference_id "ENST00000469289"; ref_gene_id "ENSG00000243485"; ref_gene_name "MIR1302-10"; cov "0.014218";
gene_id "strAD1.1"; transcript_id "strAD1.1.1"; exon_number "2"; reference_id "ENST00000469289"; ref_gene_id "ENSG00000243485"; ref_gene_name "MIR1302-10"; cov "0.072139";

Para obter os dados que desejamos, precisamos tratar transcrições e exons separadamente, já que seus atributos têm ordem diferente nos dados. Fazemos isso com awk e geramos campos diferentes nos dados de entrada, dependendo se a linha atual contém a string exon_number ou não:

$ cut -f 9 data.gtf | awk '/exon_number/ { print $2, $4, $8, $10; next } { print $2, $4, $6, $8 }'
"strAD1.1"; "strAD1.1.1"; "ENST00000469289"; "ENSG00000243485";
"strAD1.1"; "strAD1.1.1"; "ENST00000469289"; "ENSG00000243485";
"strAD1.1"; "strAD1.1.1"; "ENST00000469289"; "ENSG00000243485";

Em seguida, removemos as aspas duplas e ponto-e-vírgulas disso:

$ cut -f 9 data.gtf | awk '/exon_number/ { print $2, $4, $8, $10; next } { print $2, $4, $6, $8 }' | tr -d '";'
strAD1.1 strAD1.1.1 ENST00000469289 ENSG00000243485
strAD1.1 strAD1.1.1 ENST00000469289 ENSG00000243485
strAD1.1 strAD1.1.1 ENST00000469289 ENSG00000243485
    
por 13.05.2017 / 09:46
3

Talvez apenas:

< file cut -sd '"' -f2,4,8,10 | tr '"' ' '

Considere a entrada como uma lista de colunas separadas por " e extraindo as 2 nd , 4 th , 8 th e 10 colunas th .

Com o GNU cut , você pode substituir o | tr '"' ' ' por --output-delimiter=' ' .

Isso faz com que a suposição de que " caracteres não apareçam em outro lugar nas linhas, que esses atributos gene_id , transcript_id ... sempre apareçam e sempre nessa ordem.

Como foi observado por Kusalananda, esse não é o caso em sua amostra, onde deve ser 2,4,6,8 para a primeira linha e 2,4,8,10 para as outras.

Para fazer uma correspondência mais expressiva: que apenas a coluna delimitada por tabulação de 9 th deva ser considerada e que os nomes de atributo corretos sejam encontrados, você pode recorrer a expressões regulares como:

< file pcregrep -o1 -o2 -o3 -o4 --om-separator=' ' '(?x)
  ^(?:[^\t]*+\t){8}(?=[^\t]*? \b gene_id       \ +"([^"\t]*)")
                   (?=[^\t]*? \b transcript_id \ +"([^"\t]*)")
                   (?=[^\t]*? \b reference_id  \ +"([^"\t]*)")
                   (?=[^\t]*? \b ref_gene_id   \ +"([^"\t]*)")'

Se você não tiver uma versão pcregrep ou muito antiga para suportar -o1... , use perl :

< file perl -lne 'print "$1 $2 $3 $4" if m{
  ^(?:[^\t]*+\t){8}(?=[^\t]*? \b gene_id       \ +"([^"\t]*)")
                   (?=[^\t]*? \b transcript_id \ +"([^"\t]*)")
                   (?=[^\t]*? \b reference_id  \ +"([^"\t]*)")
                   (?=[^\t]*? \b ref_gene_id   \ +"([^"\t]*)")}x'

Esse regexp corresponde primeiro aos primeiros 8 campos ( (?:[^\t]*+\t){8} ) e, depois disso, temos quatro expressões de look-ahead ( (?=...) ), portanto, correspondemos a esses oito campos, desde que o seguinte corresponda a todos os quatro expressões de look-ahead. Cada expressão de look-ahead procura por um dos atributos e captura o valor (na (...) part). Esses valores capturados estão disponíveis em $1 , $2 , $3 , $4 .

Isso permite atributos em qualquer ordem.

Note que ele pode ser enganado por coisas como:

1 2 3 4 5 6 7 8 gene_id "transcript_id " ...

Embora seja possível abordá-lo, provavelmente não vale a pena o esforço, pois não espero que esteja ocorrendo na entrada.

Enquanto você está usando perl , você também pode fazer uma análise mais formal do campo 9 th . Algo como:

< file perl -F'\t' -lane '
  my %field;
  while ($F[8] =~ /(\w+) +"(.*?)"/g) {$field{$1}=$2}
  if (%field) {
    print join " ", @field{
      qw(gene_id transcript_id reference_id ref_gene_id
    )}
  }'

(aqui, imprimindo uma linha contanto que pelo menos um atributo seja encontrado (ao contrário de todos os atributos solicitados nas outras abordagens)).

    
por 13.05.2017 / 09:46