染色体領域を持つGFFファイルである7つの列を持つファイルがあります。 REGION = "exon"を含む行をファイルの1行に縮小したいと思います。行は、互いに重なる領域に基づいて縮小する必要があります。
REGION START END SCORE STRAND FRAME ATTRIBUTE
exon 26453 26644 . + . Transcript "XM_092971"; Name "XM_092971"
exon 26842 27020 . + . Transcript "XM_092971"; Name "XM_092971"
exon 30355 30899 . - . Transcript "XM_104663"; Name "XM_104663"
GS_TRAN 30355 34083 . - . GS_TRAN "Hs22_30444_28_1_1"; Name "Hs22_30444_28_1_1"
snp 30847 30847 . + . SNP "rs2971719"; Name "rs2971719"
exon 31012 31409 . - . Transcript "XM_104663"; Name "XM_104663"
exon 34013 34083 . - . Transcript "XM_104663"; Name "XM_104663"
exon 40932 41071 . + . Transcript "XM_092971"; Name "XM_092971"
snp 44269 44269 . + . SNP "rs2873227"; Name "rs2873227"
snp 45723 45723 . + . SNP "rs2227095"; Name "rs2227095"
exon 134031 134495 . - . Transcript "XM_086913"; Name "XM_086913"
exon 134034 134457 . - . Transcript "XM_086914"; Name "XM_086914"
上記のサンプルデータを見ると、最後の2行だけを1行にマージできます。したがって、新しい行になります。
exon 134031 134495 . - . Transcript "XM_086913"; Name "XM_086913"
他の行の終わりがその前の行よりも大きい場合、この場合はこれがEND領域になります。デフォルトでは、重複する部分がある場合は、最初に始まる領域と後で終わる領域を取ります。
これらのインスタンスには複数行があります。ここでは最後の2行しかありません。 1つの点は、ATRRIBUTE列にはほとんど同じですが、行に対して間違いなく異なる成績表名が表示されることです。
進行方法に関する提案。
更新された例:最後の2行が次の場合
exon 134031 134457 . - . Transcript "XM_086913"; Name "XM_086913"
exon 134034 134495 . - . Transcript "XM_086914"; Name "XM_086914"
その後、出力は次のようになります。
exon 134031 134495 . - . Transcript "XM_086913"; Transcript "XM_086914"
デフォルトでは、最初の行で始まり、2行目で終わります。なぜなら2~3行、それ以上の行ではなく、1行の重なる部分だけをカバーしたいからです。ここで重なる部分は2行の間にありますが、2行以上でも構いません。
更新された例(2014年3月24日)
chr1 HAVANA stop_codon 1120520 1120522 . + 0 gene_id "ENSG00000162571.9"; transcript_id "ENST00000379288.3"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "protein_coding"; transcript_status "KNOWN"; transcript_name "TTLL10-001"; level 2; tag "CCDS"; ccdsid "CCDS8.1"; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002420.2";
chr1 HAVANA UTR 1115077 1115233 . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000379288.3"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "protein_coding"; transcript_status "KNOWN"; transcript_name "TTLL10-001"; level 2; tag "CCDS"; ccdsid "CCDS8.1"; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002420.2";
chr1 HAVANA UTR 1115414 1115433 . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000379288.3"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "protein_coding"; transcript_status "KNOWN"; transcript_name "TTLL10-001"; level 2; tag "CCDS"; ccdsid "CCDS8.1"; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002420.2";
chr1 HAVANA UTR 1120520 1121244 . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000379288.3"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "protein_coding"; transcript_status "KNOWN"; transcript_name "TTLL10-001"; level 2; tag "CCDS"; ccdsid "CCDS8.1"; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002420.2";
chr1 HAVANA transcript 1115864 1119307 . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000460998.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "processed_transcript"; transcript_status "KNOWN"; transcript_name "TTLL10-004"; level 2; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002423.1";
chr1 HAVANA exon 1115864 1116240 . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000460998.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "processed_transcript"; transcript_status "KNOWN"; transcript_name "TTLL10-004"; level 2; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002423.1";
chr1 HAVANA *exon 1117121 1117195* . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000460998.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "processed_transcript"; transcript_status "KNOWN"; transcript_name "TTLL10-004"; level 2; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002423.1";
chr1 HAVANA *exon 1117150 1117826* . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000460998.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "processed_transcript"; transcript_status "KNOWN"; transcript_name "TTLL10-004"; level 2; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002423.1";
chr1 HAVANA exon 1118256 1118427 . + . gene_id "ENSG00000162571.9"; transcript_id "ENST00000460998.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "TTLL10"; transcript_type "processed_transcript"; transcript_status "KNOWN"; transcript_name "TTLL10-004"; level 2; havana_gene "OTTHUMG00000000851.3"; havana_transcript "OTTHUMT00000002423.1";
chr1 HAVANA transcript 1190648 1209229 . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA exon 1209046 1209229 . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA exon 1203113 1203372 . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA CDS 1203241 1203372 . - 0 gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA start_codon 1203370 1203372 . - 0 gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA stop_codon 1203238 1203240 . - 0 gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA exon 1198726 1198766 . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA exon 1192588 1192690 . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA exon 1192372 1192510 . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA *exon 1191425 1191505* . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
chr1 HAVANA *exon 1190648 1191470* . - . gene_id "ENSG00000160087.16"; transcript_id "ENST00000473215.1"; gene_type "protein_coding"; gene_status "KNOWN"; gene_name "UBE2J2"; transcript_type "nonsense_mediated_decay"; transcript_status "KNOWN"; transcript_name "UBE2J2-003"; level 2; havana_gene "OTTHUMG00000001911.7"; havana_transcript "OTTHUMT00000005432.2";
上半分は「+」チェーンのネストを表し、下半分は「-」チェーンのネストを表します。 「-」チェーンは面積が減少するため、重なり合う部分は最後の2行に示されているものと同じです。どちらも異なる遺伝子であるため、時には異なる遺伝子にも重複するエクソンがあるため、重複が各遺伝子に対して必要であるが、一部の投稿で読んだように、これは非常にまれです。 「gene_name」と表示される最後の列から遺伝子情報を抽出できます。
gene_name = TTLL10の2行に重複するエクソンがあるため、最終出力でマージされます。
chr1 HAVANA *exon 1117121 1117195* . + . transcript_id "ENST00000460998.1"; gene_name "TTLL10";
chr1 HAVANA *exon 1117150 1117826* . + . transcript_id "ENST00000460998.1"; gene_name "TTLL10";
gene_name = UBE2J2の2行には重複するエクソンがあります。
chr1 HAVANA *exon 1191425 1191505* . - . transcript_id "ENST00000473215.1"; gene_name "UBE2J2";
chr1 HAVANA *exon 1190648 1191470* . - . transcript_id "ENST00000473215.1"; gene_name "UBE2J2";
サンプル出力
残りの行は変更されずに保持され、上記の行は各遺伝子についてマージされます。
chr1 HAVANA *exon 1117121 1117826* . + . transcript_id "ENST00000460998.1"; gene_name "TTLL10";
chr1 HAVANA *exon 1190648 1191505* . - . transcript_id "ENST00000473215.1"; gene_name "UBE2J2";
Transcript_idが異なる場合、gene_nameは同じままですが、両方のTranscript IDが印刷されます。例えば、遺伝子の場合、転写物IDは次のように異なる。
chr1 HAVANA *exon 1191425 1191505* . - . transcript_id "ENST00000473215.1"; gene_name "UBE2J2";
chr1 HAVANA *exon 1190648 1191470* . - . transcript_id "ENST00000473215.2"; gene_name "UBE2J2";
上記のようにマージされますが、成績表の名前は2つ必要です。なぜなら、後で成績表情報を保管することが必要であり重要であると考えるからです。
chr1 HAVANA *exon 1190648 1191505* . - . transcript_id "ENST00000473215.1"; "ENST00000473215.2" gene_name "UBE2J2";
答え1
「awk」方法、
awk '
$1!="exon" { # If the first died is unequal to "exon"
if(previous)print previous # If there is a previous line then print it
print # Print the current line
previous=start=end=exon_seq="" # Set all variable to an empty string
next # Move on to the next line in the input file
}
{
if(exon_seq) { # if there is a sequence of lines with "exon in field 1
if(start<=$2 && end>=$3) # if the start value (field 2) of the previous line
# is less or equal to the current line and the end
# value of the previous line is greater than or
# equal to field 3 of the current line
next # then do nothing and read the next line
else # if there is no overlap,
print previous # then print the previous line
}
else { # if we are not already in the a sequence of
# "exon" lines, then this is the first one
exon_seq=1 # so exon_seq should become 1
}
previous=$0; start=$2; end=$3 # `start` become field2, `end` becomes field 3 and
# `previous` becomes the current record (line)
}
END{ # After all lines are processed
if(previous) print previous # If there still is a previous line, then print it
}
' file
答え2
私はこれらの複雑なタスクを解決するためにPerlを使用します。以下は部分的な解決策です。うまく機能するように調整する必要があるかもしれません。
#!/usr/bin/perl
use warnings;
use strict;
use List::Util qw{ max };
sub output {
my $previous = shift;
print join ' ', 'exon',
@{$previous}{qw(start end score strand frame attribute)};
}
$\ = "\n";
my %previous;
while (<>) {
chomp;
my ($region, $start, $end, $score, $strand, $frame, $attribute)
= split ' ', $_, 7;
if ($. == 1) {
print;
} elsif ('exon' eq $region) {
if (keys %previous
and $start < $previous{end} # Overlap.
) {
if ($end > $previous{end}) { # Not contained.
$previous{attribute} =~ s/; Name .*//;
$previous{attribute} .= '; '
. ($attribute =~ /(Transcript ".*?")/)[0];
$previous{end} = $end;
}
} else {
if (keys %previous) {
output(\%previous);
}
%previous = ( start => $start,
end => $end,
score => $score,
strand => $strand,
frame => $frame,
attribute => $attribute,
);
}
} else {
output(\%previous) if keys %previous;
%previous = ();
}
}
output(\%previous) if keys %previous;